source: nappy/tags/nappy_pre_refactor_feb2008/nappy/localRules/aircraftData.py @ 3313

Subversion URL: http://proj.badc.rl.ac.uk/svn/ndg/nappy/tags/nappy_pre_refactor_feb2008/nappy/localRules/aircraftData.py@3313
Revision 3313, 8.4 KB checked in by astephen, 12 years ago (diff)

Moved to new version as tagged backup.

RevLine 
[3313]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"""
6airCraftData.py
7===============
8
9Container module for AircraftData class that deals with
10flight data which has unusual time axes and flags in separate
11variables. This class also does sub-sampling by averaging to
121Hz frequency.
13
14
15"""
16
17# Imports from python standard library
18import os,sys
19import MA, Numeric, cdms, MV, cdtime
20
21# Imports from local package
22
23
24class AircraftData:
25    """
26    AircraftData class that deals with
27    flight data which has unusual time axes and flags in separate
28    variables. This class also does sub-sampling by averaging to
29    1Hz frequency.
30    """
31
32    def __init__(self, var, timeVar, flagVar=None, args=["av"]):
33        """
34        Parses args and calls required methods.
35        """
36        if args[0]=="av":
37            samplingMethod="averaging"
38        elif args[0]=="ss":
39            samplingMethod="selection"
40
41        flagsToUse=[0,1]
42        samplingRate=1
43        if len(args)>1:
44            samplingRate=int(args[1])
45
46            if len(args)>2:
47                flagsToUse=[int(flag[0]) for flag in args[2:]]
48
49        if samplingMethod=="averaging":
50            self.subSampledVar=self._subSampleByAveraging(var, timeVar, flagVar,
51                    samplingRate, flagsToUse)
52        elif samplingMethod=="selection":
53            self.subSampledVar=self._subSampleBySelection(var, timeVar, flagVar,
54                    samplingRate, flagsToUse)
55
56
57    def _subSampleBySelection(self, var, timeVar, flagVar, samplingRate, flagsToUse):
58        """
59        Returns a new variable which is 'var' sub-sampled by selecting the value
60        at the given samplingRate with data selected according to flagVar
61        including the flag values  specified in flagsToUse (defaults are 0 and 1).
62        """
63        if samplingRate!=1:
64            raise "Sampling rates of values other than 1Hz are not yet supported."
65        maskedArray=self._getMaskedArray(var[:,0](squeeze=1), flagVar[:,0](squeeze=1), flagsToUse)
66        # Now re-construct variable axes etc
67        newTimeAxis=self._flatten2DTimeAxis(timeVar, samplingRate)
68        newVar=self._recreateVariable(var, maskedArray, newTimeAxis, flagVar, max(flagsToUse), missingValue=maskedArray.fill_value(), sampleBy="selection")
69        return newVar
70
71       
72    def _subSampleByAveraging(self, var, timeVar, flagVar, samplingRate, flagsToUse):
73        """
74        Returns a new variable which is 'var' sub-sampled by averaging
75        at the given samplingRate with data selected according to flagVar
76        including the flag values specified in flagsToUse (defaults are 0 and 1).
77        """
78        maskedArray=self._getMaskedArray(var, flagVar, flagsToUse)
79        shape=var.shape
80        if shape[1]==1:
81            newNumArray=MV.ravel(maskedArray)
82            newArrayMask=MV.ravel(maskedArray.mask())
83            newArray=MA.masked_array(newNumArray, mask=newArrayMask, fill_value=maskedArray.fill_value())
84        else:
85            newArray=Numeric.zeros(shape[0], 'f')
86            for t0 in range(shape[0]):
87                # Set as missing if less than half are valid                 
88                t1Array=maskedArray[t0]
89
90                if samplingRate==1:
91                    # If half or more are good values then calculate the mean
92                    if t1Array.count()>=(shape[1]/2.):
93                        newArray[t0]=MA.average(t1Array)
94                    # otherwise set as missing value
95                    else:
96                        newArray[t0]=maskedArray.fill_value()
97       
98                else:
99                    raise "Averaging for non 1Hz sampling rates not yet supported!"
100
101        # Now re-construct variable axes etc
102        newTimeAxis=self._flatten2DTimeAxis(timeVar, samplingRate)
103        newVar=self._recreateVariable(var, newArray, newTimeAxis, flagVar, max(flagsToUse), missingValue=maskedArray.fill_value(), sampleBy="averaging") 
104        return newVar
105
106
107    def _flatten2DTimeAxis(self, timevar, samplingRate):
108        """
109        Returns a flattened 2D time axis.
110        """
111        if samplingRate!=1:
112            raise "Cannot yet deal with sub-sampling to non 1Hz sampling rates!"
113
114        timevalues=timevar._data
115        timeunit=timevar.units
116        newTimeValues=[]
117        rate=samplingRate
118
119        for t in timevalues:
120            for i in range(rate):
121                tvalue=t+((1./rate)*i)
122                newTimeValues.append(tvalue)
123
124        newTimeAx=cdms.createAxis(newTimeValues)
125        newTimeAx.units=timeunit
126        newTimeAx.designateTime()
127        newTimeAx.id=newTimeAx.long_name=newTimeAx.standard_name="time"
128        return newTimeAx
129
130
131    def _recreateVariable(self, var, newArray, timeAxis, flagVar, flagMax, missingValue, sampleBy="averaging"):
132        """
133        Rebuilds a variable using new array and new time axis.
134        """
135        atts=var.attributes
136        newAtts={}
137        for att,value in atts.items():
138            if type(value) in (type((1,1)), type([1,2]), type(Numeric.array([1.]))) and len(value)==1:
139                value=value[0]
140            if type(value) in (type(1), type(1.), type(long(1))):
141                newAtts[att]=value
142            elif type(value)==type(""):
143                newAtts[att]=value.strip()
144
145        # Now create the variable
146        missing=missingValue
147        data=newArray
148        newvar=cdms.createVariable(data, id=var.id, fill_value=missing, axes=[timeAxis], attributes=newAtts)
149        newvar.frequency=1
150        if sampleBy=="averaging":
151            newvar.comment="Data averaged to give a pseudo sampling frequency of 1 Hz. "
152            if flagVar:
153                newvar.comment=newvar.comment+"Values were included where the flag variable was less than %s and where both the flag value and actual data value were present (i.e. not missing values)." % (flagMax+1)
154                newvar.comment=newvar.comment+" Average data values were set to missing where more than half the values in any one second period were already missing. Missing value was also used where more than half the flags were either missing value or values of 2 or greater."
155            else:
156                newvar.comment=newvar.comment+"The flag variable was not used in creating this variable."
157                newvar.comment=newvar.comment+" Average data values were set to missing where more than half the values in any one second period were already missing."
158        elif sampleBy=="selection":
159            newvar.comment="Data sub-sampled to give a pseudo sampling frequency of 1 Hz by selecting only the per-second frequency values. "
160            if flagVar:
161                newvar.comment=newvar.comment+" Values were included where the flag variable was less than 2 and where both the flag value and actual data value were present (i.e. not missing values)."
162            else:
163                newvar.comment=newvar.comment+" The flag variable was not used in creating this variable." 
164        return newvar
165 
166
167    def _getMaskedArray(self, var, flagVar=None, flagValues=(0,1), missingValuesToTest=(-9999., -32767.)):
168        """
169        Returns a masked array that has the values only present where the flag
170        values are acceptable and the data itself is not missing.
171        """
172        if flagVar:
173            flagFillValue=-1
174            flagGTValue=max(flagValues)
175
176            flagMaskWhereUnacceptable=MA.masked_greater(flagVar, flagGTValue).mask()
177            flagMaskWhereMissing=MA.masked_equal(flagVar, flagFillValue).mask()
178            flagMask=flagMaskWhereUnacceptable+flagMaskWhereMissing
179
180        for fv in missingValuesToTest:
181            if fv in MV.ravel(var):
182                print "Setting missing value for '%s' as: %s" % (var.id, fv)
183                varFillValue=fv
184        else:
185            varFillValue=missingValuesToTest[0]
186
187        if flagVar:
188            varMask=MA.masked_array(var, mask=MA.equal(var, varFillValue), fill_value=varFillValue).mask()
189            fullmask=MA.bitwise_or(varMask, flagMask)
190            maskedArray=MA.masked_array(var, mask=fullmask, fill_value=varFillValue)
191        else:
192            maskedArray=MA.masked_array(var, mask=MA.equal(var, varFillValue), fill_value=varFillValue)
193
194        #maskedArray=cdms.createVariable(maskedArray, id=var.id, fill_value=varFillValue)
195        #maskedArray.missing_value=varFillValue
196        return maskedArray
197
Note: See TracBrowser for help on using the repository browser.