Shapely vira coordenadas lat / long

Oct 16 2020

Eu tenho um dataframe do pandas contendo um objeto de geometria Shapely. O objeto bem torneado (círculo) foi calculado a partir do seguinte SOhttps://stackoverflow.com/questions/27431528/find-the-intersection-between-two-geographical-data-points, onde tenho lat, long e radius como entrada.

import pandas as pd
from shapely.geometry import Point
from shapely.ops import transform
import pyproj
from functools import partial

WGS84 = pyproj.Proj('epsg:4326')

def latlonbuffer(lat, lon, radius_m):
    proj4str = '+proj=aeqd +lat_0=%s +lon_0=%s +x_0=0 +y_0=0' % (lat, lon)
    AEQD = pyproj.Proj(proj4str) # azimuthal equidistant
    project = partial(pyproj.transform, AEQD, WGS84)
    return transform(project, Point(0, 0).buffer(radius_m))

Por exemplo, usando a função acima, criei uma tabela como a abaixo, com uma coluna central adicional

geometry, centroid
POLYGON((26.48306 50.09625, 26.47916 50.09604, ..)), ((26.48307336330026, 50.052005610561245))

Mas quando eu converter esta tabela para geopandas e plotar no mplleaflet

import geopandas as gpd
import mplleaflet
import matplotlib.pyplot as plt

df_gpd = gpd.GeoDataFrame(df)
ax = df_gpd.plot(figsize=(20,20), color='r')

mplleaflet.display(fig=ax.figure)

ele exibe o polígono em um local diferente do esperado.

Abaixo está o que está plotado

E esta é a localização do centroide do círculo esperado

Eu descobri que embora o polígono mostre as coordenadas no formato (lat, long), quando mplleaflet está plotando, ele troca as coordenadas para (long, lat).

O objeto Polygon está mantendo as coordenadas no formato (long, lat) e mplleaflet está invertendo as coordenadas para plotar corretamente no mapa?

O resultado que estou esperando é um círculo centralizado em torno do centróide (26,48307336330026, 50,052005610561245) que está no formato (lat, long).

tldr: Eu gostaria de saber se o mplleaflet inverte as coordenadas para (long, lat) ao plotar o objeto Polygon no Geopandas, ou o objeto Polygon usando os latlonbuffer()resultados da função no formato (long, lat).

Atualização :

O problema era com o pyproj invertendo as coordenadas durante a transformação. Podemos manter o objeto Polygon no padrão GIS (long, lat) adicionando always_xy = True no parâmetro de transformação.

import pandas as pd
from shapely.geometry import Point
from shapely.ops import transform
import pyproj
from functools import partial

WGS84 = pyproj.Proj('epsg:4326')

def latlonbuffer(lat, lon, radius_m):
    proj4str = '+proj=aeqd +lat_0=%s +lon_0=%s +x_0=0 +y_0=0' % (lat, lon)
    AEQD = pyproj.Proj(proj4str) # azimuthal equidistant
    project = partial(pyproj.transform, AEQD, WGS84, always_xy=True)
    return transform(project, Point(0, 0).buffer(radius_m))

Usar isso produz

POLYGON((50.09625 26.48306, 50.09604 26.47916, ..))

que pode ser plotado diretamente usando mplleaflet sem inverter as coordenadas.

Respostas

4 snowman2 Oct 16 2020 at 10:56

Acredito que o problema seja devido a mudanças na ordem dos eixos no PROJ:

https://pyproj4.github.io/pyproj/stable/gotchas.html#axis-order-changes-in-proj-6