Mengkonversi koordinat yang diproyeksikan ke lat / lon menggunakan Python?

51

Saya menganggap ini adalah pertanyaan dasar tetapi sepertinya saya tidak dapat menemukan atau mengenali solusinya.

Situs ini kembali

Point:
X: -11705274.6374
Y: 4826473.6922

ketika Anda mencari dengan nilai kunci pertama 000090 sebagai contoh. Saya kira ini adalah referensi spasial dan saya agak mengerti apa itu.

Saya mencari instruksi atau contoh cara mengonversikan ini ke Latitude dan bujur menggunakan Python.

Vincent
sumber

Jawaban:

100

Cara paling sederhana untuk mengubah koordinat dalam Python adalah pyproj , yaitu antarmuka Python ke pustaka PROJ.4 . Faktanya:

from pyproj import Proj, transform

inProj = Proj(init='epsg:3857')
outProj = Proj(init='epsg:4326')
x1,y1 = -11705274.6374,4826473.6922
x2,y2 = transform(inProj,outProj,x1,y1)
print x2,y2

kembali -105.150271116 39.7278572773

Antonio Falciano
sumber
4
Ya. Pyproj sepanjang jalan.
sgillies
Ini bekerja untuk saya
lenhhoxung
36

Secara default situs yang Anda tautkan menggunakan Sistem Referensi Spasial EPSG 3857 (WGS84 Web Mercator). Saya menemukan informasi ini di sini .

Anda dapat menentukan Sistem Referensi Spasial lain dengan memasukkan EPSG yang diinginkan ke dalam formulir di bawah Spatial Referenceatau Anda dapat mengonversi koordinat yang dikembalikan dengan Python.

Misalnya Anda dapat menggunakan binding GDAL Python untuk mengubah titik ini dari sistem koordinat yang diproyeksikan (EPSG 3857) ke sistem koordinat geografis (EPSG 4326).

import ogr, osr

pointX = -11705274.6374 
pointY = 4826473.6922

# Spatial Reference System
inputEPSG = 3857
outputEPSG = 4326

# create a geometry from coordinates
point = ogr.Geometry(ogr.wkbPoint)
point.AddPoint(pointX, pointY)

# create coordinate transformation
inSpatialRef = osr.SpatialReference()
inSpatialRef.ImportFromEPSG(inputEPSG)

outSpatialRef = osr.SpatialReference()
outSpatialRef.ImportFromEPSG(outputEPSG)

coordTransform = osr.CoordinateTransformation(inSpatialRef, outSpatialRef)

# transform point
point.Transform(coordTransform)

# print point in EPSG 4326
print point.GetX(), point.GetY()

Ini mengembalikan untuk titik Anda koordinat -105.150271116 39.7278572773.

ustroetz
sumber
8

afalciano memiliki jawaban yang tepat tetapi ingin memasukkan varian penggunaan pyproj.

Ini tidak mengharuskan Anda tahu PROJ4 tali dan sedikit lebih cepat.

import pyproj
p = pyproj.Proj("+proj=merc +lon_0=0 +k=1 +x_0=0 +y_0=0 +a=6378137 +b=6378137 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs")
lon, lat = p(x, y, inverse=True)
print lat, lon
Marcel Wilson
sumber
2
Anda tidak memerlukan string proj4, gantikan baris kedua p = pyproj.Proj(init='epsg:3857')dan hasilnya sama.
alphabetasoup
1
Hasilnya sama tetapi terakhir saya memeriksa ini sedikit lebih cepat.
Marcel Wilson
6

Outputnya bukan sistem referensi spasial / koordinat , ini adalah sepasang koordinat. Anda perlu tahu apa referensi spasial untuk memproyeksi ulang koordinat.

Namun, itu tidak diperlukan dalam kasus ini. Cukup berikan referensi spasial keluaran yang sesuai ke layanan dan akan mengembalikan koordinat di Lon / Lat.

Ini adalah halaman dengan koordinat keluaran dalam format Lon / Lat menggunakan sistem referensi spasial geografis WGS-84 ( EPSG 4326 ).

pengguna2856
sumber
4

Mencoba kode yang disarankan oleh Marcel Wilson dan itu lebih cepat:

from pyproj import Proj, transform
import time
import pyproj


# Test 1 - 0.0006158 s
start=time.time()
inProj = Proj(init='epsg:3857')
outProj = Proj(init='epsg:4326')
x1,y1 = -11705274.6374,4826473.6922
x2,y2 = transform(inProj,outProj,x1,y1)
end=time.time()
print(y2,x2)
print('%.7f' % (end-start))

# Test 2 - 0.0000517 s --- factor 11,9
start=time.time()
x,y = -11705274.6374,4826473.6922
p = pyproj.Proj("+proj=merc +lon_0=0 +k=1 +x_0=0 +y_0=0 +a=6378137 +b=6378137 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs")
lon, lat = p(x, y, inverse=True)
end=time.time()
print(lat, lon)
print('%.7f' % (end-start))
-----------------

39.72785727727918 -105.15027111593008
0.0006158
39.72785727727918 -105.15027111593008
0.0000517
pengguna1174460
sumber
1

Saya menemukan posting ini ketika mencari cara untuk melakukan ini dalam QGIS. Seperti dijelaskan di sini , metode yang digunakan terlihat seperti ini:

def convertProjection(self,x,y,from_crs,to_crs):
    crsSrc = QgsCoordinateReferenceSystem(from_crs)
    crsDest = QgsCoordinateReferenceSystem(to_crs)
    xform = QgsCoordinateTransform(crsSrc, crsDest)
    pt = xform.transform(QgsPoint(x,y))
    return pt.x, pt.y

# Remove the "EPSG:" part
from_crs = 3857
to_crs = 4326
x = -11705274.6374    
y = 4826473.6922
lon, lat = self.convertProjection(x,y,from_crs, to_crs)
Toivo Säwén
sumber
2
Catatan, ada perubahan API yang melanggar di QGIS 3, jadi jika menggunakan v3.xAnda harus menggunakanxform = QgsCoordinateTransform(crsSrc, crsDest, QgsProject.instance())
Jonny
0

Harap perhatikan bahwa transformfungsi accept pyprojjuga arrays, yang cukup berguna ketika datang ke bingkai data

import pandas as pd
from pyproj import Proj, transform

df = pd.DataFrame({'x': [-11705274.6374]*100, 
                   'y': [4826473.6922]*100})
inProj, outProj = Proj(init='epsg:3857'), Proj(init='epsg:4326')
df['x2'], df['y2'] = transform(inProj, outProj, df['x'].tolist(), df['y'].tolist())
J. Doe
sumber
-1
import cv2
import numpy as np
def onMouse(event, x, y, flags, param):
    if event == cv2.EVENT_LBUTTONDOWN:
       # draw circle here (etc...)
       print('x = %f, y = %f'%((x*2/100),(y*2/100)))
a=float(input("enter length of original dimension in cm"))
b=float(input("enter width of original dimension in cm"))
print("origional image coordinates")
im=cv2.imread('mask1.jpg',0)
re_im=cv2.resize(im,(round(a*100/2),round(b*100/2)))
cv2.imshow('image',re_im)
cv2.setMouseCallback('image', onMouse)
cv2.waitKey(0)
cv2.destroyAllWindows()
Subhash Verma
sumber