First ideas and structures

This commit is contained in:
Karim Bahgat
2015-01-11 18:57:56 +01:00
parent 71d2c90b79
commit a8637812b7
4 changed files with 517 additions and 0 deletions
+123
View File
@@ -0,0 +1,123 @@
import json
#################
# USER FUNCTIONS
#################
# parse from text strings
# possible use module: https://github.com/rockdoc/grabbag/wiki/CRS-WKT-Parser
# also note some paramter descriptions: http://www.geoapi.org/3.0/javadoc/org/opengis/referencing/doc-files/WKT.html
# and see gdal source code: http://gis.stackexchange.com/questions/129764/how-are-esri-wkt-projections-different-from-ogc-wkt-projections
def from_epsg_code(string):
# must go online (or look up local table) to get crs details
pass
def from_esri_code(string):
# must go online (or look up local table) to get crs details
pass
def from_sr_code(string):
# must go online (or look up local table) to get crs details
pass
def from_esri_wkt(string):
# parse arguments into components
# use args to create crs
pass
def from_ogc_wkt(string):
# parse arguments into components
# use args to create crs
pass
def from_unknown_wkt(string):
# detect if ogc wkt or esri wkt
# TIPS: esri wkt datums all use "D_" before the datum name
# then load with appropriate function
pass
def from_proj4(string):
# parse arguments into components
# use args to create crs
pass
def from_ogc_urn(string):
# hmmm, seems like ogc urn could be anything incl online link, epsg, etc...
# if necessary, must go online (or lookup local table) to get details
# maybe test which of these and run their function?
pass
def from_unknown_text(text):
# detect type and load with appropriate function
if string.startswith("urn:"):
from_ogc_urn(string)
elif string.startswith("+proj="):
from_proj4(string)
elif string.startswith("PROJCS["):
from_unknown_wkt(string)
elif string.startswith("EPSG:"):
from_epsg_code(string)
elif string.startswith("ESRI:"):
from_esri_code(string)
elif string.startswith("SR-ORG:"):
from_sr_code(string)
else: raise Exception("Could not detect which type of crs")
def from_geotiff_params(**params):
pass
# convenience methods for loading from different sources
def from_url(url, format=None):
# first get string from url
# ...
# then load
if format:
# load string using specified format
pass
else:
from_unknown_text(string)
def from_file(filepath):
if filepath.endswith(".prj"):
string = open(filepath, "r").read()
from_esri_wkt(string)
elif filepath.endswith((".geojson",".json")):
crsinfo = json.load(filepath)["crs"]
if crsinfo["type"] == "name":
string = crsinfo["properties"]["name"]
from_unknown_text(string)
elif crsinfo["type"] == "link":
url = crsinfo["properties"]["name"]
type = crsinfo["properties"].get("type")
from_url(url, format=type)
else: raise Exception("invalid geojson crs type: must be either name or link")
elif filepath.endswith((".tif",".tiff",".geotiff")):
pass
# ...
+12
View File
@@ -0,0 +1,12 @@
# CRS Names
class proj_NAD87:
proj4 = ""
ogc_wkt = ""
class datum_WGS84:
proj4 = ""
ogc_wkt = ""
+77
View File
@@ -0,0 +1,77 @@
import urllib2
import re
def build_crs_table(savepath):
# create table
outfile = open(savepath, "wb")
# create fields
fields = ["codetype", "code", "proj4", "ogcwkt", "esriwkt"]
outfile.write("\t".join(fields) + "\n")
# make table from url requests
for codetype in ("epsg", "esri", "sr-org"):
print(codetype)
# collect existing proj list
print("fetching list of available codes")
codelist = []
page = 1
while True:
try:
link = 'http://spatialreference.org/ref/%s/?page=%s' %(codetype,page)
html = urllib2.urlopen(link).read()
codes = [match.groups()[0] for match in re.finditer(r'/ref/'+codetype+'/(\d+)', html) ]
if not codes: break
print("page",page)
codelist.extend(codes)
page += 1
except:
break
print("fetching string formats for each projection")
for i,code in enumerate(codelist):
# check if code exists
link = 'http://spatialreference.org/ref/%s/%s/' %(codetype,code)
urllib2.urlopen(link)
# collect each projection format in a table row
row = [codetype, code]
for resulttype in ("proj4", "ogcwkt", "esriwkt"):
try:
link = 'http://spatialreference.org/ref/%s/%s/%s/' %(codetype,code,resulttype)
result = urllib2.urlopen(link).read()
row.append(result)
except:
pass
print("projection %i of %i added" %(i,len(codelist)) )
outfile.write("\t".join(row) + "\n")
# close the file
outfile.close()
def crscode_to_string(codetype, code, format):
link = 'http://spatialreference.org/ref/%s/%s/%s/' %(codetype,code,format)
result = urllib2.urlopen(link).read()
return result
def crsstring_to_string(string, newformat):
# search string, if string is correct there should only be one correct match
link = 'http://spatialreference.org/ref/?search=%s' %string
searchresults = urllib2.urlopen(link).read()
# pick the first result
# ...regex...
# go to its url, with extension for the newformat
link = 'http://spatialreference.org/ref/%s/%s/%s/' %(codetype,code,newformat)
result = urllib2.urlopen(link).read()
return result
if __name__ == "__main__":
build_crs_table("crstable.txt")
+305
View File
@@ -0,0 +1,305 @@
#################
# CRS CLASSES
#################
# first classes for each crs element, from the proj4 common paramter listing
# https://trac.osgeo.org/proj/wiki/GenParms
# note also that in wkt the names of "PROJECTION", "DATUM", and "SPHEROID"
# matter and seem to be interpreted and computed upon, ie they carry meaning
# that is commonly understood by programs, so is not explicit in the crs specification.
# the paramters below simply modify certain aspects of the proj/datum/spheroids
# note that the names in wkt of "PROJCS" and "GEOGCS" seem to be purely
# for identifying and branding and can be changed at will.
# they do however act as shortcuts, so that proj4 can use +init=... to load
# everything automatically
# some of these names imply certain combinations of datum and spheroid and paramters.
# so in proj4 one simply needs to give that name, but in wkt one needs to spell it all out.
# +unit and +to_metre are what makes up 'UNIT["Meter",1.0]'
##+a Semimajor radius of the ellipsoid axis
class SemiMajorRadius:
proj4 = "+a"
def __init__(self, value):
pass
##+alpha ? Used with Oblique Mercator and possibly a few others
class Azimuth:
proj4 = "+alpha"
esri_wkt = "azimuth"
ogc_wkt = "azimuth"
geotiff = "AzimuthAngle"
def __init__(self, value):
pass
##+axis Axis orientation (new in 4.8.0)
class AxisOrientation:
proj4 = "+axis"
esri_wkt = None
ogc_wkt = "AXIS"
geotiff = None
def __init__(self, value):
pass
##+b Semiminor radius of the ellipsoid axis
class SemiMinorRadius:
proj4 = "+b"
def __init__(self, value):
pass
##+datum Datum name (see `proj -ld`)
class DatumName:
def __init__(self, value):
self.value = value
def to_proj4(self):
return "+datum=%s" %self.value
def to_ogc_wkt(self):
return 'DATUM[%s]' %self.value
def to_esri_wkt(self):
return self.to_ogc_wkt()
def to_geotiff(self):
pass
#return "GeogGeodeticDatum"
##+ellps Ellipsoid name (see `proj -le`)
class EllipsoidName:
def __init__(self, value):
self.value = value
def to_proj4(self):
return "+ellps=%s" %self.value
def to_ogc_wkt(self):
return 'SPHEROID[%s]' %self.value
def to_esri_wkt(self):
return self.to_ogc_wkt()
def to_geotiff(self):
pass
#return "GeogEllipsoid"
##+k Scaling factor (old name)
##+k_0 Scaling factor (new name)
class ScalingFactor:
def __init__(self, value):
self.value = value
def to_proj4(self):
return "+k_0=%s" %self.value
def to_ogc_wkt(self):
return 'PARAMETER["scale_factor", %s]' %self.value
def to_esri_wkt(self):
raise Exception("Paramater not supported by ESRI WKT")
def to_geotiff(self):
pass
#return "ScaleAtNatOrigin" # or ScaleAtCenter?
##+lat_0 Latitude of origin
class LatitudeOrigin:
def __init__(self, value):
self.value = value
def to_proj4(self):
return "+lat_0=%s" %self.value
def to_ogc_wkt(self):
return 'PARAMETER["latitude_of_origin", %s]' %self.value
def to_esri_wkt(self):
raise Exception("Paramater not supported by ESRI WKT")
def to_geotiff(self):
pass
#return "ProjCenterLat"
##+lat_1 Latitude of first standard parallel
class LatitudeFirstStndParallel:
proj4 = "+lat_1"
def __init__(self, value):
pass
##+lat_2 Latitude of second standard parallel
class LatitudeSecondStndParallel:
proj4 = "+lat_2"
def __init__(self, value):
pass
##+lat_ts Latitude of true scale
class LatitudeTrueScale:
def __init__(self, value):
self.value = value
def to_proj4(self):
return "+lat_ts=%s" %self.value
def to_ogc_wkt(self):
return 'PARAMETER["Standard_Parallel_1", %s]' %self.value
def to_esri_wkt(self):
return self.to_ogc_wkt()
def to_geotiff(self):
pass
#return "ProjStdParallel1"
##+lon_0 Central meridian
class CentralMeridian:
def __init__(self, value):
self.value = value
def to_proj4(self):
return "+lon_0=%s" %self.value
def to_ogc_wkt(self):
return 'PARAMETER["Central_Meridian", %s]' %self.value
def to_esri_wkt(self):
return self.to_ogc_wkt()
def to_geotiff(self):
pass
#return "ProjCenterLong"
##+lonc ? Longitude used with Oblique Mercator and possibly a few others
class LongitudeSpecial:
proj4 = "+lonc"
def __init__(self, value):
pass
##+lon_wrap Center longitude to use for wrapping (see below)
##+over Allow longitude output outside -180 to 180 range, disables wrapping (see below)
##+pm Alternate prime meridian (typically a city name, see below)
class PrimeMeridian:
def __init__(self, value):
self.value = value
def to_proj4(self):
return "+pm=%s" %self.value
def to_ogc_wkt(self):
return 'PRIMEM["Greenwich", %s]' %self.value
def to_esri_wkt(self):
return self.to_ogc_wkt()
##+proj Projection name (see `proj -l`)
class ProjectionName:
def __init__(self, value):
self.value = value
def to_proj4(self):
return "+proj=%s" %self.value
def to_ogc_wkt(self):
return 'PROJECTION["%s"]' %self.value
def to_esri_wkt(self):
return self.to_ogc_wkt()
##+south Denotes southern hemisphere UTM zone
##+towgs84 3 or 7 term datum transform parameters (see below)
class DatumShift:
def __init__(self, value):
self.value = value
def to_proj4(self):
return "+towgs84=%s" %"".join((str(val) for val in self.value))
def to_ogc_wkt(self):
return "TOWGS84[%s]" %"".join((str(val) for val in self.value))
def to_esri_wkt(self):
raise Exception("Paramater not supported by ESRI WKT")
##+to_meter Multiplier to convert map units to 1.0m
class MeterMultiplier:
def __init__(self, value):
self.value = value
def to_proj4(self):
return "+to_meter=%s" %self.value
def to_ogc_wkt(self):
# the stuff that comes after UNITS["meter", ... # must be combined with unittype in a unit class to make wkt
return str(self.value)
def to_esri_wkt(self):
return self.to_ogc_wkt()
##+units meters, US survey feet, etc.
class UnitType:
def __init__(self, value):
self.value = value
def to_proj4(self):
return "+units=%s" %self.value
def to_ogc_wkt(self):
# the stuff that comes after UNITS[... # must be combined with metermultiplier in a unit class to make wkt
return str(self.value)
def to_esri_wkt(self):
return self.to_ogc_wkt()
# special...
class Unit:
def __init__(self, unittype, metermultiplier):
self.unittype = unittype
self.metermultiplier = metermultiplier
def to_proj4(self):
return "+%s +%s" %(self.unittype, self.metermultiplier)
def to_ogc_wkt(self):
return 'UNIT["%s", %s]' %(self.unittype.to_ogc_wkt(), self.metermultiplier.to_ogc_wkt())
def to_esri_wkt(self):
return self.to_ogc_wkt()
##+x_0 False easting
class FalseEasting:
proj4 = "+x_0"
esri_wkt = "False_Easting"
ogc_wkt = "false_easting"
geotiff = "FalseEasting"
def __init__(self, value):
pass
##+y_0 False northing
class FalseEasting:
proj4 = "+y_0"
esri_wkt = "False_Northing"
ogc_wkt = "false_northing"
geotiff = "FalseNorthing"
def __init__(self, value):
pass
# then the final CRS object which is instantiated with all of these?
# remember to use +no_defs when outputting to proj4
# ...