Changeset 6518 for TOOLS


Ignore:
Timestamp:
06/22/23 11:45:13 (10 months ago)
Author:
omamce
Message:

O.M. - WATER_BUDGET
Evolutions pour Sechiba

Location:
TOOLS/WATER_BUDGET
Files:
2 edited

Legend:

Unmodified
Added
Removed

  • TOOLS/WATER_BUDGET/ATM_waterbudget.py

    • Property svn:keywords changed from Date Revision HeadURL Author Id to Date Revision HeadURL Author
    r6508 r6518  
    1212## 
    1313## 
    14 ## SVN information 
    15 #  $Author$ 
    16 #  $Date$ 
    17 #  $Revision$ 
    18 #  $Id$ 
    19 #  $HeadURL$ 
     14# SVN Information 
     15SVN = { 
     16    'Author'  : "$Author$", 
     17    'Date'    : "$Date$", 
     18    'Revision': "$Revision$", 
     19    'Id'      : "$Id: ATM_waterbudget.py 6508 2023-06-13 10:58:38Z omamce $", 
     20    'HeadURL' : "$HeadUrl: svn+ssh://omamce@forge.ipsl.jussieu.fr/ipsl/forge/projets/igcmg/svn/TOOLS/WATER_BUDGET/ATM_waterbudget.py $" 
     21    } 
    2022 
    2123### 
    … …  
    3436 
    3537## Creates parser for reading .ini input file 
    36 config = configparser.ConfigParser (interpolation=configparser.ExtendedInterpolation() ) 
     38config = configparser.ConfigParser ( interpolation=configparser.ExtendedInterpolation() ) 
    3739config.optionxform = str # To keep capitals 
    3840 
    3941## Experiment parameters 
    40 ATM=None ; ATM_HIS='latlon' ; SRF_HIS='latlon' ; RUN_HIS='latlon' ; ORCA=None ; NEMO=None ; OCE_relax=False ; OCE_icb=False ; Coupled=False ; Routing=None ; TarRestartPeriod_beg=None ; TarRestartPeriod_end=None ; Comment=None ; Period=None ; Title=None 
     42ATM=None ; ATM_HIS='latlon' ; SRF_HIS='latlon' ; RUN_HIS='latlon' ; ORCA=None ; NEMO=None ; OCE_relax=False 
     43OCE_icb=False ; Coupled=False ; Routing=None 
     44TarRestartPeriod_beg=None ; TarRestartPeriod_end=None ; Comment=None ; Period=None ; Title=None 
    4145 
    4246## 
    4347ARCHIVE=None ; STORAGE=None ; SCRATCHDIR=None ; R_IN=None ; rebuild=None 
    44 TmpDir=None ; RunDir=None ; FileOut=None ; dir_ATM_his=None ; dir_SRF_his=None ; dir_OCE_his=None ; dir_ICE_his=None ; FileCommon=None 
    45 file_ATM_his=None ; file_SRF_his=None ; file_RUN_his=None ; file_OCE_his=None ;  file_ICE_his=None ; file_OCE_sca=None 
    46 tar_restart_beg=None ; tar_restart_end=None ; file_ATM_beg=None ; file_ATM_end=None ; file_DYN_beg=None ; file_DYN_end=None ; file_SRF_beg=None ; file_SRF_end=None 
    47 file_RUN_beg=None ; file_RUN_end=None ; file_RUN_end=None ; file_OCE_beg=None ; file_ICE_beg=None ; file_OCE_beg=None ; file_OCE_end=None ; file_ICE_beg=None ; file_OCE_end=None ; file_ICE_end=None ; ContinueOnError=False ; ErrorCount=0 
    48 tar_restart_beg_ATM=None ; tar_restart_beg_DYN=None ; tar_restart_beg_SRF=None ; tar_restart_beg_RUN=None ; tar_restart_beg_OCE=None ; tar_restart_beg_ICE=None 
    49 tar_restart_end_ATM=None ; tar_restart_end_DYN=None ; tar_restart_end_SRF=None ; tar_restart_end_RUN=None ; tar_restart_end_OCE=None ; tar_restart_end_ICE=None 
    50 SortIco = False 
     48TmpDir=None ; RunDir=None ; FileOut=None 
     49dir_ATM_his=None ; dir_SRF_his=None ; dir_OCE_his=None ; dir_ICE_his=None 
     50FileCommon=None ; file_ATM_his=None ; file_SRF_his=None ; file_RUN_his=None 
     51file_OCE_his=None ;  file_ICE_his=None ; file_OCE_sca=None 
     52tar_restart_beg=None ; tar_restart_end=None ; file_ATM_beg=None ; file_ATM_end=None ; file_DYN_beg=None 
     53file_DYN_end=None ; file_SRF_beg=None ; file_SRF_end=None 
     54file_RUN_beg=None ; file_RUN_end=None ; file_RUN_end=None ; file_OCE_beg=None 
     55file_ICE_beg=None ; file_OCE_beg=None 
     56file_OCE_end=None ; file_ICE_beg=None ; file_OCE_end=None ; file_ICE_end=None 
     57tar_restart_beg_ATM=None ; tar_restart_beg_DYN=None ; tar_restart_beg_SRF=None 
     58tar_restart_beg_RUN=None ; tar_restart_beg_OCE=None ; tar_restart_beg_ICE=None 
     59tar_restart_end_ATM=None ; tar_restart_end_DYN=None ; tar_restart_end_SRF=None 
     60tar_restart_end_RUN=None ; tar_restart_end_OCE=None ; tar_restart_end_ICE=None 
     61ContinueOnError=False ; ErrorCount=0 ; SortIco = False 
    5162 
    5263# Read command line arguments 
    … …  
    7485            exec  ( f'{VarName} = wu.setNone ({VarName})' ) 
    7586            exec  ( f'wu.{VarName} = {VarName}' ) 
    76             print ( f'{VarName:25} set to : {locals()[VarName]:}' ) 
     87            print ( f'    {VarName:21} set to : {locals()[VarName]:}' ) 
    7788            #exec ( f'del {VarName}' ) 
    7889 
    … …  
    8192##-- Some physical constants 
    8293if wu.unDefined ( 'Ra' )           : Ra          = 6366197.7236758135 #-- Earth Radius (m) 
    83 if wu.unDefined ( 'Grav' )         : Grav        = 9.81               #-- Gravity (m^2/s 
    84 if wu.unDefined ( 'ICE_rho_ice' )  : ICE_rho_ice = 917.0              #-- Ice volumic mass (kg/m3) in LIM3 
    85 if wu.unDefined ( 'ICE_rho_sno')   : ICE_rho_sno = 330.0              #-- Snow volumic mass (kg/m3) in LIM3 
     94if wu.unDefined ( 'Grav' )         : Grav        = 9.81               #-- Gravity (m^2/s) 
     95if wu.unDefined ( 'ICE_rho_ice' )  : ICE_rho_ice = 917.0              #-- Ice volumic mass (kg/m3) in LIM3 or SI3 
     96if wu.unDefined ( 'ICE_rho_sno')   : ICE_rho_sno = 330.0              #-- Snow volumic mass (kg/m3) in LIM3 or SI3 
    8697if wu.unDefined ( 'OCE_rho_liq' )  : OCE_rho_liq = 1026.0             #-- Ocean water volumic mass (kg/m3) in NEMO 
    8798if wu.unDefined ( 'ATM_rho' )      : ATM_rho     = 1000.0             #-- Water volumic mass in atmosphere (kg/m^3) 
    … …  
    89100if wu.unDefined ( 'RUN_rho' )      : RUN_rho     = 1000.0             #-- Water volumic mass of rivers (kg/m^3) 
    90101if wu.unDefined ( 'ICE_rho_pnd' )  : ICE_rho_pnd = 1000.              #-- Water volumic mass in ice ponds (kg/m^3) 
    91 if wu.unDefined ( 'YearLength' )   : YearLength  = 365.25 * 24. * 60. * 60. #-- Year length (s) - Use only to compute drif in approximate unit  
     102if wu.unDefined ( 'YearLength' )   : YearLength  = 365.25 * 86400.    #-- Year length (s) - Use only to compute drif in approximate unit  
    92103             
    93104##-- Set libIGCM and machine dependant values 
    94105if not 'Files' in config.keys () : config['Files'] = {} 
    95106 
    96 config['Physics'] = { 'Ra':Ra, 'Grav':Grav, 'ICE_rho_ice':ICE_rho_ice, 'ICE_rho_sno':ICE_rho_sno, 'OCE_rho_liq':OCE_rho_liq, 'ATM_rho':ATM_rho, 'SRF_rho':SRF_rho, 'RUN_rho':RUN_rho} 
    97  
    98 config['Config'] = { 'ContinueOnError':ContinueOnError, 'SortIco':SortIco} 
     107config['Physics'] = { 'Ra':Ra, 'Grav':Grav, 'ICE_rho_ice':ICE_rho_ice, 'ICE_rho_sno':ICE_rho_sno, 
     108                      'OCE_rho_liq':OCE_rho_liq, 'ATM_rho':ATM_rho, 'SRF_rho':SRF_rho, 'RUN_rho':RUN_rho} 
     109 
     110config['Config'] = { 'ContinueOnError':ContinueOnError, 'SortIco':SortIco } 
    99111 
    100112## -------------------------- 
    … …  
    106118 
    107119if libIGCM : 
    108     config['libIGCM'] = {'ARCHIVE':ARCHIVE, 'STORAGE':STORAGE, 'SCRATCHDIR':SCRATCHDIR, 'R_IN':R_IN, 'rebuild':rebuild } 
     120    config['libIGCM'] = { 'ARCHIVE':ARCHIVE, 'STORAGE':STORAGE, 'SCRATCHDIR':SCRATCHDIR, 'R_IN':R_IN, 'rebuild':rebuild } 
    109121 
    110122##-- Import specific module 
    … …  
    161173if wu.unDefined ('R_FIGR'     ) : R_FIGR      = os.path.join ( STORAGE, 'IGCM_OUT', L_EXP ) 
    162174 
    163 config['libIGCM'] = { 'R_OUT':R_OUT, 'R_BUF':R_BUF, 'L_EXP':L_EXP, 'R_BUFR':R_BUFR, 'POST_DIR':POST_DIR, 'REBUILD_DIR':REBUILD_DIR, 'R_BUF_KSH':R_BUF_KSH, 'R_FIGR':R_FIGR } 
     175config['libIGCM'] = { 'R_OUT':R_OUT, 'R_BUF':R_BUF, 'L_EXP':L_EXP, 'R_BUFR':R_BUFR, 'POST_DIR':POST_DIR, 
     176                      'REBUILD_DIR':REBUILD_DIR, 'R_BUF_KSH':R_BUF_KSH, 'R_FIGR':R_FIGR } 
    164177 
    165178##- Set directory to extract files 
    … …  
    213226      
    214227echo ('\nOpen history files' ) 
    215 if wu.unDefined ( 'file_ATM_his' ) :  
    216     file_ATM_his = os.path.join (  dir_ATM_his, f'{FileCommon}_histmth.nc' ) 
     228if wu.unDefined ( 'file_ATM_his' ) : 
     229    if ATM_HIS == 'latlon' : 
     230        file_ATM_his = os.path.join ( dir_ATM_his, f'{FileCommon}_histmth.nc' ) 
     231    if ATM_HIS == 'ico' : 
     232        file_ATM_his = os.path.join ( dir_ATM_his, f'{FileCommon}_histmth_ico.nc' ) 
    217233    config['Files']['file_ATM_his'] = file_ATM_his 
    218234if wu.unDefined ( 'file_SRF_his' ) : 
    219     file_SRF_his = os.path.join (  dir_SRF_his, f'{FileCommon}_sechiba_history.nc' ) 
     235    if ATM_HIS == 'latlon' : 
     236        file_SRF_his = os.path.join ( dir_SRF_his, f'{FileCommon}_sechiba_history.nc' ) 
     237    if ATM_HIS == 'ico' : 
     238        file_SRF_his = os.path.join ( dir_SRF_his, f'{FileCommon}_sechiba_history_ico.nc' ) 
    220239    config['Files']['file_SRF_his'] = file_SRF_his 
    221 #if Routing == 'SECHIBA'  : 
    222 #     file_RUN_his = os.path.join (  dir_SRF_his, f'{FileCommon}_sechiba_history.nc' ) 
     240 
    223241if Routing == 'SIMPLE' : 
    224     if file_RUN_his == None :  
    225         file_RUN_his = os.path.join (  dir_SRF_his, f'{FileCommon}_sechiba_history.nc' ) 
     242    if file_RUN_his == None : 
     243        if ATM_HIS == 'latlon' : 
     244            file_RUN_his = os.path.join ( dir_SRF_his, f'{FileCommon}_sechiba_history.nc' ) 
     245        if ATM_HIS == 'ico' : 
     246            file_RUN_his = os.path.join ( dir_SRF_his, f'{FileCommon}_sechiba_history_ico.nc' ) 
    226247        config['Files']['file_RUN_his'] = file_RUN_his 
    227248 
    … …  
    449470 
    450471## Write the full configuration 
    451 config_out = open (FullIniFile, 'w') 
    452 config.write ( config_out ) 
    453 config_out.close () 
     472#config_out = open (FullIniFile, 'w') 
     473#config.write ( config_out ) 
     474#config_out.close () 
    454475 
    455476if ICO : 
    … …  
    483504# ATM grid with cell surfaces 
    484505if LMDZ : 
    485     #ATM_lat  = lmdz.geo2point ( d_ATM_his ['lat']         , dim1D='cell' ) 
    486     #ATM_lon  = lmdz.geo2point ( d_ATM_his ['lon']         , dim1D='cell' ) 
    487     ATM_aire = lmdz.geo2point ( d_ATM_his ['aire']     [0], cumulPoles=True, dim1D='cell' ) 
    488     ATM_fter = lmdz.geo2point ( d_ATM_his ['fract_ter'][0], dim1D='cell' ) 
    489     ATM_foce = lmdz.geo2point ( d_ATM_his ['fract_oce'][0], dim1D='cell' ) 
    490     ATM_fsic = lmdz.geo2point ( d_ATM_his ['fract_sic'][0], dim1D='cell' ) 
    491     ATM_flic = lmdz.geo2point ( d_ATM_his ['fract_lic'][0], dim1D='cell' ) 
    492     #SRF_lat  = lmdz.geo2point ( d_SRF_his ['lat']      [0], dim1D='cell' ) 
    493     #SRF_lon  = lmdz.geo2point ( d_SRF_his ['lon']      [0], dim1D='cell' ) 
    494      
    495     #??SRF_aire = lmdz.geo2point ( d_SRF_his ['Areas'] * d_SRF_his ['Contfrac'], dim1D='cell' ) 
    496     SRF_aire =  ATM_aire * lmdz.geo2point ( d_SRF_his['Contfrac'], dim1D='cell' ) 
     506    ATM_lat  = lmdz.geo2point (   d_ATM_his ['lat']+0*d_ATM_his ['lon'], dim1D='cell_latlon' ) 
     507    ATM_lon  = lmdz.geo2point ( 0*d_ATM_his ['lat']+  d_ATM_his ['lon'], dim1D='cell_latlon' ) 
     508    ATM_aire = lmdz.geo2point ( d_ATM_his ['aire']     [0], cumulPoles=True, dim1D='cell_latlon' ) 
     509    ATM_fter = lmdz.geo2point ( d_ATM_his ['fract_ter'][0], dim1D='cell_latlon' ) 
     510    ATM_foce = lmdz.geo2point ( d_ATM_his ['fract_oce'][0], dim1D='cell_latlon' ) 
     511    ATM_fsic = lmdz.geo2point ( d_ATM_his ['fract_sic'][0], dim1D='cell_latlon' ) 
     512    ATM_flic = lmdz.geo2point ( d_ATM_his ['fract_lic'][0], dim1D='cell_latlon' ) 
     513    SRF_lat  = lmdz.geo2point (   d_SRF_his ['lat']+0*d_SRF_his ['lon'], dim1D='cell_latlon' ) 
     514    SRF_lon  = lmdz.geo2point ( 0*d_SRF_his ['lat']+  d_SRF_his ['lon'], dim1D='cell_latlon' ) 
     515    SRF_aire      = lmdz.geo2point ( d_SRF_his ['Areas'] * d_SRF_his ['Contfrac'], dim1D='cell_latlonl', cumulPoles=True ) 
     516    SRF_areas     = lmdz.geo2point ( d_SRF_his ['Areas'],  dim1D='cell_latlonl', cumulPoles=True ) 
     517    SRF_contfrac  = lmdz.geo2point ( d_SRF_his ['Contfrac'], dim1D='cell_latlonl' ) 
     518 
     519    #SRF_res_aire = SRF_aire 
    497520     
    498521if ICO : 
    499522    if ATM_HIS == 'latlon' : 
    500         jpja, jpia = d_ATM_his['aire'][0].shape    
    501         file_ATM_aire = os.path.join ( R_IN, 'ATM', 'GRID', f'aire_{ATM}_to_{jpia}x{jpja}.nc' ) 
    502         config['Files']['file_ATM_aire'] = file_ATM_aire 
     523        jpja, jpia = d_ATM_his['aire'][0].shape 
     524        if wu.unDefined ( 'file_ATM_aire' ) :  
     525            file_ATM_aire = os.path.join ( R_IN, 'ATM', 'GRID', f'aire_{ATM}_to_{jpia}x{jpja}.nc' ) 
     526            config['Files']['file_ATM_aire'] = file_ATM_aire 
    503527        echo ( f'Aire sur grille reguliere : {file_ATM_aire = }' ) 
    504528        d_ATM_aire = xr.open_dataset ( file_ATM_aire, decode_times=False ).squeeze() 
    505         #ATM_lat  = d_ATM_his ['lat'] 
    506         #ATM_lon  = d_ATM_his ['lon'] 
    507         ATM_aire = lmdz.geo2point ( d_ATM_aire ['aire'].squeeze(), cumulPoles=True, dim1D='cell' ) 
    508         ATM_fter = lmdz.geo2point ( d_ATM_his ['fract_ter'][0], dim1D='cell' ) 
    509         ATM_foce = lmdz.geo2point ( d_ATM_his ['fract_oce'][0], dim1D='cell' ) 
    510         ATM_fsic = lmdz.geo2point ( d_ATM_his ['fract_sic'][0], dim1D='cell' ) 
    511         ATM_flic = lmdz.geo2point ( d_ATM_his ['fract_lic'][0], dim1D='cell' ) 
    512         #SRF_lat  = d_SRF_his ['lat'] 
    513         #SRF_lon  = d_SRF_his ['lon'] 
    514         SRF_aire = lmdz.geo2point ( d_SRF_his ['Areas'] * d_SRF_his ['Contfrac'], dim1D='cell' ) 
    515          
     529        try : 
     530            ATM_lat  = lmdz.geo2point (   d_ATM_his ['lat']+0*d_ATM_his ['lon'], dim1D='cell_latlon' ) 
     531            ATM_lon  = lmdz.geo2point ( 0*d_ATM_his ['lat']+  d_ATM_his ['lon'], dim1D='cell_latlon' ) 
     532        except : 
     533            ATM_lat  = lmdz.geo2point (   d_ATM_his ['lat_dom_out']+0*d_ATM_his ['lon_dom_out'], dim1D='cell_latlon' ) 
     534            ATM_lon  = lmdz.geo2point ( 0*d_ATM_his ['lat_dom_out']+  d_ATM_his ['lon_dom_out'], dim1D='cell_latlon' ) 
     535        ATM_aire = lmdz.geo2point ( d_ATM_aire ['aire'].squeeze(), cumulPoles=True, dim1D='cell_latlon' ) 
     536        ATM_fter = lmdz.geo2point ( d_ATM_his ['fract_ter'][0], dim1D='cell_latlon' ) 
     537        ATM_foce = lmdz.geo2point ( d_ATM_his ['fract_oce'][0], dim1D='cell_latlon' ) 
     538        ATM_fsic = lmdz.geo2point ( d_ATM_his ['fract_sic'][0], dim1D='cell_latlon' ) 
     539        ATM_flic = lmdz.geo2point ( d_ATM_his ['fract_lic'][0], dim1D='cell_latlon' ) 
     540 
     541    if SRF_HIS == 'latlon' : 
     542        try : 
     543            SRF_lat  = lmdz.geo2point (   d_SRF_his ['lat']+0*d_SRF_his ['lon'], dim1D='cell_latlon' ) 
     544            SRF_lon  = lmdz.geo2point ( 0*d_SRF_his ['lat']+  d_SRF_his ['lon'], dim1D='cell_latlon' ) 
     545        except : 
     546            try :  
     547                SRF_lat  = lmdz.geo2point (   d_SRF_his ['lat_dom_out']+0*d_SRF_his ['lon_dom_out'], dim1D='cell_latlon' ) 
     548                SRF_lon  = lmdz.geo2point ( 0*d_SRF_his ['lat_dom_out']+  d_SRF_his ['lon_dom_out'], dim1D='cell_latlon' ) 
     549            except :  
     550                SRF_lat  = lmdz.geo2point (   d_SRF_his ['lat_domain_landpoints_out']+0*d_SRF_his ['lon_domain_landpoints_out'], dim1D='cell_latlon' ) 
     551                SRF_lon  = lmdz.geo2point ( 0*d_SRF_his ['lat_domain_landpoints_out']+  d_SRF_his ['lon_domain_landpoints_out'], dim1D='cell_latlon' ) 
     552        SRF_areas    = d_SRF_his ['Areas'].values 
     553        SRF_areafrac = d_SRF_his ['AreaFrac'].values 
     554        SRF_areas    = xr.DataArray ( SRF_areas   , coords=d_SRF_his ['Contfrac'].coords, dims=d_SRF_his ['Contfrac'].dims) 
     555        SRF_areafrac = xr.DataArray ( SRF_areafrac, coords=d_SRF_his ['Contfrac'].coords, dims=d_SRF_his ['Contfrac'].dims) 
     556        SRF_contfrac = SRF_areafrac / SRF_areas 
     557        #SRF_aire      = lmdz.geo2point (SRF_areas * SRF_contfrac, dim1D='cell_latlon', cumulPoles=True ) 
     558        #SRF_aire      = lmdz.geo2point (SRF_areas, dim1D='cell_latlon', cumulPoles=True ) 
     559        #SRF_aire      = lmdz.geo2point (d_ATM_aire ['aire'].squeeze().values * SRF_contfrac, dim1D='cell_latlon', cumulPoles=True ) 
     560        SRF_areas     = lmdz.geo2point ( SRF_areas   , dim1D='cell_latlon', cumulPoles=True ) 
     561        SRF_areafrac  = lmdz.geo2point ( SRF_areafrac, dim1D='cell_latlon', cumulPoles=True ) 
     562        SRF_contfrac  = lmdz.geo2point ( SRF_contfrac, dim1D='cell_latlon' ) 
     563        SRF_aire      = SRF_areafrac 
     564 
    516565    if ATM_HIS == 'ico' : 
    517566        echo ( f'Grille icosaedre' ) 
    518567        ATM_aire =  d_ATM_his ['aire']     [0][ATM_his_keysort].squeeze() 
    519         #ATM_lat  =  d_ATM_his ['lat']         [ATM_his_keysort] 
    520         #ATM_lon  =  d_ATM_his ['lon']         [ATM_his_keysort] 
     568        ATM_lat  =  d_ATM_his ['lat']         [ATM_his_keysort] 
     569        ATM_lon  =  d_ATM_his ['lon']         [ATM_his_keysort] 
    521570        ATM_fter =  d_ATM_his ['fract_ter'][0][ATM_his_keysort] 
    522571        ATM_foce =  d_ATM_his ['fract_oce'][0][ATM_his_keysort] 
    523572        ATM_fsic =  d_ATM_his ['fract_sic'][0][ATM_his_keysort] 
    524573        ATM_flic =  d_ATM_his ['fract_lic'][0][ATM_his_keysort] 
    525         #SRF_lat  =  d_SRF_his ['lat']         [SRF_his_keysort] 
    526         #SRF_lon  =  d_SRF_his ['lon']         [SRF_his_keysort] 
    527         SRF_aire =  d_SRF_his ['Areas'][SRF_his_keysort] * d_SRF_his ['Contfrac'][SRF_his_keysort] 
     574         
     575    if SRF_HIS == 'ico' :  
     576        SRF_lat       =  d_SRF_his ['lat'][SRF_his_keysort] 
     577        SRF_lon       =  d_SRF_his ['lon'][SRF_his_keysort] 
     578        SRF_aire      =  d_SRF_his ['Areas'][SRF_his_keysort] * d_SRF_his ['Contfrac'][SRF_his_keysort] 
     579        SRF_areas     =  d_SRF_his ['Areas'][SRF_his_keysort]  
     580        SRF_contfrac  =  d_SRF_his ['Contfrac'][SRF_his_keysort] 
    528581 
    529582ATM_fsea      = ATM_foce + ATM_fsic 
    … …  
    537590 
    538591SRF_aire = SRF_aire.where ( np.abs (SRF_aire) < 1E15, 0.) 
     592 
     593## Write the full configuration 
     594config_out = open (FullIniFile, 'w') 
     595config.write ( config_out ) 
     596config_out.close () 
    539597             
    540598if ICO : 
    541599    # Area on icosahedron grid 
    542     d_DYN_aire = xr.open_dataset ( file_DYN_aire, decode_times=False).squeeze() 
     600    d_DYN_aire = xr.open_dataset ( file_DYN_aire, decode_times=False ).squeeze() 
    543601 
    544602    if SortIco :  
    … …  
    552610 
    553611    DYN_aire = d_DYN_aire['aire'][DYN_aire_keysort] 
    554     DYN_fsea = d_DYN_aire ['fract_oce'][DYN_aire_keysort] + d_DYN_aire['fract_sic'][DYN_aire_keysort] 
     612    DYN_fsea = d_DYN_aire['fract_oce'][DYN_aire_keysort] + d_DYN_aire['fract_sic'][DYN_aire_keysort] 
    555613    DYN_flnd = 1.0 - DYN_fsea 
    556614    DYN_fter = d_ATM_beg['FTER'][ATM_beg_keysort] 
    … …  
    568626 
    569627# Functions computing integrals and sum 
    570 def ATM_stock_int (stock) : 
     628def DYN_stock_int (stock) : 
    571629    '''Integrate (* surface) stock on atmosphere grid''' 
    572     ATM_stock_int  = wu.Psum ( (stock * DYN_aire).to_masked_array().ravel() )  
    573     return ATM_stock_int 
     630    DYN_stock_int  = wu.Psum ( (stock * DYN_aire).to_masked_array().ravel() )  
     631    return DYN_stock_int 
    574632 
    575633def ATM_flux_int (flux) : 
    … …  
    580638def SRF_stock_int (stock) : 
    581639    '''Integrate (* surface) stock on land grid''' 
     640    #SRF_stock_int  = wu.Ksum (  ( (stock * DYN_aire_fter).to_masked_array().ravel()) ) # Marche avec ICO et lon/lat na 
    582641    SRF_stock_int  = wu.Ksum (  ( (stock * DYN_aire_fter).to_masked_array().ravel()) ) 
    583642    return SRF_stock_int 
    … …  
    598657    return ONE_flux_int 
    599658 
    600  
    601 #if LMDZ : 
    602 #    d_ATM_beg = d_ATM_beg.assign ( coords={'lon':d_ATM_beg.lon*180./np.pi} ) 
    603  
    604659ATM_aire_sea     = ATM_aire * ATM_fsea 
    605660ATM_aire_tot     = ONE_stock_int (ATM_aire) 
    606661ATM_aire_sea_tot = ONE_stock_int (ATM_aire*ATM_fsea) 
    607662 
    608  
    609 echo ( 'Area of atmosphere/(4pi R^2) : {:12.5f}'.format(ATM_aire_tot/(Ra*Ra*4*np.pi) ) ) 
     663DYN_aire_tot     =  ONE_stock_int ( DYN_aire ) 
     664SRF_aire_tot     =  ONE_stock_int ( SRF_aire ) 
     665 
     666echo ('') 
     667echo ( f'ATM DYN : Area of atmosphere : {DYN_aire_tot:18.9e}' ) 
     668echo ( f'ATM HIS : Area of atmosphere : {ATM_aire_tot:18.9e}' ) 
     669echo ( f'ATM SRF : Area of atmosphere : {SRF_aire_tot:18.9e}' ) 
     670echo ('') 
     671echo ( 'ATM DYN : Area of atmosphere/(4pi R^2) : {:18.9f}'.format(DYN_aire_tot/(Ra*Ra*4*np.pi) ) ) 
     672echo ( 'ATM HIS : Area of atmosphere/(4pi R^2) : {:18.9f}'.format(ATM_aire_tot/(Ra*Ra*4*np.pi) ) ) 
     673echo ( 'ATM SRF : Area of atmosphere/(4pi R^2) : {:18.9f}'.format(SRF_aire_tot/(Ra*Ra*4*np.pi) ) ) 
     674echo ('') 
     675echo ( f'ATM SRF : Area of atmosphere (no contfrac): {ONE_stock_int (SRF_areas):18.9e}' ) 
     676 
    610677 
    611678if (  np.abs (ATM_aire_tot/(Ra*Ra*4*np.pi) - 1.0) > 0.01 ) : 
    … …  
    627694    DYN_psol_end = d_DYN_end['ps'][DYN_end_keysort] 
    628695if LMDZ :  
    629     DYN_psol_beg = lmdz.geo2point ( d_DYN_beg['ps'].isel(rlonv=slice(0,-1)), dim1D='cell' ) 
    630     DYN_psol_end = lmdz.geo2point ( d_DYN_end['ps'].isel(rlonv=slice(0,-1)), dim1D='cell' ) 
     696    DYN_psol_beg = lmdz.geo2point ( d_DYN_beg['ps'].isel(rlonv=slice(0,-1)), dim1D='cell_latlon' ) 
     697    DYN_psol_end = lmdz.geo2point ( d_DYN_end['ps'].isel(rlonv=slice(0,-1)), dim1D='cell_latlon' ) 
    631698     
    632699echo ( '3D Pressure at the interface layers (not scalar points)' ) 
    … …  
    635702 
    636703klevp1 = ATM_Bhyb.shape[-1] 
    637 cell        = DYN_psol_beg.shape[-1] 
    638 klev = klevp1 - 1 
     704cell   = DYN_psol_beg.shape[-1] 
     705klev   = klevp1 - 1 
    639706 
    640707echo ( 'Layer thickness (pressure)' ) 
    641708DYN_dpres_beg = xr.DataArray ( np.empty( (klev, cell )), dims = ('sigs', 'cell'), coords=(np.arange(klev), np.arange(cell) ) ) 
    642709DYN_dpres_end = xr.DataArray ( np.empty( (klev, cell )), dims = ('sigs', 'cell'), coords=(np.arange(klev), np.arange(cell) ) ) 
    643 #DYN_dpres_beg = DYN_dpres_beg 
    644 #DYN_dpres_end = DYN_dpres_end 
    645710     
    646711for k in np.arange (klevp1-1) : 
    … …  
    669734     
    670735echo ( 'Compute water content : vertical and horizontal integral' ) 
    671 DYN_mas_wat_beg = ATM_stock_int ( (DYN_dpres_beg * DYN_wat_beg).sum (dim='sigs') ) / Grav 
    672 DYN_mas_wat_end = ATM_stock_int ( (DYN_dpres_end * DYN_wat_end).sum (dim='sigs') ) / Grav 
     736DYN_mas_wat_beg = DYN_stock_int ( (DYN_dpres_beg * DYN_wat_beg).sum (dim='sigs') ) / Grav 
     737DYN_mas_wat_end = DYN_stock_int ( (DYN_dpres_end * DYN_wat_end).sum (dim='sigs') ) / Grav 
    673738 
    674739echo ( 'Variation of water content' ) 
    … …  
    683748prtFlux ( 'dMass (atm)  ', dDYN_mas_wat, 'e', True ) 
    684749 
    685 ATM_sno_beg = d_ATM_beg['SNOW01']*d_ATM_beg['FTER'] + d_ATM_beg['SNOW02']*d_ATM_beg['FLIC'] + d_ATM_beg['SNOW03']*d_ATM_beg['FOCE'] + d_ATM_beg['SNOW04']*d_ATM_beg['FSIC'] 
    686 ATM_sno_end = d_ATM_end['SNOW01']*d_ATM_end['FTER'] + d_ATM_end['SNOW02']*d_ATM_end['FLIC'] + d_ATM_end['SNOW03']*d_ATM_end['FOCE'] + d_ATM_end['SNOW04']*d_ATM_end['FSIC'] 
    687  
    688 ATM_qs_beg  = d_ATM_beg['QS01']  *d_ATM_beg['FTER'] + d_ATM_beg['QS02']  *d_ATM_beg['FLIC'] + d_ATM_beg['QS03']  *d_ATM_beg['FOCE'] + d_ATM_beg['QS04']  *d_ATM_beg['FSIC'] 
    689 ATM_qs_end  = d_ATM_end['QS01']  *d_ATM_end['FTER'] + d_ATM_end['QS02']  *d_ATM_end['FLIC'] + d_ATM_end['QS03']  *d_ATM_end['FOCE'] + d_ATM_end['QS04']  *d_ATM_end['FSIC'] 
     750ATM_sno_beg = d_ATM_beg['SNOW01']*d_ATM_beg['FTER'] + d_ATM_beg['SNOW02']*d_ATM_beg['FLIC'] + \ 
     751              d_ATM_beg['SNOW03']*d_ATM_beg['FOCE'] + d_ATM_beg['SNOW04']*d_ATM_beg['FSIC'] 
     752ATM_sno_end = d_ATM_end['SNOW01']*d_ATM_end['FTER'] + d_ATM_end['SNOW02']*d_ATM_end['FLIC'] + \ 
     753              d_ATM_end['SNOW03']*d_ATM_end['FOCE'] + d_ATM_end['SNOW04']*d_ATM_end['FSIC'] 
     754 
     755ATM_qs_beg  = d_ATM_beg['QS01']  *d_ATM_beg['FTER'] + d_ATM_beg['QS02']  *d_ATM_beg['FLIC'] + \ 
     756              d_ATM_beg['QS03']  *d_ATM_beg['FOCE'] + d_ATM_beg['QS04']  *d_ATM_beg['FSIC'] 
     757ATM_qs_end  = d_ATM_end['QS01']  *d_ATM_end['FTER'] + d_ATM_end['QS02']  *d_ATM_end['FLIC'] + \ 
     758              d_ATM_end['QS03']  *d_ATM_end['FOCE'] + d_ATM_end['QS04']  *d_ATM_end['FSIC'] 
    690759 
    691760ATM_qsol_beg = d_ATM_beg['QSOL'] 
    … …  
    707776ATM_qs04_end  = d_ATM_end['QS04'] * d_ATM_end['FSIC']   
    708777 
    709 # if LMDZ : 
    710 #     ATM_sno_beg       = lmdz.geo2point ( ATM_sno_beg    , dim1D='cell' ) 
    711 #     ATM_sno_end       = lmdz.geo2point ( ATM_sno_end    , dim1D='cell' ) 
    712 #     ATM_qs_beg        = lmdz.geo2point ( ATM_qs_beg     , dim1D='cell' ) 
    713 #     ATM_qs_end        = lmdz.geo2point ( ATM_qs_end     , dim1D='cell' ) 
    714 #     ATM_qsol_beg      = lmdz.geo2point ( ATM_qsol_beg   , dim1D='cell' ) 
    715 #     ATM_qsol_end      = lmdz.geo2point ( ATM_qsol_end   , dim1D='cell' ) 
    716 #     ATM_qs01_beg      = lmdz.geo2point ( ATM_qs01_beg   , dim1D='cell' ) 
    717 #     ATM_qs01_end      = lmdz.geo2point ( ATM_qs01_end   , dim1D='cell' ) 
    718 #     ATM_qs02_beg      = lmdz.geo_point ( ATM_qs02_beg   , dim1D='cell' ) 
    719 #     ATM_qs02_end      = lmdz.geo2point ( ATM_qs02_end   , dim1D='cell' ) 
    720 #     ATM_qs03_beg      = lmdz.geo2point ( ATM_qs03_beg   , dim1D='cell' ) 
    721 #     ATM_qs03_end      = lmdz.geo2point ( ATM_qs03_end   , dim1D='cell' ) 
    722 #     ATM_qs04_beg      = lmdz.geo2point ( ATM_qs04_beg   , dim1D='cell' ) 
    723 #     ATM_qs04_end      = lmdz.geo2point ( ATM_qs04_end   , dim1D='cell' ) 
    724 #     LIC_sno_beg       = lmdz.geo2point ( LIC_sno_beg    , dim1D='cell' ) 
    725 #     LIC_sno_end       = lmdz.geo2point ( LIC_sno_end    , dim1D='cell' ) 
    726 #     LIC_runlic0_beg   = lmdz.geo2point ( LIC_runlic0_beg, dim1D='cell' ) 
    727 #     LIC_runlic0_end   = lmdz.geo2point ( LIC_runlic0_end, dim1D='cell' ) 
    728778if ICO : 
    729779     ATM_sno_beg       = ATM_sno_beg    [ATM_beg_keysort] 
    … …  
    745795     LIC_runlic0_beg   = LIC_runlic0_beg[ATM_beg_keysort] 
    746796     LIC_runlic0_end   = LIC_runlic0_end[ATM_end_keysort] 
    747  
    748  
    749797    
    750798LIC_qs_beg = ATM_qs02_beg 
    751799LIC_qs_end = ATM_qs02_end 
    752800 
    753 ATM_mas_sno_beg     = ATM_stock_int ( ATM_sno_beg  ) 
    754 ATM_mas_sno_end     = ATM_stock_int ( ATM_sno_end  ) 
    755  
    756 ATM_mas_qs_beg      = ATM_stock_int ( ATM_qs_beg   ) 
    757 ATM_mas_qs_end      = ATM_stock_int ( ATM_qs_end   ) 
    758 ATM_mas_qsol_beg    = ATM_stock_int ( ATM_qsol_beg ) 
    759 ATM_mas_qsol_end    = ATM_stock_int ( ATM_qsol_end ) 
    760 ATM_mas_qs01_beg    = ATM_stock_int ( ATM_qs01_beg ) 
    761 ATM_mas_qs01_end    = ATM_stock_int ( ATM_qs01_end ) 
    762 ATM_mas_qs02_beg    = ATM_stock_int ( ATM_qs02_beg ) 
    763 ATM_mas_qs02_end    = ATM_stock_int ( ATM_qs02_end ) 
    764 ATM_mas_qs03_beg    = ATM_stock_int ( ATM_qs03_beg ) 
    765 ATM_mas_qs03_end    = ATM_stock_int ( ATM_qs03_end ) 
    766 ATM_mas_qs04_beg    = ATM_stock_int ( ATM_qs04_beg ) 
    767 ATM_mas_qs04_end    = ATM_stock_int ( ATM_qs04_end ) 
    768  
    769 LIC_mas_sno_beg     = ATM_stock_int ( LIC_sno_beg     ) 
    770 LIC_mas_sno_end     = ATM_stock_int ( LIC_sno_end     ) 
    771 LIC_mas_runlic0_beg = ATM_stock_int ( LIC_runlic0_beg ) 
    772 LIC_mas_runlic0_end = ATM_stock_int ( LIC_runlic0_end ) 
     801ATM_mas_sno_beg     = DYN_stock_int ( ATM_sno_beg  ) 
     802ATM_mas_sno_end     = DYN_stock_int ( ATM_sno_end  ) 
     803 
     804ATM_mas_qs_beg      = DYN_stock_int ( ATM_qs_beg   ) 
     805ATM_mas_qs_end      = DYN_stock_int ( ATM_qs_end   ) 
     806ATM_mas_qsol_beg    = DYN_stock_int ( ATM_qsol_beg ) 
     807ATM_mas_qsol_end    = DYN_stock_int ( ATM_qsol_end ) 
     808ATM_mas_qs01_beg    = DYN_stock_int ( ATM_qs01_beg ) 
     809ATM_mas_qs01_end    = DYN_stock_int ( ATM_qs01_end ) 
     810ATM_mas_qs02_beg    = DYN_stock_int ( ATM_qs02_beg ) 
     811ATM_mas_qs02_end    = DYN_stock_int ( ATM_qs02_end ) 
     812ATM_mas_qs03_beg    = DYN_stock_int ( ATM_qs03_beg ) 
     813ATM_mas_qs03_end    = DYN_stock_int ( ATM_qs03_end ) 
     814ATM_mas_qs04_beg    = DYN_stock_int ( ATM_qs04_beg ) 
     815ATM_mas_qs04_end    = DYN_stock_int ( ATM_qs04_end ) 
     816 
     817LIC_mas_sno_beg     = DYN_stock_int ( LIC_sno_beg     ) 
     818LIC_mas_sno_end     = DYN_stock_int ( LIC_sno_end     ) 
     819LIC_mas_runlic0_beg = DYN_stock_int ( LIC_runlic0_beg ) 
     820LIC_mas_runlic0_end = DYN_stock_int ( LIC_runlic0_end ) 
    773821 
    774822LIC_mas_qs_beg  = ATM_mas_qs02_beg 
    … …  
    826874 
    827875if Routing == 'SECHIBA' : 
    828         RUN_mas_wat_fast_beg   = ONE_stock_int ( d_SRF_beg ['fastres']   ) 
    829         RUN_mas_wat_slow_beg   = ONE_stock_int ( d_SRF_beg ['slowres']   ) 
    830         RUN_mas_wat_stream_beg = ONE_stock_int ( d_SRF_beg ['streamres'] ) 
    831         RUN_mas_wat_flood_beg  = ONE_stock_int ( d_SRF_beg ['floodres']  ) 
    832         RUN_mas_wat_lake_beg   = ONE_stock_int ( d_SRF_beg ['lakeres']   ) 
    833         RUN_mas_wat_pond_beg   = ONE_stock_int ( d_SRF_beg ['pondres']   ) 
    834          
    835         RUN_mas_wat_fast_end   = ONE_stock_int ( d_SRF_end ['fastres']   ) 
    836         RUN_mas_wat_slow_end   = ONE_stock_int ( d_SRF_end ['slowres']   ) 
    837         RUN_mas_wat_stream_end = ONE_stock_int ( d_SRF_end ['streamres'] ) 
    838         RUN_mas_wat_flood_end  = ONE_stock_int ( d_SRF_end ['floodres']  ) 
    839         RUN_mas_wat_lake_end   = ONE_stock_int ( d_SRF_end ['lakeres']   ) 
    840         RUN_mas_wat_pond_end   = ONE_stock_int ( d_SRF_end ['pondres']   ) 
    841  
    842 RUN_mas_wat_beg = RUN_mas_wat_fast_beg + RUN_mas_wat_slow_beg + RUN_mas_wat_stream_beg + RUN_mas_wat_flood_beg + RUN_mas_wat_lake_beg + RUN_mas_wat_pond_beg 
    843 RUN_mas_wat_end = RUN_mas_wat_fast_end + RUN_mas_wat_slow_end + RUN_mas_wat_stream_end + RUN_mas_wat_flood_end + RUN_mas_wat_lake_end + RUN_mas_wat_pond_end 
     876    RUN_mas_wat_fast_beg   = ONE_stock_int ( d_SRF_beg ['fastres']   ) 
     877    RUN_mas_wat_slow_beg   = ONE_stock_int ( d_SRF_beg ['slowres']   ) 
     878    RUN_mas_wat_stream_beg = ONE_stock_int ( d_SRF_beg ['streamres'] ) 
     879    RUN_mas_wat_flood_beg  = ONE_stock_int ( d_SRF_beg ['floodres']  ) 
     880    RUN_mas_wat_lake_beg   = ONE_stock_int ( d_SRF_beg ['lakeres']   ) 
     881    RUN_mas_wat_pond_beg   = ONE_stock_int ( d_SRF_beg ['pondres']   ) 
     882     
     883    RUN_mas_wat_fast_end   = ONE_stock_int ( d_SRF_end ['fastres']   ) 
     884    RUN_mas_wat_slow_end   = ONE_stock_int ( d_SRF_end ['slowres']   ) 
     885    RUN_mas_wat_stream_end = ONE_stock_int ( d_SRF_end ['streamres'] ) 
     886    RUN_mas_wat_flood_end  = ONE_stock_int ( d_SRF_end ['floodres']  ) 
     887    RUN_mas_wat_lake_end   = ONE_stock_int ( d_SRF_end ['lakeres']   ) 
     888    RUN_mas_wat_pond_end   = ONE_stock_int ( d_SRF_end ['pondres']   ) 
     889 
     890RUN_mas_wat_beg = RUN_mas_wat_fast_beg  + RUN_mas_wat_slow_beg + RUN_mas_wat_stream_beg + \ 
     891                  RUN_mas_wat_flood_beg + RUN_mas_wat_lake_beg + RUN_mas_wat_pond_beg 
     892RUN_mas_wat_end = RUN_mas_wat_fast_end  + RUN_mas_wat_slow_end + RUN_mas_wat_stream_end + \ 
     893                  RUN_mas_wat_flood_end + RUN_mas_wat_lake_end + RUN_mas_wat_pond_end 
    844894 
    845895dRUN_mas_wat_fast   = RUN_mas_wat_fast_end   - RUN_mas_wat_fast_beg 
    … …  
    874924echo ( f'\nWater content in routing  -- {Title} ' ) 
    875925echo ( '------------------------------------------------------------------------------------' ) 
    876 echo ( 'RUN_mas_wat_beg = {:12.6e} kg | RUN_mas_wat_end = {:12.6e} kg '.format (RUN_mas_wat_end, RUN_mas_wat_end) ) 
     926echo ( f'RUN_mas_wat_beg = {RUN_mas_wat_end:12.6e} kg | RUN_mas_wat_end = {RUN_mas_wat_end:12.6e} kg' ) 
    877927prtFlux ( 'dMass (routing) ', dRUN_mas_wat , 'e', True   ) 
    878928 
    … …  
    890940 
    891941if LMDZ : 
    892     SRF_tot_watveg_beg  = lmdz.geo2point (SRF_tot_watveg_beg , dim1D='cell') 
    893     SRF_tot_watsoil_beg = lmdz.geo2point (SRF_tot_watsoil_beg, dim1D='cell') 
    894     SRF_snow_beg        = lmdz.geo2point (SRF_snow_beg       , dim1D='cell') 
    895     SRF_lakeres_beg     = lmdz.geo2point (SRF_lakeres_beg    , dim1D='cell') 
    896     SRF_tot_watveg_end  = lmdz.geo2point (SRF_tot_watveg_end , dim1D='cell') 
    897     SRF_tot_watsoil_end = lmdz.geo2point (SRF_tot_watsoil_end, dim1D='cell') 
    898     SRF_snow_end        = lmdz.geo2point (SRF_snow_end       , dim1D='cell') 
    899     SRF_lakeres_end     = lmdz.geo2point (SRF_lakeres_end    , dim1D='cell') 
     942    SRF_tot_watveg_beg  = lmdz.geo2point (SRF_tot_watveg_beg , dim1D='cell_latlon') 
     943    SRF_tot_watsoil_beg = lmdz.geo2point (SRF_tot_watsoil_beg, dim1D='cell_latlon') 
     944    SRF_snow_beg        = lmdz.geo2point (SRF_snow_beg       , dim1D='cell_latlon') 
     945    SRF_lakeres_beg     = lmdz.geo2point (SRF_lakeres_beg    , dim1D='cell_latlon') 
     946    SRF_tot_watveg_end  = lmdz.geo2point (SRF_tot_watveg_end , dim1D='cell_latlon') 
     947    SRF_tot_watsoil_end = lmdz.geo2point (SRF_tot_watsoil_end, dim1D='cell_latlon') 
     948    SRF_snow_end        = lmdz.geo2point (SRF_snow_end       , dim1D='cell_latlon') 
     949    SRF_lakeres_end     = lmdz.geo2point (SRF_lakeres_end    , dim1D='cell_latlon') 
    900950   
    901951if ICO : 
    … …  
    942992dSRF_mas_lake     = SRF_mas_lake_end    - SRF_mas_lake_beg 
    943993 
    944  
    945994echo ( '------------------------------------------------------------------------------------' ) 
    946995echo ( f'\nSurface reservoirs  -- {Title} ') 
    … …  
    950999echo ( f'SRF_mas_lake_beg     = {SRF_mas_lake_beg   :12.6e} kg | SRF_mas_lake_end     = {SRF_mas_lake_end   :12.6e} kg ' ) 
    9511000 
    952 prtFlux ('dMass (watveg) ', dSRF_mas_watveg , 'e' , True) 
    953 prtFlux ('dMass (watsoil)', dSRF_mas_watsoil, 'e' , True) 
    954 prtFlux ('dMass (snow)   ', dSRF_mas_snow   , 'e' , True) 
    955 prtFlux ('dMass (lake)   ', dSRF_mas_lake   , 'e' , True) 
     1001prtFlux ( 'dMass (watveg) ', dSRF_mas_watveg , 'e' , True ) 
     1002prtFlux ( 'dMass (watsoil)', dSRF_mas_watsoil, 'e' , True ) 
     1003prtFlux ( 'dMass (snow)   ', dSRF_mas_snow   , 'e' , True ) 
     1004prtFlux ( 'dMass (lake)   ', dSRF_mas_lake   , 'e' , True ) 
    9561005 
    9571006SRF_mas_wat_beg = SRF_mas_watveg_beg + SRF_mas_watsoil_beg + SRF_mas_snow_beg  
    … …  
    9761025print ( 'Step 1', end='' ) 
    9771026if ATM_HIS == 'latlon' : 
    978     ATM_wbilo_oce   = lmdz.geo2point ( d_ATM_his ['wbilo_oce'], dim1D='cell' ) 
    979     ATM_wbilo_sic   = lmdz.geo2point ( d_ATM_his ['wbilo_sic'], dim1D='cell' ) 
    980     ATM_wbilo_ter   = lmdz.geo2point ( d_ATM_his ['wbilo_ter'], dim1D='cell' ) 
    981     ATM_wbilo_lic   = lmdz.geo2point ( d_ATM_his ['wbilo_lic'], dim1D='cell' ) 
    982     ATM_runofflic   = lmdz.geo2point ( d_ATM_his ['runofflic'], dim1D='cell' ) 
    983     ATM_fqcalving   = lmdz.geo2point ( d_ATM_his ['fqcalving'], dim1D='cell' ) 
    984     ATM_fqfonte     = lmdz.geo2point ( d_ATM_his ['fqfonte']  , dim1D='cell' ) 
    985     ATM_precip      = lmdz.geo2point ( d_ATM_his ['precip']   , dim1D='cell' ) 
    986     ATM_snowf       = lmdz.geo2point ( d_ATM_his ['snow']     , dim1D='cell' ) 
    987     ATM_evap        = lmdz.geo2point ( d_ATM_his ['evap']     , dim1D='cell' ) 
    988     ATM_wevap_ter   = lmdz.geo2point ( d_ATM_his ['wevap_ter'], dim1D='cell' ) 
    989     ATM_wevap_oce   = lmdz.geo2point ( d_ATM_his ['wevap_oce'], dim1D='cell' ) 
    990     ATM_wevap_lic   = lmdz.geo2point ( d_ATM_his ['wevap_lic'], dim1D='cell' ) 
    991     ATM_wevap_sic   = lmdz.geo2point ( d_ATM_his ['wevap_sic'], dim1D='cell' ) 
    992     ATM_runofflic   = lmdz.geo2point ( d_ATM_his ['runofflic'], dim1D='cell' ) 
     1027    ATM_wbilo_oce   = lmdz.geo2point ( d_ATM_his ['wbilo_oce'], dim1D='cell_latlon' ) 
     1028    ATM_wbilo_sic   = lmdz.geo2point ( d_ATM_his ['wbilo_sic'], dim1D='cell_latlon' ) 
     1029    ATM_wbilo_ter   = lmdz.geo2point ( d_ATM_his ['wbilo_ter'], dim1D='cell_latlon' ) 
     1030    ATM_wbilo_lic   = lmdz.geo2point ( d_ATM_his ['wbilo_lic'], dim1D='cell_latlon' ) 
     1031    ATM_runofflic   = lmdz.geo2point ( d_ATM_his ['runofflic'], dim1D='cell_latlon' ) 
     1032    ATM_fqcalving   = lmdz.geo2point ( d_ATM_his ['fqcalving'], dim1D='cell_latlon' ) 
     1033    ATM_fqfonte     = lmdz.geo2point ( d_ATM_his ['fqfonte']  , dim1D='cell_latlon' ) 
     1034    ATM_precip      = lmdz.geo2point ( d_ATM_his ['precip']   , dim1D='cell_latlon' ) 
     1035    ATM_snowf       = lmdz.geo2point ( d_ATM_his ['snow']     , dim1D='cell_latlon' ) 
     1036    ATM_evap        = lmdz.geo2point ( d_ATM_his ['evap']     , dim1D='cell_latlon' ) 
     1037    ATM_wevap_ter   = lmdz.geo2point ( d_ATM_his ['wevap_ter'], dim1D='cell_latlon' ) 
     1038    ATM_wevap_oce   = lmdz.geo2point ( d_ATM_his ['wevap_oce'], dim1D='cell_latlon' ) 
     1039    ATM_wevap_lic   = lmdz.geo2point ( d_ATM_his ['wevap_lic'], dim1D='cell_latlon' ) 
     1040    ATM_wevap_sic   = lmdz.geo2point ( d_ATM_his ['wevap_sic'], dim1D='cell_latlon' ) 
     1041    ATM_runofflic   = lmdz.geo2point ( d_ATM_his ['runofflic'], dim1D='cell_latlon' ) 
    9931042     
    9941043print ( ' - Step 2', end='' ) 
    … …  
    10241073print ( ' - Step 4', end='' ) 
    10251074if RUN_HIS == 'latlon' :  
    1026     RUN_coastalflow = lmdz.geo2point ( d_RUN_his ['coastalflow'], dim1D='cell' )  
    1027     RUN_riverflow   = lmdz.geo2point ( d_RUN_his ['riverflow']  , dim1D='cell' ) 
    1028     RUN_runoff      = lmdz.geo2point ( d_RUN_his ['runoff']     , dim1D='cell' ) 
    1029     RUN_drainage    = lmdz.geo2point ( d_RUN_his ['drainage']   , dim1D='cell' ) 
    1030     RUN_riversret   = lmdz.geo2point ( d_RUN_his ['riversret']  , dim1D='cell' ) 
    1031      
    1032     RUN_coastalflow_cpl = lmdz.geo2point ( d_RUN_his ['coastalflow_cpl'], dim1D='cell' )  
    1033     RUN_riverflow_cpl   = lmdz.geo2point ( d_RUN_his ['riverflow_cpl']  , dim1D='cell' ) 
     1075    RUN_coastalflow = lmdz.geo2point ( d_RUN_his ['coastalflow'], dim1D='cell_latlon' )  
     1076    RUN_riverflow   = lmdz.geo2point ( d_RUN_his ['riverflow']  , dim1D='cell_latlon' ) 
     1077    RUN_runoff      = lmdz.geo2point ( d_RUN_his ['runoff']     , dim1D='cell_latlon' ) 
     1078    RUN_drainage    = lmdz.geo2point ( d_RUN_his ['drainage']   , dim1D='cell_latlon' ) 
     1079    RUN_riversret   = lmdz.geo2point ( d_RUN_his ['riversret']  , dim1D='cell_latlon' ) 
     1080     
     1081    RUN_coastalflow_cpl = lmdz.geo2point ( d_RUN_his ['coastalflow_cpl'], dim1D='cell_latlon' )  
     1082    RUN_riverflow_cpl   = lmdz.geo2point ( d_RUN_his ['riverflow_cpl']  , dim1D='cell_latlon' ) 
    10341083 
    10351084print ( ' - Step 5', end='' ) 
    … …  
    10461095print ( ' - Step 6', end='' ) 
    10471096if SRF_HIS == 'latlon' : 
    1048     SRF_rain     = lmdz.geo2point ( d_SRF_his ['rain']  , dim1D='cell') 
    1049     SRF_evap     = lmdz.geo2point ( d_SRF_his ['evap']  , dim1D='cell') 
    1050     SRF_snowf    = lmdz.geo2point ( d_SRF_his ['snowf'] , dim1D='cell') 
    1051     SRF_subli    = lmdz.geo2point ( d_SRF_his ['subli'] , dim1D='cell') 
    1052     SRF_transpir = lmdz.geo2point ( np.sum (d_SRF_his ['transpir'], axis=1), dim1D='cell' ) 
     1097    SRF_rain     = lmdz.geo2point ( d_SRF_his ['rain']  , dim1D='cell_latlon') 
     1098    SRF_evap     = lmdz.geo2point ( d_SRF_his ['evap']  , dim1D='cell_latlon') 
     1099    SRF_snowf    = lmdz.geo2point ( d_SRF_his ['snowf'] , dim1D='cell_latlon') 
     1100    SRF_subli    = lmdz.geo2point ( d_SRF_his ['subli'] , dim1D='cell_latlon') 
     1101    SRF_transpir = lmdz.geo2point ( np.sum (d_SRF_his ['transpir'], axis=1), dim1D='cell_latlon' ) 
    10531102 
    10541103print ( '- Step 7', end='' ) 
    … …  
    11401189prtFlux ('precip        ', ATM_flx_precip     , 'f' )        
    11411190prtFlux ('snowf         ', ATM_flx_snowf      , 'f' )         
    1142 prtFlux ('evap          ', ATM_flx_evap       , 'f' )          
     1191prtFlux ('evap          ', ATM_flx_evap       , 'f' ) 
     1192prtFlux ('runoff lic    ', ATM_flx_runlic     , 'f' )    
     1193 
     1194echo ( '\n====================================================================================' ) 
     1195echo ( f'-- RUNOFF Fluxes  -- {Title} ' ) 
    11431196prtFlux ('coastalflow   ', RUN_flx_coastal    , 'f' )  
    11441197prtFlux ('riverflow     ', RUN_flx_river      , 'f' )         
    … …  
    11521205prtFlux ('river out     ', RUN_flx_output     , 'f' )    
    11531206prtFlux ('river bil     ', RUN_flx_bil        , 'f' )    
    1154 prtFlux ('runoff lic    ', ATM_flx_runlic     , 'f' )    
    11551207 
    11561208ATM_flx_budget = -ATM_flx_wbilo + ATM_flx_calving + ATM_flx_runlic #+ ATM_flx_fqfonte + RUN_flx_river 
    … …  
    12771329echo ( 'Routing error : {:12.3e} (rel)'.format ( (RUN_flx_all-dRUN_mas_wat)/dRUN_mas_wat ) ) 
    12781330 
    1279 # PRINT *,"routing water Balance ;  before : ", water_balance_before," ; after : ",water_balance_after,  &a 
    1280 # 1150              " ; delta : ", 100*(water_balance_after-water_balance_before)/(0.5*(water_balance_after+water_balance_before)),"%" 
    1281  
    12821331echo ( ' ' ) 
    12831332echo ( f'{Title = }' ) 
     1333 
     1334echo ( 'SVN Information' ) 
     1335for clef in SVN.keys () : echo ( SVN[clef] ) 

  • TOOLS/WATER_BUDGET/WaterUtils.py

    • Property svn:keywords set to Date Revision HeadURL Author
    r6508 r6518  
    1515## 
    1616## SVN information 
    17 #  $Author: omamce $ 
    18 #  $Date: 2022-12-08 10:24:05 +0100 (Thu, 08 Dec 2022) $ 
    19 #  $Revision: 6277 $ 
     17#  $Author$ 
     18#  $Date$ 
     19#  $Revision$ 
    2020#  $Id: ATM_waterbudget.py 6277 2022-12-08 09:24:05Z omamce $ 
    21 #  $HeadURL: svn+ssh://omamce@forge.ipsl.jussieu.fr/ipsl/forge/projets/igcmg/svn/TOOLS/WATER_BUDGET/ATM_waterbudget.py $ 
     21#  $HeadURL$ 
    2222 
    2323## Import system modules 
Note: See TracChangeset for help on using the changeset viewer.