#!/usr/bin/env python3 ### ### Script to check water conservation in the IPSL coupled model ### ## 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 SVN = { 'Author' : "$Author$", 'Date' : "$Date$", 'Revision': "$Revision$", 'Id' : "$Id$", 'HeadURL' : "$HeadUrl: $" } ### ### ## Import system modules import sys import os import configparser ## Import needed scientific modules import numpy as np import xarray as xr ## Import local modules import WaterUtils as wu import nemo, lmdz ## Read command line arguments ## --------------------------- print ( "Name of Python script:", sys.argv[0] ) IniFile = sys.argv[1] # Test existence of IniFile if not os.path.exists (IniFile ) : raise FileExistsError ( f"File not found : {IniFile = }" ) if 'full' in IniFile or 'ATM' in IniFile : FullIniFile = IniFile else : FullIniFile = 'ATM_' + IniFile print ("Output file : ", FullIniFile ) ## Experiment parameters ## -------------------- dpar = wu.ReadConfig ( IniFile ) ## Configure all needed parameter from existant parameters ## ------------------------------------------------------- dpar = wu.SetDatesAndFiles ( dpar ) ## Output file with water budget diagnostics ## ----------------------------------------- f_out = dpar['Files']['f_out'] ## Put dpar values in local namespace ## ---------------------------------- for Section in dpar.keys () : print ( f'\nReading [{Section}]' ) for VarName in dpar[Section].keys() : locals()[VarName] = dpar[Section][VarName] print ( f' {VarName:21} set to : {locals()[VarName]}' ) ## Debuging and timer Timer = wu.functools.partial (wu.Timer, debug=Debug, timer=Timing) ## Useful functions ## ---------------- if repr(readPrec) == "" or readPrec == float : def rprec (ptab) : '''This version does nothing rprec may be use to reduce floating precision when reading history files ''' return ptab else : def rprec (ptab) : '''Returns float with a different precision''' return ptab.astype(readPrec).astype(float) def kg2Sv (val, rho=ATM_RHO) : '''From kg to Sverdrup''' return val/dtime_sec*1.0e-6/rho def kg2myear (val, rho=ATM_RHO) : '''From kg to m/year''' return val/ATM.aire_sea_tot/rho/NbYear def var2prt (var, small=False, rho=ATM_RHO) : '''Formats value for printing''' if small : return var , kg2Sv(var, rho=rho)*1000., kg2myear(var, rho=rho)*1000 else : return var , kg2Sv(var, rho=rho) , kg2myear(var, rho=rho) def prtFlux (Desc, var, Form='F', small=False, rho=ATM_RHO, width=15) : '''Pretty print of formattd value''' if small : if Form in ['f', 'F'] : ff=" {:14.6e} kg | {:12.4f} mSv | {:12.4f} mm/year " if Form in ['e', 'E'] : ff=" {:14.6e} kg | {:12.4e} mSv | {:12.4e} mm/year " else : if Form in ['f', 'F'] : ff=" {:14.6e} kg | {:12.4f} Sv | {:12.4f} m/year " if Form in ['e', 'E'] : ff=" {:14.6e} kg | {:12.4e} Sv | {:12.4e} m/year " echo ( (' {:>{width}} = ' +ff).format (Desc, *var2prt(var, small=small, rho=rho), width=width ) ) return None def echo (string, 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 () return None d_ATM_his = xr.open_dataset ( file_ATM_his, use_cftime=True, decode_times=True, decode_cf=True ).squeeze() if SECHIBA : d_SRF_his = xr.open_dataset ( file_SRF_his, use_cftime=True, decode_times=True, decode_cf=True ).squeeze() if Routing == 'SECHIBA' : d_RUN_his = d_SRF_his if Routing == 'SIMPLE' : d_RUN_his = xr.open_dataset ( file_RUN_his, use_cftime=True, decode_times=True, decode_cf=True ).squeeze() d_OCE_his = xr.open_dataset ( file_OCE_his, use_cftime=True, decode_times=True, decode_cf=True ).squeeze() d_OCE_sca = xr.open_dataset ( file_OCE_sca, use_cftime=True, decode_times=True, decode_cf=True ).squeeze() #d_OCE_srf = xr.open_dataset ( file_OCE_srf, use_cftime=True, decode_times=True, decode_cf=True ).squeeze() d_ICE_his = xr.open_dataset ( file_ICE_his, use_cftime=True, decode_times=True, decode_cf=True ).squeeze() if NEMO == '3.6' : d_ICE_his = d_ICE_his.rename ( {'y_grid_T':'y', 'x_grid_T':'x'} ) echo ( f'{file_OCE_his = }' ) echo ( f'{file_ICE_his = }' ) echo ( f'{file_OCE_sca = }' ) echo ( f'{file_OCE_srf = }' ) ## Compute run length ## ------------------ dtime = ( d_ATM_his.time_counter_bounds.max() - d_ATM_his.time_counter_bounds.min() ) echo ('\nRun length : {:8.2f} days'.format ( (dtime/np.timedelta64(1, "D")).values ) ) dtime_sec = (dtime/np.timedelta64(1, "s")).values.item() # Convert in seconds ## Compute length of each period dtime_per = (d_ATM_his.time_counter_bounds[:,-1] - d_ATM_his.time_counter_bounds[:,0] ) echo ('\nPeriods lengths (days) : {:} days'.format ( (dtime_per/np.timedelta64(1, "D")).values ) ) dtime_per_sec = (dtime_per/np.timedelta64(1, "s")).values # In seconds dtime_per_sec = xr.DataArray (dtime_per_sec, dims=["time_counter", ], coords=[d_ATM_his.time_counter,] ) dtime_per_sec.attrs['unit'] = 's' # Number of years NbYear = dtime_sec / YEAR_LENGTH ## Write the full configuration ## ---------------------------- params_out = open (FullIniFile, 'w', encoding = 'utf-8') params = wu.dict2config ( dpar ) params.write ( params_out ) params_out.close () ## Compute aire, fractions, etc ... ## -------------------------------- if ICO : if not file_DYN_aire : file_DYN_aire = os.path.join ( R_IN, 'ATM', 'GRID', ResolAtm+'_grid.nc' ) dpar['Files']['file_DYN_aire'] = file_DYN_aire echo ( f'{file_DYN_aire = }' ) d_DYN_aire = xr.open_dataset ( file_DYN_aire, decode_times=False ).squeeze() else : d_DYN_aire = None dpar, ATM = wu.ComputeGridATM ( dpar, d_ATM_his) dpar, SRF = wu.ComputeGridSRF ( dpar, d_SRF_his) dpar, OCE = wu.ComputeGridOCE ( dpar, d_OCE_his, nperio=nperio ) def ATM_flux_int (flux) : '''Integrate (* time * surface) flux on atmosphere grid''' zATM_flux_int = wu.P1sum ( flux * dtime_per_sec * ATM.aire ) return zATM_flux_int def SRF_flux_int (flux) : '''Integrate (* time * surface) flux on land grid''' zSRF_flux_int = wu.P1sum ( flux * dtime_per_sec * SRF.aire ) return zSRF_flux_int def ONE_stock_int (stock) : '''Sum stock''' zONE_stock_int = wu.P1sum ( stock ) return zONE_stock_int def OCE_flux_int (flux) : '''Integrate flux on oce grid''' zOCE_flux_int = wu.P1sum ( flux * OCE.OCE_aire * dtime_per_sec ) return zOCE_flux_int def ONE_flux_int (flux) : '''Integrate flux on oce grid''' zOCE_flux_int = wu.P1sum ( flux * dtime_per_sec ) return zOCE_flux_int echo ( '\n====================================================================================' ) echo ( f'-- ATM Fluxes -- {Title} ' ) if ATM_HIS == 'latlon' : echo ( ' latlon case' ) ATM.wbilo_oce = lmdz.geo2point ( rprec (d_ATM_his ['wbilo_oce']), dim1D='cell' ) ATM.wbilo_sic = lmdz.geo2point ( rprec (d_ATM_his ['wbilo_sic']), dim1D='cell' ) ATM.wbilo_ter = lmdz.geo2point ( rprec (d_ATM_his ['wbilo_ter']), dim1D='cell' ) ATM.wbilo_lic = lmdz.geo2point ( rprec (d_ATM_his ['wbilo_lic']), dim1D='cell' ) ATM.runofflic = lmdz.geo2point ( rprec (d_ATM_his ['runofflic']), dim1D='cell' ) ATM.fqcalving = lmdz.geo2point ( rprec (d_ATM_his ['fqcalving']), dim1D='cell' ) ATM.fqfonte = lmdz.geo2point ( rprec (d_ATM_his ['fqfonte'] ), dim1D='cell' ) ATM.precip = lmdz.geo2point ( rprec (d_ATM_his ['precip'] ), dim1D='cell' ) ATM.snowf = lmdz.geo2point ( rprec (d_ATM_his ['snow'] ), dim1D='cell' ) ATM.evap = lmdz.geo2point ( rprec (d_ATM_his ['evap'] ), dim1D='cell' ) ATM.wevap_ter = lmdz.geo2point ( rprec (d_ATM_his ['wevap_ter']), dim1D='cell' ) ATM.wevap_oce = lmdz.geo2point ( rprec (d_ATM_his ['wevap_oce']), dim1D='cell' ) ATM.wevap_lic = lmdz.geo2point ( rprec (d_ATM_his ['wevap_lic']), dim1D='cell' ) ATM.wevap_sic = lmdz.geo2point ( rprec (d_ATM_his ['wevap_sic']), dim1D='cell' ) ATM.wrain_ter = lmdz.geo2point ( rprec (d_ATM_his ['wrain_ter']), dim1D='cell' ) ATM.wrain_oce = lmdz.geo2point ( rprec (d_ATM_his ['wrain_oce']), dim1D='cell' ) ATM.wrain_lic = lmdz.geo2point ( rprec (d_ATM_his ['wrain_lic']), dim1D='cell' ) ATM.wrain_sic = lmdz.geo2point ( rprec (d_ATM_his ['wrain_sic']), dim1D='cell' ) ATM.wsnow_ter = lmdz.geo2point ( rprec (d_ATM_his ['wsnow_ter']), dim1D='cell' ) ATM.wsnow_oce = lmdz.geo2point ( rprec (d_ATM_his ['wsnow_oce']), dim1D='cell' ) ATM.wsnow_lic = lmdz.geo2point ( rprec (d_ATM_his ['wsnow_lic']), dim1D='cell' ) ATM.wsnow_sic = lmdz.geo2point ( rprec (d_ATM_his ['wsnow_sic']), dim1D='cell' ) ATM.runofflic = lmdz.geo2point ( rprec (d_ATM_his ['runofflic']), dim1D='cell' ) echo ( f'End of LATLON case') if ATM_HIS == 'ico' : echo (' ico case') ATM.wbilo_oce = rprec (d_ATM_his ['wbilo_oce']) ATM.wbilo_sic = rprec (d_ATM_his ['wbilo_sic']) ATM.wbilo_ter = rprec (d_ATM_his ['wbilo_ter']) ATM.wbilo_lic = rprec (d_ATM_his ['wbilo_lic']) ATM.runofflic = rprec (d_ATM_his ['runofflic']) ATM.fqcalving = rprec (d_ATM_his ['fqcalving']) ATM.fqfonte = rprec (d_ATM_his ['fqfonte'] ) ATM.precip = rprec (d_ATM_his ['precip'] ) ATM.snowf = rprec (d_ATM_his ['snow'] ) ATM.evap = rprec (d_ATM_his ['evap'] ) ATM.wevap_ter = rprec (d_ATM_his ['wevap_ter']) ATM.wevap_oce = rprec (d_ATM_his ['wevap_oce']) ATM.wevap_lic = rprec (d_ATM_his ['wevap_lic']) ATM.wevap_sic = rprec (d_ATM_his ['wevap_sic']) ATM.runofflic = rprec (d_ATM_his ['runofflic']) ATM.wevap_ter = rprec (d_ATM_his ['wevap_ter']) ATM.wevap_oce = rprec (d_ATM_his ['wevap_oce']) ATM.wevap_lic = rprec (d_ATM_his ['wevap_lic']) ATM.wevap_sic = rprec (d_ATM_his ['wevap_sic']) ATM.wrain_ter = rprec (d_ATM_his ['wrain_ter']) ATM.wrain_oce = rprec (d_ATM_his ['wrain_oce']) ATM.wrain_lic = rprec (d_ATM_his ['wrain_lic']) ATM.wrain_sic = rprec (d_ATM_his ['wrain_sic']) ATM.wsnow_ter = rprec (d_ATM_his ['wsnow_ter']) ATM.wsnow_oce = rprec (d_ATM_his ['wsnow_oce']) ATM.wsnow_lic = rprec (d_ATM_his ['wsnow_lic']) ATM.wsnow_sic = rprec (d_ATM_his ['wsnow_sic']) echo ( f'End of ico case ') echo ( 'ATM wprecip_oce' ) ATM.wprecip_oce = ATM.wrain_oce + ATM.wsnow_oce ATM.wprecip_ter = ATM.wrain_ter + ATM.wsnow_ter ATM.wprecip_sic = ATM.wrain_sic + ATM.wsnow_sic ATM.wprecip_lic = ATM.wrain_lic + ATM.wsnow_lic ATM.wbilo = ATM.wbilo_oce + ATM.wbilo_sic + ATM.wbilo_ter + ATM.wbilo_lic ATM.wevap = ATM.wevap_oce + ATM.wevap_sic + ATM.wevap_ter + ATM.wevap_lic ATM.wprecip = ATM.wprecip_oce + ATM.wprecip_sic + ATM.wprecip_ter + ATM.wprecip_lic ATM.wsnow = ATM.wsnow_oce + ATM.wsnow_sic + ATM.wsnow_ter + ATM.wsnow_lic ATM.wrain = ATM.wrain_oce + ATM.wrain_sic + ATM.wrain_ter + ATM.wrain_lic ATM.wemp = ATM.wevap - ATM.wprecip ATM.emp = ATM.evap - ATM.precip ATM.wprecip_sea = ATM.wprecip_oce + ATM.wprecip_sic ATM.wsnow_sea = ATM.wsnow_oce + ATM.wsnow_sic ATM.wrain_sea = ATM.wrain_oce + ATM.wrain_sic ATM.wbilo_sea = ATM.wbilo_oce + ATM.wbilo_sic ATM.wevap_sea = ATM.wevap_sic + ATM.wevap_oce ATM.wemp_ter = ATM.wevap_ter - ATM.wprecip_ter ATM.wemp_oce = ATM.wevap_oce - ATM.wprecip_oce ATM.wemp_sic = ATM.wevap_sic - ATM.wprecip_sic ATM.wemp_lic = ATM.wevap_lic - ATM.wprecip_lic ATM.wemp_sea = ATM.wevap_sic - ATM.wprecip_oce if SECHIBA : if RUN_HIS == 'latlon' : echo ( f'RUN costalflow Grille LATLON' ) if TestInterp : echo ( f'RUN runoff TestInterp' ) SRF.RUN_runoff = lmdz.geo2point ( rprec (d_RUN_his ['runoff_contfrac_interp'] ) , dim1D='cell' ) SRF.RUN_drainage = lmdz.geo2point ( rprec (d_RUN_his ['drainage_contfrac_interp']) , dim1D='cell' ) else : echo ( f'RUN runoff' ) SRF.RUN_runoff = lmdz.geo2point ( rprec (d_RUN_his ['runoff'] ), dim1D='cell' ) SRF.RUN_drainage = lmdz.geo2point ( rprec (d_RUN_his ['drainage'] ), dim1D='cell' ) SRF.RUN_coastalflow = lmdz.geo2point ( rprec (d_RUN_his ['coastalflow'] ), dim1D='cell' ) SRF.RUN_riverflow = lmdz.geo2point ( rprec (d_RUN_his ['riverflow'] ), dim1D='cell' ) SRF.RUN_riversret = lmdz.geo2point ( rprec (d_RUN_his ['riversret'] ), dim1D='cell' ) SRF.RUN_coastalflow_cpl = lmdz.geo2point ( rprec (d_RUN_his ['coastalflow_cpl']), dim1D='cell' ) SRF.RUN_riverflow_cpl = lmdz.geo2point ( rprec (d_RUN_his ['riverflow_cpl'] ), dim1D='cell' ) if RUN_HIS == 'ico' : echo ( f'RUN costalflow Grille ICO' ) SRF.RUN_coastalflow = rprec (d_RUN_his ['coastalflow']) SRF.RUN_riverflow = rprec (d_RUN_his ['riverflow'] ) SRF.RUN_runoff = rprec (d_RUN_his ['runoff'] ) SRF.RUN_drainage = rprec (d_RUN_his ['drainage'] ) SRF.RUN_riversret = rprec (d_RUN_his ['riversret'] ) SRF.RUN_coastalflow_cpl = rprec (d_RUN_his ['coastalflow_cpl']) SRF.RUN_riverflow_cpl = rprec (d_RUN_his ['riverflow_cpl'] ) ## Correcting units of SECHIBA variables def mmd2SI ( Var ) : '''Change unit from mm/d or m^3/s to kg/s if needed''' if 'units' in Var.attrs : if Var.attrs['units'] in ['m^3/s', 'm3/s', 'm3.s-1',] : Var.values = Var.values * ATM_RHO ; Var.attrs['units'] = 'kg/s' if Var.attrs['units'] == 'mm/d' : Var.values = Var.values * ATM_RHO * (1e-3/86400.) ; Var.attrs['units'] = 'kg/s' if Var.attrs['units'] in ['m^3', 'm3', ] : Var.values = Var.values * ATM_RHO ; Var.attrs['units'] = 'kg' return Var SRF.RUN_coastalflow = mmd2SI ( SRF.RUN_coastalflow ) SRF.RUN_coastalflow_cpl = mmd2SI ( SRF.RUN_coastalflow_cpl ) SRF.RUN_drainage = mmd2SI ( SRF.RUN_drainage ) SRF.RUN_riverflow = mmd2SI ( SRF.RUN_riverflow ) SRF.RUN_riverflow_cpl = mmd2SI ( SRF.RUN_riverflow_cpl ) SRF.RUN_riversret = mmd2SI ( SRF.RUN_riversret ) SRF.RUN_runoff = mmd2SI ( SRF.RUN_runoff ) SRF.RUN_input = SRF.RUN_runoff + SRF.RUN_drainage SRF.RUN_output = SRF.RUN_coastalflow + SRF.RUN_riverflow echo ( f'ATM flw_wbilo' ) ATM.flx_wbilo = ATM_flux_int ( ATM.wbilo ) ATM.flx_wevap = ATM_flux_int ( ATM.wevap ) ATM.flx_wprecip = ATM_flux_int ( ATM.wprecip ) ATM.flx_wsnow = ATM_flux_int ( ATM.wsnow ) ATM.flx_wrain = ATM_flux_int ( ATM.wrain ) ATM.flx_wemp = ATM_flux_int ( ATM.wemp ) ATM.flx_wbilo_lic = ATM_flux_int ( ATM.wbilo_lic ) ATM.flx_wbilo_oce = ATM_flux_int ( ATM.wbilo_oce ) ATM.flx_wbilo_sea = ATM_flux_int ( ATM.wbilo_sea ) ATM.flx_wbilo_sic = ATM_flux_int ( ATM.wbilo_sic ) ATM.flx_wbilo_ter = ATM_flux_int ( ATM.wbilo_ter ) ATM.flx_calving = ATM_flux_int ( ATM.fqcalving ) ATM.flx_fqfonte = ATM_flux_int ( ATM.fqfonte ) ATM.LIC_flx_calving = ATM_flux_int ( ATM.fqcalving ) ATM.LIC_flx_fqfonte = ATM_flux_int ( ATM.fqfonte ) ATM.flx_precip = ATM_flux_int ( ATM.precip ) ATM.flx_snowf = ATM_flux_int ( ATM.snowf ) ATM.flx_evap = ATM_flux_int ( ATM.evap ) ATM.flx_runlic = ATM_flux_int ( ATM.runofflic ) ATM.flx_wrain_ter = ATM_flux_int ( ATM.wrain_ter ) ATM.flx_wrain_oce = ATM_flux_int ( ATM.wrain_oce ) ATM.flx_wrain_lic = ATM_flux_int ( ATM.wrain_lic ) ATM.flx_wrain_sic = ATM_flux_int ( ATM.wrain_sic ) ATM.flx_wrain_sea = ATM_flux_int ( ATM.wrain_sea ) ATM.flx_wsnow_ter = ATM_flux_int ( ATM.wsnow_ter ) ATM.flx_wsnow_oce = ATM_flux_int ( ATM.wsnow_oce ) ATM.flx_wsnow_lic = ATM_flux_int ( ATM.wsnow_lic ) ATM.flx_wsnow_sic = ATM_flux_int ( ATM.wsnow_sic ) ATM.flx_wsnow_sea = ATM_flux_int ( ATM.wsnow_sea ) ATM.flx_wevap_ter = ATM_flux_int ( ATM.wevap_ter ) ATM.flx_wevap_oce = ATM_flux_int ( ATM.wevap_oce ) ATM.flx_wevap_lic = ATM_flux_int ( ATM.wevap_lic ) ATM.flx_wevap_sic = ATM_flux_int ( ATM.wevap_sic ) ATM.flx_wevap_sea = ATM_flux_int ( ATM.wevap_sea ) ATM.flx_wprecip_lic = ATM_flux_int ( ATM.wprecip_lic ) ATM.flx_wprecip_oce = ATM_flux_int ( ATM.wprecip_oce ) ATM.flx_wprecip_sic = ATM_flux_int ( ATM.wprecip_sic ) ATM.flx_wprecip_ter = ATM_flux_int ( ATM.wprecip_ter ) ATM.flx_wprecip_sea = ATM_flux_int ( ATM.wprecip_sea ) ATM.flx_wemp_lic = ATM_flux_int ( ATM.wemp_lic ) ATM.flx_wemp_oce = ATM_flux_int ( ATM.wemp_oce ) ATM.flx_wemp_sic = ATM_flux_int ( ATM.wemp_sic ) ATM.flx_wemp_ter = ATM_flux_int ( ATM.wemp_ter ) ATM.flx_wemp_sea = ATM_flux_int ( ATM.wemp_sea ) ATM.flx_emp = ATM_flux_int ( ATM.emp ) if SECHIBA : SRF.RUN_flx_coastal = ONE_flux_int ( SRF.RUN_coastalflow) SRF.RUN_flx_river = ONE_flux_int ( SRF.RUN_riverflow ) SRF.RUN_flx_coastal_cpl = ONE_flux_int ( SRF.RUN_coastalflow_cpl) SRF.RUN_flx_river_cpl = ONE_flux_int ( SRF.RUN_riverflow_cpl ) SRF.RUN_flx_drainage = SRF_flux_int ( SRF.RUN_drainage ) SRF.RUN_flx_riversret = SRF_flux_int ( SRF.RUN_riversret ) SRF.RUN_flx_runoff = SRF_flux_int ( SRF.RUN_runoff ) SRF.RUN_flx_input = SRF_flux_int ( SRF.RUN_input ) SRF.RUN_flx_output = ONE_flux_int ( SRF.RUN_output ) SRF.RUN_flx_bil = ONE_flux_int ( SRF.RUN_input - SRF.RUN_output) SRF.RUN_flx_rivcoa = ONE_flux_int ( SRF.RUN_coastalflow + SRF.RUN_riverflow) prtFlux ('wbilo_oce ', ATM.flx_wbilo_oce , 'f' ) prtFlux ('wbilo_sic ', ATM.flx_wbilo_sic , 'f' ) prtFlux ('wbilo_sic+oce ', ATM.flx_wbilo_sea , 'f' ) prtFlux ('wbilo_ter ', ATM.flx_wbilo_ter , 'f' ) prtFlux ('wbilo_lic ', ATM.flx_wbilo_lic , 'f' ) prtFlux ('Sum wbilo_* ', ATM.flx_wbilo , 'f', True) prtFlux ('E-P ', ATM.flx_emp , 'f', True) prtFlux ('calving ', ATM.flx_calving , 'f' ) prtFlux ('fqfonte ', ATM.flx_fqfonte , 'f' ) prtFlux ('precip ', ATM.flx_precip , 'f' ) prtFlux ('snowf ', ATM.flx_snowf , 'f' ) prtFlux ('evap ', ATM.flx_evap , 'f' ) prtFlux ('runoff lic ', ATM.flx_runlic , 'f' ) prtFlux ('ATM.flx_wevap* ', ATM.flx_wevap , 'f' ) prtFlux ('ATM.flx_wrain* ', ATM.flx_wrain , 'f' ) prtFlux ('ATM.flx_wsnow* ', ATM.flx_wsnow , 'f' ) prtFlux ('ATM.flx_wprecip* ', ATM.flx_wprecip , 'f' ) prtFlux ('ATM.flx_wemp* ', ATM.flx_wemp , 'f', True ) prtFlux ('ERROR evap ', ATM.flx_wevap - ATM.flx_evap , 'e', True ) prtFlux ('ERROR precip ', ATM.flx_wprecip - ATM.flx_precip, 'e', True ) prtFlux ('ERROR snow ', ATM.flx_wsnow - ATM.flx_snowf , 'e', True ) prtFlux ('ERROR emp ', ATM.flx_wemp - ATM.flx_emp , 'e', True ) if SECHIBA : echo ( '\n====================================================================================' ) echo ( f'-- RUNOFF Fluxes -- {Title} ' ) prtFlux ('coastalflow ', SRF.RUN_flx_coastal , 'f' ) prtFlux ('riverflow ', SRF.RUN_flx_river , 'f' ) prtFlux ('coastal_cpl ', SRF.RUN_flx_coastal_cpl, 'f' ) prtFlux ('riverf_cpl ', SRF.RUN_flx_river_cpl , 'f' ) prtFlux ('river+coastal ', SRF.RUN_flx_rivcoa , 'f' ) prtFlux ('drainage ', SRF.RUN_flx_drainage , 'f' ) prtFlux ('riversret ', SRF.RUN_flx_riversret , 'f' ) prtFlux ('runoff ', SRF.RUN_flx_runoff , 'f' ) prtFlux ('river in ', SRF.RUN_flx_input , 'f' ) prtFlux ('river out ', SRF.RUN_flx_output , 'f' ) prtFlux ('river bil ', SRF.RUN_flx_bil , 'f' ) echo ( '\n====================================================================================' ) echo ( f'-- OCE Fluxes -- {Title} ' ) # Read variable and computes integral over space and time OCE.OCE_empmr = rprec (d_OCE_his['wfo'] ) ; OCE.OCE_mass_empmr = OCE_flux_int ( OCE.OCE_empmr ) OCE.OCE_wfob = rprec (d_OCE_his['wfob'] ) ; OCE.OCE_mass_wfob = OCE_flux_int ( OCE.OCE_wfob ) OCE.OCE_emp_oce = rprec (d_OCE_his['emp_oce'] ) ; OCE.OCE_mass_emp_oce = OCE_flux_int ( OCE.OCE_emp_oce ) OCE.OCE_emp_ice = rprec (d_OCE_his['emp_ice'] ) ; OCE.OCE_mass_emp_ice = OCE_flux_int ( OCE.OCE_emp_ice ) OCE.OCE_iceshelf = rprec (d_OCE_his['iceshelf']) ; OCE.OCE_mass_iceshelf = OCE_flux_int ( OCE.OCE_iceshelf ) OCE.OCE_calving = rprec (d_OCE_his['calving'] ) ; OCE.OCE_mass_calving = OCE_flux_int ( OCE.OCE_calving ) OCE.OCE_iceberg = rprec (d_OCE_his['iceberg'] ) ; OCE.OCE_mass_iceberg = OCE_flux_int ( OCE.OCE_iceberg ) OCE.OCE_friver = rprec (d_OCE_his['friver'] ) ; OCE.OCE_mass_friver = OCE_flux_int ( OCE.OCE_friver ) OCE.OCE_runoffs = rprec (d_OCE_his['runoffs'] ) ; OCE.OCE_mass_runoffs = OCE_flux_int ( OCE.OCE_runoffs ) if NEMO == 4.0 or NEMO == 4.2 : OCE.OCE_wfxice = rprec (d_OCE_his['vfxice']) ; OCE.OCE_mass_wfxice = OCE_flux_int ( OCE.OCE_wfxice ) OCE.OCE_wfxsnw = rprec (d_OCE_his['vfxsnw']) ; OCE.OCE_mass_wfxsnw = OCE_flux_int ( OCE.OCE_wfxsnw ) OCE.OCE_wfxsub = rprec (d_OCE_his['vfxsub']) ; OCE.OCE_mass_wfxsub = OCE_flux_int ( OCE.OCE_wfxsub ) if NEMO == 3.6 : OCE.OCE_wfxice = rprec (d_OCE_his['vfxice'])/86400.*OCE.ICE_RHO_ICE ; OCE.OCE_mass_wfxice = OCE_flux_int ( OCE.OCE_wfxice ) OCE.OCE_wfxsnw = rprec (d_OCE_his['vfxsnw'])/86400.*OCE.ICE_RHO_SNO ; OCE.OCE_mass_wfxsnw = OCE_flux_int ( OCE.OCE_wfxsnw ) OCE.OCE_wfxsub = rprec (d_OCE_his['vfxsub'])/86400.*ICE_RHO_SNO ; OCE.OCE_mass_wfxsub = OCE_flux_int ( OCE.OCE_wfxsub ) # Additional checks OCE.OCE_evap_oce = rprec (d_OCE_his['evap_ao_cea']) ; OCE.OCE_mass_evap_oce = OCE_flux_int ( OCE.OCE_evap_oce ) OCE.ICE_evap_ice = rprec (d_OCE_his['subl_ai_cea']) ; OCE.ICE_mass_evap_ice = OCE_flux_int ( OCE.ICE_evap_ice ) OCE.OCE_snow_oce = rprec (d_OCE_his['snow_ao_cea']) ; OCE.OCE_mass_snow_oce = OCE_flux_int ( OCE.OCE_snow_oce ) OCE.OCE_snow_ice = rprec (d_OCE_his['snow_ai_cea']) ; OCE.OCE_mass_snow_ice = OCE_flux_int ( OCE.OCE_snow_ice ) OCE.OCE_rain = rprec (d_OCE_his['rain'] ) ; OCE.OCE_mass_rain = OCE_flux_int ( OCE.OCE_rain ) OCE.ICE_wfxsub_err = rprec (d_ICE_his['vfxsub_err'] ) ; OCE.ICE_mass_wfxsub_err = OCE_flux_int ( OCE.ICE_wfxsub_err ) if NEMO == 4.0 or NEMO == 4.2 : OCE.ICE_wfxpnd = rprec (d_ICE_his['vfxpnd'] ) ; OCE.ICE_mass_wfxpnd = OCE_flux_int ( OCE.ICE_wfxpnd ) OCE.ICE_wfxsnw_sub = rprec (d_ICE_his['vfxsnw_sub']) ; OCE.ICE_mass_wfxsnw_sub = OCE_flux_int ( OCE.ICE_wfxsnw_sub ) OCE.ICE_wfxsnw_pre = rprec (d_ICE_his['vfxsnw_pre']) ; OCE.ICE_mass_wfxsnw_pre = OCE_flux_int ( OCE.ICE_wfxsnw_pre ) if NEMO == 3.6 : OCE.ICE_wfxpnd = 0.0 ; OCE.ICE_mass_wfxpnd = 0.0 OCE.ICE_wfxsnw_sub = rprec (d_ICE_his['vfxsub'])/86400.*ICE_RHO_SNO ; OCE.ICE_mass_wfxsnw_sub = OCE_flux_int ( OCE.ICE_wfxsnw_sub ) OCE.ICE_wfxsnw_pre = rprec (d_ICE_his['vfxspr'])/86400.*ICE_RHO_SNO ; OCE.ICE_mass_wfxsnw_pre = OCE_flux_int ( OCE.ICE_wfxsnw_pre ) OCE.OCE_wfcorr = rprec (d_OCE_his['wfcorr']) ; OCE.OCE_mass_wfcorr = OCE_flux_int ( OCE.OCE_wfcorr ) if OCE_relax : # ssr and fwb are included in emp=>empmr but not in emp_oce (outputed by sea-ice) OCE.OCE_vflx_fwb = rprec (d_OCE_his['vflx_fwb']) ; OCE.OCE_mass_vflx_fwb = OCE_flux_int ( OCE.OCE_vflx_fwb ) OCE.OCE_vflx_ssr = rprec (d_OCE_his['vflx_ssr']) ; OCE.OCE_mass_vflx_ssr = OCE_flux_int ( OCE.OCE_vflx_ssr ) else : OCE.OCE_fwb = 0.0 ; OCE.OCE_mass_fwb = 0.0 OCE.OCE_ssr = 0.0 ; OCE.OCE_mass_ssr = 0.0 if OCE_icb : OCE.OCE_berg_icb = rprec (d_OCE_his['berg_floating_melt']) ; OCE.OCE_mass_berg_icb = OCE_flux_int ( OCE.OCE_berg_icb ) OCE.OCE_calving_icb = rprec (d_OCE_his['calving_icb'] ) ; OCE.OCE_mass_calv_icb = OCE_flux_int ( OCE.OCE_calving_icb ) else : OCE.OCE_berg_icb = 0. ; OCE.OCE_mass_berg_icb = 0. OCE.OCE_calv_icb = 0. ; OCE.OCE_mass_calv_icb = 0. OCE.OCE_mass_emp = OCE.OCE_mass_emp_oce - OCE.OCE_mass_wfxice - OCE.OCE_mass_wfxsnw - OCE.ICE_mass_wfxpnd - OCE.ICE_mass_wfxsub_err OCE.OCE_mass_all = OCE.OCE_mass_emp_oce + OCE.OCE_mass_emp_ice - OCE.OCE_mass_runoffs - OCE.OCE_mass_iceshelf prtFlux ('OCE+ICE budget ', OCE.OCE_mass_all , 'e', True) prtFlux (' EMPMR ', OCE.OCE_mass_empmr , 'e', True) prtFlux (' WFOB ', OCE.OCE_mass_wfob , 'e', True) prtFlux (' EMP_OCE ', OCE.OCE_mass_emp_oce , 'e', True) prtFlux (' EMP_ICE ', OCE.OCE_mass_emp_ice , 'e', True) prtFlux (' EMP ', OCE.OCE_mass_emp , 'e', True) prtFlux (' ICEBERG ', OCE.OCE_mass_iceberg , 'e', ) prtFlux (' ICESHELF ', OCE.OCE_mass_iceshelf , 'e', True) prtFlux (' CALVING ', OCE.OCE_mass_calving , 'e', True) prtFlux (' FRIVER ', OCE.OCE_mass_friver , 'e', ) prtFlux (' RUNOFFS ', OCE.OCE_mass_runoffs , 'e', True) prtFlux (' WFXICE ', OCE.OCE_mass_wfxice , 'e', True) prtFlux (' WFXSNW ', OCE.OCE_mass_wfxsnw , 'e', True) prtFlux (' WFXSUB ', OCE.OCE_mass_wfxsub , 'e', True) prtFlux (' WFXPND ', OCE.ICE_mass_wfxpnd , 'e', True) prtFlux (' WFXSNW_SUB ', OCE.ICE_mass_wfxsnw_sub, 'e', True) prtFlux (' WFXSNW_PRE ', OCE.ICE_mass_wfxsnw_pre, 'e', True) prtFlux (' WFXSUB_ERR ', OCE.ICE_mass_wfxsub_err, 'e', True) prtFlux (' EVAP_OCE ', OCE.OCE_mass_evap_oce , 'e' ) prtFlux (' EVAP_ICE ', OCE.ICE_mass_evap_ice , 'e', True) prtFlux (' SNOW_OCE ', OCE.OCE_mass_snow_oce , 'e', True) prtFlux (' SNOW_ICE ', OCE.OCE_mass_snow_ice , 'e', True) prtFlux (' RAIN ', OCE.OCE_mass_rain , 'e' ) prtFlux (' FWB ', OCE.OCE_mass_fwb , 'e', True) prtFlux (' SSR ', OCE.OCE_mass_ssr , 'e', True) prtFlux (' WFCORR ', OCE.OCE_mass_wfcorr , 'e', True) prtFlux (' BERG_ICB ', OCE.OCE_mass_berg_icb , 'e', True) prtFlux (' CALV_ICB ', OCE.OCE_mass_calv_icb , 'e', True) echo (' ') prtFlux ( 'wbilo sea ', ATM_flux_int (ATM.wbilo_sea), 'e', ) if SECHIBA : prtFlux ( 'costalflow ', ONE_flux_int (SRF.RUN_coastalflow), 'e', ) prtFlux ( 'riverflow ', SRF.RUN_flx_river , 'e', ) prtFlux ( 'costalflow ', SRF.RUN_flx_coastal, 'e', ) prtFlux ( 'runoff ', SRF.RUN_flx_river+SRF.RUN_flx_coastal, 'e', ) #ATM.to_OCE = ATM_flux_int (ATM.wbilo_sea) - SRF.RUN_flx_river - SRF.RUN_flx_coastal - ATM.flx_calving ATM.to_OCE = ATM_flux_int (ATM.wbilo_sea) - SRF.RUN_flx_river - SRF.RUN_flx_coastal - ATM.flx_calving #OCE.OCE_from_ATM = -OCE.OCE_mass_emp_oce - OCE.OCE_mass_emp_ice + OCE.OCE_mass_runoffs + OCE.OCE_mass_iceberg + OCE.OCE_mass_calving + OCE.OCE_mass_iceshelf OCE.OCE_from_ATM = OCE.OCE_mass_all prtFlux ( 'ATM.to_OCE ', ATM.to_OCE , 'e', True ) prtFlux ( 'OCE.OCE_from_ATM', OCE.OCE_from_ATM, 'e', True ) echo ( ' ' ) echo ( f'{Title = }' ) echo ( 'SVN Information' ) for kk in SVN.keys(): print ( SVN[kk] ) ## Write the full configuration ## ---------------------------- params_out = open (FullIniFile, 'w', encoding = 'utf-8') params = wu.dict2config ( dpar ) params.write ( params_out ) params_out.close ()