crs fixes

OWS.py

@@ -11,8 +11,9 @@ import ConfigParser
 import md5
 import mapscript
 from string import Template
+from osgeo import osr
+from osgeo import ogr
 
-from owslib.wfs import WebFeatureService
 
 class OWS:
 
@@ -57,7 +58,12 @@ class OWS:
 
     def __getCapabilities(self):
         self.parsedUrl = urlparse.urlparse(self.url)
-        self.capabilities = WebFeatureService(url=self.url)
+        if self.service == "WFS":
+            from owslib.wfs import WebFeatureService
+            self.capabilities = WebFeatureService(url=self.url)
+        elif self.service == "WCS":
+            from owslib.wcs import WebCoverageService
+            self.capabilities = WebCoverageService(url=self.url,version="1.0.0")
 
     def getParams(self):
         return urlparse.parse_qs(self.qstring)
@@ -210,6 +216,39 @@ class OWS:
 
             return defFileName
 
+    def getLayerExtent(self,layer,crs=None,wkt=None):
+        """Get extent of layer in form of minx, miny, maxx,maxy
+        """
+
+        bbox = None
+        if layer.boundingBoxWGS84:
+
+            dest = osr.SpatialReference()
+            if crs:
+                dest.ImportFromEPSG(int(crs.split(":")[1]))
+            else:
+                dest.ImportFromWkt(wkt)
+            source = osr.SpatialReference()
+            source.ImportFromEPSG(4326)
+
+            bbox = []
+
+            # TODO rewrite this using ogr CoordinateTransformation
+            # http://www.gdal.org/ogr/osr_tutorial.html
+            geom = ogr.CreateGeometryFromWkt("""POINT(%s %s)""" % (layer.boundingBoxWGS84[0],layer.boundingBoxWGS84[1]),source)
+            geom.TransformTo(dest)
+            bbox.append(geom.GetX())
+            bbox.append(geom.GetY())
+
+            geom = ogr.CreateGeometryFromWkt("""POINT(%s %s)""" % (layer.boundingBoxWGS84[2],layer.boundingBoxWGS84[3]),source)
+            geom.TransformTo(dest)
+            bbox.append(geom.GetX())
+            bbox.append(geom.GetY())
+
+            logging.debug("Setting extent for layer <%s> to %s" %\
+                    (layer.id,bbox))
+        return bbox
+
 def getService():
 
     qstring = os.environ["QUERY_STRING"]

wcs/__init__.py

@@ -4,9 +4,8 @@
 from OWS import OWS
 import mapscript
 from osgeo import gdal
-from osgeo import osr
 from osgeo import ogr
-from string import Template
+from osgeo import osr
 import os
 import logging
 import re
@@ -17,7 +16,6 @@ class WCS(OWS):
     request = '' # GetCapabilies, ...
     parser = None
     version = None
-    wcsns = None
 
     def __init__(self,url=None,qstring=None,configFiles=None):
         OWS.__init__(self,url,qstring)
@@ -27,8 +25,7 @@ class WCS(OWS):
 
         lyrobj = None
 
-        self.version = self.capabilities.attrib["version"]
-        self.wcsns = self.capabilities.nsmap[self.capabilities.prefix]
+        self.version = self.capabilities.version
 
         mapobj = self.getMapObj(mapfilename)
 
@@ -36,9 +33,8 @@ class WCS(OWS):
         # for each layer from the WCS Capabilities file, transform it to
         # mapscript.layerObj
         #
-        logging.info(self.capabilities["Contents"].getchildren())
-        for layer in self.capabilities["{%s}Contents" % (self.wcsns)].getchildren():
-            name =  layer.Identifier.text
+        for name in self.capabilities.contents:
+            layer = self.capabilities.contents[name]
             logging.debug("Creating layer %s" % name)
 
             # create layer definition file
@@ -49,18 +45,31 @@ class WCS(OWS):
             # create the layer
             ds = gdal.Open(layerDefFile)
             try:
-
                 lyrobj = mapscript.layerObj(mapobj)
                 lyrobj.name = name
                 lyrobj.data = layerDefFile
