SpatialDB Advisor

Multi-Centroid Shootout · 83 days ago by Simon Greener

This blog follows on from the first “Centroid Shoot Out” and the comment left by Andy. It covers two issues. Firstly it publishes the original centroid code written by Tino Delbourgo when he worked for a company based in Hobart, Tasmania called Geometry Pty Ltd ; secondly, it shows how the PL/SQL version of the Java code handles holes in polygons, multi-part polygons and the generation of multi-point centroids.

Original Java Code.

The original Java code is:

private static void calcParaCentroid(Geometry geometry, Double centroidX, Double centroidY)
{
   Enumeration segs = ((Polygon)geometry).getAllSegments();
   int segCount = 0;
   Envelope envp = geometry.getEnvelope();
   // Go half-way along bottom edge of bounding box to get our candidate point (x,y)
   double x = (envp.getMinX()+envp.getMaxX())/2.0;
   double y = envp.getMinY();
   // Go through each line segment in turn, working out crossing points from the
   // half-way along the bottom edge of the bounding box in a line due north.
   // We only need to keep the lowest and second lowest of these points.
   double lowestCrossing = 0.0;
   double secondLowestCrossing = 0.0;
   while (segs.hasMoreElements()) {
      segCount++;
      if (debugPinP)
         System.out.println("Considering polygon segment #"<ins>segCount);
      LinearSegment seg = (LinearSegment)segs.nextElement();
      double coords[] = seg.getCoordArray();
      double x1 = coords<sup><a href="#fn89740777947bffc32601fb">0</a></sup>;
      double y1 = coords<sup><a href="#fn98879393047bffc32609ca">1</a></sup>;
      double x2 = 0.0;
      double y2 = 0.0;
      for (int i=2; i<coords.length; i</ins>=2) {
         double ycrossing = 0.0;
         x2 = coords[i];
         y2 = coords[i+1];
         if (debugPinP)
            System.out.println("   Considering line segment: ("<ins>x1</ins>","<ins>y1</ins>") to ("<ins>x2</ins>","<ins>y2</ins>")");
         if ((x1<x && x2<x) || (x1>x && x2>x) || (y1<y && y2<y)) {
            // do nothing - segment is either wholly to left, to right or below our point
         } else if (y1>=y && y2>=y) {
            // the line segment is above our point with one end to the left and the other to the right
            // - this is a definite crossing
            if ((x1<=x && x<=x2) || (x2<=x && x<=x1))
               ycrossing = y1+((y2-y1)/(x2-x1))*(x-x1);
         } else if (x1==x && x2==x) {
            // do nothing - we are comparing with a vertical line (no crossing possible)
         } else if (y1==y && y2==y) {
            // we are dealing with a horizontal line going through our point - a definite crossing
            ycrossing = y;
         } else {
            // we are dealing with a line above our point with our point being within the line's bounding box
            // - we need to project up to the line to see if it crosses
            double ytest = y1+((y2-y1)/(x2-x1))*(x-x1);
            if (ytest>y)
               ycrossing = ytest;
         }

         // If it crosses, then see if it is the crossing with the lowest or second lowest y-value
         if (ycrossing!=0.0) {
          if ( ycrossing  lowestCrossing ) {
            lowestCrossing = ycrossing;
            } else if ( ycrossing  secondLowestCrossing ) {
              secondLowestCrossing = ycrossing;
       } else if (lowestCrossing==0.0 || 
                  ycrossing <  lowestCrossing) {
              secondLowestCrossing = lowestCrossing;
              lowestCrossing = ycrossing;
            } else if (secondLowestCrossing==0.0 || 
                       ycrossing < secondLowestCrossing) {
               secondLowestCrossing = ycrossing;
            }
         }
         // Go round the loop again with next segment
         x1 = x2;
         y1 = y2;
      }
   }
   if (secondLowestCrossing==0.0 || lowestCrossing==0.0)
      throw new RuntimeException(“Para-centroid calculation failed, couldn’t find two crossing points up from centre bottom edge of bounding box”);
   centroidX = new Double(x);
   centroidY = new Double((lowestCrossing+secondLowestCrossing)/2.0);
}

You will note that the code requires some functions in the old sdoapi.jar from 9i (the api was changed with the 10g release). The function it needs from this jar file is one that “vectorizes” an sdo_geometry polygon into simple start/end 2 point linestrings (or vectors): getAllSegments().

I have looked at deploying this code inside the Oracle database JVM with PL/SQL wrappers. It is possible but because the sdoapi.jar file is deprecated I decided the only long term solution was to use the Java Topology Suite and to build a vectorisation function using it. I have not time to do this. There are also other issues such as the conversion of Sdo_Geometry/JGeometry to JTS via GeoTools (I think that the current converter does not support circular arcs).

When I converted the Java to PL/SQL, I had to build my own vectoriser (see my article on spatial pipelining). I also added in support for rectangles, circles and compound linestrings/polygons with three point circular arcs. I added the ability to choose the largest part of a multi-part geometry object (SdoGType x007). All these could be added to the Java code but have not been done. In creating the PL/SQL version I found a few centroid use cases that the original Java code didn’t implement. I have added these to the PL/SQL version (compiled and tested) and the Java version (never compiled or tested).

Multi-Part Polygons: Single Centroid

As indicated in the previous paragraph the PL/SQL implementation supports generating a centroid for multi-part polygons (SDO_GTYPE = x007). The method I use is simple: firstly, choose the largest part; secondly, generate the centroid as normal. The algorithm used gets the minimum bounding rectangle of each outer shell and selects the largest. Of course, this may still not be correct as in certain situations area would be better. Perhaps one day I will add this in.

Multi-Part Polygons: Multiple Centroids

The PL/SQL version has support for generating a multi-point sdo_geometry (SDO_GTYPE = x005) holding the centroid of each and every outer shell (EType 1003) of a multi-part polygon (Sdo_Gtype x007). The algorithm simple extracts each outer shell as a single (2003) polygon and passes it to the standard centroid function: each centroid created is appended to a multi-point sdo_geometry object. Once all parts have been processed the resultant multi-point geometry is returned.

Holes

There is little that I can say other than the code will never put a centroid into any of the holes (SDO_GTYPE = 2003) of any part of a polygon of any type.

Data and Diagram

Putting it all together. I have created a single, complex, multi-part polygon (SDO_GTYPE = 2007) object (see below). For this polygon I have generated a centroid from the standard Oracle function sdo_geom.sdo_centroid, my own geom.sdo_centroid function and finally I have generated a single multi-point geometry via my geom.sdo_multi_centroid function.

