Thursday, 5 October 2017

arcgis desktop - ArcMap error drawing postgis layer - The number of points is less than required for feature


I have one big polygon called huc and a multilinestring called streams, I am running an intersection on between the two layers to clip the streams to the huc.


here is what I have so far


drop table if exists base_layers.streams_huc_clip;
create table base_layers.streams_huc_clip as
select s.*,case
when st_within(s.shape,huc.shape) then s.shape
when st_intersects(s.shape,huc.shape) then st_intersection(s.shape,huc.shape)
end as geom2
from postgres.streams s join base_layers.huc_hl_dissolve huc on st_intersects(s.shape,huc.shape);


--pretty straightforward query in postgis, this query returns with no errors and correctly gets displayed in QGIS.


However when I try and bring it into ArcGIS I get this error and it does not draw


One or more layers failed to draw:

lucz_2017.base_layers.streams_huc_clip: The number of points is less than required for feature

so here are the following steps I tried after to do after it would not show up in ArcGIS to no avail.


dropped the original geometry so I just had one geometry column, renamed the new geom to shape


alter table base_layers.streams_huc_clip drop column shape;

alter table base_layers.streams_huc_clip rename geom2 to shape;

dropped the objectid and added new primary key, changed owners,


alter table base_layers.streams_huc_clip drop column objectid;
ALTER TABLE base_layers.streams_huc_clip ADD COLUMN id SERIAL PRIMARY KEY;
alter table base_layers.streams_huc_clip owner to base_layers;

created an index and tried updating the srid. none of this worked!


ALTER TABLE base_layers.streams_huc_clip ALTER COLUMN shape type geometry(multilinestring, 3424) using st_multi(shape);
CREATE INDEX streams ON base_layers.streams_huc_clip USING GIST (shape);


also ran this query


SELECT * FROM public.geometry_columns ORDER BY type

which returned


enter image description here


I am thinking about running an st_dump on the geometry to turn it into a linestring, I will take any suggestions on why this layer is not displaying and what I can do to make it work


I ran an st_dump and turned the geometry into a linestring but I still got the same arcmap error


drop table if exists base_layers.streams_huc_clip;
create table base_layers.streams_huc_clip as

select (st_dump(geom)).geom as geom2,t.* from(
select s.*,case
when st_within(s.shape,huc.shape) then s.shape
when st_intersects(s.shape,huc.shape) then st_intersection(s.shape,huc.shape)
end as geom
from postgres.streams s join base_layers.huc_hl_dissolve huc on st_intersects(s.shape,huc.shape))t;

ST_MakeValid() did not help either


I make edits to these tables outside the ESRI stack.



Answer




The issue is caused by the tolerance ArcMap has for determining whether two vertices are identical vs the precision of the coordinates passed to it. The general concepts are well explained in great detail here, and a more bite size form here, which also states Arc's defaults as:



  • 0.0001 map units for Tolerance (this is the cell size of the grid on which vertices are stored) - it determines the precision (no. of decimal places) at which coordinatess are stored.

  • 0.001 map units for Resolution which is the threshold below which two points are considered to be the same..


When ArcMap tries to read the geometry it maps vertices to the underlying grid, and removes repeated points, and so geometries can become invalid causing the error:


The number of points is less than required for feature


This can be solved by using PostGIS to snap vertices to an equivalent grid, and remove duplicated points before passing it to ArcMap. One way to do this is to use ST_Snap_To_Grid which snaps all vertices to an imaginary grid at a given tolerance.


This can be applied only to the table you wish to bring into ArcMap (preserving the precision of the original data from which the clip was derived) with:


UPDATE base_layers.streams_huc_clip SET shape = ST_RemoveRepeatedPoints(ST_SnapToGrid(shape, 0.001));

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...