-                lyrobj.title = layer["{%s}%s" % (self.owsns11,"Title")].text
-                lyrobj.setMetaData("wms_title",layer["{%s}%s" % (self.owsns11,"Title")].text)
-                if "{%s}%s" % (self.owsns11,"Abstract") in layer:
-                    lyrobj.setMetaData("ows_abstract",  layer["{%s}%s" % (self.owsns11,"Abstract")].text)
+                if layer.title:
+                    lyrobj.title = layer.title
+                    lyrobj.setMetaData("wms_title",layer.title)
+                #if "{%s}%s" % (self.owsns11,"Abstract") in layer:
+                #    lyrobj.setMetaData("ows_abstract",  layer["{%s}%s" % (self.owsns11,"Abstract")].text)
                 lyrobj.setMetaData("wms_srs",self.config.get("MapServer","srs"))
-                sr = osr.SpatialReference()
-                sr.ImportFromWkt(ds.GetProjection())
-                lyrobj.setProjection(sr.ExportToProj4())
+                extent = None
+                if layer.crsOptions:
+                    lyrobj.setProjection(layer.crsOptions[0])
+                    extent = self.getLayerExtent(layer,layer.crsOptions[0])
+                else:
+                    sr = osr.SpatialReference()
+                    sr.ImportFromWkt(ds.GetProjection())
+                    if sr.AutoIdentifyEPSG() == 0:
+                        epsg = "{code}:{value}".format(code=sr.GetAuthorityName("PROJCS"),value=sr.GetAuthorityCode("PROJCS"))
+                        extent = self.getLayerExtent(layer,epsg)
+                        lyrobj.setProjection(epsg)
+                    else:
+                        extent = self.getLayerExtent(layer,wkt=ds.GetProjection())
+                        lyrobj.setProjection(sr.ExportToProj4())
+                lyrobj.setMetaData("wms_extent","%s %s %s %s" % \
+                        (extent[0],extent[1],extent[2],extent[3]))
                 lyrobj.type = mapscript.MS_LAYER_RASTER
                 lyrobj.dump = mapscript.MS_TRUE
                 lyrobj.template = "foo"
@@ -70,7 +79,7 @@ class WCS(OWS):
             except Exception,e:
                 mapobj.removeLayer(lyrobj.index)
                 logging.warning("Layer %s: %s, skipping." %\
-                                (layer.Identifier.text,e.message))
+                                (name,e.message))
 
         self.saveMapfile(mapobj,mapfilename)
         return mapobj

wfs/__init__.py

@@ -9,7 +9,6 @@ import urllib
 import urlparse
 import logging
 from osgeo import ogr
-from osgeo import osr
 import os
 
 class WFS(OWS):
@@ -55,10 +54,10 @@ class WFS(OWS):
 
             if ds:
                 ogrLayer = ds.GetLayerByName(name)
-                # TODO : udelat to pomoci owslib
                 extent = self.getLayerExtent(layer,layer.crsOptions[0])
-                lyrobj.setMetaData("wms_extent","%s %s %s %s" % \
-                        (extent[0],extent[1],extent[2],extent[3]))
+                if extent:
+                    lyrobj.setMetaData("wms_extent","%s %s %s %s" % \
+                            (extent[0],extent[1],extent[2],extent[3]))
                 lyrobj.type = self._getLayerType(ogrLayer)
             else:
                 mapobj.removeLayer(mapobj.numlayers-1)
@@ -96,31 +95,6 @@ class WFS(OWS):
         if self.capabilities.identification.abstract:
             mapobj.setMetaData("wms_abstract",self.capabilities.identification.abstract)
 
-    def getLayerExtent(self,layer,crs):
-        """Get extent of layer in form of minx, miny, maxx,maxy
-        """
-
-        dest = osr.SpatialReference()
-        dest.ImportFromEPSG(int(crs.split(":")[1]))
-        source = osr.SpatialReference()
-        source.ImportFromEPSG(4326)
-
-        bbox = []
-
-        geom = ogr.CreateGeometryFromWkt("""POINT(%s %s)""" % (layer.boundingBoxWGS84[0],layer.boundingBoxWGS84[1]),source)
-        geom.TransformTo(dest)
-        bbox.append(geom.GetX())
-        bbox.append(geom.GetY())
-
-        geom = ogr.CreateGeometryFromWkt("""POINT(%s %s)""" % (layer.boundingBoxWGS84[2],layer.boundingBoxWGS84[3]),source)
-        geom.TransformTo(dest)
-        bbox.append(geom.GetX())
-        bbox.append(geom.GetY())
-
-        logging.debug("Setting extent for layer <%s> to %s" %\
-                (layer.id,bbox))
-        return bbox
-
     def _getLayerType(self,layer):
         """Returns MS layer type based on ogr.Layer.GetGeomType