source: TI05-delivery/ows_framework/branches/ows_framework-refactor/ows_common/ows_common/service/imps/wms_csmllayer.py @ 3760

"""
implementation of ILayerMapper, ILayer, IDimension, ILayerSlab interfaces, as defined in wms_iface.py

"""
import os, string
import csml
try:
    import cdms2 as cdms
except:
    import cdms
import Image
from copy import copy
from pywms.render_imp import RGBARenderer
from matplotlib import cm
import genutil
from pylons import config  #config must have tmpfilebuffer and csmlstore values
import cdtime


class CSMLLayerMapper(object):
    """
    Map keyword arguments to a collection of layers.
    Supports the retrieval of sets of layers according to arbitary
    keyword/value pairs.
    Implements  ILayerMapper 
    
    """
    
    def _getInfo(self, feature):
        ''' given a csml feature, return info about the layer/feature
        @return:   title, abstract, dimensions, units, crss '''

        try:
            title=feature.name.CONTENT
        except:
            title=''
        try:
            abstract=feature.description.CONTENT
        except:
            abstract=title
        
        units=feature.getDomainUnits()
        dimensions={}
        tmpunits=copy(units)
        tmpunits.reverse()
        domain = feature.getDomain()
        for dim in feature.getAxisLabels():
            nextdim=CSMLDimension(domain, dim, tmpunits.pop())
            if dim not in ['latitude', 'longitude']:
                dimensions[dim]=nextdim
        crs=feature.getNativeCRS()
        crss=[self._crscat.getCRS(crs).twoD]
        if 'EPSG:4326' in crss:
            crss.append('CRS:84')
            crss.append('WGS84')
        return title, abstract, dimensions, units, crss
    
    def map(self, **kwargs):
        """
        Given csml.parser.Dataset object list the names of
        all layers available.
        
        @return: A mapping of layer names to ILayer implementations.
        @raise ValueError: If no layers are available for these keywords. 
        """
        fileoruri=kwargs['fileoruri']
        
        #TODO - handle file paths/directories URIs in config
        #.xml or .csml extensions are supported:
        filename='%s/%s.csml'%(config['csmlstore'],fileoruri)
        if not os.path.exists(filename):
           filename='%s/%s.xml'%(config['csmlstore'],fileoruri)
        if not os.path.exists(filename):
            raise Exception(str('CSML File could not be found: %s')%filename)
            
        ds=csml.parser.Dataset(filename)
            
        layermap={}
        self._crscat=csml.csmllibs.csmlcrs.CRSCatalogue()
        for feature in csml.csmllibs.csmlextra.listify(ds.featureCollection.featureMembers):
            title, abstract, dimensions, units, crss=self._getInfo(feature)
            layermap[feature.id]=CSMLLayer(title,abstract, dimensions, units, crss, feature)
        if len(layermap) > 0:
            return layermap
        else:
            raise ValueError


class CSMLLayer(object):
    """
     representing a WMS layer.    Implements ILayer

    @ivar title: The layer title.  As seen in the Capabilities document.
    @ivar abstract:  Abstract as seen in the Capabilities document.
    @ivar dimensions: A dictionary of IDimension objects.
    @ivar units: A string describing the units.
    @ivar crss: A sequence of SRS/CRSs supported by this layer.

    @todo: Do we need minValue/maxValue?

    """

    def __init__(self, title, abstract, dimensions, units, crss, feature):
        self.featureInfoFormats=None #NotImplemented
        self.title=title
        self.abstract=abstract
        self.dimensions=dimensions
        self.units=units
        self.crss=crss
        self._feature=feature
        self.legendSize=(30,100)
        try:
            self.wgs84BBox = self.getBBox('EPSG:4326')
        except:
            raise ValueError("Layer must provide a bounding box in EPSG:4326 "
                             "coordinates for compatibility with WMS-1.3.0")

    def getBBox(self, crs):
        """
        @return: A 4-typle of the bounding box in the given coordinate
            reference system.
        """
        bb= self._feature.getCSMLBoundingBox().getBox()
        #convert 0 - 360 to -180, 180 as per common WMS convention
        if abs(bb[2]-bb[0]) >= 359 and abs(bb[2]-bb[0]) < 361:
            bb[0], bb[2]=-180, 180
        return bb
        #raise NotImplementedError
        
    def getSlab(self, crs, dimValues=None, renderOpts={}):
        """
        Creates a slab of the layer in a particular CRS and set of
        dimensions.

        @param crs: The coordinate reference system.
        @param dimValues: A mapping of dimension names to dimension values
            as specified in the IDimension.extent
        @param renderOpts: A generic mapping object for passing rendering
            options
        @return: An object implementing ILayerSlab
        #create netcdf for whole lat/lon for given dimValues, use to init slab
        """
        #unicode conversion
        for dimval in dimValues:
            if dimval != 'time':
                dimValues[dimval]=float(dimValues[dimval])
            else:
                #remove any trailing Zs from time string
                if dimValues[dimval] [-1:] in ['Z', 'z']:
                    dimValues[dimval]=dimValues[dimval][:-1]
        if type(self._feature) == csml.parser.GridSeriesFeature:
            randomname= csml.csmllibs.csmlextra.getRandomID() + '.nc'
            result= self._feature.subsetToGridSeries(config['tmpfilebuffer'], ncname=randomname, **dimValues)
            #for now have to read netcdf back from disk (limitiation of CSML api)
            netcdf=cdms.open(result[1])
            #and then delete the temporary file
            os.system('rm %s'%result[1])
            bbox=self.getBBox(crs)
            return CSMLLayerSlab(netcdf, self, crs, dimValues, renderOpts, bbox)
        else:
            raise NotImplementedError
        
    def getCacheKey(self, crs, dimValues=None, renderOpts={}):
        """
        Create a unique key for use in caching a slab.

        The intention here is that most of the work should be done when
        instantiating an ILayerSlab object.  These can be cached by the
        server for future use.  The server will first call getCacheKey()
        for the slab creation arguments and if the key is in it's cache
        it will use a pre-generated ILayerSlab object.

        """
        return None
        #raise NotImplementedError



