pyproj в Python возвращает неверные координаты для NAD_1927_StatePlane_Wisconsin_South_FIPS_4803

Я новичок в использовании геопространственных данных. У меня возникли проблемы с получением pyproj правильной долготы и широты из файла с другой системой координат. Я работаю с файлом, в котором границы города Милуоки находятся в координатах State Plane NAD_1927_StatePlane_Wisconsin_South_FIPS_4803, и мне нужны координаты WGS84, которые должны быть обычными долготой и широтой. Я хочу преобразовать это, потому что было бы проще работать с такими инструментами, как folium, которые я использую.

https://data.milwaukee.gov/dataset/corporate-boundary Первые координаты в файле shp есть [2516440.160073474, 440481.9301029742]

Он имеет файл prj:

PROJCS["NAD_1927_StatePlane_Wisconsin_South_FIPS_4803",
       GEOGCS["GCS_North_American_1927",
              DATUM["D_North_American_1927",
                    SPHEROID["Clarke_1866",6378206.4,294.9786982]],
              PRIMEM["Greenwich",0.0],
              UNIT["Degree",0.017453292519943295]],
       PROJECTION["Lambert_Conformal_Conic"],
       PARAMETER["False_Easting",2000000.0],
       PARAMETER["False_Northing",0.0],
       PARAMETER["Central_Meridian",-90.0],
       PARAMETER["Standard_Parallel_1",42.73333333333333],
       PARAMETER["Standard_Parallel_2",44.06666666666667],
       PARAMETER["Latitude_Of_Origin",42.0],
       UNIT["Foot_US",0.30480060960121924]]

Согласно https://spatialreference.org, это будет иметь код EPSG6948, но pyproj не распознает его, поэтому я построить пользовательскую строку. После просмотра документации и примеров я придумал код:

import pyproj
from pyproj import Transformer
cust = "+proj=lcc +ellps=clrk66 +datum=NAD27 +lat_1=42.73333333333333 +lat_2=44.06666666666667 +lat_0=42.0 +lon_0=-90.0 +x_0=2000000.0 +y_0=0.0 +units=us-ft +no_defs"
t = Transformer.from_crs(cust, "WGS84")
x1,y1 = [2516440.160073474, 440481.9301029742]
lat, lon = t.transform(x1, y1)
print(lat, ", ", lon)

Согласно сайту https://www.earthpoint.us/StatePlane.aspx этот код должен вернуть координаты: Градусы широты и долготы 43,1900328, -087,9451397

Но код выводит: 42.206953018562416 , -105.00939853775037 где-то в Вайоминге, а не в Милуоки, Висконсин. Я проверял это несколько раз и не вижу, что я делаю неправильно. Если кто-нибудь может помочь, я был бы очень признателен.


python gis geospatial pyproj coordinate-transformation
person JJF    schedule 25.03.2021    source источник

Ответы (1)


arrow_upward
0
arrow_downward

Вы можете использовать строку WKT для создания преобразователя:

>>> from pyproj import Transformer
>>> wkt = """PROJCS["NAD_1927_StatePlane_Wisconsin_South_FIPS_4803",
...        GEOGCS["GCS_North_American_1927",
...               DATUM["D_North_American_1927",
...                     SPHEROID["Clarke_1866",6378206.4,294.9786982]],
...               PRIMEM["Greenwich",0.0],
...               UNIT["Degree",0.017453292519943295]],
...        PROJECTION["Lambert_Conformal_Conic"],
...        PARAMETER["False_Easting",2000000.0],
...        PARAMETER["False_Northing",0.0],
...        PARAMETER["Central_Meridian",-90.0],
...        PARAMETER["Standard_Parallel_1",42.73333333333333],
...        PARAMETER["Standard_Parallel_2",44.06666666666667],
...        PARAMETER["Latitude_Of_Origin",42.0],
...        UNIT["Foot_US",0.30480060960121924]]"""
>>> transformer = Transformer.from_crs(wkt, "WGS84")
>>> transformer.transform(2516440.160073474, 440481.9301029742)
(43.19212933447268, -88.06334977811528)

Примечание: я думаю, что это код EPSG, который вам нужен:

<Projected CRS: EPSG:32054>
Name: NAD27 / Wisconsin South
Axis Info [cartesian]:
- X[east]: Easting (US survey foot)
- Y[north]: Northing (US survey foot)
Area of Use:
- name: United States (USA) - Wisconsin - counties of Adams; Calumet; Columbia; Crawford; Dane; Dodge; Fond Du Lac; Grant; Green; Green Lake; Iowa; Jefferson; Juneau; Kenosha; La Crosse; Lafayette; Manitowoc; Marquette; Milwaukee; Monroe; Ozaukee; Racine; Richland; Rock; Sauk; Sheboygan; Vernon; Walworth; Washington; Waukesha; Waushara; Winnebago.
- bounds: (-91.43, 42.48, -86.95, 44.33)
Coordinate Operation:
- name: Wisconsin CS27 South zone
- method: Lambert Conic Conformal (2SP)
Datum: North American Datum 1927
- Ellipsoid: Clarke 1866
- Prime Meridian: Greenwich


>>> transformer = Transformer.from_crs("EPSG:32054", "WGS84")
>>> transformer.transform(2516440.160073474, 440481.9301029742)
(43.19212933447268, -88.06334977811528)
person snowman2    schedule 26.03.2021
comment
Спасибо. Это помогло - person JJF; 27.03.2021