diff --git a/A-INSTALL b/A-INSTALL index 706b184e120f08545f02dd0d2703dd195fd41e01..3fd5df2b26cbe7ef3e56e7bee8e70989c01640d2 100644 --- a/A-INSTALL +++ b/A-INSTALL @@ -1,8 +1,8 @@ # # Version of PACKAGE MESONH "Open distribution" -# PACK-MNH-V5-5-0 -# DATE : 20/07/2021 -# VERSION : MESONH MASDEV5_5 + BUG-0 +# PACK-MNH-V5-5-1 +# DATE : 09/06/2022 +# VERSION : MESONH MASDEV5_5 + BUG-1 # # MAP # @@ -84,14 +84,14 @@ # # or directly # -# http://mesonh.aero.obs-mip.fr/mesonh/dir_open/dir_MESONH/MNH-V5-5-0.tar.gz +# http://mesonh.aero.obs-mip.fr/mesonh/dir_open/dir_MESONH/MNH-V5-5-1.tar.gz # -# Then untar the file "MNH-V5-5-0.tar.gz" where you want to. +# Then untar the file "MNH-V5-5-1.tar.gz" where you want to. # For example, in your home directory: # cd ~ -tar xvfz MNH-V5-5-0.tar.gz +tar xvfz MNH-V5-5-1.tar.gz # # Process now to the chapter to configure the MesoNH package. @@ -171,10 +171,10 @@ git config --global http.sslverify false # Finally you can clone the Meso-NH Git repository with the following command: # -git lfs clone anongit@anongit_mesonh:/gitrepos/MNH-git_open_source-lfs.git -b MNH-55-branch MNH-V5-5-0 +git lfs clone anongit@anongit_mesonh:/gitrepos/MNH-git_open_source-lfs.git -b MNH-55-branch MNH-V5-5-1 # -# that will create the MNH-V5-5-0 directory containing a clone (copy) of the +# that will create the MNH-V5-5-1 directory containing a clone (copy) of the # Meso-NH package on the remote developpement branch MNH-55-branch # # @@ -184,16 +184,16 @@ git lfs clone anongit@anongit_mesonh:/gitrepos/MNH-git_open_source-lfs.git -b MN # Once the repository is cloned, it's better for you to checkout your own branch # (by default, you are on HEAD of the MNH-55-branch development branch ). # -# To create your local branch corresponding to the V5-5-0 version, type: +# To create your local branch corresponding to the V5-5-1 version, type: # -cd MNH-V5-5-0 -git checkout -b MYB-MNH-V5-5-0 PACK-MNH-V5-5-0 +cd MNH-V5-5-1 +git checkout -b MYB-MNH-V5-5-1 PACK-MNH-V5-5-1 # -# MYB-MNH-V5-5-0 is the name of the local branch you created +# MYB-MNH-V5-5-1 is the name of the local branch you created # and -# PACK-MNH-V5-5-0 is the remote/origin tag on which it is based. +# PACK-MNH-V5-5-1 is the remote/origin tag on which it is based. # # The advantage of this way of downloading the package is that in the future # you could check/update quickly differences with the new version of the @@ -257,7 +257,7 @@ git clone anongit@anongit_mesonh:/gitrepos/MNH-DOC.git # use the "./configure" script like this # -cd ~/MNH-V5-5-0/src +cd ~/MNH-V5-5-1/src ./configure . ../conf/profile_mesonh @@ -308,7 +308,7 @@ export OPTLEVEL=O2 # Compile in O2, 4 times faster then DEBUG, but less # and then source/load the new generate file -. ../conf/profile_mesonh.LXifort.MNH-V5-5-0.MPIAUTO.O2 +. ../conf/profile_mesonh.LXifort.MNH-V5-5-1.MPIAUTO.O2 # # REM: @@ -333,7 +333,7 @@ export OPTLEVEL=O2 # Compile in O2, 4 times faster then DEBUG, but less # go to the directory "src" # -cd ~/MNH-V5-5-0/src +cd ~/MNH-V5-5-1/src # # if you have not already configured your MESONH environment @@ -562,7 +562,7 @@ make examples # cd $WORKDIR -cd MNH-V5-5-0/src +cd MNH-V5-5-1/src ./configure @@ -576,8 +576,8 @@ cd MNH-V5-5-0/src # - On JEAN-ZAY ( HPE ) the compilation is only possible in interactive : -cd MNH-V5-5-0/src -. ../conf/profile_mesonh-LXifort-R8I4-MNH-V5-5-0-MPIINTEL-O2 +cd MNH-V5-5-1/src +. ../conf/profile_mesonh-LXifort-R8I4-MNH-V5-5-1-MPIINTEL-O2 make -j16 |& tee error$XYZ make installmaster @@ -608,12 +608,12 @@ sbatch job_make_examples_BullX_occigen # - ssh irene-fr : for Intel SkyLake/KNL processors # On Intel processors the MPI use is OPENMPI/2.0.4 # the configure will generate a -# profile_mesonh-LXifort-R8I4-MNH-V5-5-0-MPIAUTO-O2 +# profile_mesonh-LXifort-R8I4-MNH-V5-5-1-MPIAUTO-O2 # # - ssh irene-amd : for AMD , processors # On AMD processors the MPI use is OPENMPI/4.02 # the configure will generate a -# profile_mesonh-LXifort-R8I4-MNH-V5-5-0-AMD-MPIAUTO-O2 +# profile_mesonh-LXifort-R8I4-MNH-V5-5-1-AMD-MPIAUTO-O2 # # - install the PACKAGE in your $CCCHOME ( default 20Go of quota ) # - Compile in interactive mode ( see IDRIS ) @@ -654,14 +654,14 @@ qsub job_make_examples_CRAY_cca # - At Meteo-France DSI on beaufix (or prolix) # # to install the whole package on your "$HOME" directory -# untar the file "MNH-V5-5-0.tar.gz" from its location : +# untar the file "MNH-V5-5-1.tar.gz" from its location : cd ~ -tar xvf $MESONH/MNH-V5-5-0.tar.gz +tar xvf $MESONH/MNH-V5-5-1.tar.gz # run the "./configure" command : -cd MNH-V5-5-0/src +cd MNH-V5-5-1/src ./configure # @@ -761,7 +761,7 @@ scandollar ## OUTPUT :: -># read default config file :: ---> CONF_DOLLAR=/home/escj/DEV64/PACK-MNH-V5-5-0/conf/post/confdollar_aeropc_default +># read default config file :: ---> CONF_DOLLAR=/home/escj/DEV64/PACK-MNH-V5-5-1/conf/post/confdollar_aeropc_default ># ># read user config file :: ---> CONFIG=confdollar ># @@ -783,7 +783,7 @@ scandollar 0* ## OUTPUT :: ># -># read default config file :: ---> CONF_DOLLAR=/home/escj/DEV64/PACK-MNH-V5-5-0/conf/post/confdollar_aeropc_default +># read default config file :: ---> CONF_DOLLAR=/home/escj/DEV64/PACK-MNH-V5-5-1/conf/post/confdollar_aeropc_default ># ># read user config file :: ---> CONFIG=confdollar ># @@ -857,22 +857,22 @@ cp -R 007_16janvier_scandollar /.../your_directory # # use this "profile_mesonh" : -. /home/rech/mnh/rmnh007/DEV/MNH-V5-5-0/conf/profile_mesonh-SX8-MNH-V5-5-0-MPIAUTO-O4 +. /home/rech/mnh/rmnh007/DEV/MNH-V5-5-1/conf/profile_mesonh-SX8-MNH-V5-5-1-MPIAUTO-O4 # And the examples are here ( link to my $WORKDIR in actually ) -/home/rech/mnh/rmnh007/DEV/MNH-V5-5-0/MY_RUN/KTEST/007_16janvier_scandollar +/home/rech/mnh/rmnh007/DEV/MNH-V5-5-1/MY_RUN/KTEST/007_16janvier_scandollar # # On vargas # --------- # use this "profile_mesonh" : -. /workgpfs/rech/mnh/rmnh007/DEV/MNH-V5-5-0/conf/profile_mesonh-AIX64-MNH-V5-5-0-MPIAUTO-O2 +. /workgpfs/rech/mnh/rmnh007/DEV/MNH-V5-5-1/conf/profile_mesonh-AIX64-MNH-V5-5-1-MPIAUTO-O2 # and examples here : -/workgpfs/rech/mnh/rmnh007/DEV/MNH-V5-5-0/MY_RUN/KTEST/007_16janvier_scandollar +/workgpfs/rech/mnh/rmnh007/DEV/MNH-V5-5-1/MY_RUN/KTEST/007_16janvier_scandollar # # - At CINES on JADE : @@ -880,11 +880,11 @@ cp -R 007_16janvier_scandollar /.../your_directory # # use -. /work/escobar/DEV/MNH-V5-5-0/conf/profile_mesonh-LXifort-MNH-V5-5-0-MPIICE-O2 +. /work/escobar/DEV/MNH-V5-5-1/conf/profile_mesonh-LXifort-MNH-V5-5-1-MPIICE-O2 # and the exemples -/work/escobar/DEV/MNH-V5-5-0/MY_RUN/KTEST/007_16janvier_scandollar +/work/escobar/DEV/MNH-V5-5-1/MY_RUN/KTEST/007_16janvier_scandollar # # - At ECMWF on cxa : @@ -892,11 +892,11 @@ cp -R 007_16janvier_scandollar /.../your_directory # # use -. /c1a/ms_perm/au5/MNH-V5-5-0/conf/profile_mesonh-AIX64-MNH-V5-5-0-MPIAUTO-O2 +. /c1a/ms_perm/au5/MNH-V5-5-1/conf/profile_mesonh-AIX64-MNH-V5-5-1-MPIAUTO-O2 # and the examples -/c1a/ms_perm/au5/MNH-V5-5-0/MY_RUN/KTEST/007_16janvier_scandollar +/c1a/ms_perm/au5/MNH-V5-5-1/MY_RUN/KTEST/007_16janvier_scandollar # diff --git a/MY_RUN/KTEST/004_Reunion/007_python/plot_004_Reunion.py b/MY_RUN/KTEST/004_Reunion/007_python/plot_004_Reunion.py index 5de8f4b897fab5ecd28427472f553f42100b6f99..99bcf6ff02731847abf571b96a571cd7062603f5 100644 --- a/MY_RUN/KTEST/004_Reunion/007_python/plot_004_Reunion.py +++ b/MY_RUN/KTEST/004_Reunion/007_python/plot_004_Reunion.py @@ -23,8 +23,7 @@ os.system('rm -f tempgraph*') LnameFiles = ['REUNI.1.00A20.004dia.nc', 'REUNI.1.00A20.004.nc'] Dvar_input = { -'f1':['ZS', 'UT', 'VT', 'WT', 'THT', 'ALT_PRESSURE','ALT_U','ALT_V','ALT_THETA','level','ZTOP', 'longitude','latitude','level_w','time', - 'RN', 'H','LE','GFLUX','HU2M','T2M','W10M','CD','CH','Z0','CE','TS','Z0H'], +'f1':['ZS', 'UT', 'VT', 'WT', 'THT', 'ALT_PRESSURE','ALT_U','ALT_V','ALT_THETA','level','ZTOP', 'longitude','latitude','level_w','time'], 'f2':['LSTHM', 'LSVM']} # Read the variables in the files diff --git a/MY_RUN/KTEST/013_Iroise_ideal_case_coupling/2_INPUT_TOY/Create_grid_and_restart_files_for_TOY.py b/MY_RUN/KTEST/013_Iroise_ideal_case_coupling/2_INPUT_TOY/Create_grid_and_restart_files_for_TOY.py index 153d1fcb664f3ff6a66fbc4749a3da7cd1b24bed..7df39f10579aade64aea10e546571cd89c156ee7 100755 --- a/MY_RUN/KTEST/013_Iroise_ideal_case_coupling/2_INPUT_TOY/Create_grid_and_restart_files_for_TOY.py +++ b/MY_RUN/KTEST/013_Iroise_ideal_case_coupling/2_INPUT_TOY/Create_grid_and_restart_files_for_TOY.py @@ -6,6 +6,8 @@ # Creating grid and restart file for toy model # Author : J. Pianezze # Date : 2015 +# Modification :Avril 2022 Joris Pianezze et Françoise Orain pour python3 +#fichier etopo2 dispo sur site SURFEX #=================================================# ################################################### # @@ -14,7 +16,7 @@ import numpy as np import scipy import matplotlib.pyplot as plt import math -from pylab import * +from pylab import * import os # curdir_path=os.path.abspath(os.curdir)+'/' @@ -24,8 +26,10 @@ curdir_path=os.path.abspath(os.curdir)+'/' #$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$ # #-- Limit of the grid (from etopo) -lat_domain=[47.0, 49.5] -lon_domain=[-6.2, -4.0] +lat_domain=[ 46.681, 50.092] +lon_domain=[-6.6348, -1.3652] +print ('londomaine[0]',lon_domain[0]) +print ('londomaine[1]',lon_domain[1]) #-- Type of forcing to create the restart file # for the toy : CNSTE or SINUS @@ -43,63 +47,63 @@ value_SINUS_LENGTH=1000. file_topo = netCDF4.Dataset('topo.nc') #------ Read variables -lon_full=file_topo.variables['lon'][:] -lat_full=file_topo.variables['lat'][:] -topo_full=file_topo.variables['topo'][:,:] +lon_full = file_topo.variables['lon'] [:] +lat_full = file_topo.variables['lat'] [:] +topo_full = file_topo.variables['topo'][:,:] -ind_min_lon=find(lon_full[:]>lon_domain[0])[0] ; print ind_min_lon -ind_max_lon=find(lon_full[:]>lon_domain[1])[0] ; print ind_max_lon -ind_min_lat=find(lat_full[:]>lat_domain[0])[0] ; print ind_min_lat -ind_max_lat=find(lat_full[:]>lat_domain[1])[0] ; print ind_max_lat +ind_min_lon = np.where(lon_full[:]>lon_domain[0])[0][0] ; print ("ind_min_lon = ",ind_min_lon) +ind_max_lon = np.where(lon_full[:]>lon_domain[1])[0][0] ; print ("ind_max_lon = ",ind_max_lon) +ind_min_lat = np.where(lat_full[:]>lat_domain[0])[0][0] ; print ("ind_min_lat = ",ind_min_lat) +ind_max_lat = np.where(lat_full[:]>lat_domain[1])[0][0] ; print ("ind_max_lat = ",ind_max_lat) -lon=lon_full[ind_min_lon:ind_max_lon] -lat=lat_full[ind_min_lat:ind_max_lat] -topo=topo_full[ind_min_lat:ind_max_lat,ind_min_lon:ind_max_lon] +lon = lon_full [ind_min_lon:ind_max_lon] +lat = lat_full [ind_min_lat:ind_max_lat] +topo = topo_full[ind_min_lat:ind_max_lat,ind_min_lon:ind_max_lon] #plt.contourf(topo) #plt.show() -nlon=np.size(lon) ; print 'nlon=', nlon -nlat=np.size(lat) ; print 'nlat=', nlat -ncorn=4 ; print 'ncorn=', ncorn +nlon=np.size(lon) ; print ('nlon =', nlon) +nlat=np.size(lat) ; print ('nlat =', nlat) +ncorn=4 ; print ('ncorn=', ncorn) -print '---- longitude/latitude' +print ('---- longitude/latitude') lon2D=np.zeros((nlat,nlon)) lat2D=np.zeros((nlat,nlon)) -for ind_lon in xrange(nlon): +for ind_lon in range(nlon): lat2D[:,ind_lon]=lat[:] -for ind_lat in xrange(nlat): +for ind_lat in range(nlat): lon2D[ind_lat,:]=lon[:] -print '---- corners longitude/latitude' +print ('---- corners longitude/latitude') clo=np.zeros((ncorn,nlat,nlon)) cla=np.zeros((ncorn,nlat,nlon)) -deltax=lon[1]-lon[0] ; print 'deltax=', deltax +deltax=lon[1]-lon[0] ; print ('deltax=', deltax) clo[0,:,:]=lon2D[:,:]+deltax/2.0 clo[1,:,:]=lon2D[:,:]-deltax/2.0 clo[2,:,:]=lon2D[:,:]-deltax/2.0 clo[3,:,:]=lon2D[:,:]+deltax/2.0 -deltay=lat[1]-lat[0] ; print 'deltay=', deltay +deltay=lat[1]-lat[0] ; print ('deltay=', deltay) cla[0,:,:]=lat2D[:,:]+deltay/2.0 cla[1,:,:]=lat2D[:,:]+deltay/2.0 cla[2,:,:]=lat2D[:,:]-deltay/2.0 cla[3,:,:]=lat2D[:,:]-deltay/2.0 -print '---- surface' +print ('---- surface') surface=np.zeros((nlat,nlon)) surface[:,:]=deltax*deltay -print '---- mask and var send by toy' +print ('---- mask and var send by toy') mask=np.zeros((nlat,nlon)) toyvarcnste=np.zeros((nlat,nlon)) toyvarsinus=np.zeros((nlat,nlon)) -for ind_lon in xrange(nlon): - for ind_lat in xrange(nlat): +for ind_lon in range(nlon): + for ind_lat in range(nlat): if topo[ind_lat,ind_lon] > 0.0 : mask[ind_lat,ind_lon]=0 toyvarcnste[ind_lat,ind_lon] = value_CNSTE @@ -111,8 +115,8 @@ for ind_lon in xrange(nlon): ################################################## -print '------------------------------------------' -print ' Creating netcdf file : grid_toy_model.nc' +print ('------------------------------------------') +print (' Creating netcdf file : grid_toy_model.nc ') grid_file=netCDF4.Dataset(curdir_path+'grid_toy_model.nc','w',format='NETCDF3_64BIT') grid_file.Description='Grid file for OASIS coupling' @@ -149,13 +153,13 @@ grid_file.variables['imask'][:,:] = mask[:,:] # --------------------------------------- grid_file.close() -print ' Closing netcdf file : grid_toy_model.nc' -print '------------------------------------------' +print (' Closing netcdf file : grid_toy_model.nc ') +print ('------------------------------------------') ################################################## ################################################## -print '------------------------------------------' -print ' Creating netcdf file : rstrt_TOY.nc' +print ('------------------------------------------') +print (' Creating netcdf file : rstrt_TOY.nc ') rstrt_file=netCDF4.Dataset(curdir_path+'rstrt_TOY.nc','w',format='NETCDF3_64BIT') rstrt_file.Description='Restart file for OASIS coupling' @@ -185,6 +189,6 @@ rstrt_file.variables['VARSIN02'][:,:] = toyvarsinus[:,:] # --------------------------------------- rstrt_file.close() -print ' Closing netcdf file : rstrt_TOY.nc' -print '-----------------------------------------' +print (' Closing netcdf file : rstrt_TOY.nc ') +print ('-----------------------------------------') ##################################################### diff --git a/MY_RUN/KTEST/013_Iroise_ideal_case_coupling/2_INPUT_TOY/Create_grid_and_restart_files_for_TOYFOJP_from_MARS.py b/MY_RUN/KTEST/013_Iroise_ideal_case_coupling/2_INPUT_TOY/Create_grid_and_restart_files_for_TOYFOJP_from_MARS.py new file mode 100644 index 0000000000000000000000000000000000000000..708ab094c3336a2ec766bac64f5fb5e416a26407 --- /dev/null +++ b/MY_RUN/KTEST/013_Iroise_ideal_case_coupling/2_INPUT_TOY/Create_grid_and_restart_files_for_TOYFOJP_from_MARS.py @@ -0,0 +1,191 @@ +#!/usr/bin/python +# -*- coding: utf-8 -*- +# +################################################### +#=================================================# +# Creating grid and restart file for toy model +# Author : J. Pianezze +# Date : 2015 +# Modification Avril 2022 :Françoise Orain et Joris Pianezze pour python3 et utiliser H0(bathy de MARC) +#=================================================# +################################################### +# +import netCDF4 +import numpy as np +import scipy +import matplotlib.pyplot as plt +import math +from pylab import * +import os +# +curdir_path=os.path.abspath(os.curdir)+'/' +# +#$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$ +# To be defined by the user +#$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$ +# +#-- Limit of the grid (from etopo) +lat_domain=[ 46.681, 50.092] +lon_domain=[-6.6348, -1.3652] +print ('londomaine[0]',lon_domain[0]) +print ('londomaine[1]',lon_domain[1]) + +#-- Type of forcing to create the restart file +# for the toy : CNSTE or SINUS + +# CNSTE +value_CNSTE=290.0 + +# SINUS +value_SINUS_COEF=0.011 +value_SINUS_LENGTH=1000. +# +#$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$ +#$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$ + +file_topo = netCDF4.Dataset('bathy_MARC.nc') + +#------ Read variables +lon_full = file_topo.variables['longitude'][:,:] +lat_full = file_topo.variables['latitude'][:,:] +topo_full = file_topo.variables['H0'][:,:] + +ind_min_lon=np.where(lon_full[0,:]>lon_domain[0])[0][0] ; print ("ind_min_lon = ",ind_min_lon) +ind_max_lon=-1 # np.where(lon_full[0,:]>lon_domain[1])[0][0] ; print ("ind_max_lon = ",ind_max_lon) +ind_min_lat=np.where(lat_full[:,0]>lat_domain[0])[0][0] ; print ("ind_min_lat = ",ind_min_lat) +ind_max_lat=np.where(lat_full[:,0]>lat_domain[1])[0][0] ; print ("ind_max_lat = ",ind_max_lat) + +print ('---- topo') +topo = topo_full[ind_min_lat:ind_max_lat,ind_min_lon:ind_max_lon] + +print ('---- longitude/latitude') +lon = lon_full[ind_min_lat:ind_max_lat,ind_min_lon:ind_max_lon] +lat = lat_full[ind_min_lat:ind_max_lat,ind_min_lon:ind_max_lon] + +nlon=np.size(lon[0,:]) ; print ('nlon =', nlon) +nlat=np.size(lat[:,0]) ; print ('nlat =', nlat) +ncorn=4 ; print ('ncorn=', ncorn) + +print ('---- corners longitude/latitude') +clo=np.zeros((ncorn,nlat,nlon)) +cla=np.zeros((ncorn,nlat,nlon)) + +deltax=lon[0,1]-lon[0,0] ; print ('deltax=', deltax) +clo[0,:,:]=lon[:,:]+deltax/2.0 +clo[1,:,:]=lon[:,:]-deltax/2.0 +clo[2,:,:]=lon[:,:]-deltax/2.0 +clo[3,:,:]=lon[:,:]+deltax/2.0 + +deltay=lat[1]-lat[0] ; print ('deltay=', deltay) +cla[0,:,:]=lat[:,:]+deltay/2.0 +cla[1,:,:]=lat[:,:]+deltay/2.0 +cla[2,:,:]=lat[:,:]-deltay/2.0 +cla[3,:,:]=lat[:,:]-deltay/2.0 + +print ('---- surface') +surface=np.zeros((nlat,nlon)) +surface[:,:]=deltax*deltay + + +print ('---- mask and var send by toy') +mask=np.zeros((nlat,nlon)) +toyvarcnste=np.zeros((nlat,nlon)) +toyvarsinus=np.zeros((nlat,nlon)) + +for ind_lon in range(nlon): + for ind_lat in range(nlat): + if topo[ind_lat,ind_lon] > 0.0 : + if ind_lon >= 132 and ind_lat >= 23 and ind_lat <= 25: + mask [ind_lat,ind_lon] = 1 + else: + mask[ind_lat,ind_lon]=0 + toyvarcnste[ind_lat,ind_lon] = value_CNSTE + toyvarsinus[ind_lat,ind_lon] = value_SINUS_COEF*math.sin(lat[ind_lat,0]*math.pi/180.0*value_SINUS_LENGTH) + else: + mask [ind_lat,ind_lon] = 1 + toyvarcnste[ind_lat,ind_lon] = value_CNSTE + toyvarsinus[ind_lat,ind_lon] = value_SINUS_COEF*math.sin(lat[ind_lat,0]*math.pi/180.0*value_SINUS_LENGTH) + +# ----------------------------------------------------------------------- +# Inverse du mask car la bathy est positive en mer et negative sur terre +# ----------------------------------------------------------------------- +mask = 1 - mask + +################################################## +print ('------------------------------------------') +print (' Creating netcdf file : grid_toy_modelMARC.nc') + +grid_file=netCDF4.Dataset(curdir_path+'grid_toy_modelMARC.nc','w',format='NETCDF3_64BIT') +grid_file.Description='Grid file for OASIS coupling' + +# ---------------------------------- +# Create the dimensions of the files +# ---------------------------------- +grid_file.createDimension ('nlon', nlon) +grid_file.createDimension ('nlat', nlat) +grid_file.createDimension ('ncorner', 4 ) + +# ---------------------------------- +# Create the variables of the files +# ---------------------------------- +varout=grid_file.createVariable('lon','d',('nlat','nlon')) +varout=grid_file.createVariable('lat','d',('nlat','nlon')) +varout=grid_file.createVariable('clo','d',('ncorner','nlat','nlon')) +varout=grid_file.createVariable('cla','d',('ncorner','nlat','nlon')) +varout=grid_file.createVariable('srf','d',('nlat','nlon')) +varout=grid_file.createVariable('imask','d',('nlat','nlon')) + +# --------------------------------------- +# Write out the data arrays into the file +# --------------------------------------- +grid_file.variables['lon'][:,:] = lon[:,:] +grid_file.variables['lat'][:,:] = lat[:,:] +grid_file.variables['clo'][:,:] = clo[:,:,:] +grid_file.variables['cla'][:,:] = cla[:,:,:] +grid_file.variables['srf'][:,:] = surface[:,:] +grid_file.variables['imask'][:,:] = mask[:,:] + +# --------------------------------------- +# close the file +# --------------------------------------- +grid_file.close() + +print (' Closing netcdf file : grid_toy_modelMARC.nc') +print ('------------------------------------------') +################################################## + +################################################## +print ('------------------------------------------') +print (' Creating netcdf file : rstrt_TOYMARC.nc') + +rstrt_file=netCDF4.Dataset(curdir_path+'rstrt_TOYMARC.nc','w',format='NETCDF3_64BIT') +rstrt_file.Description='Restart file for OASIS coupling' + +# ---------------------------------- +# Create the dimensions of the files +# ---------------------------------- +rstrt_file.createDimension ('nlon', nlon) +rstrt_file.createDimension ('nlat', nlat) + +# ---------------------------------- +# Create the variables of the files +# ---------------------------------- +varout=rstrt_file.createVariable('VARCNSTE','d',('nlat','nlon')) +varout=rstrt_file.createVariable('VARSIN01','d',('nlat','nlon')) +varout=rstrt_file.createVariable('VARSIN02','d',('nlat','nlon')) + +# --------------------------------------- +# Write out the data arrays into the file +# --------------------------------------- +rstrt_file.variables['VARCNSTE'][:,:] = toyvarcnste[:,:] +rstrt_file.variables['VARSIN01'][:,:] = toyvarsinus[:,:] +rstrt_file.variables['VARSIN02'][:,:] = toyvarsinus[:,:] + +# --------------------------------------- +# close the file +# --------------------------------------- +rstrt_file.close() + +print (' Closing netcdf file : rstrt_TOYMARC.nc ') +print ('-----------------------------------------') +##################################################### diff --git a/src/LIB/Python/Panel_Plot.py b/src/LIB/Python/Panel_Plot.py index 2072d1cf9225c27f14b79849e8b11e5e49885ce8..e4588ddde7852b734d0972de097edcbf41bf894f 100644 --- a/src/LIB/Python/Panel_Plot.py +++ b/src/LIB/Python/Panel_Plot.py @@ -9,7 +9,7 @@ MNH_LIC for details. version 1. @author: 07/2021 Quentin Rodier """ import matplotlib as mpl -#mpl.use('Agg') +# mpl.use('Agg') import matplotlib.pyplot as plt from matplotlib import cm from matplotlib.colors import ListedColormap @@ -17,9 +17,10 @@ import numpy as np import cartopy import cartopy.feature as cfeature + class PanelPlot(): - - def __init__(self, nb_l, nb_c, Lfigsize, bigtitle, titlepad=40, minmaxpad=1.03, timepad=-0.06, lateralminmaxpad=0.86, + + def __init__(self, nb_l, nb_c, Lfigsize, bigtitle, bigtitlepad=0.95, titlepad=40, minmaxpad=1.03, timepad=-0.06, lateralminmaxpad=0.86, labelcolorbarpad=6.0, colorbaraspect=20, colorbarpad=0.04, tickspad=0.8, minmaxTextSize=10, bigtitleSize=13, titleSize=12, legendSize=10, xlabelSize=11, ylabelSize=11, timeSize=11, cbTicksLabelSize=11, cbTitleSize=11, xyTicksLabelSize=10, figBoxLinewidth=1, @@ -31,6 +32,7 @@ class PanelPlot(): self.nb_c = nb_c # Panel number of rows self.nb_graph = 0 # New independent graph within the subplot + self.bigtitlepad = bigtitlepad # Panel title vertical position self.titlepad = titlepad # Title pad (vertical shift) from graph self.tickspad = tickspad # Ticks pad (between ticks and axis label) self.minmaxpad = minmaxpad # Min/Max print pad (vertical shift) @@ -56,873 +58,1007 @@ class PanelPlot(): self.xyTicksLength = xyTicksLength # Ticks length # Initialization of the panel plots - self.fig = plt.figure(figsize=(self.Lfigsize[0],self.Lfigsize[1])) + self.fig = plt.figure(figsize=(self.Lfigsize[0], self.Lfigsize[1])) self.fig.set_dpi(125) - self.fig.suptitle(self.bigtitle,fontsize=bigtitleSize) + self.fig.suptitle(self.bigtitle, fontsize=self.bigtitleSize, y=self.bigtitlepad) def save_graph(self, iplt, fig, fig_name='tempgraph'): - """ - Create a temporary png file of the panel plot which can be converted to PDF - """ - self.iplt = iplt - self.fig = fig - self.fig_name=fig_name # .png figure prefix name - - self.fig.savefig(self.fig_name+str(self.iplt)) - self.iplt+=1 - return self.iplt - - def draw_Backmap(self,drawCoastLines, ax, projo): - """ - Handle drawing of the background plot (coastlines, departements, grid lines and labels) - """ - self.drawCoastLines = drawCoastLines - self.projo = projo - - # Grid lines and labels - if 'PlateCarree' in str(projo): + """ + Create a temporary png file of the panel plot which can be converted to PDF + """ + self.iplt = iplt + self.fig = fig + self.fig_name = fig_name # .png figure prefix name + + self.fig.savefig(self.fig_name + str(self.iplt)) + self.iplt += 1 + return self.iplt + + def draw_Backmap(self, drawCoastLines, ax, projo): + """ + Handle drawing of the background plot (coastlines, departements, grid lines and labels) + """ + self.drawCoastLines = drawCoastLines + self.projo = projo + + # Grid lines and labels gl = ax.gridlines(crs=self.projo, draw_labels=True, linewidth=1, color='gray') if float(cartopy.__version__[:4]) >= 0.18: - gl.top_labels = False - gl.right_labels = False + from cartopy.mpl.ticker import (LatitudeLocator, LongitudeLocator,LongitudeFormatter, LatitudeFormatter) + gl.top_labels = False + gl.right_labels = False + gl.xlines = True + gl.ylines = True + gl.xlocator = LongitudeLocator() + gl.ylocator = LatitudeLocator() + gl.xformatter = LongitudeFormatter() + gl.yformatter = LatitudeFormatter() else: - gl.xlabels_top = False - gl.ylabels_right = False - - # Coastlines - if self.drawCoastLines and 'GeoAxes' in str(type(ax)): - ax.coastlines(resolution='10m') - - # Countries border - ax.add_feature(cfeature.BORDERS) - ax.add_feature(cfeature.LAKES, alpha=0.7) - - def addWhitecm(self, colormap_in, nb_level,whiteTop=False): + gl.xlabels_top = False + gl.ylabels_right = False + gl.xlabel_style = {'size': self.xyTicksLabelSize, 'color': 'black'} + gl.ylabel_style = {'size': self.xyTicksLabelSize, 'color': 'black'} + + # Coastlines + if self.drawCoastLines and 'GeoAxes' in str(type(ax)): + ax.coastlines(resolution='10m', color='black', linewidth=1) + + # Countries border + ax.add_feature(cfeature.BORDERS) + ax.add_feature(cfeature.LAKES, alpha=0.7) + + def addWhitecm(self, colormap_in, nb_level, whiteTop=False): """ Add a white color at the top (whiteTop=True) or bottom of the colormap w.r.t. the number of independent colors used """ color_map = cm.get_cmap(colormap_in, 256) newcolor_map = color_map(np.linspace(0, 1, 256)) - whites = np.array([1, 1, 1, 1]) #RBG code + opacity + whites = np.array([1, 1, 1, 1]) # RBG code + opacity if whiteTop: - for i in range(int(256/nb_level)): newcolor_map[-i, :] = whites + for i in range(int(256 / nb_level)): + newcolor_map[-i, :] = whites else: - for i in range(int(256/nb_level)): newcolor_map[:i, :] = whites + for i in range(int(256 / nb_level)): + newcolor_map[:i, :] = whites newcmp = ListedColormap(newcolor_map) return newcmp - def set_Title(self, ax, i, title, Lid_overlap, xlab, ylab): - """ - Handle top title of each graph - Parameters : - - titlepad : global attribute for vertical shift placement w.r.t the graph - """ - self.ax = ax - self.title = title - self.Lid_overlap = Lid_overlap - self.i = i - #self.ax[self.i].set_xlabel("test", fontweight='bold') - if not Lid_overlap: - self.ax[self.i].set_title(title, pad=self.titlepad, fontsize=self.titleSize) - else: # If graph overlap, title is concatenated - new_title = self.ax[self.i].get_title() + ' and ' + title - self.ax[self.i].set_title(new_title, pad=self.titlepad, fontsize=self.titleSize) + """ + Handle top title of each graph + Parameters : + - titlepad : global attribute for vertical shift placement w.r.t the graph + """ + self.ax = ax + self.title = title + self.Lid_overlap = Lid_overlap + self.i = i + #self.ax[self.i].set_xlabel("test", fontweight='bold') + if not Lid_overlap: + self.ax[self.i].set_title(title, pad=self.titlepad, fontsize=self.titleSize) + else: # If graph overlap, title is concatenated + new_title = self.ax[self.i].get_title() + ' and ' + title + self.ax[self.i].set_title(new_title, pad=self.titlepad, fontsize=self.titleSize) def set_xlim(self, ax, i, xlim): - """ - Handle x limits plotted if necessary - """ - self.ax = ax - self.xlim = xlim - self.i = i - - self.ax[self.i].set_xlim(xlim[0],xlim[1])#, fontsize=self.xlabelSize) - + """ + Handle x limits plotted if necessary + """ + self.ax = ax + self.xlim = xlim + self.i = i + + self.ax[self.i].set_xlim(xlim[0], xlim[1]) # , fontsize=self.xlabelSize) + def set_ylim(self, ax, i, ylim): - """ - Handle x limits plotted if necessary - """ - self.ax = ax - self.ylim = ylim - self.i = i - - self.ax[self.i].set_ylim(ylim[0],ylim[1])#, fontsize=self.ylabelSize) + """ + Handle x limits plotted if necessary + """ + self.ax = ax + self.ylim = ylim + self.i = i + + self.ax[self.i].set_ylim(ylim[0], ylim[1]) # , fontsize=self.ylabelSize) def set_XYaxislab(self, ax, i, xlab, ylab): - """ - Handle x and y axis labels - """ - self.ax = ax - self.xlab = xlab - self.ylab = ylab - self.i = i - - # This handing label is a known issue with GeoAxes of cartopy - # https://stackoverflow.com/questions/35479508/cartopy-set-xlabel-set-ylabel-not-ticklabels - # https://github.com/SciTools/cartopy/issues/1332 - if 'GeoAxes' in str(type(self.ax[self.i])): - self.ax[self.i].text(-0.11, 0.45, ylab, verticalalignment='top', horizontalalignment='left', - rotation='vertical', rotation_mode='anchor', transform=self.ax[self.i].transAxes, color='black', fontsize=self.ylabelSize) - self.ax[self.i].text(0.45, -0.06, xlab, verticalalignment='top', horizontalalignment='left', - rotation='horizontal', rotation_mode='anchor', transform=self.ax[self.i].transAxes, color='black', fontsize=self.xlabelSize) - else: - self.ax[self.i].set_xlabel(xlab, fontsize=self.xlabelSize, labelpad=0.1) - self.ax[self.i].set_ylabel(ylab, fontsize=self.ylabelSize, labelpad=0.1) - + """ + Handle x and y axis labels + """ + self.ax = ax + self.xlab = xlab + self.ylab = ylab + self.i = i + + # This handing label is a known issue with GeoAxes of cartopy + # https://stackoverflow.com/questions/35479508/cartopy-set-xlabel-set-ylabel-not-ticklabels + # https://github.com/SciTools/cartopy/issues/1332 + if 'GeoAxes' in str(type(self.ax[self.i])): + self.ax[self.i].text(-0.11, 0.45, ylab, verticalalignment='top', horizontalalignment='left', + rotation='vertical', rotation_mode='anchor', transform=self.ax[self.i].transAxes, color='black', fontsize=self.ylabelSize) + self.ax[self.i].text(0.45, -0.06, xlab, verticalalignment='top', horizontalalignment='left', + rotation='horizontal', rotation_mode='anchor', transform=self.ax[self.i].transAxes, color='black', fontsize=self.xlabelSize) + else: + self.ax[self.i].set_xlabel(xlab, fontsize=self.xlabelSize, labelpad=0.1) + self.ax[self.i].set_ylabel(ylab, fontsize=self.ylabelSize, labelpad=0.1) + def addLine(self, ax, beg_coord, end_coord, color='black', linewidth=0.2): - self.ax = ax - self.beg_coord = beg_coord - self.end_coord = end_coord - self.color = color - self.linewidth = linewidth - - x1, y1 = [self.beg_coord[0],self.end_coord[0]], [self.beg_coord[1], self.end_coord[1]] - ax.plot(x1, y1,color=self.color,linewidth=self.linewidth) + self.ax = ax + self.beg_coord = beg_coord + self.end_coord = end_coord + self.color = color + self.linewidth = linewidth + + x1, y1 = [self.beg_coord[0], self.end_coord[0]], [self.beg_coord[1], self.end_coord[1]] + ax.plot(x1, y1, color=self.color, linewidth=self.linewidth) def print_minmax(self, var, title): print(str(title) + " min = " + str(np.nanmin(var)) + " max = " + str(np.nanmax(var))) def set_minmaxText(self, ax, i, var, title, Lid_overlap, facconv): - """ - Show min and max value Text in the plot - If overlap variable, text is align to the right - TODO : handle more than 2 overlap variables - """ - self.ax = ax - self.var = var - self.i = i - self.title = title - self.facconv = facconv - - strtext = " min = " + "{:.3e}".format(np.nanmin(var*facconv)) + " max = " + "{:.3e}".format(np.nanmax(var*facconv)) - if not Lid_overlap: - self.ax[self.i].text(0.01, self.minmaxpad, strtext, verticalalignment='top', horizontalalignment='left', - transform=self.ax[self.i].transAxes, color='black', fontsize=self.minmaxTextSize) - else: - self.ax[self.i].text(self.lateralminmaxpad, self.minmaxpad, strtext, verticalalignment='top', horizontalalignment='right', - transform=self.ax[self.i].transAxes, color='black', fontsize=self.minmaxTextSize) - # Print to help choose min/max value for ticks - self.print_minmax(var*facconv, title) - + """ + Show min and max value Text in the plot + If overlap variable, text is align to the right + TODO : handle more than 2 overlap variables + """ + self.ax = ax + self.var = var + self.i = i + self.title = title + self.facconv = facconv + + strtext = " min = " + "{:.3e}".format(np.nanmin(var * facconv)) + " max = " + "{:.3e}".format(np.nanmax(var * facconv)) + if not Lid_overlap: + self.ax[self.i].text(0.01, self.minmaxpad, strtext, verticalalignment='top', horizontalalignment='left', + transform=self.ax[self.i].transAxes, color='black', fontsize=self.minmaxTextSize) + else: + self.ax[self.i].text(self.lateralminmaxpad, self.minmaxpad, strtext, verticalalignment='top', horizontalalignment='right', + transform=self.ax[self.i].transAxes, color='black', fontsize=self.minmaxTextSize) + # Print to help choose min/max value for ticks + self.print_minmax(var * facconv, title) + def showTimeText(self, ax, i, timetxt): - """ - Show time validity - """ - self.ax = ax - self.i = i - self.timetxt = timetxt - - strtext = "Time = " + timetxt - self.ax[self.i].text(0.0, self.timepad, strtext, verticalalignment='top', horizontalalignment='left', - transform=self.ax[self.i].transAxes, color='black', fontsize=self.timeSize) - - - def psectionV(self, Lxx=[], Lzz=[], Lvar=[], Lxlab=[], Lylab=[], Ltitle=[], Lminval=[], Lmaxval=[], - Lstep=[], Lstepticks=[], Lcolormap=[], Lcbarlabel=[], LcolorLine=[], - Lfacconv=[], ax=[], Lid_overlap=[], colorbar=True, orog=[], Lxlim=[], Lylim=[], Ltime=[], Lpltype=[], LaddWhite_cm=[], LwhiteTop=[]): - """ - Vertical cross section plot - Parameters : - - Lxx : List of x or y coordinate variable or time axis - - Lzz : List of z coordinates variable - - Lvar : List of variables to plot - - Lxlab : List of x-axis label - - Lylab : List of y-axis label - - Lxlim : List of x (min, max) value plotted - - Lylim : List of y (min, max) value plotted - - Ltime : List of time (validity) - - Ltitle : List of sub-title - - Lminval: List of minimum value for each colorbar - - Lmaxval: List of maximum value for each colorbar - - Lstep : List of color-steps for each colorbar - - Lstepticks : List of value of labels for each colorbar - - Lcolormap : List of colormap - - LcolorLine : List of colors for colors arg of contour (color line only) - - Lcbarlabel : List of colorbar label legend (units) - - Lfacconv : List of factors for unit conversion of each variables - - ax : List of fig.axes for ploting multiple different types of plots in a subplot panel - - Lid_overlap: List of number index of plot to overlap current variables - - Lpltype : List of types of plot 'cf' or 'c'. cf=contourf, c=contour (lines only) - - colorbar : show colorbar or not - - LaddWhite_cm : List of boolean to add white color to a colormap at the last bottom (low value) tick colorbar - - LwhiteTop : List of boolean to add the white color at the first top (high value). If false, the white is added at the bottom if Laddwhite_cm=T - - orog : Orography variable - """ - self.ax = ax - firstCall = (len(self.ax) == 0) - - # Initialize default value w.r.t to the number of plots -# D={'Lxlab':Lxlab, 'Lylab':Lylab, 'Ltitle':Ltitle,'Lminval':Lminval, 'Lmaxval':Lmaxval, + """ + Show time validity + """ + self.ax = ax + self.i = i + self.timetxt = timetxt + + strtext = "Time = " + timetxt + self.ax[self.i].text(0.0, self.timepad, strtext, verticalalignment='top', horizontalalignment='left', + transform=self.ax[self.i].transAxes, color='black', fontsize=self.timeSize) + + def psectionV(self, Lxx=[], Lzz=[], Lvar=[], Lxlab=[], Lylab=[], Ltitle=[], Lminval=[], Lmaxval=[], + Lstep=[], Lstepticks=[], Lcolormap=[], Lcbarlabel=[], LcolorLine=[], Lcbformatlabel=[], Llinewidth=[], + Lfacconv=[], ax=[], Lid_overlap=[], colorbar=True, orog=[], Lxlim=[], Lylim=[], Ltime=[], Lpltype=[], LaddWhite_cm=[], LwhiteTop=[]): + """ + Vertical cross section plot + Parameters : + - Lxx : List of x or y coordinate variable or time axis + - Lzz : List of z coordinates variable + - Lvar : List of variables to plot + - Lxlab : List of x-axis label + - Lylab : List of y-axis label + - Lxlim : List of x (min, max) value plotted + - Lylim : List of y (min, max) value plotted + - Ltime : List of time (validity) + - Ltitle : List of sub-title + - Lminval: List of minimum value for each colorbar + - Lmaxval: List of maximum value for each colorbar + - Lstep : List of color-steps for each colorbar + - Lstepticks : List of value of labels for each colorbar + - Lcolormap : List of colormap + - Llinewidth : List of lines thickness of contour + - LcolorLine : List of colors for colors arg of contour (color line only) + - Lcbarlabel : List of colorbar label legend (units) + - Lfacconv : List of factors for unit conversion of each variables + - ax : List of fig.axes for ploting multiple different types of plots in a subplot panel + - Lid_overlap: List of number index of plot to overlap current variables + - Lpltype : List of types of plot 'cf' or 'c'. cf=contourf, c=contour (lines only) + - colorbar : show colorbar or not + - LaddWhite_cm : List of boolean to add white color to a colormap at the last bottom (low value) tick colorbar + - LwhiteTop : List of boolean to add the white color at the first top (high value). If false, the white is added at the bottom if Laddwhite_cm=T + - Lcbformatlabel: List of boolean to reduce the format to exponential 1.1E+02 format colorbar label + - orog : Orography variable + """ + self.ax = ax + firstCall = (len(self.ax) == 0) + + # Initialize default value w.r.t to the number of plots +# D={'Lxlab':Lxlab, 'Lylab':Lylab, 'Ltitle':Ltitle,'Lminval':Lminval, 'Lmaxval':Lmaxval, # 'Lstep':Lstep, 'Lstepticks':Lstepticks, 'Lcolormap':Lcolormap, 'Lcbarlabel':Lcbarlabel, 'Lfacconv':Lfacconv, 'Ltime':Ltime, # 'LaddWhite_cm':LaddWhite_cm, 'Lpltype':Lpltype} # D = initialize_default_val(Lvar, D) - # Default values - if not Lfacconv: Lfacconv = [1.0]*len(Lvar) - if not Lcolormap and not LcolorLine: Lcolormap=['gist_rainbow_r']*len(Lvar) # If no color given, a cmap is given - if not Lcolormap: LcolorLine=['black']*len(Lvar) - if not Lpltype: Lpltype=['cf']*len(Lvar) - if not LaddWhite_cm: LaddWhite_cm=[False]*len(Lvar) - if not LwhiteTop: LwhiteTop=[False]*len(Lvar) - if not Lylab: Lylab = ['']*len(Lvar) - if not Lxlab: Lxlab = ['']*len(Lvar) - # Add an extra percentage of the top max value for forcing the colorbar show the true user maximum value (correct a bug) - Lmaxval = list(map(lambda x, y: x + 1E-6*y, Lmaxval, Lstep) ) #The extra value is 1E-6 times the step ticks of the colorbar - - # On all variables to plot - for i,var in enumerate(Lvar): - if firstCall: #1st call - iax = i - self.ax.append(self.fig.add_subplot(self.nb_l,self.nb_c,i+1)) - self.nb_graph+=1 - elif Lid_overlap != []: #overlapping plot - iax = Lid_overlap[i] - else: #existing ax with no overlapping (graphd appended to existing panel) - self.ax.append(self.fig.add_subplot(self.nb_l,self.nb_c,self.nb_graph+1)) - self.nb_graph+=1 - iax = len(self.ax)-1 # The ax index of the new coming plot is the length of the existant ax -1 for indices matter - - # Colors normalization - norm = mpl.colors.Normalize(vmax=Lmaxval[i], vmin=Lminval[i]) - - # Print min/max (stout and on plot) - self.set_minmaxText(self.ax, iax, var, Ltitle[i], Lid_overlap, Lfacconv[i]) - - # Print time validity - if Ltime: self.showTimeText(self.ax, iax, str(Ltime[i])) - - # Number of contours level - if not Lstep[i]: # Default value of number of steps is 20 - Lstep[i] = (Lmaxval[i] - Lminval[i])/20 - Lstepticks[i] = Lstep[i] - - levels_contour = np.arange(Lminval[i],Lmaxval[i],step=Lstep[i]) - - # Add White to colormap - if LaddWhite_cm[i] and Lcolormap: Lcolormap[i]=self.addWhitecm(Lcolormap[i], len(levels_contour), LwhiteTop[i]) - - # Plot - if Lpltype[i]=='c': # Contour - if LcolorLine: - cf = self.ax[iax].contour(Lxx[i], Lzz[i], var*Lfacconv[i], levels=levels_contour, - norm=norm, vmin=Lminval[i], vmax=Lmaxval[i], colors=LcolorLine[i], linewidths=0.1) - else: - cf = self.ax[iax].contour(Lxx[i], Lzz[i], var*Lfacconv[i], levels=levels_contour, - norm=norm, vmin=Lminval[i], vmax=Lmaxval[i], cmap=Lcolormap[i], linewidths=0.1) - else: # Contourf - cf = self.ax[iax].contourf(Lxx[i], Lzz[i], var*Lfacconv[i], levels=levels_contour, - norm=norm, vmin=Lminval[i], vmax=Lmaxval[i], cmap=Lcolormap[i]) - - # Title - self.set_Title(self.ax, iax, Ltitle[i], Lid_overlap,Lxlab[i], Lylab[i]) - - # X/Y Axis label - self.set_XYaxislab(self.ax, iax, Lxlab[i], Lylab[i]) - - # Ticks label - self.ax[iax].tick_params(axis='both', labelsize=self.xyTicksLabelSize, width=self.xyTicksWidth, length=self.xyTicksLength, pad=self.tickspad) - - # Bounding box of the plot line width - for axis in ['top','bottom','left','right']: - self.ax[iax].spines[axis].set_linewidth(self.figBoxLinewidth) - - # X/Y Axis limits value - if Lxlim: - try: - self.set_xlim(self.ax, iax, Lxlim[i]) - except: - pass - if Lylim: - try: - self.set_ylim(self.ax, iax, Lylim[i]) - except: - pass - - # Color label on contour-line - if Lpltype[i]=='c': # Contour - self.ax[iax].clabel(cf, fontsize=self.cbTicksLabelSize) - #self.ax[iax].clabel(cf, levels=np.arange(Lminval[i],Lmaxval[i],step=Lstep[i]), fontsize=self.cbTicksLabelSize) #TODO bug, levels not recognized - - #Filling area under topography - if not orog==[]: - self.ax[iax].fill_between(Lxx[i][0,:], orog, color='black', linewidth=0.2) - - # Colorbar - if colorbar: - cb=plt.colorbar(cf, ax=self.ax[iax], fraction=0.031, pad=self.colorbarpad, ticks=np.arange(Lminval[i],Lmaxval[i],Lstepticks[i]), aspect=self.colorbaraspect) - cb.ax.set_title(Lcbarlabel[i], pad=self.labelcolorbarpad, loc='left', fontsize=self.cbTitleSize) # This creates a new AxesSubplot only for the colorbar y=0 ==> location at the bottom - - return self.fig + # Default values + if not Lfacconv: + Lfacconv = [1.0] * len(Lvar) + if not Lcolormap and not LcolorLine: + Lcolormap = ['gist_rainbow_r'] * len(Lvar) # If no color given, a cmap is given + if not Lcolormap: + LcolorLine = ['black'] * len(Lvar) + if not Lpltype: + Lpltype = ['cf'] * len(Lvar) + if not LaddWhite_cm: + LaddWhite_cm = [False] * len(Lvar) + if not LwhiteTop: + LwhiteTop = [False] * len(Lvar) + if not Lylab: + Lylab = [''] * len(Lvar) + if not Lxlab: + Lxlab = [''] * len(Lvar) + if not Lcbformatlabel: + Lcbformatlabel = [False] * len(Lvar) + if not Llinewidth: + Llinewidth = [1.0] * len(Lvar) + + # Add an extra percentage of the top max value for forcing the colorbar show the true user maximum value (correct a bug) + Lmaxval = list(map(lambda x, y: x + 1E-6 * y, Lmaxval, Lstep)) # The extra value is 1E-6 times the step ticks of the colorbar + + # On all variables to plot + for i, var in enumerate(Lvar): + if firstCall: # 1st call + iax = i + self.ax.append(self.fig.add_subplot(self.nb_l, self.nb_c, i + 1)) + self.nb_graph += 1 + elif Lid_overlap != []: # overlapping plot + iax = Lid_overlap[i] + else: # existing ax with no overlapping (graphd appended to existing panel) + self.ax.append(self.fig.add_subplot(self.nb_l, self.nb_c, self.nb_graph + 1)) + self.nb_graph += 1 + iax = len(self.ax) - 1 # The ax index of the new coming plot is the length of the existant ax -1 for indices matter + + # Colors normalization + norm = mpl.colors.Normalize(vmax=Lmaxval[i], vmin=Lminval[i]) + + # Print min/max (stout and on plot) + self.set_minmaxText(self.ax, iax, var, Ltitle[i], Lid_overlap, Lfacconv[i]) + + # Print time validity + if Ltime: + self.showTimeText(self.ax, iax, str(Ltime[i])) + + # Number of contours level + if not Lstep[i]: # Default value of number of steps is 20 + Lstep[i] = (Lmaxval[i] - Lminval[i]) / 20 + Lstepticks[i] = Lstep[i] + + levels_contour = np.arange(Lminval[i], Lmaxval[i], step=Lstep[i]) + + # Add White to colormap + if LaddWhite_cm[i] and Lcolormap: + Lcolormap[i] = self.addWhitecm(Lcolormap[i], len(levels_contour), LwhiteTop[i]) + + # Plot + if Lpltype[i] == 'c': # Contour + if LcolorLine: + cf = self.ax[iax].contour(Lxx[i], Lzz[i], var * Lfacconv[i], levels=levels_contour, + norm=norm, vmin=Lminval[i], vmax=Lmaxval[i], colors=LcolorLine[i], linewidths=Llinewidth[i]) + else: + cf = self.ax[iax].contour(Lxx[i], Lzz[i], var * Lfacconv[i], levels=levels_contour, + norm=norm, vmin=Lminval[i], vmax=Lmaxval[i], cmap=Lcolormap[i], linewidths=Llinewidth[i]) + else: # Contourf + cf = self.ax[iax].contourf(Lxx[i], Lzz[i], var * Lfacconv[i], levels=levels_contour, + norm=norm, vmin=Lminval[i], vmax=Lmaxval[i], cmap=Lcolormap[i]) + + # Title + self.set_Title(self.ax, iax, Ltitle[i], Lid_overlap, Lxlab[i], Lylab[i]) + + # X/Y Axis label + self.set_XYaxislab(self.ax, iax, Lxlab[i], Lylab[i]) + + # Ticks label + self.ax[iax].tick_params( + axis='both', + labelsize=self.xyTicksLabelSize, + width=self.xyTicksWidth, + length=self.xyTicksLength, + pad=self.tickspad) + + # Bounding box of the plot line width + for axis in ['top', 'bottom', 'left', 'right']: + self.ax[iax].spines[axis].set_linewidth(self.figBoxLinewidth) + + # X/Y Axis limits value + if Lxlim: + try: + self.set_xlim(self.ax, iax, Lxlim[i]) + except BaseException: + pass + if Lylim: + try: + self.set_ylim(self.ax, iax, Lylim[i]) + except BaseException: + pass + + # Color label on contour-line + if Lpltype[i] == 'c': # Contour + self.ax[iax].clabel(cf, fontsize=self.cbTicksLabelSize) + # self.ax[iax].clabel(cf, + # levels=np.arange(Lminval[i],Lmaxval[i],step=Lstep[i]), + # fontsize=self.cbTicksLabelSize) #TODO bug, levels not recognized + + # Filling area under topography + if not orog == []: + if Lxx[i].ndim == 1: + self.ax[iax].fill_between(Lxx[i], orog, color='black', linewidth=0.2) + else: + self.ax[iax].fill_between(Lxx[i][0, :], orog, color='black', linewidth=0.2) + + # Colorbar + if colorbar: + cb = plt.colorbar( + cf, + ax=self.ax[iax], + fraction=0.031, + pad=self.colorbarpad, + ticks=np.arange( + Lminval[i], + Lmaxval[i], + Lstepticks[i]), + aspect=self.colorbaraspect) + cb.ax.tick_params(labelsize=self.cbTicksLabelSize) + # This creates a new AxesSubplot only for the colorbar y=0 ==> location at the bottom + cb.ax.set_title(Lcbarlabel[i], pad=self.labelcolorbarpad, loc='left', fontsize=self.cbTitleSize) + if Lcbformatlabel[i]: + cb.ax.set_yticklabels(["{:.1E}".format(i) for i in cb.get_ticks()]) + + return self.fig def pXY_lines(self, Lxx=[], Lyy=[], Lxlab=[], Lylab=[], Ltitle=[], Llinetype=[], Llinewidth=[], - Llinecolor=[], Llinelabel=[], LfacconvX=[], LfacconvY=[], ax=[], id_overlap=None, Lxlim=[], - Lylim=[], Ltime=[], LaxisColor=[]): - """ - XY (multiple)-lines plot - Parameters : - - Lxx : List of variables to plot or coordinates along the X axis - - Lyy : List of variables to plot or coordinates along the Y axis - - Lxlab : List of x-axis label - - Lylab : List of y-axis label - - Lxlim : List of x (min, max) value plotted - - Lylim : List of y (min, max) value plotted - - Ltime : List of time (validity) - - Ltitle : List of sub-title - - Llinewidth : List of lines thickness - - Llinetype : List of line types - - Lcolorlines: List of color lines - - Llinelabel : List of legend label lines - - LfacconvX/Y: List of factors for unit conversion of the variables/coordinates to plot on X and Y axis - - ax : List of fig.axes for ploting multiple different types of plots in a subplot panel - - Lid_overlap: List of number index of plot to overlap current variables - - LaxisColor : List of colors for multiple x-axis overlap - """ - self.ax = ax - firstCall = (len(self.ax) == 0) - # Defaults value convert to x number of variables list - if not LfacconvX: LfacconvX = [1.0]*len(Lxx) - if not LfacconvY: LfacconvY = [1.0]*len(Lxx) - if not Llinewidth: Llinewidth = [1.0]*len(Lxx) - if not Llinecolor: Llinecolor = ['blue']*len(Lxx) - if not Llinetype: Llinetype = ['-']*len(Lxx) - if not Llinelabel: Llinelabel = ['']*len(Lxx) - if not LaxisColor: LaxisColor = ['black']*len(Lxx) - if not Lylab: Lylab = ['']*len(Lxx) - if not Lxlab: Lxlab = ['']*len(Lxx) - - if firstCall: # 1st call - iax = 0 - self.ax.append(self.fig.add_subplot(self.nb_l,self.nb_c,1, label='graph axe x down')) - self.nb_graph+=1 - elif id_overlap: # overlapping plot with a different x-axis - self.ax.append(self.fig.add_subplot(self.nb_l,self.nb_c,self.nb_graph, label='graph axe x top',frame_on=False)) - iax = len(self.ax)-1 - else: # existing ax with no overlapping (graph appended to existing panel) - self.ax.append(self.fig.add_subplot(self.nb_l,self.nb_c,self.nb_graph+1)) - self.nb_graph+=1 - iax = len(self.ax)-1 # The ax index of the new coming plot is the length of the existant ax -1 for indices matter - - # On all variables to plot - for i,var in enumerate(Lxx): - # Print time validity - if Ltime: self.showTimeText(self.ax, iax, str(Ltime[i])) - - # Plot - cf = self.ax[iax].plot(Lxx[i]*LfacconvX[i], Lyy[i]*LfacconvY[i], color=Llinecolor[i], ls=Llinetype[i], - label=Llinelabel[i], linewidth=Llinewidth[i]) - # Legend - #TODO : Handling legend with overlap two axis lines in the same box. For now, placement is by hand - if not id_overlap: - self.ax[iax].legend(loc='upper right', bbox_to_anchor=(1, 0.95),fontsize=self.legendSize) - else: - self.ax[iax].legend(loc='upper right', bbox_to_anchor=(1, 0.90),fontsize=self.legendSize) + Llinecolor=[], Llinelabel=[], LfacconvX=[], LfacconvY=[], ax=[], id_overlap=None, Lxlim=[], + Lylim=[], Ltime=[], LaxisColor=[], LlocLegend=[]): + """ + XY (multiple)-lines plot + Parameters : + - Lxx : List of variables to plot or coordinates along the X axis + - Lyy : List of variables to plot or coordinates along the Y axis + - Lxlab : List of x-axis label + - Lylab : List of y-axis label + - Lxlim : List of x (min, max) value plotted + - Lylim : List of y (min, max) value plotted + - Ltime : List of time (validity) + - Ltitle : List of sub-title + - Llinewidth : List of lines thickness + - Llinetype : List of line types + - Lcolorlines: List of color lines + - Llinelabel : List of legend label lines + - LfacconvX/Y: List of factors for unit conversion of the variables/coordinates to plot on X and Y axis + - ax : List of fig.axes for ploting multiple different types of plots in a subplot panel + - Lid_overlap: List of number index of plot to overlap current variables + - LaxisColor : List of colors for multiple x-axis overlap + - LlocLegend : List of localisation of the legend : 'best', 'upper left', 'upper right', 'lower left', 'lower right', + 'upper center', 'lower center', 'center left', 'center right', 'center' + """ + self.ax = ax + firstCall = (len(self.ax) == 0) + # Defaults value convert to x number of variables list + if not LfacconvX: + LfacconvX = [1.0] * len(Lxx) + if not LfacconvY: + LfacconvY = [1.0] * len(Lxx) + if not Llinewidth: + Llinewidth = [1.0] * len(Lxx) + if not Llinecolor: + Llinecolor = ['blue'] * len(Lxx) + if not Llinetype: + Llinetype = ['-'] * len(Lxx) + if not Llinelabel: + Llinelabel = [''] * len(Lxx) + if not LaxisColor: + LaxisColor = ['black'] * len(Lxx) + if not Lylab: + Lylab = [''] * len(Lxx) + if not Lxlab: + Lxlab = [''] * len(Lxx) + if not LlocLegend: + LlocLegend = ['upper right'] * len(Lxx) - # Title - if Ltitle: self.set_Title(self.ax, iax, Ltitle[i], id_overlap,Lxlab[i], Lylab[i]) - - # Ticks label - self.ax[iax].tick_params(axis='both', labelsize=self.xyTicksLabelSize, width=self.xyTicksWidth, length=self.xyTicksLength, pad=self.tickspad) - - # X/Y Axis label - if id_overlap: - self.ax[iax].xaxis.tick_top() - self.ax[iax].xaxis.set_label_position('top') - self.ax[iax].set_xlabel(Lxlab[i]) - self.ax[iax].set_ylabel(Lylab[i]) - else: - self.set_XYaxislab(self.ax, iax, Lxlab[i], Lylab[i]) - - self.ax[iax].tick_params(axis='x', colors=LaxisColor[i]) - - # X/Y Axis limits value - if Lxlim: - try: - self.set_xlim(self.ax, iax, Lxlim[i]) - except: - pass - if Lylim: - try: - self.set_ylim(self.ax, iax, Lylim[i]) - except: - pass - return self.fig - - def psectionH(self, lon=[],lat=[], Lvar=[], Lcarte=[], Llevel=[], Lxlab=[], Lylab=[], Ltitle=[], Lminval=[], Lmaxval=[], - Lstep=[], Lstepticks=[], Lcolormap=[], LcolorLine=[], Lcbarlabel=[], Lproj=[], Lfacconv=[], coastLines=True, ax=[], - Lid_overlap=[], colorbar=True, Ltime=[], LaddWhite_cm=[], LwhiteTop=[], Lpltype=[], Lcbformatlabel=[]): - """ - Horizontal cross section plot - Parameters : - - lon : longitude 2D array - - lat : latitude 2D array - - Lvar : List of variables to plot - - Lcarte : Zooming [lonmin, lonmax, latmin, latmax] - - Llevel : List of k-level value for the section plot (ignored if variable is already 2D) - - Lxlab : List of x-axis label - - Lylab : List of y-axis label - - Ltitle : List of sub-title - - Ltime : List of time (validity) - - Lminval: List of minimum value for each colorbar - - Lmaxval: List of maximum value for each colorbar - - Lstep : List of color-steps for each colorbar - - Lstepticks : List of value of labels for each colorbar - - Lcolormap : List of colormap - - LcolorLine : List of colors for colors arg of contour (color line only) - - Lcbarlabel : List of colorbar label legend (units) - - Lproj : List of ccrs cartopy projection ([] for cartesian coordinates) - - Lfacconv : List of factors for unit conversion of each variables - - coastLines : Boolean to plot coast lines and grid lines - - ax : List of fig.axes for ploting multiple different types of plots in a subplot panel - - Lid_overlap : List of number index of plot to overlap current variables - - colorbar : show colorbar or not - - Lpltype : List of types of plot 'cf' or 'c'. cf=contourf, c=contour (lines only) - - LaddWhite_cm : List of boolean to add white color to a colormap at the first (low value) tick colorbar - - LwhiteTop : List of boolean to add the white color at the first top (high value). If false, the white is added at the bottom if Laddwhite_cm=T - - Lcbformatlabel: List of boolean to reduce the format to exponential 1.1E+02 format colorbar label - """ - self.ax = ax - firstCall = (len(self.ax) == 0) - - # Initialize default value w.r.t to the number of plots -# D={'lon':lon, 'lat':lat, 'Lcarte':Lcarte, 'Llevel':Llevel, 'Lxlab':Lxlab, 'Lylab':Lylab, 'Ltitle':Ltitle,'Lminval':Lminval, 'Lmaxval':Lmaxval, + if firstCall: # 1st call + iax = 0 + self.ax.append(self.fig.add_subplot(self.nb_l, self.nb_c, 1, label='graph axe x down')) + self.nb_graph += 1 + elif id_overlap: # overlapping plot with a different x-axis + self.ax.append(self.fig.add_subplot(self.nb_l, self.nb_c, self.nb_graph, label='graph axe x top', frame_on=False)) + iax = len(self.ax) - 1 + else: # existing ax with no overlapping (graph appended to existing panel) + self.ax.append(self.fig.add_subplot(self.nb_l, self.nb_c, self.nb_graph + 1)) + self.nb_graph += 1 + iax = len(self.ax) - 1 # The ax index of the new coming plot is the length of the existant ax -1 for indices matter + + # On all variables to plot + for i, var in enumerate(Lxx): + # Print time validity + if Ltime: + self.showTimeText(self.ax, iax, str(Ltime[i])) + + # Plot + cf = self.ax[iax].plot(Lxx[i] * LfacconvX[i], Lyy[i] * LfacconvY[i], color=Llinecolor[i], ls=Llinetype[i], + label=Llinelabel[i], linewidth=Llinewidth[i]) + # Legend + # TODO : Handling legend with overlap two axis lines in the same box. For now, placement is by hand + if not id_overlap: + self.ax[iax].legend(loc=LlocLegend[i], bbox_to_anchor=(1, 0.95), fontsize=self.legendSize) + else: + self.ax[iax].legend(loc=LlocLegend[i], bbox_to_anchor=(1, 0.90), fontsize=self.legendSize) + + # Title + if Ltitle: + self.set_Title(self.ax, iax, Ltitle[i], id_overlap, Lxlab[i], Lylab[i]) + + # Ticks label + self.ax[iax].tick_params( + axis='both', + labelsize=self.xyTicksLabelSize, + width=self.xyTicksWidth, + length=self.xyTicksLength, + pad=self.tickspad) + + # X/Y Axis label + if id_overlap: + self.ax[iax].xaxis.tick_top() + self.ax[iax].xaxis.set_label_position('top') + self.ax[iax].set_xlabel(Lxlab[i]) + self.ax[iax].set_ylabel(Lylab[i]) + else: + self.set_XYaxislab(self.ax, iax, Lxlab[i], Lylab[i]) + + self.ax[iax].tick_params(axis='x', colors=LaxisColor[i]) + + # X/Y Axis limits value + if Lxlim: + try: + self.set_xlim(self.ax, iax, Lxlim[i]) + except BaseException: + pass + if Lylim: + try: + self.set_ylim(self.ax, iax, Lylim[i]) + except BaseException: + pass + return self.fig + + def psectionH(self, lon=[], lat=[], Lvar=[], Lcarte=[], Llevel=[], Lxlab=[], Lylab=[], Ltitle=[], Lminval=[], Lmaxval=[], + Lstep=[], Lstepticks=[], Lcolormap=[], LcolorLine=[], Lcbarlabel=[], Lproj=[], Lfacconv=[], coastLines=True, ax=[], + Lid_overlap=[], colorbar=True, Ltime=[], LaddWhite_cm=[], LwhiteTop=[], Lpltype=[], Lcbformatlabel=[], Llinewidth=[]): + """ + Horizontal cross section plot + Parameters : + - lon : longitude 2D array + - lat : latitude 2D array + - Lvar : List of variables to plot + - Lcarte : Zooming [lonmin, lonmax, latmin, latmax] + - Llevel : List of k-level value for the section plot (ignored if variable is already 2D) + - Lxlab : List of x-axis label + - Lylab : List of y-axis label + - Ltitle : List of sub-title + - Ltime : List of time (validity) + - Lminval: List of minimum value for each colorbar + - Lmaxval: List of maximum value for each colorbar + - Lstep : List of color-steps for each colorbar + - Lstepticks : List of value of labels for each colorbar + - Llinewidth : List of lines thickness of contour + - Lcolormap : List of colormap + - LcolorLine : List of colors for colors arg of contour (color line only) + - Lcbarlabel : List of colorbar label legend (units) + - Lproj : List of ccrs cartopy projection ([] for cartesian coordinates) + - Lfacconv : List of factors for unit conversion of each variables + - coastLines : Boolean to plot coast lines and grid lines + - ax : List of fig.axes for ploting multiple different types of plots in a subplot panel + - Lid_overlap : List of number index of plot to overlap current variables + - colorbar : show colorbar or not + - Lpltype : List of types of plot 'cf' or 'c'. cf=contourf, c=contour (lines only) + - LaddWhite_cm : List of boolean to add white color to a colormap at the first (low value) tick colorbar + - LwhiteTop : List of boolean to add the white color at the first top (high value). If false, the white is added at the bottom if Laddwhite_cm=T + - Lcbformatlabel: List of boolean to reduce the format to exponential 1.1E+02 format colorbar label + """ + self.ax = ax + firstCall = (len(self.ax) == 0) + + # Initialize default value w.r.t to the number of plots +# D={'lon':lon, 'lat':lat, 'Lcarte':Lcarte, 'Llevel':Llevel, 'Lxlab':Lxlab, 'Lylab':Lylab, 'Ltitle':Ltitle,'Lminval':Lminval, 'Lmaxval':Lmaxval, # 'Lstep':Lstep, 'Lstepticks':Lstepticks, 'Lcolormap':Lcolormap, 'Lcbarlabel':Lcbarlabel, 'Lproj':Lproj, 'Lfacconv':Lfacconv, 'Ltime':Ltime, # 'LaddWhite_cm':LaddWhite_cm, 'Lpltype':Lpltype} # D = initialize_default_val(Lvar, D) - - # If all plots are not using conversion factor, convert it to List - if not Lfacconv: Lfacconv = [1.0]*len(Lvar) - if not Lylab: Lylab = ['']*len(Lvar) - if not Lxlab: Lxlab = ['']*len(Lvar) - if not Lcolormap and not LcolorLine: Lcolormap=['gist_rainbow_r']*len(Lvar) # If no color given, a cmap is given - if not Lcolormap: LcolorLine=['black']*len(Lvar) - if not Lpltype: Lpltype=['cf']*len(Lvar) - if not LaddWhite_cm: LaddWhite_cm=[False]*len(Lvar) - if not LwhiteTop: LwhiteTop=[False]*len(Lvar) - if not Lcbformatlabel: Lcbformatlabel=[False]*len(Lvar) - # Add an extra percentage of the top max value for forcing the colorbar show the true user maximum value (correct a bug) - if Lstep: Lmaxval = list(map(lambda x, y: x + 1E-6*y, Lmaxval, Lstep) ) #The extra value is 1E-6 times the step ticks of the colorbar - - # This following must be after the extra percentage top value addition - if not Lstep: Lstep=[None]*len(Lvar) - if not Lstepticks: Lstepticks=[None]*len(Lvar) - - # On all variables to plot - for i,var in enumerate(Lvar): - if firstCall: #1st call - iax = i - if Lproj: - self.ax.append(self.fig.add_subplot(self.nb_l,self.nb_c,i+1, projection = Lproj[i])) - else: # Cartesian coordinates - self.ax.append(self.fig.add_subplot(self.nb_l,self.nb_c,i+1)) - self.nb_graph+=1 - elif Lid_overlap != []: #overlapping plot - iax = Lid_overlap[i] - else: #existing ax with no overlapping (graphd appended to existing panel) - if Lproj: - self.ax.append(self.fig.add_subplot(self.nb_l,self.nb_c,self.nb_graph+1, projection = Lproj[i])) - else: # Cartesian coordinates - self.ax.append(self.fig.add_subplot(self.nb_l,self.nb_c,self.nb_graph+1)) - self.nb_graph+=1 - iax = len(self.ax)-1 # The ax index of the new coming plot is the length of the existant ax -1 for indices matter - - #Colors normalization - norm = mpl.colors.Normalize(vmax=Lmaxval[i], vmin=Lminval[i]) - - # Zooming - if len(Lcarte) == 4: #zoom - self.ax[iax].set_xlim(Lcarte[0], Lcarte[1]) - self.ax[iax].set_ylim(Lcarte[2], Lcarte[3]) - # Variable to plot w.r.t dimensions - if var.ndim==2: - vartoPlot = var[:,:] - else: - vartoPlot = var[Llevel[i],:,:] - - # Print min/max (stout and on plot) - self.set_minmaxText(self.ax, iax, vartoPlot, Ltitle[i], Lid_overlap, Lfacconv[i]) - - # Print time validity - if Ltime: self.showTimeText(self.ax, iax, str(Ltime[i])) - - # Number of contours level - if not Lstep[i]: # Default value of number of steps is 20 - Lstep[i] = (Lmaxval[i] - Lminval[i])/20 - Lstepticks[i] = Lstep[i] - - levels_contour = np.arange(Lminval[i],Lmaxval[i],step=Lstep[i]) - - # Add White to colormap - if LaddWhite_cm[i] and Lcolormap: Lcolormap[i]=self.addWhitecm(Lcolormap[i], len(levels_contour), LwhiteTop[i]) - - # Plot - if Lproj: - if Lpltype[i]=='c': # Contour - if LcolorLine: - cf = self.ax[iax].contour(lon[i], lat[i], vartoPlot*Lfacconv[i], levels=levels_contour,transform=Lproj[i], - norm=norm, vmin=Lminval[i], vmax=Lmaxval[i], colors=LcolorLine[i]) - else: - cf = self.ax[iax].contour(lon[i], lat[i], vartoPlot*Lfacconv[i], levels=levels_contour,transform=Lproj[i], - norm=norm, vmin=Lminval[i], vmax=Lmaxval[i], cmap=Lcolormap[i]) - else: - cf = self.ax[iax].contourf(lon[i], lat[i], vartoPlot*Lfacconv[i], levels=levels_contour,transform=Lproj[i], - norm=norm, vmin=Lminval[i], vmax=Lmaxval[i], cmap=Lcolormap[i]) - else: # Cartesian coordinates - if Lpltype[i]=='c': # Contour - cf = self.ax[iax].contour(lon[i], lat[i], vartoPlot*Lfacconv[i], levels=levels_contour, - norm=norm, vmin=Lminval[i], vmax=Lmaxval[i], colors=LcolorLine[i]) - else: - cf = self.ax[iax].contourf(lon[i], lat[i], vartoPlot*Lfacconv[i], levels=levels_contour, - norm=norm, vmin=Lminval[i], vmax=Lmaxval[i], cmap=Lcolormap[i]) - # Title - self.set_Title(self.ax, iax, Ltitle[i], Lid_overlap,Lxlab[i], Lylab[i]) - - # Coastlines / Grid lines and labels - if Lproj: self.draw_Backmap(coastLines, self.ax[iax], Lproj[i]) - - # X/Y Axis - self.set_XYaxislab(self.ax, iax, Lxlab[i], Lylab[i]) - # Ticks label - self.ax[iax].tick_params(axis='both', labelsize=self.xyTicksLabelSize, width=self.xyTicksWidth, length=self.xyTicksLength, pad=self.tickspad) - - # Color label on contour-line - if Lpltype[i]=='c': # Contour - if 'GeoAxes' in str(type(self.ax[self.i])): # cartopy does not like the levels arguments in clabel, known issue - self.ax[iax].clabel(cf, fontsize=self.cbTicksLabelSize) - else: - self.ax[iax].clabel(cf, levels=np.arange(Lminval[i],Lmaxval[i],step=Lstep[i]), fontsize=self.cbTicksLabelSize) - - # Colorbar - if colorbar: - cb=plt.colorbar(cf, ax=self.ax[iax], fraction=0.031, pad=self.colorbarpad, ticks=np.arange(Lminval[i],Lmaxval[i],Lstepticks[i]), aspect=self.colorbaraspect) - cb.ax.set_title(Lcbarlabel[i], pad = self.labelcolorbarpad, loc='left', fontsize=self.cbTitleSize) #This creates a new AxesSubplot only for the colorbar y=0 ==> location at the bottom - if Lcbformatlabel[i]: cb.ax.set_yticklabels(["{:.1E}".format(i) for i in cb.get_ticks()]) - - return self.fig - - def pvector(self, Lxx=[], Lyy=[], Lvar1=[], Lvar2=[], Lcarte=[], Llevel=[], Lxlab=[], Lylab=[], - Ltitle=[], Lwidth=[], Larrowstep=[], Lcolor=[], Llegendval=[], Llegendlabel=[], + + # If all plots are not using conversion factor, convert it to List + if not Lfacconv: + Lfacconv = [1.0] * len(Lvar) + if not Lylab: + Lylab = [''] * len(Lvar) + if not Lxlab: + Lxlab = [''] * len(Lvar) + if not Lcolormap and not LcolorLine: + Lcolormap = ['gist_rainbow_r'] * len(Lvar) # If no color given, a cmap is given + if not Lcolormap: + LcolorLine = ['black'] * len(Lvar) + if not Lpltype: + Lpltype = ['cf'] * len(Lvar) + if not LaddWhite_cm: + LaddWhite_cm = [False] * len(Lvar) + if not LwhiteTop: + LwhiteTop = [False] * len(Lvar) + if not Lcbformatlabel: + Lcbformatlabel = [False] * len(Lvar) + if not Llinewidth: + Llinewidth = [1.0] * len(Lvar) + + # Add an extra percentage of the top max value for forcing the colorbar show the true user maximum value (correct a bug) + if Lstep: + Lmaxval = list(map(lambda x, y: x + 1E-6 * y, Lmaxval, Lstep)) # The extra value is 1E-6 times the step ticks of the colorbar + + # This following must be after the extra percentage top value addition + if not Lstep: + Lstep = [None] * len(Lvar) + if not Lstepticks: + Lstepticks = [None] * len(Lvar) + + # On all variables to plot + for i, var in enumerate(Lvar): + if firstCall: # 1st call + iax = i + if Lproj: + self.ax.append(self.fig.add_subplot(self.nb_l, self.nb_c, i + 1, projection=Lproj[i])) + else: # Cartesian coordinates + self.ax.append(self.fig.add_subplot(self.nb_l, self.nb_c, i + 1)) + self.nb_graph += 1 + elif Lid_overlap != []: # overlapping plot + iax = Lid_overlap[i] + else: # existing ax with no overlapping (graphd appended to existing panel) + if Lproj: + self.ax.append(self.fig.add_subplot(self.nb_l, self.nb_c, self.nb_graph + 1, projection=Lproj[i])) + else: # Cartesian coordinates + self.ax.append(self.fig.add_subplot(self.nb_l, self.nb_c, self.nb_graph + 1)) + self.nb_graph += 1 + iax = len(self.ax) - 1 # The ax index of the new coming plot is the length of the existant ax -1 for indices matter + + # Colors normalization + norm = mpl.colors.Normalize(vmax=Lmaxval[i], vmin=Lminval[i]) + + # Zooming + if len(Lcarte) == 4: # zoom + self.ax[iax].set_xlim(Lcarte[0], Lcarte[1]) + self.ax[iax].set_ylim(Lcarte[2], Lcarte[3]) + # Variable to plot w.r.t dimensions + if var.ndim == 2: + vartoPlot = var[:, :] + else: + vartoPlot = var[Llevel[i], :, :] + + # Print min/max (stout and on plot) + self.set_minmaxText(self.ax, iax, vartoPlot, Ltitle[i], Lid_overlap, Lfacconv[i]) + + # Print time validity + if Ltime: + self.showTimeText(self.ax, iax, str(Ltime[i])) + + # Number of contours level + if not Lstep[i]: # Default value of number of steps is 20 + Lstep[i] = (Lmaxval[i] - Lminval[i]) / 20 + Lstepticks[i] = Lstep[i] + + levels_contour = np.arange(Lminval[i], Lmaxval[i], step=Lstep[i]) + + # Add White to colormap + if LaddWhite_cm[i] and Lcolormap: + Lcolormap[i] = self.addWhitecm(Lcolormap[i], len(levels_contour), LwhiteTop[i]) + + # Plot + if Lproj: + if Lpltype[i] == 'c': # Contour + if LcolorLine: + cf = self.ax[iax].contour(lon[i], lat[i], vartoPlot * Lfacconv[i], levels=levels_contour, transform=Lproj[i], + norm=norm, vmin=Lminval[i], vmax=Lmaxval[i], colors=LcolorLine[i],linewidths=Llinewidth[i]) + else: + cf = self.ax[iax].contour(lon[i], lat[i], vartoPlot * Lfacconv[i], levels=levels_contour, transform=Lproj[i], + norm=norm, vmin=Lminval[i], vmax=Lmaxval[i], cmap=Lcolormap[i],linewidths=Llinewidth[i]) + else: + cf = self.ax[iax].contourf(lon[i], lat[i], vartoPlot * Lfacconv[i], levels=levels_contour, transform=Lproj[i], + norm=norm, vmin=Lminval[i], vmax=Lmaxval[i], cmap=Lcolormap[i]) + else: # Cartesian coordinates + if Lpltype[i] == 'c': # Contour + cf = self.ax[iax].contour(lon[i], lat[i], vartoPlot * Lfacconv[i], levels=levels_contour, + norm=norm, vmin=Lminval[i], vmax=Lmaxval[i], colors=LcolorLine[i],linewidths=Llinewidth[i]) + else: + cf = self.ax[iax].contourf(lon[i], lat[i], vartoPlot * Lfacconv[i], levels=levels_contour, + norm=norm, vmin=Lminval[i], vmax=Lmaxval[i], cmap=Lcolormap[i]) + # Title + self.set_Title(self.ax, iax, Ltitle[i], Lid_overlap, Lxlab[i], Lylab[i]) + + # Coastlines / Grid lines and labels + if Lproj: + self.draw_Backmap(coastLines, self.ax[iax], Lproj[i]) + + # X/Y Axis + self.set_XYaxislab(self.ax, iax, Lxlab[i], Lylab[i]) + # Ticks label + self.ax[iax].tick_params( + axis='both', + labelsize=self.xyTicksLabelSize, + width=self.xyTicksWidth, + length=self.xyTicksLength, + pad=self.tickspad) + + # Color label on contour-line + if Lpltype[i] == 'c': # Contour + if 'GeoAxes' in str(type(self.ax[self.i])): # cartopy does not like the levels arguments in clabel, known issue + self.ax[iax].clabel(cf, fontsize=self.cbTicksLabelSize) + else: + self.ax[iax].clabel(cf, levels=np.arange(Lminval[i], Lmaxval[i], step=Lstep[i]), fontsize=self.cbTicksLabelSize) + + # Colorbar + if colorbar: + cb = plt.colorbar( + cf, + ax=self.ax[iax], + fraction=0.031, + pad=self.colorbarpad, + ticks=np.arange( + Lminval[i], + Lmaxval[i], + Lstepticks[i]), + aspect=self.colorbaraspect) + cb.ax.tick_params(labelsize=self.cbTicksLabelSize) + # This creates a new AxesSubplot only for the colorbar y=0 ==> location at the bottom + cb.ax.set_title(Lcbarlabel[i], pad=self.labelcolorbarpad, loc='left', fontsize=self.cbTitleSize) + if Lcbformatlabel[i]: + cb.ax.set_yticklabels(["{:.1E}".format(i) for i in cb.get_ticks()]) + + return self.fig + + def pvector(self, Lxx=[], Lyy=[], Lvar1=[], Lvar2=[], Lcarte=[], Llevel=[], Lxlab=[], Lylab=[], + Ltitle=[], Lwidth=[], Larrowstep=[], Lcolor=[], Llegendval=[], Llegendlabel=[], Lproj=[], Lfacconv=[], ax=[], coastLines=True, Lid_overlap=[], Ltime=[], Lscale=[], Lylim=[], Lxlim=[]): - """ - Vectors - Parameters : - - Lxx : List of x or y coordinate variable (lat or ni or nm) - - Lyy : List of y coordinates variable (lon or level) - - Lvar1 : List of wind-component along x/y or oblic axis (3D for hor. section, 2D for vertical section) - - Lvar2 : List of wind-component along y-axis : v-component for horizontal section / w-component for vertical section - - Lcarte : Zooming [lonmin, lonmax, latmin, latmax] - - Llevel : List of k-level value for the horizontal section plot (ignored if variable is already 2D) - - Lxlab : List of x-axis label - - Lylab : List of y-axis label - - Lxlim : List of x (min, max) value plotted - - Lylim : List of y (min, max) value plotted - - Ltitle : List of sub-titles - - Ltime : List of time (validity) - - Lwidth : List of thickness of the arrows - - Lscale : List of scale for the length of the arrows (high value <=> small length) - - Larrowstep : List of sub-sample (frequency) if too much arrows - - Lcolor : List of colors for the arrows (default: black) - - Llegendval : List of value for the legend of the default arrow - - Llegendlabel : List of labels for the legend of the default arrow - - Lproj : List of ccrs cartopy projection - - Lfacconv : List of factors for unit conversion of each variables - - coastLines : Boolean to plot coast lines and grid lines - - ax : List of fig.axes for ploting multiple different types of plots in a subplot panel - - Lid_overlap : List of number index of plot to overlap current variables - """ - self.ax = ax - firstCall = (len(self.ax) == 0) - # If all plots are not using conversion factor, convert it to List - if not Lfacconv: Lfacconv = [1.0]*len(Lvar1) - if not Lcolor: Lcolor = ['black']*len(Lvar1) - if not Lscale : Lscale = [None]*len(Lvar1) - if not Lylab: Lylab = ['']*len(Lvar1) - if not Lxlab: Lxlab = ['']*len(Lvar1) - # On all variables to plot - for i,var1 in enumerate(Lvar1): - if firstCall: #1st call - iax = i - if Lproj: - self.ax.append(self.fig.add_subplot(self.nb_l,self.nb_c,i+1, projection = Lproj[i])) - else: - self.ax.append(self.fig.add_subplot(self.nb_l,self.nb_c,i+1)) - self.nb_graph+=1 - elif Lid_overlap != []: #overlapping plot - iax = Lid_overlap[i] - else: #existing ax with no overlapping (graphd appended to existing panel) - if Lproj: - self.ax.append(self.fig.add_subplot(self.nb_l,self.nb_c,self.nb_graph+1, projection = Lproj[i])) - else: - self.ax.append(self.fig.add_subplot(self.nb_l,self.nb_c,self.nb_graph+1)) - self.nb_graph+=1 - iax = len(self.ax)-1 # The ax index of the new coming plot is the length of the existant ax -1 for indices matter - - # Zooming - if len(Lcarte) == 4: #zoom - self.ax[iax].set_xlim(Lcarte[0], Lcarte[1]) - self.ax[iax].set_ylim(Lcarte[2], Lcarte[3]) - - # Variable to plot w.r.t dimensions - if var1.ndim==2: - vartoPlot1 = var1[:,:] - vartoPlot2 = Lvar2[i][:,:] - else: # Variable is 3D : only for horizontal section - vartoPlot1 = var1[Llevel[i],:,:] - vartoPlot2 = Lvar2[i][Llevel[i],:,:] - - # Print min/max val to help choose colorbar steps - self.set_minmaxText(self.ax, iax, np.sqrt(vartoPlot1**2 + vartoPlot2**2), Ltitle[i], Lid_overlap, Lfacconv[i]) - - # Print time validity - if Ltime: self.showTimeText(self.ax, iax, str(Ltime[i])) - - # Plot - axeX = Lxx[i] - axeY = Lyy[i] - if Lxx[i].ndim == 2: - cf = self.ax[iax].quiver(axeX[::Larrowstep[i],::Larrowstep[i]], axeY[::Larrowstep[i],::Larrowstep[i]], - vartoPlot1[::Larrowstep[i],::Larrowstep[i]], vartoPlot2[::Larrowstep[i],::Larrowstep[i]], - width=Lwidth[i] ,angles='uv', color=Lcolor[i],scale=Lscale[i]) - else: - cf = self.ax[iax].quiver(axeX[::Larrowstep[i]], axeY[::Larrowstep[i]], - vartoPlot1[::Larrowstep[i],::Larrowstep[i]], vartoPlot2[::Larrowstep[i],::Larrowstep[i]], - width=Lwidth[i] ,angles='uv', color=Lcolor[i], scale=Lscale[i]) - # Title - self.set_Title(self.ax, iax, Ltitle[i], Lid_overlap,Lxlab[i], Lylab[i]) - - # X/Y Axis Label - self.set_XYaxislab(self.ax, iax, Lxlab[i], Lylab[i]) - - # X/Y Axis limits value - if Lxlim: - try: - self.set_xlim(self.ax, iax, Lxlim[i]) - except: - pass - if Lylim: - try: - self.set_ylim(self.ax, iax, Lylim[i]) - except: - pass - - # Coastlines: - if Lproj: self.draw_Backmap(coastLines, self.ax[iax], Lproj[i]) - - # Arrow legend key - qk = self.ax[iax].quiverkey(cf, 1.0, -0.05, Llegendval[i], str(Llegendval[i]) + Llegendlabel[i], labelpos='E', color='black') - - return self.fig + """ + Vectors + Parameters : + - Lxx : List of x or y coordinate variable (lat or ni or nm) + - Lyy : List of y coordinates variable (lon or level) + - Lvar1 : List of wind-component along x/y or oblic axis (3D for hor. section, 2D for vertical section) + - Lvar2 : List of wind-component along y-axis : v-component for horizontal section / w-component for vertical section + - Lcarte : Zooming [lonmin, lonmax, latmin, latmax] + - Llevel : List of k-level value for the horizontal section plot (ignored if variable is already 2D) + - Lxlab : List of x-axis label + - Lylab : List of y-axis label + - Lxlim : List of x (min, max) value plotted + - Lylim : List of y (min, max) value plotted + - Ltitle : List of sub-titles + - Ltime : List of time (validity) + - Lwidth : List of thickness of the arrows + - Lscale : List of scale for the length of the arrows (high value <=> small length) + - Larrowstep : List of sub-sample (frequency) if too much arrows + - Lcolor : List of colors for the arrows (default: black) + - Llegendval : List of value for the legend of the default arrow + - Llegendlabel : List of labels for the legend of the default arrow + - Lproj : List of ccrs cartopy projection + - Lfacconv : List of factors for unit conversion of each variables + - coastLines : Boolean to plot coast lines and grid lines + - ax : List of fig.axes for ploting multiple different types of plots in a subplot panel + - Lid_overlap : List of number index of plot to overlap current variables + """ + self.ax = ax + firstCall = (len(self.ax) == 0) + # If all plots are not using conversion factor, convert it to List + if not Lfacconv: + Lfacconv = [1.0] * len(Lvar1) + if not Lcolor: + Lcolor = ['black'] * len(Lvar1) + if not Lscale: + Lscale = [None] * len(Lvar1) + if not Lylab: + Lylab = [''] * len(Lvar1) + if not Lxlab: + Lxlab = [''] * len(Lvar1) + # On all variables to plot + for i, var1 in enumerate(Lvar1): + if firstCall: # 1st call + iax = i + if Lproj: + self.ax.append(self.fig.add_subplot(self.nb_l, self.nb_c, i + 1, projection=Lproj[i])) + else: + self.ax.append(self.fig.add_subplot(self.nb_l, self.nb_c, i + 1)) + self.nb_graph += 1 + elif Lid_overlap != []: # overlapping plot + iax = Lid_overlap[i] + else: # existing ax with no overlapping (graphd appended to existing panel) + if Lproj: + self.ax.append(self.fig.add_subplot(self.nb_l, self.nb_c, self.nb_graph + 1, projection=Lproj[i])) + else: + self.ax.append(self.fig.add_subplot(self.nb_l, self.nb_c, self.nb_graph + 1)) + self.nb_graph += 1 + iax = len(self.ax) - 1 # The ax index of the new coming plot is the length of the existant ax -1 for indices matter + + # Zooming + if len(Lcarte) == 4: # zoom + self.ax[iax].set_xlim(Lcarte[0], Lcarte[1]) + self.ax[iax].set_ylim(Lcarte[2], Lcarte[3]) + + # Variable to plot w.r.t dimensions + if var1.ndim == 2: + vartoPlot1 = var1[:, :] + vartoPlot2 = Lvar2[i][:, :] + else: # Variable is 3D : only for horizontal section + vartoPlot1 = var1[Llevel[i], :, :] + vartoPlot2 = Lvar2[i][Llevel[i], :, :] + + # Print min/max val to help choose colorbar steps + self.set_minmaxText(self.ax, iax, np.sqrt(vartoPlot1**2 + vartoPlot2**2), Ltitle[i], Lid_overlap, Lfacconv[i]) + + # Print time validity + if Ltime: + self.showTimeText(self.ax, iax, str(Ltime[i])) + + # Plot + axeX = Lxx[i] + axeY = Lyy[i] + if Lxx[i].ndim == 2: + cf = self.ax[iax].quiver(axeX[::Larrowstep[i], ::Larrowstep[i]], axeY[::Larrowstep[i], ::Larrowstep[i]], + vartoPlot1[::Larrowstep[i], ::Larrowstep[i]], vartoPlot2[::Larrowstep[i], ::Larrowstep[i]], + width=Lwidth[i], angles='uv', color=Lcolor[i], scale=Lscale[i]) + else: + cf = self.ax[iax].quiver(axeX[::Larrowstep[i]], axeY[::Larrowstep[i]], + vartoPlot1[::Larrowstep[i], ::Larrowstep[i]], vartoPlot2[::Larrowstep[i], ::Larrowstep[i]], + width=Lwidth[i], angles='uv', color=Lcolor[i], scale=Lscale[i]) + # Title + self.set_Title(self.ax, iax, Ltitle[i], Lid_overlap, Lxlab[i], Lylab[i]) + + # X/Y Axis Label + self.set_XYaxislab(self.ax, iax, Lxlab[i], Lylab[i]) + + # X/Y Axis limits value + if Lxlim: + try: + self.set_xlim(self.ax, iax, Lxlim[i]) + except BaseException: + pass + if Lylim: + try: + self.set_ylim(self.ax, iax, Lylim[i]) + except BaseException: + pass + + # Coastlines: + if Lproj: + self.draw_Backmap(coastLines, self.ax[iax], Lproj[i]) + + # Arrow legend key + qk = self.ax[iax].quiverkey(cf, 1.0, -0.05, Llegendval[i], str(Llegendval[i]) + Llegendlabel[i], + labelpos='E', color='black', fontproperties={'size': self.legendSize}) + + return self.fig def pstreamline(self, Lxx=[], Lyy=[], Lvar1=[], Lvar2=[], Lcarte=[], Llevel=[], Lxlab=[], Lylab=[], Llinewidth=[], Ldensity=[], - Ltitle=[], Lcolor=[], Lproj=[], Lfacconv=[], ax=[], coastLines=True, Lid_overlap=[], Ltime=[], - Lylim=[], Lxlim=[]): - """ - Wind stream lines - Parameters : - - Lxx : List of x or y coordinate variable (lat or ni or nm) - - Lyy : List of y coordinates variable (lon or level) - - Lvar1 : List of wind-component along x/y or oblic axis (3D for hor. section, 2D for vertical section) - - Lvar2 : List of wind-component along y-axis : v-component for horizontal section / w-component for vertical section - - Lcarte : Zooming [lonmin, lonmax, latmin, latmax] - - Llevel : List of k-level value for the horizontal section plot (ignored if variable is already 2D) - - Lxlab : List of x-axis label - - Lylab : List of y-axis label - - Lxlim : List of x (min, max) value plotted - - Lylim : List of y (min, max) value plotted - - Ltitle : List of sub-titles - - Ltime : List of time (validity) - - Llinewidth : List of lines thickness - - Ldensity : List of density that control the closeness of streamlines - - Lcolor : List of colors for the streamline (default: black) - - Lproj : List of ccrs cartopy projection - - Lfacconv : List of factors for unit conversion of each variables - - coastLines : Boolean to plot coast lines and grid lines - - ax : List of fig.axes for ploting multiple different types of plots in a subplot panel - - Lid_overlap : List of number index of plot to overlap current variables - """ - self.ax = ax - firstCall = (len(self.ax) == 0) - # If all plots are not using conversion factor, convert it to List - if not Lfacconv: Lfacconv = [1.0]*len(Lvar1) - if not Lcolor: Lcolor = ['black']*len(Lvar1) - if not Lylab: Lylab = ['']*len(Lvar1) - if not Lxlab: Lxlab = ['']*len(Lvar1) - if not Llinewidth: Llinewidth = [1.0]*len(Lvar1) - if not Ldensity: Ldensity = [1.0]*len(Lvar1) - - # On all variables to plot - for i,var1 in enumerate(Lvar1): - if firstCall: #1st call - iax = i - if Lproj: - self.ax.append(self.fig.add_subplot(self.nb_l,self.nb_c,i+1, projection = Lproj[i])) - else: - self.ax.append(self.fig.add_subplot(self.nb_l,self.nb_c,i+1)) - self.nb_graph+=1 - elif Lid_overlap != []: #overlapping plot - iax = Lid_overlap[i] - else: #existing ax with no overlapping (graphd appended to existing panel) - if Lproj: - self.ax.append(self.fig.add_subplot(self.nb_l,self.nb_c,self.nb_graph+1, projection = Lproj[i])) - else: - self.ax.append(self.fig.add_subplot(self.nb_l,self.nb_c,self.nb_graph+1)) - self.nb_graph+=1 - iax = len(self.ax)-1 # The ax index of the new coming plot is the length of the existant ax -1 for indices matter - - # Zooming - if len(Lcarte) == 4: #zoom - self.ax[iax].set_xlim(Lcarte[0], Lcarte[1]) - self.ax[iax].set_ylim(Lcarte[2], Lcarte[3]) - - # Variable to plot w.r.t dimensions - if var1.ndim==2: - vartoPlot1 = var1[:,:] - vartoPlot2 = Lvar2[i][:,:] - else: # Variable is 3D : only for horizontal section - vartoPlot1 = var1[Llevel[i],:,:] - vartoPlot2 = Lvar2[i][Llevel[i],:,:] - - # Print min/max val to help choose steps - self.set_minmaxText(self.ax, iax, np.sqrt(vartoPlot1**2 + vartoPlot2**2), Ltitle[i], Lid_overlap, Lfacconv[i]) - - # Print time validity - if Ltime: self.showTimeText(self.ax, iax, str(Ltime[i])) + Ltitle=[], Lcolor=[], Lproj=[], Lfacconv=[], ax=[], coastLines=True, Lid_overlap=[], Ltime=[], + Lylim=[], Lxlim=[]): + """ + Wind stream lines + Parameters : + - Lxx : List of x or y coordinate variable (lat or ni or nm) + - Lyy : List of y coordinates variable (lon or level) + - Lvar1 : List of wind-component along x/y or oblic axis (3D for hor. section, 2D for vertical section) + - Lvar2 : List of wind-component along y-axis : v-component for horizontal section / w-component for vertical section + - Lcarte : Zooming [lonmin, lonmax, latmin, latmax] + - Llevel : List of k-level value for the horizontal section plot (ignored if variable is already 2D) + - Lxlab : List of x-axis label + - Lylab : List of y-axis label + - Lxlim : List of x (min, max) value plotted + - Lylim : List of y (min, max) value plotted + - Ltitle : List of sub-titles + - Ltime : List of time (validity) + - Llinewidth : List of lines thickness + - Ldensity : List of density that control the closeness of streamlines + - Lcolor : List of colors for the streamline (default: black) + - Lproj : List of ccrs cartopy projection + - Lfacconv : List of factors for unit conversion of each variables + - coastLines : Boolean to plot coast lines and grid lines + - ax : List of fig.axes for ploting multiple different types of plots in a subplot panel + - Lid_overlap : List of number index of plot to overlap current variables + """ + self.ax = ax + firstCall = (len(self.ax) == 0) + # If all plots are not using conversion factor, convert it to List + if not Lfacconv: + Lfacconv = [1.0] * len(Lvar1) + if not Lcolor: + Lcolor = ['black'] * len(Lvar1) + if not Lylab: + Lylab = [''] * len(Lvar1) + if not Lxlab: + Lxlab = [''] * len(Lvar1) + if not Llinewidth: + Llinewidth = [1.0] * len(Lvar1) + if not Ldensity: + Ldensity = [1.0] * len(Lvar1) - # Plot - cf = self.ax[iax].streamplot(Lxx[i], Lyy[i], vartoPlot1, vartoPlot2, density=Ldensity[i], linewidth=Llinewidth[i], color=Lcolor[i]) - - # Title - self.set_Title(self.ax, iax, Ltitle[i], Lid_overlap,Lxlab[i], Lylab[i]) + # On all variables to plot + for i, var1 in enumerate(Lvar1): + if firstCall: # 1st call + iax = i + if Lproj: + self.ax.append(self.fig.add_subplot(self.nb_l, self.nb_c, i + 1, projection=Lproj[i])) + else: + self.ax.append(self.fig.add_subplot(self.nb_l, self.nb_c, i + 1)) + self.nb_graph += 1 + elif Lid_overlap != []: # overlapping plot + iax = Lid_overlap[i] + else: # existing ax with no overlapping (graphd appended to existing panel) + if Lproj: + self.ax.append(self.fig.add_subplot(self.nb_l, self.nb_c, self.nb_graph + 1, projection=Lproj[i])) + else: + self.ax.append(self.fig.add_subplot(self.nb_l, self.nb_c, self.nb_graph + 1)) + self.nb_graph += 1 + iax = len(self.ax) - 1 # The ax index of the new coming plot is the length of the existant ax -1 for indices matter + + # Zooming + if len(Lcarte) == 4: # zoom + self.ax[iax].set_xlim(Lcarte[0], Lcarte[1]) + self.ax[iax].set_ylim(Lcarte[2], Lcarte[3]) + + # Variable to plot w.r.t dimensions + if var1.ndim == 2: + vartoPlot1 = var1[:, :] + vartoPlot2 = Lvar2[i][:, :] + else: # Variable is 3D : only for horizontal section + vartoPlot1 = var1[Llevel[i], :, :] + vartoPlot2 = Lvar2[i][Llevel[i], :, :] + + # Print min/max val to help choose steps + self.set_minmaxText(self.ax, iax, np.sqrt(vartoPlot1**2 + vartoPlot2**2), Ltitle[i], Lid_overlap, Lfacconv[i]) + + # Print time validity + if Ltime: + self.showTimeText(self.ax, iax, str(Ltime[i])) + + # Plot + cf = self.ax[iax].streamplot(Lxx[i], Lyy[i], vartoPlot1, vartoPlot2, density=Ldensity[i], linewidth=Llinewidth[i], color=Lcolor[i]) + + # Title + self.set_Title(self.ax, iax, Ltitle[i], Lid_overlap, Lxlab[i], Lylab[i]) + + # X/Y Axis Label + self.set_XYaxislab(self.ax, iax, Lxlab[i], Lylab[i]) + + # X/Y Axis limits value + if Lxlim: + try: + self.set_xlim(self.ax, iax, Lxlim[i]) + except BaseException: + pass + if Lylim: + try: + self.set_ylim(self.ax, iax, Lylim[i]) + except BaseException: + pass + + # Coastlines: + if Lproj: + self.draw_Backmap(coastLines, self.ax[iax], Lproj[i]) + + return self.fig - # X/Y Axis Label - self.set_XYaxislab(self.ax, iax, Lxlab[i], Lylab[i]) - - # X/Y Axis limits value - if Lxlim: - try: - self.set_xlim(self.ax, iax, Lxlim[i]) - except: - pass - if Lylim: - try: - self.set_ylim(self.ax, iax, Lylim[i]) - except: - pass - - # Coastlines: - if Lproj: self.draw_Backmap(coastLines, self.ax[iax], Lproj[i]) - - return self.fig - def pXY_bar(self, Lbins=[], Lvar=[], Lxlab=[], Lylab=[], Ltitle=[], Lcolor=[], Lwidth=[], - Llinecolor=[], Llinewidth=[], Lfacconv=[], ax=[], id_overlap=None, Lxlim=[], - Lylim=[], Ltime=[], LaxisColor=[]): - """ - XY Histogram - Parameters : - - Lbins : List of bins - - Lvar : List of the value for each bin - - Lxlab : List of x-axis label - - Lylab : List of y-axis label - - Ltitle : List of sub-title - - Lcolor : List of color (or sequence of colors for each value) of the bars - - Lwidth : List of width of the bars - - Llinecolor : List of line color of the bar edges - - Llinewidth : List of lines thickness of the bar edges - - Lfacconv : List of factors for unit conversion of each variables - - Lxlim : List of x (min, max) value plotted - - Lylim : List of y (min, max) value plotted - - Ltime : List of time (validity) - - ax : List of fig.axes for ploting multiple different types of plots in a subplot panel - - Lid_overlap: List of number index of plot to overlap current variables - - LaxisColor : List of colors for multiple x-axis overlap - """ - self.ax = ax - firstCall = (len(self.ax) == 0) - # Defaults value convert to x number of variables list - if not Lfacconv: Lfacconv = [1.0]*len(Lvar) - if not LaxisColor: LaxisColor = ['black']*len(Lvar) - if not Lylab: Lylab = ['']*len(Lvar) - if not Lxlab: Lxlab = ['']*len(Lvar) - if not Lcolor: Lcolor = ['black']*len(Lvar) - if not Lwidth: Lwidth = [1]*len(Lvar) - if not Llinecolor: Llinecolor = ['black']*len(Lvar) - if not Llinewidth: Llinewidth = [0]*len(Lvar) - - # On all variables to plot - for i,var in enumerate(Lvar): - if firstCall: # 1st call - iax = i - self.ax.append(self.fig.add_subplot(self.nb_l,self.nb_c,i+1, label='graph axe x down')) - self.nb_graph+=1 - elif id_overlap: # overlapping plot with a different x-axis - self.ax.append(self.fig.add_subplot(self.nb_l,self.nb_c,self.nb_graph, label='graph axe x top',frame_on=False)) - iax = len(self.ax)-1 - else: # existing ax with no overlapping (graph appended to existing panel) - self.ax.append(self.fig.add_subplot(self.nb_l,self.nb_c,self.nb_graph+1)) - self.nb_graph+=1 - print(len(self.ax)) - iax = len(self.ax)-1 # The ax index of the new coming plot is the length of the existant ax -1 for indices matter - - # Print time validity - if Ltime: self.showTimeText(self.ax, iax, str(Ltime[i])) - - # Bins by labels - labels = np.array([str(L) for L in Lbins[i]]) - - # Plot - cf = self.ax[iax].bar(labels, var*Lfacconv[i], width=Lwidth[i], color=Lcolor[i], linewidth=Llinewidth[i], edgecolor=Llinecolor[i]) - - # Legend - #TODO : Handling legend with overlap two axis lines in the same box. For now, placement is by hand - if not id_overlap: - self.ax[iax].legend(loc='upper right', bbox_to_anchor=(1, 0.95)) - else: - self.ax[iax].legend(loc='upper right', bbox_to_anchor=(1, 0.90)) - - # Title - if Ltitle: self.set_Title(self.ax, iax, Ltitle[i], id_overlap,Lxlab[i], Lylab[i]) - - # X/Y Axis label - if id_overlap: - self.ax[iax].xaxis.tick_top() - self.ax[iax].xaxis.set_label_position('top') - self.ax[iax].set_xlabel(Lxlab[i]) - self.ax[iax].set_ylabel(Lylab[i]) - else: - self.set_XYaxislab(self.ax, iax, Lxlab[i], Lylab[i]) - - self.ax[iax].tick_params(axis='x', colors=LaxisColor[i]) - - # X/Y Axis limits value - if Lxlim: - try: - self.set_xlim(self.ax, iax, Lxlim[i]) - except: - pass - if Lylim: - try: - self.set_ylim(self.ax, iax, Lylim[i]) - except: - pass - - return self.fig -#def initialize_default_val(Lvar, Dparam): + Llinecolor=[], Llinewidth=[], Lfacconv=[], ax=[], id_overlap=None, Lxlim=[], + Lylim=[], Ltime=[], LaxisColor=[], LlocLegend=[]): + """ + XY Histogram + Parameters : + - Lbins : List of bins + - Lvar : List of the value for each bin + - Lxlab : List of x-axis label + - Lylab : List of y-axis label + - Ltitle : List of sub-title + - Lcolor : List of color (or sequence of colors for each value) of the bars + - Lwidth : List of width of the bars + - Llinecolor : List of line color of the bar edges + - Llinewidth : List of lines thickness of the bar edges + - Lfacconv : List of factors for unit conversion of each variables + - Lxlim : List of x (min, max) value plotted + - Lylim : List of y (min, max) value plotted + - Ltime : List of time (validity) + - ax : List of fig.axes for ploting multiple different types of plots in a subplot panel + - Lid_overlap: List of number index of plot to overlap current variables + - LaxisColor : List of colors for multiple x-axis overlap + - LlocLegend : List of localisation of the legend : 'best', 'upper left', 'upper right', 'lower left', 'lower right', + 'upper center', 'lower center', 'center left', 'center right', 'center' + """ + self.ax = ax + firstCall = (len(self.ax) == 0) + # Defaults value convert to x number of variables list + if not Lfacconv: + Lfacconv = [1.0] * len(Lvar) + if not LaxisColor: + LaxisColor = ['black'] * len(Lvar) + if not Lylab: + Lylab = [''] * len(Lvar) + if not Lxlab: + Lxlab = [''] * len(Lvar) + if not Lcolor: + Lcolor = ['black'] * len(Lvar) + if not Lwidth: + Lwidth = [1] * len(Lvar) + if not Llinecolor: + Llinecolor = ['black'] * len(Lvar) + if not Llinewidth: + Llinewidth = [0] * len(Lvar) + if not LlocLegend: + LlocLegend = ['upper right'] * len(Lvar) + + # On all variables to plot + for i, var in enumerate(Lvar): + if firstCall: # 1st call + iax = i + self.ax.append(self.fig.add_subplot(self.nb_l, self.nb_c, i + 1, label='graph axe x down')) + self.nb_graph += 1 + elif id_overlap: # overlapping plot with a different x-axis + self.ax.append(self.fig.add_subplot(self.nb_l, self.nb_c, self.nb_graph, label='graph axe x top', frame_on=False)) + iax = len(self.ax) - 1 + else: # existing ax with no overlapping (graph appended to existing panel) + self.ax.append(self.fig.add_subplot(self.nb_l, self.nb_c, self.nb_graph + 1)) + self.nb_graph += 1 + print(len(self.ax)) + iax = len(self.ax) - 1 # The ax index of the new coming plot is the length of the existant ax -1 for indices matter + + # Print time validity + if Ltime: + self.showTimeText(self.ax, iax, str(Ltime[i])) + + # Bins by labels + labels = np.array([str(L) for L in Lbins[i]]) + + # Plot + cf = self.ax[iax].bar(labels, var * Lfacconv[i], width=Lwidth[i], color=Lcolor[i], linewidth=Llinewidth[i], edgecolor=Llinecolor[i]) + + # Legend + # TODO : Handling legend with overlap two axis lines in the same box. For now, placement is by hand + if not id_overlap: + self.ax[iax].legend(loc=LlocLegend[i], bbox_to_anchor=(1, 0.95)) + else: + self.ax[iax].legend(loc=LlocLegend[i], bbox_to_anchor=(1, 0.90)) + + # Title + if Ltitle: + self.set_Title(self.ax, iax, Ltitle[i], id_overlap, Lxlab[i], Lylab[i]) + + # X/Y Axis label + if id_overlap: + self.ax[iax].xaxis.tick_top() + self.ax[iax].xaxis.set_label_position('top') + self.ax[iax].set_xlabel(Lxlab[i]) + self.ax[iax].set_ylabel(Lylab[i]) + else: + self.set_XYaxislab(self.ax, iax, Lxlab[i], Lylab[i]) + + self.ax[iax].tick_params(axis='x', colors=LaxisColor[i]) + + # X/Y Axis limits value + if Lxlim: + try: + self.set_xlim(self.ax, iax, Lxlim[i]) + except BaseException: + pass + if Lylim: + try: + self.set_ylim(self.ax, iax, Lylim[i]) + except BaseException: + pass + + return self.fig +# def initialize_default_val(Lvar, Dparam): # #TODO : initialize default value for all parameter of all type of graphs # #Returns All the parameters given in Dparam where : # "- If no value is found (empty list []) : return the default value (if exist) * number of graph # #- If ONE value only is found : return the value copied x times the number of graph # #- If the list is complete, nothing is done # #CURRENT PROBLEM -# #The returned value do not change the true referenced variable given as argument +# #The returned value do not change the true referenced variable given as argument # # Number of graphs # l = len(Lvar) -# +# # # Add an extra percentage of the top max value for forcing the colorbar show the true user maximum value (correct a bug) # #if Dparam['Lstep'] and Dparam['Lmaxval']: Dparam['Lmaxval'] = list(map(lambda x, y: x + 1E-6*y, Dparam['Lmaxval'], Dparam['Lstep']) ) #The extra value is 1E-6 times the step ticks of the colorbar # print(Dparam.items()) -# +# # # Read on which parameters initialize the default values # for args_t in list(Dparam.items()): # Test on all arguments present, if they are empty list, default values apply for each graph # args = list(args_t) @@ -936,16 +1072,16 @@ class PanelPlot(): # elif args[0] == 'Lstepticks' and not args[1]: args[1] = [None]*l # Dparam[args[0]] = args[1] # print(args) -# +# # # Check if there is no value for a parameter -## for args_t in list(Dparam.items()): +# for args_t in list(Dparam.items()): ## args = list(args_t) -## if args[1] +# if args[1] # # Number of contours level -## for i in range(l): -## if 'Lstepticks' in arguments.args and 'Lmaxval' in arguments.args and 'Lminval' in arguments.args: -## Lstep[i] = (Lmaxval[i] - Lminval[i])/20 # Default value of number of steps is 20 -## Lstepticks[i] = Lstep[i] # Default value is stepticks are the same as steps values +# for i in range(l): +# if 'Lstepticks' in arguments.args and 'Lmaxval' in arguments.args and 'Lminval' in arguments.args: +# Lstep[i] = (Lmaxval[i] - Lminval[i])/20 # Default value of number of steps is 20 +# Lstepticks[i] = Lstep[i] # Default value is stepticks are the same as steps values +# # -# # return Dparam diff --git a/src/LIB/Python/misc_functions.py b/src/LIB/Python/misc_functions.py index efc0916e04a60a7ee3cd2a276d1bb1eba3de0182..53577de6ee71fbc12524c19af6daff7657bd5fe9 100644 --- a/src/LIB/Python/misc_functions.py +++ b/src/LIB/Python/misc_functions.py @@ -13,207 +13,214 @@ from scipy.interpolate import RectBivariateSpline import numpy as np import math + def convert_date(datesince, time_in_sec): - return str(time_in_sec) + datesince[:33] + return str(time_in_sec) + datesince[:33] + class mean_operator(): - def MYM(self,var): + def MYM(self, var): ny = var.shape[1] out = copy.deepcopy(var) - for j in range(ny-1): - out[:,j,:] = (var[:,j,:] + var[:,j+1,:])*0.5 + for j in range(ny - 1): + out[:, j, :] = (var[:, j, :] + var[:, j + 1, :]) * 0.5 return out - - def MXM(self,var): + + def MXM(self, var): nx = var.shape[2] out = copy.deepcopy(var) - for i in range(nx-1): - out[:,:,i] = (var[:,:,i] + var[:,:,i+1])*0.5 + for i in range(nx - 1): + out[:, :, i] = (var[:, :, i] + var[:, :, i + 1]) * 0.5 return out - - def MZM(self,var): + + def MZM(self, var): nz = var.shape[0] out = copy.deepcopy(var) - for k in range(nz-1): - out[k,:,:] = (var[k,:,:] + var[k+1,:,:])*0.5 + for k in range(nz - 1): + out[k, :, :] = (var[k, :, :] + var[k + 1, :, :]) * 0.5 return out + def windvec_verti_proj(u, v, level, angle): """Compute the projected horizontal wind vector on an axis with a given angle w.r.t. the x/ni axes (West-East) - + Parameters ---------- u : array 3D U-wind component - + v : array 3D V-wind component level : array 1D level dimension array - + angle : float angle (radian) of the new axe w.r.t the x/ni axes (West-East). angle = 0 for (z,x) sections, angle=pi/2 for (z,y) sections - + Returns ------- - + projected_wind : array 3D a 3D wind component projected on the axe to be used with Panel_Plot.pvector as Lvar1 - """ + """ projected_wind = copy.deepcopy(u) for k in range(len(level)): - projected_wind[k,:,:] = u[k,:,:]*math.cos(angle) + v[k,:,:]*math.sin(angle) + projected_wind[k, :, :] = u[k, :, :] * math.cos(angle) + v[k, :, :] * math.sin(angle) return projected_wind + def oblique_proj(var, ni, nj, lvl, i_beg, j_beg, i_end, j_end): """Compute an oblique projection of a variable w.r.t. its axes - + Parameters ---------- var : array 3D or 2D the variable to project (e.g. THT, ZS) - + ni : array 1D 1D x-axis of the 3D dimension - + nj : array 1D 1D y-axis of the 3D dimension - + level : array 1D 1D z-axe of the 3D dimension - + i_beg, j_beg : int coordinate of the begin point of the new axe - + i_end, j_end : int - coordinate of the end point of the new axe - + coordinate of the end point of the new axe + Returns ------- angle_proj : float the angle (radian) of the new axe w.r.t the x/ni axes (West-East) - + out_var : array 2D or 1D a 2D (z,m) or 1D (m) variable projected on the oblique axe - + axe_m : array 1D a 1D m new axe (distance from the beggining point) - + """ - dist_seg=np.sqrt((i_end-i_beg)**2.0 + (j_end-j_beg)**2.0) # Distance de la section oblique m - if var.ndim ==3: - out_var = np.zeros((len(lvl),int(dist_seg)+1)) # Initialisation du nouveau champs projeté dans la coupe (z,m) + dist_seg = np.sqrt((i_end - i_beg)**2.0 + (j_end - j_beg)**2.0) # Distance de la section oblique m + if var.ndim == 3: + out_var = np.zeros((len(lvl), int(dist_seg) + 1)) # Initialisation du nouveau champs projeté dans la coupe (z,m) else: # 2D - out_var = np.zeros(int(dist_seg)+1) # Initialisation du nouveau champs projeté dans la coupe (m) + out_var = np.zeros(int(dist_seg) + 1) # Initialisation du nouveau champs projeté dans la coupe (m) - axe_m = np.zeros(int(dist_seg)+1) # Axe des abscisses qui sera tracé selon la coupe - axe_m_coord = [] # Coordonnées x,y des points qui composent l'axe - axe_m_coord.append( (ni[i_beg],nj[j_beg]) ) # Le premier point est celui donné par l'utilisateur + axe_m = np.zeros(int(dist_seg) + 1) # Axe des abscisses qui sera tracé selon la coupe + axe_m_coord = [] # Coordonnées x,y des points qui composent l'axe + axe_m_coord.append((ni[i_beg], nj[j_beg])) # Le premier point est celui donné par l'utilisateur - for m in range(int(dist_seg)): # Discrétisation selon distance de la coupe / int(distance_de_la_coupe) - axe_m_coord.append( (axe_m_coord[0][0] + (ni[i_end]-ni[i_beg])/(int(dist_seg))*(m+1), - axe_m_coord[0][1] + (nj[j_end]-nj[j_beg])/(int(dist_seg))*(m+1) )) - axe_m[m+1] = np.sqrt((ni[i_beg]-axe_m_coord[m+1][0])**2 + (nj[j_beg]-axe_m_coord[m+1][1])**2) + for m in range(int(dist_seg)): # Discrétisation selon distance de la coupe / int(distance_de_la_coupe) + axe_m_coord.append((axe_m_coord[0][0] + (ni[i_end] - ni[i_beg]) / (int(dist_seg)) * (m + 1), + axe_m_coord[0][1] + (nj[j_end] - nj[j_beg]) / (int(dist_seg)) * (m + 1))) + axe_m[m + 1] = np.sqrt((ni[i_beg] - axe_m_coord[m + 1][0])**2 + (nj[j_beg] - axe_m_coord[m + 1][1])**2) - if var.ndim ==3: # 3D variables to project + if var.ndim == 3: # 3D variables to project for k in range(len(lvl)): - a=RectBivariateSpline(ni, nj,var[k,:,:],kx=1,ky=1) # Interpolation par niveau à l'ordre 1 pour éviter des valeurs négatives de champs strictement > 0 - for m in range(int(dist_seg)+1): - out_var[k,m] = a.ev(axe_m_coord[m][0],axe_m_coord[m][1]) # La fonction ev de RectBivariate retourne la valeur la plus proche du point considéré + a = RectBivariateSpline(nj, ni, var[k, :, :], kx=1, ky=1) + for m in range(int(dist_seg) + 1): + # La fonction ev de RectBivariate retourne la valeur la plus proche du point considéré + out_var[k, m] = a.ev(axe_m_coord[m][1], axe_m_coord[m][0]) else: # 2D variables to project - a=RectBivariateSpline(ni, nj,var[:,:],kx=1,ky=1) - for m in range(int(dist_seg)+1): - out_var[m] = a.ev(axe_m_coord[m][0],axe_m_coord[m][1]) - - angle_proj = math.acos((ni[i_end]-ni[i_beg])/axe_m[-1]) + a = RectBivariateSpline(nj, ni, var[:, :], kx=1, ky=1) + for m in range(int(dist_seg) + 1): + out_var[m] = a.ev(axe_m_coord[m][1], axe_m_coord[m][0]) + + angle_proj = math.acos((ni[i_end] - ni[i_beg]) / axe_m[-1]) return angle_proj, out_var, axe_m + def comp_altitude1DVar(oneVar2D, orography, ztop, level, n_xory): """Compute and returns an altitude and x or y grid mesh variable in 2D following the topography in 1D To be used with 2D simulations - + Parameters ---------- oneVar2D : array 2D a 2D array (e.g. UT, THT) - + orography : array 1D 1D orography (ZS) - + ztop : real scalar of the top height of the model (ZTOP) - + level : array 1D 1D level variable (level or level_w) - - n_xory : array 1D + + n_xory : array 1D 1D directionnal grid variable (ni_u, nj_u, ni_v or nj_v) - + Returns ------- altitude a 2D altitude variable with topography taken into account - + n_xory_2D a 2D directionnal variable duplicated from n_xory """ n_xory_2D = copy.deepcopy(oneVar2D) - altitude = copy.deepcopy(oneVar2D) - + altitude = copy.deepcopy(oneVar2D) + for k in range(len(level)): - n_xory_2D[k,:] = n_xory + n_xory_2D[k, :] = n_xory for j in range(len(n_xory)): for k in range(len(level)): - altitude[k,j] = orography[j] + level[k]*((ztop-orography[j])/ztop) + altitude[k, j] = orography[j] + level[k] * ((ztop - orography[j]) / ztop) return altitude, n_xory_2D + def comp_altitude2DVar(oneVar3D, orography, ztop, level, n_y, n_x): """Compute and returns an altitude and x or y grid mesh variable in 3D following the topography in 2D To be used with 3D simulations - + Parameters ---------- oneVar3D : array 3D a 3D array (e.g. UT, THT) - + orography : array 2D 2D orography (ZS) - + ztop : real scalar of the top height of the model (ZTOP) - + level : array 1D 1D level variable (level or level_w) - - n_x : array 1D + + n_x : array 1D 1D directionnal grid variable along i (ni_u, or ni_v) - - n_y : array 1D + + n_y : array 1D 1D directionnal grid variable along j (nj_u, or nj_v) - + Returns ------- altitude a 3D altitude variable with topography taken into account - + n_x3D a 3D directionnal variable duplicated from n_x - + n_y3D a 3D directionnal variable duplicated from n_y """ n_x3D = copy.deepcopy(oneVar3D) n_y3D = copy.deepcopy(oneVar3D) - altitude = copy.deepcopy(oneVar3D) + altitude = copy.deepcopy(oneVar3D) for i in range(len(level)): - n_y3D[i,:] = n_y - n_x3D[i,:] = n_x + n_y3D[i, :] = n_y + n_x3D[i, :] = n_x for i in range(oneVar3D.shape[2]): for j in range(oneVar3D.shape[1]): - if ztop==0: - altitude[:,i,j] = level[:] + if ztop == 0: + altitude[:, i, j] = level[:] else: for k in range(len(level)): - altitude[k,j,i] = orography[j,i] + level[k]*((ztop-orography[j,i])/ztop) + altitude[k, j, i] = orography[j, i] + level[k] * ((ztop - orography[j, i]) / ztop) return altitude, n_x3D, n_y3D diff --git a/src/LIB/Python/read_MNHfile.py b/src/LIB/Python/read_MNHfile.py index 38072becffe3b1d8c066c6f370683b4988cf7ccf..7e240f65c91f1abbb80576f99655cdad016f080d 100644 --- a/src/LIB/Python/read_MNHfile.py +++ b/src/LIB/Python/read_MNHfile.py @@ -11,33 +11,67 @@ MNH_LIC for details. version 1. import netCDF4 as nc import numpy as np -def read_withEPY(LnameFiles,Dvar_input, Dvar_output={}, path='.'): + +def read_withEPY(LnameFiles, Dvar_input, Dvar_output={}, path='.'): + """Read a netCDF4 Meso-NH file, LFI, FA or GRIB 1/2 file with the Epygram Meteo-France library + Parameters + ---------- + LnameFiles : list of str + list of Meso-NH netCDF4 file (diachronic or synchronous) + + Dvar_input : Dict{'keyFile' : [field_identifier]} + where + 'keyFile' is a str corresponding to a key for the file number in LnameFiles (by order) + [field_identifier] is + - GRIB1 [indicatorOfParameter, paramId, indicatorOfTypeOfLevel, level, casual name] for GRIB1 + - GRIB2 [discipline, parameterCategory, typeOfFirstFixedSurface, parameterNumber, level, casual name] for GRIB2 + - MNH netcdf ('variable string', level) for netcdf MNH; set level=0 for 2D variables + - LFI ('variable string', level) for LFI; set level=0 for 2D variables + - FA 'variable string' + + path : str or list + path(s) of the files to read + + Returns + ------- + Dvar : Dict + Dvar[ifile]['var_name'] an epygram list of resources + + """ import epygram epygram.init_env() - for i,keyFiles in enumerate(Dvar_input.keys()): + for i, keyFiles in enumerate(Dvar_input.keys()): + if isinstance(path, list): + ipath = path[i] + else: + ipath = path print('Reading file ' + keyFiles) - theFile = epygram.formats.resource(LnameFiles[i],'r') - Dvar_output[keyFiles] = {} #initialize dic for each files - for var in Dvar_input[keyFiles]: #For each files + theFile = epygram.formats.resource(ipath + LnameFiles[i], 'r') + Dvar_output[keyFiles] = {} # initialize dic for each files + for var in Dvar_input[keyFiles]: # For each files # Read variables if(theFile.format == 'FA'): Dvar_output[keyFiles][var] = theFile.readfield(var) elif(theFile.format == 'LFI'): - if(var[1]==None or var[1]==0): # 2D Field + if(var[1] is None or var[1] == 0): # 2D Field Dvar_output[keyFiles][var[0]] = theFile.readfield(var) - else: # 3D Field - Dvar_output[keyFiles][var[0]+str(var[1])] = theFile.readfield(var).getlevel(k=var[1]) + else: # 3D Field + Dvar_output[keyFiles][var[0] + str(var[1])] = theFile.readfield(var).getlevel(k=var[1]) elif(theFile.format == 'netCDFMNH'): - if(var[1]==None or var[1]==0): # 2D Field + if(var[1] is None or var[1] == 0): # 2D Field Dvar_output[keyFiles][var[0]] = theFile.readfield(var[0]) else: - Dvar_output[keyFiles][var[0]+str(var[1])] = theFile.readfield(var[0]).getlevel(k=var[1]) - elif(theFile.format == 'GRIB'): - if len(var)==6: # GRIB2 - Dvar_output[keyFiles][var[5]] = theFile.readfield({'discipline': var[0], 'parameterCategory': var[1], 'typeOfFirstFixedSurface': var[2],'parameterNumber': var[3], 'level': var[4]}) - elif len(var)==5: # GRIB1 - Dvar_output[keyFiles][var[4]] = theFile.readfield({'indicatorOfParameter': var[0], 'paramId': var[1], 'indicatorOfTypeOfLevel': var[2], 'level': var[3]}) - else: epygramError("GRIB format error. GRIB1 expects 4 values : [indicatorOfParameter, paramId, indicatorOfTypeOfLevel, level, 'casual name'], GRIB2 expects 5 values [discipline, parameterCategory, typeOfFirstFixedSurface, parameterNumber, level, casual name]") + Dvar_output[keyFiles][var[0] + str(var[1])] = theFile.readfield(var[0]).getlevel(k=var[1]) + elif(theFile.format == 'GRIB'): + if len(var) == 6: # GRIB2 + Dvar_output[keyFiles][var[5]] = theFile.readfield( + {'discipline': var[0], 'parameterCategory': var[1], 'typeOfFirstFixedSurface': var[2], 'parameterNumber': var[3], 'level': var[4]}) + elif len(var) == 5: # GRIB1 + Dvar_output[keyFiles][var[4]] = theFile.readfield( + {'indicatorOfParameter': var[0], 'paramId': var[1], 'indicatorOfTypeOfLevel': var[2], 'level': var[3]}) + else: + epygramError( + "GRIB format error. GRIB1 expects 4 values : [indicatorOfParameter, paramId, indicatorOfTypeOfLevel, level, 'casual name'], GRIB2 expects 5 values [discipline, parameterCategory, typeOfFirstFixedSurface, parameterNumber, level, casual name]") else: raise epygramError("Unknown format file, please use FA, LFI, GRIB or MNH NetCDF") theFile.close() @@ -45,22 +79,23 @@ def read_withEPY(LnameFiles,Dvar_input, Dvar_output={}, path='.'): # Transform spectral data to physics space (for AROME and ARPEGE) for f in Dvar_output: for var in Dvar_output[f]: - try: + try: if(Dvar_output[f][var].spectral): Dvar_output[f][var].sp2gp() - except: + except BaseException: break return Dvar_output + def read_netcdf(LnameFiles, Dvar_input, path='.', get_data_only=True, del_empty_dim=True, removeHALO=True): """Read a netCDF4 Meso-NH file For each file, call functions to read diachronic or synchronous file - + Parameters ---------- LnameFiles : list of str list of Meso-NH netCDF4 file (diachronic or synchronous) - + Dvar_input : Dict{'keyFile' : 'var_name',('group_name','var_name')} where 'keyFile' is a str corresponding to a key for the file number in LnameFiles (by order) @@ -69,211 +104,225 @@ def read_netcdf(LnameFiles, Dvar_input, path='.', get_data_only=True, del_empty_ e.g. : {'f1':['ZS', 'WT','ni', 'level'], 'f2':[('/LES_budgets/Cartesian/Not_time_averaged/Not_normalized/cart/',MEAN_TH'),('/Budgets/RI','AVEF')] } - - path : str - unique path of the files - + + path : str or list + path(s) of the files to read + get_data_only : bool, default: True if True, the function returns Dvar as masked_array (only data) if False, the function returns Dvar as netCDF4._netCDF4.Variable - + del_empty_dim : bool, default: True if get_data_only=True and del_empty_dim=True, returns Dvar as an array without dimensions with size 1 and 0 e.g. : an array of dimensions (time_budget, cart_level, cart_nj, cart_ni) with shape (180,1,50,1) is returned (180,50) - + removeHALO : bool, default: True - if True, remove first and last (NHALO=1) point [1:-1] if get_data_only=True on each + if True, remove first and last (NHALO=1) point [1:-1] if get_data_only=True on each level, level_w, ni, ni_u, ni_v, nj, nj_u, nj_v dimensions - + Returns ------- - Dvar : Dict + Dvar : Dict Dvar[ifile]['var_name'] if the group contains only one variable Dvar[ifile][('group_name','var_name')] if the group contains more than one variable """ Dvar = {} - for i,keyFiles in enumerate(Dvar_input.keys()): + for i, keyFiles in enumerate(Dvar_input.keys()): + if isinstance(path, list): + ipath = path[i] + else: + ipath = path print('Reading file ' + keyFiles) - print(path + LnameFiles[i]) - theFile = nc.Dataset(path + LnameFiles[i],'r') + print(ipath + LnameFiles[i]) + theFile = nc.Dataset(ipath + LnameFiles[i], 'r') Dvar[keyFiles] = {} - if '000' in LnameFiles[i][-6:-3]: + if '000' in LnameFiles[i][-6:-3]: if theFile['MASDEV'][0] <= 54: raise TypeError('The python lib is available for MNH >= 5.5') else: Dvar[keyFiles] = read_TIMESfiles_55(theFile, Dvar_input[keyFiles], Dvar[keyFiles], get_data_only, del_empty_dim, removeHALO) else: - Dvar[keyFiles]= read_BACKUPfile(theFile, Dvar_input[keyFiles], Dvar[keyFiles], get_data_only, del_empty_dim, removeHALO) - #theFile.close() + Dvar[keyFiles] = read_BACKUPfile(theFile, Dvar_input[keyFiles], Dvar[keyFiles], get_data_only, del_empty_dim, removeHALO) + # theFile.close() return Dvar + def read_var(theFile, Dvar, var_name, get_data_only, del_empty_dim=True, removeHALO=True): """Read a netCDF4 variable - + Parameters ---------- theFile : netCDF4._netCDF4.Dataset a Meso-NH diachronic netCDF4 file - + var_name : str a Meso-NH netCDF4 variable name - + get_data_only : bool, default: True if True, the function returns Dvar as masked_array (only data) if False, the function returns Dvar as netCDF4._netCDF4.Variable - + del_empty_dim : bool, default: True if get_data_only=True and del_empty_dim=True, returns Dvar as an array without dimensions with size 1 and 0 e.g. : an array of dimensions (time_budget, cart_level, cart_nj, cart_ni) with shape (180,1,50,1) is returned (180,50) - + removeHALO : bool, default: True - if True, remove first and last (NHALO=1) point [1:-1] if get_data_only=True on each + if True, remove first and last (NHALO=1) point [1:-1] if get_data_only=True on each level, level_w, ni, ni_u, ni_v, nj, nj_u, nj_v dimensions - + Returns ------- - Dvar : Dict + Dvar : Dict Dvar['var_name'] if the group contains only one variable Dvar[('group_name','var_name')] if the group contains more than one variable """ try: var_dim = theFile.variables[var_name].ndim var_dim_name = theFile.variables[var_name].dimensions - except: - raise KeyError("Group and variable name not found in the file, check the group/variable name with ncdump -h YourMNHFile.000.nc. You asked for variable : " + var_name) + except BaseException: + raise KeyError( + "Group and variable name not found in the file, check the group/variable name with ncdump -h YourMNHFile.000.nc. You asked for variable : " + + var_name) if not get_data_only: Dvar[var_name] = theFile.variables[var_name] - else: + else: if var_dim == 0: Dvar[var_name] = theFile.variables[var_name][0].data elif var_dim == 1: Dvar[var_name] = theFile.variables[var_name][:] elif var_dim == 2: - Dvar[var_name] = theFile.variables[var_name][:,:] + Dvar[var_name] = theFile.variables[var_name][:, :] elif var_dim == 3: - Dvar[var_name] = theFile.variables[var_name][:,:,:] + Dvar[var_name] = theFile.variables[var_name][:, :, :] elif var_dim == 4: - Dvar[var_name] = theFile.variables[var_name][:,:,:,:] + Dvar[var_name] = theFile.variables[var_name][:, :, :, :] elif var_dim == 5: - Dvar[var_name] = theFile.variables[var_name][:,:,:,:,:] + Dvar[var_name] = theFile.variables[var_name][:, :, :, :, :] elif var_dim == 6: - Dvar[var_name] = theFile.variables[var_name][:,:,:,:,:,:] + Dvar[var_name] = theFile.variables[var_name][:, :, :, :, :, :] elif var_dim == 7: - Dvar[var_name] = theFile.variables[var_name][:,:,:,:,:,:,:] + Dvar[var_name] = theFile.variables[var_name][:, :, :, :, :, :, :] if removeHALO: for i in range(8): try: - if var_dim_name[i]=='level' or var_dim_name[i]=='level_w' or \ - var_dim_name[i]=='ni' or var_dim_name[i]=='ni_u' or var_dim_name[i]=='ni_v' or \ - var_dim_name[i]=='nj' or var_dim_name[i]=='nj_u' or var_dim_name[i]=='nj_v': + if var_dim_name[i] == 'level' or var_dim_name[i] == 'level_w' or \ + var_dim_name[i] == 'ni' or var_dim_name[i] == 'ni_u' or var_dim_name[i] == 'ni_v' or \ + var_dim_name[i] == 'nj' or var_dim_name[i] == 'nj_u' or var_dim_name[i] == 'nj_v': if var_dim != 0: - Dvar[var_name] = removetheHALO(i+1, Dvar[var_name]) - except: + Dvar[var_name] = removetheHALO(i + 1, Dvar[var_name]) + except BaseException: break if del_empty_dim: - Ldimtosqueeze=[] + Ldimtosqueeze = [] var_shape = theFile.variables[var_name].shape for i in range(8): try: - if var_shape[i]==1: Ldimtosqueeze.append(i) + if var_shape[i] == 1: + Ldimtosqueeze.append(i) except IndexError: break - Ldimtosqueeze=tuple(Ldimtosqueeze) + Ldimtosqueeze = tuple(Ldimtosqueeze) Dvar[var_name] = np.squeeze(Dvar[var_name], axis=Ldimtosqueeze) return Dvar -def read_from_group(theFile, Dvar, group_name, var_name, get_data_only, del_empty_dim=True,removeHALO=True): - """Read a variable from a netCDF4 group - + +def read_from_group(theFile, Dvar, group_name, var_name, get_data_only, del_empty_dim=True, removeHALO=True): + """Read a variable from a netCDF4 group + Parameters ---------- theFile : netCDF4._netCDF4.Dataset a Meso-NH diachronic netCDF4 file - + group_name : str a Meso-NH netCDF4 group name - + var_name : str a Meso-NH netCDF4 variable name - + get_data_only : bool, default: True if True, the function returns Dvar as masked_array (only data) if False, the function returns Dvar as netCDF4._netCDF4.Variable - + del_empty_dim : bool, default: True if get_data_only=True and del_empty_dim=True, returns Dvar as an array without dimensions with size 1 and 0 e.g. : an array of dimensions (time_budget, cart_level, cart_nj, cart_ni) with shape (180,1,50,1) is returned (180,50) - + removeHALO : bool, default: True - if True, remove first and last (NHALO=1) point [1:-1] if get_data_only=True on each + if True, remove first and last (NHALO=1) point [1:-1] if get_data_only=True on each level, level_w, ni, ni_u, ni_v, nj, nj_u, nj_v dimensions - + Returns ------- - Dvar : Dict + Dvar : Dict Dvar['var_name'] if the group contains only one variable Dvar[('group_name','var_name')] if the group contains more than one variable """ try: var_dim = theFile[group_name].variables[var_name].ndim var_dim_name = theFile[group_name].variables[var_name].dimensions - except: - raise KeyError("Group and variable name not found in the file, check the group/variable name with ncdump -h YourMNHFile.000.nc. You asked for group/variable : " + group_name + var_name) - + except BaseException: + raise KeyError( + "Group and variable name not found in the file, check the group/variable name with ncdump -h YourMNHFile.000.nc. You asked for group/variable : " + + group_name + + var_name) + if not get_data_only: - Dvar[(group_name,var_name)] = theFile[group_name].variables[var_name] + Dvar[(group_name, var_name)] = theFile[group_name].variables[var_name] else: if var_dim == 0: - Dvar[(group_name,var_name)] = theFile[group_name].variables[var_name][0].data + Dvar[(group_name, var_name)] = theFile[group_name].variables[var_name][0].data if var_dim == 1: - Dvar[(group_name,var_name)] = theFile[group_name].variables[var_name][:] + Dvar[(group_name, var_name)] = theFile[group_name].variables[var_name][:] elif var_dim == 2: - Dvar[(group_name,var_name)] = theFile[group_name].variables[var_name][:,:] + Dvar[(group_name, var_name)] = theFile[group_name].variables[var_name][:, :] elif var_dim == 3: - Dvar[(group_name,var_name)] = theFile[group_name].variables[var_name][:,:,:] + Dvar[(group_name, var_name)] = theFile[group_name].variables[var_name][:, :, :] elif var_dim == 4: - Dvar[(group_name,var_name)] = theFile[group_name].variables[var_name][:,:,:,:] + Dvar[(group_name, var_name)] = theFile[group_name].variables[var_name][:, :, :, :] elif var_dim == 5: - Dvar[(group_name,var_name)] = theFile[group_name].variables[var_name][:,:,:,:,:] + Dvar[(group_name, var_name)] = theFile[group_name].variables[var_name][:, :, :, :, :] elif var_dim == 6: - Dvar[(group_name,var_name)] = theFile[group_name].variables[var_name][:,:,:,:,:,:] + Dvar[(group_name, var_name)] = theFile[group_name].variables[var_name][:, :, :, :, :, :] elif var_dim == 7: - Dvar[(group_name,var_name)] = theFile[group_name].variables[var_name][:,:,:,:,:,:,:] + Dvar[(group_name, var_name)] = theFile[group_name].variables[var_name][:, :, :, :, :, :, :] if removeHALO: for i in range(8): try: - if var_dim_name[i]=='level' or var_dim_name[i]=='level_w' or \ - var_dim_name[i]=='ni' or var_dim_name[i]=='ni_u' or var_dim_name[i]=='ni_v' or \ - var_dim_name[i]=='nj' or var_dim_name[i]=='nj_u' or var_dim_name[i]=='nj_v': + if var_dim_name[i] == 'level' or var_dim_name[i] == 'level_w' or \ + var_dim_name[i] == 'ni' or var_dim_name[i] == 'ni_u' or var_dim_name[i] == 'ni_v' or \ + var_dim_name[i] == 'nj' or var_dim_name[i] == 'nj_u' or var_dim_name[i] == 'nj_v': if var_dim != 0: - Dvar[(group_name,var_name)] = removetheHALO(i+1, Dvar[(group_name,var_name)]) - except: + Dvar[(group_name, var_name)] = removetheHALO(i + 1, Dvar[(group_name, var_name)]) + except BaseException: break if del_empty_dim: - Ldimtosqueeze=[] - var_shape = Dvar[(group_name,var_name)].shape + Ldimtosqueeze = [] + var_shape = Dvar[(group_name, var_name)].shape for i in range(8): try: - if var_shape[i]==1: Ldimtosqueeze.append(i) + if var_shape[i] == 1: + Ldimtosqueeze.append(i) except IndexError: break - Ldimtosqueeze=tuple(Ldimtosqueeze) - Dvar[(group_name,var_name)] = np.squeeze(Dvar[(group_name,var_name)], axis=Ldimtosqueeze) - + Ldimtosqueeze = tuple(Ldimtosqueeze) + Dvar[(group_name, var_name)] = np.squeeze(Dvar[(group_name, var_name)], axis=Ldimtosqueeze) + # LES budget, ZTSERIES needs to be transposed to use psection functions without specifying .T each time if 'LES_budget' in group_name or 'ZTSERIES' in group_name or 'XTSERIES' in group_name: - Dvar[(group_name,var_name)] = Dvar[(group_name,var_name)].T + Dvar[(group_name, var_name)] = Dvar[(group_name, var_name)].T return Dvar + def read_BACKUPfile(theFile, Dvar_input, Dvar, get_data_only, del_empty_dim=True, removeHALO=True): """Read variables from Meso-NH MASDEV >= 5.5.0 synchronous file For all variables in Dvar_input of one file, call functions to read the variable of the group+variable - + Parameters ---------- theFile : netCDF4._netCDF4.Dataset a Meso-NH diachronic netCDF4 file - + Dvar_input : Dict{'var_name',('group_name','var_name')} with 'var_name' is the exact str of the netCDF4 variable name @@ -281,35 +330,40 @@ def read_BACKUPfile(theFile, Dvar_input, Dvar, get_data_only, del_empty_dim=True e.g. : {'f1':['ZS', 'WT','ni', 'level'], 'f2':[('/LES_budgets/Cartesian/Not_time_averaged/Not_normalized/cart/',MEAN_TH'),('/Budgets/RI','AVEF')] } - + get_data_only: bool, default: True if True, the function returns Dvar as masked_array (only data) if False, the function returns Dvar as netCDF4._netCDF4.Variable - + del_empty_dim: bool, default: True if get_data_only=True and del_empty_dim=True, returns Dvar as masked_array without dimensions with size 1 and 0 e.g. : an array of dimensions (time_budget, cart_level, cart_nj, cart_ni) with shape (180,1,50,1) is returned (180,50) - + Returns ------- - Dvar : Dict + Dvar : Dict Dvar['var_name'] if the group contains only one variable Dvar[('group_name','var_name')] if the group contains more than one variable - """ - # Reading date since beginning of the model run - Dvar['time'] = theFile.variables['time'][0] - Dvar['date'] = nc.num2date(Dvar['time'],units=theFile.variables['time'].units, calendar = theFile.variables['time'].calendar) + """ + # Reading date since beginning of the model run if theFile is a Meso-NH netCDF file + if 'MNH_cleanly_closed' in theFile.ncattrs(): + try: + Dvar['time'] = theFile.variables['time'][0] + Dvar['date'] = nc.num2date(Dvar['time'], units=theFile.variables['time'].units, calendar=theFile.variables['time'].calendar) + except BaseException: + pass for var in Dvar_input: - if type(var) == tuple: + if isinstance(var, tuple): Dvar = read_from_group(theFile, Dvar, var[0], var[1], get_data_only, del_empty_dim, removeHALO) else: Dvar = read_var(theFile, Dvar, var, get_data_only, del_empty_dim, removeHALO) # For all variables except scalars, change Fill_Value to NaN if get_data_only: - Dvar[var]= np.where(Dvar[var] != -99999.0, Dvar[var], np.nan) - Dvar[var]= np.where(Dvar[var] != 999.0, Dvar[var], np.nan) + Dvar[var] = np.where(Dvar[var] != -99999.0, Dvar[var], np.nan) + Dvar[var] = np.where(Dvar[var] != 999.0, Dvar[var], np.nan) return Dvar + def read_TIMESfiles_55(theFile, Dvar_input, Dvar, get_data_only, del_empty_dim=True, removeHALO=True): """Read variables from Meso-NH MASDEV >= 5.5.0 diachronic file For all variables in Dvar_input of one file, call functions to read the variable of the group+variable @@ -318,7 +372,7 @@ def read_TIMESfiles_55(theFile, Dvar_input, Dvar, get_data_only, del_empty_dim=T ---------- theFile : netCDF4._netCDF4.Dataset a Meso-NH diachronic netCDF4 file - + Dvar_input : Dict{'var_name',('group_name','var_name')} with 'var_name' is the exact str of the netCDF4 variable name @@ -326,57 +380,58 @@ def read_TIMESfiles_55(theFile, Dvar_input, Dvar, get_data_only, del_empty_dim=T e.g. : {'f1':['ZS', 'WT','ni', 'level'], 'f2':[('/LES_budgets/Cartesian/Not_time_averaged/Not_normalized/cart/',MEAN_TH'),('/Budgets/RI','AVEF')] } - + get_data_only: bool, default: True if True, the function returns Dvar as masked_array (only data) if False, the function returns Dvar as netCDF4._netCDF4.Variable - + del_empty_dim: bool, default: True if get_data_only=True and del_empty_dim=True, returns Dvar as masked_array without dimensions with size 1 and 0 e.g. : an array of dimensions (time_budget, cart_level, cart_nj, cart_ni) with shape (180,1,50,1) is returned (180,50) - + Returns ------- - Dvar : Dict + Dvar : Dict Dvar[ifile]['var_name'] if the group contains only one variable Dvar[ifile][('group_name','var_name')] if the group contains more than one variable - """ + """ for var in Dvar_input: - if type(var) == tuple: + if isinstance(var, tuple): Dvar = read_from_group(theFile, Dvar, var[0], var[1], get_data_only, del_empty_dim, removeHALO) else: Dvar = read_var(theFile, Dvar, var, get_data_only, del_empty_dim, removeHALO) return Dvar + def removetheHALO(idim, var): """Remove a NHALO=1 point [1:-1] at a given dimension idim of a variable var TODO: with NHALO=3 with WENO5 TODO: with the use of get_data_only=False (NetCDF4 output) - + Parameters ---------- idim: int the dimension over which remove the first and last point - + var: array a Meso-NH netCDF4 variable name - + Returns ------- - var : array - """ + var : array + """ if idim == 1: var = var[1:-1] elif idim == 2: - var = var[:,1:-1] + var = var[:, 1:-1] elif idim == 3: - var = var[:,:,1:-1] + var = var[:, :, 1:-1] elif idim == 4: - var = var[:,:,:,1:-1] + var = var[:, :, :, 1:-1] elif idim == 5: - var = var[:,:,:,:,1:-1] + var = var[:, :, :, :, 1:-1] elif idim == 6: - var = var[:,:,:,:,:,1:-1] + var = var[:, :, :, :, :, 1:-1] elif idim == 7: - var = var[:,:,:,:,:,:,1:-1] + var = var[:, :, :, :, :, :, 1:-1] return var diff --git a/src/LIB/SURCOUCHE/src/modd_field.f90 b/src/LIB/SURCOUCHE/src/modd_field.f90 index b81b59f1f0354c2041f9908b0d3c7fb7230e0f3d..cf3e014068caca76b99fbcba63b96680894a39d6 100644 --- a/src/LIB/SURCOUCHE/src/modd_field.f90 +++ b/src/LIB/SURCOUCHE/src/modd_field.f90 @@ -12,6 +12,7 @@ ! P. Wautelet 27/01/2020: create the tfield_metadata_base abstract datatype ! P. Wautelet 14/09/2020: add ndimlist field to tfield_metadata_base ! P. Wautelet 10/11/2020: new data structures for netCDF dimensions +! P. Wautelet 08/10/2021: add 2 new dimensions: LW_bands (NMNHDIM_NLWB) and SW_bands (NMNHDIM_NSWB) !----------------------------------------------------------------- module modd_field @@ -42,49 +43,52 @@ integer, parameter :: NMNHDIM_TIME = 9 integer, parameter :: NMNHDIM_ONE = 10 -integer, parameter :: NMNHDIM_LASTDIM_NODIACHRO = 10 ! Index of the last defined dimension for non-diachronic files +integer, parameter :: NMNHDIM_NSWB = 11 +integer, parameter :: NMNHDIM_NLWB = 12 -integer, parameter :: NMNHDIM_COMPLEX = 11 +integer, parameter :: NMNHDIM_LASTDIM_NODIACHRO = 12 ! Index of the last defined dimension for non-diachronic files -integer, parameter :: NMNHDIM_BUDGET_CART_NI = 12 -integer, parameter :: NMNHDIM_BUDGET_CART_NJ = 13 -integer, parameter :: NMNHDIM_BUDGET_CART_NI_U = 14 -integer, parameter :: NMNHDIM_BUDGET_CART_NJ_U = 15 -integer, parameter :: NMNHDIM_BUDGET_CART_NI_V = 16 -integer, parameter :: NMNHDIM_BUDGET_CART_NJ_V = 17 -integer, parameter :: NMNHDIM_BUDGET_CART_LEVEL = 18 -integer, parameter :: NMNHDIM_BUDGET_CART_LEVEL_W = 19 +integer, parameter :: NMNHDIM_COMPLEX = 13 -integer, parameter :: NMNHDIM_BUDGET_MASK_LEVEL = 20 -integer, parameter :: NMNHDIM_BUDGET_MASK_LEVEL_W = 21 -integer, parameter :: NMNHDIM_BUDGET_MASK_NBUMASK = 22 +integer, parameter :: NMNHDIM_BUDGET_CART_NI = 14 +integer, parameter :: NMNHDIM_BUDGET_CART_NJ = 15 +integer, parameter :: NMNHDIM_BUDGET_CART_NI_U = 16 +integer, parameter :: NMNHDIM_BUDGET_CART_NJ_U = 17 +integer, parameter :: NMNHDIM_BUDGET_CART_NI_V = 18 +integer, parameter :: NMNHDIM_BUDGET_CART_NJ_V = 19 +integer, parameter :: NMNHDIM_BUDGET_CART_LEVEL = 20 +integer, parameter :: NMNHDIM_BUDGET_CART_LEVEL_W = 21 -integer, parameter :: NMNHDIM_BUDGET_TIME = 23 +integer, parameter :: NMNHDIM_BUDGET_MASK_LEVEL = 22 +integer, parameter :: NMNHDIM_BUDGET_MASK_LEVEL_W = 23 +integer, parameter :: NMNHDIM_BUDGET_MASK_NBUMASK = 24 -integer, parameter :: NMNHDIM_BUDGET_LES_TIME = 24 -integer, parameter :: NMNHDIM_BUDGET_LES_AVG_TIME = 25 -integer, parameter :: NMNHDIM_BUDGET_LES_LEVEL = 26 -integer, parameter :: NMNHDIM_BUDGET_LES_SV = 27 -integer, parameter :: NMNHDIM_BUDGET_LES_PDF = 28 +integer, parameter :: NMNHDIM_BUDGET_TIME = 25 + +integer, parameter :: NMNHDIM_BUDGET_LES_TIME = 26 +integer, parameter :: NMNHDIM_BUDGET_LES_AVG_TIME = 27 +integer, parameter :: NMNHDIM_BUDGET_LES_LEVEL = 28 +integer, parameter :: NMNHDIM_BUDGET_LES_SV = 29 +integer, parameter :: NMNHDIM_BUDGET_LES_PDF = 30 integer, parameter :: NMNHDIM_BUDGET_LES_MASK = 100 ! This is not a true dimension -integer, parameter :: NMNHDIM_SPECTRA_2PTS_NI = 29 -integer, parameter :: NMNHDIM_SPECTRA_2PTS_NJ = 30 -integer, parameter :: NMNHDIM_SPECTRA_SPEC_NI = 31 -integer, parameter :: NMNHDIM_SPECTRA_SPEC_NJ = 32 -integer, parameter :: NMNHDIM_SPECTRA_LEVEL = 33 +integer, parameter :: NMNHDIM_SPECTRA_2PTS_NI = 31 +integer, parameter :: NMNHDIM_SPECTRA_2PTS_NJ = 32 +integer, parameter :: NMNHDIM_SPECTRA_SPEC_NI = 33 +integer, parameter :: NMNHDIM_SPECTRA_SPEC_NJ = 34 +integer, parameter :: NMNHDIM_SPECTRA_LEVEL = 35 -integer, parameter :: NMNHDIM_SERIES_LEVEL = 34 -integer, parameter :: NMNHDIM_SERIES_LEVEL_W = 35 -integer, parameter :: NMNHDIM_SERIES_TIME = 36 ! Time dimension for time series +integer, parameter :: NMNHDIM_SERIES_LEVEL = 36 +integer, parameter :: NMNHDIM_SERIES_LEVEL_W = 37 +integer, parameter :: NMNHDIM_SERIES_TIME = 38 ! Time dimension for time series -integer, parameter :: NMNHDIM_FLYER_TIME = 37 ! Time dimension for aircraft/balloon (dimension local to each flyer) -integer, parameter :: NMNHDIM_PROFILER_TIME = 38 ! Time dimension for profilers -integer, parameter :: NMNHDIM_STATION_TIME = 39 ! Time dimension for stations +integer, parameter :: NMNHDIM_FLYER_TIME = 39 ! Time dimension for aircraft/balloon (dimension local to each flyer) +integer, parameter :: NMNHDIM_PROFILER_TIME = 40 ! Time dimension for profilers +integer, parameter :: NMNHDIM_STATION_TIME = 41 ! Time dimension for stations -integer, parameter :: NMNHDIM_PAIR = 40 ! For values coming by pair (ie boundaries) +integer, parameter :: NMNHDIM_PAIR = 42 ! For values coming by pair (ie boundaries) -integer, parameter :: NMNHDIM_LASTDIM_DIACHRO = 40 ! Index of the last defined dimension for diachronic files +integer, parameter :: NMNHDIM_LASTDIM_DIACHRO = 42 ! Index of the last defined dimension for diachronic files integer, parameter :: NMNHDIM_BUDGET_NGROUPS = 101 ! This is not a true dimension integer, parameter :: NMNHDIM_FLYER_PROC = 102 ! This is not a true dimension diff --git a/src/LIB/SURCOUCHE/src/mode_field.f90 b/src/LIB/SURCOUCHE/src/mode_field.f90 index 648ac8b7546ff473b262245e84e3f186abd6fa6c..27ab3f729991c45afae925bbc6824681d8d470ad 100644 --- a/src/LIB/SURCOUCHE/src/mode_field.f90 +++ b/src/LIB/SURCOUCHE/src/mode_field.f90 @@ -3988,6 +3988,7 @@ IF (CPROGRAM == 'MESONH') THEN TFIELDLIST(IID)%TFIELD_X3D(KFROM)%DATA => XRRS_CLD(:,:,:,CONF_MODEL(KFROM)%IDX_RHT) END IF CALL FIND_FIELD_ID_FROM_MNHNAME('CLDFR',IID,IRESP); TFIELDLIST(IID)%TFIELD_X3D(KFROM)%DATA => XCLDFR + CALL FIND_FIELD_ID_FROM_MNHNAME('ICEFR',IID,IRESP); TFIELDLIST(IID)%TFIELD_X3D(KFROM)%DATA => XICEFR CALL FIND_FIELD_ID_FROM_MNHNAME('CIT', IID,IRESP); TFIELDLIST(IID)%TFIELD_X3D(KFROM)%DATA => XCIT CALL FIND_FIELD_ID_FROM_MNHNAME('RAINFR',IID,IRESP); TFIELDLIST(IID)%TFIELD_X3D(KFROM)%DATA => XRAINFR ! @@ -4317,6 +4318,7 @@ IF (CPROGRAM == 'MESONH') THEN CALL FIND_FIELD_ID_FROM_MNHNAME('WS_PRES',IID,IRESP); XRWS_PRES => TFIELDLIST(IID)%TFIELD_X3D(KTO)%DATA CALL FIND_FIELD_ID_FROM_MNHNAME('THS_CLD',IID,IRESP); XRTHS_CLD => TFIELDLIST(IID)%TFIELD_X3D(KTO)%DATA CALL FIND_FIELD_ID_FROM_MNHNAME('CLDFR', IID,IRESP); XCLDFR => TFIELDLIST(IID)%TFIELD_X3D(KTO)%DATA + CALL FIND_FIELD_ID_FROM_MNHNAME('ICEFR', IID,IRESP); XICEFR => TFIELDLIST(IID)%TFIELD_X3D(KTO)%DATA CALL FIND_FIELD_ID_FROM_MNHNAME('CIT', IID,IRESP); XCIT => TFIELDLIST(IID)%TFIELD_X3D(KTO)%DATA CALL FIND_FIELD_ID_FROM_MNHNAME('RAINFR', IID,IRESP); XRAINFR => TFIELDLIST(IID)%TFIELD_X3D(KTO)%DATA END IF diff --git a/src/LIB/SURCOUCHE/src/mode_io_field_write.f90 b/src/LIB/SURCOUCHE/src/mode_io_field_write.f90 index b48a37b2649d4d8ae06964e24a736db43439cbaf..9e95c6e518de2c22d3c373e9d8994cfe72e12a11 100644 --- a/src/LIB/SURCOUCHE/src/mode_io_field_write.f90 +++ b/src/LIB/SURCOUCHE/src/mode_io_field_write.f90 @@ -1,4 +1,4 @@ -!MNH_LIC Copyright 1994-2021 CNRS, Meteo-France and Universite Paul Sabatier +!MNH_LIC Copyright 1994-2022 CNRS, Meteo-France and Universite Paul Sabatier !MNH_LIC This is part of the Meso-NH software governed by the CeCILL-C licence !MNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt !MNH_LIC for details. version 1. @@ -1079,6 +1079,8 @@ end subroutine IO_Ndimlist_reduce iresp = 0 iresp_lfi = 0 iresp_nc4 = 0 + iresp_tmp_lfi = 0 + iresp_tmp_nc4 = 0 GALLOC = .FALSE. GALLOC_ll = .FALSE. IHEXTOT = 2*JPHEXT+1 diff --git a/src/LIB/SURCOUCHE/src/mode_io_file.f90 b/src/LIB/SURCOUCHE/src/mode_io_file.f90 index e269accf0320ac26d31879e4ea155bc61ecc2d17..6ed3a03c255e7690965e101065a6892fee7c1d7c 100644 --- a/src/LIB/SURCOUCHE/src/mode_io_file.f90 +++ b/src/LIB/SURCOUCHE/src/mode_io_file.f90 @@ -1,4 +1,4 @@ -!MNH_LIC Copyright 1994-2021 CNRS, Meteo-France and Universite Paul Sabatier +!MNH_LIC Copyright 1994-2022 CNRS, Meteo-France and Universite Paul Sabatier !MNH_LIC This is part of the Meso-NH software governed by the CeCILL-C licence !MNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt !MNH_LIC for details. version 1. @@ -38,6 +38,7 @@ ! P. Wautelet 12/03/2019: simplify opening of IO split files ! P. Wautelet 05/09/2019: disable IO_Coordvar_write_nc4 for Z-split files ! P. Wautelet 01/10/2020: bugfix: add missing initializations for IRESP +! P. Wautelet 19/08/2022: bugfix: IO_File_check_format_exist: broadcast cformat if changed !----------------------------------------------------------------- module mode_io_file @@ -666,9 +667,11 @@ end subroutine IO_Transfer_list_addto subroutine IO_File_check_format_exist( tpfile ) +use modd_mpif type(tfiledata), intent(inout) :: tpfile ! File structure +integer :: ierr logical :: gexist_lfi, gexist_nc4 @@ -720,6 +723,9 @@ IF (TPFILE%LMASTER) THEN end if MODE END IF +if ( tpfile%cmode == 'READ' ) & + call MPI_BCAST( tpfile%cformat, Len( tpfile%cformat ), MPI_CHARACTER, tpfile%nmaster_rank - 1, tpfile%nmpicomm, ierr ) + end subroutine IO_File_check_format_exist diff --git a/src/LIB/SURCOUCHE/src/mode_io_manage_struct.f90 b/src/LIB/SURCOUCHE/src/mode_io_manage_struct.f90 index f7954781fd46754679b6f25aeb2c6fc018a19b6b..60cf44de04d57ca5056621b06c3fda899d42bdb5 100644 --- a/src/LIB/SURCOUCHE/src/mode_io_manage_struct.f90 +++ b/src/LIB/SURCOUCHE/src/mode_io_manage_struct.f90 @@ -1,4 +1,4 @@ -!MNH_LIC Copyright 2016-2021 CNRS, Meteo-France and Universite Paul Sabatier +!MNH_LIC Copyright 2016-2022 CNRS, Meteo-France and Universite Paul Sabatier !MNH_LIC This is part of the Meso-NH software governed by the CeCILL-C licence !MNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt !MNH_LIC for details. version 1. @@ -19,6 +19,8 @@ ! P. Wautelet 11/02/2020: bugfix: TDADFILE was wrongly constructed for output files ! S. Donnier 28/02/2020: type STREAM needed for use of ECOCLIMAP SG ! P. Wautelet 08/01/2021: allow output files with empty variable list (useful if IO_Field_user_write is not empty) +! P. Wautelet 18/03/2022: minor bugfix in ISTEP_MAX computation + adapt diagnostics messages + (change verbosity level and remove some unnecessary warnings) !----------------------------------------------------------------- MODULE MODE_IO_MANAGE_STRUCT ! @@ -86,8 +88,8 @@ END IF DO IMI = 1, NMODEL IBAK_NUMB = 0 IOUT_NUMB = 0 - ISTEP_MAX = NINT(XSEGLEN/DYN_MODEL(IMI)%XTSTEP)+1 - IF (IMI == 1) ISTEP_MAX = ISTEP_MAX - KSUP + !Reduce XSEGLEN by time added to XSEGLEN for 1st domain (see set_grid subroutine) + ISTEP_MAX = NINT( ( XSEGLEN - KSUP * DYN_MODEL(1)%XTSTEP ) / DYN_MODEL(IMI)%XTSTEP ) + 1 ! !* Insert regular backups/outputs into XBAK_TIME/XOUT_TIME arrays ! @@ -468,7 +470,7 @@ SUBROUTINE FIND_REMOVE_DUPLICATES(KNUMB,KSTEPS) DO JOUT = 1,KNUMB DO JKLOOP = JOUT+1,KNUMB IF ( KSTEPS(JKLOOP) == KSTEPS(JOUT) .AND. KSTEPS(JKLOOP) > 0 ) THEN - CALL PRINT_MSG(NVERB_INFO,'IO','FIND_REMOVE_DUPLICATES','found duplicated backup/output step (removed extra one)') + CALL PRINT_MSG(NVERB_DEBUG,'IO','FIND_REMOVE_DUPLICATES','found duplicated backup/output step (removed extra one)') KSTEPS(JKLOOP) = NNEGUNDEF END IF END DO @@ -484,7 +486,8 @@ SUBROUTINE FIND_REMOVE_OUTOFTIMERANGE(KNUMB,KSTEPS) ! DO JOUT = 1,KNUMB IF ( KSTEPS(JOUT) < 1 .OR. KSTEPS(JOUT) > ISTEP_MAX ) THEN - CALL PRINT_MSG(NVERB_WARNING,'IO','FIND_REMOVE_OUTOFTIMERANGE','found backup/output step outside of time range') + IF ( KSTEPS(JOUT) /= NNEGUNDEF ) & + CALL PRINT_MSG(NVERB_WARNING,'IO','FIND_REMOVE_OUTOFTIMERANGE','found backup/output step outside of time range') KSTEPS(JOUT) = NNEGUNDEF END IF END DO diff --git a/src/LIB/SURCOUCHE/src/mode_io_tools_nc4.f90 b/src/LIB/SURCOUCHE/src/mode_io_tools_nc4.f90 index 477c389dbf6ffe6691029ff11a627e431b2e0a84..b436b0d18cf7763e6d5cdfc213050d36ec6cfec5 100644 --- a/src/LIB/SURCOUCHE/src/mode_io_tools_nc4.f90 +++ b/src/LIB/SURCOUCHE/src/mode_io_tools_nc4.f90 @@ -1,4 +1,4 @@ -!MNH_LIC Copyright 1994-2021 CNRS, Meteo-France and Universite Paul Sabatier +!MNH_LIC Copyright 1994-2022 CNRS, Meteo-France and Universite Paul Sabatier !MNH_LIC This is part of the Meso-NH software governed by the CeCILL-C licence !MNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt !MNH_LIC for details. version 1. @@ -19,6 +19,7 @@ ! P. Wautelet 18/03/2021: workaround for an intel compiler bug ! P. Wautelet 04/05/2021: improve IO_Vdims_fill_nc4 if l2d and lpack ! P. Wautelet 27/05/2021: improve IO_Mnhname_clean to autocorrect names to be CF compliant +! P. Wautelet 08/10/2021: add 2 new dimensions: LW_bands and SW_bands !----------------------------------------------------------------- #ifdef MNH_IOCDF4 module mode_io_tools_nc4 @@ -251,10 +252,10 @@ USE MODD_CONF, ONLY: CPROGRAM, l2d, lpack USE MODD_CONF_n, ONLY: CSTORAGE_TYPE USE MODD_DIM_n, ONLY: NIMAX_ll, NJMAX_ll, NKMAX use modd_dyn, only: xseglen -use modd_dyn_n, only: xtstep +use modd_dyn_n, only: dyn_model use modd_field, only: NMNHDIM_NI, NMNHDIM_NJ, NMNHDIM_NI_U, NMNHDIM_NJ_U, NMNHDIM_NI_V, NMNHDIM_NJ_V, & NMNHDIM_LEVEL, NMNHDIM_LEVEL_W, NMNHDIM_TIME, & - NMNHDIM_ONE, NMNHDIM_COMPLEX, & + NMNHDIM_ONE, NMNHDIM_NSWB, NMNHDIM_NLWB, NMNHDIM_COMPLEX, & NMNHDIM_BUDGET_CART_NI, NMNHDIM_BUDGET_CART_NJ, NMNHDIM_BUDGET_CART_NI_U, & NMNHDIM_BUDGET_CART_NJ_U, NMNHDIM_BUDGET_CART_NI_V, NMNHDIM_BUDGET_CART_NJ_V, & NMNHDIM_BUDGET_CART_LEVEL, NMNHDIM_BUDGET_CART_LEVEL_W, & @@ -274,7 +275,9 @@ use modd_les, only: lles_pdf, nles_k, npdf, nspectra_k, xles_temp_mean use modd_les_n, only: nles_times, nspectra_ni, nspectra_nj use modd_nsv, only: nsv USE MODD_PARAMETERS_ll, ONLY: JPHEXT, JPVEXT +use modd_param_n, only: crad use modd_profiler_n, only: numbprofiler, tprofiler +use modd_radiations_n, only: nlwb_mnh, nswb_mnh use modd_series, only: lseries use modd_series_n, only: nsnbstept use modd_station_n, only: numbstat, tstation @@ -336,6 +339,11 @@ if ( tpfile%ctype == 'MNHDIACHRONIC' .or. ( lpack .and. l2d ) ) then call IO_Add_dim_nc4( tpfile, NMNHDIM_ONE, 'one', 1 ) end if +if ( tpfile%ctype /= 'MNHDIACHRONIC' .and. crad /= 'NONE' ) then + call IO_Add_dim_nc4( tpfile, NMNHDIM_NSWB, 'SW_bands', nswb_mnh ) !number of SW bands practically used + call IO_Add_dim_nc4( tpfile, NMNHDIM_NLWB, 'LW_bands', nlwb_mnh ) !number of LW bands practically used +end if + !Write dimensions used in diachronic files if ( tpfile%ctype == 'MNHDIACHRONIC' ) then !Dimension of size 2 used for NMNHDIM_COMPLEX @@ -413,13 +421,13 @@ if ( tpfile%ctype == 'MNHDIACHRONIC' ) then !Dimension for the number of profiler times if ( numbprofiler > 0 ) then - iprof = Int ( ( xseglen - xtstep ) / tprofiler%step ) + 1 + iprof = Nint ( ( xseglen - dyn_model(1)%xtstep ) / tprofiler%step ) + 1 call IO_Add_dim_nc4( tpfile, NMNHDIM_PROFILER_TIME, 'time_profiler', iprof ) end if !Dimension for the number of station times if ( numbstat > 0 ) then - istation = Int ( ( xseglen - xtstep ) / tstation%step ) + 1 + istation = Nint ( ( xseglen - dyn_model(1)%xtstep ) / tstation%step ) + 1 call IO_Add_dim_nc4( tpfile, NMNHDIM_STATION_TIME, 'time_station', istation ) end if diff --git a/src/LIB/SURCOUCHE/src/mode_io_write_nc4.f90 b/src/LIB/SURCOUCHE/src/mode_io_write_nc4.f90 index 53053069c936e71c6df416594babe000686782a7..10ed5d33ba1b1993d9e3baad7c9c3a87f3f4a1c2 100644 --- a/src/LIB/SURCOUCHE/src/mode_io_write_nc4.f90 +++ b/src/LIB/SURCOUCHE/src/mode_io_write_nc4.f90 @@ -1,4 +1,4 @@ -!MNH_LIC Copyright 1994-2021 CNRS, Meteo-France and Universite Paul Sabatier +!MNH_LIC Copyright 1994-2022 CNRS, Meteo-France and Universite Paul Sabatier !MNH_LIC This is part of the Meso-NH software governed by the CeCILL-C licence !MNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt !MNH_LIC for details. version 1. @@ -28,6 +28,8 @@ ! P. Wautelet 14/01/2021: add IO_Field_write_nc4_N4, IO_Field_partial_write_nc4_N2, ! IO_Field_partial_write_nc4_N3 and IO_Field_partial_write_nc4_N4 subroutines ! P. Wautelet 30/03/2021: budgets: LES cartesian subdomain limits are defined in the physical domain +! P. Wautelet 22/03/2022: correct time_les_avg and time_les_avg_bounds coordinates +! P. Wautelet 21/06/2022: bugfix: time_budget was not computed correctly (tdtexp -> tdtseg) !----------------------------------------------------------------- #ifdef MNH_IOCDF4 module mode_io_write_nc4 @@ -86,6 +88,8 @@ use modd_parameters, only: jphext use mode_tools_ll, only: Get_globaldims_ll +use NETCDF, only: NF90_FLOAT + type(tfiledata), intent(in) :: tpfile type(tfielddata), intent(in) :: tpfield integer, intent(in) :: knblocks @@ -112,7 +116,11 @@ call IO_Mnhname_clean( tpfield%cmnhname, yvarname ) istatus = NF90_INQ_VARID( incid, yvarname, ivarid ) if ( istatus /= NF90_NOERR ) then - istatus = NF90_DEF_VAR( incid, yvarname, MNHREAL_NF90, ivarid) + if ( tpfile%lncreduce_float_precision ) then + istatus = NF90_DEF_VAR( incid, yvarname, NF90_FLOAT, ivarid ) + else + istatus = NF90_DEF_VAR( incid, yvarname, MNHREAL_NF90, ivarid ) + end if if ( tpfield%ndims /= 3 ) call Print_msg( NVERB_FATAL, 'IO', 'IO_Field_header_split_write_nc4', & trim( tpfile%cname )//': '//trim( yvarname )//': NDIMS should be 3' ) @@ -1444,9 +1452,9 @@ use modd_grid, only: xlatori, xlonori use modd_grid_n, only: lsleve, xxhat, xyhat, xzhat use modd_les, only: cles_level_type, cspectra_level_type, nlesn_iinf, nlesn_isup, nlesn_jinf, nlesn_jsup, & nles_k, nles_levels, nspectra_k, nspectra_levels, & - xles_altitudes, xles_temp_mean_start, xles_temp_mean_end, xles_temp_mean_step, & - xspectra_altitudes -use modd_les_n, only: nles_times, nspectra_ni, nspectra_nj, tles_dates + xles_altitudes, xspectra_altitudes +use modd_les_n, only: nles_dtcount, nles_mean_end, nles_mean_start, nles_mean_step, nles_mean_times, & + nles_times, nspectra_ni, nspectra_nj, tles_dates, xles_times use modd_netcdf, only: tdimnc use modd_parameters, only: jphext, JPVEXT use modd_profiler_n, only: numbprofiler, tprofiler @@ -1460,7 +1468,6 @@ use modd_type_date, only: date_time use mode_field, only: Find_field_id_from_mnhname use mode_gridproj, only: Sm_latlon use mode_nest_ll, only: Get_model_number_ll, Go_tomodel_ll -use mode_time, only: tdtexp type(tfiledata), intent(in) :: tpfile character(len=*), optional, intent(in) :: hprogram_orig !To emulate a file coming from this program @@ -1470,14 +1477,13 @@ character(len=:), allocatable :: yprogram integer :: iiu, iju, iku integer :: id, iid, iresp integer :: imi -integer :: iavg integer :: ji integer :: jt +integer :: jtb, jte integer(kind=cdfint) :: incid logical :: gchangemodel logical :: gdealloc logical, pointer :: gsleve -real :: zles_temp_mean_start, zles_temp_mean_end real, dimension(:), pointer :: zxhat, zyhat, zzhat real, dimension(:), allocatable :: zxhatm, zyhatm, zzhatm !Coordinates at mass points in the transformed space real, dimension(:), allocatable :: zles_levels @@ -1713,20 +1719,20 @@ if ( tpfile%lmaster ) then Allocate( tzdates_bound(2, nbutotwrite) ) do jt = 1, nbutotwrite - tzdates(jt)%nyear = tdtexp%nyear - tzdates(jt)%nmonth = tdtexp%nmonth - tzdates(jt)%nday = tdtexp%nday - tzdates(jt)%xtime = tdtexp%xtime + nbustep * ( ( jt - 1 ) + 0.5 ) * xtstep - - tzdates_bound(1, jt)%nyear = tdtexp%nyear - tzdates_bound(1, jt)%nmonth = tdtexp%nmonth - tzdates_bound(1, jt)%nday = tdtexp%nday - tzdates_bound(1, jt)%xtime = tdtexp%xtime + nbustep * ( jt - 1 ) * xtstep - - tzdates_bound(2, jt)%nyear = tdtexp%nyear - tzdates_bound(2, jt)%nmonth = tdtexp%nmonth - tzdates_bound(2, jt)%nday = tdtexp%nday - tzdates_bound(2, jt)%xtime = tdtexp%xtime + nbustep * jt * xtstep + tzdates(jt)%nyear = tdtseg%nyear + tzdates(jt)%nmonth = tdtseg%nmonth + tzdates(jt)%nday = tdtseg%nday + tzdates(jt)%xtime = tdtseg%xtime + nbustep * ( ( jt - 1 ) + 0.5 ) * xtstep + + tzdates_bound(1, jt)%nyear = tdtseg%nyear + tzdates_bound(1, jt)%nmonth = tdtseg%nmonth + tzdates_bound(1, jt)%nday = tdtseg%nday + tzdates_bound(1, jt)%xtime = tdtseg%xtime + nbustep * ( jt - 1 ) * xtstep + + tzdates_bound(2, jt)%nyear = tdtseg%nyear + tzdates_bound(2, jt)%nmonth = tdtseg%nmonth + tzdates_bound(2, jt)%nday = tdtseg%nday + tzdates_bound(2, jt)%xtime = tdtseg%xtime + nbustep * jt * xtstep end do call Write_time_coord( tpfile%tncdims%tdims(NMNHDIM_BUDGET_TIME), 'time axis for budgets', tzdates, tzdates_bound ) @@ -1740,31 +1746,34 @@ if ( tpfile%lmaster ) then call Write_time_coord( tpfile%tncdims%tdims(NMNHDIM_BUDGET_LES_TIME), 'time axis for LES budgets', tles_dates ) !Coordinates for the number of LES budget time averages - iavg = int( xles_temp_mean_end - 1.e-10 - xles_temp_mean_start ) / xles_temp_mean_step + 1 !Condition also on nles_times to not create this coordinate when not used (no time average if nles_times=0) - if ( nles_times > 0 .and. iavg > 0 ) then - Allocate( tzdates(iavg) ) - Allocate( tzdates_bound(2, iavg) ) - - do jt = 1, iavg - zles_temp_mean_start = xles_temp_mean_start + ( jt - 1 ) * xles_temp_mean_step - zles_temp_mean_end = Min( xles_temp_mean_end, zles_temp_mean_start + xles_temp_mean_step ) + if ( nles_times > 0 .and. nles_mean_times > 0 ) then + Allocate( tzdates(nles_mean_times) ) + Allocate( tzdates_bound(2, nles_mean_times) ) + + do jt = 1, nles_mean_times + jtb = ( nles_mean_start + ( jt - 1 ) * nles_mean_step ) / nles_dtcount + jte = MIN( jtb + nles_mean_step / nles_dtcount, nles_mean_end / nles_dtcount, nles_times ) + ! jtb could be 0 if nles_mean_start is smaller than the first LES measurement + ! For example, it occurs if xles_temp_mean_start is smaller than xles_temp_sampling (if xles_temp_mean_start=0.) + ! Do this correction only after computation of jte + if ( jtb < 1 ) jtb = 1 tzdates(jt)%nyear = tdtseg%nyear tzdates(jt)%nmonth = tdtseg%nmonth tzdates(jt)%nday = tdtseg%nday - tzdates(jt)%xtime = tdtseg%xtime + ( zles_temp_mean_start + zles_temp_mean_end ) / 2. + tzdates(jt)%xtime = tdtseg%xtime + ( xles_times(jtb) + xles_times(jte) ) / 2. !Not necessary: call Datetime_correctdate( tzdates(jt ) ) tzdates_bound(1, jt)%nyear = tdtseg%nyear tzdates_bound(1, jt)%nmonth = tdtseg%nmonth tzdates_bound(1, jt)%nday = tdtseg%nday - tzdates_bound(1, jt)%xtime = tdtseg%xtime + zles_temp_mean_start + tzdates_bound(1, jt)%xtime = tdtseg%xtime + xles_times(jtb) tzdates_bound(2, jt)%nyear = tdtseg%nyear tzdates_bound(2, jt)%nmonth = tdtseg%nmonth tzdates_bound(2, jt)%nday = tdtseg%nday - tzdates_bound(2, jt)%xtime = tdtseg%xtime + zles_temp_mean_end + tzdates_bound(2, jt)%xtime = tdtseg%xtime + xles_times(jte) end do call Write_time_coord( tpfile%tncdims%tdims(NMNHDIM_BUDGET_LES_AVG_TIME), 'time axis for LES budget time averages', & tzdates, tzdates_bound ) diff --git a/src/LIB/eccodes-2.25.0-Source.tar.gz b/src/LIB/eccodes-2.25.0-Source.tar.gz new file mode 100644 index 0000000000000000000000000000000000000000..bce1fa13227e80b843b193de77d21ecbd1e69f39 --- /dev/null +++ b/src/LIB/eccodes-2.25.0-Source.tar.gz @@ -0,0 +1,3 @@ +version https://git-lfs.github.com/spec/v1 +oid sha256:8975131aac54d406e5457706fd4e6ba46a8cc9c7dd817a41f2aa64ce1193c04e +size 12439952 diff --git a/src/LIB/megan.tar.gz b/src/LIB/megan.tar.gz deleted file mode 100644 index c9f186360fd3c305db25d9c72af24cca87988203..0000000000000000000000000000000000000000 --- a/src/LIB/megan.tar.gz +++ /dev/null @@ -1,3 +0,0 @@ -version https://git-lfs.github.com/spec/v1 -oid sha256:0d9386c64cbbdf2ffe32f1c52e51e2606869e136a95a65e21f8081961f8e48d3 -size 42865 diff --git a/src/MNH/aircraft_balloon.f90 b/src/MNH/aircraft_balloon.f90 index 153d76fe61b3174ea15d8e7546a8acf6cc4ddc74..ef6648cd6972e9caf618a7fa101eca6d44f82d66 100644 --- a/src/MNH/aircraft_balloon.f90 +++ b/src/MNH/aircraft_balloon.f90 @@ -34,7 +34,7 @@ REAL, DIMENSION(:,:), INTENT(IN) :: PTS ! surface temperature ! ++ OC REAL, DIMENSION(:,:,:), INTENT(IN) :: PRHODREF ! dry air density of the reference state REAL, DIMENSION(:,:,:), INTENT(IN) :: PCIT ! pristine ice concentration -REAL, DIMENSION(:,:),INTENT(IN) :: PSEA +REAL, DIMENSION(:,:),OPTIONAL,INTENT(IN) :: PSEA ! -- OC ! !------------------------------------------------------------------------------- diff --git a/src/MNH/aircraft_balloon_evol.f90 b/src/MNH/aircraft_balloon_evol.f90 index 1b4fcc77b3c459438d089866c5d19898daadb3cb..452113c95439e64ed6de56aec3ca1d3d1aa3bebe 100644 --- a/src/MNH/aircraft_balloon_evol.f90 +++ b/src/MNH/aircraft_balloon_evol.f90 @@ -37,7 +37,7 @@ REAL, DIMENSION(:,:,:), INTENT(IN) :: PRHODREF ! dry air density of the re REAL, DIMENSION(:,:,:), INTENT(IN) :: PCIT ! pristine ice concentration ! TYPE(FLYER), INTENT(INOUT) :: TPFLYER! balloon/aircraft -REAL, DIMENSION(:,:), INTENT(IN) :: PSEA +REAL, DIMENSION(:,:),OPTIONAL,INTENT(IN) :: PSEA ! !------------------------------------------------------------------------------- ! diff --git a/src/MNH/ch_update_jvalues.f90 b/src/MNH/ch_update_jvalues.f90 index f62d1d96eb657846610f5e8988627f9c417ee1c7..4c4754b7bf10fb1fa07d2c078b8687863e6a2fee 100644 --- a/src/MNH/ch_update_jvalues.f90 +++ b/src/MNH/ch_update_jvalues.f90 @@ -85,7 +85,9 @@ END MODULE MODI_CH_UPDATE_JVALUES !! ------------- !! Original 05/03/97 !! 05/03/05 P. Tulet (CNRM/GMEI) Update for Arome -! P. Wautelet 10/04/2019: replace ABORT and STOP calls by Print_msg +!! P. Wautelet 10/04/2019: replace ABORT and STOP calls by Print_msg +!! M. Leriche 01/06/2022: correct J by the presence of clouds +!! if nrr >= 2 and not nrr = 2 (meaning only cloud and water vapor) !! !!------------------------------------------------------------------------------ !! @@ -273,7 +275,7 @@ ELSE IF (.NOT.ALLOCATED(ZFCLD)) ALLOCATE(ZFCLD(1:IIU,1:IJU,1:IKU,JPJVMAX)) ! ZAZ3D(:,:,:) = PZZ(:,:,:) - IF (SIZE(PRT,4) == 2 ) THEN + IF (SIZE(PRT,4) >= 2 ) THEN ZLWC3D(:,:,:) = PRT(:,:,:,2) * PRHODREF(:,:,:) ELSE ZLWC3D(:,:,:) = 0. diff --git a/src/MNH/contrav.f90 b/src/MNH/contrav.f90 index 50e81468802efef8d210ae1736c96011e3bbbede..c890435f523160d94c0da938e847bf36c21b1da6 100644 --- a/src/MNH/contrav.f90 +++ b/src/MNH/contrav.f90 @@ -201,8 +201,20 @@ END IF ! ------------------------------------ ! IF (LFLAT) THEN - PRWCT(:,:,:) = PRWT(:,:,:) / PDZZ(:,:,:) - RETURN + ! + PRWCT(:,:,:) = PRWT(:,:,:) / PDZZ(:,:,:) + ! + IF (KADV_ORDER == 4 ) THEN + NULLIFY(TZFIELD_U) + NULLIFY(TZFIELD_V) + CALL ADD3DFIELD_ll( TZFIELD_U, PRUCT, 'CONTRAV::PRUCT' ) + CALL ADD3DFIELD_ll( TZFIELD_V, PRVCT, 'CONTRAV::PRVCT' ) + CALL UPDATE_HALO_ll(TZFIELD_U,IINFO_ll) + CALL UPDATE_HALO_ll(TZFIELD_V,IINFO_ll) + ENDIF + ! + RETURN + ! END IF ! !* 3. Compute the vertical contravariant components (general case) diff --git a/src/MNH/drag_veg.f90 b/src/MNH/drag_veg.f90 index 453bf94c2ae5ed4b865a02fa87e2180003eaf228..037369dc1d16a30cb4cd77e858131494727cd736 100644 --- a/src/MNH/drag_veg.f90 +++ b/src/MNH/drag_veg.f90 @@ -73,12 +73,14 @@ SUBROUTINE DRAG_VEG(PTSTEP,PUT,PVT,PTKET,ODEPOTREE, PVDEPOTREE, & ! P. Wautelet 28/01/2020: use the new data structures and subroutines for budgets for U ! C. Lac 02/2020: correction missing condition for budget on RC and SV ! P. Wautelet 04/02/2021: budgets: bugfixes for LDRAGTREE if LIMA + small optimisations and verifications +! R. Schoetter 04/2022: bug add update halo for vegetation drag variables !!--------------------------------------------------------------- ! ! !* 0. DECLARATIONS ! ------------ ! +USE MODD_ARGSLIST_ll, ONLY: LIST_ll use modd_budget, only: lbudget_u, lbudget_v, lbudget_rc, lbudget_sv, lbudget_tke, & NBUDGET_U, NBUDGET_V, NBUDGET_RC, NBUDGET_SV1, NBUDGET_TKE, & tbudgets @@ -96,6 +98,7 @@ USE MODD_VEG_n use mode_budget, only: Budget_store_init, Budget_store_end use mode_msg +USE MODE_ll USE MODI_MNHGET_SURF_PARAM_n USE MODI_SHUMAN @@ -125,8 +128,10 @@ REAL, DIMENSION(:,:,:,:), INTENT(INOUT) :: PSVS ! !* 0.2 Declarations of local variables : ! -INTEGER :: IIU,IJU,IKU,IKV ! array size along the k direction -INTEGER :: JI, JJ, JK ! loop index +INTEGER :: IIU,IJU,IKU,IKV ! array size along the k direction +INTEGER :: JI, JJ, JK ! loop index +INTEGER :: IINFO_ll +TYPE(LIST_ll), POINTER :: TZFIELDS_ll ! list of fields to exchange ! ! REAL, DIMENSION(SIZE(PUT,1),SIZE(PUT,2),SIZE(PUT,3)) :: & @@ -186,6 +191,16 @@ WHERE ( ZLAI (:,:) > (XUNDEF-1.) ) ZLAI (:,:) = 0.0 ZUT_SCAL(:,:,:) = MXF(PUT(:,:,:)) ZVT_SCAL(:,:,:) = MYF(PVT(:,:,:)) ZTKET(:,:,:) = PTKET(:,:,:) +! +! Update halo +! +NULLIFY(TZFIELDS_ll) +CALL ADD3DFIELD_ll( TZFIELDS_ll, ZUT_SCAL, 'DRAG_VEG::ZUT_SCAL') +CALL ADD3DFIELD_ll( TZFIELDS_ll, ZVT_SCAL, 'DRAG_VEG::ZVT_SCAL') +CALL ADD3DFIELD_ll( TZFIELDS_ll, ZTKET , 'DRAG_VEG::ZTKET' ) +CALL UPDATE_HALO_ll(TZFIELDS_ll,IINFO_ll) +CALL CLEANLIST_ll(TZFIELDS_ll) +! !------------------------------------------------------------------------------- ! !* 1. Computations of wind tendency due to canopy drag @@ -243,6 +258,15 @@ ENDDO ! To exclude the first vertical level already dealt in rain_ice or rain_c2r2_khko GDEP(:,:,2) = .FALSE. ! +! Update halo +! +NULLIFY(TZFIELDS_ll) +CALL ADD3DFIELD_ll( TZFIELDS_ll, ZCDRAG , 'DRAG_VEG::ZCDRAG') +CALL ADD3DFIELD_ll( TZFIELDS_ll, ZDENSITY, 'DRAG_VEG::ZDENSITY') +CALL UPDATE_HALO_ll(TZFIELDS_ll,IINFO_ll) +CALL CLEANLIST_ll(TZFIELDS_ll) +! +! !* 1.2 Drag force by wall surfaces ! --------------------------- ! diff --git a/src/MNH/ground_paramn.f90 b/src/MNH/ground_paramn.f90 index 876a976d4e3db0888894c4289d51aa448205b619..c438ead7fe900e3a1748ab41a9ca0656ccee065b 100644 --- a/src/MNH/ground_paramn.f90 +++ b/src/MNH/ground_paramn.f90 @@ -111,6 +111,7 @@ END MODULE MODI_GROUND_PARAM_n !! (V. Vionnet) 18/07/2017 add coupling for blowing snow module !! (Bielli S.) 02/2019 Sea salt : significant sea wave height influences salt emission; 5 salt modes ! P. Wautelet 20/05/2019: add name argument to ADDnFIELD_ll + new ADD4DFIELD_ll subroutine +! P. Wautelet 09/02/2022: bugfix: add missing XCURRENT_LEI computation !------------------------------------------------------------------------------- ! !* 0. DECLARATIONS @@ -122,7 +123,8 @@ USE MODI_GET_HALO USE MODI_MNH_OASIS_RECV USE MODI_MNH_OASIS_SEND USE MODD_SFX_OASIS, ONLY : LOASIS -USE MODD_DYN, ONLY : XSEGLEN +USE MODD_DYN, ONLY : XSEGLEN +USE MODD_DYN_n, ONLY : DYN_MODEL #endif ! USE MODD_LUNIT_n, ONLY: TLUOUT @@ -336,7 +338,8 @@ REAL, DIMENSION(:), ALLOCATABLE :: ZP_PET_B_COEF REAL, DIMENSION(:), ALLOCATABLE :: ZP_PEQ_B_COEF REAL, DIMENSION(:), ALLOCATABLE :: ZP_RN ! net radiation (W/m2) REAL, DIMENSION(:), ALLOCATABLE :: ZP_H ! sensible heat flux (W/m2) -REAL, DIMENSION(:), ALLOCATABLE :: ZP_LE ! latent heat flux (W/m2) +REAL, DIMENSION(:), ALLOCATABLE :: ZP_LE ! Total latent heat flux (W/m2) +REAL, DIMENSION(:), ALLOCATABLE :: ZP_LEI ! Solid Latent heat flux (W/m2) REAL, DIMENSION(:), ALLOCATABLE :: ZP_GFLUX ! ground flux (W/m2) REAL, DIMENSION(:), ALLOCATABLE :: ZP_T2M ! Air temperature at 2 meters (K) REAL, DIMENSION(:), ALLOCATABLE :: ZP_Q2M ! Air humidity at 2 meters (kg/kg) @@ -574,7 +577,7 @@ CALL DATETIME_DISTANCE(TDTSEG,TDTCUR,ZTIMEC) #ifdef CPLOASIS IF (LOASIS) THEN IF ( MOD(ZTIMEC,1.0) .LE. 1E-2 .OR. (1.0 - MOD(ZTIMEC,1.0)) .LE. 1E-2 ) THEN - IF ( NINT(ZTIMEC-(XSEGLEN-XTSTEP)) .LT. 0 ) THEN + IF ( NINT(ZTIMEC-(XSEGLEN-DYN_MODEL(1)%XTSTEP)) .LT. 0 ) THEN WRITE(ILUOUT,*) '----------------------------' WRITE(ILUOUT,*) ' Reception des champs avec OASIS' WRITE(ILUOUT,*) 'NINT(ZTIMEC)=', NINT(ZTIMEC) @@ -605,7 +608,7 @@ CALL COUPLING_SURF_ATM_n(YSURF_CUR,'MESONH', 'E',ZTIMEC, #ifdef CPLOASIS IF (LOASIS) THEN IF ( MOD(ZTIMEC,1.0) .LE. 1E-2 .OR. (1.0 - MOD(ZTIMEC,1.0)) .LE. 1E-2 ) THEN - IF (NINT(ZTIMEC-(XSEGLEN-XTSTEP)) .LT. 0) THEN + IF (NINT(ZTIMEC-(XSEGLEN-DYN_MODEL(1)%XTSTEP)) .LT. 0) THEN WRITE(ILUOUT,*) '----------------------------' WRITE(ILUOUT,*) ' Envoi des champs avec OASIS' WRITE(ILUOUT,*) 'NINT(ZTIMEC)=', NINT(ZTIMEC) @@ -618,9 +621,9 @@ END IF ! IF (CPROGRAM=='DIAG ' .OR. LDIAG_IN_RUN) THEN CALL DIAG_SURF_ATM_n(YSURF_CUR,'MESONH') - CALL MNHGET_SURF_PARAM_n(PRN=ZP_RN,PH=ZP_H,PLE=ZP_LE,PGFLUX=ZP_GFLUX, & - PT2M=ZP_T2M,PQ2M=ZP_Q2M,PHU2M=ZP_HU2M, & - PZON10M=ZP_ZON10M,PMER10M=ZP_MER10M ) + CALL MNHGET_SURF_PARAM_n( PRN = ZP_RN, PH = ZP_H, PLE = ZP_LE, PLEI = ZP_LEI, & + PGFLUX = ZP_GFLUX, PT2M = ZP_T2M, PQ2M = ZP_Q2M, PHU2M = ZP_HU2M, & + PZON10M = ZP_ZON10M, PMER10M = ZP_MER10M ) END IF ! ! Transform 1D output fields into 2D: @@ -843,6 +846,7 @@ ALLOCATE(ZP_QSURF (KDIM1D)) ALLOCATE(ZP_RN (KDIM1D)) ALLOCATE(ZP_H (KDIM1D)) ALLOCATE(ZP_LE (KDIM1D)) +ALLOCATE(ZP_LEI (KDIM1D)) ALLOCATE(ZP_GFLUX (KDIM1D)) ALLOCATE(ZP_T2M (KDIM1D)) ALLOCATE(ZP_Q2M (KDIM1D)) @@ -966,6 +970,7 @@ IF (LDIAG_IN_RUN) THEN XCURRENT_RN (IIB:IIE,IJB:IJE) = RESHAPE(ZP_RN(:), ISHAPE_2) XCURRENT_H (IIB:IIE,IJB:IJE) = RESHAPE(ZP_H (:), ISHAPE_2) XCURRENT_LE (IIB:IIE,IJB:IJE) = RESHAPE(ZP_LE(:), ISHAPE_2) + XCURRENT_LEI (IIB:IIE,IJB:IJE) = RESHAPE(ZP_LEI(:), ISHAPE_2) XCURRENT_GFLUX (IIB:IIE,IJB:IJE) = RESHAPE(ZP_GFLUX(:), ISHAPE_2) XCURRENT_T2M (IIB:IIE,IJB:IJE) = RESHAPE(ZP_T2M(:), ISHAPE_2) XCURRENT_Q2M (IIB:IIE,IJB:IJE) = RESHAPE(ZP_Q2M(:), ISHAPE_2) @@ -1014,6 +1019,7 @@ DEALLOCATE(ZP_EMIS ) DEALLOCATE(ZP_RN ) DEALLOCATE(ZP_H ) DEALLOCATE(ZP_LE ) +DEALLOCATE(ZP_LEI ) DEALLOCATE(ZP_GFLUX ) DEALLOCATE(ZP_T2M ) DEALLOCATE(ZP_Q2M ) diff --git a/src/MNH/ice4_tendencies.f90 b/src/MNH/ice4_tendencies.f90 index 2832464bf64ed83229d3779531ef93c206db078d..a80e037aab72dc9260224f09ad5d4453957573b3 100644 --- a/src/MNH/ice4_tendencies.f90 +++ b/src/MNH/ice4_tendencies.f90 @@ -467,9 +467,12 @@ IF(KSIZE>0) THEN DO JL=1,KSIZE ZRHT3D (K1(JL), K2(JL), K3(JL)) = ZRHT(JL) ENDDO + CALL ICE4_RAINFR_VERT(KIB, KIE, KIT, KJB, KJE, KJT, KKB, KKE, KKT, KKL, PRAINFR(:,:,:), & + &ZRRT3D(:,:,:), ZRST3D(:,:,:), ZRGT3D(:,:,:), ZRHT3D(:,:,:)) + ELSE + CALL ICE4_RAINFR_VERT(KIB, KIE, KIT, KJB, KJE, KJT, KKB, KKE, KKT, KKL, PRAINFR(:,:,:), & + &ZRRT3D(:,:,:), ZRST3D(:,:,:), ZRGT3D(:,:,:)) ENDIF - CALL ICE4_RAINFR_VERT(KIB, KIE, KIT, KJB, KJE, KJT, KKB, KKE, KKT, KKL, PRAINFR(:,:,:), ZRRT3D(:,:,:), & - &ZRST3D(:,:,:), ZRGT3D(:,:,:), ZRHT3D(:,:,:)) DO JL=1,KSIZE ZRF(JL)=PRAINFR(K1(JL), K2(JL), K3(JL)) END DO diff --git a/src/MNH/ini_aircraft_balloon.f90 b/src/MNH/ini_aircraft_balloon.f90 index 0b22d34031958913fa591c450709ede432d98301..cdb990bb73ae06e3440e85858b68161bc9806c92 100644 --- a/src/MNH/ini_aircraft_balloon.f90 +++ b/src/MNH/ini_aircraft_balloon.f90 @@ -1,4 +1,4 @@ -!MNH_LIC Copyright 2000-2020 CNRS, Meteo-France and Universite Paul Sabatier +!MNH_LIC Copyright 2000-2022 CNRS, Meteo-France and Universite Paul Sabatier !MNH_LIC This is part of the Meso-NH software governed by the CeCILL-C licence !MNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt !MNH_LIC for details. version 1. @@ -360,7 +360,7 @@ IF (IMI /= TPFLYER%NMODEL .AND. .NOT. (IMI==1 .AND. TPFLYER%NMODEL==0) ) RETURN IF ( CPROGRAM == 'DIAG ' ) THEN ISTORE = INT ( NTIME_AIRCRAFT_BALLOON / TPFLYER%STEP ) + 1 ELSE - ISTORE = INT ( (PSEGLEN-XTSTEP) / TPFLYER%STEP ) + 1 + ISTORE = NINT ( ( PSEGLEN - DYN_MODEL(1)%XTSTEP ) / TPFLYER%STEP ) + 1 ENDIF ! IF (TPFLYER%NMODEL == 0) ISTORE=0 diff --git a/src/MNH/ini_cturb.f90 b/src/MNH/ini_cturb.f90 index 245dfa0632b533d2855e56b9fc3459a9e51a48cd..ffe97a66684a46723380e1b3b03f9066c0099113 100644 --- a/src/MNH/ini_cturb.f90 +++ b/src/MNH/ini_cturb.f90 @@ -217,7 +217,7 @@ XCPR5= XCPR2 ! 3. MINIMUM VALUES ! -------------- ! -XTKEMIN=0.01 +XTKEMIN=0.01 ! This value is replaced by XKEMIN in &NAM_TURBn ! !XLINI=10. ! BL mixing length XLINI=0.1 ! BL mixing length diff --git a/src/MNH/ini_lesn.f90 b/src/MNH/ini_lesn.f90 index 406ee30a2feac649f7b708aa132f538ed52288b0..ff089c8c480444dd13b575e408294467283ff8c4 100644 --- a/src/MNH/ini_lesn.f90 +++ b/src/MNH/ini_lesn.f90 @@ -1,4 +1,4 @@ -!MNH_LIC Copyright 2000-2021 CNRS, Meteo-France and Universite Paul Sabatier +!MNH_LIC Copyright 2000-2022 CNRS, Meteo-France and Universite Paul Sabatier !MNH_LIC This is part of the Meso-NH software governed by the CeCILL-C licence !MNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt !MNH_LIC for details. version 1. @@ -40,6 +40,7 @@ ! P. Wautelet 04/01/2021: bugfix: nles_k was used instead of nspectra_k for a loop index ! P. Wautelet 30/03/2021: budgets: LES cartesian subdomain limits are defined in the physical domain ! P. Wautelet 09/07/2021: bugfix: altitude levels are on the correct grid position (mass point) +! P. Wautelet 22/03/2022: LES averaging periods are more reliable (compute with integers instead of reals) ! -------------------------------------------------------------------------- ! !* 0. DECLARATIONS @@ -321,7 +322,7 @@ XLES_TEMP_SAMPLING = XTSTEP * NLES_DTCOUNT ! ---------------------------------------- ! ! -NLES_TIMES = ( INT( (XSEGLEN-XTSTEP+1.E-6) / XTSTEP ) ) / NLES_DTCOUNT +NLES_TIMES = ( NINT( ( XSEGLEN - DYN_MODEL(1)%XTSTEP ) / XTSTEP ) ) / NLES_DTCOUNT ! !* 3.5 current LES time counter ! ------------------------ @@ -342,6 +343,46 @@ IF (NLES_TIMES==0) THEN RETURN END IF ! +!* 3.8 Averaging +! --------- +IF ( XLES_TEMP_MEAN_END == XUNDEF & + .OR. XLES_TEMP_MEAN_START == XUNDEF & + .OR. XLES_TEMP_MEAN_STEP == XUNDEF ) THEN + !No LES temporal averaging + NLES_MEAN_TIMES = 0 + NLES_MEAN_STEP = NNEGUNDEF + NLES_MEAN_START = NNEGUNDEF + NLES_MEAN_END = NNEGUNDEF +ELSE + !LES temporal averaging is enabled + !Ensure that XLES_TEMP_MEAN_END is not after segment end + XLES_TEMP_MEAN_END = MIN( XLES_TEMP_MEAN_END, XSEGLEN - DYN_MODEL(1)%XTSTEP ) + + NLES_MEAN_START = NINT( XLES_TEMP_MEAN_START / XTSTEP ) + + IF ( MODULO( NLES_MEAN_START, NLES_DTCOUNT ) /= 0 ) THEN + CMNHMSG(1) = 'XLES_TEMP_MEAN_START is not a multiple of XLES_TEMP_SAMPLING' + CMNHMSG(2) = 'LES averaging periods could be wrong' + CALL Print_msg( NVERB_WARNING, 'IO', 'INI_LES_n' ) + END IF + + NLES_MEAN_END = NINT( XLES_TEMP_MEAN_END / XTSTEP ) + + NLES_MEAN_STEP = NINT( XLES_TEMP_MEAN_STEP / XTSTEP ) + + IF ( NLES_MEAN_STEP < NLES_DTCOUNT ) & + CALL Print_msg( NVERB_ERROR, 'IO', 'INI_LES_n', 'XLES_TEMP_MEAN_STEP < XLES_TEMP_SAMPLING not allowed' ) + + IF ( MODULO( NLES_MEAN_STEP, NLES_DTCOUNT ) /= 0 ) THEN + CMNHMSG(1) = 'XLES_TEMP_MEAN_STEP is not a multiple of XLES_TEMP_SAMPLING' + CMNHMSG(2) = 'LES averaging periods could be wrong' + CALL Print_msg( NVERB_WARNING, 'IO', 'INI_LES_n' ) + END IF + + NLES_MEAN_TIMES = ( NLES_MEAN_END - NLES_MEAN_START ) / NLES_MEAN_STEP + !Add 1 averaging period if the last one is incomplete (for example: start=0., end=10., step=3.) + IF ( MODULO( NLES_MEAN_END - NLES_MEAN_START, NLES_MEAN_STEP ) > 0 ) NLES_MEAN_TIMES = NLES_MEAN_TIMES + 1 +END IF !------------------------------------------------------------------------------- ! !* 4. Number of vertical levels for local diagnostics diff --git a/src/MNH/ini_modeln.f90 b/src/MNH/ini_modeln.f90 index 911e543fc9d4be37e4f44207574d2b64831cd270..42447401a055a7354fec871864ee0e837084522b 100644 --- a/src/MNH/ini_modeln.f90 +++ b/src/MNH/ini_modeln.f90 @@ -1799,7 +1799,7 @@ gles = lles_mean .or. lles_resolved .or. lles_subgrid .or. lles_updraft & .or. lles_downdraft .or. lles_spectra !Called if budgets are enabled via NAM_BUDGET !or if LES budgets are enabled via NAM_LES (condition on kmi==1 to call it max once) -if ( ( cbutype /= "NONE" .and. nbumod == kmi ) .or. ( gles .and. kmi == 1 ) .or. LCHECK ) THEN +if ( ( cbutype /= "NONE" .and. nbumod == kmi ) .or. ( ( gles .or. lcheck ) .and. kmi == 1 ) ) THEN call Budget_preallocate() end if diff --git a/src/MNH/ini_posprofilern.f90 b/src/MNH/ini_posprofilern.f90 index 08b25fba13d2dfdc78cdaf669b31ee81a1bdd62f..00279c8a81488360ebcc1e2ce4bb677094e55361 100644 --- a/src/MNH/ini_posprofilern.f90 +++ b/src/MNH/ini_posprofilern.f90 @@ -1,4 +1,4 @@ -!MNH_LIC Copyright 1994-2021 CNRS, Meteo-France and Universite Paul Sabatier +!MNH_LIC Copyright 1994-2022 CNRS, Meteo-France and Universite Paul Sabatier !MNH_LIC This is part of the Meso-NH software governed by the CeCILL-C licence !MNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt !MNH_LIC for details. version 1. @@ -156,7 +156,7 @@ SUBROUTINE ALLOCATE_PROFILER_n(TPROFILER) ! TYPE(PROFILER), INTENT(INOUT) :: TPROFILER ! -ISTORE = INT ( (PSEGLEN-XTSTEP) / TPROFILER%STEP ) + 1 +ISTORE = NINT( ( PSEGLEN - DYN_MODEL(1)%XTSTEP ) / TPROFILER%STEP ) + 1 ! allocate( tprofiler%tpdates( istore ) ) ALLOCATE(TPROFILER%ERROR (NUMBPROFILER)) diff --git a/src/MNH/ini_segn.f90 b/src/MNH/ini_segn.f90 index 590efa55c1c37baaab8c9a2ccf4e5eef8c4b2ad4..fe1fe7a55f86b4a023a8096e9fc7522616d2588b 100644 --- a/src/MNH/ini_segn.f90 +++ b/src/MNH/ini_segn.f90 @@ -1,4 +1,4 @@ -!MNH_LIC Copyright 1994-2021 CNRS, Meteo-France and Universite Paul Sabatier +!MNH_LIC Copyright 1994-2022 CNRS, Meteo-France and Universite Paul Sabatier !MNH_LIC This is part of the Meso-NH software governed by the CeCILL-C licence !MNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt !MNH_LIC for details. version 1. @@ -455,7 +455,7 @@ CALL READ_EXSEG_n(KMI,TZFILE_DES,YCONF,GFLAT,GUSERV,GUSERC, & GLNOX_EXPLICIT, & GCONDSAMP,GBLOWSNOW, IRIMX,IRIMY,ISV, & YTURB,YTOM,GRMC01,YRAD,YDCONV,YSCONV,YCLOUD,YELEC,YEQNSYS, & - PTSTEP_ALL,CSTORAGE_TYPE,CINIFILEPGD_n ) + PTSTEP_ALL,CINIFILEPGD_n ) ! if ( cprogram == 'MESONH' .and. kmi == 1 ) then !Do this only once call Fieldlist_nmodel_resize(NMODEL) diff --git a/src/MNH/ini_surfstationn.f90 b/src/MNH/ini_surfstationn.f90 index f53ee35d11e73fd6873e8e7abf40d90a227579af..52c9df797e2008d15cd424e3f9e5c14140ab6036 100644 --- a/src/MNH/ini_surfstationn.f90 +++ b/src/MNH/ini_surfstationn.f90 @@ -1,4 +1,4 @@ -!MNH_LIC Copyright 1994-2021 CNRS, Meteo-France and Universite Paul Sabatier +!MNH_LIC Copyright 1994-2022 CNRS, Meteo-France and Universite Paul Sabatier !MNH_LIC This is part of the Meso-NH software governed by the CeCILL-C licence !MNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt !MNH_LIC for details. version 1. @@ -171,16 +171,12 @@ if ( tstation%step < xtstep ) then tstation%step = xtstep end if -IF (KMI==1) THEN - ISTORE = NINT ( (PSEGLEN-XTSTEP) / TSTATION%STEP ) + 1 -ELSE - ISTORE = NINT ( (PSEGLEN-XTSTEP * NDTRATIO(KMI)) / TSTATION%STEP ) + 1 -END IF +ISTORE = NINT ( ( PSEGLEN - DYN_MODEL(1)%XTSTEP ) / TSTATION%STEP ) + 1 allocate( tstation%tpdates( istore ) ) ALLOCATE(TSTATION%ERROR (NUMBSTAT)) -!ALLOCATE(TSTATION%X (NUMBSTAT)) -!ALLOCATE(TSTATION%Y (NUMBSTAT)) +ALLOCATE(TSTATION%X (NUMBSTAT)) +ALLOCATE(TSTATION%Y (NUMBSTAT)) ALLOCATE(TSTATION%SV (ISTORE,NUMBSTAT,KSV)) ALLOCATE(TSTATION%TSRAD (ISTORE,NUMBSTAT)) ALLOCATE(TSTATION%ZS (NUMBSTAT)) diff --git a/src/MNH/initial_guess.f90 b/src/MNH/initial_guess.f90 index 004d6be13168a80d05113e16940d8d9283283f72..cba1058bc1b554620e674215843be01b2218f9ae 100644 --- a/src/MNH/initial_guess.f90 +++ b/src/MNH/initial_guess.f90 @@ -153,6 +153,7 @@ use modd_budget, only: lbudget_u, lbudget_v, lbudget_w, lbudget_th, lbudg NBUDGET_RR, NBUDGET_RI, NBUDGET_RS, NBUDGET_RG, NBUDGET_RH, NBUDGET_SV1, & lbu_beg, lbu_enable, tbudgets USE MODD_CONF +use modd_conf_n, only: luserv USE MODD_GRID_n use mode_budget, only: Budget_store_init, Budget_store_end @@ -213,11 +214,14 @@ IF (SIZE(PTKET,1) /= 0) THEN PRTKES(:,:,:) = PTKET(:,:,:) * ZINVTSTEP * PRHODJ(:,:,:) END IF ! -! Case with KRR moist variables -DO JRR=1,KRR - PRRS(:,:,:,JRR) = PRT(:,:,:,JRR) * ZINVTSTEP * PRHODJ(:,:,:) -END DO -CALL MPPDB_CHECK3DM("initial_guess:PRRS/RT/RHO",PRECISION,PRRS(:,:,:,1) , PRT(:,:,:,1) , PRHODJ) +! Case with KRR moist variables +! +IF (LUSERV) THEN + DO JRR=1,KRR + PRRS(:,:,:,JRR) = PRT(:,:,:,JRR) * ZINVTSTEP * PRHODJ(:,:,:) + END DO + CALL MPPDB_CHECK3DM("initial_guess:PRRS/RT/RHO",PRECISION,PRRS(:,:,:,1) , PRT(:,:,:,1) , PRHODJ) +ENDIF ! ! *** passive tracers ! diff --git a/src/MNH/isocom.f b/src/MNH/isocom.f index d92460fe1ac19f237be0ac2bd6c7a9fa2c79deab..42481f83a3f489686c70fa8a696f9ffe8966ff04 100644 --- a/src/MNH/isocom.f +++ b/src/MNH/isocom.f @@ -1,4 +1,4 @@ -CMNH_LIC Copyright 1996-2021 CNRS, Meteo-France and Universite Paul Sabatier +CMNH_LIC Copyright 1996-2022 CNRS, Meteo-France and Universite Paul Sabatier CMNH_LIC This is part of the Meso-NH software governed by the CeCILL-C licence CMNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt CMNH_LIC for details. version 1. diff --git a/src/MNH/mesonh.f90 b/src/MNH/mesonh.f90 index ac5f3cdef04db6c7ffa837bbbd1f849b8974e9a8..dd28504bb919d89831b3b9b1c56e8cc63828331d 100644 --- a/src/MNH/mesonh.f90 +++ b/src/MNH/mesonh.f90 @@ -134,7 +134,7 @@ INTEGER :: IINFO_ll ! return code of // routines ! Switch to model 1 variables #ifndef CPLOASIS CALL MPPDB_INIT() -CALL MPPDB_STOP_DEBUG() +!CALL MPPDB_STOP_DEBUG() #endif ! CALL GOTO_MODEL(1,ONOFIELDLIST=.TRUE.) diff --git a/src/MNH/modd_lesn.f90 b/src/MNH/modd_lesn.f90 index ac78ef503a2edbb805192fd505df5d174686b7ad..2d2009b060cd4742732796d94958bd006eb20704 100644 --- a/src/MNH/modd_lesn.f90 +++ b/src/MNH/modd_lesn.f90 @@ -1,4 +1,4 @@ -!MNH_LIC Copyright 1995-2020 CNRS, Meteo-France and Universite Paul Sabatier +!MNH_LIC Copyright 1995-2022 CNRS, Meteo-France and Universite Paul Sabatier !MNH_LIC This is part of the Meso-NH software governed by the CeCILL-C licence !MNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt !MNH_LIC for details. version 1. @@ -42,6 +42,7 @@ ! P. Wautelet 08/02/2019: add missing NULL association for pointers ! C. Lac 02/2019: add rain fraction as a LES diagnostic ! P. Wautelet 13/09/2019: budget: simplify and modernize date/time management +! P. Wautelet 22/03/2022: LES averaging periods are more reliable (compute with integers instead of reals) !------------------------------------------------------------------------------- ! !* 0. DECLARATIONS @@ -61,6 +62,10 @@ TYPE LES_t INTEGER :: NLES_TIMES ! number of LES computations in time INTEGER :: NLES_DTCOUNT ! number of time steps between two LES comp. INTEGER :: NLES_TCOUNT ! current time counter for LES comp. + INTEGER :: NLES_MEAN_TIMES ! Number of LES averaging periods + INTEGER :: NLES_MEAN_STEP ! number of time steps between two LES average comp. + INTEGER :: NLES_MEAN_START ! First time step number taken into account for LES averaging + INTEGER :: NLES_MEAN_END ! Last time step number taken into account for LES averaging INTEGER :: NSPECTRA_NI ! number of wave lengths in I direction INTEGER :: NSPECTRA_NJ ! number of wave lengths in J direction ! @@ -661,6 +666,10 @@ TYPE(LES_t), DIMENSION(JPMODELMAX), TARGET, SAVE :: LES_MODEL INTEGER, POINTER :: NLES_TIMES=>NULL() INTEGER, POINTER :: NLES_DTCOUNT=>NULL() INTEGER, POINTER :: NLES_TCOUNT=>NULL() +INTEGER, POINTER :: NLES_MEAN_TIMES => NULL() +INTEGER, POINTER :: NLES_MEAN_STEP => NULL() +INTEGER, POINTER :: NLES_MEAN_START => NULL() +INTEGER, POINTER :: NLES_MEAN_END => NULL() INTEGER, POINTER :: NSPECTRA_NI=>NULL() INTEGER, POINTER :: NSPECTRA_NJ=>NULL() type(date_time), dimension(:), pointer :: tles_dates => null() @@ -1512,6 +1521,10 @@ LES_MODEL(KFROM)%XLES_RADEFF=>XLES_RADEFF NLES_TIMES=>LES_MODEL(KTO)%NLES_TIMES NLES_DTCOUNT=>LES_MODEL(KTO)%NLES_DTCOUNT NLES_TCOUNT=>LES_MODEL(KTO)%NLES_TCOUNT +NLES_MEAN_TIMES => LES_MODEL(KTO)%NLES_MEAN_TIMES +NLES_MEAN_STEP => LES_MODEL(KTO)%NLES_MEAN_STEP +NLES_MEAN_START => LES_MODEL(KTO)%NLES_MEAN_START +NLES_MEAN_END => LES_MODEL(KTO)%NLES_MEAN_END NSPECTRA_NI=>LES_MODEL(KTO)%NSPECTRA_NI NSPECTRA_NJ=>LES_MODEL(KTO)%NSPECTRA_NJ tles_dates=>les_model(kto)%tles_dates diff --git a/src/MNH/mode_fscatter.f90 b/src/MNH/mode_fscatter.f90 index 6c5bdbb670eba3e73137a8e01323136edfa1b4f9..3e66e3a4394b708b49acc9dbe6c7c4d65a37ef1c 100644 --- a/src/MNH/mode_fscatter.f90 +++ b/src/MNH/mode_fscatter.f90 @@ -153,6 +153,12 @@ CONTAINS ! ! Modification (C.Lac) 04/2014 : exclude very small values of x ! + IF (X <= 1.E-07) THEN + QEXT = 0. + QSCA = 0. + QBACK = 0. + RETURN + ELSE ! ----------------------------------------------------------- ! del is the inner sphere convergence criterion ! ----------------------------------------------------------- @@ -169,13 +175,7 @@ CONTAINS qext = 0.0 xback = (0.0,0.0) n=1 - IF (x <= 1.E-07) THEN - QEXT = 0. - QSCA = 0. - QBACK = 0. - RETURN - ELSE - do while(n<=nstop) + do while(n<=nstop) ! DO n=1,nstop DX = 1.0/(n/x-dx) - n/x DY = 1.0/(n/y-dy) - n/y diff --git a/src/MNH/mode_les_diachro.f90 b/src/MNH/mode_les_diachro.f90 index f56b2ce0f607015dde207609bd257cbdcf1050da..7d8d698722653f0bb5255fcbe99e88a33992bf0e 100644 --- a/src/MNH/mode_les_diachro.f90 +++ b/src/MNH/mode_les_diachro.f90 @@ -1,4 +1,4 @@ -!MNH_LIC Copyright 1994-2021 CNRS, Meteo-France and Universite Paul Sabatier +!MNH_LIC Copyright 1994-2022 CNRS, Meteo-France and Universite Paul Sabatier !MNH_LIC This is part of the Meso-NH software governed by the CeCILL-C licence !MNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt !MNH_LIC for details. version 1. @@ -12,6 +12,7 @@ ! P. Wautelet 10/2020: restructure subroutines to use tfield_metadata_base type ! P. Wautelet 03/03/2021: budgets: add tbudiachrometadata type (useful to pass more information to Write_diachro) ! P. Wautelet 11/03/2021: budgets: remove ptrajx/y/z optional dummy arguments of Write_diachro +! P. Wautelet 22/03/2022: LES averaging periods are more reliable (compute with integers instead of reals) !----------------------------------------------------------------- !####################### MODULE MODE_LES_DIACHRO @@ -524,8 +525,7 @@ END SUBROUTINE LES_TIME_AVG subroutine Les_time_avg_4D( pwork4, tpdates, kresp ) !######################################################## -use modd_les, only: nles_current_times, xles_temp_mean_start, xles_temp_mean_end, xles_temp_mean_step -use modd_parameters, only: XUNDEF +use modd_les_n, only: nles_dtcount, nles_mean_start, nles_mean_end, nles_mean_step, nles_mean_times, nles_times use modd_time, only: tdtseg use modd_type_date, only: date_time @@ -537,60 +537,33 @@ real, dimension(:,:,:,:), allocatable, intent(inout) :: pwork4 ! cont type(date_time), dimension(:), allocatable, intent(inout) :: tpdates integer, intent(out) :: kresp ! return code (0 is ok) !------------------------------------------------------------------------------ -integer :: jt ! time counter -integer :: itime ! nb of avg. points -integer :: iavg ! nb of avg. periods integer :: javg ! loop counter on avg. periods integer :: jk ! vertical loop counter integer :: jp ! process loop counter integer :: jsv ! scalar loop counter -integer :: jx ! first spatial or spectral coordinate loop counter -integer :: jy ! second spatial or spectral coordinate loop counter integer :: jtb, jte -real :: zles_temp_mean_start ! initial and end times -real :: zles_temp_mean_end ! of one avergaing preiod real, dimension(:,:,:,:), allocatable :: zwork4 ! contains averaged physical field !------------------------------------------------------------------------------ -if ( xles_temp_mean_end == XUNDEF & - .or. xles_temp_mean_start == XUNDEF & - .or. xles_temp_mean_step == XUNDEF ) then - kresp = -1 - return -end if - -iavg = Int( xles_temp_mean_end - 1.e-10 - xles_temp_mean_start ) / xles_temp_mean_step + 1 -if ( iavg <= 0 ) then +if ( nles_mean_times == 0 ) then kresp = -1 return end if Deallocate( tpdates ) -Allocate( tpdates(iavg) ) -Allocate( zwork4(Size( pwork4, 1 ), iavg, Size( pwork4, 3 ), Size( pwork4, 4 )) ) +Allocate( tpdates(nles_mean_times) ) +Allocate( zwork4(Size( pwork4, 1 ), nles_mean_times, Size( pwork4, 3 ), Size( pwork4, 4 )) ) zwork4(:, :, :, :) = 0. -do javg = 1, iavg - zles_temp_mean_start = xles_temp_mean_start + (javg - 1) * xles_temp_mean_step - zles_temp_mean_end = Min( xles_temp_mean_end, zles_temp_mean_start + xles_temp_mean_step ) - - jtb = -1 - jte = -2 - do jt = 1, nles_current_times - if ( xles_times(jt) >= zles_temp_mean_start ) then - jtb = jt - exit - end if - end do - do jt = jtb, nles_current_times - if ( xles_times(jt) <= zles_temp_mean_end ) then - jte = jt - else - exit - end if - end do +do javg = 1, nles_mean_times + jtb = ( nles_mean_start + ( javg - 1 ) * nles_mean_step ) / nles_dtcount + jte = MIN( jtb + nles_mean_step / nles_dtcount, nles_mean_end / nles_dtcount, nles_times ) + ! jtb could be 0 if nles_mean_start is smaller than the first LES measurement + ! For example, it occurs if xles_temp_mean_start is smaller than xles_temp_sampling (if xles_temp_mean_start=0.) + ! Do this correction only after computation of jte + if ( jtb < 1 ) jtb = 1 do jp = 1, Size( pwork4, 4 ) do jsv = 1, Size( pwork4, 3 ) @@ -603,12 +576,12 @@ do javg = 1, iavg tpdates(javg)%nyear = tdtseg%nyear tpdates(javg)%nmonth = tdtseg%nmonth tpdates(javg)%nday = tdtseg%nday - tpdates(javg)%xtime = tdtseg%xtime + ( zles_temp_mean_start + zles_temp_mean_end ) / 2. + tpdates(javg)%xtime = tdtseg%xtime + ( xles_times(jtb) + xles_times(jte) ) / 2. call Datetime_correctdate( tpdates(javg) ) end do Deallocate( pwork4 ) -Allocate( pwork4(Size( zwork4, 1 ), iavg, Size( zwork4, 3 ), Size( zwork4, 4 )) ) +Allocate( pwork4(Size( zwork4, 1 ), nles_mean_times, Size( zwork4, 3 ), Size( zwork4, 4 )) ) pwork4 = zwork4 diff --git a/src/MNH/modeln.f90 b/src/MNH/modeln.f90 index 8672e488dd9daf956d657273bf10b4622df20f5a..1e183f9a4a6ff16e98568d5084adc43b94316561 100644 --- a/src/MNH/modeln.f90 +++ b/src/MNH/modeln.f90 @@ -1906,6 +1906,7 @@ IF (CCLOUD /= 'NONE' .AND. CELEC == 'NONE') THEN XHLC_HRC, XHLC_HCF, XHLI_HRI, XHLI_HCF, & ZSEA, ZTOWN ) DEALLOCATE(ZTOWN) + DEALLOCATE(ZSEA) ELSE CALL RESOLVED_CLOUD ( CCLOUD, CACTCCN, CSCONV, CMF_CLOUD, NRR, NSPLITR, & NSPLITG, IMI, KTCOUNT, & @@ -1998,6 +1999,7 @@ IF (CELEC /= 'NONE' .AND. (CCLOUD(1:3) == 'ICE')) THEN XINPRS, XINPRG, XINPRH, & ZSEA, ZTOWN ) DEALLOCATE(ZTOWN) + DEALLOCATE(ZSEA) ELSE CALL RESOLVED_ELEC_n (CCLOUD, CSCONV, CMF_CLOUD, & NRR, NSPLITR, IMI, KTCOUNT, OEXIT, & @@ -2105,12 +2107,23 @@ XT_STEP_SWA = XT_STEP_SWA + ZTIME2 - ZTIME1 - XTIME_BU_PROCESS ! ZTIME1 = ZTIME2 ! -IF (LFLYER) & - CALL AIRCRAFT_BALLOON(XTSTEP, & +IF (LFLYER) THEN + IF (CSURF=='EXTE') THEN + ALLOCATE(ZSEA(IIU,IJU)) + ZSEA(:,:) = 0. + CALL MNHGET_SURF_PARAM_n (PSEA=ZSEA(:,:)) + CALL AIRCRAFT_BALLOON(XTSTEP, & XXHAT, XYHAT, XZZ, XMAP, XLONORI, XLATORI, & XUT, XVT, XWT, XPABST, XTHT, XRT, XSVT, XTKET, XTSRAD, & XRHODREF,XCIT,PSEA=ZSEA(:,:)) - + DEALLOCATE(ZSEA) + ELSE + CALL AIRCRAFT_BALLOON(XTSTEP, & + XXHAT, XYHAT, XZZ, XMAP, XLONORI, XLATORI, & + XUT, XVT, XWT, XPABST, XTHT, XRT, XSVT, XTKET, XTSRAD, & + XRHODREF,XCIT) + END IF +END IF !------------------------------------------------------------------------------- ! @@ -2127,11 +2140,23 @@ IF (LSTATION) & !* 24.3 PROFILER (observation diagnostic) ! --------------------------------- ! -IF (LPROFILER) & - CALL PROFILER_n(XTSTEP, & +IF (LPROFILER) THEN + IF (CSURF=='EXTE') THEN + ALLOCATE(ZSEA(IIU,IJU)) + ZSEA(:,:) = 0. + CALL MNHGET_SURF_PARAM_n (PSEA=ZSEA(:,:)) + CALL PROFILER_n(XTSTEP, & XXHAT, XYHAT, XZZ,XRHODREF, & XUT, XVT, XWT, XTHT, XRT, XSVT, XTKET, XTSRAD, XPABST, & XAER, MAX(XCLDFR,XICEFR), XCIT,PSEA=ZSEA(:,:)) + DEALLOCATE(ZSEA) + ELSE + CALL PROFILER_n(XTSTEP, & + XXHAT, XYHAT, XZZ,XRHODREF, & + XUT, XVT, XWT, XTHT, XRT, XSVT, XTKET, XTSRAD, XPABST, & + XAER, MAX(XCLDFR,XICEFR), XCIT) + END IF +END IF ! IF (ALLOCATED(ZSEA)) DEALLOCATE (ZSEA) ! diff --git a/src/MNH/modn_param_ice.f90 b/src/MNH/modn_param_ice.f90 index 085f74d7381939312012e03d5f1cf082779e233d..5f669897704cbd6fc0adee895d96c432b4b68327 100644 --- a/src/MNH/modn_param_ice.f90 +++ b/src/MNH/modn_param_ice.f90 @@ -26,6 +26,7 @@ NAMELIST/NAM_PARAM_ICE/LWARM,LSEDIC,LCONVHG,CPRISTINE_ICE,CSEDIM,LDEPOSC,XVDEPOS LEVLIMIT,LNULLWETG,LWETGPOST,LNULLWETH,LWETHPOST, & CSNOWRIMING,XFRACM90,NMAXITER,XMRSTEP,XTSTEP_TS, & LADJ_BEFORE, LADJ_AFTER, CFRAC_ICE_ADJUST, LCRFLIMIT, & - XSPLIT_MAXCFL, CFRAC_ICE_SHALLOW_MF, LSEDIM_AFTER, LSNOW_T + XSPLIT_MAXCFL, CFRAC_ICE_SHALLOW_MF, LSEDIM_AFTER, & + CSUBG_RC_RR_ACCR, CSUBG_RR_EVAP, CSUBG_PR_PDF, LSNOW_T ! END MODULE MODN_PARAM_ICE diff --git a/src/MNH/profilern.f90 b/src/MNH/profilern.f90 index c9638ee4b3b035314f7605bb3670a2ef4830dbb8..1e307b9eaa5f86571552630438a476abab2cc7c7 100644 --- a/src/MNH/profilern.f90 +++ b/src/MNH/profilern.f90 @@ -1,4 +1,4 @@ -!MNH_LIC Copyright 2002-2021 CNRS, Meteo-France and Universite Paul Sabatier +!MNH_LIC Copyright 2002-2022 CNRS, Meteo-France and Universite Paul Sabatier !MNH_LIC This is part of the Meso-NH software governed by the CeCILL-C licence !MNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt !MNH_LIC for details. version 1. @@ -31,7 +31,7 @@ REAL, DIMENSION(:,:,:), INTENT(IN) :: PP ! pressure REAL, DIMENSION(:,:,:,:), INTENT(IN) :: PAER ! aerosol extinction REAL, DIMENSION(:,:,:), INTENT(IN) :: PCLDFR ! cloud fraction REAL, DIMENSION(:,:,:), INTENT(IN) :: PCIT ! ice concentration -REAL, DIMENSION(:,:), INTENT(IN) :: PSEA ! for radar +REAL, DIMENSION(:,:), OPTIONAL, INTENT(IN) :: PSEA ! for radar ! !------------------------------------------------------------------------------- ! @@ -85,6 +85,9 @@ END MODULE MODI_PROFILER_n !! Philippe Wautelet: 05/2016-04/2018: new data structures and calls for I/O ! P. Wautelet 13/09/2019: budget: simplify and modernize date/time management ! M. Taufour 05/07/2021: modify RARE for hydrometeors containing ice and add bright band calculation for RARE +! P. Wautelet 09/02/2022: add message when some variables not computed +! + bugfix: put values in variables in this case +! + move some operations outside a do loop ! -------------------------------------------------------------------------- ! !* 0. DECLARATIONS @@ -97,12 +100,13 @@ USE MODD_GRID USE MODD_SUB_PROFILER_n USE MODD_NSV USE MODD_PARAMETERS -USE MODD_PARAM_n, ONLY: CCLOUD, CRAD +USE MODD_PARAM_n, ONLY: CCLOUD, CRAD, CSURF USE MODD_PROFILER_n USE MODD_TIME, only: tdtexp USE MODD_TIME_n, only: tdtcur ! USE MODE_ll +USE MODE_MSG ! USE MODI_GPS_ZENITH_GRID USE MODI_LIDAR @@ -159,7 +163,7 @@ REAL, DIMENSION(:,:,:), INTENT(IN) :: PP ! pressure REAL, DIMENSION(:,:,:,:), INTENT(IN) :: PAER ! aerosol extinction REAL, DIMENSION(:,:,:), INTENT(IN) :: PCLDFR ! cloud fraction REAL, DIMENSION(:,:,:), INTENT(IN) :: PCIT ! ice concentration -REAL, DIMENSION(:,:), INTENT(IN) :: PSEA ! for radar +REAL, DIMENSION(:,:), OPTIONAL, INTENT(IN) :: PSEA ! for radar ! !------------------------------------------------------------------------------- ! @@ -518,6 +522,14 @@ IF (GSTORE) THEN TPROFILER%ZTD(IN,I)= ZZTD_PROFILER TPROFILER%ZWD(IN,I)= ZZWD_PROFILER TPROFILER%ZHD(IN,I)= ZZHD_PROFILER + ELSE + CMNHMSG(1) = 'altitude of profiler ' // TRIM( TPROFILER%NAME(I) ) // ' is too far from orography' + CMNHMSG(2) = 'some variables are therefore not computed (IWV, ZTD, ZWD, ZHD)' + CALL PRINT_MSG( NVERB_WARNING, 'GEN', 'PROFILER_n' ) + TPROFILER%IWV(IN,I)= XUNDEF + TPROFILER%ZTD(IN,I)= XUNDEF + TPROFILER%ZWD(IN,I)= XUNDEF + TPROFILER%ZHD(IN,I)= XUNDEF END IF TPROFILER%ZON (IN,:,I) = ZU_PROFILER(:) * COS(ZGAM) + ZV_PROFILER(:) * SIN(ZGAM) TPROFILER%MER (IN,:,I) = - ZU_PROFILER(:) * SIN(ZGAM) + ZV_PROFILER(:) * COS(ZGAM) @@ -551,10 +563,15 @@ IF (GSTORE) THEN DO JLOOP=3,6 ZRZ(:,JLOOP)=PROFILER_INTERP(PR(:,:,:,JLOOP)) END DO - DO JK=1,IKU - ZRZ(JK,2)=PROFILER_INTERP_2D(PR(:,:,JK,2)*PSEA(:,:)) ! becomes cloud mixing ratio over sea - ZRZ(JK,7)=PROFILER_INTERP_2D(PR(:,:,JK,2)*(1.-PSEA(:,:))) ! becomes cloud mixing ratio over land - END DO + IF (CSURF=="EXTE") THEN + DO JK=1,IKU + ZRZ(JK,2)=PROFILER_INTERP_2D(PR(:,:,JK,2)*PSEA(:,:)) ! becomes cloud mixing ratio over sea + ZRZ(JK,7)=PROFILER_INTERP_2D(PR(:,:,JK,2)*(1.-PSEA(:,:))) ! becomes cloud mixing ratio over land + END DO + ELSE + ZRZ(:,2)=PROFILER_INTERP(PR(:,:,:,2)) + ZRZ(:,7)=0. + END IF ALLOCATE(ZAELOC(IKU)) ! ZAELOC(:)=0. @@ -897,10 +914,7 @@ IF (GSTORE) THEN CALL DISTRIBUTE_PROFILER(TPROFILER%CRARE(IN,JK,I)) CALL DISTRIBUTE_PROFILER(TPROFILER%CRARE_ATT(IN,JK,I)) CALL DISTRIBUTE_PROFILER(TPROFILER%CIZ(IN,JK,I)) - CALL DISTRIBUTE_PROFILER(TPROFILER%IWV(IN,I)) - CALL DISTRIBUTE_PROFILER(TPROFILER%ZTD(IN,I)) - CALL DISTRIBUTE_PROFILER(TPROFILER%ZHD(IN,I)) - CALL DISTRIBUTE_PROFILER(TPROFILER%ZWD(IN,I)) + ! IF (LDIAG_IN_RUN) CALL DISTRIBUTE_PROFILER(TPROFILER%TKE_DISS(IN,JK,I)) ! @@ -912,6 +926,11 @@ IF (GSTORE) THEN END DO IF (SIZE(PTKE)>0) CALL DISTRIBUTE_PROFILER(TPROFILER%TKE (IN,JK,I)) ENDDO + + CALL DISTRIBUTE_PROFILER(TPROFILER%IWV(IN,I)) + CALL DISTRIBUTE_PROFILER(TPROFILER%ZTD(IN,I)) + CALL DISTRIBUTE_PROFILER(TPROFILER%ZHD(IN,I)) + CALL DISTRIBUTE_PROFILER(TPROFILER%ZWD(IN,I)) ENDDO ! END IF diff --git a/src/MNH/read_all_data_grib_case.f90 b/src/MNH/read_all_data_grib_case.f90 index 841ae62b0c834902a5fec42fcbbac2c947e6d8c3..04a8ff5e8b5556152d2a61599a05cbe5ed0c25ef 100644 --- a/src/MNH/read_all_data_grib_case.f90 +++ b/src/MNH/read_all_data_grib_case.f90 @@ -1480,6 +1480,7 @@ END IF ! After this projection, the file is simil ! IF (HFILE(1:3)=='ATM') THEN +ISTARTLEVEL = 1 ALLOCATE (XU_LS(IIU,IJU,INLEVEL)) ALLOCATE (XV_LS(IIU,IJU,INLEVEL)) ALLOCATE (ZTU_LS(INO)) diff --git a/src/MNH/read_chem_data_netcdf_case.f90 b/src/MNH/read_chem_data_netcdf_case.f90 index 547a55b535791a53149ef6e20f232706ddfbff71..3ecc3f594c36347e046515c7e9a5995490a3ad02 100644 --- a/src/MNH/read_chem_data_netcdf_case.f90 +++ b/src/MNH/read_chem_data_netcdf_case.f90 @@ -87,6 +87,7 @@ END MODULE MODI_READ_CHEM_DATA_NETCDF_CASE ! P. Wautelet 10/04/2019: replace ABORT and STOP calls by Print_msg ! P. Wautelet 18/09/2019: correct support of 64bit integers (MNH_INT=8) ! P. Wautelet 09/03/2021: move some chemistry initializations to ini_nsv +! P. Wautelet 01/04/2022: add error if CDUMMY1 not set correctly !------------------------------------------------------------------------------- ! !* 0. DECLARATIONS @@ -395,6 +396,9 @@ ELSEIF (CDUMMY1=="18") THEN itimeindex=3 ELSEIF ((CDUMMY1=="24").OR.(CDUMMY1=="00")) THEN itimeindex=4 +ELSE + call Print_msg( NVERB_ERROR, 'GEN', 'READ_CHEM_DATA_NETCDF_CASE', 'CDUMMY1 is not set correctly (or not set at all)' ) + itimeindex=1 ENDIF istart3d(4) = itimeindex ! @@ -423,7 +427,7 @@ enddo istatus = nf90_get_var(incid, ips_varid, ZPSMOZ(:,:), start=istart2d, count=icount2d) if (istatus /= nf90_noerr) call handle_err(istatus) - + !------------------------------------------------------------------------ !* 3 Interpolation of MOZART variable !--------------------------------------------------------------------- @@ -580,7 +584,7 @@ DO JI = 1,IMOZ !for every MNH species existing in MOZ1.nam ENDIF ENDDO ! JNCHEM - DO JNAER = NSV_AERBEG, NSV_AEREND + DO JNAER = NSV_AERBEG, NSV_AEREND IF (trim(CAERONAMES(JNAER-NSV_AERBEG+1))==trim(YSPCMNH(JI))) THEN !MNH mechanism species IF (ISPCMOZ(JI)==1) THEN istatus = nf90_inq_varid(incid, trim(YCHANGE(JI,1)), ind_netcdf) @@ -661,7 +665,7 @@ DO JI = 1,IMOZ !for every MNH species existing in MOZ1.nam ENDIF ENDDO ! JNAER ENDDO ! JIDO JNCHEM = NSV_CHEMBEG, NSV_CHEMEND !loop on all MNH species -DEALLOCATE(YSPCMNH) +DEALLOCATE(YSPCMNH) DEALLOCATE(TZSTOC) DEALLOCATE(ISPCMOZ) DEALLOCATE(ZCOEFMOZART) diff --git a/src/MNH/read_exsegn.f90 b/src/MNH/read_exsegn.f90 index f1cb9ef1f6be28f28c87f8d7edf332c8a7f74e13..29b1f57918a3ee48bdc268c7f4e69d3cab8d3ba8 100644 --- a/src/MNH/read_exsegn.f90 +++ b/src/MNH/read_exsegn.f90 @@ -21,7 +21,7 @@ INTERFACE OCONDSAMP,OBLOWSNOW, & KRIMX,KRIMY, KSV_USER, & HTURB,HTOM,ORMC01,HRAD,HDCONV,HSCONV,HCLOUD,HELEC, & - HEQNSYS,PTSTEP_ALL,HSTORAGE_TYPE,HINIFILEPGD ) + HEQNSYS,PTSTEP_ALL,HINIFILEPGD ) ! USE MODD_IO, ONLY: TFILEDATA ! @@ -71,7 +71,6 @@ CHARACTER (LEN=4), INTENT(IN) :: HCLOUD ! Kind of microphysical scheme CHARACTER (LEN=4), INTENT(IN) :: HELEC ! Kind of electrical scheme CHARACTER (LEN=*), INTENT(IN) :: HEQNSYS! type of equations' system REAL,DIMENSION(:), INTENT(INOUT):: PTSTEP_ALL ! Time STEP of ALL models -CHARACTER (LEN=*), INTENT(IN) :: HSTORAGE_TYPE ! type of initial file CHARACTER (LEN=*), INTENT(IN) :: HINIFILEPGD ! name of PGD file ! END SUBROUTINE READ_EXSEG_n @@ -94,7 +93,7 @@ END MODULE MODI_READ_EXSEG_n OCONDSAMP, OBLOWSNOW, & KRIMX,KRIMY, KSV_USER, & HTURB,HTOM,ORMC01,HRAD,HDCONV,HSCONV,HCLOUD,HELEC, & - HEQNSYS,PTSTEP_ALL,HSTORAGE_TYPE,HINIFILEPGD ) + HEQNSYS,PTSTEP_ALL,HINIFILEPGD ) ! ######################################################################### ! !!**** *READ_EXSEG_n * - routine to read the descriptor file EXSEG @@ -302,7 +301,8 @@ END MODULE MODI_READ_EXSEG_n ! P. Wautelet 09/03/2021: move some chemistry initializations to ini_nsv ! P. Wautelet 10/03/2021: move scalar variable name initializations to ini_nsv ! R. Honnert 23/04/2021: add ADAP mixing length and delete HRIO and BOUT from CMF_UPDRAFT -! S. Riette 11/05/2021 HighLow cloud +! S. Riette 11/05/2021: HighLow cloud +! P. Wautelet 24/06/2022: remove check on CSTORAGE_TYPE for restart of ForeFire variables !------------------------------------------------------------------------------ ! !* 0. DECLARATIONS @@ -313,7 +313,6 @@ USE MODD_CH_AEROSOL USE MODD_CH_M9_n, ONLY : NEQ USE MODD_CONDSAMP USE MODD_CONF -USE MODD_CONF_n, ONLY: CSTORAGE_TYPE USE MODD_CONFZ ! USE MODD_DRAG_n USE MODD_DUST @@ -453,7 +452,6 @@ CHARACTER (LEN=4), INTENT(IN) :: HCLOUD ! Kind of microphysical scheme CHARACTER (LEN=4), INTENT(IN) :: HELEC ! Kind of electrical scheme CHARACTER (LEN=*), INTENT(IN) :: HEQNSYS! type of equations' system REAL,DIMENSION(:), INTENT(INOUT):: PTSTEP_ALL ! Time STEP of ALL models -CHARACTER (LEN=*), INTENT(IN) :: HSTORAGE_TYPE ! type of initial file CHARACTER (LEN=*), INTENT(IN) :: HINIFILEPGD ! name of PGD file ! !* 0.2 declarations of local variables @@ -1761,7 +1759,6 @@ END IF IF (CCLOUD == 'LIMA') THEN IF (HCLOUD == 'LIMA') THEN CGETSVT(NSV_LIMA_BEG:NSV_LIMA_END)='READ' -!!JPP IF(HSTORAGE_TYPE=='TT') CGETSVT(NSV_LIMA_BEG:NSV_LIMA_END)='INIT' ELSE WRITE(UNIT=ILUOUT,FMT=9001) KMI WRITE(UNIT=ILUOUT,FMT='("THERE IS NO SCALAR VARIABLES FOR LIMA & @@ -2056,9 +2053,6 @@ END IF IF (LFOREFIRE) THEN IF (OFOREFIRE) THEN CGETSVT(NSV_FFBEG:NSV_FFEND)='READ' - IF(HSTORAGE_TYPE=='TT') THEN - CGETSVT(NSV_FFBEG:NSV_FFEND)='INIT' - END IF ELSE WRITE(UNIT=ILUOUT,FMT=9001) KMI WRITE(UNIT=ILUOUT,FMT='("THERE IS NO FOREFIRE SCALAR VARIABLES IN INITIAL FMFILE",/,& diff --git a/src/MNH/read_surf_mnh.f90 b/src/MNH/read_surf_mnh.f90 index f051ed7393a5eb3042200c278fb0c6af3ebcbf0a..30b2a93712df54151a1f2dd8927e3e2de626d8d9 100644 --- a/src/MNH/read_surf_mnh.f90 +++ b/src/MNH/read_surf_mnh.f90 @@ -1,4 +1,4 @@ -!MNH_LIC Copyright 2003-2020 CNRS, Meteo-France and Universite Paul Sabatier +!MNH_LIC Copyright 2003-2022 CNRS, Meteo-France and Universite Paul Sabatier !MNH_LIC This is part of the Meso-NH software governed by the CeCILL-C licence !MNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt !MNH_LIC for details. version 1. @@ -1023,7 +1023,7 @@ IF (.NOT. GCOVER_PACKED) THEN TZFIELD%CSTDNAME = '' TZFIELD%CLONGNAME = TRIM(YREC) TZFIELD%CUNITS = '' - TZFIELD%CDIR = YDIR1 + TZFIELD%CDIR = YDIR TZFIELD%CCOMMENT = 'X_Y_'//TRIM(YREC) TZFIELD%NGRID = 4 TZFIELD%NTYPE = TYPEREAL diff --git a/src/MNH/stationn.f90 b/src/MNH/stationn.f90 index 4de047e716b3d60dda687a4abe41342d8f6dabcb..8774f316af1407a1f08b151fd469aa351461ae21 100644 --- a/src/MNH/stationn.f90 +++ b/src/MNH/stationn.f90 @@ -1,4 +1,4 @@ -!MNH_LIC Copyright 2002-2020 CNRS, Meteo-France and Universite Paul Sabatier +!MNH_LIC Copyright 2002-2022 CNRS, Meteo-France and Universite Paul Sabatier !MNH_LIC This is part of the Meso-NH software governed by the CeCILL-C licence !MNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt !MNH_LIC for details. version 1. @@ -80,6 +80,7 @@ END MODULE MODI_STATION_n ! P. Wautelet 05/2016-04/2018: new data structures and calls for I/O ! P. Wautelet 13/09/2019: budget: simplify and modernize date/time management ! R. Schoetter 11/2019: use LCARTESIAN instead of LSTATLAT for multiproc in cartesian +! P. Wautelet 09/05/2022: bugfix: use correct indices for U and V interpolation ! ! -------------------------------------------------------------------------- ! @@ -491,7 +492,7 @@ ELSEIF (L1D) THEN JI=2 JJ=2 ELSE - JI=II(I) + JI=IU(I) JJ=IJ(I) END IF ! @@ -521,7 +522,7 @@ ELSEIF (L1D) THEN JJ=2 ELSE JI=II(I) - JJ=IJ(I) + JJ=IV(I) END IF ! IF ((JI .GT. 0).AND. (JI .LT. SIZE(PA,1)) .AND. & diff --git a/src/MNH/write_lfifm1_for_diag_supp.f90 b/src/MNH/write_lfifm1_for_diag_supp.f90 index 0175e2c28a7b797264a2578763d00077218aadcd..c067e62189f73bf435dff6b1160824cd1c3edc74 100644 --- a/src/MNH/write_lfifm1_for_diag_supp.f90 +++ b/src/MNH/write_lfifm1_for_diag_supp.f90 @@ -847,15 +847,15 @@ END IF ! IF (NRTTOVINFO(1,1) /= NUNDEF) THEN ! PRINT*,'YOU ASK FOR BRIGHTNESS TEMPERATURE COMPUTED BY THE RTTOV CODE' -#ifdef MNH_RTTOV_8 +#if defined(MNH_RTTOV_8) CALL CALL_RTTOV8(NDLON, NFLEV, NSTATM, XEMIS(:,:,1), XTSRAD, XSTATM, XTHT, XRT, & XPABST, XZZ, XMFCONV, MAX(XCLDFR,XICEFR), XUT(:,:,IKB), XVT(:,:,IKB), & LUSERI, NRTTOVINFO, TPFILE ) -#elif MNH_RTTOV_11 +#elif defined(MNH_RTTOV_11) CALL CALL_RTTOV11(NDLON, NFLEV, XEMIS(:,:,1), XTSRAD, XTHT, XRT, & XPABST, XZZ, XMFCONV, MAX(XCLDFR,XICEFR), XUT(:,:,IKB), XVT(:,:,IKB), & LUSERI, NRTTOVINFO, TPFILE ) -#elif MNH_RTTOV_13 +#elif defined(MNH_RTTOV_13) CALL CALL_RTTOV13(NDLON, NFLEV, XEMIS(:,:,1), XTSRAD, XTHT, XRT, & XPABST, XZZ, XMFCONV, MAX(XCLDFR,XICEFR), XUT(:,:,IKB), XVT(:,:,IKB), & LUSERI, NRTTOVINFO, TPFILE ) diff --git a/src/Makefile b/src/Makefile index 03c5d4492c4f7bc0cba2db78ed82a3ad46a152ae..d996ff973445adca08d33056dbc64e8425d7c3c4 100644 --- a/src/Makefile +++ b/src/Makefile @@ -321,7 +321,7 @@ eccodes_lib : $(ECCODES_MOD) $(ECCODES_MOD) : - [ ! -d $(DIR_ECCODES_BUILD) ] && mkdir -p $(DIR_ECCODES_BUILD) cd ${DIR_ECCODES_BUILD} && \ - cmake ${DIR_ECCODES_SRC} -DCMAKE_INSTALL_PREFIX=${DIR_ECCODES_INSTALL} -DBUILD_SHARED_LIBS=OFF \ + AEC_PATH=$(CDF_PATH) cmake ${DIR_ECCODES_SRC} -DCMAKE_INSTALL_PREFIX=${DIR_ECCODES_INSTALL} -DBUILD_SHARED_LIBS=OFF \ -DENABLE_NETCDF=OFF -DENABLE_JPG=OFF -DENABLE_PYTHON=OFF -DENABLE_EXAMPLES=OFF \ -DCMAKE_Fortran_COMPILER=$(FC) -DCMAKE_C_COMPILER=$(CC) \ -DCMAKE_Fortran_FLAGS=$(ECCODES_FFLAGS) -DCMAKE_C_FLAGS=$(ECCODES_CFLAGS) && \ diff --git a/src/SURFEX/mode_read_grib.F90 b/src/SURFEX/mode_read_grib.F90 index f92149d19a20c6d00eeb29570069aecf84e819ce..50ba15ed1cff423436c155bcdab4843e26376a81 100644 --- a/src/SURFEX/mode_read_grib.F90 +++ b/src/SURFEX/mode_read_grib.F90 @@ -446,6 +446,7 @@ INTEGER :: ILEV ! level INTEGER :: INUM_ZS,ISIZE,ICOUNT,JLOOP,IPARAM,IGRIB,IPARAM2 CHARACTER(LEN=24) :: YLTYPELU REAL(KIND=JPRB) :: ZHOOK_HANDLE +CHARACTER(LEN=50) :: CNAME ! name of the parameter (for ARPEGE GRIB2 converted with Epygram 1.4.8) !------------------------------------------------------------------- IF (LHOOK) CALL DR_HOOK('MODE_READ_GRIB:READ_GRIB_LAND_MASK',0,ZHOOK_HANDLE) WRITE (KLUOUT,'(A)') 'MODE_READ_GRIB:READ_GRIB_LAND_MASK: | Reading land mask from ',HINMODEL @@ -458,7 +459,8 @@ SELECT CASE (HINMODEL) CASE ('ARPEGE','ALADIN','MOCAGE') IF(HINMODEL=='ARPEGE' .AND. NGRIB_VERSION==2) THEN ILTYPE=1 - CALL READ_GRIB(HGRIB,HINMODEL,KLUOUT,0,IRET,PMASK,KLTYPE=ILTYPE) + CNAME='Land-sea mask' + CALL READ_GRIB(HGRIB,HINMODEL,KLUOUT,0,IRET,PMASK,KLTYPE=ILTYPE,HNAME=CNAME) ELSE CALL READ_GRIB(HGRIB,HINMODEL,KLUOUT,81,IRET,PMASK) ENDIF diff --git a/src/SURFEX/pgd_grid.F90 b/src/SURFEX/pgd_grid.F90 index 5a7158fc109dd922e185beef9c6b51261a77ee8e..9c17d153a3428715d0cb1593b01b253b90c94fb8 100644 --- a/src/SURFEX/pgd_grid.F90 +++ b/src/SURFEX/pgd_grid.F90 @@ -316,7 +316,7 @@ ALLOCATE(UG%XJPDIR (NL)) CALL MPPDB_CHECK_SURFEX2D(UG%G%XLAT,"PGD_GRID after LATLON_GRID:XLAT",PRECISION,ILUOUT) CALL MPPDB_CHECK_SURFEX2D(UG%G%XLON,"PGD_GRID after LATLON_GRID:XLON",PRECISION,ILUOUT) CALL MPPDB_CHECK_SURFEX2D(UG%G%XMESH_SIZE,"PGD_GRID after LATLON_GRID:XMESH_SIZE",PRECISION,ILUOUT) -#endif! +#endif !------------------------------------------------------------------------------ ! !* 7. Average grid length (in degrees) diff --git a/src/configure b/src/configure index bf9ce3128cccb7730631c1ec4bb936abd15aeaf5..348092bcedb6e841ea16a41b21f4488f56c579ad 100755 --- a/src/configure +++ b/src/configure @@ -9,7 +9,7 @@ if [ "x$XYZ" = "x" ] then # export VERSION_MASTER=${VERSION_MASTER:-MNH-V5-5} -export VERSION_BUG=${VERSION_BUG:-0} +export VERSION_BUG=${VERSION_BUG:-1} export VERSION_XYZ=${VERSION_XYZ:-${VERSION_MASTER}-${VERSION_BUG}${VER_OASIS:+-${VER_OASIS}}} export VERSION_DATE=${VERSION_DATE:-"19/03/2021"} export VERSION_LIBAEC=${VERSION_LIBAEC:-"0.3.4"} @@ -551,12 +551,6 @@ if [ "x${VER_OASIS}" == "xOASISAUTO" ] ; then ( cd $LOCAL/src/LIB ; [ ! -d oasis3-${VERSION_OASIS} ] && tar xvfz oasis3-${VERSION_OASIS}.tar.gz ; [ ! -d toy_${VERSION_TOY} ] && tar xvfz toy_${VERSION_TOY}.tar.gz ) fi # -# Install MEGAN if MNH_MEGAN=1 -# -if [ "x${MNH_MEGAN}" == "x1" ] ; then -( cd $LOCAL/src/LIB ; [ ! -d MEGAN ] && tar xvfz megan.tar.gz ) -fi -# # Install GRIBAPI or ecCodes # if [ "x${MNH_GRIBAPI}" == "xyes" ] ; then diff --git a/src/job_make_examples_BG b/src/job_make_examples_BG index 8f18e061117c3a2ae65983ee927bb242926e7304..aaaec9e56bafa748c1bc75dc8519aa4537114d63 100755 --- a/src/job_make_examples_BG +++ b/src/job_make_examples_BG @@ -18,7 +18,7 @@ set -x cd $LOADL_STEP_INITDIR -. ../conf/profile_mesonh-BG-R8I4-MNH-V5-5-0-MPIAUTO-O2 +. ../conf/profile_mesonh-BG-R8I4-MNH-V5-5-1-MPIAUTO-O2 #001_2Drelief 002_3Drelief 003_KW78 004_Reunion 007_16janvier diff --git a/src/job_make_examples_BGQ b/src/job_make_examples_BGQ index a6785fd9ff187b157b4567db85f44ebd08e982f0..644bd8a961b8e6be89bcba4556425fbf19d450b8 100755 --- a/src/job_make_examples_BGQ +++ b/src/job_make_examples_BGQ @@ -18,7 +18,7 @@ cd $LOADL_STEP_INITDIR -. ../conf/profile_mesonh-BGQ-R8I4-MNH-V5-5-0-MPIAUTO-O2NAN +. ../conf/profile_mesonh-BGQ-R8I4-MNH-V5-5-1-MPIAUTO-O2NAN set -x diff --git a/src/job_make_examples_BullX b/src/job_make_examples_BullX index 2a29d889ee1533ecd009169d6427f4c2dada3b18..8297eb5b9142c7353ad0cf22fb9b36bd01bccfcc 100755 --- a/src/job_make_examples_BullX +++ b/src/job_make_examples_BullX @@ -19,7 +19,7 @@ set -e hostname # Echo des commandes -. ../conf/profile_mesonh-LXifort-R8I4-MNH-V5-5-0-MPIINTEL-O3 +. ../conf/profile_mesonh-LXifort-R8I4-MNH-V5-5-1-MPIINTEL-O3 export MONORUN="Mpirun -np 1 " export MPIRUN="Mpirun -np 2 " export POSTRUN="time " diff --git a/src/job_make_examples_BullX_belenos b/src/job_make_examples_BullX_belenos index 4aa4f1129fdbcf671a4d475cc998234ddc22f0d1..d27459f20cf6b2e38cc1d910cb9b99a5fb1500a9 100755 --- a/src/job_make_examples_BullX_belenos +++ b/src/job_make_examples_BullX_belenos @@ -18,7 +18,7 @@ set -e hostname # Echo des commandes -. ../conf/profile_mesonh-LXifort-R8I4-MNH-V5-5-0-MPIAUTO-O3 +. ../conf/profile_mesonh-LXifort-R8I4-MNH-V5-5-1-MPIAUTO-O3 export MONORUN="Mpirun -np 1 " export MPIRUN="Mpirun -np 2 " export POSTRUN="echo " diff --git a/src/job_make_examples_BullX_irene b/src/job_make_examples_BullX_irene index b59dd5bf59e0a9365cab6f715c89903575d63ad1..852288dec893574c1a5280ca57c41f0d1f72628e 100755 --- a/src/job_make_examples_BullX_irene +++ b/src/job_make_examples_BullX_irene @@ -21,7 +21,7 @@ set +x # Nom de la machine hostname -. ../conf/profile_mesonh-LXifort-R8I4-MNH-V5-5-0-MPIAUTO-O2 +. ../conf/profile_mesonh-LXifort-R8I4-MNH-V5-5-1-MPIAUTO-O2 set -x diff --git a/src/job_make_examples_BullX_irene_AMD b/src/job_make_examples_BullX_irene_AMD index e2e9c6d7d98f1322c881ae0c0138fb4670e840e4..f1afe22f74ef63b55919b251975f95c229161dd1 100755 --- a/src/job_make_examples_BullX_irene_AMD +++ b/src/job_make_examples_BullX_irene_AMD @@ -20,7 +20,7 @@ set +x # Nom de la machine hostname -. ../conf/profile_mesonh-LXifort-R8I4-MNH-V5-5-0-AMD-MPIAUTO-O2 +. ../conf/profile_mesonh-LXifort-R8I4-MNH-V5-5-1-AMD-MPIAUTO-O2 set -x diff --git a/src/job_make_examples_BullX_occigen b/src/job_make_examples_BullX_occigen index 1a1cd8bc9c58bdc8ed8584db9caa25f46d608018..63aa7239118836bc18772147f73f65b2f5fc2766 100755 --- a/src/job_make_examples_BullX_occigen +++ b/src/job_make_examples_BullX_occigen @@ -18,7 +18,7 @@ set -x # Nom de la machine hostname -. ../conf/profile_mesonh-LXifort-R8I4-MNH-V5-5-0-MPIINTEL-O2 +. ../conf/profile_mesonh-LXifort-R8I4-MNH-V5-5-1-MPIINTEL-O2 export MONORUN="Mpirun -prepend-rank -np 1 " export MPIRUN="Mpirun -prepend-rank -np 4 " export POSTRUN="echo " diff --git a/src/job_make_examples_BullX_olympe b/src/job_make_examples_BullX_olympe index 27f0e8c970502733f890ecd81837aea4f68acb54..ba3fa83fb515d22365976e3b088520c02d7232d8 100755 --- a/src/job_make_examples_BullX_olympe +++ b/src/job_make_examples_BullX_olympe @@ -17,7 +17,7 @@ set -x # Nom de la machine hostname -. ../conf/profile_mesonh-LXifort-R8I4-MNH-V5-5-0-MPIINTEL-O2 +. ../conf/profile_mesonh-LXifort-R8I4-MNH-V5-5-1-MPIINTEL-O2 export MONORUN="mpirun -prepend-rank -np 1 " export MPIRUN="mpirun -prepend-rank -np 4 " export POSTRUN="echo " diff --git a/src/job_make_examples_CRAY_cca b/src/job_make_examples_CRAY_cca index 786947f40e6fdbe4689be5aae55636b36e17ec61..9296c3799e1856c0e2cad8127331efbf5fe377f3 100755 --- a/src/job_make_examples_CRAY_cca +++ b/src/job_make_examples_CRAY_cca @@ -28,7 +28,7 @@ cd ${PBS_O_WORKDIR} ARCH=LXifort #ARCH=LXcray -. ../conf/profile_mesonh-${ARCH}-R8I4-MNH-V5-5-0-MPICRAY-O2 +. ../conf/profile_mesonh-${ARCH}-R8I4-MNH-V5-5-1-MPICRAY-O2 export MONORUN="aprun -n 1 " diff --git a/src/job_make_examples_HPE_jeanzay b/src/job_make_examples_HPE_jeanzay index bd03e566f6090cdf224e080e65c11f8bc549f233..1687368fd3b871ff8123d957104070bedf5ef1e5 100755 --- a/src/job_make_examples_HPE_jeanzay +++ b/src/job_make_examples_HPE_jeanzay @@ -19,7 +19,7 @@ set -x # Nom de la machine hostname -. ../conf/profile_mesonh-LXifort-R8I4-MNH-V5-5-0-MPIINTEL-O2 +. ../conf/profile_mesonh-LXifort-R8I4-MNH-V5-5-1-MPIINTEL-O2 export MONORUN="Exec srun -l -n 1 --export=ALL numabind_core_slurm" export MPIRUN="Exec srun -l -n 4 --export=ALL numabind_core_slurm" export POSTRUN="echo " diff --git a/src/job_make_examples_IBM_ada b/src/job_make_examples_IBM_ada index a1b1b51f54543532d5eff7d9f115d950a6c83c80..044996e27a3a490d66d3a4d5c5b55a80c64f8237 100755 --- a/src/job_make_examples_IBM_ada +++ b/src/job_make_examples_IBM_ada @@ -19,7 +19,7 @@ cd $LOADL_STEP_INITDIR -. ../conf/profile_mesonh-LXifort-R8I4-MNH-V5-5-0-MPIINTEL-O2 +. ../conf/profile_mesonh-LXifort-R8I4-MNH-V5-5-1-MPIINTEL-O2 # Pour avoir l'echo des commandes set -x diff --git a/src/job_make_examples_IBM_sp6_vargas b/src/job_make_examples_IBM_sp6_vargas index c8d65b89b0dd0089e8822d27a3264a4cc6a877d9..3d087d03d42416d69f26f1db42778e4fe3c12a40 100755 --- a/src/job_make_examples_IBM_sp6_vargas +++ b/src/job_make_examples_IBM_sp6_vargas @@ -24,7 +24,7 @@ set -x cd $LOADL_STEP_INITDIR -. ../conf/profile_mesonh-AIX64-R8I4-MNH-V5-5-0-MPIAUTO-O2 +. ../conf/profile_mesonh-AIX64-R8I4-MNH-V5-5-1-MPIAUTO-O2 #001_2Drelief 002_3Drelief 003_KW78 004_Reunion 007_16janvier diff --git a/src/job_make_examples_NEC_SX8 b/src/job_make_examples_NEC_SX8 index d33de70361b74d52899c0f8d9866686c76d1d115..97abd2caf8716875d30ca6deb9569687cfcb6b85 100755 --- a/src/job_make_examples_NEC_SX8 +++ b/src/job_make_examples_NEC_SX8 @@ -18,7 +18,7 @@ hostname [ -d $PBS_O_WORKDIR ] && cd $PBS_O_WORKDIR # -. ../conf/profile_mesonh-SX8-R8I4-MNH-V5-5-0-MPIAUTO-O4 +. ../conf/profile_mesonh-SX8-R8I4-MNH-V5-5-1-MPIAUTO-O4 export MONORUN="Mpirun -np 1 " export MPIRUN="Mpirun -np 2 " diff --git a/src/job_make_examples_SX8 b/src/job_make_examples_SX8 index fdd3553013dfc4ba5e40ff30945f81ebcf8a4acd..0dad22c62c62e840b29cbf2508c19cdba11ab613 100755 --- a/src/job_make_examples_SX8 +++ b/src/job_make_examples_SX8 @@ -19,7 +19,7 @@ hostname [ -d $PBS_O_WORKDIR ] && cd $PBS_O_WORKDIR # -. ../conf/profile_mesonh-SX8-R8I4-MNH-V5-5-0-MPIAUTO-O2 +. ../conf/profile_mesonh-SX8-R8I4-MNH-V5-5-1-MPIAUTO-O2 export MONORUN="Mpirun -np 1 " export MPIRUN="Mpirun -np 2 " diff --git a/src/job_make_examples_cxa b/src/job_make_examples_cxa index ef8507f6bb94d9d0ac2964f8b808c7dd99325022..458d780da0a7349845cc1edf7bd74bc4a8ef29da 100755 --- a/src/job_make_examples_cxa +++ b/src/job_make_examples_cxa @@ -34,7 +34,7 @@ echo SHELL=$SHELL cd $LOADL_STEP_INITDIR -. ../conf/profile_mesonh-AIX64-R8I4-MNH-V5-5-0-MPIAUTO-O2 +. ../conf/profile_mesonh-AIX64-R8I4-MNH-V5-5-1-MPIAUTO-O2 ulimit -c 0 # pas de core diff --git a/src/job_make_mesonh_BG b/src/job_make_mesonh_BG index ce08d470e8b49270a701f2b31ffd9461a2932efb..a2b996576ae8ae5c912d612f694364b205cf801f 100755 --- a/src/job_make_mesonh_BG +++ b/src/job_make_mesonh_BG @@ -18,7 +18,7 @@ set -x cd $LOADL_STEP_INITDIR -. ../conf/profile_mesonh-BG-R8I4-MNH-V5-5-0-MPIAUTO-O2 +. ../conf/profile_mesonh-BG-R8I4-MNH-V5-5-1-MPIAUTO-O2 #time gmake time gmake -r -j8 diff --git a/src/job_make_mesonh_BGQ b/src/job_make_mesonh_BGQ index a8e332a69e6e87ee81f4b79577fe18034eb014b2..dd82b6dbae59418381421bfceb994a718e89ee9c 100755 --- a/src/job_make_mesonh_BGQ +++ b/src/job_make_mesonh_BGQ @@ -20,7 +20,7 @@ set -x cd $LOADL_STEP_INITDIR -. ../conf/profile_mesonh-BGQ-R8I4-MNH-V5-5-0-MPIAUTO-O2NAN +. ../conf/profile_mesonh-BGQ-R8I4-MNH-V5-5-1-MPIAUTO-O2NAN time gmake -j8 time gmake -j8 diff --git a/src/job_make_mesonh_BullX b/src/job_make_mesonh_BullX index 967c115e1dcee2800b1b01ed89733a9e924bf30b..ece2b8a1ed3103115efb4174268153a8d71a68ee 100755 --- a/src/job_make_mesonh_BullX +++ b/src/job_make_mesonh_BullX @@ -19,7 +19,7 @@ set -x # On va lancer la compilation dans le répertoire de lancement du job pwd -. ../conf/profile_mesonh-LXifort-R8I4-MNH-V5-5-0-MPIINTEL-O3 +. ../conf/profile_mesonh-LXifort-R8I4-MNH-V5-5-1-MPIINTEL-O3 time gmake -j 4 time gmake -j 1 installmaster diff --git a/src/job_make_mesonh_BullX_belenos b/src/job_make_mesonh_BullX_belenos index 5fe6d8a062e5c331d3e66c16aa044ee52ed5bbf6..c9cd95a18d19cf125dea392c635a717f39efdaca 100755 --- a/src/job_make_mesonh_BullX_belenos +++ b/src/job_make_mesonh_BullX_belenos @@ -16,7 +16,7 @@ set -x # On va lancer la compilation dans le répertoire de lancement du job pwd -. ../conf/profile_mesonh-LXifort-R8I4-MNH-V5-5-0-MPIAUTO-O3 +. ../conf/profile_mesonh-LXifort-R8I4-MNH-V5-5-1-MPIAUTO-O3 time make -j 16 time make -j 1 installmaster diff --git a/src/job_make_mesonh_CRAY_cca b/src/job_make_mesonh_CRAY_cca index 56dba703cf90a3b45a985cc4b68bc64551331a92..5951a84069c96d98b124eeb606d0bff9e25f085b 100755 --- a/src/job_make_mesonh_CRAY_cca +++ b/src/job_make_mesonh_CRAY_cca @@ -21,7 +21,7 @@ pwd ARCH=LXifort #ARCH=LXcray -. ../conf/profile_mesonh-${ARCH}-R8I4-MNH-V5-5-0-MPICRAY-O2 +. ../conf/profile_mesonh-${ARCH}-R8I4-MNH-V5-5-1-MPICRAY-O2 time gmake -j 4 2>&1 | tee sortie_compile_${ARCH}.$$ time gmake -j 4 2>&1 | tee sortie_compile_${ARCH}2.$$ diff --git a/src/job_make_mesonh_HPE_jeanzay b/src/job_make_mesonh_HPE_jeanzay index 0dfa8598dda30dfaef3ef7268e2e02d184cfc0d3..417802b7da16d50ff373515b5f9a1c706a772f80 100755 --- a/src/job_make_mesonh_HPE_jeanzay +++ b/src/job_make_mesonh_HPE_jeanzay @@ -14,7 +14,7 @@ set -x # On va lancer la compilation dans le répertoire de lancement du job pwd -. ../conf/profile_mesonh-LXifort-R8I4-MNH-V5-5-0-MPIINTEL-O2 +. ../conf/profile_mesonh-LXifort-R8I4-MNH-V5-5-1-MPIINTEL-O2 time gmake -j 16 time gmake -j 1 installmaster diff --git a/src/job_make_mesonh_IBM_ada b/src/job_make_mesonh_IBM_ada index bcdfcf110cdeb34e78dfcb55c06801c7ab1b2cd5..b64aa08d990e865e5005051c78cf95391251c5df 100755 --- a/src/job_make_mesonh_IBM_ada +++ b/src/job_make_mesonh_IBM_ada @@ -16,7 +16,7 @@ cd $LOADL_STEP_INITDIR -. ../conf/profile_mesonh-LXifort-R8I4-MNH-V5-5-0-MPIINTEL-O2 +. ../conf/profile_mesonh-LXifort-R8I4-MNH-V5-5-1-MPIINTEL-O2 # Pour avoir l'echo des commandes set -x diff --git a/src/job_make_mesonh_IBM_sp6_vargas b/src/job_make_mesonh_IBM_sp6_vargas index b5f606d6381285d09b9b62a66944331ad329b758..7f01ee8b03a80c8e169cd75e40df81ffaa4ecf9c 100755 --- a/src/job_make_mesonh_IBM_sp6_vargas +++ b/src/job_make_mesonh_IBM_sp6_vargas @@ -24,7 +24,7 @@ set -x cd $LOADL_STEP_INITDIR -. ../conf/profile_mesonh-AIX64-R8I4-MNH-V5-5-0-MPIAUTO-O2 +. ../conf/profile_mesonh-AIX64-R8I4-MNH-V5-5-1-MPIAUTO-O2 time gmake -r -j8 time gmake installmaster diff --git a/src/job_make_mesonh_MFSX8 b/src/job_make_mesonh_MFSX8 index 811e5aa92ab70293537a9d11cbf048c80f1d0089..326b7b9dda718c6acfc2dcb08479db736bbd2905 100644 --- a/src/job_make_mesonh_MFSX8 +++ b/src/job_make_mesonh_MFSX8 @@ -12,7 +12,7 @@ set -x # On va lancer la compilation dans le répertoire de lancement du job [ ${PBS_O_WORKDIR} ] && cd ${PBS_O_WORKDIR} -. ../conf/profile_mesonh-SX8-R8I4-MNH-V5-5-0-MPIAUTO-O4 +. ../conf/profile_mesonh-SX8-R8I4-MNH-V5-5-1-MPIAUTO-O4 time gmake -j 4 ########## compile on four processors to speedup the compilation time gmake -j 1 installmaster diff --git a/src/job_make_mesonh_NEC_SX8 b/src/job_make_mesonh_NEC_SX8 index 75d68db4ade1024e25a21842fcad41eccb63f7e1..9160501349021c3267d15ab219416e9bb3535f8b 100755 --- a/src/job_make_mesonh_NEC_SX8 +++ b/src/job_make_mesonh_NEC_SX8 @@ -11,7 +11,7 @@ set -x # On va lancer la compilation dans le répertoire de lancement du job [ $PBS_O_WORKDIR ] && cd $PBS_O_WORKDIR -. ../conf/profile_mesonh-SX8-R8I4-MNH-V5-5-0-MPIAUTO-O4 +. ../conf/profile_mesonh-SX8-R8I4-MNH-V5-5-1-MPIAUTO-O4 time gmake -j 4 time gmake -j 4 # some time problem with first pass in parallel compilation diff --git a/src/job_make_mesonh_cxa b/src/job_make_mesonh_cxa index d531fb9a75484db4fa6dd5c88d509c2188454d9e..f9f3c1c573e52bb9cea9f32b465a2549a01daa85 100755 --- a/src/job_make_mesonh_cxa +++ b/src/job_make_mesonh_cxa @@ -27,7 +27,7 @@ set -x cd $LOADL_STEP_INITDIR -. ../conf/profile_mesonh-AIX64-R8I4-MNH-V5-5-0-MPIAUTO-O2 +. ../conf/profile_mesonh-AIX64-R8I4-MNH-V5-5-1-MPIAUTO-O2 time gmake -r -j1 time gmake installmaster diff --git a/src/job_make_mesonh_user_BullX b/src/job_make_mesonh_user_BullX index 0b982175890de1f9249f8f361a7b5e2492a5ce66..3251ac554a66c6e37b31ff2313e08e4474acadad 100755 --- a/src/job_make_mesonh_user_BullX +++ b/src/job_make_mesonh_user_BullX @@ -19,7 +19,7 @@ export VER_USER= ########## Your own USER Directory set -x # On va lancer la compilation dans le répertoire de lancement du job -. ../conf/profile_mesonh-LXifort-R8I4-MNH-V5-5-0-${VER_USER}-MPIINTEL-O3 +. ../conf/profile_mesonh-LXifort-R8I4-MNH-V5-5-1-${VER_USER}-MPIINTEL-O3 time gmake user time gmake -j 1 installuser diff --git a/src/job_make_mesonh_user_BullX_belenos b/src/job_make_mesonh_user_BullX_belenos index fddbb21e59d239bc8be88461fc41d3f3fba6d353..ad8037af1faa94c6567f53d96d3e044f48519259 100755 --- a/src/job_make_mesonh_user_BullX_belenos +++ b/src/job_make_mesonh_user_BullX_belenos @@ -18,7 +18,7 @@ set -x # On va lancer la compilation dans le répertoire de lancement du job pwd -. ../conf/profile_mesonh-LXifort-R8I4-MNH-V5-5-0-${VER_USER}-MPIAUTO-O3 +. ../conf/profile_mesonh-LXifort-R8I4-MNH-V5-5-1-${VER_USER}-MPIAUTO-O3 time make user -j 2 time make -j 1 installuser diff --git a/src/job_make_mesonh_user_MFSX8 b/src/job_make_mesonh_user_MFSX8 index e6a87c197a2c5157aa4a5cc82bf8a50189b1e157..0804a1caf9aa323d9bb0112b4c6815c89f52ce9c 100644 --- a/src/job_make_mesonh_user_MFSX8 +++ b/src/job_make_mesonh_user_MFSX8 @@ -14,7 +14,7 @@ set -x [ ${PBS_O_WORKDIR} ] && cd ${PBS_O_WORKDIR} -. ../conf/profile_mesonh-SX8-R8I4-MNH-V5-5-0-${VER_USER}-MPIAUTO-O4 +. ../conf/profile_mesonh-SX8-R8I4-MNH-V5-5-1-${VER_USER}-MPIAUTO-O4 time gmake user time gmake -j 1 installuser