codesys@XE> CREATE TABLE MultiCentroidShootOut AS
SELECT 
SDO_GEOMETRY(2007, 28355, NULL, 
SDO_ELEM_INFO_ARRAY(
 1, 1003, 1, 
25, 2003, 1, 
35, 2003, 1, 
45, 2003, 1, 
63, 1003, 1, 
87, 1003, 1, 
99, 1003, 1), 
SDO_ORDINATE_ARRAY(
17.0158371, -492.10068, 
60.6568627, -602.4868, 
186.445701, -612.75528, 
378.979638, -712.87293, 
740.943439, -712.87293, 
669.064103, -515.20475, 
835.926848, -284.16403, 
879.567873, -140.40535, 
645.96003, 16.188914, 
299.398944, 13.6217949, 
76.0595777, -171.21078, 
17.0158371, -492.10068, 
386.680995, -576.81561, 
396.949472, -366.31184, 
602.319005, -489.53356, 
566.379336, -651.26207, 
386.680995, -576.81561, 
91.4622926, -474.13084, 
168.475867, -399.68439, 
307.100302, -481.8322, 
181.311463, -558.84578, 
91.4622926, -474.13084, 
176.177225, -209.71757, 
299.398944, -86.495852, 
561.245098, -30.019231, 
692.168175, -94.19721, 
789.718703, -171.21078, 
633.124434, -340.64065, 
376.412519, -255.92572, 
214.684012, -361.1776, 
176.177225, -209.71757, 
69, 9.5, 
206, 86.5, 
397, 185.5, 
553, 189.5, 
698, 143.5, 
920, -7.5, 
853, 105.5, 
704, 259.5, 
537, 307.5, 
403, 271.5, 
183, 134.5, 
69, 9.5, 
412.352187, -158.37519, 
468.828808, -225.12029, 
592.050528, -237.95588, 
661.362745, -202.01621, 
568.946456, -89.062971, 
412.352187, -158.37519, 
886, -424.5, 
999, -357.5, 
1153, -208.5, 
1201, -41.5, 
1165, 92.5, 
1028, 312.5, 
903, 426.5, 
980, 289.5, 
1079, 98.5, 
1083, -57.5, 
1037, -202.5, 
886, -424.5)) as geom
from dual;

codesys@XE> select sdo_geom.sdo_centroid(geom,0.05) from MultiCentroidShootout;
SDO_GEOM.SDO_CENTROID(GEOM,0.05)(SDO_GTYPE, SDO_SRID, SDO_POINT(X, Y, Z), SDO_ELEM_INFO, SDO_ORDINATES)
-----------------------------------------------------------------------------------------------------------------------------------
SDO_GEOMETRY(2001, 28355, SDO_POINT_TYPE(545.679984, -227.04292, NULL), NULL, NULL)

codesys@XE> select geom.sdo_centroid(geom,0.05) from MultiCentroidShootout;
GEOM.SDO_CENTROID(GEOM,0.05)(SDO_GTYPE, SDO_SRID, SDO_POINT(X, Y, Z), SDO_ELEM_INFO, SDO_ORDINATES)
-----------------------------------------------------------------------------------------------------------------------------------
SDO_GEOMETRY(2001, 28355, SDO_POINT_TYPE(448.3, -657.6, NULL), NULL, NULL)

codesys@XE> select geom.sdo_multi_centroid(geom,0.05) from MultiCentroidShootout;
GEOM.SDO_MULTI_CENTROID(GEOM,0.05)(SDO_GTYPE, SDO_SRID, SDO_POINT(X, Y, Z), SDO_ELEM_INFO, SDO_ORDINATES)
-----------------------------------------------------------------------------------------------------------------------------------
SDO_GEOMETRY(2005, 28355, NULL, SDO_ELEM_INFO_ARRAY(1, 1, 4), SDO_ORDINATE_ARRAY(448.3, -657.6, 494.5, 242.1, 536.9, -167.8, 1043.5, -248.2))

Hope all this helps in understanding how the centroid functions in my PL/SQL packages work.

Comment

Spatial Pipelining · 94 days ago by Simon Greener

This technical blog article describes the benefits in using pipelined functions in Oracle to manipulate sdo_geometry objects. In doing so it will describe a number of functions available in the COGO and GEOM packages in the PL/SQL packages I make available for a free download from this site.

What I will do is introduce a realistic “business need” and then show how to construct a non-pipelined function that can be used in implementing that need. I will also create a pipelined version of the function and “compare and contrast” the two approaches in terms of memory use and performance.

Business Need

Imagine a company has a large Oracle Spatial database in which are stored land parcel (land record) polygons. The company would like to be able to display the bearings and distances (metes and bounds) of each boundary of each polygon dynamically by not relying on a second layer of sdo_geometry linestrings. This is displayed pictorially as follows.

Business Requirement: Dynamic bearings and distances

How can we achieve this?

Steps to create table

First let’s start by creating a table and populating it.

CREATE TABLE land_parcel (
 gid  INTEGER,
 geom SDO_GEOMETRY 
);

INSERT INTO land_parcel 
VALUES(1,
       SDO_GEOMETRY(2003, NULL, NULL, 
                    SDO_ELEM_INFO_ARRAY(1,1003,1), 
                    SDO_ORDINATE_ARRAY(100,0,400,0,400,150,250,100,250,200,400,150,400,300,100,300,100,0)));

INSERT INTO user_sdo_geom_metadata 
VALUES('LAND_PARCEL','GEOM',SDO_DIM_ARRAY(SDO_DIM_ELEMENT('X',0,1000,0.05), SDO_DIM_ELEMENT('Y',0,1000,0.05)),NULL);
COMMIT;

Vector Elements

Now we need a method for accessing the individual lines that make up our sdo_geometry polygon. For this we will need a function that takes the polygon and splits it up into its constituent vectors (where a vector is defined as being a linestring made up of only a starting and ending vertex).

First we create a vector data structure as follows:

CREATE OR REPLACE TYPE Coord2DType AS OBJECT (
  x NUMBER,
  y NUMBER 
);

CREATE OR REPLACE TYPE Vector2DType AS OBJECT (
  startCoord Coord2DType,
  endCoord   Coord2DType 
);

Once defined we can now create a set of vectors (to hold all those in any one land parcel polygon) as follows:

CREATE OR REPLACE TYPE Vector2DSetType AS TABLE OF Vector2DType;

Finally, having created our data structures we can create our function. The one I created in my GEOM package is GetVector2D. Its essentials (for a full implementation see my packages) are:

FUNCTION GetVector2D ( p_geometry IN mdsys.sdo_geometry)
  RETURN CODESYS.Vector2DSetType DETERMINISTIC;
    vectors Vector2DSetType := Vector2DSetType();
BEGIN
  ...
  WHILE v_partToProcess Loop
    ...
    IF v_vertex = 1 THEN
      vectors.EXTEND;
      v_vector := vectors.LAST;
      vectors(v_vector) := Vector2DType(Coord2DType(-1,1),Coord2DType(-1,1));
      vectors(v_vector).startCoord.x := v_coord.x;
      vectors(v_vector).startCoord.y := v_coord.y;
    ELSE
      vectors(v_vector).endCoord.x := v_coord.x;
      vectors(v_vector).endCoord.y := v_coord.y;
      vectors.EXTEND;
      v_vector := vectors.LAST;
      vectors(v_vector) := Vector2DType(Coord2DType(-1,1),Coord2DType(-1,1));
      vectors(v_vector).startCoord.x := v_coord.x;
      vectors(v_vector).startCoord.y := v_coord.y;
    END IF;
    ...
  END LOOP;
  ...
  RETURN vectors;
END;

Finally, we will need two functions that, given a single vector, can return a bearing and distance. I won’t go into the details of how to do this, all I will do is point out that, in my COGO package, are two functions:

CREATE OR REPLACE PACKAGE COGO
AS
    FUNCTION Bearing( dE1 in number,
                      dN1 in number,
                      dE2 in number,
                      dN2 in number)
    RETURN NUMBER DETERMINISTIC;
    FUNCTION Distance( dE1 in number,
                       dN1 in number,
                       dE2 in number,
                       dN2 in number)
    RETURN NUMBER DETERMINISTIC;
...
END COGO;

(In the CONSTANTS package is a definition of PI which we will also use.)

View

