source: nappy/trunk/cdms2na.py @ 343

Subversion URL: http://proj.badc.rl.ac.uk/svn/ndg/nappy/trunk/cdms2na.py@343
Revision 343, 13.3 KB checked in by astephen, 15 years ago (diff)

Latest version with new files in test directory.

  • Property svn:eol-style set to native
  • Property svn:executable set to *
  • Property svn:keywords set to Author Date Id Revision
Line 
1#!/usr/bin/env python
2
3helpMessage="""
4
5cdms2na.py
6==========
7
8Converts a CdmsFile object into one a NASA Ames file.
9
10Usage
11=====
12
13    cdms2na.py <infilename> <outfilename> [<options>]
14
15Where
16-----
17
18    infilename  - name of input file (NetCDF).
19    outfilename - name of output file (NASA Ames).
20    options     - list of options [Not yet implemented].
21   
22"""
23
24# Imports from python standard library
25import sys, os, time, string, fileinput
26sys.path.append("..")
27
28# Import cdat modules
29if sys.platform.find("win")>-1:
30    pass
31else:
32    import cdms, cdtime
33    # Define cdms variables
34    cdms.setAutoBounds("off")
35
36# Import nappy modules
37import cdmsMap
38import version
39import general
40from naCore import *
41from naError import *
42
43def compareAxes(ax1, ax2):
44    """Takes 2 cmds axis objects returning 1 if they are essentially
45    the same and 0 if not."""
46    for axtype in ("time", "level", "latitude", "longitude"):
47        if cdms.axisMatches(ax1, axtype) and not cdms.axisMatches(ax2, axtype):
48            return 0
49
50    # Check ids
51    if ax1.id!=ax2.id: return 0
52    # Check lengths
53    if len(ax1)!=len(ax2): return 0
54    # Check values
55    if ax1._data_!=ax2._data_: return 0
56    # Check units
57    if ax1.units!=ax2.units: return 0
58    # OK, I think they are the same axis!
59    return 1
60
61def compareVariables(var1, var2):
62    """Compares two cdms variables to see if they are defined on identical
63    axes. Returns 1 if true, 0 if not."""
64    try:
65        for i in range(len(var1.getAxisList())):
66            ax1=var1.getAxis(i)
67            ax2=var2.getAxis(i)       
68            if compareAxes(ax1, ax2)==0:
69                return 0
70    except:
71        return 0
72    return 1
73
74def isAuxAndVar(avar, var):
75    """Compares two cdms variables and returns true if the first dimension of the
76    main variable is the only dimension the auxiliary variable is defined against."""
77    if len(avar.shape)>1:
78        return 0
79    auxax=avar.getAxis(0)
80    varax=var.getAxis(0)
81    return compareAxes(auxax, varax)
82
83def arrayToList(array, inlist):
84
85    """ Takes an n-dimensional Numeric array and converts it to an
86    n-dimensional list object.
87    """
88    dimlist=array.shape
89    if len(dimlist[1:])>0:
90            for i in range(dimlist[0]):
91                arrayToList(inlist[i], array[i])
92    else:
93            for i in range(dimlist[0]):
94                inlist.append(array[i])
95    return inlist
96
97def listOfListsCreator(inlist, dimlist):
98    """
99    Creates a list of lists with dimensions defined in dimlist (a list of integers).
100    """
101    if len(dimlist[1:])>0:
102        for i in range(dimlist[0]):
103            inlist.append([])
104            listOfListsCreator(inlist[i], array, dimlist[1:])
105    return inlist
106
107def getBestName(var): 
108    """Returns the most appropriate variable name for a NASA Ames header"""
109    if hasattr(var, "id"): name=var.id
110    if hasattr(var, "name"): name=var.name
111    if hasattr(var, "long_name"): name=var.long_name
112    if hasattr(var, "standard_name"): name=var.standard_name
113    if hasattr(var, "units"): name="%s (%s)" % (name, var.units)
114    return name
115
116def getMissingValue(var):
117    if hasattr(var, "missing_value"): return var.missing_value
118    if hasattr(var, "fill_value"): return var.fill_value
119    if hasattr(var, "_FillValue"): return var._FillValue
120    return 1.E20
121
122
123def fixHeaderLength(file):
124    lines=open(file).readlines()[:200]
125    count=0
126    for line in lines:
127        count=count+1
128        if line=="***Data Section begins on the next line***\n":
129            break
130    # Now replace 1st line NLHEAD with count
131    firstline=lines[0]
132    newline="%s    %s" % (count, firstline.split()[-1])
133    for line in fileinput.input(file, inplace=1):
134        if line==firstline:
135          # write the changed line
136          sys.stdout.write(newline+"\n")
137        else:
138          # write the line unchanged
139          sys.stdout.write(line)
140
141
142def cdms2na(infilename, outfilename, variable=None, ffi="automatic", rules=None):
143    cdmsfile=cdms.open(infilename)
144    globals=cdmsfile.attributes
145    if not variable:
146        vars=[]
147        for var in cdmsfile.variables:
148            vars.append(cdmsfile(var))   
149    else:
150        vars=[cdmsfile(variable)]
151    builder=CdmsToNABuilder(vars, globals)
152    #print builder.naDict["X"]
153    general.openNAFile(outfilename, 'w', builder.naDict)
154    nlhead=fixHeaderLength(outfilename)   
155    print "%s written successfully." % outfilename
156
157
158class CdmsToNABuilder(NACore):
159
160    def __init__(self, vars, global_attributes={}):
161        self.naDict={}
162        self.vars=vars
163        self.globals=global_attributes
164        (vars, auxVars)=self.analyseVariables()
165        self.naDict["NLHEAD"]="-999"
166        self.defineNAVars(vars)
167        self.defineNAAuxVars(auxVars)
168        self.defineNAGlobals()
169        self.defineNAComments()
170        self.defineGeneralHeader()
171        # Quick fudge
172        if self.naDict["FFI"]==1001: self.naDict["X"]=self.naDict["X"][0]
173
174    def analyseVariables(self):
175        vars4NA=[self.vars[0]]
176        auxVars4NA=[]
177        shape=self.vars[0].shape
178        ndims=len(shape)
179        self.naDict["NIV"]=ndims
180
181        # Work out which File Format Index is appropriate
182        if ndims in (2,3,4):
183            self.naDict["FFI"]=10+(ndims*1000)
184        elif ndims>4:
185            raise "Cannot write variables defined against greater than 4 axes in NASA Ames format."
186        else:
187            if len(auxVars4NA)>0 or (self.naDict.has_key("NAUXV") and self.naDict["NAUXV"]>0):
188                self.naDict["FFI"]=1010
189            else:
190                self.naDict["FFI"]=1001
191
192        axes=self.vars[0].getAxisList()
193       
194        # Get other variable info
195        for var in self.vars[1:]:
196            if len(var.shape)!=ndims or var.shape!=shape: 
197                # Could it be an auxiliary variable
198                if len(var.shape)!=1:  continue
199                caxis=var.getAxis(0)
200                if compareAxes(axes[0], caxis)==0:  continue
201                # I think it is an auxiliary variable
202                auxVars4NA.append(var) 
203            else:
204                caxes=var.getAxisList()
205                for i in range(ndims):           
206                    if compareAxes(axes[i], caxes[i])==0:
207                        continue
208                # OK, I think they are compatible
209                vars4NA.append(var)
210
211        return (vars4NA, auxVars4NA)
212
213    def defineNAVars(self, vars):
214        self.naDict["NV"]=len(vars)
215        self.naDict["VNAME"]=[]
216        self.naDict["VMISS"]=[]
217        self.naDict["VSCAL"]=[]
218        self.naDict["V"]=[]
219        for var in vars:
220            name=getBestName(var)
221            self.naDict["VNAME"].append(name)
222            miss=getMissingValue(var)
223            if type(miss)!=float:  miss=miss[0]
224            self.naDict["VMISS"].append(miss)
225            #print self.naDict["VMISS"]
226            self.naDict["VSCAL"].append(1)
227            # AND THE ARRAY
228            # Populate the variable list 
229            ######## NOTE - might not have to do this #####
230            ######## It  might handle writing from a Numeric array ########
231            self.naDict["V"].append(var._data)
232            #listOfListsCreator(inlist, var.shape)
233            #arrayToList(var, inlist)
234
235            if not self.naDict.has_key("X"):
236                self.naDict["NXDEF"]=[]
237                self.naDict["NX"]=[]
238                # Create independent variable information
239                self.ax0=var.getAxis(0)
240                self.naDict["X"]=[list(self.ax0._data_)]
241                self.naDict["XNAME"]=[getBestName(self.ax0)]
242                if len(self.ax0)==1:
243                    self.naDict["DX"]=[0]
244                else:
245                    incr=self.ax0[1]-self.ax0[0]
246                    for i in range(1, len(self.ax0)):
247                        if (self.ax0[i]-self.ax0[i-1])!=incr:
248                            self.naDict["DX"]=[0]
249                            break
250                    else:
251                        self.naDict["DX"]=[incr]
252                # Now sort the rest of the axes
253                for axis in var.getAxisList()[1:]:
254                    self.getAxisDefinition(axis)
255
256    def defineNAAuxVars(self, auxVars):
257        self.naDict["NAUXV"]=len(auxVars)
258        self.naDict["ANAME"]=[]
259        self.naDict["AMISS"]=[]
260        self.naDict["ASCAL"]=[]
261        self.naDict["A"]=[]
262        for var in auxVars:
263            name=getBestName(var)
264            self.naDict["ANAME"].append(name)
265            miss=getMissingValue(var)
266            if type(miss)!=float:  miss=miss[0]
267            self.naDict["SMISS"].append(miss)
268            self.naDict["ASCAL"].append(1)
269            # AND THE ARRAY
270            # Populate the variable list 
271            ######## NOTE - might not have to do this #####
272            ######## It  might handle writing from a Numeric array ########
273            self.naDict["A"].append(var._data)
274            #listOfListsCreator(inlist, var.shape)
275            #arrayToList(var, inlist)     
276
277    def getAxisDefinition(self, axis):
278        self.naDict["NX"].append(len(axis))
279        self.naDict["XNAME"].append(getBestName(axis))
280        # If only one item in axis values
281        if len(axis)==1:
282            self.naDict["DX"].append(0)
283            self.naDict["NXDEF"].append(len(axis))
284            self.naDict["X"].append(list(axis._data_))       
285            return
286   
287        incr=axis[1]-axis[0]
288        for i in range(1, len(axis)):
289            if (axis[i]-axis[i-1])!=incr:
290                self.naDict["DX"].append(0)
291                self.naDict["NXDEF"].append(len(axis))
292                self.naDict["X"].append(list(axis._data_))
293                break
294        else:
295            self.naDict["DX"].append(incr)
296            self.naDict["NXDEF"].append(3)
297            self.naDict["X"].append(axis[:3])
298        return
299
300    def defineNAGlobals(self):
301        globalmap=cdmsMap.cdmsMap["toNA"]
302        self.extra_comments=[]
303        for key in self.globals.keys():
304            if type(self.globals[key]) not in (str, float, int): continue
305            if key in globalmap.keys():
306                if key=="history":
307                    timestring=time.strftime("%Y-%m-%d %H:%M:%S", time.localtime(time.time()))
308                    self.history="History:\n\t%s - Converted to NASA Ames format using nappy-%s.\n\t%s" % (timestring, version.version, self.globals[key])
309                    self.history=self.history.split("\n") 
310                   
311                elif key=="institution":
312                    self.naDict["ONAME"]=self.globals[key]
313                    self.naDict["ORG"]=self.globals[key]
314                elif key=="comment":
315                    self.ncom=[self.globals[key]]
316                elif key in ("Conventions", "references"):
317                    pass
318                else:
319                    self.naDict[globalmap[key]]=self.globals[key]
320            else:
321                self.extra_comments.append("%s:   %s" % (key, self.globals[key]))
322        return
323
324
325    def defineNAComments(self, normal_comments=[], special_comments=[]):
326        NCOM=["***NASA Ames Normal Comments follow***"]
327        if hasattr(self, "ncom"):  normal_comments=self.ncom+normal_comments
328        for ncom in normal_comments:
329            NCOM.append(ncom)
330        NCOM.append("")
331        NCOM.append("Additional Global Attributes defined in Cdmsfile and not translated elsewhere:")
332
333        if hasattr(self, "extra_comments"):
334            for excom in self.extra_comments:
335                NCOM.append(excom)
336
337        if hasattr(self, "history"):
338            for h in self.history:
339                NCOM.append(h)
340
341        NCOM.append("")
342        NCOM.append("***Data Section begins on the next line***")
343
344        SCOM=["***NASA Ames Special Comments follow***"]
345        for scom in special_comments:
346            SCOM.append(scom)
347
348        SCOM.append("Additional Variable Attributes defined in Cdmsfile and not translated elsewhere:")
349        used_var_atts=("name", "long_name", "standard_name", "id", "missing_value", "fill_value", "units")
350        for var in self.vars:
351            for scom in var.attributes.keys():
352                name=getBestName(var)
353                if type(scom) in (str, float, int) and scom not in used_var_atts:
354                    SCOM.append("\t%s: %s - %s" % (name, scom, var.attributes[scom]))
355
356        self.naDict["NCOM"]=NCOM
357        self.naDict["NNCOML"]=len(NCOM)
358        self.naDict["SCOM"]=SCOM
359        self.naDict["NSCOML"]=len(SCOM)
360        return
361
362    def defineGeneralHeader(self, header_items={}):
363        """Defines known header items and overwrites any with header_items key/value pairs."""
364        self.naDict["DATE"]=time.strftime("%Y %m %d", time.localtime(time.time())).split()
365        if self.ax0.isTime():
366            # Get first date in list
367            (unit, start_date)=re.match("(\w+)\s+?since\s+?(\d+-\d+-\d+)", self.ax0.units).groups()
368           
369            comptime=cdtime.s2c(start_date)
370            first_day=comptime.add(self.naDict["X"][0][0], getattr(cdtime, unit.capitalize()))
371            self.naDict["RDATE"]=string.replace(str(first_day).split(" ")[0], "-", " ").split()
372        else: 
373            self.naDict["RDATE"]=self.naDict["DATE"]
374        self.naDict["IVOL"]=1
375        self.naDict["NVOL"]=1
376        for key in header_items.keys():
377             self.naDict[key]=header_items[key]
378        return
379               
380
381
382if __name__=="__main__":
383
384    args=sys.argv[1:]
385    if len(args)<2:
386        print helpMessage
387        print "Incorrect number of arguments used."
388        sys.exit()
389    elif len(args)==3:
390        cdms2na(args[0], args[1], args[2]) 
391    else:
392        cdms2na(args[0], args[1]) 
Note: See TracBrowser for help on using the repository browser.