source: cows/trunk/cows/service/imps/csmlbackend/wcs_csmllayer.py @ 4695

Subversion URL: http://proj.badc.rl.ac.uk/svn/ndg/cows/trunk/cows/service/imps/csmlbackend/wcs_csmllayer.py@4695
Revision 4695, 5.8 KB checked in by domlowe, 12 years ago (diff)

Merging WCS branch back into cows trunk

Line 
1"""
2implementation of ILayerMapper, IwmsLayer, IwmsDimension, ILayerSlab interfaces, as defined in wms_iface.py & wxs_iface.py
3
4"""
5import os, string
6import csml
7try:
8    import cdms2 as cdms
9except:
10    import cdms
11from pylons import config, request, session
12from copy import copy
13
14import logging
15log = logging.getLogger(__name__)
16
17from cows.model.wcs import WcsDatasetSummary
18from cows.service.imps.csmlbackend.csmlcommon import CSMLLayerMapper, CSMLConnector, extractToNetCDF
19
20
21class CSMLwcsCoverageMapper(CSMLLayerMapper): 
22    """
23    Map keyword arguments to a collection of layers (coverages)
24    Supports the retrieval of coverages according to arbitrary
25    keyword/value pairs. 
26    """
27    def __init__(self):
28        super(CSMLwcsCoverageMapper, self).__init__()
29       
30   
31    def getInfo(self, feature):
32        ''' given a csml feature, return info about the layer/feature
33        @return:   title, abstract, units, crss '''
34
35        title, abstract = super(CSMLwcsCoverageMapper, self).getInfo(feature)
36        units=feature.getDomainUnits()
37        tmpunits=copy(units)
38        tmpunits.reverse()
39        domain = feature.getDomain()
40        tax= feature.getTimeAxis()
41        timelimits=[domain[tax][0],domain[tax][len(domain[tax])-1]] 
42        crs=feature.getNativeCRS()
43        crss=[self._crscat.getCRS(crs).twoD]
44        if 'EPSG:4326' in crss:
45            crss.append('CRS:84')
46            crss.append('WGS84')
47
48        return title, abstract, timelimits, units, crss
49           
50    def getCoverageDescription(self):
51        pass
52       
53     
54    def map(self, **kwargs):
55        """
56        Given csml.parser.Dataset object list the names of
57        all layers available.
58       
59        @return: A mapping of layer names to ILayer implementations.
60        @raise ValueError: If no layers are available for these keywords.
61        """
62        fileoruri=kwargs['fileoruri']
63        if fileoruri in self.layermapcache.keys():
64            #we've accessed this layer map before, get it from the cache dictionary
65            return self.layermapcache[fileoruri]
66         
67        ds = self.connector.getCsmlDoc(fileoruri)
68        layermap={}
69        self._crscat=csml.csmllibs.csmlcrs.CRSCatalogue()
70        for feature in csml.csmllibs.csmlextra.listify(ds.featureCollection.featureMembers):
71            title, abstract, timelimits, units, crss=self.getInfo(feature)
72            layermap[feature.id]=CSMLCoverage([title],[abstract], timelimits, units, crss, feature)
73        if len(layermap) > 0:
74            self.layermapcache[fileoruri]=layermap
75            return layermap
76        else:
77            raise ValueError
78
79
80
81
82class CSMLCoverage(): #TODO: define ICoverage
83    """ represents a WCS Coverage. Implements ICoverage """
84   
85    def __init__(self, title, abstract, timelimits, units, crss, feature):
86        self.title=title
87        self.abstract=abstract
88        self.description='TO DO - coverage description'
89        self.timeLimits=timelimits
90        self.units=units
91        self.crss=crss
92        self._feature=feature
93        self.id=feature.id
94        bb = self._feature.getCSMLBoundingBox().getBox()
95        #convert 0 - 360 to -180, 180 as per common WXS convention
96        if abs(bb[2]-bb[0]) >= 359 and abs(bb[2]-bb[0]) < 361:
97            bb[0], bb[2]=-180, 180
98        self.wgs84BBox = bb
99        self.featureInfoFormats = ['text/html']
100   
101    def getBBox(self, crs):
102        """
103        @return: A 4-typle of the bounding box in the given coordinate
104            reference system.
105        """
106        #TODO: make this generic
107        return self.wgs84BBox
108   
109    def getCvg(self, bbox, time=None, crs=None, response_crs=None):
110        """
111        Creates a subset of the layer in a particular CRS and set of
112        dimensions.
113        #TODO: this is all out of synch
114        @param crs: The coordinate reference system.
115        @param dimValues: A mapping of dimension names to dimension values
116            as specified in the IDimension.extent
117        @param renderOpts: A generic mapping object for passing rendering
118            options
119        @return: An object implementing ILayerSlab
120        #create netcdf for whole lat/lon for given dimValues, use to init slab
121        """
122
123        log.debug('WCS: getSubset(%s, %s, %s)' % (bbox, time, crs))
124        #unpack the Boundingbox.
125       
126       
127        ############################# from old WCS stack #################
128        boundingbox=bbox
129        lon=self._feature.getLongitudeAxis()
130        lat=self._feature.getLatitudeAxis()
131        t=self._feature.getTimeAxis()
132        if None in [lon, lat, t]:
133            #TODO need to return a suitable wcs error.
134            print 'warning, could not get correct axis info'
135            #best guess!
136            if t is None:
137                t='time'
138            if lon is None:
139                lon = 'longitude'
140            if lat is None:
141                lat = 'latitude'
142       
143        #create selection dictionary:
144        sel={}
145        sel[lat]=(boundingbox[1], boundingbox[3])
146        sel[lon]=(boundingbox[0], boundingbox[2])
147        if time is not None:
148            if  type(time) is unicode:
149                sel[t]=str(time)
150            else:
151                sel[t]=time
152        #z is the 4th axis (eg height or pressure).
153        #NOTE, need to decide whether bounding box is of form: x1,y1,z1,x2,y2,z2 or x1,y1,x2,y2,z1,z2
154        #currently the latter is implemented.
155           
156        if len(boundingbox)  == 6:
157            for ax in self._feature.getAxisLabels():
158                if ax not in [lat, lon, t]:
159                    #must be Z 
160                    z=str(ax)
161                    sel[z]=(boundingbox[4], boundingbox[5])
162           
163           
164        axisNames=self._feature.getAxisLabels()
165        ##################################################################
166       
167        filename = extractToNetCDF(self._feature, sel)
168        return filename
169       
Note: See TracBrowser for help on using the repository browser.