This content has been marked as final. Show 8 replies
Well, I think this might have to do with interpreting the invalid geometry. For some contrast I headed over to a certain Open Source competitor and tried
And got back a multipolygon equivalent of your unwanted results:
SELECT ST_AsText( ST_MakeValid( ST_GeomFromText('POLYGON ((0.0 1.0, 4.0 1.0, 4.0 4.0, 1.0 4.0, 1.0 0.0, 2.0 0.0, 2.0 3.0, 3.0 3.0, 3.0 2.0, 0.0 2.0, 0.0 1.0))') ) ) AS shape
So my thought is both Oracle Spatial and PostGIS interpret the correction of the polygon similarly and thus perhaps its about the math and not so much about what we visually expect.
MULTIPOLYGON(((1 2,1 4,4 4,4 1,2 1,2 2,1 2),(2 2,3 2,3 3,2 3,2 2)),((1 1,0 1,0 2,1 2,1 1)),((2 1,2 0,1 0,1 1,2 1)))
Can anyone chime in with a explanation of how the code flattens out all the intersections into rings? The "good" hole is clearly clockwise and thus remains a hole, but the lower left unwanted hole has edges that look counter-clockwise on paper to me but the logic must see things differently. This is an interesting question. If you submit this to support, please tell us what the response is. I am not sure there are any rules to interpreting a bad polygon that you can cite are being broken here. But clearly PostGIS and Oracle have a common logic they are following.
I had a play around as well, even trying self-unions (and other things), and I can't get the overlapping area to remain.
I thought it was just me, but I am glad to see you had similar issues. Yes, it might be some interpretation that SFS demands on implementors. I don't know.
I could get what Markus wants but I had to:
1. Smash the polygon apart into its constituent linestrings (as 2 vertex vectors).
2. Intersect the linestring to correctly node all intersections;
3. Then pump the result into a polygon builder to get all the polygons back which could then be unioned.
The short and simple answer to your question is: No there isn't. This polygon is intersecting itself in to many points, thus making it impossible for oracle to decide which bit is which. In general you will see that when two intersections occur that the area between the two intersections is seen as a hole, instead of an area.
The only way to obtain the geometry you want is to split the valid hole into the geometries you want, then add them together again. AFAIK there is no existing software that can do that out of the box (and I've worked with/for the biggest vendors in the GIS-software industry), but this will be the easiest done using something that gives you visual tools. It might be possible to write some custom code that would be able to deal with this kind of thing, but that will be difficult too: the amount of possible permutations of this kind of self-intersection/overlapping is simply so large that it is nigh impossble to account for all possibilities. The sample geometry you drew is nice and straight edged, but that is something that rarely happens in real life :-)
Out of curiosity: what kind of software/sensor generates this kind of geometries?
HTH (but I don't think it does ;-) ),
I agree with Stefan when he said about custom code:
It might be possible to write some custom code that would be able to deal with this kind of thing, but that will be difficult too: the amount of possible permutations of this kind of self-intersection/overlapping is simply so large that it is nigh impossble to account for all possibilities.Here is some custom code for this use case but it uses a lot of non-Oracle code that I have put together over the years. And, in the end, it generates all the constituent polygons (including the one that went missing) but you still have to put it together because the inner ring becomes a polygon too.
with geom as ( select SDO_GEOMETRY(2003, NULL, NULL, SDO_ELEM_INFO_ARRAY(1, 1003, 1), SDO_ORDINATE_ARRAY(0,1, 4,1, 4,4, 1,4, 1,0, 2,0, 2,3, 3,3, 3,2, 0,2, 0,1)) as geom from dual ) select t.geom, SDO_GEOM.relate(t.geom,'DETERMINE',c.geom,0.005) as rel from TABLE(ST_EXTRACTRINGS( (select codesys.SC4O.ST_PolygonBuilder(g.gArray,3) AS POLYGONS from (select CAST(COLLECT(T_Vector(va.vector_id, va.element_id, va.subelement_id, va.startCoord, va.endCoord).AsSdoGeometry(null)) as mdsys.sdo_geometry_array) as gArray from (select sdo_geom.sdo_union(SDO_UTIL.POLYGONTOLINE(a.geom),SDO_UTIL.POLYGONTOLINE(a.geom),0.05) as uGeom from geom a ) f, table(ST_Vectorize(f.uGeom)) va ) g ))) t, geom c; -- Result -- GEOM REL ---------------------------------------------------------------------------------------------------------- ------------------- SDO_GEOMETRY(2003,NULL,NULL,SDO_ELEM_INFO_ARRAY(1,1003,1),SDO_ORDINATE_ARRAY(1,1,1,0,2,0,2,1,1,1)) COVEREDBY SDO_GEOMETRY(2003,NULL,NULL,SDO_ELEM_INFO_ARRAY(1,1003,1),SDO_ORDINATE_ARRAY(2,1,2,2,1,2,1,1,2,1)) TOUCH SDO_GEOMETRY(2003,NULL,NULL,SDO_ELEM_INFO_ARRAY(1,1003,1),SDO_ORDINATE_ARRAY(2,1,4,1,4,4,1,4,1,2,2,2,2,1)) OVERLAPBDYINTERSECT SDO_GEOMETRY(2003,NULL,NULL,SDO_ELEM_INFO_ARRAY(1,1003,1),SDO_ORDINATE_ARRAY(2,2,3,2,3,3,2,3,2,2)) TOUCH SDO_GEOMETRY(2003,NULL,NULL,SDO_ELEM_INFO_ARRAY(1,1003,1),SDO_ORDINATE_ARRAY(1,2,0,2,0,1,1,1,1,2)) COVEREDBY SDO_GEOMETRY(2003,NULL,NULL,SDO_ELEM_INFO_ARRAY(1,1003,1),SDO_ORDINATE_ARRAY(2,2,3,2,3,3,2,3,2,2)) TOUCH 6 rows selected
Hi all, thanks for the quick response.
@Simon: Thanks for the code, but what are those special non-Oracle functions? Where can I find them?
@Paul: I also think this has something to do with the ordering (clockwise/anti-clockwise). The calculation of a polygon area is really simple but only works correct if its in the correct order and if it is not self intersecting.
@Stefan: Its a really simple scenario: Our software allows users to select areas in the map. This areas are used for monitoring of vehicles (entering, leaving this area). And now think about a city ring.
The non-Oracle functions are available from my website.
SC4O is Spatial Companion For Oracle.
ST_EXTRACTRINGS is in GEOM package (not named as such) but here is the source code:
create or replace FUNCTION ST_ExtractRings(p_geom in mdsys.sdo_geometry) RETURN T_Geometries PIPELINED IS c_i_empty_geometry Constant pls_integer := -20120; c_s_empty_geometry Constant VarChar2(100) := 'Input geometry must not be null or empty'; v_num_rings pls_integer; v_ring mdsys.sdo_geometry; BEGIN If ( p_geom is null or p_geom.get_gtype() not in (3,7) ) Then raise_application_error(c_i_empty_geometry,c_s_empty_geometry,true); End If; <<process_all_elements>> FOR v_elem_no IN 1..MDSYS.SDO_UTIL.GETNUMELEM(p_geom) LOOP v_num_rings := ST_GetNumRings( MDSYS.SDO_UTIL.EXTRACT(p_geom,v_elem_no,0),0); If ( v_num_rings = 1 ) Then PIPE ROW (T_Geometry(MDSYS.SDO_UTIL.EXTRACT(p_geom,v_elem_no,1))); Else <<process_all_rings>> FOR v_ring_no IN 1..v_num_rings LOOP v_ring := MDSYS.SDO_UTIL.EXTRACT(p_geom,v_elem_no,v_ring_no); IF ( v_ring is not null ) Then PIPE ROW(T_Geometry(v_ring)); End If; END LOOP process_all_rings; End If; END LOOP process_all_elements; RETURN; END ST_ExtractRings; / -- With ST_GetNumRings being. create or replace Function ST_GetNumRings(p_geom in mdsys.sdo_geometry, p_ring_type in integer default 0 /* 0 = ALL; 1 = OUTER; 2 = INNER */ ) Return Number Deterministic Is c_i_empty_geometry Constant pls_integer := -20120; c_s_empty_geometry Constant VarChar2(100) := 'Input geometry must not be null or empty'; c_i_ring_type Constant pls_integer := -20120; c_s_ring_type Constant VarChar2(100) := 'p_ring_type must be one of 0(ALL),1(OUTER),2(INNER) only.'; v_elements pls_integer := 0; v_ring_count pls_integer := 0; v_etype pls_integer; Begin If ( ST_isEmpty(p_geom)<>'f' ) Then raise_application_error(c_i_empty_geometry,c_s_empty_geometry,true); End If; If ( p_ring_type is null OR p_ring_type not in (0,1,2) ) Then raise_application_error(c_i_ring_type,c_s_ring_type); End If; v_elements := ( ( p_geom.sdo_elem_info.COUNT / 3 ) - 1 ); <<element_extraction>> FOR v_i IN 0 .. v_elements LOOP v_etype := p_geom.sdo_elem_info(v_i * 3 + 2); If ( ( v_etype in (1003,1005,2003,2005) and 0 = p_ring_type ) OR ( v_etype in (1003,1005) and 1 = p_ring_type ) OR ( v_etype in (2003,2005) and 2 = p_ring_type ) ) Then v_ring_count := v_ring_count + 1; End If; END LOOP element_extraction; Return v_ring_count; End ST_GetNumRings;
Hi Simon, thank you very much for your help. I will have a look at it.
The SC4O.ST_PolygonBuilder is written in Java and is a front-end over JTS.
If you download SC4O the source code is included.
The vectorisation can be done in Java.
The relate can also be done using JTS - see SC4O.ST_Relate - as could the decomposition of the polygon into its rings.
If you need some guidance email me directly simon at spatialdbadvisor dot com