Extrapolando una línea en PostGIS

19

Estoy tratando de extrapolar desde un segmento de línea para encontrar un punto en la línea pero un tercio del camino 'atrás', es decir, tratando de encontrar un punto new, puntos dados Ay Babajo:

ingrese la descripción de la imagen aquí

Dada una línea, puedo interpolarla para obtener una posición en cualquier porcentaje en particular:

=# select st_line_interpolate_point(
   st_makeline('0101000020E6100000300DC347C49418C03EE8D9ACFAA44A40', 
               '0101000020E6100000FB743C66A03218C0CDCCCCCCCC7C4A40'), 
   0.333);
0101000020E6100000ED45B41D537718C069C6A2E9EC984A40

Intenté ingresar un número negativo para encontrar un punto a lo largo de la línea en la dirección opuesta, pero eso falla ya que el argumento de interpolación tiene que estar en el rango [0, 1]

Pensé en escalar primero la línea, pero eso no usa el centro de la línea como el origen, por lo que es inútil para mis propósitos.

EoghanM
fuente

Respuestas:

21

Otra forma en que he resuelto un problema similar anteriormente es dividirlo en los siguientes pasos.

-- get the points A and B given a line L
A := ST_STARTPOINT(L);
B := ST_ENDPOINT(L);

-- get the bearing from point B --> A
azimuth := ST_AZIMUTH(B,A);

-- get the length of the line A --> B
length := ST_DISTANCE(A,B);
newlength := length + (length * (1/3));   -- increase the line length by 1/3

-- create a new point 1/3 as far away from A as B is from A
newpoint := ST_TRANSLATE(A, sin(azimuth) * newlength, cos(azimuth) * newlength);

EDITAR: se arregló la asignación de nueva longitud para que sea 1 1/3 de la longitud, en lugar de 1/3 y se cambió A y B para que coincida con el diagrama.

Jayden
fuente
¡Tenemos un ganador! Mucho más comprensible.
EoghanM
Esto es genial
Nathan W
Gracias. Originalmente tenía este fragmento de un trabajo que estaba haciendo interpolando manualmente entre líneas de contorno, resultó que era una pérdida de tiempo lo que estaba tratando de hacer con los contornos, ¡pero me alegra que este fragmento haya ayudado a alguien más! :)
Jayden
6

Lo hemos resuelto con:

F = 1.3333
st_affine(A, F, 0, 
             0, F, 
            (F-1)*-st_x(st_line_interpolate_point(st_makeline(A, B), 0.5)), 
            (F-1)*-st_y(st_line_interpolate_point(st_makeline(A, B), 0.5))
          )

Explicación:

(2-d) Escale el punto de inicio en un factor de 1.3333, tomando el punto medio del segmento de línea como el origen de la escala.

¡Saca el papel cuadriculado!

http://en.wikipedia.org/wiki/Affine_transformation

EoghanM
fuente
2

Escribí una función para esto:

CREATE OR REPLACE FUNCTION st_extend (
    geom geometry,
    head_rate double precision,
    head_constant double precision,
    tail_rate double precision,
    tail_constant double precision)
  RETURNS geometry AS
$BODY$
-- Extends a linestring.
-- First segment get extended by length * head_rate + head_constant.
-- Last segment get extended by length * tail_rate + tail_constant.
--
-- References:
-- http://blog.cleverelephant.ca/2015/02/breaking-linestring-into-segments.html
-- /gis//a/104451/44921
-- /gis//a/16701/44921
WITH segment_parts AS (
SELECT
(pt).path[1]-1 as segment_num
,
CASE
WHEN
  (nth_value((pt).path, 2) OVER ()) = (pt).path
AND
  (last_value((pt).path) OVER ()) = (pt).path
THEN
  3
WHEN
  (nth_value((pt).path, 2) OVER ()) = (pt).path
THEN
  1
WHEN
  (last_value((pt).path) OVER ()) = (pt).path
THEN
  2
ELSE
  0
END AS segment_flag
,
(pt).geom AS a
,
lag((pt).geom, 1, NULL) OVER () AS b
FROM ST_DumpPoints($1) pt
)
,
extended_segment_parts
AS
(
SELECT
  *
  ,
  ST_Azimuth(a,b) AS az1
  ,
  ST_Azimuth(b,a) AS az2
  ,
  ST_Distance(a,b) AS len
FROM
segment_parts
where b IS NOT NULL
)
,
expanded_segment_parts
AS
(
SELECT
  segment_num
  ,
  CASE
  WHEN
    bool(segment_flag & 2)
  THEN
    ST_Translate(b, sin(az2) * (len*tail_rate+tail_constant), cos(az2) * (len*tail_rate+tail_constant))
  ELSE
    a
  END
  AS a
  ,
  CASE
  WHEN
    bool(segment_flag & 1)
  THEN
    ST_Translate(a, sin(az1) * (len*head_rate+head_constant), cos(az1) * (len*head_rate+head_constant))
  ELSE
    b
  END
  AS b
FROM extended_segment_parts
)
,
expanded_segment_lines
AS
(
SELECT
  segment_num
  ,
  ST_MakeLine(a, b) as geom
FROM
expanded_segment_parts
)
SELECT
  ST_LineMerge(ST_Collect(geom ORDER BY segment_num)) AS geom
FROM expanded_segment_lines
;
$BODY$
LANGUAGE sql;

Uso:

SELECT st_extend(
st_makeline(
  '0101000020E6100000300DC347C49418C03EE8D9ACFAA44A40', 
  '0101000020E6100000FB743C66A03218C0CDCCCCCCCC7C4A40'
),
1.333::double precision,
0::double precision,
1::double precision,
0::double precision
);

Tenga en cuenta que esto produce la cadena lineal más larga pero no el punto final.

Código en GitHub Gist (si votas aquí también agradecería una estrella)

Descripción de los parámetros (en caso de que los haya perdido en el comentario de la función sql):

  • La longitud del primer segmento será original_length * head_rate + head_constant.
  • Si desea que se duplique, la tasa de carga es 2, la constante es 0.
  • En Hungría normalmente usamos la proyección EOV que se basa en medidores. Entonces, si quiero agregar 2 metros al final de la línea, configuro tail: rate en 1 y tail_constant en 2.
SzieberthAdam
fuente
Esto funciona muy bien ¿Puedes agregar alguna información sobre head_rate, head_constant, tail_rate y tail_constant? No se explican aquí ni en su GitHub. Supongo que la tasa de cabecera = factor de escala para la extensión de línea después del punto final y la tasa de cola = factor de escala para la extensión de línea antes del punto de inicio. ¿Cómo funcionan las constantes? ¿Qué efecto tienen?
jbalk
Está en el comentario. La longitud del primer segmento será original_length * head_rate + head_constant. Si desea que se duplique, la tasa de carga es 2, la constante es 0. En Hungría, normalmente utilizamos la proyección EOV, que se basa en medidores. Entonces, si quiero agregar 2 metros al final de la línea, configuro tail: rate a 1 y tail_constant a 2.
SzieberthAdam
¡Gracias! Y gracias por compartir esta función. Funciona perfectamente y se ejecuta rápidamente.
jbalk