Skip to content

Pyproj

Python interface to PROJ (cartographic projections and coordinate transformations library)

Initializing CRS

The CRS class can be initialized in many different ways. Here are some examples of initialization. 

from pyproj import CRS
crs = CRS.from_epsg(4326)
crs = CRS.from_string("epsg:4326")
crs = CRS.from_proj4("+proj=latlon")
crs = CRS.from_user_input(4326)

Converting projected coordinates

transform을 사용한 방법: 

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

Proj4 인자를 사용한 방법: 

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

두 좌표간 거리 계산 방법

You can use the library pyproj, that allows you to change the coordinate system. I use EPSG3035 that is convenient for Europe, but depending on where your points are you might have to use a different one. 

import pyproj
import numpy as np

lon1 = 2
lat1 = 48

lon2 = 2
lat2 = 49

x1, y1 = pyproj.transform(wgs84, epsg3035, lon1, lat1)
print(x1, y1)

x2, y2 = pyproj.transform(wgs84, epsg3035, lon2, lat2)
print(x2, y2)

# a Pythagore's theorem is sufficient to compute an approximate distance
distance_m = np.sqrt((x2-x1)**2 + (y2-y1)**2)
print(distance_m)

The error will be typically of a few percents for points in Europe, depending on where they are and how far they are. Have a nice day !

Geod를 사용한 방법 

import pyproj

geod = pyproj.Geod(ellps='WGS84')

for city, coord in cities.items():
    lat0, lon0 = london_coord
    lat1, lon1 = coord

    azimuth1, azimuth2, distance = geod.inv(lon0, lat0, lon1, lat1)
    print(city, distance)
    print('    azimuth', azimuth1, azimuth2)

그 밖의 방법

Geod class

특정 거리만큼 이동시 위치 구하는 방법:

  1. 두 경위도를 사용한 방위각(1->2, 2->1)과 거리 획득.
  2. 특정 위치(경위도)에서 특정 방위각으로 meter 거리만큼 이동했을 때의 결과 위치. 
import pyproj

meter = 30
lons1 = 127.30594507347541
lats1 = 36.63872883265555
lons2 = 127.306055839949
lats2 = 36.638730247186125

geodetic = pyproj.geod.Geod(ellps='WGS84')
az12, az21, dist = geodetic.inv(lons1, lats1, lons2, lats2)  # 1
endlon, endlat, backaz = geodetic.fwd(lons1, lats1, az12, meter)  # 2

See also

Favorite site