Traduzione EPSG: 4326 coordinate lon lat in xy utilizzando pyproj

Aug 19 2020

Sto cercando di elaborare i dettagli dell'utilizzo di pyproj per tradurre le coordinate latlon dal feed USGS Earthquake in coordinate xy corrispondenti a quelle posizioni su una mappa base di una determinata larghezza, altezza e origine in alto a sinistra. Il mio codice finora utilizza ingenuamente PIL / pillow per tradurre posizioni da long lat a xy e disegnare cerchi senza tenere conto della proiezione EPSG: 4326 della mappa di base:

from PIL import Image
from PIL import ImageDraw

# open basemap image file
basemap = Image.open(basemap_path).convert('RGBA)

# resize to desired map size
basemap.thumbnail(width, height, Image.LANCZOS)

# get proportional height
width_bmp, height_bmp = basemap.size

# create background frame and paste basemap on it
img=Image.new('RGB',(width, height), color = '#000000')
img.paste(basemap, (0,0), basemap)
draw = ImageDraw.Draw(img, 'RGBA')

width_scale = width/360
height_scale = height_bmp/180

# usgs data has been parsed into a list
for quake in earthquake_list:
    lon = float(quake["longitude"])
    lat = float(quake["latitude"])
    mag = float(quake["mag"])

    # want to use pyproj to translate coordinates here instead of the following
    cx = (lon + 180) * width_scale
    cy = (90-lat) * height_scale
    r  = scaleRadius(mag)

    # draw earthquake circles
    draw.ellipse((cx-r, cy-r, cx+r, cy+r), fill = colormap(mag))

draw = ImageDraw.Draw(img)
img.save(filepath, quality=100)

Ho esaminato la documentazione e ho la sensazione generale che avrei usato qualcosa del tipo:

cx, cy = pyproj.transform("EPSG:4326", "xy", lon, lat)

con alcune dimensioni di larghezza e altezza. L'esecuzione della riga sopra dà l'errore:

TypeError: p1 must be a Proj class

Aggiorna ottengo valori per:

p = Proj(proj = 'longlat', ellps='WGS84')
cx,cy = p(lon, lat)

In questo caso un input (lon, lat) di (-110, 39)fornisce l'output 1.923, 0.688che non è corretto. Ciò di cui ho bisogno è una funzione di proiezione che prenda la lunghezza, la latitudine, la larghezza e l'altezza dell'immagine e restituisca le coordinate x, y. Lo scriverei da solo, ma non riesco a trovare un'equazione adatta per quella conversione e speravo di trovare questa capacità in proj.

AGGIORNAMENTO Ecco la mappa di base che sto utilizzando.

Ed ecco il disegno risultante. Suppongo che i punti centrali per i cerchi che segnano i terremoti siano disattivati, perché sto usando coordinate cartesiane per lon, lat, ma la mappa di base è stata creata con una proiezione e quindi deve avere una distorsione da x, y coordinate cartesiane.

Inoltre, mi rendo conto che ho bisogno di mettere la mia leggenda su una scala logaritmica. Non ci sono ancora arrivato.

Risposte

5 FSimardGIS Aug 22 2020 at 04:12

Quando sovrappongo la tua mappa di base su una proiezione equirettangolare, si adatta perfettamente. Quindi questo significa che la tua mappa di base utilizza solo la semplice "proiezione" di tracciare le longitudini come x e le latitudini come y. Il che significa che l'unica cosa di cui devi preoccuparti è controllare la dimensione dell'immagine, cosa che fai già nelle tue formule.

cx = (lon + 180) * width_scale
cy = (90-lat) * height_scale

Questi funzioneranno correttamente per le tue esigenze, al pixel più vicino. In effetti, ho verificato i terremoti su USGS in quella data ei tuoi punti sono tutti nella posizione giusta!