source: TI03-DataExtractor/trunk/pydxs/CDMSDataHandler.py @ 1715

Subversion URL: http://proj.badc.rl.ac.uk/svn/ndg/TI03-DataExtractor/trunk/pydxs/CDMSDataHandler.py@1715
Revision 1715, 17.4 KB checked in by astephen, 13 years ago (diff)

Merged with titania version.

Line 
1#   Copyright (C) 2004 CCLRC & NERC( Natural Environment Research Council ).
2#   This software may be distributed under the terms of the
3#   Q Public License, version 1.0 or later. http://ndg.nerc.ac.uk/public_docs/QPublic_license.txt
4
5"""
6CDMSDataHandler.py
7==================
8
9CDMSDataHandler module for the dx package.
10
11This module holds the CDMSDataHandler class that is used
12to hold and access information about datasets held in CDMS-type
13formats visible to the dx package.
14
15"""
16
17# Import required modules
18import os
19import cdms
20import re
21
22# Import global variables
23from serverConfig import *
24from common import *
25from DXDMLHandler import DXDMLHandler
26from DXErrors import *
27
28
29class CDMSDataHandler:
30    """
31    A group of methods to connect to a dataset group or
32    dataset to extract information about the contents.
33    """
34 
35    def __init__(self, datasetURI=None):
36        """
37        Set up instance variables.
38        """
39        self.DXDML=DXDMLHandler()
40        self.file=datasetURI
41        if self.file: self._openDataFile(datasetURI=self.file)
42   
43   
44    def _openDataFile(self, datasetGroup=None, dataset=None, datasetURI=None):
45        """
46        Opens a file and allocates to file handle called: self.file.
47        """
48        if datasetURI:
49            cdmlfile=datasetURI
50        else:
51            for item in self.DXDML.getDatasetsAndDatasetURIs(datasetGroup):
52                if item[0]==dataset:
53                    cdmlfile=item[1]
54       
55        try:
56            self.file=cdms.open(cdmlfile)
57        except IOError, error:
58            raise DXDataIOError, error
59
60
61    def _getVariable(self, varname):
62        """
63        Gets variable metadata object from a data file.
64        """
65        try:
66            rtvalue=self.file[varname]
67        except:
68            raise DXOptionHandlingError, "Cannot find variable %s in file %s" % (varname, self.file.id)
69        return rtvalue
70       
71
72    def _getBestName(self, v, vcount=0):
73        """
74        Returns the best name for a cdms variable.
75        """
76        if not hasattr(v, "standard_name"): 
77            if not hasattr(v, "long_name"):
78                if not hasattr(v, "title"):
79                    if not hasattr(v, "name"):
80                        if hasattr(v, "id"):
81                            name=v.id
82                        else:
83                            vcount=vcount+1
84                            name="unnamed_var_%s" % vcount
85                    else:
86                       name=v.name
87                else:
88                    name=v.title
89            else:
90                name=v.long_name
91        else:
92            name=v.standard_name
93
94        if hasattr(v, "cell_methods"):
95            name="%s [%s]" % (name, v.cell_methods)
96
97        return name
98
99
100    def getVariables(self, datasetGroup=None, dataset=None, datasetURI=None):
101        """
102        Returns a list of variables for the given dataset
103        group/dataset combination or datasetURI. The variable name used is selected
104        hierarchically depending on the available attributes. Each returned item in
105        the list includes a [<long_name>, <id>].
106        """ 
107        self._openDataFile(datasetGroup, dataset, datasetURI)
108        vars=self.file.listvariables()
109        rtvars=[]
110        vcount=0
111        for var in vars:
112           v=self.file[var]
113           name=self._getBestName(v, vcount)
114               
115           # Fix name to remove leading asterisks and lower case surface start.
116           name=name.replace("_", " ")
117           if name[:2]=="**": name=(name[2:]).strip()
118           if name[:7]=="surface": name=(name[7:]).strip()
119           # Remove variables they are actually bounds on axes or coefficients in formulae
120           if v.id not in ("bounds_longitude", "bounds_latitude", "bounds_level", "bounds_time", "p0"):
121               rtvars.append([name, v.id]) 
122
123        return rtvars
124
125
126    def getDomain(self, datasetGroup=None, dataset=None, variable=None, datasetURI=None):
127        """
128        Returns the full domain listing for a variable returning:
129       
130        [knownAxisString, id, longName, units, listType, unusedItem,
131        listValue-1, listValue-2, ..., listValue-n]
132       
133        For example:
134       
135        ["time", "time", "Time", "hours since 1999-09-09 00:00:00", "start end interval",
136        "", 0, 3, 6]
137       
138        This listType represents 6-hourly time steps of 0,1,2,3 past the base time.
139       
140        listType can also take the value "full list" where all values in the list are provided,
141        or "start end" where only the first and last value are given.
142        """ 
143        self._openDataFile(datasetGroup, dataset, datasetURI)
144        var=self._getVariable(variable)
145        varList=self.file.listvariables()
146        varList.sort()
147        varIndex=varList.index(var.id)
148        rtlist=[]
149       
150        axcount=0
151        for axis in var.getAxisList():
152            #axisIndexString="axis%s.%s" % ("77777", axcount) #(varIndex+1, axcount)
153            units=None
154            axisSpecificInfo=""
155
156            if axis.isTime():
157                knownAxis="time"
158                (start, end, (intervalValue, intervalUnits))=self.getTemporalDomain(datasetGroup, dataset, variable, datasetURI)
159                arrayValues=[start, end, intervalValue]
160                if hasattr(axis, "calendar") and axis.calendar=="360_day":
161                    axisSpecificInfo=axis.calendar
162                listType="start end interval"
163                units=intervalUnits
164            elif axis.isLevel():
165                knownAxis="level"
166                arrayValues=axis[:]
167                listType="full list"
168            elif axis.isLatitude():
169                knownAxis="latitude"
170                arrayValues=[axis[0], axis[-1]]
171                arrayValues.sort()
172                arrayValues.reverse()
173                listType="start end" 
174            elif axis.isLongitude():
175                knownAxis="longitude"
176                arrayValues=[axis[0], axis[-1]]
177                arrayValues.sort()     
178                listType="start end"                                   
179            else:
180                # For any axis not known as above
181                knownAxis=""
182                if len(axis[:])>200:
183                    arrayValues=[axis[0], axis[-1]]
184                    listType="start end"                   
185                else:
186                    arrayValues=axis[:]
187                    listType="full list"
188           
189            id=axis.id
190            longName=self._getBestName(axis).title()
191            if not units:  units=getattr(axis, "units", "")
192
193            if type(arrayValues)!=type([2,3]):
194                arrayValues=list(arrayValues)
195
196            # Need two values minimum so re-use if only one but really don't want to send at all
197            if len(arrayValues)==1: 
198                arrayValues=[arrayValues[0]]*2
199           
200            rtlist.append([knownAxis, id, longName, units, listType, axisSpecificInfo]+arrayValues)
201            axcount=axcount+1
202           
203        return rtlist
204       
205
206    def getHorizontalDomain(self, datasetGroup=None, dataset=None, variable=None, datasetURI=None):
207        """
208        Returns the horizontal domain as (northernExtent, westernExtent, southernExtent, easternExtent).
209        """
210        self._openDataFile(datasetGroup, dataset, datasetURI)
211        var=self._getVariable(variable)
212        lat=list(var.getLatitude()[:])
213        if lat[-1]<lat[0]: lat.reverse()
214        (southernExtent, northernExtent)=(lat[0], lat[-1])
215        lon=var.getLongitude()[:]
216        (westernExtent, easternExtent)=(lon[0], lon[-1])
217        return (northernExtent, westernExtent, southernExtent, easternExtent)
218
219
220    def getVerticalSpatialDomain(self, datasetGroup=None, dataset=None, variable=None, datasetURI=None):
221        """
222        Returns the vertical domain as a tuple containing
223        a list of levels (or "Single level" string) and the units.
224        """
225        self._openDataFile(datasetGroup, dataset, datasetURI)
226        var=self._getVariable(variable)
227        try:
228            levels=var.getLevel()
229            vertical_units=levels.units
230            vertical_domain=tuple(map(lambda x: float(x), levels[:]))
231     
232        except:
233            vertical_domain=("Single level",)
234            vertical_units=None
235        return (vertical_domain, vertical_units)
236
237
238    def getTemporalDomain(self, datasetGroup=None, dataset=None, variable=None, datasetURI=None):
239        """
240        Returns the temporal domain as a tuple of (start time, end time,
241        (interval value, interval units)).
242        """
243        self._openDataFile(datasetGroup, dataset, datasetURI)
244        var=self._getVariable(variable)
245        time_axis=var.getTime()
246        #time_keys=("year", "month", "day", "hour", "minute", "second")
247        start_time=str(time_axis.asComponentTime()[0]).replace(" ", "T")
248        end_time=str(time_axis.asComponentTime()[-1]).replace(" ", "T")
249        """start_time=[]
250        end_time=[]
251        for key in time_keys:
252            start_time.append(getattr(start_time_components, key))
253            end_time.append(getattr(end_time_components, key))""" 
254
255        time_units=time_axis.units.split()[0]
256        if time_units[-1]=="s":  time_units=time_units[:-1]
257        if len(time_axis)>1:
258            interval_value=abs(time_axis[1]-time_axis[0])
259        else:
260            interval_value=1  # To stop loops breaking later!
261        return (start_time, end_time, (interval_value, time_units)) 
262
263
264    def getSelectedTimeSteps(self, datasetURI, variable, axisSelectionDict):
265        """
266        Returns a list of time step strings based on the selection.
267        """     
268        self._openDataFile(datasetURI=datasetURI)
269        var=self._getVariable(variable)
270       
271        timeAxis=None
272        axcount=0
273        for axis in var.getAxisList():
274            if axis.isTime():
275                timeAxisIndex=axcount
276                timeAxis=axis
277            axcount=axcount+1
278
279        if timeAxis==None:
280            return []
281       
282        startDateTime=None
283        for key in axisSelectionDict.keys():
284            #print key, axisSelectionDict[key]
285            axisIndex=int(key.split(".")[-1])-1
286            if axisIndex==timeAxisIndex:
287                (startDateTime, endDateTime)=axisSelectionDict[key][:2]
288       
289        if startDateTime==None:
290            return [str(tst) for tst in timeAxis.asComponentTime()]
291       
292        startDateTime=startDateTime.replace("T", " ")
293        items=startDateTime.split(":")
294        startDateTime=":".join(items[:-1])+":"+("%f" % float(items[-1]))
295        endDateTime=endDateTime.replace("T", " ")
296        items=endDateTime.split(":")
297        endDateTime=":".join(items[:-1])+":"+("%f" % float(items[-1]))
298       
299        timeSteps=timeAxis.asComponentTime()
300        selectedTimes=[]
301
302        for timeStep in timeSteps:
303            ts=timeStep
304            timeStep="%.4d-%.2d-%.2d %.2d:%.2d:%f" % (ts.year, ts.month, ts.day, \
305                                            ts.hour, ts.minute, ts.second)
306            #print str(timeStep), startDateTime
307            ts=str(timeStep)
308            if ts>endDateTime:
309                break
310            elif ts<startDateTime:
311                continue
312            else:
313                selectedTimes.append(ts)
314       
315        if selectedTimes==[]:
316            raise DXOptionHandlingError, "All selected time steps for '%s' are out of range, please go back and re-select." % variable
317               
318        return selectedTimes
319
320
321    def getSelectedVariableArrayDetails(self, datasetURI, variable, axisSelectionDict):
322        """
323        Returns a tuple representing the (array shape, grid shape, size)
324        of the selected subset of a variable. Grid shape can be None if both latitude
325        and longitude axes are not present.
326        """     
327        self._openDataFile(datasetURI=datasetURI)
328        var=self._getVariable(variable)
329        varType=var.typecode()
330               
331        timeAxisIndex=None
332        axcount=0
333        for axis in var.getAxisList():
334            if axis.isTime():
335                timeAxisIndex=axcount
336            axcount=axcount+1
337         
338        startDateTime=None   
339       
340        axesCounted=[]
341        axLens=[]
342        latLength=None
343        lonLength=None
344        size=1
345
346        seenAxes=[]
347        for key in axisSelectionDict.keys():
348            #print key, axisSelectionDict[key]
349            # Set the axisIndex as value minus one because we index from 1 in the dx lists
350            axisIndex=int(key.split(".")[-1])-1
351            (low, high)=axisSelectionDict[key][:2]
352            axis=var.getAxis(axisIndex)
353           
354            axisID=axis.id
355            if axisID in seenAxes: 
356                continue
357            else:
358                seenAxes.append(axisID)
359
360            if axisIndex==timeAxisIndex:
361                axlen=len(self.getSelectedTimeSteps(datasetURI, variable, {key:(low, high)}))           
362            elif low==high: 
363                axlen=1
364            else:
365                axlen=len(getValuesInRange(low, high, axis[:]))
366                if axlen==0:
367                    if (axis.isLongitude() or axis.isLatitude()) and low==high:
368                        print "Lat and lon can be axis length zero because we'll nudge to nearest if only one value given."
369                        axlen=1
370                    else:
371                        raise DXOptionHandlingError, "All selected '%s' axis values for '%s' are out of range, please go back and re-select." % (axis.id, variable)
372       
373                if axis.isLatitude():
374                    latLength=axlen
375                elif axis.isLongitude():
376                    lonLength=axlen
377                       
378            size=size*axlen
379            print "Axlen and size:", axlen, size
380            axesCounted.append(axisIndex)
381            axLens.append(axlen)
382       
383#       print axesCounted
384#       print axLens
385       
386        axcount=0
387        arrayShape=[]
388        for axis in var.getAxisList():
389            print "Sizing:", size
390            if axcount not in axesCounted:
391                size=size*len(axis)
392                arrayShape.append(len(axis))
393            else:
394                #print axcount
395                arrayShape.append(axLens[axesCounted.index(axcount)])   
396            if axis.isLatitude() and latLength==None:
397                latLength=len(axis)
398            elif axis.isLongitude() and lonLength==None:
399                lonLength=len(axis)   
400            axcount=axcount+1
401       
402        # Now work out gridShape if appropriate
403        if latLength and lonLength:
404            gridShape=(latLength, lonLength)
405        else:
406            gridShape=None
407       
408        if varType=="f":
409            size=size*4.
410        elif varType=="d":
411            size=size*8
412        elif varType=="i":
413            size=size
414
415        print "SUPERSIZE:", size
416
417       
418        return (tuple(arrayShape), gridShape, size)
419
420
421    def getSelectedVariableSubsetSize(self, datasetURI, varID, axisSelectionDict):
422        """
423        Returns the size in bytes of the selected subset of a variable.
424        """
425        return self.getSelectedVariableArrayDetails(datasetURI, varID, axisSelectionDict)[2]
426
427       
428    def readVariableSubsetIntoMemory(self, datasetURI, variable, axisSelectionDict, timeStep=None):
429        """
430        Reads the variable with ID 'variable' into memory from file
431        'datasetURI' - sub-setting across all axes indicated in 'axisSelectionDict'.
432        If 'timeStep' is provided then override the time selection in 'axisSelectionDict'
433        with the 'timeStep' given.
434        """
435        self._openDataFile(datasetURI=datasetURI)
436        var=self._getVariable(variable)
437
438        axisList=var.getAxisList()
439           
440        selectionStrings=[]
441        for key in axisSelectionDict.keys():
442            #print key, axisSelectionDict[key]
443            axisIndex=int(key.split(".")[-1])-1
444            axis=axisList[axisIndex]   
445            id=axis.id
446            print id
447           
448            # deal with time differently
449            if axis.isTime():
450                if timeStep!=None:     
451                    timeStep=timeStep.replace("T", " ")                 
452                    selectionStrings.append('time="%s"' % timeStep)
453                else:
454                    selector=axisSelectionDict[key][:2]
455                    selector=(selector[0].replace("T", " "), selector[1].replace("T", " "))
456                    selectionStrings.append('%s=%s' % (axis.id, str(selector)))         
457            elif axis.isLatitude():
458                (low, high)=axisSelectionDict[key][:2]
459                if low==high:
460                    loworig=low
461                    (low, high, nudgeMessage)=nudgeSingleValuesToAxisValues(low, axis[:], "Latitude")   
462                    if loworig!=low:   print "Nudged latitudes to nearest points..."
463                selectionStrings.append('%s=(%s, %s)' % (axis.id, low, high))   
464            elif axis.isLongitude():
465                (low, high)=axisSelectionDict[key][:2]
466                if low==high:
467                    loworig=low
468                    (low, high, nudgeMessage)=nudgeSingleValuesToAxisValues(low, axis[:], "Longitude") 
469                    if loworig!=low:   print "Nudged latitudes to nearest points..."
470                selectionStrings.append('%s=(%s, %s)' % (axis.id, low, high))           
471            else:
472                selector=axisSelectionDict[key][:2]
473                s=str(selector)
474                if s[0]=="[" and s[-1]=="]":
475                    s="("+s[1:-1]+")"
476                selectionStrings.append('%s=%s' % (axis.id, s))
477       
478        fullSelectionString=", ".join(selectionStrings)
479        print fullSelectionString
480        variableData=eval("self.file('%s', %s)" % (variable, fullSelectionString))
481        return variableData         
482       
483
484    def getCFGlobalAttributes(self, datafile):
485        """
486        Gets any CF metadta global attributes that are available
487        from the source dataset/file.
488        """
489        # Make sure data file is open
490        if self.file==None: self._openDataFile(datasetURI=datafile)
491        gatts={}
492
493        for gatt in CF_METADATA_GLOBAL_ATTRIBUTE_KEYS:
494            if hasattr(self.file, gatt):
495                gatts[gatt]=self.file.__getattr__(gatt)
496       
497        return gatts
498
499
500if __name__=="__main__":       
501    a=CDMSDataHandler()
502    print a.getDomain(variable="vo", datasetURI="file:/badc/ecmwf-e40/metadata/e40-1.0-am-195709-200208.xml")
503    print a.getVariables(datasetGroup='Test Data Group 1', dataset='Test Dataset 1')
504    print a.getVariables(datasetGroup='Test Data Group 3', datasetURI='file:/usr/local/test/dxs/testdata/testdata3.xml')   
505    print a.getDomain('Test Data Group 3', 'Test Dataset 3', "var3")
506    print a.getHorizontalDomain('Test Data Group 3', 'Test Dataset 3', "var3")
507    print a.getVerticalSpatialDomain('Test Data Group 3', 'Test Dataset 3', "var3")
508    print a.getTemporalDomain('Test Data Group 3', 'Test Dataset 3', "var3")
509    print a.getCFGlobalAttributes("file:/usr/local/test/dxs/testdata/testdata3.xml")
510    print a.getSelectedVariableArrayDetails('file:/usr/local/test/dxs/testdata/testdata1.xml', "pqn", 
511                     {"axis_1.1.1.1":("1999-01-01T00:00:00", "1999-01-01T06:00:00"),
512                     "axis_1.1.1.2":(10,20), "axis_1.1.1.3":(0,90)})
513    print a.getSelectedVariableArrayDetails('file:/usr/local/test/dxs/testdata/testdata1.xml', "pqn", 
514                     {"axis_1.1.1.1":("1999-01-01T00:00:00", "1999-01-01T06:00:00"),
515                     "axis_1.1.1.3":(30,-30)})
516    print a.getSelectedVariableSubsetSize("/usr/local/test/dxs/testdata/testdata1.xml", "pqn",
517                     {'axis_1.1.1.1':('1999-01-01T00:00:0.000000', '1999-01-01T00:00:0.000000'), 
518                      'axis_1.1.1.3': [0.0, 355.0], 'axis_1.1.1.2': [-90.0, 90.0]})
519                     
520    print a.readVariableSubsetIntoMemory("/usr/local/test/dxs/testdata/testdata1.xml", "pqn",
521                     {'axis_1.1.1.1':('1999-01-01T00:00:0.000000', '1999-01-01T00:00:0.000000'), 
522                      'axis_1.1.1.3': [0.0, 35.0], 'axis_1.1.1.2': [40.0, 90.0]})
Note: See TracBrowser for help on using the repository browser.