pggeodb.nancy.inra.fr/db_cefs - db_cefs on pggeodb.nancy.inra.fr

Function: st_concavehull

Schema

public

Owner

postgres

Descriptions

args: geomA, target_percent, allow_holes=false - The concave hull of a geometry represents a possibly concave geometry that encloses all geometries within the set. You can think of it as shrink wrapping.

Options

Definition

CREATE OR REPLACE FUNCTION public.st_concavehull (
 param_geom public.geometry,
 param_pctconvex double precision,
 param_allow_holes boolean = false
)
RETURNS public.geometry AS
$span$
DECLARE
var_convhull geometry := ST_ConvexHull(param_geom);
var_param_geom geometry := param_geom;
var_initarea float := ST_Area(var_convhull);
var_newarea float := var_initarea;
var_div integer := 6;
var_tempgeom geometry;
var_tempgeom2 geometry;
var_cent geometry;
var_geoms geometry[4];
var_enline geometry;
var_resultgeom geometry;
var_atempgeoms geometry[];
var_buf float := 1;
BEGIN
-- We start with convex hull as our base
var_resultgeom := var_convhull;

IF param_pctconvex = 1 THEN
return var_resultgeom;
ELSIF ST_GeometryType(var_param_geom) = 'ST_Polygon' THEN -- it is as concave as it is going to get
IF param_allow_holes THEN -- leave the holes
RETURN var_param_geom;
ELSE -- remove the holes
var_resultgeom := ST_MakePolygon(ST_ExteriorRing(var_param_geom));
RETURN var_resultgeom;
END IF;
END IF;
IF ST_Dimension(var_resultgeom) > 1 AND param_pctconvex BETWEEN 0 and 0.98 THEN
-- get linestring that forms envelope of geometry
var_enline := ST_Boundary(ST_Envelope(var_param_geom));
var_buf := ST_Length(var_enline)/1000.0;
IF ST_GeometryType(var_param_geom) = 'ST_MultiPoint' AND ST_NumGeometries(var_param_geom) BETWEEN 4 and 200 THEN
-- we make polygons out of points since they are easier to cave in.
-- Note we limit to between 4 and 200 points because this process is slow and gets quadratically slow
var_buf := sqrt(ST_Area(var_convhull)*0.8/(ST_NumGeometries(var_param_geom)*ST_NumGeometries(var_param_geom)));
var_atempgeoms := ARRAY(SELECT geom FROM ST_DumpPoints(var_param_geom));
-- 5 and 10 and just fudge factors
var_tempgeom := ST_Union(ARRAY(SELECT geom
FROM (
-- fuse near neighbors together
SELECT DISTINCT ON (i) i,  ST_Distance(var_atempgeoms[i],var_atempgeoms[j]), ST_Buffer(ST_MakeLine(var_atempgeoms[i], var_atempgeoms[j]) , var_buf*5, 'quad_segs=3') As geom
FROM generate_series(1,array_upper(var_atempgeoms, 1)) As i
INNER JOIN generate_series(1,array_upper(var_atempgeoms, 1)) As j
ON (
NOT ST_Intersects(var_atempgeoms[i],var_atempgeoms[j])
AND ST_DWithin(var_atempgeoms[i],var_atempgeoms[j], var_buf*10)
)
UNION ALL
-- catch the ones with no near neighbors
SELECT i, 0, ST_Buffer(var_atempgeoms[i] , var_buf*10, 'quad_segs=3') As geom
FROM generate_series(1,array_upper(var_atempgeoms, 1)) As i
LEFT JOIN generate_series(ceiling(array_upper(var_atempgeoms,1)/2)::integer,array_upper(var_atempgeoms, 1)) As j
ON (
NOT ST_Intersects(var_atempgeoms[i],var_atempgeoms[j])
AND ST_DWithin(var_atempgeoms[i],var_atempgeoms[j], var_buf*10)
)
WHERE j IS NULL
ORDER BY 1, 2
) As foo ) );
IF ST_IsValid(var_tempgeom) AND ST_GeometryType(var_tempgeom) = 'ST_Polygon' THEN
var_tempgeom := ST_ForceSFS(ST_Intersection(var_tempgeom, var_convhull));
IF param_allow_holes THEN
var_param_geom := var_tempgeom;
ELSE
var_param_geom := ST_MakePolygon(ST_ExteriorRing(var_tempgeom));
END IF;
return var_param_geom;
ELSIF ST_IsValid(var_tempgeom) THEN
var_param_geom := ST_ForceSFS(ST_Intersection(var_tempgeom, var_convhull));
END IF;
END IF;