Now we can construct a view that will take a land parcel (or set of land parcels) and return the vectors that compose it. From these vectors we will create sdo_geometry linestrings and also columns containing the bearings and distances computed from those vectors.

CREATE OR REPLACE VIEW metes_and_bounds
AS
SELECT rownum AS gid,
       codesys.Cogo.DD2DMS(
              codesys.Cogo.Bearing(startx,starty,endx,endy)
              *
              (180/codesys.Constants.PI) )
              AS bearing,
       ROUND(codesys.Cogo.Distance(startx,starty,endx,endy),2)
              AS distance,
       MDSYS.sdo_geometry(2002,NULL,NULL,
              MDSYS.SDO_ELEM_INFO_ARRAY(1,2,1),
              MDSYS.SDO_ORDINATE_ARRAY(startx,startY,endX,endY))
              AS geometry
  FROM ( SELECT DISTINCT c.StartCoord.X AS startX,
                         c.StartCoord.Y AS startY,
                         c.EndCoord.X AS endX,
                         c.EndCoord.Y AS endY
           FROM land_parcel a,
                TABLE(codesys.Geom.GetVector2D(a.geom)) c
       );

Note that all the “heavy lifting” is done by the GetVector2D function that returns a set of vectors into the TABLE disaggregator.

Finally, if the land_parcel table held a lot of records, we should really consider using a fast refreshable materialized view instead and create a spatial index over its geom column. This is because we cannot create a function based index over the above view.

Pipelining and Pipelined Functions

OK, so we can now implement our business requirement. Why do we need to go further?

We should always be on the lookout for performance improvements even in these days of fast computer hardware. This is because, in any organisation, database servers are shared across multiple databases and accessed by multiple applications and users. Anything we can do to reduce our application “footprint” helps all who use the resource. Pipelining is one method of improving the performance and memory use of a function.

Firstly, what is a pipelined function?

Tom Kyte says in this article that:

Pipelined functions are simply code that you can pretend is a database table.
Pipelined functions give you the (amazing, to me) ability to use

    SELECT * FROM <PLSQL_FUNCTION>;

...
A pipelined function needs to return a collection type …

In fact it is

SELECT * FROM TABLE(<PLSQL_FUNCTION>)

In the metes_and_bounds view above we have used a function that returns a collection type (Vector2DSetType) inside a TABLE function. This is allowed, but this is still not a pipelined function!

Let’s turn it into a pipelined function and then I can explain the difference.

FUNCTION GetVector2D ( p_geometry IN mdsys.sdo_geometry)
  RETURN CODESYS.Vector2DSetType PIPELINED;
    v_vector         codesys.Vector2DType := codesys.Vector2DType(
                             codesys.Coord2DType(codesys.Constants.c_MinVal,codesys.Constants.c_MinVal),
                             codesys.Coord2DType(codesys.Constants.c_MinVal,codesys.Constants.c_MinVal));
BEGIN
  ...
  WHILE v_partToProcess Loop
    ...
    IF v_vertex = 1 THEN
      v_vector.startCoord.x := v_coord.x;
      v_vector.startCoord.y := v_coord.y;
    ELSE
      v_vector.endCoord.x := v_coord.x;
      v_vector.endCoord.y := v_coord.y;
      PIPE ROW(v_vector);
      v_vector := vectors.LAST;
      v_vector.startCoord.x := v_coord.x;
      v_vector.startCoord.y := v_coord.y;
    END IF;
    ...
  END LOOP;
  ...
  RETURN;
END;

The three points of difference are:

The definition of the function uses the PIPELINED keyword rather than the DETERMINISTIC keyword.
Instead of allocating memory as in the use of EXTEND on the private collection variable vectors (EXEND will allocate memory) in the pipelined function a single vector is simply pushed into the pipeline as it is constructed.
In the pipelined function no variable is returned using the final RETURN statement; in the non-pipelined function the local collection variable, vectors, which holds all vector objects constructed during function execution are return as one.

So, fairly obviously, a function is pipelined if it uses the PIPELINED keyword and pushes what it creates into the pipe via the PIPE ROW statement as they are created.

The big thing to notice is that very little memory is created or used by the pipelined function whereas the non-pipelined function has to allocated memory to hold all the vector objects until the function’s end. Because the non-piplined function waits until the end before it can return its results, the calling SELECT statement itself if forced to wait before it can do any processing: not so with the pipelined function. The Oracle 10gR2 help confirms this:

Rows from a collection returned by a table function can also be pipelined, that is, iteratively returned as they are produced instead of in a batch after all processing of the table function’s input is completed.

That help also outlines the benefits of pipelining:

Streaming, pipelining, and parallel execution of table functions can improve performance:

By enabling multithreaded, concurrent execution of table functions
By eliminating intermediate staging between processes
By improving query response time: With non-pipelined table functions, the entire collection returned by a table function must be constructed and returned to the server before the query can return a single result row. Pipelining enables rows to be returned iteratively, as they are produced. This also reduces the memory that a table function requires, as the object cache does not need to materialize the entire collection.
By iteratively providing result rows from the collection returned by a table function as the rows are produced instead of waiting until the entire collection is staged in tables or memory and then returning the entire collectio

Demonstration of Performance Improvements

With the land_parcel table above we don’t have enough objects to demonstrate the performance improvements the pipelined function has over the non-pipelined function. So, I used some customer data (I have permission for this).

SELECT count(*)
  FROM parcel;

  COUNT(*)
----------
     57453

CREATE TABLE pipelined_version
AS
SELECT rownum AS gid,
       MDSYS.SDO_GEOMETRY(2002,NULL,NULL,
       MDSYS.SDO_ELEM_INFO_ARRAY(1,2,1),
       MDSYS.SDO_ORDINATE_ARRAY(startx,startY,endX,endY))
       AS geometry
  FROM ( SELECT DISTINCT 
                c.StartCoord.X AS startX,
                c.StartCoord.Y AS startY,
                c.EndCoord.X AS endX,
                c.EndCoord.Y AS endY
           FROM ( SELECT geometry
                    FROM SP_PARCEL
                ) a,
                TABLE(codesys.Geom.Get{Piped}Vector2D(a.geometry)) c
       );

Note that I run this statement twice with the braces {} removed.

SELECT COUNT(*)
  FROM PIPELINED_VERSION;

  COUNT(*)
----------
    763916

Performance numbers were:

Function	TimeInSeconds
GetVector2D	02:18.13
GetPipedVector2D	00:47.90

