How to split a linestring every 1000 points in postgis

Question

To assume compatibility with my survey software, I would like to lighten geometries of my linestrings by cutting them every 1000 points without creating or interpolating points.

I transformed a function found on the net (resample) but it cuts only the first 1000 points...

Does someone knows how to loop it correctly?

CREATE OR REPLACE FUNCTION simplify_npoints(inGeom geometry, maxPoints integer) RETURNS geometry AS $$

DECLARE nPoints integer; outGeom geometry; fraction float; points geometry[];

BEGIN

nPoints := ST_NPoints(inGeom); outGeom := inGeom;

IF nPoints > maxPoints THEN FOR i IN 1..maxPoints LOOP points := array_append(points, ST_PointN(inGeom,i)); END LOOP; outGeom := ST_MakeLine(points);

END IF;

RETURN outGeom;

END; $$ LANGUAGE plpgsql IMMUTABLE;

I'll risk a silly question: how many points do your linestrings actually have? SELECT ST_NPoints(linestr.the_geom) — Micha
– Micha, Commented May 3, 2013 at 16:01
I've got a set 144 simple linestrings with 8000 to 75000 points. — aurel_nc
– aurel_nc, Commented May 5, 2013 at 21:21

Jakub Kania · Accepted Answer · 2013-05-07 20:54:53Z

You haven't stated what output you want so an array of geometry should do:

CREATE OR REPLACE FUNCTION simplify_npoints(inGeom geometry, maxPoints integer) RETURNS geometry[] AS $$ DECLARE outGeom geometry[]=ARRAY[]::geometry[]; points geometry[]=ARRAY[]::geometry[]; counter integer:=0; BEGIN IF maxPoints=1 THEN RAISE EXCEPTION 'maxPoints must be >1. Split not possible'; END IF; FOR i IN 1..ST_NPoints(inGeom) LOOP counter:=counter+1; points := array_append(points, ST_PointN(inGeom,i)); IF counter=maxPoints THEN outgeom:=array_append(outGeom,ST_MakeLine(points)); points := ARRAY[ST_PointN(inGeom,i)]; counter:=1; END IF; END LOOP; IF counter>1 THEN outgeom:=array_append(outGeom,ST_MakeLine(points)); END IF; RETURN outGeom; END; $$ LANGUAGE plpgsql IMMUTABLE;

Thank you so much! Initially, geometry array scared me a little. But after some research, they are treated as simple arrays[]. In the end I created a test tables by the following formula: DROP TABLE cut_line; CREATE TABLE cut_line AS SELECT gid , type , unnest(s.geom2) as geom FROM (SELECT gid , type , simplify_npoints(geom,1000)as geom2 FROM long_line WHERE gid = '35852') as s; And it perfectly works! — aurel_nc
– aurel_nc, Commented May 6, 2013 at 4:59

aurel_nc · Accepted Answer · 2013-05-06 23:06:34Z

I realized that the ending trunk of the linestrings was missing. We were only making new line when counter matched the maxpoints value, leaving the rest.
I add a new counter "totalpoints" and a new case system to reach the end of the linestring:

CREATE OR REPLACE FUNCTION simplify_npoints(ingeom geometry, maxpoints integer) RETURNS geometry[] AS $BODY$ DECLARE outGeom geometry[]=ARRAY[]::geometry[]; points geometry[]=ARRAY[]::geometry[]; counter integer:=0; countmax integer:=0; totalpoints integer:=0; BEGIN IF maxPoints=1 THEN RAISE EXCEPTION 'maxPoints must be >1. Split not possible'; END IF; totalpoints:=ST_NPoints(inGeom); FOR i IN 1..totalpoints LOOP countmax:=countmax+1; counter:=counter+1; points := array_append(points, ST_PointN(inGeom,i)); CASE WHEN counter=maxPoints THEN outgeom:=array_append(outGeom,ST_MakeLine(points)); points := ARRAY[ST_PointN(inGeom,i)]; counter:=1; WHEN countmax=totalpoints THEN outgeom:=array_append(outGeom,ST_MakeLine(points)); countmax:=0; ELSE END CASE; END LOOP; RETURN outGeom; END; $BODY$ LANGUAGE plpgsql IMMUTABLE

You are correct, I missed that. But the check shouldn't be in the loop, checking after the loop ends is enough. I updated my answer. — Jakub Kania
– Jakub Kania, Commented May 7, 2013 at 20:59

Stack Exchange Network

How to split a linestring every 1000 points in postgis

2 Answers 2

Hot Network Questions

How to split a linestring every 1000 points in postgis

2 Answers 2

Related

Hot Network Questions