Changeset 5434 for qesdi/data_processing


Ignore:
Timestamp:
30/06/09 16:03:18 (10 years ago)
Author:
dgrant
Message:

Can now read in precip data and write functionaly correct netCDF. Still
need to make this conform to CF

Location:
qesdi/data_processing
Files:
1 added
3 edited

Legend:

Unmodified
Added
Removed
  • qesdi/data_processing/COORDDataTable.py

    r5419 r5434  
     1import numpy 
    12 
    23class COORDDataTable: 
    34    """Class to store a lookup tree (dict of dicts) of data by latlon""" 
    45     
    5     def __init__(self): 
     6    def __init__(self,emptyVal): 
    67        self.lookupTable={} 
    7          
     8        self.emptyVal = emptyVal 
     9        self.prepopulateTable() 
     10 
     11    def prepopulateTable(self): 
     12        f32 = numpy.float32 
     13        increments = [f32(0.083), f32(0.25), f32(0.417), f32(0.583), f32(0.75), f32(0.917)] 
     14        self.lats = [] 
     15        increments.reverse() 
     16        for wholeLat in numpy.arange(-89,1,1,f32): 
     17            self.lats.extend(["%.3f" % (wholeLat - inc) for inc in increments]) 
     18        increments.reverse() 
     19        for wholeLat in numpy.arange(0,90,1,f32): 
     20            self.lats.extend(["%.3f" % (wholeLat + inc) for inc in increments]) 
     21        self.lons = [] 
     22        increments.reverse() 
     23        for wholeLon in numpy.arange(-179,1): 
     24            self.lons.extend(["%.3f" % (wholeLon - inc) for inc in increments]) 
     25        increments.reverse() 
     26        for wholeLon in numpy.arange(0,180): 
     27            self.lons.extend(["%.3f" % (wholeLon + inc) for inc in increments]) 
     28        for lat in self.lats: 
     29            for lon in self.lons: 
     30                self.lookupTable.setdefault(lat,dict())[lon]=self.emptyVal 
     31 
     32 
    833    def insert(self,dataPt, lat, lon): 
    9         self.lookupTable.setdefault(lat,dict())[lon]=dataPt 
    10          
     34        self.lookupTable[lat][lon]=dataPt 
     35 
    1136    def lookup(self, lat, lon): 
    12             return self.lookupTable[lat][lon] 
     37        return self.lookupTable[lat][lon] 
    1338 
    1439    def getLats(self): 
    15         return sorted(self.lookupTable.keys()) 
     40        return self.lats 
    1641 
    1742    def getLons(self): 
    18         return sorted(self.lookupTable.values()[0].keys()) 
     43        return self.lons 
    1944 
    2045    def getLatSize(self): 
    21         return len(self.lookupTable) 
     46        return len(self.lats) 
    2247 
    2348    def getLonSize(self): 
    24         return len(self.lookupTable.values()[0]) 
     49        return len(self.lons) 
     50 
  • qesdi/data_processing/TMCReader.py

    r5419 r5434  
    22from COORDDataTable import COORDDataTable 
    33 
    4 class TMCReader: 
     4class TMCReader(object): 
    55    """A class to read the ascii text from the precipetation data file""" 
    66     
    7     def __init__(self): 
    8         self.coordTable = COORDDataTable()        
     7    def __init__(self,emptyVal): 
     8        self.coordTable = COORDDataTable(emptyVal)        
    99 
    1010    def readData(self, path): 
  • qesdi/data_processing/convertPre.py

    r5419 r5434  
    55 
    66from TMCReader import TMCReader 
     7from TMCWriter import TMCWriter 
    78from netCDF4 import Dataset 
     9 
     10 
     11f32 = numpy.float32 
    812 
    913class PrecipetationReader(TMCReader): 
    1014    """A class to read the ascii text from the precipetation data file""" 
     15 
     16 
     17    def __init__(self):  
     18        emptyVal = PrecipetationDataPt([-1.e20,numpy.float32(-1.e20),numpy.float32(-1.e20),-1.e20,numpy.float32(-1.e20),numpy.float32(-1.e20),numpy.float32(-1.e20),numpy.float32(-1.e20),numpy.float32(-1.e20),numpy.float32(-1.e20),numpy.float32(-1.e20),numpy.float32(-1.e20),numpy.float32(-1.e20),numpy.float32(-1.e20),numpy.float32(-1.e20),numpy.float32(-1.e20),numpy.float32(-1.e20),numpy.float32(-1.e20),numpy.float32(-1.e20),numpy.float32(-1.e20),numpy.float32(-1.e20),numpy.float32(-1.e20),numpy.float32(-1.e20),numpy.float32(-1.e20),numpy.float32(-1.e20),numpy.float32(-1.e20)]) 
     19        super(PrecipetationReader,self).__init__(emptyVal) 
    1120 
    1221    def storeLine(self, line): 
     
    1423        self.coordTable.insert(precipDataPt, precipDataPt.lat, precipDataPt.lon) 
    1524 
    16  
    17  
    1825class PrecipetationDataPt: 
    1926    """Class to store data for a specific point""" 
    2027 
    2128    def __init__(self,preVals): 
    22         self.lat=float(preVals[0]) 
    23         self.lon=float(preVals[1]) 
    24         self.jan=float(preVals[2]) 
    25         self.feb=float(preVals[3]) 
    26         self.mar=float(preVals[4]) 
    27         self.apr=float(preVals[5]) 
    28         self.may=float(preVals[6]) 
    29         self.jun=float(preVals[7]) 
    30         self.jul=float(preVals[8]) 
    31         self.aug=float(preVals[9]) 
    32         self.sep=float(preVals[10]) 
    33         self.oct=float(preVals[11]) 
    34         self.nov=float(preVals[12]) 
    35         self.dec=float(preVals[13]) 
     29        self.lat=preVals[0] 
     30        self.lon=preVals[1] 
     31        self.jan=f32(preVals[2]) 
     32        self.feb=f32(preVals[3]) 
     33        self.mar=f32(preVals[4]) 
     34        self.apr=f32(preVals[5]) 
     35        self.may=f32(preVals[6]) 
     36        self.jun=f32(preVals[7]) 
     37        self.jul=f32(preVals[8]) 
     38        self.aug=f32(preVals[9]) 
     39        self.sep=f32(preVals[10]) 
     40        self.oct=f32(preVals[11]) 
     41        self.nov=f32(preVals[12]) 
     42        self.dec=f32(preVals[13]) 
     43        self.jancv=f32(preVals[14]) 
     44        self.febcv=f32(preVals[15]) 
     45        self.marcv=f32(preVals[16]) 
     46        self.aprcv=f32(preVals[17]) 
     47        self.maycv=f32(preVals[18]) 
     48        self.juncv=f32(preVals[19]) 
     49        self.julcv=f32(preVals[20]) 
     50        self.augcv=f32(preVals[21]) 
     51        self.sepcv=f32(preVals[22]) 
     52        self.octcv=f32(preVals[23]) 
     53        self.novcv=f32(preVals[24]) 
     54        self.deccv=f32(preVals[25]) 
    3655 
    37 class PrecipetationWriter: 
     56 
     57class PrecipetationWriter(TMCWriter): 
    3858    """Class to read precipetation data""" 
    3959 
    40     def __init__(self,lookupTable): 
    41         self.dataTable = lookupTable 
    42  
    4360    def writeData(self, path): 
    44         self.rootgrp = Dataset(path, 'w',format='NETCDF4') 
     61        self.rootgrp = Dataset(path, 'w',format='NETCDF3_CLASSIC') 
    4562        self.createDimensions() 
    4663        self.createVariables() 
    47         preVar = self.rootgrp.createVariable('Precipetation','f',('time','lat','lon')) 
     64        self.preVar = self.rootgrp.createVariable('precipetation_means','f4',('time','lat','lon',)) 
     65        self.cvVar = self.rootgrp.createVariable('precipetation_cvs','f4',('time','lat','lon',)) 
    4866        self.setGlobalAttributes() 
    4967        self.setVariableAttributes() 
    50         self.writeVariables() 
     68        self.preVar._fillValue = numpy.float32(-1.e20) 
     69        self.preVar.missing_value = numpy.float32(-1.e20) 
     70        self.preVar.unit = 'mm/month' 
     71        self.cvVar._fillValue = numpy.float32(-1.e20) 
     72        self.cvVar.missing_value = numpy.float32(-1.e20) 
     73        self.cvVar.unit = 'percent' 
     74        self.writeVariables()        
     75        self.preVar[:] = self.getPrecipetationData() 
     76        self.cvVar[:] = self.getPrecipetationCVsData() 
    5177        self.rootgrp.close() 
    5278 
    53     def createDimensions(self): 
    54         self.rootgrp.createDimension('time', 12) 
    55         self.rootgrp.createDimension('nv', 2) 
    56         self.rootgrp.createDimension('n',self.lookupTable.getLength) 
     79    def getPrecipetationData(self): 
     80        print 'Writing precipetation data' 
     81        self.precips = numpy.zeros((12,self.dataTable.getLatSize(),self.dataTable.getLonSize())) 
     82        i = 0 
     83        j = 0 
     84        print self.precips.shape 
     85        for lt in self.dataTable.getLats(): 
     86            j = 0 
     87            for ln in self.dataTable.getLons(): 
     88                precipetationDataPt = self.dataTable.lookup(lt,ln) 
     89                self.precips[0][i][j] = precipetationDataPt.jan 
     90                self.precips[1][i][j] = precipetationDataPt.feb 
     91                self.precips[2][i][j] = precipetationDataPt.mar 
     92                self.precips[3][i][j] = precipetationDataPt.apr 
     93                self.precips[4][i][j] = precipetationDataPt.may 
     94                self.precips[5][i][j] = precipetationDataPt.jun 
     95                self.precips[6][i][j] = precipetationDataPt.jul 
     96                self.precips[7][i][j] = precipetationDataPt.aug 
     97                self.precips[8][i][j] = precipetationDataPt.sep 
     98                self.precips[9][i][j] = precipetationDataPt.oct 
     99                self.precips[10][i][j] = precipetationDataPt.nov 
     100                self.precips[11][i][j] = precipetationDataPt.dec 
     101                j = j + 1 
     102            i = i + 1 
     103        return self.precips 
     104     
     105    def getPrecipetationCVsData(self): 
     106        print 'Writing precipetation data' 
     107        self.precips = numpy.zeros((12,self.dataTable.getLatSize(),self.dataTable.getLonSize())) 
     108        i = 0 
     109        j = 0 
     110        print self.precips.shape 
     111        for lt in self.dataTable.getLats(): 
     112            j = 0 
     113            for ln in self.dataTable.getLons(): 
     114                precipetationDataPt = self.dataTable.lookup(lt,ln) 
     115                self.precips[0][i][j] = precipetationDataPt.jancv 
     116                self.precips[1][i][j] = precipetationDataPt.febcv 
     117                self.precips[2][i][j] = precipetationDataPt.marcv 
     118                self.precips[3][i][j] = precipetationDataPt.aprcv 
     119                self.precips[4][i][j] = precipetationDataPt.maycv 
     120                self.precips[5][i][j] = precipetationDataPt.juncv 
     121                self.precips[6][i][j] = precipetationDataPt.julcv 
     122                self.precips[7][i][j] = precipetationDataPt.augcv 
     123                self.precips[8][i][j] = precipetationDataPt.sepcv 
     124                self.precips[9][i][j] = precipetationDataPt.octcv 
     125                self.precips[10][i][j] = precipetationDataPt.novcv 
     126                self.precips[11][i][j] = precipetationDataPt.deccv 
     127                j = j + 1 
     128            i = i + 1 
     129        return self.precips 
    57130 
    58     def createVariables(self): 
    59         self.timeVar = self.rootgrp.createVariable('time','d', ('time',)) 
    60         self.latVar = self.rootgrp.createVariable('lat', 'f4', ('lat',)) 
    61         self.lonVar = self.rootgrp.createVariable('lon', 'f4', ('lon',)) 
    62         self.boundsVar = self.rootgrp.createVariable('climatology_bounds','d',('time','nv',)) 
    63         self.preVar = self.rootgrp.createVariable('precipetation','f4',('time','lat','lon',)) 
     131def testReader(): 
     132    precipetationReader = PrecipetationReader() 
     133    testTable = precipetationReader.getLookupTable() 
     134    print testTable.getLats() 
     135    precipetationReader.storeLine('10.083 10.083 1 2 3 4 5 6 7 8 9 10 11 12') 
     136    precipetationReader.storeLine('11.083 10.083 1 2 3 4 5 6 7 8 9 10 11 12') 
     137    precipetationReader.storeLine('12.083 10.083 1 2 3 4 5 6 7 8 9 10 11 12') 
     138    precipetationReader.storeLine('13.083 10.083 1 2 3 4 5 6 7 8 9 10 11 12') 
     139    precipetationReader.storeLine('14.083 10.083 1 2 3 4 5 6 7 8 9 10 11 12') 
     140    precipetationReader.storeLine('15.083 10.083 1 2 3 4 5 6 7 8 9 10 11 12') 
     141    precipetationReader.storeLine('16.083 10.083 1 2 3 4 5 6 7 8 9 10 11 12') 
     142    precipetationReader.storeLine('17.083 10.083 1 2 3 4 5 6 7 8 9 10 11 12') 
     143    precipetationReader.storeLine('18.083 10.083 1 2 3 4 5 6 7 8 9 10 11 12') 
     144    precipetationReader.storeLine('19.083 10.083 1 2 3 4 5 6 7 8 9 10 11 12') 
     145    print testTable.lookup(f32(11.083), f32(10.083)).jan 
     146    print testTable.lookup(f32(11.083), f32(11.083)).jan 
    64147 
    65     def setGlobalAttributes(self): 
    66         self.rootgrp.history = 'Created ' + time.ctime(time.time()) 
     148def testWriter(): 
     149    table = COORDDataTable(10)  
    67150 
    68     def setVariableAttributes(self): 
    69         self.latVar.units = 'degrees north' 
    70         self.lonVar.units ='degrees east' 
    71         self.timeVar.units = 'days since 1961-1-1' 
    72         self.timeVar.climatology = 'climatology bounds' 
    73  
    74     def writeVariables(self): 
    75         self.latVar[:] = self.dataTable.getLats() 
    76         self.lonVar[:] = self.dataTable.getLons() 
    77         self.timeVar[:] = [15,45,75,105,135,165,195,225,255,285,315,345] 
    78         self.preVar[:] = self.getPrecipetationData() 
    79  
    80     def getPrecipetationData(self): 
    81         for lt in self.latVar: 
    82             for ln in self.lonVar: 
    83                 print lt + ln 
    84                 precipetationDataPt = self.dataTable.lookup(lt,ln) 
    85                 self.preVar[0][lt][ln] = precipetationDataPt.jan 
    86                 self.preVar[1][lt][ln] = precipetationDataPt.feb 
    87                 self.preVar[2][lt][ln] = precipetationDataPt.mar 
    88                 self.preVar[3][lt][ln] = precipetationDataPt.apr 
    89                 self.preVar[4][lt][ln] = precipetationDataPt.may 
    90                 self.preVar[5][lt][ln] = precipetationDataPt.jun 
    91                 self.preVar[6][lt][ln] = precipetationDataPt.jul 
    92                 self.preVar[7][lt][ln] = precipetationDataPt.aug 
    93                 self.preVar[8][lt][ln] = precipetationDataPt.sep 
    94                 self.preVar[9][lt][ln] = precipetationDataPt.oct 
    95                 self.preVar[10][lt][ln] = precipetationDataPt.nov 
    96                 self.preVar[11][lt][ln] = precipetationDataPt.dec 
    97  
     151def test(): 
     152    testReader() 
     153    testWriter() 
    98154 
    99155def main(): 
     
    101157    precipetationReader = PrecipetationReader() 
    102158    precipetationReader.readData(sys.argv[1]) 
    103     print precipetationReader.getLookupTable().lookup(-59.083,  -26.583).jan 
    104     print 'lat len ' + str(len(precipetationReader.getLookupTable().getLats())) 
    105     print 'lon len ' + str(len(precipetationReader.getLookupTable().getLons())) 
    106159    precipetationWriter = PrecipetationWriter(precipetationReader.getLookupTable()) 
    107160    precipetationWriter.writeData(sys.argv[2]) 
    108161 
    109162if __name__ == "__main__": 
     163    #test() 
    110164    main() 
Note: See TracChangeset for help on using the changeset viewer.