In summary, pipelining improved the performance of the operation by 287% ( ( 1 / ( 48 / 138 ) * 100 ). Now that is quite an improvement even for this small amount of data. I have used this GetPipedVector2D function on some huge spatial datasets and have seen such substantial performance improvements I now use them in preference to previous techniques (I have left the previous implementations in my PL/SQL packages as occasionally they are useful).

I hope this little article helped whet your appetite for pipelined table functions in Oracle Spatial.

Comment

Using Oracle's SDO_NN Operator - Some examples · 103 days ago by Simon Greener

This little article was occasioned by someone emailing me and asking:

Actually I do have a topic for you – clustering. Before we moved to Oracle Spatial, we had a request for a client to: “Show us all the stores that are close to each other”. I’m somewhat embarrassed to say that we did this by coordinate match.

What would be nice is to understand how to implement the query: Show me all stores within x distance of each other.

Well, there is nothing to be embarrassed about because you did this before you had Oracle Spatial (actually all the functionality used below is available to Locator users – thanks Oracle)! What I will do is outline the Oracle functionality that can do what the person wanted and give a practical example I used recently to solve a particular problem.

First off the function we need in Oracle is the SDO_NN spatial operator. (My preference for everything I do in Oracle is to try and push as much processing into a SQL statement before launching into PL/SQL ie programming.) The Oracle documentation on SDO_NN is quite thorough so I recommend that you read it in consort with this article: For example, the only thing I will say say about the SDO_NN_DISTANCE ancillary operator (used in the first example) is that it returns the actual distance Oracle computed between the searched for object and the search object supplied to the SDO_NN operator. The examples below will clarify its use. For a detailed discussion on the operators read the documentation. For now, I will limit myself to a few comments and examples.

SDO_BATCH_SIZE vs SDO_NUM_RES

One of the first things to understand is use of the sdo_batch_size and sdo_num_res parameters. (Again the documentation is quite thorough on these parameters.)

sdo_num_res=<N>

simply returns exactly N nearest objects to the search object. So, in the example that follows the nearest object only is returned.

gis@XE> select id,
  2            storetype,
  3            sdo_nn_distance(1)
  4       from store s
  5      where sdo_nn(s.geom,
  6                   mdsys.sdo_geometry(2003,NULL,NULL,
  7                         mdsys.sdo_elem_info_array(1,1003,3),
  8                         mdsys.sdo_ordinate_array(380004,5100003,390000,5160000)),
  9*                  'sdo_num_res=1',1) = 'TRUE'
gis@XE> /

ID STORETYPE  SDO_NN_DISTANCE(1)
---------- ---------- ------------------
       271 SHOE                215195.12

(The script to create the STORE table can be accessed here.)

Note that the store is a SHOE store. Let’s “widen” the search a little (note the use of the ORDER BY clause):

gis@XE> select id,
  2         storetype,
  3         sdo_nn_distance(1)
  4    from store s
  5   where sdo_nn(s.geom,
  6                mdsys.sdo_geometry(2003,NULL,NULL,
  7                      mdsys.sdo_elem_info_array(1,1003,3),
  8                      mdsys.sdo_ordinate_array(380004,5100003,390000,5160000)),
  9                'sdo_num_res=5',1) = 'TRUE'
 10*  order by 3
gis@XE> /

        ID STORETYPE  SDO_NN_DISTANCE(1)
---------- ---------- ------------------
       271 SHOE                215195.12
       100 COMPUTER           215228.103
        60 SHOE                215307.04
       494 VIDEO              215341.161
       493 FURNITURE          215347.298

So, what if we only want the nearest three SHOE stores to our search point? Can we just add the “storetype = ‘SHOE’”?

gis>XE>  select id,
  2         storetype,
  3         sdo_nn_distance(1)
  4    from store s
  5   where sdo_nn(s.geom,
  6                mdsys.sdo_geometry(2003,NULL,NULL,
  7                      mdsys.sdo_elem_info_array(1,1003,3),
  8                      mdsys.sdo_ordinate_array(380004,5100003,390000,5160000)),
  9                'sdo_num_res=5',1) = 'TRUE'
 10     and s.storetype = 'SHOE'
 11*  order by 3
gis@XE> /

        ID STORETYPE  SDO_NN_DISTANCE(1)
---------- ---------- ------------------
       271 SHOE                215195.12
        60 SHOE                215307.04

No this will not do what we want. Why? Because the sdo_num_res parameter returns the objects based on their geometry not their attributes: the 5 stores returned by SDO_NN with sdo_num_res=5 are the 5 nearest. No more are returned. These 5 objects are then passed to the added predicate resulting in only 2 of the 5 passing the test. Can we move the predicate around? No because the SDO_NN will still only return the same 5 objects to the SELECT statement. We have to find a way to increase the set that can be returned by SDO_NN. Why not just increase the sdo_num_res to, say, 10?

gis@XE> select id,
  2         storetype,
  3         sdo_nn_distance(1)
  4    from store s
  5   where sdo_nn(s.geom,
  6                mdsys.sdo_geometry(2003,NULL,NULL,
  7                      mdsys.sdo_elem_info_array(1,1003,3),
  8                      mdsys.sdo_ordinate_array(380004,5100003,390000,5160000)),
  9                'sdo_num_res=10',1) = 'TRUE'
 10     and s.storetype = 'SHOE'
 11*  order by 3
gis@XE> /

        ID STORETYPE  SDO_NN_DISTANCE(1)
---------- ---------- ------------------
       271 SHOE                215195.12
        60 SHOE                215307.04
       380 SHOE                215909.35

Yes, it works, but it is subject to one guessing a value for sdo_num_res that will work in all cases. For example, sdo_num_res=10 does not work for COMPUTER stores:

gis@XE>  select id,
  2         storetype,
  3         sdo_nn_distance(1)
  4    from store s
  5   where sdo_nn(s.geom,
  6                mdsys.sdo_geometry(2003,NULL,NULL,
  7                      mdsys.sdo_elem_info_array(1,1003,3),
  8                      mdsys.sdo_ordinate_array(380004,5100003,390000,5160000)),
  9                'sdo_num_res=10',1) = 'TRUE'
 10     and s.storetype = 'COMPUTER'
 11*  order by 3
gis@XE> /

        ID STORETYPE  SDO_NN_DISTANCE(1)
---------- ---------- ------------------
       100 COMPUTER           215228.103
       465 COMPUTER            215717.65

See, we only got two! Sure, we would increase the sdo_num_res batch size again but we are playing a guessing game.

This is why the sdo_batch_size parameter exists. It exists, as the documentation says,”If any geometries in the table might be nearer than the geometries specified in the WHERE clause”. We know SHOE stores in the table are nearer than our COMPUTER stores but it is COMPUTER stores that we want.

gis@XE>  select id,
  2         storetype,
  3         sdo_nn_distance(1)
  4    from store s
  5   where sdo_nn(s.geom,
  6                mdsys.sdo_geometry(2003,NULL,NULL,
  7                      mdsys.sdo_elem_info_array(1,1003,3),
  8                      mdsys.sdo_ordinate_array(380004,5100003,390000,5160000)),
  9                'sdo_batch_size=0',1) = 'TRUE'
 10     and s.storetype = 'COMPUTER'
 11     and rownum < 4
 12*  order by 3
gis@XE> /

        ID STORETYPE  SDO_NN_DISTANCE(1)
---------- ---------- ------------------
       100 COMPUTER           215228.103
       465 COMPUTER            215717.65
         1 COMPUTER           216857.982

A correct result. I will leave you to read up on sdo_batch_size parameter value setting in the documentation.

How do we find the 3 nearest SHOE stores to every SHOE store in the STORE table? As follows?

gis@XE> select /*+ ORDERED USE_NL(s,s2)*/
  2         s.id,
  3         s2.id as nearestStoreId,
  4         sdo_nn_distance(1) as distance
  5    from store s,
  6         store s2
  7   where s.storetype = 'SHOE'
  8     and sdo_nn(s2.geom,
  9                s.geom,
 10                'sdo_batch_size=10',1) = 'TRUE'
 11     and s2.storetype = 'SHOE'
 12     and s2.id <> s.id
 13     and rownum < 4
 14   order by 1,3
 15  /

        ID NEARESTSTOREID   DISTANCE
---------- -------------- ----------
       394             39 5631.82065
       394            413  6387.1917
       394            462 6927.07224

Anyone spot the error? The rownum < 4 is applied to the final rowset and not to each s.store and its 3 nearest neighbours!

What we need to do is order the query such that the rownum test is applied to the neighbours of each and every store. One way to do this is via a correlated subquery as follows:

gis@XE> select /*+ ORDERED USE_NL(s,s2)*/
  2         s.id,
  3         s2.id as nearestStoreId,
  4         mdsys.sdo_geom.sdo_distance(s.geom,s2.geom,0.05) as distance
  5    from store s,
  6         store s2
  7   where s.storetype = 'SHOE'
  8     AND s2.id in (select id
  9                from store s3
 10               where sdo_nn(s3.geom,s.geom,'sdo_batch_size=10',1) = 'TRUE'
 11                      and s3.storetype = 'SHOE'
 12                      and s3.id <> s.id
 13                      and rownum < 4)
 14*  order by 1,2
gis@XE> /

        ID NEARESTSTOREID   DISTANCE
---------- -------------- ----------
....
       490             82 7933.95789
       490            233 7704.35301
       490            480 5899.48249
       495            159 4850.15973
       495            243 5073.66333
       495            318 7008.11984
       499             41 2635.73385
       499             74 4104.74433
       499            430 4327.98326

A correct result.

Note: Because of limitations with returned values from correlated sub-queries, we cannot access any sdo_nn_distance() values so we simply use the sdo_geom.sdo_distance function to get out required distance. If I can come up with a better query (or perhaps my readers may have a better approach) I will amend this blog posting.

Quality of Returned Distance

What I want to turn to the types of distances SDO_NN calculates. It is an approximation to the actual geometric distance done so that the speed of the SDO_NN operator remains fast (which it is) through mainly RTree processing or is it computed by looking at the actual geometries? All I will do is conduct a very simple test using a point and a line.

In the Codesys.ProjLine2D table there is a simple, straight line geometry with linetype of ‘STRAIGHTVERTEX’ as follows:

INSERT INTO ProjLine2D VALUES(
'STRAIGHTVERTEX', 
SDO_GEOMETRY(2002, NULL, NULL, SDO_ELEM_INFO_ARRAY(1, 2, 1), SDO_ORDINATE_ARRAY( 380000.0, 5100000.0, 380000.0, 5160000.0)));

You will note that it is a horizontal line of length 60,000 meters. If I search this feature using two points to show that the SDO_NN distance calculation is quite accurate.

gis@XE> select sdo_nn_distance(1)
  2    from codesys.projline2d p
  3    where p.linetype = 'STRAIGHTVERTEX'
  4    and sdo_nn(p.geom,mdsys.sdo_geometry(2001,NULL,mdsys.sdo_point_type(380010.0, 5160000.0,NULL),NULL,NULL),
  5*  'sdo_num_res=1',1) = 'TRUE'
gis@XE> /

SDO_NN_DISTANCE(1)
------------------
                10

gis@XE> select sdo_nn_distance(1)
  2    from codesys.projline2d p
  3    where p.linetype = 'STRAIGHTVERTEX'
  4    and sdo_nn(p.geom,mdsys.sdo_geometry(2001,NULL,mdsys.sdo_point_type(380010.0, 5155550.0,NULL),NULL,NULL),
  5*  'sdo_num_res=1',1) = 'TRUE'
gis@XE> /

SDO_NN_DISTANCE(1)
------------------
                10

So, SDO_NN must look at the geometries of the objects it deals with to calculate its distances and not rely on MBRs and object vertices.

Hints

You will have noted that that all the queries above that involve a single table require no hints (as expected). However, the query against the two STORE tables above used the ORDERED and USE_NL (USE_Nested_Loops) hints. Why is this? Occasionally, queries just won’t run. For example, here is a query that does not work for the data I ship with my free PL/SQL packages on the copy of XE (10gr) I am running:

gis@XE> select p.id,
  2         l.linetype,
  3         sdo_nn_distance(1)
  4    from codesys.projpoint2d sample (0.1) p,
  5         codesys.projline2d l
  6   where sdo_nn(l.geom,p.geom,'sdo_num_res=5',1) = 'TRUE'
  7     and l.linetype like '<span>VERTEX</span>'
  8  /
select p.id,
*
ERROR at line 1:
ORA-13249: SDO_NN cannot be evaluated without using index
ORA-06512: at "MDSYS.MD", line 1723
ORA-06512: at "MDSYS.MDERR", line 17
ORA-06512: at "MDSYS.PRVT_IDX", line 9

I have note that this can happen where a non-spatial attribute used in a predicate has its own index (and sometimes not). This issue is documented in the Spatial Operators chapter of the Oracle Spatial documentation and also in the excellent book Pro Oracle Spatial. The recommended solution is to include hints. Sometimes adding the /*+ORDERED*/ hint on its own can work but in other cases not. The ones recommended in the book involve knowing about the name of the RTRee and ordinary indexes over the CODESYS.PROJLINE2D table (/*+INDEX*/). Similarly for the documentation though it does suggested the LEADING hint in cases involving a join between the two tables in the FROM clause used in the SDO_NN operator. I have found that, at times, the ORDERED, INDEX and NO_INDEX hints are not enough: I have had to specify the LEADING or USE_NL (USE_Nested_Loops) hints as in the following examples:

gis@XE> select /*+ORDERED USE_NL(p,l) */
  2         p.id,
  3         l.linetype,
  4         sdo_nn_distance(1)
  5    from codesys.projpoint2d sample (0.1) p,
  6         codesys.projline2d l
  7   where sdo_nn(l.geom,p.geom,'sdo_num_res=5',1) = 'TRUE'
  8     and l.linetype like '<span>VERTEX</span>'
  9  /

        ID LINETYPE                                 SDO_NN_DISTANCE(1)
---------- ---------------------------------------- ------------------
       393 VERTEXARC                                        156738.398
       393 VERTEX                                           241507.531
       393 STRAIGHTVERTEX                                   249322.965
       393 45DEGREEVERTEX                                   300109.767

You should read the documentation and be prepared to experiment with different hints.

A Real Example

Finally, I thought I would leave you with a real example derived from something I had to do for a customer recently. What we will do is find the nearest vertex-to-vertex segment (greater than 5 meters in length) of a land parcel polygon to a road centreline with the same name. Here is a picture of some data and a road centreline. Note that the distance to the road centreline will be the same for a parcel side boundary where it joins a front boundary as both share the same (corner) vertex. To get around this we will use the mid-point of each segment.

Before sdo_nn/min-point processing

The identifiers of the hightlighted land parcels are: oid in ( 26388, 26386, 26387, 26392, 26390, 26391, 26389 )

Here is a query that achieves what we want to do.

gis@XE> select /*+ ORDERED USE_NL(a,s)*/
  2           a.oid,
  3             sdo_geom.sdo_length(sdo_geometry(2002,8307,null,
  4                                              sdo_elem_info_array(1,2,1),
  5                                              sdo_ordinate_array(b.startcoord.x,b.startcoord.y,b.endcoord.x,b.endcoord.y)),
  6                                    0.05) length
  7          from land_parcel a,
  8               table(codesys.geom.getpipedvector2d(a.geom)) b,
  9               road_centreline s
 10          where a.oid in ( 26388, 26386, 26387, 26392, 26390, 26391, 26389 )
 11            and sdo_nn(s.geom,
 12                       sdo_geometry(2001,8307, sdo_point_type((b.startcoord.x+b.endcoord.x)/2,(b.startcoord.y+b.endcoord.y)/2,null),
 13                                    NULL,NULL),
 14                       'sdo_num_res=1') = 'TRUE'
 15            and s.street_name = a.street_name
 16            and sdo_geom.sdo_length(sdo_geometry(2002,8307,null,
 17                                                 sdo_elem_info_array(1,2,1),
 18                                                 sdo_ordinate_array(b.startcoord.x,b.startcoord.y,b.endcoord.x,b.endcoord.y)),
 19                                    0.05) > 5
 20  order by 1,2
 21  /

       OID     LENGTH
---------- ----------
     26386  16.598107
     26386 35.6314592
     26386 35.6419791
     26387 15.5455833
     26387 17.5226582
     26387 32.3342913
     26387 35.6419791
     26388 11.9649317
     26388 32.3342913
     26388 37.0525111
     26389 13.0049383
     26389 36.2842051
     26389 37.0525111
     26390 16.9069198
     26390 36.2708158
     26390 36.2842051
     26391 18.5583617
     26391 19.6084577
     26391 36.2663495
     26391 36.6570959
     26392  18.900781
     26392 36.2663495
     26392 36.2708158

23 rows selected.

But we want the nearest vector for each OID. The SQL is:

gis@XE> select /*+ ORDERED USE_NL(a,s) <strong>/ 2 a.oid, 3 min(sdo_geom.sdo_length(sdo_geometry(2002,8307,null, 4 sdo_elem_info_array(1,2,1), 5 sdo_ordinate_array(b.startcoord.x,b.startcoord.y,b.endcoord.x,b.endcoord.y)), 6 0.05) ) as length 7 from land_parcel a,cad_poly 8 table(codesys.geom.getpipedvector2d(a.geom)) b, 9 road_centreline s 10 where a.oid in ( 26388, 26386, 26387, 26392, 26390, 26391, 26389 ) 11 and sdo_nn(s.geom, 12 sdo_geometry(2001,8307, sdo_point_type((b.startcoord.x+b.endcoord.x)/2,(b.startcoord.y+b.endcoord.y)/2,null), 13 NULL,NULL), 14 'sdo_num_res=1') = 'TRUE' 15 and s.street_name = a.street_name 16 and sdo_geom.sdo_length(sdo_geometry(2002,8307,null, 17 sdo_elem_info_array(1,2,1), 18 sdo_ordinate_array(b.startcoord.x,b.startcoord.y,b.endcoord.x,b.endcoord.y)), 19 0.05) > 5 20 group by a.oid 21 / <code id="source_code"> OID LENGTH ---------- ---------- 26388 11.9649317 26392 18.900781 26391 18.5583617 26387 15.5455833 26390 16.9069198 26386 16.598107 26389 13.0049383

7 rows selected.

Finally, while we could use the MIN function with an appropriate GROUP BY clause to return the vector with the minimum distance to the road centreline, we cannot easily return the midpoint of the nearest vector of each OID using these operators. To do so requires some tricky SQL using a rank/partition analytic to find the first vector whose distance to the road centreline is the minimum of all vectors that make up any one land parcel (rank = 1). Here is the final result.

gis@XE> select oid,cldist,midpoint,parcel_vector 2 from (select oid,cldist,midpoint,parcel_vector, 3 RANK() OVER (PARTITION BY oid ORDER BY cldist) as oidrank 4 from ( select /</strong>+ ORDERED USE_NL(a,s)*/ 5 a.oid, 6 sdo_geom.sdo_distance(s.geom, 7 sdo_geometry(2001,8307, 8 sdo_point_type((b.startcoord.x+b.endcoord.x)/2,(b.startcoord.y+b.endcoord.y)/2,null), 9 NULL,NULL), 10 0.05) as cldist, 11 sdo_geometry(2002,8307,null, 12 sdo_elem_info_array(1,2,1), 13 sdo_ordinate_array(b.startcoord.x,b.startcoord.y,b.endcoord.x,b.endcoord.y)) as parcel_vector, 14 sdo_geometry(2001,8307, 15 sdo_point_type((b.startcoord.x+b.endcoord.x)/2,(b.startcoord.y+b.endcoord.y)/2,null), 16 NULL,NULL) as midpoint 17 from land_parcel a, 18 table(codesys.geom.getpipedvector2d(a.geom)) b, 19 road_centreline s 20 where a.oid in ( 26388, 26386, 26387, 26392, 26390, 26391, 26389 ) 21 and sdo_nn(s.geom, 22 sdo_geometry(2001,8307,sdo_point_type((b.startcoord.x+b.endcoord.x)/2,(b.startcoord.y+b.endcoord.y)/2,null),NULL,NULL), 23 'sdo_num_res=1') = 'TRUE' 24 and s.street_name = a.street_name 25 and sdo_geom.sdo_length(sdo_geometry(2002,8307,null, 26 sdo_elem_info_array(1,2,1), 27 sdo_ordinate_array(b.startcoord.x,b.startcoord.y,b.endcoord.x,b.endcoord.y)), 28 0.05) > 5 29 ) 30 ) 31 where oidrank = 1 32 / OID CLDIST ---------- ---------- MIDPOINT(SDO_GTYPE, SDO_SRID, SDO_POINT(X, Y, Z), SDO_ELEM_INFO, SDO_ORDINATES) ------------------------------------------------------------------------------- PARCEL_VECTOR(SDO_GTYPE, SDO_SRID, SDO_POINT(X, Y, Z), SDO_ELEM_INFO, SDO_ORDINATES) ------------------------------------------------------------------------------------ 26386 5.72710973 SDO_GEOMETRY(2001, 4030, SDO_POINT_TYPE(144.673904, -37.862795, NULL), NULL, NULL) SDO_GEOMETRY(2002, 4030, NULL, SDO_ELEM_INFO_ARRAY(1, 2, 1), SDO_ORDINATE_ARRAY(144.673898, -37.862869, 144.673909, -37.86272)) 26387 10.1054417 SDO_GEOMETRY(2001, 4030, SDO_POINT_TYPE(144.67386, -37.862933, NULL), NULL, NULL) SDO_GEOMETRY(2002, 4030, NULL, SDO_ELEM_INFO_ARRAY(1, 2, 1), SDO_ORDINATE_ARRAY(144.673822, -37.862996, 144.673898, -37.862869)) 26388 16.3701483 SDO_GEOMETRY(2001, 4030, SDO_POINT_TYPE(144.673864, -37.863038, NULL), NULL, NULL) SDO_GEOMETRY(2002, 4030, NULL, SDO_ELEM_INFO_ARRAY(1, 2, 1), SDO_ORDINATE_ARRAY(144.673905, -37.863081, 144.673822, -37.862996)) 26389 11.580563 SDO_GEOMETRY(2001, 4030, SDO_POINT_TYPE(144.673974, -37.863058, NULL), NULL, NULL) SDO_GEOMETRY(2002, 4030, NULL, SDO_ELEM_INFO_ARRAY(1, 2, 1), SDO_ORDINATE_ARRAY(144.674042, -37.863036, 144.673905, -37.863081)) 26390 6.43017772 SDO_GEOMETRY(2001, 4030, SDO_POINT_TYPE(144.674138, -37.863039, NULL), NULL, NULL) SDO_GEOMETRY(2002, 4030, NULL, SDO_ELEM_INFO_ARRAY(1, 2, 1), SDO_ORDINATE_ARRAY(144.674234, -37.863043, 144.674042, -37.863036)) 26391 5.70850051 SDO_GEOMETRY(2001, 4030, SDO_POINT_TYPE(144.674554, -37.863055, NULL), NULL, NULL) SDO_GEOMETRY(2002, 4030, NULL, SDO_ELEM_INFO_ARRAY(1, 2, 1), SDO_ORDINATE_ARRAY(144.674659, -37.863059, 144.674449, -37.863051)) 26392 6.07655561 SDO_GEOMETRY(2001, 4030, SDO_POINT_TYPE(144.674341, -37.863047, NULL), NULL, NULL) SDO_GEOMETRY(2002, 4030, NULL, SDO_ELEM_INFO_ARRAY(1, 2, 1), SDO_ORDINATE_ARRAY(144.674449, -37.863051, 144.674234, -37.863043))