IF ST_GeometryType(var_param_geom) = 'ST_Polygon' THEN
IF NOT param_allow_holes THEN
var_param_geom := ST_MakePolygon(ST_ExteriorRing(var_param_geom));
END IF;
return var_param_geom;
END IF;
           var_cent := ST_Centroid(var_param_geom);
           IF (ST_XMax(var_enline) - ST_XMin(var_enline) ) > var_buf AND (ST_YMax(var_enline) - ST_YMin(var_enline) ) > var_buf THEN
                   IF ST_Dwithin(ST_Centroid(var_convhull) , ST_Centroid(ST_Envelope(var_param_geom)), var_buf/2) THEN
               -- If the geometric dimension is > 1 and the object is symettric (cutting at centroid will not work -- offset a bit)
                       var_cent := ST_Translate(var_cent, (ST_XMax(var_enline) - ST_XMin(var_enline))/1000,  (ST_YMAX(var_enline) - ST_YMin(var_enline))/1000);
                   ELSE
                       -- uses closest point on geometry to centroid. I can't explain why we are doing this
                       var_cent := ST_ClosestPoint(var_param_geom,var_cent);
                   END IF;
                   IF ST_DWithin(var_cent, var_enline,var_buf) THEN
                       var_cent := ST_centroid(ST_Envelope(var_param_geom));
                   END IF;
                   -- break envelope into 4 triangles about the centroid of the geometry and returned the clipped geometry in each quadrant
                   FOR i in 1 .. 4 LOOP
                      var_geoms[i] := ST_MakePolygon(ST_MakeLine(ARRAY[ST_PointN(var_enline,i), ST_PointN(var_enline,i+1), var_cent, ST_PointN(var_enline,i)]));
                      var_geoms[i] := ST_ForceSFS(ST_Intersection(var_param_geom, ST_Buffer(var_geoms[i],var_buf)));
                      IF ST_IsValid(var_geoms[i]) THEN
                           
                      ELSE
                           var_geoms[i] := ST_BuildArea(ST_MakeLine(ARRAY[ST_PointN(var_enline,i), ST_PointN(var_enline,i+1), var_cent, ST_PointN(var_enline,i)]));
                      END IF;
                   END LOOP;
                   var_tempgeom := ST_Union(ARRAY[ST_ConvexHull(var_geoms[1]), ST_ConvexHull(var_geoms[2]) , ST_ConvexHull(var_geoms[3]), ST_ConvexHull(var_geoms[4])]);
                   --RAISE NOTICE 'Curr vex % ', ST_AsText(var_tempgeom);
                   IF ST_Area(var_tempgeom) <= var_newarea AND ST_IsValid(var_tempgeom)  THEN --AND ST_GeometryType(var_tempgeom) ILIKE '%Polygon'
                       
                       var_tempgeom := ST_Buffer(ST_ConcaveHull(var_geoms[1],least(param_pctconvex + param_pctconvex/var_div),true),var_buf, 'quad_segs=2');
                       FOR i IN 1 .. 4 LOOP
                           var_geoms[i] := ST_Buffer(ST_ConcaveHull(var_geoms[i],least(param_pctconvex + param_pctconvex/var_div),true), var_buf, 'quad_segs=2');
                           IF ST_IsValid(var_geoms[i]) Then
                               var_tempgeom := ST_Union(var_tempgeom, var_geoms[i]);
                           ELSE
                               RAISE NOTICE 'Not valid % %', i, ST_AsText(var_tempgeom);
                               var_tempgeom := ST_Union(var_tempgeom, ST_ConvexHull(var_geoms[i]));
                           END IF;
                       END LOOP;

                       --RAISE NOTICE 'Curr concave % ', ST_AsText(var_tempgeom);
                       IF ST_IsValid(var_tempgeom) THEN
                           var_resultgeom := var_tempgeom;
                       END IF;
                       var_newarea := ST_Area(var_resultgeom);
                   ELSIF ST_IsValid(var_tempgeom) THEN
                       var_resultgeom := var_tempgeom;
                   END IF;

                   IF ST_NumGeometries(var_resultgeom) > 1  THEN
                       var_tempgeom := _ST_ConcaveHull(var_resultgeom);
                       IF ST_IsValid(var_tempgeom) AND ST_GeometryType(var_tempgeom) ILIKE 'ST_Polygon' THEN
                           var_resultgeom := var_tempgeom;
                       ELSE
                           var_resultgeom := ST_Buffer(var_tempgeom,var_buf, 'quad_segs=2');
                       END IF;
                   END IF;
                   IF param_allow_holes = false THEN
                   -- only keep exterior ring since we do not want holes
                       var_resultgeom := ST_MakePolygon(ST_ExteriorRing(var_resultgeom));
                   END IF;
               ELSE
                   var_resultgeom := ST_Buffer(var_resultgeom,var_buf);
               END IF;
               var_resultgeom := ST_ForceSFS(ST_Intersection(var_resultgeom, ST_ConvexHull(var_param_geom)));
           ELSE
               -- dimensions are too small to cut
               var_resultgeom := _ST_ConcaveHull(var_param_geom);
           END IF;
           RETURN var_resultgeom;
END;
$span$
LANGUAGE 'plpgsql'
IMMUTABLE
RETURNS NULL ON NULL INPUT
SECURITY INVOKER
COST 100;

COMMENT ON FUNCTION public.st_concavehull(param_geom public.geometry, param_pctconvex double precision, param_allow_holes boolean)
IS 'args: geomA, target_percent, allow_holes=false - The concave hull of a geometry represents a possibly concave geometry that encloses all geometries within the set. You can think of it as shrink wrapping.';

This file was generated with SQL Manager for PostgreSQL (www.pgsqlmanager.com) at 13/03/2014 13:23