Friday, 24 January 2020

Why can PostGIS create LINESTRING from the same values of coordinates of the POINT type?


So, I needed to create linear objects from point objects.





  1. I split the objects into points with the following query:

    CREATE TABLE roads_arc_var3 AS SELECT (ST_DumpPoints(geom)).geom FROM roads_arc_var1;




  2. I created the following table, so that it created two identical fields with the same values of coordinates:


    CREATE TABLE roads_arc_var4 AS SELECT a.geom AS a, b.geom AS b FROM roads_arc_var3 AS a, roads_arc_var3 AS b WHERE a.geom = b.geom




  3. I added the field "geoline" to the table:


    SELECT AddGeometryColumn('roads_arc_var4', 'geoline', 4326, 'LINESTRING', 2);





  4. I changed the values of the initial and final coordinates in some records in order to create lines from them, and I expected that in the 'geoline' field line values will appear only in those records in which the values of the coordinates in the starting point are not equal to the values of the coordinates end point, the rest will remain empty


    UPDATE roads_arc_var4 SET geoline = ST_MakeLine(a, b);



  5. But, the picture was different, see figure


enter image description here


What is the explanation for the PostGIS developers of this situation, because it's about 2D?


Or it is a problem with rounding off the values of coordinates, but in this case the user needs to understand what is its length, etc.




Answer



Many of the PostGIS functions may return topologically invalid geometries. It is often more practical to do so than to throw an error or stop the process.


SELECT ST_AsText(
ST_MakeLine(
ST_GeomFromText('POINT (1 1)'),
ST_GeomFromText('POINT (1 1)')));

Result: "LINESTRING(1 1,1 1)"


SELECT ST_IsValid(
ST_MakeLine(

ST_GeomFromText('POINT (1 1)'),
ST_GeomFromText('POINT (1 1)')));

Result: "f"


SELECT ST_IsValidReason(
ST_MakeLine(
ST_GeomFromText('POINT (1 1)'),
ST_GeomFromText('POINT (1 1)')));

Result: "Too few points in geometry component[1 1]"



Fortunately it is not hard to add some logic into your SQL. You can for example add WHERE NOT ST_Equals(a.geom,b.geom) or then check the result with ST_IsValid and save only the valid values.


No comments:

Post a Comment

arcpy - Changing output name when exporting data driven pages to JPG?

Is there a way to save the output JPG, changing the output file name to the page name, instead of page number? I mean changing the script fo...