Guess my EPSG: Difference between revisions
From Wildsong
Jump to navigationJump to search
Brian Wilson (talk | contribs) m Created page with 'Let's say someone has a shapefile and they want to know the appropriate EPSG code for it. How do I do that?? So far I can see that most of the files have a descriptive string fo…' |
Brian Wilson (talk | contribs) mNo edit summary |
||
Line 12: | Line 12: | ||
That's often all I need. | That's often all I need. | ||
Here is the code that works for me... | |||
<pre> | |||
#!/usr/bin/python | |||
import sys,os | |||
import fnmatch | |||
lov = [ | |||
('"NAD_1983_StatePlane_California_II_FIPS_0402"', 'UNIT["Meter', 26942), | |||
('"NAD_1983_StatePlane_California_II_FIPS_0402"', 'UNIT["Foot', 2226), | |||
('NAD_1983_California_Albers', 'UNIT[Meter,', 3310), | |||
('"NAD83 / California Albers"', '', 3310), | |||
('"NAD_1927_StatePlane_California_I_FIPS_0401"', 'UNIT["Foot_US"', 26741), | |||
('PROJCS["Custom",GEOGCS["GCS_North_American_1927",DATUM["D_North_American_1927",SPHEROID["Clarke_1866",6378206.4,294.9786982]],PRIMEM["Greenwich",0],UNIT["Degree",0.0174532925199433]],PROJECTION["Albers"],PARAMETER["False_Easting",0],PARAMETER["False_Northing",-4000000],PARAMETER["Central_Meridian",-120],PARAMETER["Standard_Parallel_1",34],PARAMETER["Standard_Parallel_2",40.5]', 'UNIT["Meter",1]', 3309), | |||
('"Old_Hawaiian_StatePlane_Hawaii_3_FIPS_5103"', 'UNIT["Foot_US"', 3563), | |||
('"NAD_1983_StatePlane_Nevada_West_FIPS_2703_Feet"', '', 102309), | |||
('"NAD_1983_HARN_Lambert_Conformal_Conic"', 'UNIT["User_Defined_Unit",0.3048', 2913), | |||
('"NAD_1983_HARN_StatePlane_Oregon_North_FIPS_3601_Feet_Intl"', '', 2913), | |||
('"NAD_1983_HARN_StatePlane_Oregon_North_FIPS_3601"', 'UNIT["Foot"', 2913), | |||
('"NAD_1983_HARN_StatePlane_Oregon_North_FIPS_3601_IntlFeet"', '', 2913), | |||
('"NAD_1983_HARN_StatePlane_Oregon_South_FIPS_3602_Feet_Intl"', 'UNIT["Foot', 2914), | |||
('"NAD_1983_HARN_StatePlane_Oregon_South_FIPS_3602_AltUnits"', 'UNIT["User_Defined_Unit",0.3047999902464]', 2914), | |||
# DOGAMI uses this and it's sort of a guess | |||
('PROJCS["_MI_0",GEOGCS["NAD83",DATUM["D_NAD83",SPHEROID["Geodetic_Reference_System_of_1980",6378137,298.2572221009113]],PRIMEM["Greenwich",0],UNIT["Degree",0.017453292519943295]],PROJECTION["Lambert_Conformal_Conic"],PARAMETER["standard_parallel_1",43],PARAMETER["standard_parallel_2",45.5],PARAMETER["latitude_of_origin",41.75],PARAMETER["central_meridian",-120.5],PARAMETER["false_easting",1312335.958],PARAMETER["false_northing",0]', 'UNIT["Foot (International)",0.3048]]', 3644), | |||
('"NAD_1983_Oregon_Statewide_Lambert"', 'UNIT["Meter', 2991), | |||
('"NAD_1983_Oregon_Statewide_Lambert_Feet_Intl"', '', 2992), | |||
('"OGICdefaultAG82"', 'UNIT["Foot"', 2992), | |||
('"NAD_1983_HARN_Lambert_Conformal_Conic"', 'UNIT["Foot"', 2994), | |||
('"NAD_1927_StatePlane_Oregon_North_FIPS_3601"', 'UNIT["Foot_US', 32026), | |||
('"NAD_1927_Oregon_North"', 'UNIT["Foot_US', 32026), | |||
('"NAD_1983_StatePlane_Washington_South_FIPS_4602_Feet"', 'UNIT["Foot_US', 2286), | |||
('"NAD_1983_HARN_StatePlane_Washington_South_FIPS_4602_Feet"', 'UNIT["Foot_US', 2927), | |||
('"NAD_1983_HARN_StatePlane_Washington_South_FIPS_4602"', 'UNIT["Foot_US', 2927), | |||
('"NAD_1983_StatePlane_Washington_North_FIPS_4601_Feet"', '', 2285), | |||
('NAD_1983_UTM_Zone_11N', 'UNIT["Meter"', 26911), | |||
('"WGS_1984_Web_Mercator_Auxiliary_Sphere"', 'UNIT["Meter"', 3857), | |||
('"WGS_1984_Web_Mercator"', 'UNIT["Meter"', 3857), | |||
('GEOGCS["GCS_North_American_1983_HARN",DATUM["D_North_American_1983_HARN",SPHEROID["GRS_1980",6378137.0,298.257222101]],PRIMEM["Greenwich",0.0]', 'UNIT["Degree",0.0174532925199433]]', 4152), | |||
('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]]', '', 4269), | |||
('GEOGCS["GCS_WGS_1984",DATUM["D_WGS_1984",SPHEROID["WGS_1984"', 'UNIT["Degree"', 4326), | |||
('GEOGCS["GCS_North_American_1983",DATUM["D_North_American_1983",SPHEROID["GRS_1980",6378137,298.257222101]],PRIMEM["Greenwich",0],UNIT["Degree",0.017453292519943295]]', '', 4326), | |||
] | |||
def find_epsg(prj): | |||
""" Point me at a PRJ string and I will attempt to tell you the EPSG code as an integer. """ | |||
try: | |||
for tup in lov: | |||
(desc, units, epsg) = tup | |||
if prj.find(desc) != -1 and (len(units)==0 or prj.find(units) != -1): | |||
return epsg | |||
except: | |||
print tup | |||
exit(0) | |||
return None # huh? | |||
######################################################################### | |||
if __name__ == '__main__': | |||
matches = [] | |||
for root, dirnames, filenames in os.walk('/g/GISData'): | |||
for filename in fnmatch.filter(filenames, '*.prj'): | |||
matches.append(os.path.join(root, filename)) | |||
for prjfile in matches: | |||
fp = open(prjfile, 'r') | |||
line = fp.readline() | |||
fp.close() | |||
epsg = find_epsg(line) | |||
if epsg == None: | |||
print prjfile, line | |||
exit(1) | |||
else: | |||
print prjfile, epsg | |||
exit(0) | |||
</pre> | |||
Latest revision as of 23:53, 19 July 2013
Let's say someone has a shapefile and they want to know the appropriate EPSG code for it. How do I do that??
So far I can see that most of the files have a descriptive string for example
"NAD_1983_StatePlane_California_II_FIPS_0402"
and they also have units
UNIT["Meter",1.0]
That's often all I need.
Here is the code that works for me...
#!/usr/bin/python import sys,os import fnmatch lov = [ ('"NAD_1983_StatePlane_California_II_FIPS_0402"', 'UNIT["Meter', 26942), ('"NAD_1983_StatePlane_California_II_FIPS_0402"', 'UNIT["Foot', 2226), ('NAD_1983_California_Albers', 'UNIT[Meter,', 3310), ('"NAD83 / California Albers"', '', 3310), ('"NAD_1927_StatePlane_California_I_FIPS_0401"', 'UNIT["Foot_US"', 26741), ('PROJCS["Custom",GEOGCS["GCS_North_American_1927",DATUM["D_North_American_1927",SPHEROID["Clarke_1866",6378206.4,294.9786982]],PRIMEM["Greenwich",0],UNIT["Degree",0.0174532925199433]],PROJECTION["Albers"],PARAMETER["False_Easting",0],PARAMETER["False_Northing",-4000000],PARAMETER["Central_Meridian",-120],PARAMETER["Standard_Parallel_1",34],PARAMETER["Standard_Parallel_2",40.5]', 'UNIT["Meter",1]', 3309), ('"Old_Hawaiian_StatePlane_Hawaii_3_FIPS_5103"', 'UNIT["Foot_US"', 3563), ('"NAD_1983_StatePlane_Nevada_West_FIPS_2703_Feet"', '', 102309), ('"NAD_1983_HARN_Lambert_Conformal_Conic"', 'UNIT["User_Defined_Unit",0.3048', 2913), ('"NAD_1983_HARN_StatePlane_Oregon_North_FIPS_3601_Feet_Intl"', '', 2913), ('"NAD_1983_HARN_StatePlane_Oregon_North_FIPS_3601"', 'UNIT["Foot"', 2913), ('"NAD_1983_HARN_StatePlane_Oregon_North_FIPS_3601_IntlFeet"', '', 2913), ('"NAD_1983_HARN_StatePlane_Oregon_South_FIPS_3602_Feet_Intl"', 'UNIT["Foot', 2914), ('"NAD_1983_HARN_StatePlane_Oregon_South_FIPS_3602_AltUnits"', 'UNIT["User_Defined_Unit",0.3047999902464]', 2914), # DOGAMI uses this and it's sort of a guess ('PROJCS["_MI_0",GEOGCS["NAD83",DATUM["D_NAD83",SPHEROID["Geodetic_Reference_System_of_1980",6378137,298.2572221009113]],PRIMEM["Greenwich",0],UNIT["Degree",0.017453292519943295]],PROJECTION["Lambert_Conformal_Conic"],PARAMETER["standard_parallel_1",43],PARAMETER["standard_parallel_2",45.5],PARAMETER["latitude_of_origin",41.75],PARAMETER["central_meridian",-120.5],PARAMETER["false_easting",1312335.958],PARAMETER["false_northing",0]', 'UNIT["Foot (International)",0.3048]]', 3644), ('"NAD_1983_Oregon_Statewide_Lambert"', 'UNIT["Meter', 2991), ('"NAD_1983_Oregon_Statewide_Lambert_Feet_Intl"', '', 2992), ('"OGICdefaultAG82"', 'UNIT["Foot"', 2992), ('"NAD_1983_HARN_Lambert_Conformal_Conic"', 'UNIT["Foot"', 2994), ('"NAD_1927_StatePlane_Oregon_North_FIPS_3601"', 'UNIT["Foot_US', 32026), ('"NAD_1927_Oregon_North"', 'UNIT["Foot_US', 32026), ('"NAD_1983_StatePlane_Washington_South_FIPS_4602_Feet"', 'UNIT["Foot_US', 2286), ('"NAD_1983_HARN_StatePlane_Washington_South_FIPS_4602_Feet"', 'UNIT["Foot_US', 2927), ('"NAD_1983_HARN_StatePlane_Washington_South_FIPS_4602"', 'UNIT["Foot_US', 2927), ('"NAD_1983_StatePlane_Washington_North_FIPS_4601_Feet"', '', 2285), ('NAD_1983_UTM_Zone_11N', 'UNIT["Meter"', 26911), ('"WGS_1984_Web_Mercator_Auxiliary_Sphere"', 'UNIT["Meter"', 3857), ('"WGS_1984_Web_Mercator"', 'UNIT["Meter"', 3857), ('GEOGCS["GCS_North_American_1983_HARN",DATUM["D_North_American_1983_HARN",SPHEROID["GRS_1980",6378137.0,298.257222101]],PRIMEM["Greenwich",0.0]', 'UNIT["Degree",0.0174532925199433]]', 4152), ('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]]', '', 4269), ('GEOGCS["GCS_WGS_1984",DATUM["D_WGS_1984",SPHEROID["WGS_1984"', 'UNIT["Degree"', 4326), ('GEOGCS["GCS_North_American_1983",DATUM["D_North_American_1983",SPHEROID["GRS_1980",6378137,298.257222101]],PRIMEM["Greenwich",0],UNIT["Degree",0.017453292519943295]]', '', 4326), ] def find_epsg(prj): """ Point me at a PRJ string and I will attempt to tell you the EPSG code as an integer. """ try: for tup in lov: (desc, units, epsg) = tup if prj.find(desc) != -1 and (len(units)==0 or prj.find(units) != -1): return epsg except: print tup exit(0) return None # huh? ######################################################################### if __name__ == '__main__': matches = [] for root, dirnames, filenames in os.walk('/g/GISData'): for filename in fnmatch.filter(filenames, '*.prj'): matches.append(os.path.join(root, filename)) for prjfile in matches: fp = open(prjfile, 'r') line = fp.readline() fp.close() epsg = find_epsg(line) if epsg == None: print prjfile, line exit(1) else: print prjfile, epsg exit(0)