7 rows selected.

Mapping these midpoints (large gray circles) produces the following.

After sdo_nn/min-point processing

If you reached the end of this article, thanks for persevering with the material. I hope you got something out of it.

Comment

Converting distances and units of measure in Oracle Locator · 105 days ago by Simon Greener

Ever wanted to know what a decimal degree was in meters, nautical miles or feet?

It is something I often need to do in Oracle so I decided to do something about creating a function that would do this. This involved me having to hack my way into some of the mdsys coordinate system and distance units tables but I managed to come up with something that I have integrated in to the GEOM package of the free PL/SQL code available for download from my website.

The first table that we need to look at is the mdsys.SDO_DIST_UNITS table. This table is described as follows:

codesys@XE> desc mdsys.sdo_dist_units Name Null? Type ----------------------------------------------------------------------- -------- ------------ SDO_UNIT VARCHAR2(80) UNIT_NAME NOT NULL VARCHAR2(80) CONVERSION_FACTOR NUMBER

Let’s have a look at some of the entries of this table (I ignore NULL sdo_unit names in this article and in my code):

codesys@XE> select substr(sdo_unit,1,25) as unit, 2 substr(unit_name,1,30) as name, 3 conversion_factor 4 from mdsys.sdo_dist_units 5* where sdo_unit is not null codesys@XE> / UNIT NAME CONVERSION_FACTOR ------------------------- ------------------------------ ----------------- M Meter 1 METER Meter 1 KM Kilometer 1000 KILOMETER Kilometer 1000 CM Centimeter .01 CENTIMETER Centimeter .01 MM Millemeter .001 MILLIMETER Millemeter .001 MILE Mile 1609.344 NAUT_MILE Nautical Mile 1852 SURVEY_FOOT U.S. Foot .30480061 FOOT Foot (International) .3048 INCH Inch .0254 YARD Yard .9144 CHAIN Chain 20.1168 ROD Rod 5.0292 LINK Link .201166195 MOD_USFT Modified American Foot .304812253 CL_FT Clarke's Foot .304797265 IND_FT Indian Foot .304799518 LINK_BEN Link (Benoit) .201167651 LINK_SRS Link (Sears) .201167651 CHN_BEN Chain (Benoit) 20.1167825 CHN_SRS Chain (Sears) 20.1167651 IND_YARD Yard (Indian) .914398554 SRS_YARD Yard (Sears) .914398415 FATHOM Fathom 1.8288 British foot (1936) British foot (1936) .304800749

