Top 5 Recent Articles
ARTICLES CATEGORIES
- Algorithms (22)
- All (399)
- Biography (1)
- Blog (44)
- Business Requirements (1)
- Commentary (1)
- Conversion (2)
- Customers (2)
- Data Models (1)
- Education (2)
- GeoRaptor (13)
- GPS (1)
- Image Processing (2)
- Import Export (8)
- Licensing (2)
- LiDAR (1)
- Linear Referencing (4)
- Manifold GIS (3)
- Mapping (1)
- MySQL Spatial (7)
- Networking and Routing (including Optimization) (5)
- Open Source (18)
- Oracle Spatial and Locator (194)
- Partitioning (1)
- PostGIS (36)
- Projections (1)
- Published Articles (1)
- qGIS (1)
- Recommendations (1)
- Services (1)
- Software Change Log (1)
- Source Code (37)
- Space Curves (9)
- Spatial Database Functions (109)
- Spatial DB comparison (1)
- Spatial XML Processing (11)
- SQL Server Spatial (92)
- Standards (3)
- Stored Procedure (17)
- Tessellation or Gridding (10)
- Tools (2)
- Topological Relationships (1)
- Training (2)
- Triangulation (2)
Useful Package of Wrapper Functions for Sdo_Util.AffineTransforms
—————————————————————————————-
UPDATE (21 Mar 2018): Wrapper functions have been added to the AFFINE package to provide easier access to the SDO_UTIL.AFFINETRANSFORMS function. The AFFINE package now supports three different implementations of most of its functions:
- Original PL/SQL functions
- ST_AFFINE based on SYS.UTL_NLA_ARRAY_DBL
- SDO_UTIL.AFFINETRANSFORMS
Improved testing of the AFFINE functions has also been implemented.
The AFFINE package is capable of being installed INDEPENDENTLY of the other packages such as GEOM.
—————————————————————————————-
UPDATE (18 Jul 2012): I’ve now included my ST_Affine function which wraps Oracle’s matrix packages eg SYS.UTL_NLA_ARRAY_DBL. The implementation makes it easier to compare outputs with those from PostGIS’s ST_Affine.
Before 11g, I created a set of affine functions in my free GEOM package that enabled a user to manipulate an sdo_geometry object as follows:
- Rotate
- Move (translate)
- Scale
These functions do not use matrix algebra.
In a later version (10g onwards), I used the SYS.UTL_NLA_ARRAY_DBL package to expose an Affine function which aligns closely with the PostGIS ST_Affine functions:
FUNCTION ST_Affine(p_geom IN mdsys.sdo_geometry, p_a NUMBER, p_b NUMBER, p_c NUMBER, p_d NUMBER, p_e NUMBER, p_f NUMBER, p_g NUMBER, p_h NUMBER, p_i NUMBER, p_xoff NUMBER, p_yoff NUMBER, p_zoff NUMBER) RETURN mdsys.sdo_geometry IS -- Transformation matrix is represented by: -- / a b c xoff \ -- | d e f yoff | -- | g h i zoff | -- \ 0 0 0 1 / -- -- For 2D only need to supply: a, b, d, e, xoff, yoff ...
However, the Oracle Spatial team, sensibly, decided to release their own functionality for affine transformations in the Sdo_Util.AffineTransforms function.
The problem with this function is that it a veritable kitchen sink of capability which, coupled with some pretty limited error handling, is very hard to use.
So, in this blog I propose a PL/SQL package that wraps this function so that it is more useable.
Here it is.
CREATE OR REPLACE package affine AUTHID CURRENT_USER AS FUNCTION PI RETURN NUMBER deterministic; FUNCTION ST_Rotate(p_geometry IN mdsys.sdo_geometry, p_angle_rad IN NUMBER) RETURN mdsys.sdo_geometry deterministic; FUNCTION ST_Rotate(p_geometry IN mdsys.sdo_geometry, p_angle_rad IN NUMBER, p_dir IN pls_integer, p_rotate_point IN mdsys.sdo_geometry, p_line1 IN mdsys.sdo_geometry) RETURN mdsys.sdo_geometry deterministic; FUNCTION ST_Rotate(p_geometry IN mdsys.sdo_geometry, p_angle_rad IN NUMBER, p_rotate_x IN NUMBER, p_rotate_y IN NUMBER) RETURN mdsys.sdo_geometry deterministic; FUNCTION ST_Rotate(p_geometry IN mdsys.sdo_geometry, p_angle_rad IN NUMBER, p_rotate_point IN mdsys.sdo_geometry) RETURN mdsys.sdo_geometry deterministic; FUNCTION ST_Scale(p_geometry IN mdsys.sdo_geometry, p_sx IN NUMBER, p_sy IN NUMBER, p_sz IN NUMBER, p_scale_pt IN mdsys.sdo_geometry) RETURN mdsys.sdo_geometry deterministic; FUNCTION ST_Scale(p_geometry IN mdsys.sdo_geometry, p_sx IN NUMBER, p_sy IN NUMBER, p_sz IN NUMBER) RETURN mdsys.sdo_geometry deterministic; FUNCTION ST_Scale(p_geometry IN mdsys.sdo_geometry, p_sx IN NUMBER, p_sy IN NUMBER) RETURN mdsys.sdo_geometry deterministic; FUNCTION ST_Translate(p_geometry IN mdsys.sdo_geometry, p_tx IN NUMBER, p_ty IN NUMBER, p_tz IN NUMBER) RETURN mdsys.sdo_geometry deterministic; FUNCTION ST_Translate(p_geometry IN mdsys.sdo_geometry, p_tx IN NUMBER, p_ty IN NUMBER) RETURN mdsys.sdo_geometry deterministic; FUNCTION ST_RotateTranslateScale(p_geometry IN mdsys.sdo_geometry, p_angle_rad IN NUMBER, p_rs_point IN mdsys.sdo_geometry, p_sx IN NUMBER, p_sy IN NUMBER, p_sz IN NUMBER, p_tx IN NUMBER, p_ty IN NUMBER, p_tz IN NUMBER) RETURN mdsys.sdo_geometry deterministic; /* ---------------------------------------------------------------------------------------- * @function : Affine * @precis : Applies a 3d affine transformation to the geometry to do things like translate, rotate, scale in one step. * @version : 1.0 * @description: Applies a 3d affine transformation to the geometry to do things like translate, rotate, scale in one step. * To apply a 2D affine transformation only supply a, b, d, e, xoff, yoff * @usage : Function Affine ( p_geom IN MDSYS.SDO_GEOMETRY, * @param : p_geom : MDSYS.SDO_GEOMETRY : The shape to rotate. * @param : a, b, c, d, e, f, g, h, i, xoff, yoff, zoff : * Represent the transformation matrix * / a b c xoff \ * | d e f yoff | * | g h i zoff | * \ 0 0 0 1 / * and the vertices are transformed as follows: * x' = a*x + b*y + c*z + xoff * y' = d*x + e*y + f*z + yoff * z' = g*x + h*y + i*z + zoff * @requires : SDO_UTIL.GetVertices Function * SYS.UTL_NLA Package * SYS.UTL_NLA_ARRAY_DBL Type * SYS.UTL_NLA_ARRAY_INT Type * @return : newGeom : MDSYS.SDO_GEOMETRY : Transformed input geometry. * @note : Cartesian arithmetic only * : Not for Oracle XE. Only 10g and above. * @history : Simon Greener, SpatialDB Advisor - Feb 2009 - Original coding. * @copyright : Simon Greener, 2011, 2012 * @license : Licensed under a Creative Commons Attribution-Share Alike 2.5 Australia License. * http://creativecommons.org/licenses/by-sa/2.5/au/. * Any bugs or improvements to be supplied back to Simon Greener **/ FUNCTION ST_Affine( p_geometry IN mdsys.sdo_geometry, p_a IN NUMBER, p_b IN NUMBER, p_c IN NUMBER, p_d IN NUMBER, p_e IN NUMBER, p_f IN NUMBER, p_g IN NUMBER, p_h IN NUMBER, p_i IN NUMBER, p_xoff IN NUMBER, p_yoff IN NUMBER, p_zoff IN NUMBER) RETURN mdsys.sdo_geometry deterministic; END Affine; / SHOW errors CREATE OR REPLACE package body affine AS c_PI CONSTANT NUMBER(16,14) := 3.14159265358979; c_i_null_geometry CONSTANT INTEGER := -20120; c_s_null_geometry CONSTANT VARCHAR2(100) := 'Input geometry must not be null'; c_i_null_parameter CONSTANT INTEGER := -20127; c_s_null_parameter CONSTANT VARCHAR2(100) := 'Input parameters must not be null'; FUNCTION PI RETURN NUMBER AS BEGIN RETURN c_PI; END PI; FUNCTION ST_Rotate(p_geometry IN mdsys.sdo_geometry, p_angle_rad IN NUMBER, p_dir IN pls_integer, p_rotate_point IN mdsys.sdo_geometry, p_line1 IN mdsys.sdo_geometry) RETURN mdsys.sdo_geometry AS NULL_GEOMETRY EXCEPTION; BEGIN IF ( p_geometry IS NULL ) THEN raise NULL_GEOMETRY; END IF; -- dir - Rotation parameter for x(0), y(1), or z(2)-axis roll. -- You cannot set p_dir => 0, 1 or 2. Only -1, -2, -3. They don't see to affect the result. -- For two-dimensional geometries, rotation uses the p1 and angle values. IF ( p_geometry.get_dims() = 2 ) THEN IF ( p_angle_rad IS NULL AND p_rotate_point IS NULL ) THEN raise_application_error(-20001,'For 2D geometry rotation, p_angle_rad and p_rotate_point must not be null',TRUE); END IF; ELSE -- For three-dimensional geometries, rotation uses either: -- 1. the angle and dir values or -- 2. the angle and line1 values. IF ( p_angle_rad IS NULL ) THEN raise_application_error(-20001,'For 3D geometry rotation, p_angle_rad must not be null',TRUE); END IF; IF ( p_dir IS NULL AND p_line1 IS NULL ) THEN raise_application_error(-20001,'For 3D geometry rotation, both p_dir and p_line1 cannot be null',TRUE); END IF; END IF; RETURN SDO_UTIL.AffineTransforms ( geometry => p_geometry, rotation => 'TRUE', p1 => p_rotate_point, angle => p_angle_rad, dir => p_dir, line1 => p_line1, translation => 'FALSE', tx => 0.0, ty => 0.0, tz => 0.0, scaling => 'FALSE', psc1 => NULL, sx => 0.0, sy => 0.0, sz => 0.0, shearing => 'FALSE', shxy => 0.0, shyx => 0.0, shxz => 0.0, shzx => 0.0, shyz => 0.0, shzy => 0.0, reflection => 'FALSE', pref => NULL, lineR => NULL, dirR => -1, planeR => 'FALSE', n => NULL, bigD => NULL ); EXCEPTION WHEN NULL_GEOMETRY THEN raise_application_error(c_i_null_geometry,c_s_null_geometry,TRUE); RETURN p_geometry; END ST_Rotate; FUNCTION ST_Rotate(p_geometry IN mdsys.sdo_geometry, p_angle_rad IN NUMBER, p_rotate_x IN NUMBER, p_rotate_y IN NUMBER) RETURN mdsys.sdo_geometry AS NULL_GEOMETRY EXCEPTION; BEGIN IF ( p_geometry IS NULL ) THEN raise NULL_GEOMETRY; END IF; -- dir - Rotation parameter for x(0), y(1), or z(2)-axis roll. -- For two-dimensional geometries, rotation uses the p1 and angle values. IF ( p_geometry.get_dims() = 2 ) THEN IF ( p_angle_rad IS NULL AND p_rotate_x IS NULL AND p_rotate_y IS NULL ) THEN raise_application_error(-20001,'For 2D geometry rotation, p_angle_rad, p_rotate_x and p_rotate_y must not be null.',TRUE); END IF; ELSE raise_application_error(-20001,'This version of ST_Rotate only supports 2D geometry rotation.',TRUE); END IF; RETURN ST_Rotate(p_geometry => p_geometry, p_angle_rad => p_angle_rad, p_dir => -1, p_rotate_point => mdsys.sdo_geometry(2001,p_geometry.sdo_srid,mdsys.sdo_point_type(p_rotate_x,p_rotate_y,NULL),NULL,NULL), p_line1 => NULL); EXCEPTION WHEN NULL_GEOMETRY THEN raise_application_error(c_i_null_geometry,c_s_null_geometry,TRUE); RETURN p_geometry; END ST_Rotate; FUNCTION ST_Rotate(p_geometry IN mdsys.sdo_geometry, p_angle_rad IN NUMBER, p_rotate_point IN mdsys.sdo_geometry) RETURN mdsys.sdo_geometry AS NULL_GEOMETRY EXCEPTION; BEGIN IF ( p_geometry IS NULL ) THEN raise NULL_GEOMETRY; END IF; -- dir - Rotation parameter for x(0), y(1), or z(2)-axis roll. -- For two-dimensional geometries, rotation uses the p1 and angle values. IF ( p_geometry.get_dims() = 2 ) THEN IF ( p_angle_rad IS NULL AND p_rotate_point IS NULL ) THEN raise_application_error(-20001,'For 2D geometry rotation, p_angle_rad and p_rotate_point must not be null.',TRUE); END IF; ELSE raise_application_error(-20001,'This version of ST_Rotate only supports 2D geometry rotation.',TRUE); END IF; RETURN ST_Rotate(p_geometry => p_geometry, p_angle_rad => p_angle_rad, p_dir => -1, p_rotate_point => p_rotate_point, p_line1 => NULL); EXCEPTION WHEN NULL_GEOMETRY THEN raise_application_error(c_i_null_geometry,c_s_null_geometry,TRUE); RETURN p_geometry; END ST_Rotate; FUNCTION ST_Rotate(p_geometry IN mdsys.sdo_geometry, p_angle_rad IN NUMBER) RETURN mdsys.sdo_geometry AS NULL_GEOMETRY EXCEPTION; BEGIN IF ( p_geometry IS NULL ) THEN raise NULL_GEOMETRY; END IF; -- dir - Rotation parameter for x(0), y(1), or z(2)-axis roll. -- For two-dimensional geometries, rotation uses the p1 and angle values. IF ( p_geometry.get_dims() = 2 ) THEN IF ( p_angle_rad IS NULL ) THEN raise_application_error(-20001,'For 2D geometry rotation, p_angle_rad must not be null.',TRUE); END IF; ELSE raise_application_error(-20001,'This version of ST_Rotate only supports 2D geometry rotation.',TRUE); END IF; RETURN ST_Rotate(p_geometry => p_geometry, p_angle_rad => p_angle_rad, p_dir => -1, p_rotate_point => mdsys.sdo_geometry(2001,p_geometry.sdo_srid,mdsys.sdo_point_type(0.0,0.0,0.0),NULL,NULL), p_line1 => NULL); EXCEPTION WHEN NULL_GEOMETRY THEN raise_application_error(c_i_null_geometry,c_s_null_geometry,TRUE); RETURN p_geometry; END ST_Rotate; /** =================================================================== **/ /** Scales the geometry to a new size by multiplying the ordinates with the parameters: * ST_Scale(geom, Xfactor, Yfactor, Zfactor, scale Point). * ST_Scale(geom, Xfactor, Yfactor, Zfactor). * ST_Scale(geom, Xfactor, Yfactor). * **/ FUNCTION ST_Scale(p_geometry IN mdsys.sdo_geometry, p_sx IN NUMBER, p_sy IN NUMBER, p_sz IN NUMBER, p_scale_pt IN mdsys.sdo_geometry) RETURN mdsys.sdo_geometry AS -- psc1 is Point on the input geometry about which to perform the scaling v_dims pls_integer; v_gtype pls_integer; v_vertex mdsys.vertex_type; v_psc1 mdsys.sdo_geometry; v_sx NUMBER := CASE WHEN p_sx IS NULL THEN 0.0 ELSE p_sx END; v_sy NUMBER := CASE WHEN p_sy IS NULL THEN 0.0 ELSE p_sy END; v_sz NUMBER := CASE WHEN p_sz IS NULL THEN 0.0 ELSE p_sz END; NULL_GEOMETRY EXCEPTION; BEGIN IF ( p_geometry IS NULL ) THEN raise NULL_GEOMETRY; END IF; -- Point on the input geometry about which to perform the scaling. -- This geometry should be either a zero point (with 0,0 or 0,0,0 -- ordinates for scaling about the origin) or a nonzero point (with ordinates -- for scaling about a point other than the origin). v_dims := p_geometry.get_dims(); v_gtype := (v_dims * 1000) + 1; IF ( p_scale_pt IS NULL ) THEN v_psc1 := mdsys.sdo_geometry(v_gtype,p_geometry.sdo_srid,mdsys.sdo_point_type(0.0,0.0,CASE WHEN v_gtype=3001 THEN 0.0 ELSE NULL END),NULL,NULL); ELSE v_vertex := sdo_util.getVertices(p_geometry)(1); v_psc1 := mdsys.sdo_geometry(v_gtype,p_geometry.sdo_srid,mdsys.sdo_point_type(v_vertex.x,v_vertex.y,v_vertex.z),NULL,NULL); END IF; RETURN SDO_UTIL.AffineTransforms ( geometry => p_geometry, scaling => 'TRUE', psc1 => v_psc1, sx => v_sx, sy => v_sy, sz => v_sz, rotation => 'FALSE', p1 => NULL, angle => 0.0, dir => -1, line1 => NULL, translation => 'FALSE', tx => 0.0, ty => 0.0, tz => 0.0, shearing => 'FALSE', shxy => 0.0, shyx => 0.0, shxz => 0.0, shzx => 0.0, shyz => 0.0, shzy => 0.0, reflection => 'FALSE', pref => NULL, lineR => NULL, dirR => -1, planeR => 'FALSE', n => NULL, bigD => NULL ); EXCEPTION WHEN NULL_GEOMETRY THEN raise_application_error(c_i_null_geometry,c_s_null_geometry,TRUE); RETURN p_geometry; END ST_Scale; FUNCTION ST_Scale(p_geometry IN mdsys.sdo_geometry, p_sx IN NUMBER, p_sy IN NUMBER, p_sz IN NUMBER) RETURN mdsys.sdo_geometry AS NULL_GEOMETRY EXCEPTION; BEGIN IF ( p_geometry IS NULL ) THEN raise NULL_GEOMETRY; END IF; RETURN ST_Scale(p_geometry, p_sx, p_sy, p_sz, NULL); EXCEPTION WHEN NULL_GEOMETRY THEN raise_application_error(c_i_null_geometry,c_s_null_geometry,TRUE); RETURN p_geometry; END ST_Scale; FUNCTION ST_Scale(p_geometry IN mdsys.sdo_geometry, p_sx IN NUMBER, p_sy IN NUMBER) RETURN mdsys.sdo_geometry AS NULL_GEOMETRY EXCEPTION; BEGIN IF ( p_geometry IS NULL ) THEN raise NULL_GEOMETRY; END IF; RETURN ST_Scale(p_geometry, p_sx, p_sy, 0.0, NULL); EXCEPTION WHEN NULL_GEOMETRY THEN raise_application_error(c_i_null_geometry,c_s_null_geometry,TRUE); RETURN p_geometry; END ST_Scale; /** =================================================================== **/ /** Translates the geometry to a new location using the numeric parameters as offsets. * ST_Translate(geom, X, Y) or * ST_Translate(geom, X, Y, Z) */ FUNCTION ST_Translate(p_geometry IN mdsys.sdo_geometry, p_tx IN NUMBER, p_ty IN NUMBER, p_tz IN NUMBER) RETURN mdsys.sdo_geometry AS v_tx NUMBER := CASE WHEN p_tx IS NULL THEN 0.0 ELSE p_tx END; v_ty NUMBER := CASE WHEN p_ty IS NULL THEN 0.0 ELSE p_ty END; v_tz NUMBER := CASE WHEN p_tz IS NULL THEN 0.0 ELSE p_tz END; NULL_GEOMETRY EXCEPTION; BEGIN IF ( p_geometry IS NULL ) THEN raise NULL_GEOMETRY; END IF; IF ( p_geometry.get_dims() = 2 ) THEN v_tz := 0.0; END IF; RETURN SDO_UTIL.AffineTransforms ( geometry => p_geometry, translation => 'TRUE', tx => v_tx, ty => v_ty, tz => v_tz, scaling => 'FALSE', psc1 => NULL, sx => 0.0, sy => 0.0, sz => 0.0, rotation => 'FALSE', p1 => NULL, angle => 0.0, dir => -1, line1 => NULL, shearing => 'FALSE', shxy => 0.0, shyx => 0.0, shxz => 0.0, shzx => 0.0, shyz => 0.0, shzy => 0.0, reflection => 'FALSE', pref => NULL, lineR => NULL, dirR => -1, planeR => 'FALSE', n => NULL, bigD => NULL ); EXCEPTION WHEN NULL_GEOMETRY THEN raise_application_error(c_i_null_geometry,c_s_null_geometry,TRUE); RETURN p_geometry; END ST_Translate; FUNCTION ST_Translate(p_geometry IN mdsys.sdo_geometry, p_tx IN NUMBER, p_ty IN NUMBER) RETURN mdsys.sdo_geometry AS NULL_GEOMETRY EXCEPTION; BEGIN IF ( p_geometry IS NULL ) THEN raise NULL_GEOMETRY; END IF; RETURN AFFINE.ST_Translate(p_geometry, p_tx, p_ty, 0.0); EXCEPTION WHEN NULL_GEOMETRY THEN raise_application_error(c_i_null_geometry,c_s_null_geometry,TRUE); RETURN p_geometry; END ST_Translate; FUNCTION ST_RotateTranslateScale(p_geometry IN mdsys.sdo_geometry, p_angle_rad IN NUMBER, p_rs_point IN mdsys.sdo_geometry, p_sx IN NUMBER, p_sy IN NUMBER, p_sz IN NUMBER, p_tx IN NUMBER, p_ty IN NUMBER, p_tz IN NUMBER) RETURN mdsys.sdo_geometry AS v_dims pls_integer; v_gtype pls_integer; v_psc1 mdsys.sdo_geometry; v_vertex mdsys.vertex_type; v_sx NUMBER := CASE WHEN p_sx IS NULL THEN 0.0 ELSE p_sx END; v_sy NUMBER := CASE WHEN p_sy IS NULL THEN 0.0 ELSE p_sy END; v_sz NUMBER := CASE WHEN p_sz IS NULL THEN 0.0 ELSE p_sz END; v_tx NUMBER := CASE WHEN p_tx IS NULL THEN 0.0 ELSE p_tx END; v_ty NUMBER := CASE WHEN p_ty IS NULL THEN 0.0 ELSE p_ty END; v_tz NUMBER := CASE WHEN p_tz IS NULL THEN 0.0 ELSE p_tz END; NULL_GEOMETRY EXCEPTION; BEGIN IF ( p_geometry IS NULL ) THEN raise NULL_GEOMETRY; END IF; IF ( p_geometry.get_dims() = 2 ) THEN v_tz := 0.0; END IF; -- dir - Rotation parameter for x(0), y(1), or z(2)-axis roll. -- You cannot set p_dir => 0, 1 or 2. Only -1, -2, -3. They don't see to affect the result. -- For two-dimensional geometries, rotation uses the p1 and angle values. IF ( p_geometry.get_dims() = 2 ) THEN IF ( p_angle_rad IS NULL AND p_rs_point IS NULL ) THEN raise_application_error(-20001,'For 2D geometry rotation, p_angle_rad and p_rs_point must not be null.',TRUE); END IF; ELSE raise_application_error(-20001,'This function does not support 3D geometry rotation - Use other functions.',TRUE); END IF; v_dims := p_geometry.get_dims(); v_gtype := (v_dims * 1000) + 1; IF ( p_rs_point IS NULL ) THEN v_psc1 := mdsys.sdo_geometry(v_gtype,p_geometry.sdo_srid,mdsys.sdo_point_type(0.0,0.0,CASE WHEN v_gtype=3001 THEN 0.0 ELSE NULL END),NULL,NULL); ELSE v_vertex := sdo_util.getVertices(p_rs_point)(1); v_psc1 := mdsys.sdo_geometry(v_gtype,p_geometry.sdo_srid,mdsys.sdo_point_type(v_vertex.x,v_vertex.y,CASE WHEN v_gtype=3001 THEN v_vertex.z ELSE NULL END),NULL,NULL); END IF; RETURN SDO_UTIL.AffineTransforms ( geometry => p_geometry, rotation => 'TRUE', p1 => p_rs_point, angle => p_angle_rad, dir => -1, line1 => NULL, translation => 'TRUE', tx => v_tx, ty => v_ty, tz => v_tz, scaling => 'TRUE', psc1 => v_psc1, sx => v_sx, sy => v_sy, sz => v_sz, shearing => 'FALSE', shxy => 0.0, shyx => 0.0, shxz => 0.0, shzx => 0.0, shyz => 0.0, shzy => 0.0, reflection => 'FALSE', pref => NULL, lineR => NULL, dirR => -1, planeR => 'FALSE', n => NULL, bigD => NULL ); EXCEPTION WHEN NULL_GEOMETRY THEN raise_application_error(c_i_null_geometry,c_s_null_geometry,TRUE); RETURN p_geometry; END ST_RotateTranslateScale; /** =================================================================== **/ FUNCTION ST_Affine( p_geometry IN mdsys.sdo_geometry, p_a IN NUMBER, p_b IN NUMBER, p_c IN NUMBER, p_d IN NUMBER, p_e IN NUMBER, p_f IN NUMBER, p_g IN NUMBER, p_h IN NUMBER, p_i IN NUMBER, p_xoff IN NUMBER, p_yoff IN NUMBER, p_zoff IN NUMBER) RETURN mdsys.sdo_geometry IS -- Transformation matrix is represented by: -- / a b c xoff \ -- | d e f yoff | -- | g h i zoff | -- \ 0 0 0 1 / -- -- For 2D only need to supply: a, b, d, e, xoff, yoff v_A SYS.UTL_NLA_ARRAY_DBL := SYS.UTL_NLA_ARRAY_DBL( p_a, p_d, NVL(p_g,0), 0, p_b, p_e, NVL(p_h,0), 0, NVL(p_c,0), NVL(p_f,0), NVL(p_i,1), 0, p_xoff, p_yoff, NVL(p_zoff,0), 1 ); v_C SYS.UTL_NLA_ARRAY_DBL; -- Coordinates to be transformed v_ipiv SYS.utl_nla_array_int := SYS.utl_nla_array_int(0,0,0,0); -- Geometry variables v_dims PLS_Integer; v_measure_dim PLS_Integer; v_ord PLS_Integer; v_sdo_point mdsys.sdo_point_type := NULL; v_trans_point mdsys.sdo_point_type; v_ordinates mdsys.sdo_ordinate_array := NULL; -- EXCEPTION HANDLING NULL_GEOMETRY EXCEPTION; NULL_PARAMETER EXCEPTION; -- Cursor over vertices CURSOR c_coordinates( p_geometry IN mdsys.sdo_geometry) IS SELECT v.* FROM TABLE(mdsys.sdo_util.GetVertices(p_geometry)) v; FUNCTION TransformPoint(p_x IN NUMBER, p_y IN NUMBER, p_z IN NUMBER) RETURN mdsys.sdo_point_type IS v_info INTEGER; v_point mdsys.sdo_point_type := mdsys.sdo_point_type(p_x,p_y,p_z); BEGIN v_C := SYS.UTL_NLA_ARRAY_DBL(p_x, p_y, CASE WHEN p_z IS NULL THEN 0 ELSE p_z END, -- Supply 0 instead of NULL as this will crash LAPACK_GESV 0); -- Vertices are transformed as follows: -- x' = a*x + b*y + c*z + xoff -- y' = d*x + e*y + f*z + yoff -- z' = g*x + h*y + i*z + zoff -- SYS.UTL_NLA.LAPACK_GESV ( n => 4, -- A number of rows and columns nrhs => 1, -- B number of columns a => v_A, -- matrix A lda => 4, -- max(1, n) ipiv => v_ipiv, -- pivot indices (set to zeros) b => v_C, -- matrix Result ldb => 4, -- ldb >= max(1,n) info => v_info, -- operation status (0=sucess) pack => 'C' -- how the matrices are stored -- (C=column-wise) ); IF (v_info = 0) THEN v_point.x := v_C(1); v_point.y := v_C(2); v_point.z := CASE WHEN p_z IS NULL THEN NULL ELSE v_C(3) END; -- Return correct value only if one supplied. ELSE raise_application_error( -20001, 'Matrix transformation by LAPACK_GESV failed with error ' || v_info, FALSE ); END IF; RETURN v_point; END TransformPoint; BEGIN IF ( p_geometry IS NULL ) THEN raise NULL_GEOMETRY; END IF; IF ( p_a IS NULL OR p_b IS NULL OR p_d IS NULL OR p_e IS NULL OR p_xoff IS NULL OR p_yoff IS NULL ) THEN raise NULL_PARAMETER; END IF; v_dims := TRUNC(p_geometry.sdo_gtype/1000,0); v_measure_dim := MOD(trunc(p_geometry.sdo_gtype/100),10); -- Transform any point in the geometry object IF ( p_geometry.sdo_point IS NOT NULL ) THEN v_sdo_point := TransformPoint(p_geometry.sdo_point.x, p_geometry.sdo_point.y, p_geometry.sdo_point.z); END IF; IF ( p_geometry.sdo_ordinates IS NOT NULL ) THEN v_ordinates := NEW mdsys.sdo_ordinate_array(1); v_ordinates.DELETE; v_ordinates.EXTEND(p_geometry.sdo_ordinates.COUNT); v_ord := 1; -- Loop around coordinates and apply matrix to them. <<for_all_coords>> FOR coord IN c_coordinates( p_geometry ) loop v_trans_point := TransformPoint(coord.x, coord.y, CASE WHEN v_measure_dim=3 THEN NULL ELSE coord.z END); v_ordinates(v_ord) := v_trans_point.x; v_ord := v_ord + 1; v_ordinates(v_ord) := v_trans_point.y; v_ord := v_ord + 1; IF ( v_dims >= 3 ) THEN v_ordinates(v_ord) := v_trans_point.z; v_ord := v_ord + 1; END IF; IF ( v_dims >= 4 ) THEN v_ordinates(v_ord) := coord.w; v_ord := v_ord + 1; END IF; END LOOP for_all_coords; END IF; RETURN mdsys.sdo_geometry(p_geometry.sdo_gtype, p_geometry.sdo_srid, v_sdo_point, p_geometry.sdo_elem_info, v_ordinates ); EXCEPTION WHEN NULL_GEOMETRY THEN raise_application_error(c_i_null_geometry,c_s_null_geometry,TRUE); RETURN p_geometry; WHEN NULL_PARAMETER THEN raise_application_error( c_i_null_parameter,c_s_null_parameter,TRUE ); RETURN p_geometry; END ST_Affine; END Affine; / SHOW errors
To test this package I created the following:
WITH testGeom AS ( SELECT mdsys.sdo_geometry(2002,NULL,NULL,SDO_ELEM_INFO_ARRAY(1,2,1),SDO_ORDINATE_ARRAY(2.0,2.0,2.0,4.0,8.0,4.0,12.0,4.0,12.0,10.0,8.0,10.0,5.0,14.0)) AS geom, SDO_GEOMETRY(3002,NULL,NULL,SDO_ELEM_INFO_ARRAY(1,2,1),SDO_ORDINATE_ARRAY(2.0,2.0,0.0, 2.0,4.0,22.22, 8.0,4.0,37.04, 12.0,4.0,59.26, 12.0,10.0,74.07, 8.0,10.0,92.59, 5.0,14.0,100.0)) AS geom3D, mdsys.sdo_geometry(2001,NULL,sdo_point_type(2.0,2.0,NULL),NULL,NULL) AS rotatePoint FROM dual ) SELECT Affine.ST_Rotate(a.geom, cogo.radians(45), -1, rotatePoint, NULL ) AS geom FROM testGeom a UNION ALL SELECT Affine.ST_Rotate(a.geom, cogo.radians(90)) AS geom FROM testGeom a UNION ALL SELECT Affine.ST_Rotate(a.geom, cogo.radians(135), 2.0, 2.0) AS geom FROM testGeom a UNION ALL SELECT Affine.ST_Rotate(a.geom, cogo.radians(225), rotatePoint) AS geom FROM testGeom a UNION ALL SELECT Affine.ST_Scale(a.geom,2,2,0) AS geom FROM testGeom a UNION ALL SELECT Affine.ST_Scale(a.geom,3,3) AS geom FROM testGeom a UNION ALL SELECT Affine.ST_Translate(a.geom,10,10) AS geom FROM testGeom a UNION ALL SELECT Affine.ST_Translate(a.geom,10,10,10) AS geom FROM testGeom a UNION ALL SELECT Affine.ST_RotateTranslateScale(a.geom,cogo.radians(45),rotatePoint,2,2,0,10,10,0) AS geom FROM testGeom a; -- Results -- GEOM ------------------------------------------------------------------------------- SDO_GEOMETRY(2002,NULL,NULL,SDO_ELEM_INFO_ARRAY(1,2,1),SDO_ORDINATE_ARRAY(2.0,2.0, 0.585786437626906,3.4142135623731, 4.8284271247462,7.65685424949238, 7.65685424949239,10.4852813742386, 3.41421356237311,14.7279220613579, 0.585786437626913,11.8994949366117, -4.36396103067892,12.6066017177982)) SDO_GEOMETRY(2002,NULL,NULL,SDO_ELEM_INFO_ARRAY(1,2,1),SDO_ORDINATE_ARRAY(-2.0,2.0, -4.0,2.00000000000001, -3.99999999999999,8.00000000000001, -3.99999999999998,12.0, -9.99999999999998,12.0, -9.99999999999999,8.00000000000002, -14.0,5.00000000000002)) SDO_GEOMETRY(2002,NULL,NULL,SDO_ELEM_INFO_ARRAY(1,2,1),SDO_ORDINATE_ARRAY(2.0,2.0, 0.585786437626901,0.585786437626908, -3.65685424949237,4.8284271247462, -6.48528137423856,7.6568542494924, -10.7279220613579,3.41421356237312, -7.89949493661167,0.585786437626927, -8.60660171779823,-4.3639610306789)) SDO_GEOMETRY(2002,NULL,NULL,SDO_ELEM_INFO_ARRAY(1,2,1),SDO_ORDINATE_ARRAY(2.0,2.0, 3.41421356237309,0.585786437626899, -0.828427124746214,-3.65685424949237, -3.65685424949242,-6.48528137423855, 0.585786437626852,-10.7279220613578, 3.41421356237305,-7.89949493661167, 8.36396103067888,-8.60660171779824)) SDO_GEOMETRY(2002,NULL,NULL,SDO_ELEM_INFO_ARRAY(1,2,1),SDO_ORDINATE_ARRAY(4.0,4.0, 4.0,8.0, 16.0,8.0, 24.0,8.0, 24.0,20.0, 16.0,20.0, 10.0,28.0)) SDO_GEOMETRY(2002,NULL,NULL,SDO_ELEM_INFO_ARRAY(1,2,1),SDO_ORDINATE_ARRAY(6.0,6.0, 6.0,12.0, 24.0,12.0, 36.0,12.0, 36.0,30.0, 24.0,30.0, 15.0,42.0)) SDO_GEOMETRY(2002,NULL,NULL,SDO_ELEM_INFO_ARRAY(1,2,1),SDO_ORDINATE_ARRAY(12.0,12.0, 12.0,14.0, 18.0,14.0, 22.0,14.0, 22.0,20.0, 18.0,20.0, 15.0,24.0)) SDO_GEOMETRY(2002,NULL,NULL,SDO_ELEM_INFO_ARRAY(1,2,1),SDO_ORDINATE_ARRAY(12.0,12.0, 12.0,14.0, 18.0,14.0, 22.0,14.0, 22.0,20.0, 18.0,20.0, 15.0,24.0)) SDO_GEOMETRY(2002,NULL,NULL,SDO_ELEM_INFO_ARRAY(1,2,1),SDO_ORDINATE_ARRAY(12.0,12.0, 9.17157287525381,14.8284271247462, 17.6568542494924,23.3137084989848, 23.3137084989848,28.9705627484771, 14.8284271247462,37.4558441227157, 9.17157287525383,31.7989898732233, -0.727922061357837,33.2132034355964))
I have compared the output of some of these to those in my GEOM package and they appear to give the correct results.
Now, comparing my ST_Affine function against this example from PostGIS 2.0, I get the same result.
PostGIS (edited)
--Rotate a 3d line 180 degrees about the z axis.
SELECT ST_AsEWKT(ST_Affine(the_geom, cos(pi()), -sin(pi()), 0, sin(pi()), cos(pi()), 0, 0, 0, 1, 0, 0, 0)) AS using_affine
FROM (SELECT ST_GeomFromEWKT('LINESTRING(1 2 3, 1 4 3)') AS the_geom) AS foo;
using_affine
-----------------------------
LINESTRING(-1 -2 3,-1 -4 3)
(1 ROW)
My ST_Affine Function
-- Documented test for PostGIS -- WITH geom3d AS ( SELECT sdo_geometry(3002,NULL,NULL,sdo_elem_info_array(1,2,1),sdo_ordinate_array(1,2,3,1,4,3)) AS the_geom FROM dual ) SELECT AFFINE.ST_Affine(p_geometry=>the_geom, p_a=>cos(CONSTANTS.pi()), p_b=>0-sin(CONSTANTS.pi()), p_c=>0, p_d=>sin(CONSTANTS.pi()), p_e=>cos(CONSTANTS.pi()), p_f=>0, p_g=>0, p_h=>0, p_i=>1, p_xoff=>0, p_yoff=>0, p_zoff=>0) AS using_affine FROM geom3d a; -- Results -- USING_AFFINE ----------------- SDO_GEOMETRY(3002,NULL,NULL,SDO_ELEM_INFO_ARRAY(1,2,1),SDO_ORDINATE_ARRAY(-0.999999999999994,-2.0,3.0, -0.999999999999987,-4.0,3.0))
As you can see they are the same.
I think this AFFINE package is a start to something more helpful than the current one-size-fits-all approach.
Any fixes, corrections, improvements, email me and I will modify the package.
I hope this is of help to somone.