[GIS] Correctly Transforming (Unprojecting?) from Lambert Conformal Conic

coordinate systemgeoprocessingshapefile

I started out by posting https://stackoverflow.com/questions/6963991/how-can-i-transform-the-coordinates-of-a-shapefile

In short, I found a shapefile here, and I want to be able to convert that to Latitude/Longitude (Geodetic?) coordinates that I can store in my PostGIS database. I am having trouble working with the Ruby libraries, but that is a separate issue. First, I need to correctly understand the projection.

In the linked zip file, there is a .prj file, which contains the following (It is only one line, but I broke it up for readability):

PROJCS["NAD_1983_StatePlane_California_III_FIPS_0403_Feet",
GEOGCS["GCS_North_American_1983",DATUM["D_North_American_1983",
SPHEROID["GRS_1980",6378137.0,298.257222101]],PRIMEM["Greenwich",0.0],UNIT["Degree",0.0174532925199433]],
PROJECTION["Lambert_Conformal_Conic"],PARAMETER["False_Easting",6561666.666666666],PARAMETER["False_Northing",1640416.666666667],
PARAMETER["Central_Meridian",-120.5],PARAMETER["Standard_Parallel_1",37.06666666666667],
PARAMETER["Standard_Parallel_2",38.43333333333333],PARAMETER["Latitude_Of_Origin",36.5],
UNIT["Foot_US",0.3048006096012192]]

I was unable to find tools that translate .prj files into a projection string, so I used this as a reference and came up with:

+proj=lcc +datum=NAD83 +ellps=GRS80 +lat_1=37.06666666666667 +lat_2=38.43333333333333 +lat_0=36.5 +lon_0=-120.5 +x_0=6561666.666666666 +y_0=1640416.666666667 +units=us-ft +pm=greenwich +no_defs

I happen to know that the points in this shapefile should correspond to locations in San Francisco, so I could test my work:

cs2cs -f "%.8f" +proj=lcc +datum=NAD83 +ellps=GRS80 +lat_1=37.06666666666667 +lat_2=38.43333333333333 +lat_0=36.5 +lon_0=-120.5 +x_0=6561666.666666666 +y_0=1640416.666666667 +units=us-ft +pm=greenwich +no_defs +to +proj=lonlat +datum=WGS84 +ellps=WGS84
6011287.4999795845 2100857.2499904726
-164.18117482   17.72087079 0.00000000

That's about 1000km west of Hawaii, so clearly that's not right. I tried many different things, but the closest I got was by stripping the +units=us-ft off.

cs2cs -f "%.8f" +proj=lcc +datum=NAD83 +ellps=GRS80 +lat_1=37.06666666666667 +lat_2=38.43333333333333 +lat_0=36.5 +lon_0=-120.5 +x_0=6561666.666666666 +y_0=1640416.666666667 +pm=greenwich +no_defs +to +proj=lonlat +datum=WGS84 +ellps=WGS84
6011287.4999795845 2100857.2499904726
-126.98864867   40.47509388 0.00000000

That's about 500 km NW of San Francisco.

What am I missing here?
If I import the shapefile into the online ArcGIS Explorer, it understands it perfectly.

Best Answer

I believe that the correct string for this projection is: 

+proj=lcc +lat_1=37.06666666666667 +lat_2=38.43333333333333 +lat_0=36.5 +lon_0=-120.5 +x_0=2000000 +y_0=500000.0000000002 +ellps=GRS80 +datum=NAD83 +to_meter=0.3048006096012192 +no_defs

See: http://spatialreference.org/ref/esri/102643/

Related Question