From a8637812b79322524416ef15d5c300f1e9bb5304 Mon Sep 17 00:00:00 2001
From: Karim Bahgat <Karim.bahgat.norway@gmail.com>
Date: Sun, 11 Jan 2015 18:57:56 +0100
Subject: [PATCH] First ideas and structures

---
 io.py     | 123 ++++++++++++++++++++++
 names.py  |  12 +++
 online.py |  77 ++++++++++++++
 params.py | 305 ++++++++++++++++++++++++++++++++++++++++++++++++++++++
 4 files changed, 517 insertions(+)
 create mode 100644 io.py
 create mode 100644 names.py
 create mode 100644 online.py
 create mode 100644 params.py

diff --git a/io.py b/io.py
new file mode 100644
index 0000000..892d3d0
--- /dev/null
+++ b/io.py
@@ -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
+        # ...
+
+
+
+
+                                      
+        
diff --git a/names.py b/names.py
new file mode 100644
index 0000000..bb430c7
--- /dev/null
+++ b/names.py
@@ -0,0 +1,12 @@
+
+# CRS Names
+
+class proj_NAD87:
+    proj4 = ""
+    ogc_wkt = ""
+
+class datum_WGS84:
+    proj4 = ""
+    ogc_wkt = ""
+    
+
diff --git a/online.py b/online.py
new file mode 100644
index 0000000..43f914f
--- /dev/null
+++ b/online.py
@@ -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")
+
+
diff --git a/params.py b/params.py
new file mode 100644
index 0000000..b0136ca
--- /dev/null
+++ b/params.py
@@ -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
+# ...
+
+
+
+