- Timestamp:
- 02/10/08 12:08:51 (13 years ago)
- Location:
- cows/trunk/cows/service/imps
- Files:
-
- 2 edited
- 1 moved
Legend:
- Unmodified
- Added
- Removed
-
cows/trunk/cows/service/imps/csmlbackend/wms_csmllayer.py
r4269 r4270 1 """ 2 implementation of ILayerMapper, ILayer, IDimension, ILayerSlab interfaces, as defined in wms_iface.py 3 4 """ 5 import os, string 6 import csml 7 try: 8 import cdms2 as cdms 9 except: 10 import cdms 11 import Image 12 from copy import copy 13 from pywms.render_imp import RGBARenderer 14 from matplotlib import cm 15 import genutil 16 from pylons import config, request, session #config must have tmpfilebuffer and csmlstore values 17 import cdtime 18 import logging 19 log = logging.getLogger(__name__) 20 import ConfigParser 21 22 from cows.service.wxs_iface import ILayerMapper 23 from cows.service.wms_iface import IwmsLayer, IwmsDimension, IwmsLayerSlab 24 25 try: 26 from ndgUtils import ndgObject, ndgRetrieve 27 except: 28 log.warning("ndgUtils library could not be loaded, files in the eXist database won't be available, although you should still be able to access files from the csmlstore directory referenced in the ini file.") 29 30 31 32 33 class CSMLLayerMapper(ILayerMapper): 34 """ 35 Map keyword arguments to a collection of layers. 36 Supports the retrieval of sets of layers according to arbitary 37 keyword/value pairs. 38 Implements ILayerMapper 39 40 """ 41 def __init__(self): 42 self.layermapcache={} 43 44 def _getInfo(self, feature): 45 ''' given a csml feature, return info about the layer/feature 46 @return: title, abstract, dimensions, units, crss ''' 47 48 try: 49 title=feature.name.CONTENT 50 except: 51 title='' 52 try: 53 abstract=feature.description.CONTENT 54 except: 55 abstract=title 56 57 units=feature.getDomainUnits() 58 dimensions={} 59 tmpunits=copy(units) 60 tmpunits.reverse() 61 domain = feature.getDomain() 62 for dim in feature.getAxisLabels(): 63 nextdim=CSMLDimension(domain, dim, tmpunits.pop()) 64 if dim not in ['latitude', 'longitude']: 65 dimensions[dim]=nextdim 66 crs=feature.getNativeCRS() 67 crss=[self._crscat.getCRS(crs).twoD] 68 if 'EPSG:4326' in crss: 69 crss.append('CRS:84') 70 crss.append('WGS84') 71 return title, abstract, dimensions, units, crss 72 73 74 75 def _get_csml_doc(self, fileoruri): 76 """ 77 Gets a csml document from file or exist when passed a file name or uri 78 79 Note, access control is not implemented on file access, only on exist access. 80 81 """ 82 if string.find(fileoruri,'__NDG-A0__') == -1: 83 #it's a local file not an identifier, read from local csmlstore 84 file=fileoruri 85 csml_dir = config['csmlstore'] 86 path = os.path.join(csml_dir, file) 87 if os.path.exists(path+'.csml'): 88 f = path+'.csml' 89 elif os.path.exists(path+'.xml'): 90 f = path +'.xml' 91 else: 92 raise ValueError("Cannot find CSML file %s" % file) 93 94 d = csml.parser.Dataset() 95 d.parse(f) 96 97 else: 98 #it's an NDG identifier, get the document from exist. 99 uri=fileoruri 100 uriN=ndgObject(uri) 101 conflocation=config['ndgconfig'] 102 cf=ConfigParser.ConfigParser() 103 cf.read(conflocation) 104 status,x=ndgRetrieve(uriN, cf, discovery=0) 105 d=csml.parser.Dataset() 106 d.parse(x) 107 return d 108 109 def map(self, **kwargs): 110 """ 111 Given csml.parser.Dataset object list the names of 112 all layers available. 113 114 @return: A mapping of layer names to ILayer implementations. 115 @raise ValueError: If no layers are available for these keywords. 116 """ 117 fileoruri=kwargs['fileoruri'] 118 if fileoruri in self.layermapcache.keys(): 119 #we've accessed this layer map before, get it from the cache dictionary 120 return self.layermapcache[fileoruri] 121 122 ds = self._get_csml_doc(fileoruri) 123 layermap={} 124 self._crscat=csml.csmllibs.csmlcrs.CRSCatalogue() 125 for feature in csml.csmllibs.csmlextra.listify(ds.featureCollection.featureMembers): 126 title, abstract, dimensions, units, crss=self._getInfo(feature) 127 layermap[feature.id]=CSMLLayer(title,abstract, dimensions, units, crss, feature) 128 if len(layermap) > 0: 129 self.layermapcache[fileoruri]=layermap 130 return layermap 131 else: 132 raise ValueError 133 134 135 class CSMLLayer(IwmsLayer): 136 """ 137 representing a WMS layer. Implements ILayer 138 139 @ivar title: The layer title. As seen in the Capabilities document. 140 @ivar abstract: Abstract as seen in the Capabilities document. 141 @ivar dimensions: A dictionary of IDimension objects. 142 @ivar units: A string describing the units. 143 @ivar crss: A sequence of SRS/CRSs supported by this layer. 144 145 @todo: Do we need minValue/maxValue? 146 147 """ 148 149 def __init__(self, title, abstract, dimensions, units, crss, feature): 150 self.featureInfoFormats=None #NotImplemented 151 self.title=title 152 self.abstract=abstract 153 self.dimensions=dimensions 154 self.units=units 155 self.crss=crss 156 self._feature=feature 157 self.legendSize=(30,100) 158 bb= self._feature.getCSMLBoundingBox().getBox() 159 #convert 0 - 360 to -180, 180 as per common WMS convention 160 if abs(bb[2]-bb[0]) >= 359 and abs(bb[2]-bb[0]) < 361: 161 bb[0], bb[2]=-180, 180 162 self.wgs84BBox = bb 163 self.featureInfoFormats = ['text/html'] 164 try: 165 self.wgs84BBox = self.getBBox('EPSG:4326') 166 except: 167 raise ValueError("Layer must provide a bounding box in EPSG:4326 " 168 "coordinates for compatibility with WMS-1.3.0") 169 self.featureinfofilecache={} #used for caching netcdf file in getFeatureInfo 170 171 def getBBox(self, crs): 172 """ 173 @return: A 4-typle of the bounding box in the given coordinate 174 reference system. 175 """ 176 #bb= self._feature.getCSMLBoundingBox().getBox() 177 #convert 0 - 360 to -180, 180 as per common WMS convention 178 #if abs(bb[2]-bb[0]) >= 359 and abs(bb[2]-bb[0]) < 361: 179 # bb[0], bb[2]=-180, 180 180 #self.wgs84BBox = bb 181 return self.wgs84BBox 182 #raise NotImplementedError 183 184 def getSlab(self, crs, dimValues=None, renderOpts={}): 185 """ 186 Creates a slab of the layer in a particular CRS and set of 187 dimensions. 188 189 @param crs: The coordinate reference system. 190 @param dimValues: A mapping of dimension names to dimension values 191 as specified in the IDimension.extent 192 @param renderOpts: A generic mapping object for passing rendering 193 options 194 @return: An object implementing ILayerSlab 195 #create netcdf for whole lat/lon for given dimValues, use to init slab 196 """ 197 198 log.debug('getSlab(%s, %s)' % (crs, dimValues)) 199 200 #unicode conversion 201 for dimval in dimValues: 202 if dimval != 'time': 203 dimValues[dimval]=float(dimValues[dimval]) 204 else: 205 #remove any trailing Zs from time string 206 if dimValues[dimval] [-1:] in ['Z', 'z']: 207 dimValues[dimval]=dimValues[dimval][:-1] 208 if type(self._feature) == csml.parser.GridSeriesFeature: 209 randomname= csml.csmllibs.csmlextra.getRandomID() + '.nc' 210 result= self._feature.subsetToGridSeries(config['tmpfilebuffer'], ncname=randomname, **dimValues) 211 #for now have to read netcdf back from 212 #disk (limitiation of CSML api) 213 netcdf=cdms.open(result[1]) 214 #and then delete the temporary file 215 os.system('rm %s'%result[1]) 216 bbox=self.getBBox(crs) 217 return CSMLLayerSlab(netcdf, self, crs, dimValues, renderOpts, bbox) 218 else: 219 raise NotImplementedError 220 221 def getCacheKey(self, crs, dimValues=None, renderOpts={}): 222 """ 223 Create a unique key for use in caching a slab. 224 225 The intention here is that most of the work should be done when 226 instantiating an ILayerSlab object. These can be cached by the 227 server for future use. The server will first call getCacheKey() 228 for the slab creation arguments and if the key is in it's cache 229 it will use a pre-generated ILayerSlab object. 230 231 """ 232 dimList = list(dimValues.items()) 233 dimList.sort() 234 return '%s:%s' % (crs, dimList) 235 236 def getFeatureInfo(self, format, crs, point, dimValues): 237 """ 238 Return a response string descibing the feature at a given 239 point in a given CRS. 240 241 Currently only "html" is supported as output format 242 243 @param format: One of self.featureInfoFormats. Defines which 244 format the response will be in. 245 @param crs: One of self.crss 246 @param point: a tuple (x, y) in the supplied crs of the point 247 being selected. 248 @param dimValues: A mapping of dimension names to dimension values. 249 @return: A string containing the response. 250 251 """ 252 253 #cached netcdf is indexed by a tuple of feature id and dimvalues - i.e. typically a particular time and Z value for that feature. 254 #look in dictionary for cached copy, and if so use that as the file object. 255 dictindex=str((self._feature.id, dimValues)) 256 if dictindex in self.featureinfofilecache: 257 log.debug('calling cache') 258 f=self.featureinfofilecache[dictindex] 259 else: #else, use the csml api to subset the feature afresh 260 log.debug('not calling cache') 261 randomname= csml.csmllibs.csmlextra.getRandomID() + '.nc' 262 result= self._feature.subsetToGridSeries(config['tmpfilebuffer'], ncname=randomname, **dimValues) 263 #for now have to read netcdf back from disk (limitation of CSML api) 264 f=cdms.open(result[1]) 265 #append to cache: 266 self.featureinfofilecache[dictindex]=f 267 #and then delete the temporary file 268 os.system('rm %s'%result[1]) 269 270 netcdf = f(self.title) #netcdf here is a cdms transient variable 271 272 273 #Now grab the netCDF object for the point specified. 274 #The reason for the 'cob' option is so that if the grid the data 275 #is defined on does not have a grid point at the point specified, 276 #we should still get the nearest location 277 278 t_point = netcdf(latitude=(point[1], point[1], 'cob'), longitude=(point[0], point[0], 'cob')) 279 #now get the value recorded at this location 280 value = t_point.getValue().tolist() 281 log.debug(value) 282 log.debug(t_point.fill_value()) 283 #and the fill_value too 284 fill_value = t_point.fill_value() 285 #value is actually embedded in a multi dimensional list, 286 #so we need to extract the actual value from the list 287 while type(value) is list: 288 value = value[0] 289 290 #now check if the value is actually the fill_value rather than 291 #a value recorded at the point specified 292 log.debug('%s %s' % (value, fill_value)) 293 if (2*fill_value) == value: 294 value = "No value found at position: "+str(point[1])+", "+str(point[0]) 295 else: 296 value = "Value found at position: "+str(point[1])+", "+str(point[0])+" is: "+str(value) 297 # finally return the value 298 return value 299 300 301 class CSMLDimension(IwmsDimension): 302 """ 303 implements IDimension 304 @ivar units: The units string. 305 @ivar extent: Sequence of extent values. 306 307 """ 308 309 def __init__(self, domain, dimname, unit): 310 self.units = unit 311 self.extent = [] 312 #patch to handle current limitations of multiple time dimension scanning in csml. 313 if string.lower(self.units)[:10] in ['days_since', 'seconds_si', 'minutes_si', 'hours_sinc','months _sin', 'years_sinc']: 314 if type(domain[dimname][0]) is not str : 315 tunits=self.units.replace('_', ' ') 316 for val in domain[dimname]: 317 csmltime= csml.csmllibs.csmltime.UDtimeToCSMLtime(cdtime.reltime(float(val), tunits).tocomp()) 318 self.extent.append(csmltime) 319 else: 320 for val in domain[dimname]: 321 self.extent.append(str(val)) 322 else: 323 for val in domain[dimname]: 324 self.extent.append(str(val)) 325 #for time axis replace units with iso string 326 if dimname == 'time': 327 self.units='ISO8601' 328 329 class CSMLLayerSlab(IwmsLayerSlab): 330 """ 331 Implements LayerSlab 332 Represents a particular horizontal slice of a WMS layer. 333 334 ILayerSlab objects are designed to be convenient to cache. 335 They should be pickleable to enable memcached support in the future. 336 337 @ivar layer: The source ILayer instance. 338 @ivar crs: The coordinate reference system. 339 @ivar dimValues: A mapping of dimension values of this view. 340 @ivar renderOpts: The renderOpts used to create this view. 341 @ivar bbox: The bounding box as a 4-tuple. 342 """ 343 344 def __init__(self, netcdf, layer, crs, dimValues, renderOpts, bbox): 345 self._netcdf=netcdf 346 self.layer = layer 347 self.crs = crs 348 self.dimValues = dimValues 349 self.renderOpts=renderOpts 350 self.bbox=bbox 351 352 #set colour map for ALL images from this slab 353 v=self._netcdf(layer.title) 354 tvar=v(squeeze=1) 355 #get the min and max values to use for the colourmapping. 356 #change the fill values to ensure they aren't picked up as false max/mins. 357 tvar.missing_value=999999999 358 value=tvar.getValue() 359 self.minval=min(min(l) for l in value) 360 tvar.missing_value=-0.000000001 361 value=tvar.getValue() 362 self.maxval=max(max(l) for l in value) 363 364 365 def getImage(self, bbox, width, height): 366 """ 367 Create an image of a sub-bbox of a given size. 368 369 @ivar bbox: A bbox 4-tuple. 370 @ivar width: width in pixels.` 371 @ivar height: height in pixels. 372 @return: A PIL Image object. 373 374 """ 375 376 log.debug('getImage(%s, %s, %s)' % (bbox, width, height)) 377 378 379 cmap=eval(config['colourmap']) # renderOpts is hook for colourmap, for now use config 380 grid=Grid(self.layer, self._netcdf, bbox, width, height) 381 #how to handle varmin,varmax? ..read array? 382 #minval, maxval=genutil.minmax(grid.value) 383 #minval=min(min(l) for l in grid.value) 384 #maxval=max(max(l) for l in grid.value) 385 minval=self.minval 386 maxval=self.maxval 387 renderer=RGBARenderer(minval, maxval) 388 return renderer.renderGrid(grid, bbox, width, height, cmap) 389 390 class Grid(object): 391 """A class encapsulating a simple regularly spaced, rectilinear 392 grid. This is the only type of grid pywms is expected to 393 understand and adaptors should be provided to connect to 394 underlying implementations such as cdms or csml. 395 396 @cvar crs: Coordinate reference system 397 398 @ivar x0: coordinate of the lower left corner. 399 @ivar y0: coordinate of the lower left corner. 400 @ivar dx: The x grid spacing. 401 @ivar dy: The y grid spacing. 402 @ivar nx: The number of x grid points. 403 @ivar ny: The number of y grid points. 404 @ivar value: A masked array of the grid values. 405 @ivar ix: The dimension index of longidude in value 406 @ivar iy: The dimension index of latitude in value 407 @ivar long_name: The name of the field. 408 @ivar units: The units of the field. 409 """ 410 def __init__(self, layer, netcdf, bbox, width, height): 411 #we know the axes are called latitude and longitude as the CSML code has written it: 412 v=netcdf(layer.title) 413 tvar=v(latitude=(bbox[1], bbox[3]), longitude=(bbox[0],bbox[2]),squeeze=1) 414 order=tvar.getOrder() 415 #array of data 416 self.value=tvar.getValue() 417 if order == 'xy': 418 self.ix=0 419 self.iy=1 420 else: 421 self.ix=1 422 self.iy=0 423 lat = tvar.getLatitude() 424 lon = tvar.getLongitude() 425 self.x0=lon[0] 426 self.y0=lat[0] 427 self.dx=abs(lon[0]-lon[1]) 428 self.dy=abs(lat[0]-lat[1]) 429 430 431 432 433 434 self.nx=len(lon) 435 self.ny=len(lat) 436 self.long_name=tvar.id #TODO, get long name from feature 437 self.units=tvar.units -
cows/trunk/cows/service/imps/wms_csmllayer.py
r4266 r4270 1 """2 implementation of ILayerMapper, ILayer, IDimension, ILayerSlab interfaces, as defined in wms_iface.py3 4 """5 import os, string6 import csml7 try:8 import cdms2 as cdms9 except:10 import cdms11 import Image12 from copy import copy13 from pywms.render_imp import RGBARenderer14 from matplotlib import cm15 import genutil16 from pylons import config, request, session #config must have tmpfilebuffer and csmlstore values17 import cdtime18 import logging19 log = logging.getLogger(__name__)20 import ConfigParser21 22 from cows.service.wxs_iface import ILayerMapper23 from cows.service.wms_iface import IwmsLayer, IwmsDimension, IwmsLayerSlab24 25 try:26 from ndgUtils import ndgObject, ndgRetrieve27 except:28 log.warning("ndgUtils library could not be loaded, files in the eXist database won't be available, although you should still be able to access files from the csmlstore directory referenced in the ini file.")29 30 31 32 33 class CSMLLayerMapper(ILayerMapper):34 """35 Map keyword arguments to a collection of layers.36 Supports the retrieval of sets of layers according to arbitary37 keyword/value pairs.38 Implements ILayerMapper39 40 """41 def __init__(self):42 self.layermapcache={}43 44 def _getInfo(self, feature):45 ''' given a csml feature, return info about the layer/feature46 @return: title, abstract, dimensions, units, crss '''47 48 try:49 title=feature.name.CONTENT50 except:51 title=''52 try:53 abstract=feature.description.CONTENT54 except:55 abstract=title56 57 units=feature.getDomainUnits()58 dimensions={}59 tmpunits=copy(units)60 tmpunits.reverse()61 domain = feature.getDomain()62 for dim in feature.getAxisLabels():63 nextdim=CSMLDimension(domain, dim, tmpunits.pop())64 if dim not in ['latitude', 'longitude']:65 dimensions[dim]=nextdim66 crs=feature.getNativeCRS()67 crss=[self._crscat.getCRS(crs).twoD]68 if 'EPSG:4326' in crss:69 crss.append('CRS:84')70 crss.append('WGS84')71 return title, abstract, dimensions, units, crss72 73 74 75 def _get_csml_doc(self, fileoruri):76 """77 Gets a csml document from file or exist when passed a file name or uri78 79 Note, access control is not implemented on file access, only on exist access.80 81 """82 if string.find(fileoruri,'__NDG-A0__') == -1:83 #it's a local file not an identifier, read from local csmlstore84 file=fileoruri85 csml_dir = config['csmlstore']86 path = os.path.join(csml_dir, file)87 if os.path.exists(path+'.csml'):88 f = path+'.csml'89 elif os.path.exists(path+'.xml'):90 f = path +'.xml'91 else:92 raise ValueError("Cannot find CSML file %s" % file)93 94 d = csml.parser.Dataset()95 d.parse(f)96 97 else:98 #it's an NDG identifier, get the document from exist.99 uri=fileoruri100 uriN=ndgObject(uri)101 conflocation=config['ndgconfig']102 cf=ConfigParser.ConfigParser()103 cf.read(conflocation)104 status,x=ndgRetrieve(uriN, cf, discovery=0)105 d=csml.parser.Dataset()106 d.parse(x)107 return d108 109 def map(self, **kwargs):110 """111 Given csml.parser.Dataset object list the names of112 all layers available.113 114 @return: A mapping of layer names to ILayer implementations.115 @raise ValueError: If no layers are available for these keywords.116 """117 fileoruri=kwargs['fileoruri']118 if fileoruri in self.layermapcache.keys():119 #we've accessed this layer map before, get it from the cache dictionary120 return self.layermapcache[fileoruri]121 122 ds = self._get_csml_doc(fileoruri)123 layermap={}124 self._crscat=csml.csmllibs.csmlcrs.CRSCatalogue()125 for feature in csml.csmllibs.csmlextra.listify(ds.featureCollection.featureMembers):126 title, abstract, dimensions, units, crss=self._getInfo(feature)127 layermap[feature.id]=CSMLLayer(title,abstract, dimensions, units, crss, feature)128 if len(layermap) > 0:129 self.layermapcache[fileoruri]=layermap130 return layermap131 else:132 raise ValueError133 134 135 class CSMLLayer(IwmsLayer):136 """137 representing a WMS layer. Implements ILayer138 139 @ivar title: The layer title. As seen in the Capabilities document.140 @ivar abstract: Abstract as seen in the Capabilities document.141 @ivar dimensions: A dictionary of IDimension objects.142 @ivar units: A string describing the units.143 @ivar crss: A sequence of SRS/CRSs supported by this layer.144 145 @todo: Do we need minValue/maxValue?146 147 """148 149 def __init__(self, title, abstract, dimensions, units, crss, feature):150 self.featureInfoFormats=None #NotImplemented151 self.title=title152 self.abstract=abstract153 self.dimensions=dimensions154 self.units=units155 self.crss=crss156 self._feature=feature157 self.legendSize=(30,100)158 bb= self._feature.getCSMLBoundingBox().getBox()159 #convert 0 - 360 to -180, 180 as per common WMS convention160 if abs(bb[2]-bb[0]) >= 359 and abs(bb[2]-bb[0]) < 361:161 bb[0], bb[2]=-180, 180162 self.wgs84BBox = bb163 self.featureInfoFormats = ['text/html']164 try:165 self.wgs84BBox = self.getBBox('EPSG:4326')166 except:167 raise ValueError("Layer must provide a bounding box in EPSG:4326 "168 "coordinates for compatibility with WMS-1.3.0")169 self.featureinfofilecache={} #used for caching netcdf file in getFeatureInfo170 171 def getBBox(self, crs):172 """173 @return: A 4-typle of the bounding box in the given coordinate174 reference system.175 """176 #bb= self._feature.getCSMLBoundingBox().getBox()177 #convert 0 - 360 to -180, 180 as per common WMS convention178 #if abs(bb[2]-bb[0]) >= 359 and abs(bb[2]-bb[0]) < 361:179 # bb[0], bb[2]=-180, 180180 #self.wgs84BBox = bb181 return self.wgs84BBox182 #raise NotImplementedError183 184 def getSlab(self, crs, dimValues=None, renderOpts={}):185 """186 Creates a slab of the layer in a particular CRS and set of187 dimensions.188 189 @param crs: The coordinate reference system.190 @param dimValues: A mapping of dimension names to dimension values191 as specified in the IDimension.extent192 @param renderOpts: A generic mapping object for passing rendering193 options194 @return: An object implementing ILayerSlab195 #create netcdf for whole lat/lon for given dimValues, use to init slab196 """197 198 log.debug('getSlab(%s, %s)' % (crs, dimValues))199 200 #unicode conversion201 for dimval in dimValues:202 if dimval != 'time':203 dimValues[dimval]=float(dimValues[dimval])204 else:205 #remove any trailing Zs from time string206 if dimValues[dimval] [-1:] in ['Z', 'z']:207 dimValues[dimval]=dimValues[dimval][:-1]208 if type(self._feature) == csml.parser.GridSeriesFeature:209 randomname= csml.csmllibs.csmlextra.getRandomID() + '.nc'210 result= self._feature.subsetToGridSeries(config['tmpfilebuffer'], ncname=randomname, **dimValues)211 #for now have to read netcdf back from212 #disk (limitiation of CSML api)213 netcdf=cdms.open(result[1])214 #and then delete the temporary file215 os.system('rm %s'%result[1])216 bbox=self.getBBox(crs)217 return CSMLLayerSlab(netcdf, self, crs, dimValues, renderOpts, bbox)218 else:219 raise NotImplementedError220 221 def getCacheKey(self, crs, dimValues=None, renderOpts={}):222 """223 Create a unique key for use in caching a slab.224 225 The intention here is that most of the work should be done when226 instantiating an ILayerSlab object. These can be cached by the227 server for future use. The server will first call getCacheKey()228 for the slab creation arguments and if the key is in it's cache229 it will use a pre-generated ILayerSlab object.230 231 """232 dimList = list(dimValues.items())233 dimList.sort()234 return '%s:%s' % (crs, dimList)235 236 def getFeatureInfo(self, format, crs, point, dimValues):237 """238 Return a response string descibing the feature at a given239 point in a given CRS.240 241 Currently only "html" is supported as output format242 243 @param format: One of self.featureInfoFormats. Defines which244 format the response will be in.245 @param crs: One of self.crss246 @param point: a tuple (x, y) in the supplied crs of the point247 being selected.248 @param dimValues: A mapping of dimension names to dimension values.249 @return: A string containing the response.250 251 """252 253 #cached netcdf is indexed by a tuple of feature id and dimvalues - i.e. typically a particular time and Z value for that feature.254 #look in dictionary for cached copy, and if so use that as the file object.255 dictindex=str((self._feature.id, dimValues))256 if dictindex in self.featureinfofilecache:257 log.debug('calling cache')258 f=self.featureinfofilecache[dictindex]259 else: #else, use the csml api to subset the feature afresh260 log.debug('not calling cache')261 randomname= csml.csmllibs.csmlextra.getRandomID() + '.nc'262 result= self._feature.subsetToGridSeries(config['tmpfilebuffer'], ncname=randomname, **dimValues)263 #for now have to read netcdf back from disk (limitation of CSML api)264 f=cdms.open(result[1])265 #append to cache:266 self.featureinfofilecache[dictindex]=f267 #and then delete the temporary file268 os.system('rm %s'%result[1])269 270 netcdf = f(self.title) #netcdf here is a cdms transient variable271 272 273 #Now grab the netCDF object for the point specified.274 #The reason for the 'cob' option is so that if the grid the data275 #is defined on does not have a grid point at the point specified,276 #we should still get the nearest location277 278 t_point = netcdf(latitude=(point[1], point[1], 'cob'), longitude=(point[0], point[0], 'cob'))279 #now get the value recorded at this location280 value = t_point.getValue().tolist()281 log.debug(value)282 log.debug(t_point.fill_value())283 #and the fill_value too284 fill_value = t_point.fill_value()285 #value is actually embedded in a multi dimensional list,286 #so we need to extract the actual value from the list287 while type(value) is list:288 value = value[0]289 290 #now check if the value is actually the fill_value rather than291 #a value recorded at the point specified292 log.debug('%s %s' % (value, fill_value))293 if (2*fill_value) == value:294 value = "No value found at position: "+str(point[1])+", "+str(point[0])295 else:296 value = "Value found at position: "+str(point[1])+", "+str(point[0])+" is: "+str(value)297 # finally return the value298 return value299 300 301 class CSMLDimension(IwmsDimension):302 """303 implements IDimension304 @ivar units: The units string.305 @ivar extent: Sequence of extent values.306 307 """308 309 def __init__(self, domain, dimname, unit):310 self.units = unit311 self.extent = []312 #patch to handle current limitations of multiple time dimension scanning in csml.313 if string.lower(self.units)[:10] in ['days_since', 'seconds_si', 'minutes_si', 'hours_sinc','months _sin', 'years_sinc']:314 if type(domain[dimname][0]) is not str :315 tunits=self.units.replace('_', ' ')316 for val in domain[dimname]:317 csmltime= csml.csmllibs.csmltime.UDtimeToCSMLtime(cdtime.reltime(float(val), tunits).tocomp())318 self.extent.append(csmltime)319 else:320 for val in domain[dimname]:321 self.extent.append(str(val))322 else:323 for val in domain[dimname]:324 self.extent.append(str(val))325 #for time axis replace units with iso string326 if dimname == 'time':327 self.units='ISO8601'328 329 class CSMLLayerSlab(IwmsLayerSlab):330 """331 Implements LayerSlab332 Represents a particular horizontal slice of a WMS layer.333 334 ILayerSlab objects are designed to be convenient to cache.335 They should be pickleable to enable memcached support in the future.336 337 @ivar layer: The source ILayer instance.338 @ivar crs: The coordinate reference system.339 @ivar dimValues: A mapping of dimension values of this view.340 @ivar renderOpts: The renderOpts used to create this view.341 @ivar bbox: The bounding box as a 4-tuple.342 """343 344 def __init__(self, netcdf, layer, crs, dimValues, renderOpts, bbox):345 self._netcdf=netcdf346 self.layer = layer347 self.crs = crs348 self.dimValues = dimValues349 self.renderOpts=renderOpts350 self.bbox=bbox351 352 #set colour map for ALL images from this slab353 v=self._netcdf(layer.title)354 tvar=v(squeeze=1)355 #get the min and max values to use for the colourmapping.356 #change the fill values to ensure they aren't picked up as false max/mins.357 tvar.missing_value=999999999358 value=tvar.getValue()359 self.minval=min(min(l) for l in value)360 tvar.missing_value=-0.000000001361 value=tvar.getValue()362 self.maxval=max(max(l) for l in value)363 364 365 def getImage(self, bbox, width, height):366 """367 Create an image of a sub-bbox of a given size.368 369 @ivar bbox: A bbox 4-tuple.370 @ivar width: width in pixels.`371 @ivar height: height in pixels.372 @return: A PIL Image object.373 374 """375 376 log.debug('getImage(%s, %s, %s)' % (bbox, width, height))377 378 379 cmap=eval(config['colourmap']) # renderOpts is hook for colourmap, for now use config380 grid=Grid(self.layer, self._netcdf, bbox, width, height)381 #how to handle varmin,varmax? ..read array?382 #minval, maxval=genutil.minmax(grid.value)383 #minval=min(min(l) for l in grid.value)384 #maxval=max(max(l) for l in grid.value)385 minval=self.minval386 maxval=self.maxval387 renderer=RGBARenderer(minval, maxval)388 return renderer.renderGrid(grid, bbox, width, height, cmap)389 390 class Grid(object):391 """A class encapsulating a simple regularly spaced, rectilinear392 grid. This is the only type of grid pywms is expected to393 understand and adaptors should be provided to connect to394 underlying implementations such as cdms or csml.395 396 @cvar crs: Coordinate reference system397 398 @ivar x0: coordinate of the lower left corner.399 @ivar y0: coordinate of the lower left corner.400 @ivar dx: The x grid spacing.401 @ivar dy: The y grid spacing.402 @ivar nx: The number of x grid points.403 @ivar ny: The number of y grid points.404 @ivar value: A masked array of the grid values.405 @ivar ix: The dimension index of longidude in value406 @ivar iy: The dimension index of latitude in value407 @ivar long_name: The name of the field.408 @ivar units: The units of the field.409 """410 def __init__(self, layer, netcdf, bbox, width, height):411 #we know the axes are called latitude and longitude as the CSML code has written it:412 v=netcdf(layer.title)413 tvar=v(latitude=(bbox[1], bbox[3]), longitude=(bbox[0],bbox[2]),squeeze=1)414 order=tvar.getOrder()415 #array of data416 self.value=tvar.getValue()417 if order == 'xy':418 self.ix=0419 self.iy=1420 else:421 self.ix=1422 self.iy=0423 lat = tvar.getLatitude()424 lon = tvar.getLongitude()425 self.x0=lon[0]426 self.y0=lat[0]427 self.dx=abs(lon[0]-lon[1])428 self.dy=abs(lat[0]-lat[1])429 430 431 432 433 434 self.nx=len(lon)435 self.ny=len(lat)436 self.long_name=tvar.id #TODO, get long name from feature437 self.units=tvar.units
Note: See TracChangeset
for help on using the changeset viewer.