source: TOOLS/MOSAIX/check_conserv.py @ 4085

Last change on this file since 4085 was 4085, checked in by omamce, 4 years ago

OM.M. : add keywords to check_conserv.py

  • Property svn:keywords set to Date Revision HeadURL Author
File size: 7.8 KB
Line 
1### ===========================================================================
2###
3### Checks integrals
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##
16import cdms2
17import numpy as np
18import glob
19import sys, argparse, textwrap
20import nemo
21
22# Creating a parser
23parser = argparse.ArgumentParser (
24    description = """Compute integral of matching coupling fields on atmosphere and ocean side.""",
25    epilog='-------- This is the end of the help message --------')
26
27# Adding arguments
28parser.add_argument ( '--file'     , help='Input file, e.g. dia_tICO40_to_teORCA1.2_1storder_false.nc', type=str, default='dia_tico_to_torc_1storder_false.nc' )
29parser.add_argument ( '--grids'    , help='grids file  ', type=str, default='grids.nc' )
30parser.add_argument ( '--areas'    , help='areas file  ', type=str, default='areas.nc' )
31parser.add_argument ( '--masks'    , help='masks file  ', type=str, default='masks.nc' )
32parser.add_argument ( '--fracs'    , help='fracs file  ', type=str, default=None )
33parser.add_argument ( '--guess'    , help='If set, build grids, areas and masks file name from the dia file name', action="store_true", default=False )
34parser.add_argument ( '--src_int'  , help='Source field is integrated'     , action="store_true" )
35parser.add_argument ( '--dst_int'  , help='Destination field is integrated', action="store_true" )
36parser.add_argument ( '--src_nofrac' , help='Do not apply frac on source'     , action="store_true" )
37parser.add_argument ( '--dst_nofrac' , help='Do not apply frac on destination', action="store_true" )
38
39# Parse command line
40myargs = parser.parse_args()
41
42#
43name_input = myargs.file
44
45# Determination des modeles
46Elements = myargs.file.split('_')
47
48#
49Name_src = Elements[1][1:] ; Name_dst = Elements[3][1:]
50grid_src = Elements[1][0]  ; grid_dst = Elements[3][0]
51
52Grid_src = grid_src.upper() ; Grid_dst = grid_dst.upper()
53if Grid_src in ( 'A', 'C') : Grid_src = 'T'
54if Grid_dst in ( 'A', 'C') : Grid_dst = 'T'
55
56# src side
57if Name_src.count('ico') != 0 or Name_src.count('lmd') != 0 : name_src = Name_src ; ShortNames = True
58else :
59    if Name_src.count('ICO') : name_src = 'ico' ; ShortNames = False
60    if Name_src.count('LMD') : name_src = 'lmd' ; ShortNames = False
61
62if Name_dst.count('orc') != 0    : name_src = Name_src ; ShortNames = True
63else :
64    if Name_src.count('ORC')!= 0 : name_src = 'orc'    ; ShortNames = False
65       
66# dst side
67if Name_dst.count('ico')!= 0 or Name_dst.count('lmd') != 0 : name_dst = Name_dst ; ShortNames = True
68else :
69    if Name_dst.count('ICO') : name_dst = 'ico' ; ShortNames = False
70    if Name_dst.count('LMD') : name_dst = 'lmd' ; ShortNames = False
71
72if Name_dst.count('orc') != 0 : name_dst = Name_dst ; ShortNames = True
73else :
74    if Name_dst.count('ORC') != 0 : name_dst = 'orc' ; ShortNames = False
75   
76if Name_src.count('ORC') : CplModel = Name_src + 'x' + Name_dst
77if Name_dst.count('ORC') : CplModel = Name_dst + 'x' + Name_src
78
79print ('CplModel : ' + CplModel )
80print ('Input file: ' + name_input )
81     
82# Coordonnees et masques
83if myargs.fracs == None :
84    if Name_src.count('ORC') : name_frc = Name_dst + '_grid_maskFrom_' + Name_src + '.nc'
85    if Name_dst.count('ORC') : name_frc = Name_src + '_grid_maskFrom_' + Name_dst + '.nc'
86else:
87    name_frc = myargs.fracs
88
89if myargs.guess :
90    n_grids =  'grids_' + CplModel + '.nc'
91    n_masks =  'masks_' + CplModel + '.nc'
92    n_areas =  'areas_' + CplModel + '.nc'
93else :
94    n_grids = myargs.grids
95    n_masks = myargs.masks
96    n_areas = myargs.areas
97
98print ('Opening grids file: ' + n_grids )
99f_grids = cdms2.open ( n_grids )
100print ('Opening masks file: ' + n_masks )
101f_masks = cdms2.open ( n_masks )
102print ('Opening areas file: ' + n_areas )
103f_areas = cdms2.open ( n_areas )
104print ('Opening fracs file: ' + name_frc )   
105f_frac  = cdms2.open ( name_frc     )
106
107#
108msk_src = np.float64(1.0) - f_masks ( grid_src + name_src + '.msk' )
109lon_src = f_grids ( grid_src + name_src + '.lon' )
110lat_src = f_grids ( grid_src + name_src + '.lat' )
111if myargs.dst_int : srf_src = np.ones ( (msk_src.shape), dtype=np.float64 )
112else              : srf_src = f_areas ( grid_src + name_src + '.srf' )
113
114#
115msk_dst = np.float64(1.0) - f_masks ( grid_dst + name_dst + '.msk' )
116lon_dst = f_grids ( grid_dst + name_dst + '.lon' )
117lat_dst = f_grids ( grid_dst + name_dst + '.lat' )
118if myargs.dst_int : srf_dst = np.ones ( (msk_dst.shape), dtype=np.float64)
119else              : srf_dst = f_areas ( grid_dst + name_dst + '.srf' )
120
121#
122if name_src in ( 'ico', 'lmd' ) :
123    if myargs.src_nofrac : frc_src = np.ones ( (msk_src.shape), dtype=np.float64)
124    else                 : frc_src = f_frac ( 'OceFrac' )
125if name_src == 'orc'     : frc_src = np.ones ( (msk_src.shape), dtype=np.float64)
126if name_dst in ( 'ico', 'lmd' ) :
127    if myargs.dst_nofrac : frc_dst = np.ones ( (msk_dst.shape), dtype=np.float64)
128    else                 : frc_dst = f_frac ( 'OceFrac' )
129if name_dst == 'orc'     : frc_dst = np.ones ( (msk_dst.shape), dtype=np.float64)
130
131# Periodicity
132nperio_src = 0
133nperio_dst = 0
134
135if Name_src == 'ORCA2.3'  : nperio_src = 4
136if Name_src == 'eORCA1.2' : nperio_src = 6
137if Name_src == 'eORCA025' : nperio_src = 6
138
139if Name_dst == 'ORCA2.3'  : nperio_dst = 4
140if Name_dst == 'eORCA1.2' : nperio_dst = 6
141if Name_dst == 'eORCA025' : nperio_dst = 6
142
143# NEMO periodicity.
144if Name_src in ('ORCA2.3', 'eORCA1.2', 'eORCA025') : msk_src = nemo.lbc_mask ( msk_src, nperio=nperio_src, cd_type=Grid_src )
145if Name_dst in ('ORCA2.3', 'eORCA1.2', 'eORCA025') : msk_dst = nemo.lbc_mask ( msk_dst, nperio=nperio_dst, cd_type=Grid_dst )
146 
147#frc_atm = f_frac ('OceFrac')[0,:] * msk_atm
148#frc_oce = np.float64(msk_oce)
149#srf_atm = f_areas ('tico.srf') * frc_atm
150#srf_oce = f_areas ('torc.srf') * frc_oce
151
152area_src = np.sum (srf_src * msk_src * frc_src)
153area_dst = np.sum (srf_dst * msk_dst * frc_dst)
154
155print ( Name_src )
156print 'mask: {:12.3} {:12.3} {:12.3}'.format( np.min(msk_src), np.max(msk_src), np.sum(msk_src) )
157print 'frac: {:12.3} {:12.3} {:12.3}'.format( np.min(frc_src), np.max(frc_src), np.sum(frc_src) )
158print 'area: {:12.3} {:12.3} {:12.3}'.format( np.min(srf_src), np.max(srf_src), np.sum(srf_src) )
159
160print ( Name_dst )
161print 'mask: {:12.3} {:12.3} {:12.3}'.format( np.min(msk_dst), np.max(msk_dst), np.sum(msk_dst) )
162print 'frac: {:12.3} {:12.3} {:12.3}'.format( np.min(frc_dst), np.max(frc_dst), np.sum(frc_dst) )
163print 'area: {:12.3} {:12.3} {:12.3}'.format( np.min(srf_dst), np.max(srf_dst), np.sum(srf_dst) )
164print (' ')
165print ( "Surfaces   : %14.6e %14.6e %11.3e"%( area_src, area_dst, (area_src-area_dst)/(area_src+area_dst)*0.5 ) )
166print (' ')
167
168##
169f_input = cdms2.open ( name_input )
170# Loop over fields
171for num in np.arange (1, 7) :
172    name_src = 'field{:02}_src'.format(num)
173    name_dst = 'field{:02}_dst'.format(num)
174    print ( '{:02}'.format(num) + ':' + name_src + ':' + name_dst)
175
176    # Reading field
177    v_src = f_input ( name_src )
178    v_dst = f_input ( name_dst )
179
180    # Compute integrals
181    sum_src = np.sum (v_src * srf_src * msk_src * frc_src)
182    sum_dst = np.sum (v_dst * srf_dst * msk_dst * frc_dst)
183   
184    # Compute average
185    mean_src = sum_src / area_src
186    mean_dst = sum_dst / area_dst
187
188    #
189    print ( "Integrals : %14.6e %14.6e %11.3e"%( sum_src , sum_dst , (sum_src - sum_dst)/(sum_src + sum_dst)*0.5 ) )
190    print ( "Averages  : %14.6e %14.6e %11.3e"%( mean_src, mean_dst, (mean_src-mean_dst)/(mean_src+mean_dst)*0.5 ) )
191    print ( " " )
192#=======================================
Note: See TracBrowser for help on using the repository browser.