#!/usr/bin/env python3 ''' Script to check water conservation in the IPSL coupled model Here, some common utilities to all scripts Warning, to install, configure, run, use any of included software or to read the associated documentation you'll need at least one (1) brain in a reasonably working order. Lack of this implement will void any warranties (either express or implied). Authors assumes no responsability for errors, omissions, data loss, or any other consequences caused directly or indirectly by the usage of his software by incorrectly or partially configured personal SVN information $Author$ $Date$ $Revision$ $Id: $ $HeadURL$ ''' import os, sys, types import functools import time import configparser, re import numpy as np import xarray as xr import libIGCM_sys, libIGCM_date as ldate import nemo import lmdz def ReadConfig ( ini_file, default_ini_file='defaults.ini' ) : ''' Reads experiment parameters Reads config file to set all defaults parameters Reads config file to set all experiment parameters Reads config.card files if specified in Returns a dictionary with all parameters ''' ## Read file with defaults values ## ------------------------------ params = configparser.ConfigParser ( interpolation=configparser.ExtendedInterpolation() ) params.optionxform = str # To keep capitals params.read (default_ini_file) print ( f'\nConfig file readed : {default_ini_file} ' ) ## Read Experiment config file ## ---------------------------- exp_params = configparser.ConfigParser ( interpolation=configparser.ExtendedInterpolation() ) exp_params.optionxform = str # To keep capitals exp_params.read (ini_file) print ( f'\nConfig file readed : {ini_file} ' ) ## Reading config.card if possible ## ------------------------------- if 'Experiment' in params.keys () : ## Read Experiment on parameter file if possible if 'ConfigCard' in exp_params['Experiment'].keys () : ConfigCard = exp_params['Experiment']['ConfigCard'] print ( f'{ConfigCard=}' ) if ConfigCard : ## Read config card if it exists # Text existence of ConfigCard if os.path.exists ( ConfigCard ) : print ( f'Reading Config Card : {ConfigCard}' ) ## Creates parser for reading .ini input file config = configparser.ConfigParser (interpolation=configparser.ExtendedInterpolation() ) config.optionxform = str # To keep capitals config.read (ConfigCard) for VarName in ['JobName', 'ExperimentName', 'SpaceName', 'LongName', 'ModelName', 'TagName', 'ORCA_version', 'ExpType', 'PeriodLength', 'ResolAtm', 'ResolOce', 'CalendarType' ] : if VarName in config['UserChoices'].keys() : exp_params['Experiment'][VarName] = config['UserChoices'][VarName] if 'Post' in config.sections() : if 'PackFrequency' in config['Post'] : PackFrequency = config['Post']['PackFrequency'] exp_params['Experiment']['PackFrequency'] = PackFrequency else : raise FileExistsError ( f"File not found : {ConfigCard = }" ) ## Reading config file ## ------------------- for Section in exp_params.keys () : params[Section].update ( exp_params[Section] ) print ( f'\nConfig file readed : {ini_file} ' ) return config2dict (params) def SetDatesAndFiles ( pdict ) : ''' From readed experiment parameters, set all needed parameters ''' ## Set libIGCM and machine dependant values ## ---------------------------------------- pdict['Experiment']['ICO'] = 'ICO' in pdict['Experiment']['ResolAtm'] #pdict['Experiment']['LMDZ'] = 'LMD' in pdict['Experiment']['ResolAtm'] pdict['Experiment']['LMDZ'] = not pdict['Experiment']['ICO'] mm = libIGCM_sys.config ( JobName = pdict['Experiment']['JobName'], ExperimentName = pdict['Experiment']['ExperimentName'], TagName = pdict['Experiment']['TagName'], SpaceName = pdict['Experiment']['SpaceName'], User = pdict['Experiment']['User'], Group = pdict['Experiment']['Group'], TmpDir = pdict['Files']['TmpDir'], ARCHIVE = pdict['libIGCM']['ARCHIVE'], SCRATCHDIR = pdict['libIGCM']['SCRATCHDIR'], STORAGE = pdict['libIGCM']['STORAGE'], R_IN = pdict['libIGCM']['R_IN'], R_OUT = pdict['libIGCM']['R_OUT'], R_FIG = pdict['libIGCM']['R_FIG'], R_SAVE = pdict['libIGCM']['R_SAVE'], R_FIGR = pdict['libIGCM']['R_FIGR'], R_BUFR = pdict['libIGCM']['R_BUFR'], R_BUF_KSH = pdict['libIGCM']['R_BUF_KSH'], POST_DIR = pdict['libIGCM']['POST_DIR'], L_EXP = pdict['libIGCM']['L_EXP'] ) pdict['Files']['TmpDir'] = mm['TmpDir'] pdict['libIGCM'].update (mm) ## Complete configuration ## ---------------------- if not pdict['Experiment']['NEMO'] and pdict['Experiment']['ORCA_version'] : if '1.2' in pdict['Experiment']['ORCA_version'] : pdict['Experiment']['NEMO'] = 3.6 if '4.0' in pdict['Experiment']['ORCA_version'] : pdict['Experiment']['NEMO'] = 4.0 if '4.2' in pdict['Experiment']['ORCA_version'] : pdict['Experiment']['NEMO'] = 4.2 if pdict['Experiment']['ATM_HIS'] : if not pdict['Experiment']['SRF_HIS'] : pdict['Experiment']['SRF_HIS'] = pdict['Experiment']['ATM_HIS'] if not pdict['Experiment']['RUN_HIS'] : pdict['Experiment']['RUN_HIS'] = pdict['Experiment']['ATM_HIS'] ## ## Reading prec if not pdict['Config']['readPrec'] : pdict['Config']['readPrec'] = np.float64 else : if pdict['Config']['readPrec'] in ["float", "float64", "r8", "double", "" ] : pdict['Config']['readPrec'] = float if pdict['Config']['readPrec'] in [ "float32", "r4", "single", "" ] : pdict['Config']['readPrec'] = np.float32 if pdict['Config']['readPrec'] in [ "float16", "r2", "half" , "" ] : pdict['Config']['readPrec'] = np.float16 ## Defines begining and end of experiment ## -------------------------------------- if not pdict['Experiment']['DateBegin'] : pdict['Experiment']['DateBegin'] = f"{pdict['Experiment']['YearBegin']}0101" else : YearBegin, MonthBegin, DayBegin = ldate.GetYearMonthDay ( pdict['Experiment']['DateBegin'] ) DateBegin = FormatToGregorian (pdict['Experiment']['DateBegin']) pdict['Experiment']['YearBegin'] = YearBegin if not pdict['Experiment']['DateEnd'] : pdict['Experiment']['DateEnd'] = f"{pdict['Experiment']['YearEnd']}1231" else : YearEnd, MonthEnd, DayEnd = ldate.GetYearMonthDay ( pdict['Experiment']['DateEnd'] ) DateEnd = FormatToGregorian (pdict['Experiment']['DateEnd']) pdict['Experiment']['DateEnd'] = DateEnd if not pdict['Experiment']['PackFrequency'] : PackFrequencypdict['Experiment']['PackFrequency'] = f"{YearEnd - YearBegin + 1}YE" ## Set directory to extract files ## ------------------------------ if not pdict['Files']['FileDir'] : pdict['Files']['FileDir'] = os.path.join ( pdict['Files']['TmpDir'], f"WATER_{pdict['Experiment']['JobName']}" ) if not os.path.isdir ( pdict['Files']['FileDir'] ) : os.makedirs ( pdict['Files']['FileDir'] ) FileOut = str2value (pdict['Files']['FileOut']) if not FileOut : FileOut = f"ATM_waterbudget_{pdict['Experiment']['JobName']}_{pdict['Experiment']['YearBegin']}_{pdict['Experiment']['YearEnd']}" if pdict['Experiment']['ICO'] : if pdict['Experiment']['ATM_HIS'] == 'latlon' : FileOut = f'{FileOut}_LATLON' if pdict['Experiment']['ATM_HIS'] == 'ico' : FileOut = f'{FileOut}_ICO' if pdict['Config']['readPrec'] == np.float32 : FileOut = f'{FileOut}_float32' FileOut = f'{FileOut}.out' pdict['Files']['FileOut'] = FileOut f_out = open ( FileOut, mode = 'w', encoding="utf-8" ) pdict['Files']['f_out'] = f_out echo (' ', f_out) echo ( f"JobName : {pdict['Experiment']['JobName']}" , f_out=f_out ) echo ( f"Comment : {pdict['Experiment']['Comment']}" , f_out=f_out ) echo ( f"TmpDir : {pdict['Files']['TmpDir']}" , f_out=f_out ) echo ( f"FileDir : {pdict['Files']['FileDir']}" , f_out=f_out ) echo ( f"\nDealing with {pdict['libIGCM']['L_EXP']}", f_out ) ## Creates model output directory names ## ------------------------------------ if not pdict['Files']['FreqDir'] : if pdict['Experiment']['Freq'] == "MO" : pdict['Files']['FreqDir'] = os.path.join ('Output' , 'MO' ) if pdict['Experiment']['Freq'] == "SE" : pdict['Files']['FreqDir'] = os.path.join ('Analyse', 'SE' ) if not pdict['Files']['dir_ATM_his'] : pdict['Files']['dir_ATM_his'] = os.path.join ( pdict['libIGCM']['R_SAVE'], "ATM", pdict['Files']['FreqDir'] ) if not pdict['Files']['dir_SRF_his'] : pdict['Files']['dir_SRF_his'] = os.path.join ( pdict['libIGCM']['R_SAVE'], "SRF", pdict['Files']['FreqDir'] ) if not pdict['Files']['dir_OCE_his'] : pdict['Files']['dir_OCE_his'] = os.path.join ( pdict['libIGCM']['R_SAVE'], "OCE", pdict['Files']['FreqDir'] ) if not pdict['Files']['dir_ICE_his'] : pdict['Files']['dir_ICE_his'] = os.path.join ( pdict['libIGCM']['R_SAVE'], "ICE", pdict['Files']['FreqDir'] ) echo ( 'The analysis relies on files from the following model output directories : ', f_out ) echo ( f"{pdict['Files']['dir_ATM_his'] = }" , f_out ) echo ( f"{pdict['Files']['dir_SRF_his'] = }" , f_out ) echo ( f"{pdict['Files']['dir_OCE_his'] = }" , f_out ) echo ( f"{pdict['Files']['dir_ICE_his'] = }" , f_out ) ## Creates files names ## -------------------- if not pdict['Experiment']['Period'] : if pdict['Experiment']['Freq'] == 'MO' : pdict['Experiment']['Period'] = f"{pdict['Experiment']['DateBegin']}_{pdict['Experiment']['DateEnd']}_1M" if pdict['Experiment']['Freq'] == 'SE' : pdict['Experiment']['Period'] = f"SE_{pdict['Files']['DateBegin']}_{pdict['Files']['DateEnd']}_1M" echo ( f"Period : {pdict['Experiment']['Period']}", f_out ) if not pdict['Files']['FileCommon'] : pdict['Files']['FileCommon'] = f"{pdict['Experiment']['JobName']}_{pdict['Experiment']['Period']}" if not pdict['Files']['Title'] : pdict['Files']['Title'] = f"{pdict['Experiment']['JobName']} : {pdict['Experiment']['Freq']} : {pdict['Experiment']['DateBegin']} - {pdict['Experiment']['DateEnd']}" echo ('\nOpen history files', f_out ) if not pdict['Files']['file_ATM_his'] : if pdict['Experiment']['ATM_HIS'] == 'latlon' : pdict['Files']['file_ATM_his'] = os.path.join ( pdict['Files']['dir_ATM_his'], f"{pdict['Files']['FileCommon']}_histmth.nc" ) if pdict['Experiment']['ATM_HIS'] == 'ico' : pdict['Files']['file_ATM_his'] = os.path.join ( pdict['Files']['dir_ATM_his'], f"{pdict['Files']['FileCommon']}_histmth_ico.nc" ) if not pdict['Files']['file_SRF_his'] : if pdict['Experiment']['ATM_HIS'] == 'latlon' : pdict['Files']['file_SRF_his'] = os.path.join ( pdict['Files']['dir_SRF_his'], f"{pdict['Files']['FileCommon']}_sechiba_history.nc" ) if pdict['Experiment']['ATM_HIS'] == 'ico' : pdict['Files']['file_SRF_his'] = os.path.join ( pdict['Files']['dir_SRF_his'], f"{pdict['Files']['FileCommon']}_sechiba_history_ico.nc" ) if pdict['Experiment']['SECHIBA'] and pdict['Experiment']['Routing'] == 'SIMPLE' : if not pdict['Files']['file_RUN_his'] : if pdict['Experiment']['RUN_HIS'] == 'latlon' : pdict['Files']['file_RUN_his'] = os.path.join ( pdict['Files']['dir_SRF_his'], f"{pdict['Files']['FileCommon']}_sechiba_history.nc" ) if pdict['Experiment']['SRF_HIS'] == 'ico' : pdict['Files']['file_RUN_his'] = os.path.join ( pdict['Files']['dir_SRF_his'], f"{pdict['Files']['FileCommon']}_sechiba_history_ico.nc" ) echo ( f"{pdict['Files']['file_ATM_his'] = }", f_out ) if pdict['Experiment']['SECHIBA'] : echo ( f"{pdict['Files']['file_SRF_his'] = }", f_out ) if pdict['Experiment']['Routing'] == 'SIMPLE' : echo ( f"{pdict['Files']['file_RUN_his'] = }", f_out ) if not pdict['Files']['file_OCE_his'] : pdict['Files']['file_OCE_his'] = os.path.join ( pdict['Files']['dir_OCE_his'], f"{pdict['Files']['FileCommon']}_grid_T.nc" ) if not pdict['Files']['file_OCE_sca'] : pdict['Files']['file_OCE_sca'] = os.path.join ( pdict['Files']['dir_OCE_his'], f"{pdict['Files']['FileCommon']}_scalar.nc" ) if not pdict['Files']['file_OCE_srf'] : pdict['Files']['file_OCE_srf'] = os.path.join ( pdict['Files']['dir_OCE_his'], f"{pdict['Files']['FileCommon']}_sbc.nc" ) if not pdict['Files']['file_ICE_his'] : pdict['Files']['file_ICE_his'] = os.path.join ( pdict['Files']['dir_ICE_his'], f"{pdict['Files']['FileCommon']}_icemod.nc" ) return pdict def SetRestartNames ( pdict, f_out) : ''' Defines dates for restart files Define names of tar files with restart Returns a dictionary ''' if not pdict['Files']['TarRestartDate_beg'] : pdict['Files']['TarRestartDate_beg'] = ldate.SubOneDayToDate ( pdict['Experiment']['DateBegin'], pdict['Experiment']['CalendarType'] ) if not pdict['Files']['TarRestartDate_end'] : pdict['Files']['TarRestartDate_end'] = ldate.ConvertFormatToGregorian ( pdict['Experiment']['DateEnd'] ) if not pdict['Files']['TarRestartPeriod_beg'] : pdict['Files']['TarRestartPeriod_beg_DateEnd'] = pdict['Files']['TarRestartDate_beg'] pdict['Files']['TarRestartPeriod_beg_DateBeg'] = ldate.DateAddPeriod ( pdict['Files']['TarRestartPeriod_beg_DateEnd'], f"-{pdict['Experiment']['PackFrequency']}") pdict['Files']['TarRestartPeriod_beg_DateBeg'] = ldate.AddOneDayToDate ( pdict['Files']['TarRestartPeriod_beg_DateBeg'], pdict['Experiment']['CalendarType'] ) pdict['Files']['TarRestartPeriod_beg'] = f"{pdict['Files']['TarRestartPeriod_beg_DateBeg']}_{pdict['Files']['TarRestartPeriod_beg_DateEnd']}" echo ( f"Tar period for initial restart : {pdict['Files']['TarRestartPeriod_beg']}", f_out) if not pdict['Files']['TarRestartPeriod_end'] : pdict['Files']['TarRestartPeriod_end_DateEnd'] = pdict['Files']['TarRestartDate_end'] pdict['Files']['TarRestartPeriod_end_DateBeg'] = ldate.DateAddPeriod ( pdict['Files']['TarRestartPeriod_end_DateEnd'], f"-{pdict['Experiment']['PackFrequency']}" ) pdict['Files']['TarRestartPeriod_end_DateBeg'] = ldate.AddOneDayToDate ( pdict['Files']['TarRestartPeriod_end_DateBeg'], pdict['Experiment']['CalendarType'] ) pdict['Files']['TarRestartPeriod_end'] = f"{pdict['Files']['TarRestartPeriod_end_DateBeg']}_{pdict['Files']['TarRestartPeriod_end_DateEnd']}" echo ( f"Tar period for final restart : {pdict['Files']['TarRestartPeriod_end']}", f_out ) echo ( f"Restart dates - Start : {pdict['Files']['TarRestartPeriod_beg']} / End : {pdict['Files']['TarRestartPeriod_end']}", f_out) if not pdict['Files']['tar_restart_beg'] : pdict['Files']['tar_restart_beg'] = os.path.join ( pdict['libIGCM']['R_SAVE'], 'RESTART', f"{pdict['Experiment']['JobName']}_{pdict['Files']['TarRestartPeriod_beg']}_restart.tar" ) if not pdict['Files']['tar_restart_end'] : pdict['Files']['tar_restart_end'] = os.path.join ( pdict['libIGCM']['R_SAVE'], 'RESTART', f"{pdict['Experiment']['JobName']}_{pdict['Files']['TarRestartPeriod_end']}_restart.tar" ) echo ( f"{ pdict['Files']['tar_restart_beg'] = }", f_out ) echo ( f"{ pdict['Files']['tar_restart_end'] = }", f_out ) ##-- Names of tar files with restarts if not pdict['Files']['tar_restart_beg_ATM'] : pdict['Files']['tar_restart_beg_ATM'] = pdict['Files']['tar_restart_beg'] if not pdict['Files']['tar_restart_beg_DYN'] : pdict['Files']['tar_restart_beg_DYN'] = pdict['Files']['tar_restart_beg'] if not pdict['Files']['tar_restart_beg_RUN'] : pdict['Files']['tar_restart_beg_RUN'] = pdict['Files']['tar_restart_beg'] if not pdict['Files']['tar_restart_beg_OCE'] : pdict['Files']['tar_restart_beg_OCE'] = pdict['Files']['tar_restart_beg'] if not pdict['Files']['tar_restart_beg_ICE'] : pdict['Files']['tar_restart_beg_ICE'] = pdict['Files']['tar_restart_beg'] if not pdict['Files']['tar_restart_end_ATM'] : pdict['Files']['tar_restart_end_ATM'] = pdict['Files']['tar_restart_end'] if not pdict['Files']['tar_restart_end_DYN'] : pdict['Files']['tar_restart_end_DYN'] = pdict['Files']['tar_restart_end'] if not pdict['Files']['tar_restart_end_RUN'] : pdict['Files']['tar_restart_end_RUN'] = pdict['Files']['tar_restart_end'] if not pdict['Files']['tar_restart_end_OCE'] : pdict['Files']['tar_restart_end_OCE'] = pdict['Files']['tar_restart_end'] if not pdict['Files']['tar_restart_end_ICE'] : pdict['Files']['tar_restart_end_ICE'] = pdict['Files']['tar_restart_end'] if not pdict['Files']['file_DYN_beg'] : if pdict['Experiment']['LMDZ'] : pdict['Files']['file_DYN_beg'] = f"{pdict['Files']['FileDir']}/ATM_{pdict['Experiment']['JobName']}_{pdict['Files']['TarRestartDate_beg']}_restart.nc" if pdict['Experiment']['ICO'] : pdict['Files']['file_DYN_beg'] = f"{pdict['Files']['FileDir']}/ICO_{pdict['Experiment']['JobName']}_{pdict['Files']['TarRestartDate_beg']}_restart.nc" if not pdict['Files']['file_DYN_end'] : if pdict['Experiment']['LMDZ'] : pdict['Files']['file_DYN_end'] = f"{pdict['Files']['FileDir']}/ATM_{pdict['Experiment']['JobName']}_{pdict['Files']['TarRestartDate_end']}_restart.nc" if pdict['Experiment']['ICO'] : pdict['Files']['file_DYN_end'] = f"{pdict['Files']['FileDir']}/ICO_{pdict['Experiment']['JobName']}_{pdict['Files']['TarRestartDate_end']}_restart.nc" if pdict['Experiment']['SECHIBA'] : if not pdict['Files']['file_SRF_beg'] : pdict['Files']['file_SRF_beg'] = f"{pdict['Files']['FileDir']}/SRF_{pdict['Experiment']['JobName']}_{pdict['Files']['TarRestartDate_beg']}_sechiba_rest.nc" if not pdict['Files']['file_SRF_end'] : pdict['Files']['file_SRF_end'] = f"{pdict['Files']['FileDir']}/SRF_{pdict['Experiment']['JobName']}_{pdict['Files']['TarRestartDate_end']}_sechiba_rest.nc" if pdict['Experiment']['SECHIBA'] : if not pdict['Files']['tar_restart_beg_SRF'] : pdict['Files']['tar_restart_beg_SRF'] = pdict['Files']['tar_restart_beg'] if not pdict['Files']['tar_restart_end_SRF'] : pdict['Files']['tar_restart_end_SRF'] = pdict['Files']['tar_restart_end'] if not pdict['Files']['file_ATM_beg'] : pdict['Files']['file_ATM_beg'] = f"{pdict['Files']['FileDir']}/ATM_{pdict['Experiment']['JobName']}_{pdict['Files']['TarRestartDate_beg']}_restartphy.nc" if not pdict['Files']['file_ATM_end'] : pdict['Files']['file_ATM_end'] = f"{pdict['Files']['FileDir']}/ATM_{pdict['Experiment']['JobName']}_{pdict['Files']['TarRestartDate_end']}_restartphy.nc" if pdict['Experiment']['SECHIBA'] : echo ( f"{pdict['Files']['file_SRF_beg'] = }", f_out) echo ( f"{pdict['Files']['file_SRF_end'] = }", f_out) if pdict['Experiment']['ICO'] : if not pdict['Files']['file_DYN_aire'] : pdict['Files']['file_DYN_aire'] = os.path.join ( pdict['libIGCM']['R_IN'], 'ATM', 'GRID', f"{pdict['Experiment']['ResolAtm']}_grid.nc" ) if pdict['Experiment']['SECHIBA'] and pdict['Experiment']['Routing'] == 'SIMPLE' : if not pdict['Files']['file_RUN_beg'] : pdict['Files']['file_RUN_beg'] = f"{ pdict['Files']['FileDir']}/SRF_{pdict['Experiment']['JobName']}_{pdict['Files']['TarRestartDate_beg']}_routing_restart.nc" if not pdict['Files']['file_RUN_end'] : pdict['Files']['file_RUN_end'] = f"{ pdict['Files']['FileDir']}/SRF_{pdict['Experiment']['JobName']}_{pdict['Files']['TarRestartDate_end']}_routing_restart.nc" if not pdict['Files']['tar_restart_beg_OCE'] : pdict['Files']['tar_restart_beg_OCE'] = pdict['Files']['tar_restart_beg'] if not pdict['Files']['tar_restart_beg_ICE'] : pdict['Files']['tar_restart_beg_ICE'] = pdict['Files']['tar_restart_beg'] if not pdict['Files']['tar_restart_end_OCE'] : pdict['Files']['tar_restart_end_OCE'] = pdict['Files']['tar_restart_end'] if not pdict['Files']['tar_restart_end_ICE'] : pdict['Files']['tar_restart_end_ICE'] = pdict['Files']['tar_restart_end'] if not pdict['Files']['file_OCE_beg'] : pdict['Files']['file_OCE_beg'] = f"{pdict['Files']['FileDir']}/OCE_{pdict['Experiment']['JobName']}_{pdict['Files']['TarRestartDate_beg']}_restart.nc" if not pdict['Files']['file_OCE_end'] : pdict['Files']['file_OCE_end'] = f"{pdict['Files']['FileDir']}/OCE_{pdict['Experiment']['JobName']}_{pdict['Files']['TarRestartDate_end']}_restart.nc" if not pdict['Files']['file_ICE_beg'] : pdict['Files']['file_ICE_beg'] = f"{pdict['Files']['FileDir']}/ICE_{pdict['Experiment']['JobName']}_{pdict['Files']['TarRestartDate_beg']}_restart_icemod.nc" if not pdict['Files']['file_ICE_end'] : pdict['Files']['file_ICE_end'] = f"{pdict['Files']['FileDir']}/ICE_{pdict['Experiment']['JobName']}_{pdict['Files']['TarRestartDate_end']}_restart_icemod.nc" return pdict def ComputeGridATM ( dpar, d_ATM_his ) : readPrec = dpar['Config']['readPrec'] if repr(readPrec) in [ "" , float ] : def rprec (ptab) : return ptab else : def rprec (ptab) : return ptab.astype(readPrec).astype(float) ATM = types.SimpleNamespace () # ATM grid with cell surfaces if dpar['Experiment']['LMDZ'] : #echo ('ATM grid with cell surfaces : LMDZ', f_out) ATM.lat = lmdz.geo2point ( rprec (d_ATM_his ['lat'])+0*rprec (d_ATM_his ['lon']), dim1d='cell' ) ATM.lon = lmdz.geo2point ( 0*rprec (d_ATM_his ['lat'])+ rprec (d_ATM_his ['lon']), dim1d='cell' ) ATM.aire = lmdz.geo2point ( rprec (d_ATM_his ['aire'] [0]), cumul_poles=True, dim1d='cell' ) ATM.fter = lmdz.geo2point ( rprec (d_ATM_his ['fract_ter'][0]), dim1d='cell' ) ATM.foce = lmdz.geo2point ( rprec (d_ATM_his ['fract_oce'][0]), dim1d='cell' ) ATM.fsic = lmdz.geo2point ( rprec (d_ATM_his ['fract_sic'][0]), dim1d='cell' ) ATM.flic = lmdz.geo2point ( rprec (d_ATM_his ['fract_lic'][0]), dim1d='cell' ) if dpar['Experiment']['ICO'] : if dpar['Experiment']['ATM_HIS'] == 'latlon' : #echo ( 'ATM areas and fractions on LATLON grid' ) if 'lat_dom_out' in d_ATM_his.variables : ATM.lat = lmdz.geo2point ( rprec (d_ATM_his ['lat_dom_out'])+0*rprec (d_ATM_his ['lon_dom_out']), dim1d='cell' ) ATM.lon = lmdz.geo2point ( 0*rprec (d_ATM_his ['lat_dom_out'])+ rprec (d_ATM_his ['lon_dom_out']), dim1d='cell' ) else : ATM.lat = lmdz.geo2point ( rprec (d_ATM_his ['lat'])+0*rprec (d_ATM_his ['lon']), dim1d='cell' ) ATM.lon = lmdz.geo2point ( 0*rprec (d_ATM_his ['lat'])+ rprec (d_ATM_his ['lon']), dim1d='cell' ) ATM.aire = lmdz.geo2point ( rprec (d_ATM_his ['aire'][0]).squeeze(), cumul_poles=True, dim1d='cell' ) ATM.fter = lmdz.geo2point ( rprec (d_ATM_his ['fract_ter'][0]), dim1d='cell' ) ATM.foce = lmdz.geo2point ( rprec (d_ATM_his ['fract_oce'][0]), dim1d='cell' ) ATM.fsic = lmdz.geo2point ( rprec (d_ATM_his ['fract_sic'][0]), dim1d='cell' ) ATM.flic = lmdz.geo2point ( rprec (d_ATM_his ['fract_lic'][0]), dim1d='cell' ) if dpar['Experiment']['ATM_HIS'] == 'ico' : #echo ( 'ATM areas and fractions on ICO grid' ) ATM.aire = rprec (d_ATM_his ['aire'] [0]).squeeze() ATM.lat = rprec (d_ATM_his ['lat'] ) ATM.lon = rprec (d_ATM_his ['lon'] ) ATM.fter = rprec (d_ATM_his ['fract_ter'][0]) ATM.foce = rprec (d_ATM_his ['fract_oce'][0]) ATM.fsic = rprec (d_ATM_his ['fract_sic'][0]) ATM.flic = rprec (d_ATM_his ['fract_lic'][0]) ATM.fsea = ATM.foce + ATM.fsic ATM.flnd = ATM.fter + ATM.flic ATM.aire_fter = ATM.aire * ATM.fter ATM.aire_flic = ATM.aire * ATM.flic ATM.aire_fsic = ATM.aire * ATM.fsic ATM.aire_foce = ATM.aire * ATM.foce ATM.aire_flnd = ATM.aire * ATM.flnd ATM.aire_fsea = ATM.aire * ATM.fsea ATM.aire_sea = ATM.aire * ATM.fsea ATM.aire_tot = P1sum ( ATM.aire ) ATM.aire_sea_tot = P1sum ( ATM.aire_fsea ) ATM.aire_ter_tot = P1sum ( ATM.aire_fter ) ATM.aire_lic_tot = P1sum ( ATM.aire_flic ) return dpar, ATM def ComputeGridSRF ( dpar, d_SRF_his ) : readPrec = dpar['Config']['readPrec'] if repr(readPrec) in [ "" , float ] : def rprec (ptab) : return ptab else : def rprec (ptab) : return ptab.astype(readPrec).astype(float) SRF = types.SimpleNamespace () if dpar['Experiment']['SECHIBA'] : if dpar['Experiment']['LMDZ'] : SRF.lat = lmdz.geo2point ( rprec (d_SRF_his ['lat'])+0*rprec (d_SRF_his ['lon']), dim1d='cell' ) SRF.lon = lmdz.geo2point ( 0*rprec (d_SRF_his ['lat'])+ rprec (d_SRF_his ['lon']), dim1d='cell' ) SRF.aire = lmdz.geo2point ( rprec (d_SRF_his ['Areas']) * rprec (d_SRF_his ['Contfrac']), dim1d='cell', cumul_poles=True ) SRF.areas = lmdz.geo2point ( rprec (d_SRF_his ['Areas']) , dim1d='cell', cumul_poles=True ) SRF.contfrac = lmdz.geo2point ( rprec (d_SRF_his ['Contfrac']), dim1d='cell' ) if dpar['Experiment']['ICO'] : if dpar['Experiment']['SRF_HIS'] == 'latlon' : #echo ( 'SRF areas and fractions on LATLON grid' ) if 'lat_domain_landpoints_out' in d_SRF_his : SRF.lat = lmdz.geo2point ( rprec (d_SRF_his ['lat_domain_landpoints_out'])+0*rprec (d_SRF_his ['lon_domain_landpoints_out']), dim1d='cell' ) SRF.lon = lmdz.geo2point ( 0*rprec (d_SRF_his ['lat_domain_landpoints_out'])+ rprec (d_SRF_his ['lon_domain_landpoints_out']), dim1d='cell' ) else : if 'lat_domain_landpoints_out' in d_SRF_his : SRF.lat = lmdz.geo2point ( rprec (d_SRF_his ['lat_dom_out'])+0*rprec (d_SRF_his ['lon_dom_out']), dim1d='cell' ) SRF.lon = lmdz.geo2point ( 0*rprec (d_SRF_his ['lat_dom_out'])+ rprec (d_SRF_his ['lon_dom_out']), dim1d='cell' ) else : SRF.lat = lmdz.geo2point ( rprec (d_SRF_his ['lat'])+0*rprec (d_SRF_his ['lon']), dim1d='cell' ) SRF.lon = lmdz.geo2point ( 0*rprec (d_SRF_his ['lat'])+ rprec (d_SRF_his ['lon']), dim1d='cell' ) SRF.areas = lmdz.geo2point ( rprec (d_SRF_his ['Areas'] ) , dim1d='cell', cumul_poles=True ) SRF.areafrac = lmdz.geo2point ( rprec (d_SRF_his ['AreaFrac']) , dim1d='cell', cumul_poles=True ) SRF.contfrac = lmdz.geo2point ( rprec (d_SRF_his ['Contfrac']) , dim1d='cell', cumul_poles=True ) SRF.aire = SRF.areafrac if dpar['Experiment']['SRF_HIS'] == 'ico' : #echo ( 'SRF areas and fractions on ICO grid' ) SRF.lat = rprec (d_SRF_his ['lat'] ) SRF.lon = rprec (d_SRF_his ['lon'] ) SRF.areas = rprec (d_SRF_his ['Areas'] ) SRF.contfrac = rprec (d_SRF_his ['Contfrac']) SRF.aire = SRF.areas * SRF.contfrac SRF.aire_tot = P1sum ( SRF.aire ) return dpar, SRF def ComputeGridDYN ( dpar, ATM, d_DYN_aire, d_ATM_beg ) : readPrec = dpar['Config']['readPrec'] if repr(readPrec) in [ "" , float ] : def rprec (ptab) : return ptab else : def rprec (ptab) : return ptab.astype(readPrec).astype(float) DYN = types.SimpleNamespace () if dpar['Experiment']['ICO'] : if dpar['Config']['SortIco'] : # Creation d'une clef de tri pour le fichier aire DYN.aire_keysort = np.lexsort ( (d_DYN_aire['lat'], d_DYN_aire['lon']) ) else : DYN.aire_keysort = np.arange ( len ( d_DYN_aire['lat'] ) ) DYN.lat = d_DYN_aire['lat'] DYN.lon = d_DYN_aire['lon'] DYN.aire = d_DYN_aire['aire'] DYN.fsea = d_DYN_aire['fract_oce'] + d_DYN_aire['fract_sic'] DYN.flnd = 1.0 - DYN.fsea DYN.fter = d_ATM_beg['FTER'] DYN.flic = d_ATM_beg['FLIC'] DYN.foce = d_ATM_beg['FOCE'] DYN.aire_fter = DYN.aire * DYN.fter if dpar['Experiment']['ATM_HIS'] == 'ico' : DYN.aire = ATM.aire DYN.foce = ATM.foce DYN.fsic = ATM.fsic DYN.flic = ATM.flic DYN.fter = ATM.fter DYN.fsea = ATM.fsea DYN.flnd = ATM.flnd DYN.aire_fter = ATM.aire_fter DYN.aire_flic = ATM.aire_flic DYN.aire_fsic = ATM.aire_fsic DYN.aire_foce = ATM.aire_foce DYN.aire_flnd = ATM.aire_flnd DYN.aire_fsea = ATM.aire_fsea if dpar['Experiment']['LMDZ'] : # Area on lon/lat grid DYN.aire = ATM.aire DYN.fsea = ATM.fsea DYN.flnd = ATM.flnd DYN.fter = rprec (d_ATM_beg['FTER']) DYN.flic = rprec (d_ATM_beg['FLIC']) DYN.aire_fter = DYN.aire * DYN.fter DYN.aire_tot = P1sum ( DYN.aire ) return dpar, DYN def ComputeGridOCE ( dpar, d_OCE_his, nperio=None ) : OCE = types.SimpleNamespace () # Get mask and surfaces sos = d_OCE_his ['sos'][0].squeeze() OCE.OCE_msk = nemo.lbc_mask ( xr.where ( sos>0., 1., 0. ), cd_type = 'T', nperio=nperio, sval = 0. ) so = sos = d_OCE_his ['sos'][0].squeeze() OCE.msk3 = nemo.lbc_mask ( xr.where ( so>0., 1., 0. ), cd_type = 'T', nperio=nperio, sval = 0. ) # lbc_mask removes the duplicate points (periodicity and north fold) OCE.OCE_aire = nemo.lbc_mask ( d_OCE_his ['area'] * OCE.OCE_msk, cd_type = 'T', nperio=nperio, sval = 0.0 ) OCE.ICE_aire = OCE.OCE_aire OCE.OCE_aire_tot = P1sum ( OCE.OCE_aire ) OCE.ICE_aire_tot = P1sum ( OCE.ICE_aire ) return dpar, OCE def config2dict ( pconf ) : ''' Convert a config parser object into a dictionary ''' pdict = {} for section in pconf.keys () : pdict[section] = {} for key in pconf[section].keys() : zz = str2value ( pconf[section][key] ) pdict[section].update ( {key:zz} ) return pdict def dict2config ( pdict ): ''' Convert a dictionary into a configparser object All values are converted to strings The dictionary must have two levels (see configparser) ''' zconf = configparser.ConfigParser ( interpolation=configparser.ExtendedInterpolation() ) zconf.optionxform = str # To keep capitals zconf.read_dict(dict2str(pdict)) return zconf def updateConf ( pconf, pdict ): ''' Update a config parser with a dictionary All values are converted to strings The dictionary must have two levels (see configparser) ''' pconf.read_dict (dict1str(pdict)) #for section in pdict.keys() : # if section not in pconf.keys() : pconf[section] = {} # for key in pdict[section].keys() : # pconf[section][key] = str(pdict[section][key]) return pconf def dict2str (pdict) : ''' Convert all values of a dictionary to strings ''' zout = {} for section in pdict.keys() : zout[section] = {} for key in pdict[section].keys() : zout[section][key] = str(pdict[section][key]) return zout def echo (string, f_out, end='\n') : '''Function to print to stdout *and* output file''' print ( str(string), end=end ) sys.stdout.flush () f_out.write ( str(string) + end ) f_out.flush () def str2value ( pz ) : ''' Tries to convert a string into integer, real, boolean or None ''' zz = pz zz = setBool (zz) zz = setNum (zz) zz = setNone (zz) return zz def setBool (chars) : '''Convert specific char string in boolean if possible''' z_set_bool = chars if isinstance (chars, str) : for key in configparser.ConfigParser.BOOLEAN_STATES.keys () : if chars.lower() == key : z_set_bool = configparser.ConfigParser.BOOLEAN_STATES[key] return z_set_bool def setNum (chars) : '''Convert specific char string in integer or real if possible''' zset_num = chars if isinstance (chars, str) : realnum = re.compile ("^[-+]?[0-9]*\.?[0-9]+(e[-+]?[0-9]+)?$") is_number = realnum.match(chars.strip()) != None is_int = chars.strip().isdigit() if is_number : if is_int : zset_num = int (chars) else : zset_num = float (chars) return zset_num def setNone (chars) : '''Convert specific char string to None if possible''' zset_none = chars if isinstance (chars, str) : if chars in ['None', 'NONE', 'none'] : zset_none = None return zset_none def unDefined (char) : '''True if a variable is not defined, or set to None''' if char in globals () : if globals()[char] : zun_defined = False else : zun_defined = True else : zun_defined = True return zun_defined def Defined (char) : '''True if a variable is defined and not equal to None''' if char in globals () : if globals()[char] is None : zdefined = False else : zdefined = True else : zdefined = False return zdefined def Ksum (tab) : '''Kahan summation algorithm, also known as compensated summation https://en.wikipedia.org/wiki/Kahan_summation_algorithm ''' # Prepare the accumulator. ksum = 0.0 # A running compensation for lost low-order bits. comp = 0.0 for xx in np.where ( np.isnan(tab), 0., tab ) : # comp is zero the first time around. yy = xx - comp # Alas, sum is big, y small, so low-order digits of y are lost. tt = ksum + yy # (tt - Ksum) cancels the high-order part of y; subtracting y recovers negative (low part of yy) comp = (tt - ksum) - yy # Algebraically, comp should always be zero. Beware overly-aggressive optimizing compilers ! ksum = tt # Next time around, the lost low part will be added to y in a fresh attempt. return ksum def Ssum (tab) : '''Precision summation by sorting by absolute values''' ssum = np.sum ( tab[np.argsort(np.abs(tab))] ) return ssum def KSsum (tab) : '''Precision summation by sorting by absolute value, and applying Kahan summation''' kssum = Ksum ( tab[np.argsort(np.abs(tab))] ) return kssum # Choosing algorithm for precise summation Psum = KSsum def P1sum (ptab) : return Psum ( ptab.to_masked_array().ravel() ) class Timer : '''Measures execution time of a function''' def __str__ (self): return str (self.__class__) def __name__ (self): return self.__class__.__name__ def __init__ (self, func, debug=False, timer=True) : functools.update_wrapper (self, func) self.func = func self.num_calls = 0 self.cumul_time = 0. self.debug = debug self.timer = timer def __call__ (self, *args, **kwargs) : self.num_calls += 1 if self.debug : print ( f'>-- Calling {self.__name__} --------------------------------------------------' ) args_repr = [f"{repr (a)}" for a in args] kwargs_repr = [f"{k}={v!r}" for k, v in kwargs.items ()] signature = ", ".join (args_repr + kwargs_repr) print ( f"Signature : ({signature}") start_time = time.perf_counter () values = self.func (*args, **kwargs) end_time = time.perf_counter () run_time = end_time - start_time self.cumul_time += run_time if self.timer : print ( f"--> Calling {self.__name__!r} : {run_time:.3f} secs / cumul {self.cumul_time:.3f} # {self.num_calls:2d} calls") if self.debug : print ( '<------------------------------------------------------------' ) return values