source: TI02-CSML/branches/csml-cdms2/csmllibs/csmldataiface.py @ 3627

Subversion URL: http://proj.badc.rl.ac.uk/svn/ndg/TI02-CSML/branches/csml-cdms2/csmllibs/csmldataiface.py@3627
Revision 3627, 29.5 KB checked in by spascoe, 12 years ago (diff)

This branch contains CSML converted to use cdat_lite-5.

  • convertcdms was run on the source
  • the MA.set_print_limit call was changed to the numpy equivilent
  • The tests were changed to account for the existence of numpy scalar types.

All tests, except the two known to fail, pass on i686 ubuntu.

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
Line 
1#!/usr/bin/env python
2# Adapted for numpy/ma/cdms2 by convertcdms.py
3
4'''
5**************************************************************************************
6csmldataIface.py
7contains classes for interfacing with various files
8currently supports cdunif (NetCDF, PP, Grib(untested)) And Nappy (NASAAmes)
9use by instantiating the factory class: DataInterface
10last updated 28 January 2008
11Dominic Lowe, BADC
12**************************************************************************************
13'''
14
15import pdb
16import cdms2 as cdms 
17import numpy.oldnumeric as Numeric
18
19try:
20    import nappy 
21except ImportError:
22    print "NASAAmes interface not available. CSML will still work, but won't support NASA Ames files"
23import string
24import sys
25import csml.csmllibs.csmltime
26# This is required to prevent Numeric arrays being truncated when printed.
27import numpy.oldnumeric.ma as MA
28#MA.set_print_limit(0)
29import numpy; numpy.set_printoptions(threshold=0)
30
31# Generic mathematical expression solver, required by the raw and image
32# interfaces
33import numpy.oldnumeric as NumericTransform
34
35try:
36    # Part of the PIL. Required for ImageFileInterface:
37    import Image
38except ImportError:
39    print "Python Image Library not available. CSML will still work, but won't support image files"
40
41try:
42    # Required for RawFileInterface:
43    import struct
44except ImportError:
45    print 'could not import struct'
46    pass
47
48
49class DataInterface(object):
50        '''DataInterface is the entry point for all data types. Use DataInterface and setInterfaceType to instantiate the correct
51        subclass for data''' 
52        def __init__(self):
53                self.iface ='None'
54               
55        def setInterfaceType(self,interfaceType):
56                '''returns approprate data interface
57                @param interfaceType:     should correspond to available datainterface subclasses
58                @return:    a subclass of DataInterface e.g. a ImageFileInterface'''
59                self.iface =interfaceType
60                if self.iface == 'nappy':
61                        return NappyInterface()
62                elif self.iface == 'cdunif':
63                        return cdunifInterface()
64                elif self.iface == 'raw':
65                        return RawFileInterface()
66                elif self.iface == 'image':
67                        return ImageFileInterface()
68               
69               
70        def getUnknownInterfaceType(self, filename):
71                '''if the interface type is not known at the time of instantiation, then use
72                this method to examine the file and return the correct interface (if it exists).
73                @param filename:      name of data file'''
74                fileExtension = str(filename)[-3:]
75                if fileExtension == '.nc':
76                        return cdunifInterface()
77                if fileExtension == 'qxf':
78                        return cdunifInterface()
79                elif fileExtension == '.pp':
80                        return cdunifInterface()
81                elif fileExtension == 'ctl':
82                        return cdunifInterface()
83                elif fileExtension == 'xml':
84                        return cdmlInterface()
85                elif fileExtensions == 'png' or fileExtension == 'gif':
86                        return ImageFileInterface()
87                else:
88                        try:
89                                nappy.readFFI(filename) in [1001,1010,1020,2010,2110,2160,2310,3010,4010]
90                                return NappyInterface()
91                        except:
92                                print "Could not establish file type"
93                                print "See csmldataiface.py"
94                       
95       
96                                                       
97class AbstractDI(object):               
98        '''Abstract data interface class
99        does nothing, but contains templates for methods required for a data interface class
100        individual interfaces (e.g NappyInterface) should override these methods. This is the minimum
101        set required for data access and subsetting.
102        Note that if scanning of data (for CSML creation) is required more methods may be needed - see cdunif interface'''
103
104        def __init__(self):
105                ''' When subclassed: unique 'type'  and 'prefix' must be defined'''
106                self.extractType=''
107                self.extractPrefix = ''
108                                       
109        def openFile(self, filename):
110           ''' When subclassed: opens file, sets self.file'''
111           raise NotImplementedError 
112               
113        def closeFile(self):
114           '''closes file, probably needs to be overwritten by subclass'''
115           try:
116               self.file.close()
117           except:
118               raise NotImplementedError
119         
120        def setAxis(self,axis):
121               '''set' the name of the current axis , must be overwritten by subclass
122           this may just involve a stub (see NASAAmes interface) or may involve
123           calling a real set method of the underlying api (see cdunif Interface)'''
124               raise NotImplementedError 
125                       
126        def getDataForAxis(self):
127               ''' return all data for axis, must be overwritten by subclass '''
128               raise NotImplementedError       
129       
130        def setVariable(self,varname):
131             '''As for setAxis, 'set' the name of the current axis , must be overwritten by subclass'''
132             raise NotImplementedError
133       
134        def getDataForVar(self):
135               '''return all data for variable, must be overwritten by subclass'''
136               raise NotImplementedError       
137       
138        def getSubsetOfDataForVar(self,**kwargs):
139            '''return subset of data for variable, must be overwritten by subclass'''
140            raise NotImplementedError
141       
142       
143class NappyInterface(AbstractDI):       
144        ''' Data Interface for Nappy (NASA Ames Processing in Python)'''
145   
146        def __init__(self):
147                self.extractType='NASAAmesExtract'
148                self.extractPrefix = '_naextract_'
149                                         
150        def openFile(self, filename):
151                ''' open file'''
152                self.file=nappy.openNAFile(filename)
153
154        def getListOfAxes(self):
155                ''' gets list of axes in file'''
156                axes=self.file.XNAME
157                #print 'before units stripped' + str(axes)
158                axes=self.stripunits(axes)
159                #print 'after units stripped' + str(axes)
160                return axes
161
162        def setAxis(self,axis):
163                '''sets selected axis in Data Interface'''
164                axes = self.getListOfAxes()
165                self.axisstub=axes.index(axis)
166
167        def getAxisAttribute(self, att):
168                ''' not implemented yet -should return axis attribute'''
169                #need to do something here...? maybe
170                attValue=None
171                return attValue
172
173        def getTimeUnits(self):
174                ''' returns units of time axis'''
175                axes = self.getListOfAxes()
176                for axis in axes:
177                        if string.find(string.upper(axis),'SECONDS SINCE') != -1:
178                                #found possible time axis.
179                                if axis[-3:]=='UTC':
180                                    units =string.lower(axis[:-4]) #hack!
181                                    units=units.replace('/','-') #try and clean it up
182                                else:
183                                    units=string.lower(axis)
184                                break
185                        elif string.find(string.upper(axis),'HOURS SINCE') != -1:
186                                #found possible time axis.
187                                units =(str(axis))
188                                break
189                        elif string.find(string.upper(axis),'DAYS SINCE') != -1:
190                                #found possible time axis.
191                                units =(str(axis))
192                                break
193                       
194                #revisit with udunits python library?
195                return units
196
197
198        def getDataForAxis(self):           
199            ''' gets data for an axis'''
200            if self.file.X == None:
201                    #print 'reading data....'
202                    self.file.readData()
203            for x in self.file.X:
204                pass
205                #print type(x)
206           
207           
208            d=Numeric.array(self.file.X)
209           
210            if type(self.file.X[1])==list:
211            #if len(self.file.X) > 0:
212                    data = self.file.X[self.axisstub]
213            else:
214                    data =self.file.X
215            #print data
216            data=Numeric.array(data)
217            return data
218
219        def getSizeOfAxis(self,axis):
220                ''' gets size of an axis'''
221                #check this function is okay
222                #is time always the first dimension in NA??
223                axes = self.getListOfAxes()
224                axisPosition=axes.index(axis)
225                #print "axis position" + str( axisPosition)
226                #print "NX" + str(self.file.NX)
227                try :
228                        axisSize=self.file.NX[axisPosition-1]
229                except:
230                        axisSize ='Unknown axis size'
231                return axisSize
232
233        def getListofVariables(self):
234                '''gets list of variables in the file'''
235                variableList=self.stripunits(self.file.VNAME)
236                return variableList
237
238        def setVariable(self,varname):
239                ''' sets the specified variable to be the selected variable'''
240                vlist=self.getListofVariables()
241                self.varstub=vlist.index(varname)
242
243        def getVariableAxes(self):
244                ''' Get axis names for a  variable '''
245                #hmm, now with Nasa Ames the Axis will be the same for all variables.
246                #so just use the getListOfAxes function again
247                #I think... check this!
248                varAxesList=self.getListOfAxes()
249                return varAxesList
250
251        def getVariableAttribute(self,attName):
252                ''' get an attribute for a variable'''
253                if attName =='units':
254                        #strip the units (attribute) from the variable
255                        unitslist=self.getUnits(self.file.VNAME)
256                        attribValue = unitslist[self.varstub]
257                        try:
258                                attribValue = unitslist[self.varstub]
259                        except:
260                                attribValue = 'unknown'
261                else:
262                        attribValue = 'unknown'
263                return attribValue
264
265        def getDataForVar(self):
266            ''' get data for a  variable'''
267            #NOTE TO SELF:
268            #Review this function (and in fact all of nasa ames data interface...)
269            if self.file.V == None:
270                    #print 'reading data....'
271                    self.file.readData()
272
273            try:
274                if type(self.file.V[1])==list:
275                    data = self.file.V[self.varstub]
276            #else:
277            #   data =self.file.X
278            #   print data
279                return data
280            except:
281                data = self.file.X
282                # print data
283                return data
284
285        def getArraySizeOfVar(self):
286            '''iterates through all dimensions in variable to get array size i.e a 3x3x3 array would have a size of 27'''
287
288            dimlist=self.file.NX
289            varsize =1
290            for item in dimlist:
291                    varsize = varsize * item
292                    #print "VARSISZE" + str(varsize)
293            return varsize
294
295        def getShapeOfVar(self):
296            '''this returns a list with the shape of the variable'''
297            varShape = []
298            for item in self.file.NX:
299                varShape.append(item)
300            return varShape
301
302        def getLowLimits(self):
303            ''' gets low limits of the data'''
304            lowlims = ""
305            for i in range (0, len(self.file.NX)):
306                    #for now, assume low limit is always of form 1 1 1 ..
307                    lowlims =lowlims + str(1)  +' '
308            return lowlims
309
310        def getHighLimits(self):
311            ''' gets high limits of the data'''
312            highlims = ""
313            for i in range (0, len(self.file.NX)):
314                    dimValue = self.file.NX[i]
315                    highlims =highlims  + str(dimValue) +' '
316            return highlims
317
318
319        def stripunits(self,listtostrip):
320            ''' strips units of measure from list
321            eg ['Universal time (hours)', 'Altitude (km)', 'Latitude (degrees)', 'Longitude (degrees)']
322            becomes ['Universal time', 'Altitude', 'Latitude', 'Longitude']'''
323            cleanlist = []
324            for item in listtostrip:
325                    openbracket=string.find(item,'(')
326                    if openbracket != -1:
327                            #if brackets exist, strip units.
328                            item=item[:openbracket-1]
329                    cleanlist.append(item)
330            return cleanlist
331
332        def getUnits(self,listwithunits):
333            '''gets units from list
334            eg ['Universal time (hours)', 'Altitude (km)', 'Latitude (degrees)', 'Longitude (degrees)']
335            becomes ['hours', 'km', 'degrees', 'degrees']'''
336            unitlist=[]
337            for item in listwithunits:
338                    openbracket=string.find(item,'(')
339                    item = item[openbracket+1:-1]
340                    unitlist.append(item)
341            return unitlist
342
343        def getTimes(self):
344            '''This function attempts to determine the time axis and read the time data
345            it may well not manage it.'''
346            axes = self.getListOfAxes()
347            for axis in axes:
348                if string.find(string.upper(axis),'TIME') != -1:
349                        #found possible time axis.
350                        self.setAxis(axis)
351                        times=self.getDataForAxis()
352                        break
353                elif string.find(string.upper(axis),'SECONDS SINCE') != -1:
354                        #found possible time axis.
355                        self.setAxis(axis)
356                        times=self.getDataForAxis()
357                        break
358                elif string.find(string.upper(axis),'HOURS SINCE') != -1:
359                        #found possible time axis.
360                        self.setAxis(axis)
361                        times=self.getDataForAxis()
362                        break
363                elif string.find(string.upper(axis),'DAYS SINCE') != -1:
364                        #found possible time axis.
365                        self.setAxis(axis)
366                        times=self.getDataForAxis()
367                        break
368            times=Numeric.array(times)
369            return times
370
371
372
373class cdunifInterface(AbstractDI):
374    '''Data Interface for cdunif (netcdf & pp formats & grib (not tested with grib)'''
375
376    def __init__(self):
377        '''these are just temporary values until we can determine whether the
378        file is netcdf pp or grib'''
379        self.extractType='cdunifExtract'
380        self.extractPrefix = '_cdunifextract_'
381   
382    def openFile(self, filename):
383        ''' open file'''
384        self._filename=filename
385        self.file=cdms.open(filename)
386       
387        #now we have the file name can properly determine extractType/Prefix
388        fileExtension = str(filename)[-3:]
389        if fileExtension == '.nc':
390            self.extractType = 'NetCDFExtract'
391            self.extractPrefix = '_ncextract_'
392        elif fileExtension == '.qxf':
393            self.extractType = 'NetCDFExtract'
394            self.extractPrefix = '_ncextract_'
395        elif fileExtension == '.pp':
396            self.extractType  = 'PPExtract'
397            self.extractPrefix = '_ppextract_'
398        elif fileExtension == 'ctl':
399            self.extractType = 'GRIBExtract'
400            self.extractPrefix = '_gribextract_'
401        elif fileExtension == 'xml': 
402            self.extractType = 'NetCDFExtract'  #okay this isn't true, but ok for testing
403            self.extractPrefix = '_ncextract__' 
404   
405    def getListOfAxes(self):
406        ''' returns list of dimensions/axes'''
407        axes=self.file.dimensions.keys()
408        return axes
409
410    def getSizeOfAxis(self,axis):
411        ''' returns size of axis'''
412        axisSize=self.file.dimensions[axis]
413        return axisSize
414
415    def getListofVariables(self):
416        ''' returns list of variables'''
417        variableList=self.file.variables.keys()
418
419        # Hack to test if removing climatology_bounds fixes pywms bug
420        if 'climatology_bounds' in variableList:
421            variableList.remove('climatology_bounds')
422
423        return variableList
424
425    def setAxis(self,axis):
426        ''' set axis'''
427        if not hasattr(self, 'file'):
428            raise Exception, 'Could not open/find underlying file: %s  -  If you are the system maintainer check file paths and permissions'%self._filename
429        self.axisobj=self.file.getAxis(axis)
430
431           
432
433    def getAxisAttribute(self, att):
434        ''' get attribute of set axis'''
435        attValue=self.axisobj.attributes[att]
436        return attValue
437   
438    def getTimeUnits(self):
439        ''' get time units - only works if set axis is time...'''
440        #this does the same as getAxisAttribute, but is a separate function as different formats handle time differently.
441        return self.getAxisAttribute('units')
442
443    def getDataForAxis(self):
444        ''' get data for entire axis'''
445        data = self.axisobj.getValue()
446        return data
447
448    def setVariable(self,varname):
449        ''' sets variable '''
450        self.varobj=self.file.variables[varname]
451
452    def getVariableAxes(self):
453        ''' gets axes for set variable'''
454        varAxesList=self.varobj.getAxisIds()
455        return varAxesList
456
457    def getVariableAttribute(self,attName):
458        ''' gets requested attribute for set variable'''
459        if attName == '_FillValue':
460            try:
461                attribValue=self.varobj._FillValue
462                attribValue=attribValue.item()
463            except:
464                try:
465                    attribValue=self.varobj.missing_value
466                    attribValue=attribValue.item()
467                except:
468                    attribValue = None 
469        else:
470            attribValue = self.varobj.getattribute(attName)
471        return attribValue
472
473    def getDataForVar(self):
474        ''' get data for set variable'''
475        data = self.varobj.getValue()
476        return data
477   
478    def _fixLongitudeRequest(self, **kwargs):
479        ''' fixes problem of longitude requests taking different forms'''
480        lonkey='longitude'
481        if lonkey in kwargs.keys():   #this test needs to be much more robust...!
482            if  type(kwargs[lonkey]) is tuple:
483                newkws=[]
484                for val in kwargs[lonkey]:
485                    newval=((val + 180) % 360) - 180
486                    newkws.append(newval)
487               
488                 #lon = ((t_lon + 180) % 360) - 180
489                #lonmin=kwargs[lonkey][0] - (kwargs[lonkey][0]/180)*360.0
490                #lonmax=kwargs[lonkey][1] - (kwargs[lonkey][1]/180)*360.0
491               
492                #if newkws[0] > newkws[1]:
493                    #flip=[newkws[1], newkws[0]]
494                    #newkws=flip
495                kwargs[lonkey]=tuple(newkws)
496        return kwargs
497               
498    def _componentTimes(self, times):
499        ''' handles time formats'''
500        if type(times) is tuple:
501            comptimes=[]
502            for t in times:
503                compt=csml.csmllibs.csmltime.getCDtime(t)
504                comptimes.append(compt)
505            result=tuple(comptimes)
506        else: #just a single time
507            result=csml.csmllibs.csmltime.getCDtime(times)
508        return result
509       
510   
511    def getSubsetOfDataForVar(self, **kwargs):     
512        ''' Main subsetting method. Kwargs contains dictionary of requested subset.'''
513        #put any slicing indices aside for later and use names       
514        try:
515            upper=kwargs['upper']
516            del kwargs['upper']
517        except:
518            pass
519        try:
520            lower=kwargs['lower']
521            del kwargs['lower']
522        except:
523            pass
524        if 'time' in kwargs:
525            kwargs['time']=self._componentTimes(kwargs['time'])
526       
527        #kwargs=self._fixLongitudeRequest(**kwargs)
528        #return self.file(self.varobj.id,**kwargs)
529        try:
530            #takes keyword args defining subset eg
531            #subset=getSubsetOfDataForVar(latitude=(0.,10.0), longitude=(90, 100.0))
532            subset=None
533            lonkey='longitude'
534            if lonkey in kwargs.keys():   #this test needs to be much more robust...!
535                if  type(kwargs[lonkey]) is tuple:
536                    if kwargs[lonkey][0] > kwargs[lonkey][1]:
537                        #subsetting greenwich meridian around 0
538                        lonMin = kwargs[lonkey][0]
539                        lonMax =kwargs[lonkey][1]
540                        kwargs[lonkey]=(0.1, lonMax)
541                        subset1=self.file(self.varobj.id,**kwargs)
542                        kwargs[lonkey]=(lonMin,360)
543                        try:
544                            subset2=self.file(self.varobj.id,**kwargs)
545                            longitudeAxis=subset1.getAxisIndex(lonkey)
546                            #concatenate arrays along longitude             
547                            subset = cdms.MV2.concatenate([subset2,subset1],axis=longitudeAxis)
548                        except:
549                            subset=subset1
550            if type(subset) is not cdms.tvariable.TransientVariable:
551                subset=self.file(self.varobj.id,**kwargs)
552        except:             
553                return None
554                #try and slice with indices instead
555        return subset
556   
557   
558    def getArraySizeOfVar(self):
559        '''iterates through all dimensions in variable to get array size i.e a 3x3x3 array would have a size of 27'''
560        var = self.varobj
561        size = var.shape
562        varsize = 1
563        for item in size:
564            varsize = item *varsize
565        return varsize
566
567    def getShapeOfVar(self):
568        ''' get shape of set variable '''
569        varShape = []
570        for item in self.varobj.shape:
571            varShape.append(item)
572        return varShape
573
574    def getLowLimits(self):
575        ''' get gml low limits of set variable'''
576        dimNames = self.varobj.getAxisIds()
577        lowlims = ""
578        for i in range (1, len(dimNames)):
579            #for now, assume low limit is always of form 1 1 1 ..
580            lowlims =lowlims + str(1)  +' '
581        return lowlims
582
583    def getHighLimits(self):
584        ''' get gml high limits of set variable'''
585        dimNames = self.varobj.getAxisIds()
586        highlims = ""
587        for i in range (1, len(dimNames)):
588            dimValue = self.file.dimensions[dimNames[i]]
589            highlims =highlims  + str(dimValue) +' '
590        return highlims
591       
592   
593class cdmlInterface(cdunifInterface):
594    ''' this is more  or less the cdunif interface but a few methods have been overwritten to read CDML'''
595    def __init__(self):
596        #this all needs to be revisited in csml v2.
597        self.extractType='cdmlExtract'
598        self.extractPrefix = '_cdmlextract_'
599       
600    def getListOfAxes(self):
601        ''' get list of axis in file'''
602        axes=self.file.axes.keys() 
603        return axes
604
605    def getSizeOfAxis(self,axis):
606        ''' get size of set axis'''
607        axisSize=self.file.axes[axis].length
608        return axisSize
609
610
611class RawFileInterface(AbstractDI):
612    ''' interface for raw file types'''
613    def __init__(self):
614        self.extractType   = 'rawExtract'
615        self.extractPrefix = '_rawextract_'
616 
617 
618    def openFile(self, filename):
619        self.file = open(filename, "rb")
620
621    def closeFile(self):
622            self.file.close()
623
624
625    # Read the data from the raw file into a multidimensional Numeric array
626    def readFile(self, **kwargs):
627        # Determine the numeric type:
628        if 'numericType' in kwargs:
629            try:
630                numericTypeCode = {
631                    'uint8':Numeric.UInt8,
632                    'uint16':Numeric.UInt16,
633                    'uint32':Numeric.UInt32,
634                    'int8':Numeric.Int8,
635                    'int16':Numeric.Int16,
636                    'int32':Numeric.Int32,
637                    'float':Numeric.Float32,
638                    'double':Numeric.Float64
639                }[kwargs['numericType']]
640            except KeyError:
641                raise TypeError("Invalid numericType: " + str(kwargs['numericType']))
642        else:
643            raise KeyError("numericType is mandatory.")
644       
645        # Read the file into a numpy array:
646        self.data = Numeric.fromstring(self.file.read(), numericTypeCode)
647        # If necessary, swap the byte order:
648        if 'endianess' in kwargs:
649            if ((kwargs['endianness'] == 'little' and not Numeric.LittleEndian) or (kwargs['endianness'] == 'big' and Numeric.LittleEndian)):
650                self.data = self.data.byteswap()   
651        # Declare the shape of the array:
652        dimensions = map(int,kwargs['dimensions'])
653        dimensions.reverse()
654        self.data.shape = tuple(dimensions)
655        # If numericTransform or fillValue were provided, store them as
656        # attributes.
657        if 'numericTransform' in kwargs:
658            self.numericTransform = NumericTransform.infixExpression(kwargs['numericTransform'])
659        if 'fillValue' in kwargs:
660            self.fillValue = kwargs['fillValue']
661   
662   # Return the fill value, if set, and transform if necessary:
663    def getFillValue(self):
664      # Both fillValue and numericTransform attributes may or may
665      # not exist...
666        try:
667            return self.numericTransform.solve( n = float(self.fillValue) )
668        except AttributeError:
669            try:
670                return self.fillValue
671            except AttributeError:
672                return None
673
674
675   # This is a just a special case of getSubsetOfDataForVar. To avoid
676   # duplication of code, just subset the entire array. (getSubset.. is
677   # optimised for this case)
678    def getDataForVar(self):
679        return self.getSubsetOfDataForVar(lower = (0,)*len(self.data.shape),
680                                        upper = self.data.shape)
681
682
683   # Subset the n-dimensional file based on array indices. Accepts parameters
684   # in the form of e.g.
685   #
686   # getSubsetOfDataForVar( lower=(0,0), upper=(1,2) )
687   #
688   # Note: The rank of the required subset, must exactly match the
689   # rank of the original data: len(lower) == len(upper) == rank of file
690   #
691    def getSubsetOfDataForVar(self, **kwargs):
692        # Assume subset parameters are passed as: lower=(0,0) upper=(512,512)
693        if 'lower' not in kwargs or 'upper' not in kwargs:
694            # Have not specified any subset parameters that we recognise, so raise
695            # an exception:
696            raise NotImplementedError("Only supports subsetting with lower+upper array indices")
697        elif not len(kwargs['lower']) == len(kwargs['upper']) == len(self.data.shape):
698            # Rank of subset is not the same as rank of full data array. so raise
699            # an exception:
700            raise NotImplementedError("Only supports subsets of same rank as full dataset")
701        elif Numeric.sometrue(Numeric.greater(kwargs['upper'], self.data.shape)):
702            # Requested upper bound of subset is beyond the size of the the full
703            # data array, so raise an exception
704            raise IndexError("Subset out of range")
705        elif Numeric.sometrue(Numeric.less( kwargs['upper'], Numeric.zeros(len(self.data.shape)))):
706            # Requested lower bound of subset is beyond the size of the the full
707            # data array, so raise an exception
708            raise IndexError("Subset out of range")
709        elif Numeric.sometrue(Numeric.less_equal(kwargs['upper'], kwargs['lower'])):
710            # lower bound <= upper_bound for at least one dimension, so raise an
711            # exception
712            raise IndexError("Upper bound less than lower bound")
713        elif tuple(kwargs['lower']) == (0,)*len(self.data.shape) and tuple(kwargs['upper']) == self.data.shape:
714            # Special case of requested subset == entire data file.
715            subset = self.data
716        else:
717            # We are okay to subset.
718   
719            # I cant see any nice (and speedy) way of subsetting a Numeric
720            # array of arbitrary rank without the use of eval. By doing it
721            # this way, we can make use of the (possible) acceleration provided
722            # by Numeric/NumPy.
723            slices = []
724            for i in range(len(self.data.shape)):
725                lower = int(kwargs['lower'][i])
726                upper = int(kwargs['upper'][i]) 
727                slices.append(str(lower)+':'+str(upper))
728            subset = eval("self.data["+','.join(slices)+"]")
729   
730        # Attempt to perform the numericTransform on the data array, if we get
731        # AttributeError, it most likely means that numericTransform was not
732        # specified, so return the data as-is
733        try:
734            return self.numericTransform.solve( n = subset )
735        except AttributeError:
736            return subset.copy()
737
738class ImageFileInterface(RawFileInterface):
739    '''Interface for reading data from image files. Requires PIL Image module.'''
740    def __init__(self):
741        self.extractType   = 'imageExtract'
742        self.extractPrefix = '_imageextract_'
743   
744    def image2array(self,im):
745        #Adapted from code by Fredrik Lundh, http://www.pythonware.com
746        #http://effbot.org/zone/pil-numpy.htm
747            if im.mode not in ("L", "F"):
748                raise ValueError, "can only convert single-layer images"
749            if im.mode == "L":
750                a = Numeric.fromstring(im.tostring(), Numeric.UnsignedInt8)
751            else:
752                a = Numeric.fromstring(im.tostring(), Numeric.Float32)
753            a.shape = im.size[1], im.size[0]
754            return a
755
756    def openFile(self, filename):
757        self.file = Image.open(filename)
758
759    def closeFile(self):
760        self.file = None #...Image does not seem to have a close() method.
761
762    def readFile(self, **kwargs):
763        # Convert the image to a Numeric array
764     
765        self.data=self.image2array(self.file)
766        #slower method:
767        #self.data = Numeric.array(self.file.getdata())
768
769        if 'numericTransform' in kwargs:
770            # numericTransform was specified, so compile the expression:
771            self.numericTransform = NumericTransform.infixExpression(kwargs['numericTransform'])
772        if 'fillValue' in kwargs:
773            self.fillValue = kwargs['fillValue']
774
775class HDF5Interface(AbstractDI):
776    '''Data Interface for HDF5'''
777    def __init__(self):
778        self.extractType='hdfExtract'
779        self.extractPrefix='_hdfextract_'
780   
781    def openFile(self, filename):
782        '''some code to open the file'''
783        pass
784
785    def setAxis(self,axis):
786        '''some code to set an axis to be queried, may not need to do much, depending on your format'''
787        pass
788       
789    def getDataForAxis(self):
790        '''some code to return the values for an axis'''
791        return data
792
793    def setVariable(self,varname):
794        '''some code to set a variable to be queried, may not need to do much, depending on your format'''
795        pass
796
797    def getDataForVar(self):
798        '''some code to return all values for a variable'''
799        return data
800
801    def getSubsetOfDataForVar(self, **kwargs):
802        '''takes keyword args defining subset eg
803        subset=getSubsetOfDataForVar(latitude=(0.,10.0), longitude=(90, 100.0), ...)
804        and returns a subset of data for tha variable '''
805        return data
806
807    def closeFile(self):
808        '''some code to close the file'''
809        pass
Note: See TracBrowser for help on using the repository browser.