From this we can see that CONVERSION_FACTOR expresses how long a unit (1) of measure (eg FOOT) is terms of meters (eg .3048)

OK, so now we know how to convert between any two units in this table. Here is a function that will do it.

create or replace Function Convert_Unit( p_from_unit in varchar2, p_value in number, p_to_unit in varchar2 ) return number Is v_from_conversion_factor number; v_to_conversion_factor number; Begin If ( p_value is null or p_from_unit is null or p_to_unit is null ) Then raise_application_error( codesys.Constants.c_i_null_parameter, codesys.Constants.c_s_null_parameter,False ); End If; -- Check if p_from_unit exists by getting the necessary conversion factor to meters BEGIN -- Note that the conversion_factor is a conversion factor between v_from_unit and 1 metre. SELECT conversion_factor INTO v_from_conversion_factor FROM mdsys.sdo_dist_units WHERE sdo_unit = UPPER(p_from_unit) AND ROWNUM = 1; EXCEPTION WHEN NO_DATA_FOUND THEN raise_application_error( codesys.Constants.c_i_invalid_unit, codesys.Constants.c_s_invalid_unit || ' ' || p_from_unit); END; -- Check if p_to_unit exists by getting the necessary conversion factor to meters BEGIN -- Note that the conversion_factor is a conversion factor between v_to_unit and 1 metre. SELECT conversion_factor INTO v_to_conversion_factor FROM mdsys.sdo_dist_units WHERE sdo_unit = UPPER(p_to_unit) AND ROWNUM = 1; EXCEPTION WHEN NO_DATA_FOUND THEN raise_application_error( codesys.Constants.c_i_invalid_unit, codesys.Constants.c_s_invalid_unit || ' ' || p_to_unit); END; -- Do the computation RETURN ( p_value * v_from_conversion_factor ) / v_to_conversion_factor; End Convert_Unit;