class CSMLDimension(object):
    """
    implements IDimension
    @ivar units: The units string.
    @ivar extent: Sequence of extent values.

    """
    
    def __init__(self, domain, dimname, unit):
        self.units = unit
        self.extent = []
        #patch to handle current limitations of multiple time dimension scanning in csml. 
        if string.lower(self.units)[:10] in ['days_since', 'seconds_si', 'minutes_si', 'hours_sinc','months _sin', 'years_sinc']:
            if type(domain[dimname][0]) is not str   :
                tunits=self.units.replace('_', ' ')
                for val in domain[dimname]:
                    csmltime= csml.csmllibs.csmltime.UDtimeToCSMLtime(cdtime.reltime(float(val), tunits).tocomp())
                    self.extent.append(csmltime)
            else:
                for val in domain[dimname]:
                    self.extent.append(str(val))
        else:
            for val in domain[dimname]:
                self.extent.append(str(val))
        #for time axis replace units with iso string
        if dimname == 'time':
            self.units='ISO8601'
       
class CSMLLayerSlab(object):
    """
    Implements LayerSlab
    Represents a particular horizontal slice of a WMS layer.

    ILayerSlab objects are designed to be convenient to cache.
    They should be pickleable to enable memcached support in the future.

    @ivar layer: The source ILayer instance.
    @ivar crs: The coordinate reference system.
    @ivar dimValues: A mapping of dimension values of this view.
    @ivar renderOpts: The renderOpts used to create this view.
    @ivar bbox: The bounding box as a 4-tuple.
    """
    
    def __init__(self, netcdf, layer, crs, dimValues, renderOpts, bbox):
        self._netcdf=netcdf
        self.layer = layer
        self.crs = crs
        self.dimValues = dimValues
        self.renderOpts=renderOpts
        self.bbox=bbox
        
        
    def getImage(self, bbox, width, height):
        """
        Create an image of a sub-bbox of a given size.

        @ivar bbox: A bbox 4-tuple.
        @ivar width: width in pixels.`  
        @ivar height: height in pixels.
        @return: A PIL Image object.

        """
        cmap=eval(config['colourmap']) # renderOpts is hook for colourmap, for now use config
        grid=Grid(self.layer, self._netcdf, bbox, width, height)
        #how to handle varmin,varmax? ..read array?
        #minval, maxval=genutil.minmax(grid.value)
        minval=min(min(l) for l in grid.value)
        maxval=max(max(l) for l in grid.value)
        renderer=RGBARenderer(minval, maxval)         
        return renderer.renderGrid(grid, bbox, width, height, cmap)
    
class Grid(object):
    """A class encapsulating a simple regularly spaced, rectilinear
    grid.  This is the only type of grid pywms is expected to
    understand and adaptors should be provided to connect to
    underlying implementations such as cdms or csml.

    @cvar crs: Coordinate reference system

    @ivar x0: coordinate of the lower left corner.
    @ivar y0: coordinate of the lower left corner.
    @ivar dx: The x grid spacing.
    @ivar dy: The y grid spacing.
    @ivar nx: The number of x grid points.
    @ivar ny: The number of y grid points.
    @ivar value: A masked array of the grid values.
    @ivar ix: The dimension index of longidude in value
    @ivar iy: The dimension index of latitude in value
    @ivar long_name: The name of the field.
    @ivar units: The units of the field.
    """
    def __init__(self, layer, netcdf, bbox, width, height):
        #we know the axes are called latitude and longitude as the CSML code has written it:

        v=netcdf(layer.title)
        tvar=v(latitude=(bbox[1], bbox[3]), longitude=(bbox[0],bbox[2]),squeeze=1)
        order=tvar.getOrder()
        #array of data
        self.value=tvar.getValue()
        #order of axes
        if order == 'xy':
            self.ix=0
            self.iy=1
        else:
            self.ix=1
            self.iy=0
        lat = tvar.getLatitude()
        lon = tvar.getLongitude()
        self.x0=lon[0]
        self.y0=lat[0]
        self.dx=abs(lon[0]-lon[1])
        self.dy=abs(lat[0]-lat[1])
        self.nx=len(lon)
        self.ny=len(lat)
        self.long_name=tvar.id  #TODO, get long name from feature
        self.units=tvar.units