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)
ST_VertexN / ST_PointN – Extracting a specific point from any geometry
I have published a number of times and implementation of functions what would extract a specific point from any geometry object. These appear in my LINEAR PL/SQL package and the T_GEOMETRY object that was published as part of the PACT book I wrote with Dr Siva Ravada.
When working for a customer recently I had extracted the ST_VertexN, ST_StartVertex and ST_EndVertex functions from the source of the PACKT book for use by a process that needed to extract a point.
After having deployed it, I accidently deleted what I had done. For specific reasons I won’t go in to I could not access the original source, so I re-wrote the main ST_VertexN function from scratch. I publish that source here for no other reason than as a reflection on the way that algorithms or implementations can change.
FUNCTION ST_VertexN(p_geom IN mdsys.sdo_geometry, p_vertex IN INTEGER) RETURN mdsys.sdo_geometry IS c_i_invalid_vertex Constant pls_integer := -20123; c_s_invalid_vertex Constant VarChar2(100) := 'Vertex position (*POSN*) is invalid.'; v_vertex pls_integer := NVL(p_vertex,-1); -- Index into coordinates array v_vertices mdsys.vertex_set_type; v_num_vertices pls_integer; v_dims pls_integer; BEGIN IF (p_geom IS NULL) THEN NULL; END IF; v_dims := p_geom.get_dims(); v_vertices := mdsys.sdo_util.getVertices(p_geom); v_num_vertices := v_vertices.COUNT; -- Check / Set vertex request -- IF ( v_vertex < 0 ) THEN v_vertex := v_num_vertices - (ABS(v_vertex)-1); END IF; IF ( v_vertex = 0 OR ABS(v_vertex) > v_num_vertices ) THEN raise_application_error(c_i_invalid_vertex, REPLACE(c_s_invalid_vertex,'*POSN*',to_char(v_vertex)),TRUE); END IF; RETURN CASE WHEN (p_geom.sdo_point IS NOT NULL AND p_geom.sdo_ordinates IS NULL AND v_vertex IN (-1,1) ) THEN p_geom; ELSE SDO_GEOMETRY((v_dims*1000)+1, p_geom.sdo_srid, CASE WHEN ( v_dims <= 3 ) THEN mdsys.sdo_point_type(v_vertices(v_vertex).x, v_vertices(v_vertex).y, CASE WHEN (v_dims=3) THEN v_vertices(v_vertex).y ELSE NULL END) ELSE NULL END, CASE WHEN ( v_dims <= 3 ) THEN NULL ELSE mdsys.sdo_elem_info_array(1,1,1) END, CASE WHEN ( v_dims <= 3 ) THEN NULL ELSE mdsys.sdo_ordinate_array(v_vertices(v_vertex).x,v_vertices(v_vertex).y,v_vertices(v_vertex).z,v_vertices(v_vertex).w) END ); END ST_VertexN; FUNCTION ST_StartVertex(p_geom IN mdsys.sdo_geometry) RETURN mdsys.sdo_geometry IS BEGIN RETURN ST_VertexN(p_geom,1); END ST_StartVertex; FUNCTION ST_EndVertex(p_geom IN mdsys.sdo_geometry) RETURN mdsys.sdo_geometry IS BEGIN RETURN ST_VertexN(p_geom,-1); END ST_EndVertex;
The code had come from the PACT book. Now, I would not normally use the mdsys.vertex_set_type in PL/SQL. Why? Because it involves:
- An addition pass over the whole array;
- The creation of a new memory object that holds a copy of the original sdo_ordinates.
What if the person had only asked for the first point? It would be faster just to access the sdo_ordinate_array by converting vertex references to ordinate references.
I needed to get this written quickly to I decided to use the mdsys.vertex_set_type approach as it involved no conversion of the input p_vertex number to ordinates. The main saving, in the end, was that in this implementation I went straight from checking the input p_vertex value to returning the sdo_geometry where previously I had used IF statements to check if the geometry was a point, whether it was encoded in sdo_point_type or sdo_ordinates etc. This checking is done via case statements.
Actually, I think the new implementation looks fresher and cleaner and easier to debug.
Update: I had a bit of time, so I modified the ST_VertexN code so that it no longer needed mdsys.vertex_set_type processing directly against the sdo_ordinate_array. I decided to call this function ST_PointN as follows:
CREATE OR REPLACE FUNCTION ST_PointN(p_geom IN mdsys.sdo_geometry, p_Point IN INTEGER) RETURN mdsys.sdo_geometry IS c_i_invalid_Point Constant pls_integer := -20123; c_s_invalid_Point Constant VarChar2(100) := 'Point position (*POSN*) is invalid.'; v_Point pls_integer := NVL(p_Point,-1); -- Index into coordinates array v_num_points pls_integer; v_dims pls_integer; v_ord pls_integer; BEGIN IF (p_geom IS NULL) THEN NULL; END IF; v_dims := p_geom.get_dims(); v_num_points := CASE WHEN p_geom.sdo_ordinates IS NULL THEN 1 ELSE p_geom.sdo_ordinates.COUNT / v_dims END; -- Check / Set Point request -- IF ( v_Point < 0 ) THEN v_Point := v_num_points - (ABS(v_Point)-1); END IF; IF ( v_Point = 0 OR ABS(v_Point) > v_num_points ) THEN raise_application_error(c_i_invalid_Point, REPLACE(c_s_invalid_Point,'*POSN*',to_char(v_Point)),TRUE); END IF; v_ord := ( v_point - 1 ) * v_dims + 1; RETURN CASE WHEN (p_geom.sdo_point IS NOT NULL AND p_geom.sdo_ordinates IS NULL AND v_Point IN (-1,1) ) THEN p_geom ELSE SDO_GEOMETRY((v_dims*1000)+1, p_geom.sdo_srid, CASE WHEN ( v_dims <= 3 ) THEN mdsys.sdo_point_type(p_geom.SDO_ORDINATES(v_ord), p_geom.SDO_ORDINATES(v_ord+1), CASE WHEN (v_dims=3) THEN p_geom.SDO_ORDINATES(v_ord+2) ELSE NULL END) ELSE NULL END, CASE WHEN ( v_dims <= 3 ) THEN NULL ELSE mdsys.sdo_elem_info_array(1,1,1) END, CASE WHEN ( v_dims <= 3 ) THEN NULL ELSE mdsys.sdo_ordinate_array( p_geom.SDO_ORDINATES(v_ord), p_geom.SDO_ORDINATES(v_ord+1), p_geom.SDO_ORDINATES(v_ord+2), p_geom.SDO_ORDINATES(v_ord+3)) END ) END; END ST_PointN; / SHOW errors -- Now create start point wrapper -- CREATE OR REPLACE FUNCTION ST_StartPoint(p_geom IN mdsys.sdo_geometry) RETURN mdsys.sdo_geometry IS BEGIN RETURN ST_PointN(p_geom,1); END ST_StartPoint; / SHOW errors -- now create end point wrapper -- CREATE OR REPLACE FUNCTION ST_EndPoint(p_geom IN mdsys.sdo_geometry) RETURN mdsys.sdo_geometry IS BEGIN RETURN ST_PointN(p_geom,-1); END ST_EndPoint; / SHOW errors
Some tests.
WITH geom AS ( SELECT SDO_GEOMETRY(2002, 8307, NULL, SDO_ELEM_INFO_ARRAY(1, 2, 1), SDO_ORDINATE_ARRAY(59.8833333, 25.8666667, 59.975, 26.1116667, 60.1333333, 26.55, 60.205, 27.3, 60.2166667, 27.4666667, 60.5494444, 27.8016667, 60.555, 27.8369444)) AS line FROM dual ) SELECT st_startPoint(line) AS point FROM geom UNION ALL SELECT st_endPoint(line) AS point FROM geom; -- POINT ---------------------------------------------------------------------------- SDO_GEOMETRY(2001,8307,SDO_POINT_TYPE(59.8833333,25.8666667,NULL),NULL,NULL) SDO_GEOMETRY(2001,8307,SDO_POINT_TYPE(60.555,27.8369444,NULL),NULL,NULL) -- -- WITH geom AS ( SELECT SDO_GEOMETRY(2002, 8307, NULL, SDO_ELEM_INFO_ARRAY(1, 2, 1), SDO_ORDINATE_ARRAY(59.8833333, 25.8666667, 59.975, 26.1116667, 60.1333333, 26.55, 60.205, 27.3, 60.2166667, 27.4666667, 60.5494444, 27.8016667, 60.555, 27.8369444)) AS line FROM dual ) SELECT level AS pointN, st_PointN(line,level) AS point FROM geom a CONNECT BY level <= sdo_util.getNumVertices(a.line); -- -- POINTN POINT ------ --------------------------------------------------------------------------- 1 SDO_GEOMETRY(2001,8307,SDO_POINT_TYPE(59.8833333,25.8666667,NULL),NULL,NULL) 2 SDO_GEOMETRY(2001,8307,SDO_POINT_TYPE(59.975,26.1116667,NULL),NULL,NULL) 3 SDO_GEOMETRY(2001,8307,SDO_POINT_TYPE(60.1333333,26.55,NULL),NULL,NULL) 4 SDO_GEOMETRY(2001,8307,SDO_POINT_TYPE(60.205,27.3,NULL),NULL,NULL) 5 SDO_GEOMETRY(2001,8307,SDO_POINT_TYPE(60.2166667,27.4666667,NULL),NULL,NULL) 6 SDO_GEOMETRY(2001,8307,SDO_POINT_TYPE(60.5494444,27.8016667,NULL),NULL,NULL) 7 SDO_GEOMETRY(2001,8307,SDO_POINT_TYPE(60.555,27.8369444,NULL),NULL,NULL) -- 7 ROWS selected -- -- WITH geom AS ( SELECT SDO_GEOMETRY(2001, 8307, sdo_point_type(59.8833333, 25.8666667, NULL),NULL,NULL) AS point FROM dual ) SELECT st_startPoint(point) AS sp, st_endPoint(point) AS ep FROM geom ; -- POINT ---------------------------------------------------------------------------- SDO_GEOMETRY(2001,8307,SDO_POINT_TYPE(59.8833333,25.8666667,NULL),NULL,NULL) SDO_GEOMETRY(2001,8307,SDO_POINT_TYPE(59.8833333,25.8666667,NULL),NULL,NULL)
I hope this is useful to someone out there.