source: qesdi/data_processing/convertPre.py @ 5434

Subversion URL: http://proj.badc.rl.ac.uk/svn/ndg/qesdi/data_processing/convertPre.py@5434
Revision 5434, 7.0 KB checked in by dgrant, 12 years ago (diff)

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

Line 
1import netCDF4
2import sys
3import numpy
4import time
5
6from TMCReader import TMCReader
7from TMCWriter import TMCWriter
8from netCDF4 import Dataset
9
10
11f32 = numpy.float32
12
13class PrecipetationReader(TMCReader):
14    """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)
20
21    def storeLine(self, line):
22        precipDataPt = PrecipetationDataPt(line.split())
23        self.coordTable.insert(precipDataPt, precipDataPt.lat, precipDataPt.lon)
24
25class PrecipetationDataPt:
26    """Class to store data for a specific point"""
27
28    def __init__(self,preVals):
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])
55
56
57class PrecipetationWriter(TMCWriter):
58    """Class to read precipetation data"""
59
60    def writeData(self, path):
61        self.rootgrp = Dataset(path, 'w',format='NETCDF3_CLASSIC')
62        self.createDimensions()
63        self.createVariables()
64        self.preVar = self.rootgrp.createVariable('precipetation_means','f4',('time','lat','lon',))
65        self.cvVar = self.rootgrp.createVariable('precipetation_cvs','f4',('time','lat','lon',))
66        self.setGlobalAttributes()
67        self.setVariableAttributes()
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()
77        self.rootgrp.close()
78
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
130
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
147
148def testWriter():
149    table = COORDDataTable(10) 
150
151def test():
152    testReader()
153    testWriter()
154
155def main():
156    print time.ctime(time.time())
157    precipetationReader = PrecipetationReader()
158    precipetationReader.readData(sys.argv[1])
159    precipetationWriter = PrecipetationWriter(precipetationReader.getLookupTable())
160    precipetationWriter.writeData(sys.argv[2])
161
162if __name__ == "__main__":
163    #test()
164    main()
Note: See TracBrowser for help on using the repository browser.