Échantillonnage ponctuel le long d'un littoral enveloppant les poteaux avec PostGIS

Je travaille sur une tâche qui m'oblige à obtenir des points d'échantillonnage tous les 1000 km le long des côtes et j'ai rencontré un problème avec l'Antarctique. D'après ce que je peux dire, cela semble être un problème avec l'utilisation de la géométrie dans les fonctions, alors que je pense vraiment que la géographie devrait être utilisée pour cette opération.

En utilisant la fonction de cette question très similaire , je suis en mesure de produire un résultat qui ressemble à ceci: .

Comme vous pouvez le voir, ST_AddMeasure()et ST_LocateAlong()ne semblez pas traiter la géométrie de manière sphérique, ce qui entraîne de nombreux points qui se trouvent sur le pôle Sud. Un point a même été ajouté sur le clip le long de la ligne de date (côté gauche). Selon la documentation de ces deux fonctions, seule la géométrie peut être utilisée.

Le code utilisé pour générer le polygone et les points peut être trouvé ici , mais c'est le SQL utilisé pour générer les points:

CREATE TABLE atest AS WITH line AS 
  (SELECT
      id,
      ST_ExteriorRing((ST_Dump(geom)).geom) AS geom
    FROM line_sample_test),
linemeasure AS
    (SELECT
        ST_AddMeasure(line.geom, 0, (ST_Length(line.geom))::int) AS linem,
    generate_series(0, (ST_Length(line.geom))::int, 10) AS i
FROM line),

geometries AS (
    SELECT
        i,
        ST_LocateAlong(linem, i) AS geom 
    FROM linemeasure)

SELECT
    * from geometries;

Comment générer des points tous les 1000 km le long de ce littoral?

Comme suggéré dans l'un des commentaires, je voudrais d'abord transformer la géométrie d'entrée en projection stéréographique polaire.

De plus, vous voudrez l'utiliser ST_Bufferdessus (avec une quantité de 0) pour vous débarrasser de la ligne de coupe résultante.

Ainsi, vous obtiendrez le résultat souhaité:

-- ST_Transform(geom,3031) reprojects to south polar stereographic,
-- in meters.  ST_Buffer(...) doesn't change the shape, but removes
-- the cut line to the pole (at 180 degrees).
WITH line AS (
    SELECT ST_ExteriorRing(
        ST_Buffer(ST_Transform(geom, 3031), 0)
    ) AS geom
    FROM line_sample_test
),

-- This just generates a table of numbers.  In this case, from 0
-- to the geometry length, counting by 1,000,000 (1000 km).
linemeasure AS (
    SELECT generate_series(0, ST_Length(geom)::int, 1000000) AS i
    FROM line
),

-- Convert those values to a fraction of the overall length (for
-- use as input to ST_LineInterpolatePoint)
linefraction AS (
    SELECT i / ST_Length(geom) AS fraction
    FROM line, linemeasure
),

-- Do the interpolation
geometries AS (
    SELECT ST_LineInterpolatePoint(l.geom, lf.fraction) AS geom
    FROM linefraction lf, line l
),

-- Convert back to EPSG:4326 (i.e. lat/lon coords)
geometries_4326 AS (
    SELECT ST_Transform(geom, 4326) AS geom FROM geometries
)
SELECT * FROM geometries_4326

Notez que cette requête suppose qu'il n'y a qu'une seule ligne dans le line_sample_testtableau, ajustez-la si nécessaire pour vos données d'entrée réelles.