And some examples on how to use this function (also a part of the GEOM package):

codesys@XE> select convert_unit('CHAIN',1,'LINK') from dual; CONVERT_UNIT('CHAIN',1,'LINK') ------------------------------ 100.000897

But what if we have data coded to a SRID and want to convert from its unit of measure to one of those in the mdsys.sdo_dist_units table? For example, you will note that there is no sdo_unit for ‘Decimal Degrees’ in which longitude/latitude data is expressed. We need to look at the definition of a SRID to find the conversion information we need. This is held in the MDSYS.CS_SRS table.

codesys@XE> desc mdsys.cs_srs Name Null? Type ----------------------------------------------------------------------- -------- ------------------ CS_NAME VARCHAR2(80) SRID NOT NULL NUMBER(38) AUTH_SRID NUMBER(38) AUTH_NAME VARCHAR2(256) WKTEXT VARCHAR2(2046) CS_BOUNDS MDSYS.SDO_GEOMETRY

For example, let’s look at the projection information for WGS84 (SRID = 8307). I have formatted the output for readability.

codesys@XE> select wktext 2 from mdsys.cs_srs 3 where srid = 8307; WKTEXT ----------------------------------------------------------------------------------------------------------------------------------- GEOGCS [ "Longitude / Latitude (WGS 84)", DATUM ["WGS 84", SPHEROID ["WGS 84", 6378137, 298.257223563] ], PRIMEM [ "Greenwich", 0.000000 ], UNIT ["Decimal Degree", 0.01745329251994330] ]

OK, we can see that the unit of measure for this geographic coordinate system is our “Decimal Degrees”! But how do we access it?

In my downloadable PL/SQL code there is a string tokenizer which I can use as follows:

For a projected coordinate system eg SRID 2964 “NAD27 / Alaska Albers” the tokens would be:

codesys@XE> SELECT rownum as id, 2 substr(trim(both ' ' from replace(b.column_value,'"')),1,40) as token 3 FROM mdsys.cs_srs a, 4 TABLE(codesys.Tokenizer(a.wktext,',[]')) b 5* WHERE srid = 2964 codesys@XE> / ID TOKEN ---------- ---------------------------------------- 1 PROJCS 2 NAD27 / Alaska Albers 3 GEOGCS 4 NAD27 5 DATUM 6 North American Datum 1927 (EPSG ID 6267) 7 SPHEROID 8 Clarke 1866 (EPSG ID 7008) 9 6378206.4 10 294.978698213905820761610537123195175418 11 PRIMEM 12 Greenwich 13 0.000000 14 UNIT 15 Decimal Degree 16 0.01745329251994328 17 PROJECTION 18 Alaska Albers (EPSG OP 15020) 19 UNIT 20 U.S. Foot 21 .304800609601219202438404876809753619507

