Creé una secuencia de comandos pequeña e ingenua que convierte LineStrings de entrada a CompoundCurves en función de algunas heurísticas.
Que hace:
- Corta esquinas afiladas para crear resultados visualmente más atractivos que los datos originales.
- Utiliza plpgsql. No se requieren extensiones adicionales.
- Acepta un "factor de suavizado" opcional entre 0 y 100 además de una geometría.
Lo que no hace:
- Procesos MultiLineStrings. Para cualquier otro tipo de geometría, simplemente devuelve la entrada.
- Utiliza los valores Z y M. Simplemente los deja caer. Use esto solo para fines cartográficos 2D.
- Crea resultados matemáticamente correctos. Los resultados están lejos de ser correctos e incluso pueden ser visualmente poco estéticos en algunos casos (por ejemplo, esquinas afiladas). No lo probé a fondo. ¡Siempre revise los resultados!
- Corre rápido. Estoy seguro de que puede reescribirse en una forma mucho más óptima.
- Real suavizado. Hay algoritmos mucho mejores (p. Ej., Chaiken o los mencionados en la pregunta) para usar para un suavizado real. Esta respuesta se dirige a personas como yo que buscan un enfoque puro de PostGIS creando automáticamente algún tipo de líneas curvas a partir de datos reales.
La secuencia de comandos:
CREATE OR REPLACE FUNCTION CreateCurve(geom geometry, percent int DEFAULT 40)
RETURNS geometry AS
$$
DECLARE
result text;
p0 geometry;
p1 geometry;
p2 geometry;
intp geometry;
tempp geometry;
geomtype text := ST_GeometryType(geom);
factor double precision := percent::double precision / 200;
i integer;
BEGIN
IF percent < 0 OR percent > 100 THEN
RAISE EXCEPTION 'Smoothing factor must be between 0 and 100';
END IF;
IF geomtype != 'ST_LineString' OR factor = 0 THEN
RETURN geom;
END IF;
result := 'COMPOUNDCURVE((';
p0 := ST_PointN(geom, 1);
IF ST_NPoints(geom) = 2 THEN
p1:= ST_PointN(geom, 2);
result := result || ST_X(p0) || ' ' || ST_Y(p0) || ',' || ST_X(p1) || ' ' || ST_Y(p1) || '))';
ELSE
FOR i IN 2..(ST_NPoints(geom) - 1) LOOP
p1 := ST_PointN(geom, i);
p2 := ST_PointN(geom, i + 1);
result := result || ST_X(p0) || ' ' || ST_Y(p0) || ',';
tempp := ST_Line_Interpolate_Point(ST_MakeLine(p1, p0), factor);
p0 := ST_Line_Interpolate_Point(ST_MakeLine(p1, p2), factor);
intp := ST_Line_Interpolate_Point(
ST_MakeLine(
ST_Line_Interpolate_Point(ST_MakeLine(p0, p1), 0.5),
ST_Line_Interpolate_Point(ST_MakeLine(tempp, p1), 0.5)
), 0.5);
result := result || ST_X(tempp) || ' ' || ST_Y(tempp) || '),CIRCULARSTRING(' || ST_X(tempp) || ' ' || ST_Y(tempp) || ',' || ST_X(intp) || ' ' ||
ST_Y(intp) || ',' || ST_X(p0) || ' ' || ST_Y(p0) || '),(';
END LOOP;
result := result || ST_X(p0) || ' ' || ST_Y(p0) || ',' || ST_X(p2) || ' ' || ST_Y(p2) || '))';
END IF;
RETURN ST_SetSRID(result::geometry, ST_SRID(geom));
END;
$$
LANGUAGE 'plpgsql' IMMUTABLE;
Como devuelve curvas en un tipo de geometría, si desea utilizarlo en un SIG como QGIS, debe ajustarlo en funciones PostGIS convirtiéndolos. La sintaxis de uso prevista es:
SELECT ST_AsText(ST_CurveToLine(CreateCurve(geom))) AS geom FROM linestringtable;
Esto sigue siendo un problema abierto en PostGIS (y otras herramientas SIG) como se indica en el libro "PostGIS en acción" en el capítulo 2.2.6 "Geometrías curvas".
Aquí hay algunas referencias a algoritmos y código:
fuente
Puede intentar convertir sus cadenas lineales en curvas con ST_LineToCurve y luego volver a las cadenas lineales con ST_CurveToLine .
Puede establecer el número de segmentos por cuarto de círculo que desee en ST_CurveToLine.
fuente