# Converting state plane coordinates to latitude/longitude

I'm trying to convert pairs of x,y coordinates to latitude/longitude using python. The projection is "NAD_1983_HARN_StatePlane_Oregon_North_FIPS_3601_Feet_Intl" and a sample x,y pair is 7625659, 714714. So far, I have found this question and have tried using the code in the first solution:

``````from pyproj import Proj, transform

inProj = Proj(init='epsg:????')
outProj = Proj(init='epsg:????')
x1,y1 = 7625659,714714
x2,y2 = transform(inProj,outProj,x1,y1)
print x2,y2
``````

However, I'm having trouble finding an epsg code for the location. I came across this but I can't seem to find the epsg code for that area. All I can see is the sr-org code. Also, it seems strange that the location on the map is not in Oregon. I was able to convert it (the output looked reasonable) using the tool at http://beta.ngs.noaa.gov/gtkweb/ using datum=NAD83, projection for conversion=SPC, Norting=714714, Easting=7625659, Units=International feet and SPC Zone=OR N-3601. How should I convert this projection to latitude and longitude?

Edit: I tried using the tool at http://prj2epsg.org/search which gave me 2913 (I also tried 3646 per BradHards suggestion, it came up the same) for inProj and 4326 (BradHards) for outProj, but it's giving me a location (-62.0472510133 32.7243167795) which is in the middle of the ocean, not in Portland. The lat/lon from http://beta.ngs.noaa.gov/gtkweb/ is (45.605018016, -122.752592697) which is in Portland (Sorry, I can't comment on the site...).

• Have a look at this, spatialreference.org/ref/epsg/4152 – Fezter Nov 16 '16 at 3:02
• The outProj is probably meant to be `epsg:4326` if you want WGS-84 longitude / latitude. It might be something else if you want a different ellipsoid. Perhaps `epsg:3646` for the input. – BradHards Nov 16 '16 at 3:08

## 1 Answer

I needed to add `preserve_units = True` and the epsg codes provided by BradHards.

``````from pyproj import Proj, transform

inProj = Proj(init='epsg:3646', preserve_units = True)
outProj = Proj(init='epsg:4326')
x1,y1 = 7625659,714714
x2,y2 = transform(inProj,outProj,x1,y1)
print x2,y2
``````

This returns lat/lon `(-122.752652268, 45.6050168495)`. Strangely enough http://beta.ngs.noaa.gov/gtkweb/ gives me lon/lat `(-122.752592697, 45.605018016)` which is a couple of feet over. I'm going to assume that's a round off error somewhere and it shouldn't make too much of a difference.