So, all we have to do is iterate over this list to extract the conversion unit for a coordinate system. In this case the second parameter of the SPHEROID entry (line 7 of the tokens for srid 8307). Since both (sdo_dist_unit and wktext) conversion units are expressed relative to meters we now have the ability to convert a distance expressed in the units of measure of a coodinate system to any unit of measure in the mdsys.sdo_dist_units table via a simple equation:

new_value = ( value x srid_conversion_factor ) / unit_conversion_factor

One other thing, for geographic coordinate systems (first parameter = GEOCS and not PROJCS) the value associated with the “Decimal Degrees” UNIT must be multiplied by the radius of the earth. Thus the equation would be:

new_value = ( value x srid_conversion_factor x radius_of_earth ) / unit_conversion_factor

So, we can now construct a function that will do this conversion (please excuse the length of this):

create or replace Function Convert_Distance( p_srid in number, p_value in number, p_unit in varchar2 := 'Meter' ) Return number Is v_unit varchar2(1000) := UPPER(p_unit); v_unit_conversion_factor number; v_srid_conversion_factor number; v_radius_of_earth number := 6378137; -- Default v_length number; v_srid mdsys.cs_srs.SRID%TYPE; v_token_id number; v_token varchar2(4000); v_geocs boolean; cursor c_cs_tokens(p_srid in number) Is select rownum as id, substr(trim(both ' ' from replace(b.column_value,'"')),1,40) as token from mdsys.cs_srs a, table(codesys.Tokenizer(a.wktext,',[]')) b where srid = p_srid; Begin If ( p_srid is null ) Then -- Normally Oracle assumes a NULL srid is planar but -- this could be planar feet, or meters etc so throw an error raise_application_error( codesys.Constants.c_i_null_srid, codesys.Constants.c_s_null_srid,False ); End If; If ( p_value is null ) Then raise_application_error( codesys.Constants.c_i_null_parameter, codesys.Constants.c_s_null_parameter,False ); End If; -- Check if p_unit exists by getting the necessary conversion factor to meters BEGIN -- Note that the conversion_factor is a conversion factor between v_unit and 1 metre. SELECT conversion_factor INTO v_unit_conversion_factor FROM mdsys.sdo_dist_units WHERE sdo_unit = v_unit AND ROWNUM = 1; EXCEPTION WHEN NO_DATA_FOUND THEN raise_application_error( codesys.Constants.c_i_invalid_unit, codesys.Constants.c_s_invalid_unit || v_unit); END; -- Check if SRID exists BEGIN SELECT srid INTO v_srid FROM mdsys.cs_srs WHERE srid = p_srid; EXCEPTION WHEN NO_DATA_FOUND THEN raise_application_error( codesys.Constants.c_i_invalid_srid, codesys.Constants.c_s_invalid_srid || p_srid); END; -- We need to get the conversion factor to meters and earth's radius for the supplied SRID. -- This can only be gotten by getting the WKTEXT in mdsys.cs_srs, breaking it into tokens, -- and finding the right ones: -- SPHEROID + 2 tokens = Radius -- Last UNIT + 1 = conversion unit -- Last UNIT + 2 = conversion unit value FOR rec IN c_cs_tokens(p_srid) LOOP If ( rec.id = 1 ) Then v_geocs := case rec.token when 'GEOGCS' then true else false end; ElsIf ( rec.token = 'SPHEROID' ) Then v_token := rec.token; v_token_id := rec.id + 2; ElsIf ( rec.token = 'UNIT' ) Then v_token := rec.token; v_token_id := rec.id + 2; End If; If ( rec.id = v_token_id ) Then If ( v_token = 'SPHEROID' ) Then v_radius_of_earth := to_number(rec.token); ElsIf ( v_token = 'UNIT' ) Then v_srid_conversion_factor := to_number(rec.token); End If; End If; END LOOP; If ( v_geocs ) Then v_srid_conversion_factor := v_srid_conversion_factor * v_radius_of_earth; End If; -- OK, now we have a conversion factor from p_unit to meters -- and a conversion factor for the units to meters -- The returned value is: p_value * v_srid_conversion_factor (to get value in meters) / v_unit_conversion_factor (to convert from meters to the unit) -- return ( p_value * v_srid_conversion_factor ) / v_unit_conversion_factor; End Convert_Distance;

Some examples:

codesys@XE> select Convert_Distance(8311,1,'Meter') as meters_per_degree, 2 1 / Convert_Distance(8311,1,'Meter') as degrees_per_metre 3 from dual; METERS_PER_DEGREE DEGREES_PER_METRE ----------------- ----------------- 111319.491 8.9832E-06 codesys@XE> select Convert_Distance(8311,1,'Foot') as feet_per_degree, 2 1 / Convert_Distance(8311,1,'Foot') as degrees_per_foot 3 from dual; FEET_PER_DEGREE DEGREES_PER_FOOT --------------- ---------------- 365221.426 2.7381E-06 codesys@XE> select Convert_Distance(2964,1,'Meter') as feet_per_metre, 2 1 / Convert_Distance(2964,1,'Meter') as metres_per_foot 3 from dual; FEET_PER_METRE METRES_PER_FOOT -------------- --------------- .30480061 3.28083333

These functions have been integrated into my GEOM PL/SQL package. These functions are useful in the context of my packages because I have recently added the ability to convert special elements of an sdo_geometry (eg rectangles, circles and circular arcs) to vertex-to-vertex connected segments. These require an arc2chord value which is expressed in dataset units. So, for geographic data this is decimal degrees not meters (unlike Oracle itself). Similarly, the tolerance parameter of the sdo_centroid function is similarly expressed in dataset units.

If anyone finds any errors in my work, please drop me a line.

Comment

>>Tips and Tricks

SpatialDB Advisor

Converting distances and units of measure in Oracle Locator

Generating multi-points from single point records in Oracle Spatial

Multi-Centroid Shootout

Oracle Locator vs Oracle Spatial: A Reflection on Oracle Licensing of the SDO_GEOM Package

Oracle Spatial Centroid Shootout

Oracle Spatial Forum - Melbourne April 2007

Oracle Spatial Mapping and Map Rendering Performance Tips

Oracle's SQL/MM Compliant Types

Spatial Pipelining

Split Sdo_Geometry Linestring at a known point

The significance of sdo_lb/sdo_ub in USER_SDO_GEOM_METDATA: Do I need it?

Tip #2: layer_gtypes for spatial indexes

Tip #3: Generating sample (test) point data

Tip #4: Loading Shapefiles (SHP) into Oracle Spatial

Tip #5: FAST REFRESHing of Oracle Materialized Views containing Sdo_Geometry columns

Tip #6: Unpacking USER_SDO_GEOM_METADATA's DIMINFO structure

Tip #7: Forcing Sdo_Geometry to contain only points, lines or areas

Tips and Tricks

Trick #1: Australian MGA/AMG Zone Calculation from geographic (longitude/latitude) data

Trick #2: Object Tables of Sdo_Geometry

Using Oracle's SDO_NN Operator - Some examples

Multi-Centroid Shootout · 83 days ago by Simon Greener

Spatial Pipelining · 94 days ago by Simon Greener

Using Oracle's SDO_NN Operator - Some examples · 103 days ago by Simon Greener

Converting distances and units of measure in Oracle Locator · 105 days ago by Simon Greener