source: TOOLS/CPLRESTART/FillOceRestart.py @ 3741

Last change on this file since 3741 was 3741, checked in by omamce, 6 years ago

O.M. : improve documentation
  • Property svn:keywords set to Author Date Revision Id HeadURL
File size: 6.9 KB
Line 
1### ===========================================================================
2###
3### Fill a coupler restart, ocean side
4###
5### ===========================================================================
6##
7##  Warning, to install, configure, run, use any of Olivier Marti's
8##  software or to read the associated documentation you'll need at least
9##  one (1) brain in a reasonably working order. Lack of this implement
10##  will void any warranties (either express or implied).
11##  O. Marti assumes no responsability for errors, omissions,
12##  data loss, or any other consequences caused directly or indirectly by
13##  the usage of his software by incorrectly or partially configured
14##  personal.
15##
16## SVN information
17__Author__   = "$Author$"
18__Date__     = "$Date$"
19__Revision__ = "$Revision$"
20__Id__       = "$Id$"
21__HeadURL    = "$HeadURL$"
22##
23import netCDF4
24import numpy.ma as np
25from scipy import ndimage
26import nemo
27import shutil, getopt, sys
28
29def MyFill (InputData=None) :
30    """
31    Replace the value of masked 'InputData' cells by the value of the nearest valid data cell
32    From https://stackoverflow.com/questions/5551286/filling-gaps-in-a-numpy-array
33   
34    Input:  InputData : numpy.ma array of any dimension.
35
36    Output: Return a filled array.
37    """   
38
39    Invalid = np.where ( InputData[:,:].mask, True, False)
40    Indices = ndimage.distance_transform_bf(Invalid, metric='euclidean', sampling=None, return_distances=False,
41                                            return_indices=True, distances=None )
42    FilledData = InputData[tuple(Indices)]
43    return FilledData
44
45def usage () :
46    __help__ = """%(prog)s usage :
47python %(prog)s [-d] -i <sstoce file> [-n <perio>] [-v variables] [-o <output file>]
48 -d           | --debug                : debug
49 -i <file>    | --input=<file>         : input file  (default: none)
50 -o <file>    | --output=<file>        : output file (default : build a name form input file)
51 -r           | --replace              : replace input file by new file with filled variables
52 -v <varlist> | --variable=<variables> : list of variable to fill (defautlt: all variable in file)
53 -x           | --exclude              : fills all variable in files, except those given in -v|--variable
54 -n <perio>   | --perio=<perio> : periodicity type (default: try to guess)
55    1, 4, 6 : Cyclic on i dimension (generaly longitudes)
56    2       : Obsolete (was symmetric condition at southern boundary ?)
57    3, 4    : North fold T-point pivot (legacy ORCA2)
58    5, 6    : North fold F-point pivot (ORCA1, ORCA025, ORCA2 with new grid for paleo)
59    If <perio> is not specified, %(prog)s will try to guess it from the grid dimensions
60example :
61    python %(prog)s -n 4 -i sstoce_ORCA2.3.nc
62    python %(prog)s -n 4 -i sstoce_ORCA2.3.nc -v O_SSTSST,OIceFrc
63    """
64    print ( __help__ % { 'prog':sys.argv[0] } )
65
66## Default input parameters
67InFile      = None
68OuFile      = None
69replace     = False
70ListVarName = None
71nperio      = None
72Debug       = False
73ListExclude = None ; Exclude=False
74
75## Command line options
76try:
77    myopts, myargs = getopt.getopt ( sys.argv[1:], 'i:o:rv:xp:hd',
78                                         [ 'input=', 'output=', 'replace', 'exclude', 'variable=', 'variables=', 'perio=', 'debug', 'help' ] )
79except getopt.GetoptError as cmdle :
80    print ( "Command line error : "+str(cmdle)+"\n" )
81    usage ()
82    sys.exit(1)
83
84for myopt, myval in myopts :
85    if   myopt in [ '-h', '--help'    ] : usage () ; sys.exit (0) ;
86    elif myopt in [ '-d', '--debug'   ] : Debug    = True  ;
87    elif myopt in [ '-i', '--input'   ] : InFile   = myval ;
88    elif myopt in [ '-o', '--output'  ] :
89        if replace :
90            print ("Error : you can not specify both -r|--replace and -o|--output" )
91            usage ()
92            sys.exit(1)
93        else :
94            if Debug : print ("Out file set to "+OuFile)
95            OuFile   = myval ;
96    elif myopt in [ '-r', '--replace' ] :
97        if OuFile != None :
98            print ("Error : you can not specify both -r|--replace and -o|--output" )
99            usage ()
100            sys.exit(1)
101        else :
102            if Debug : print ("Out file set to input file")
103            OuFile   = InFile ;
104    elif myopt in [ '-p', '--perio'     ] : nperio = int(myval) ;
105    elif myopt in [ '-v', '--variable', '--variables'  ] :
106        if Exclude : ListExclude = myval.split(',')         
107        else :       ListVarName = myval.split(',')
108    elif myopt in [ '-x', '--exclude'  ] :
109        Exclude = True
110        if ListVarName != None :
111            ListExclude = ListVarName
112            ListVarName  = None
113 
114if OuFile == None :
115    OuFile = InFile.replace ( ".nc", "_filled.nc" )
116    print ("Creates output file name: "+OuFile)
117
118# Copy the input file if needed
119if OuFile != InFile : shutil.copyfile ( InFile, OuFile )
120
121print ("Output file: " + OuFile )
122
123# Open file
124OuFile = netCDF4.Dataset ( OuFile , "r+" )
125
126# Try to guess periodicity type
127if nperio == None :
128    print ( "Trying to guess nperio parameter" )
129    jpoi = OuFile.dimensions["x"].size
130    jpoj = OuFile.dimensions["y"].size
131    print ("Grid dimensions: ("+str(jpoj)+", "+str(jpoi)+")")
132    if   'ORCA2' in InFile :
133        print ("ORCA 2 grid found from file name: nperio may vary for this configuration")
134        print ("Choosen nperio=4")
135    elif 'ORCA1' in InFile :
136        if 'eORCA1' in InFile : 
137            print ("eORCA 1 grid found from file name, nperio=6")
138        else :
139            print ("ORCA 1 grid found from file name, nperio=6")
140    elif (jpoj, jpoi) == (149, 182) :
141        print ("ORCA 2 grid found from dimension: nperio may vary for this configuration")
142        print ("Choosen nperio=4")
143        nperio = 4
144    elif (jpoj, jpoi) == (332, 292) :
145        nperio = 6
146        print ("ORCA1 grid found from dimensions, nperio=6" )
147    elif (jpoj, jpoi) == (332, 362) :
148        nperio = 6
149        print ("eORCA1 grid found from dimensions, nperio=6" )
150       
151if nperio == None :
152    print ("%(prog)s couldn't guess the periodicity type of your file")
153    print ("Please specify -p|--perio")
154    usage ()
155    sys.exit(1)
156   
157# Get variables from file is need
158if ListVarName == None : ListVarName = list ( OuFile.variables.keys() )
159
160# Exclude some var if needed
161if Exclude : 
162    for Var in ListExclude : 
163        if Var in ListVarName : ListVarName.remove(Var)
164
165# Loop on variables
166for VarName in ListVarName :
167    Var    = OuFile.variables[VarName]
168    if 'mask' in dir(Var[...]) :
169        print ( "Working on " + VarName )
170        NewVar = MyFill ( InputData = Var[...] )
171        NewVar = nemo.lbc (NewVar, nperio=nperio, cd_type='T', psgn=1.0)
172        OuFile.variables[VarName][:,:] = NewVar[:,:]
173    else :
174        print ( VarName + " is not masked" )
175       
176# Close file : writes update variables.   
177OuFile.close()
178
179## ===========================================================================
180##
181##                               That's all folk's !!!
182##
183## ===========================================================================
Note: See TracBrowser for help on using the repository browser.