Index: /TOOLS/CMIP6_FORCING/SCENARIOS/create_O3_hybrid_CNRM_Tier2.sh
===================================================================
--- /TOOLS/CMIP6_FORCING/SCENARIOS/create_O3_hybrid_CNRM_Tier2.sh	(revision 4226)
+++ /TOOLS/CMIP6_FORCING/SCENARIOS/create_O3_hybrid_CNRM_Tier2.sh	(revision 4226)
@@ -0,0 +1,319 @@
+# Script gÃ©nÃ©rant les fichiers ozone hybride CNRM pour les ssp du Tier 2 (en l'absence de l'ozone officiel UReading)
+# -----------------------------------
+# Th. Lurton - janvier 2019
+# -----------------------------------
+# Ce script est Ã  exÃ©cuter sur Ciclad
+# inclus dans un job.
+# -----------------------------------
+
+# Chemins Ã  personnaliser
+ici=/data/tlurton/CMIP6/O3_Tier2/
+repout=${ici}Tier2_OUT/
+
+# Chemin des donnÃ©es fournies par le CNRM
+repCNRM=/data/tlurton/CMIP6/OzoneCNRM/
+
+# Grande boucle sur les annÃ©es
+# On a 2100 du CNRM, mais pas 2100 d'UReading, ce qui nous oblige Ã  ne pas traiter la derniÃšre annÃ©e. Au final on dupliquera 2100 = 2099.
+for year in {2015..2099}
+do
+
+echo "========================================================================="
+echo "================================= "${year}" ================================="
+echo "========================================================================="
+
+# PHASE 1
+# INTERPOLATION VERTICALE NIVEAUX CNRM --> NIVEAUX UREADING
+
+# --> ici boucle sur tous les ssp
+# On n'aura pas besoin de ssp585
+
+for ssp in "ssp119" "ssp126" "ssp245" "ssp370" "ssp434" "ssp460" "ssp534-over"
+do
+
+echo "Processing year "${year}" for "${ssp}"..."
+
+fileCNRM=${repCNRM}o3_HOMEPmon_CNRM-ESM2-1_${ssp}_r1i1p1f2_gr_${year}01-${year}12.nc
+
+filetmp1=${ici}${ssp}_${year}_CNRM_hPa.nc
+filetmp2=${ici}${ssp}_${year}_CNRM_hPa_interpURea.nc
+
+# Attention !
+# Les niveaux plev des fichiers CNRM sont en Pa
+# Les niveaux plev des fichiers UReading sont en hPa
+
+# Niveaux CNRM entre 1000 hPa 3 Pa :
+# [Pa] 100000,92500,85000,70000,60000,50000,40000,30000,25000,20000,17000,15000,13000,11500,10000,9000,8000,7000,5000,3000,2000,1500,1000,700,500,300,200,150,100,70,50,40,30,20,15,10,7,5,3
+# [hPa] 1000,925,850,700,600,500,400,300,250,200,170,150,130,115,100,90,80,70,50,30,20,15,10,7,5,3,2,1.5,1,.7,.5,.4,.3,.2,.15,.1,.07,.05,.03
+
+# Niveaux Reading, liste totale
+# [hPa] 1000,925,850,800,780,750,700,650,600,500,450,400,350,300,285,250,200,170,150,130,115,100,90,80,70,60,50,40,35,30,25,20,15,10,7,5,4,3,2,1.5,1,0.7,0.5,0.4,0.3,0.2,0.15,0.1,0.07,0.05,0.04,0.03,0.02,0.015,0.01,0.007,0.005,0.004,0.003,0.002,0.0015,0.001,0.0008,0.0005,0.0003,0.0001
+
+# Niveaux Reading entre 1000 hPa et 3 Pa
+# [hPa] 1000,925,850,800,780,750,700,650,600,500,450,400,350,300,285,250,200,170,150,130,115,100,90,80,70,60,50,40,35,30,25,20,15,10,7,5,4,3,2,1.5,1,0.7,0.5,0.4,0.3,0.2,0.15,0.1,0.07,0.05,0.04,0.03
+
+# Conversion des plev CNRM vers hPa, en utilisant le fichier "CNRM_axis_en_hPa" : overwrite de la dimension verticale
+echo "... Converting CNRM pressure levels towards hPa..."
+cdo setzaxis,CNRM_axis_en_hPa ${fileCNRM} ${filetmp1}
+
+# Projection sur les niveaux Reading
+echo "... Interpolating CNRM file onto UReading pressure levels..."
+cdo intlevel,1000,925,850,800,780,750,700,650,600,500,450,400,350,300,285,250,200,170,150,130,115,100,90,80,70,60,50,40,35,30,25,20,15,10,7,5,4,3,2,1.5,1,0.7,0.5,0.4,0.3,0.2,0.15,0.1,0.07,0.05,0.04,0.03 ${filetmp1} ${filetmp2}
+
+rm -f ${filetmp1}
+
+# --> Fin des deux boucles sur les ssp
+done
+
+
+# PHASE 2
+# "RÃGLE DE TROIS"
+
+# --> Boucle sur les ssp du Tier 2 uniquement : ssp119, ssp434, ssp460 et ssp534-over
+
+for sspT2 in "ssp119" "ssp434" "ssp460" "ssp534-over"
+do
+
+# choix des ssp-Tier1 associÃ©s
+if [ ${sspT2} = "ssp119" ]
+then
+sspT1a=ssp126
+sspT1b=ssp126
+fi
+if [ ${sspT2} = "ssp434" ]
+then
+sspT1a=ssp126
+sspT1b=ssp245
+fi
+if [ ${sspT2} = "ssp460" ]
+then
+sspT1a=ssp245
+sspT1b=ssp370
+fi
+if [ ${sspT2} = "ssp534-over"Â ]
+then
+sspT1a=ssp126
+sspT1b=ssp245
+fi
+
+fileCNRMT1a=${ici}${sspT1a}_${year}_CNRM_hPa_interpURea.nc
+fileCNRMT1b=${ici}${sspT1b}_${year}_CNRM_hPa_interpURea.nc
+fileCNRMT2=${ici}${sspT2}_${year}_CNRM_hPa_interpURea.nc
+
+if [ ${year} -le 2049 ]
+then
+fileURead1aori=/prodigfs/project/input4MIPs/CMIP6/ScenarioMIP/UReading/UReading-CCMI-${sspT1a}-1-0/atmos/mon/vmro3/gn/v20181101/vmro3_input4MIPs_ozone_ScenarioMIP_UReading-CCMI-${sspT1a}-1-0_gn_201501-204912.nc
+fileURead1bori=/prodigfs/project/input4MIPs/CMIP6/ScenarioMIP/UReading/UReading-CCMI-${sspT1b}-1-0/atmos/mon/vmro3/gn/v20181101/vmro3_input4MIPs_ozone_ScenarioMIP_UReading-CCMI-${sspT1b}-1-0_gn_201501-204912.nc
+else
+fileURead1aori=/prodigfs/project/input4MIPs/CMIP6/ScenarioMIP/UReading/UReading-CCMI-${sspT1a}-1-0/atmos/mon/vmro3/gn/v20181101/vmro3_input4MIPs_ozone_ScenarioMIP_UReading-CCMI-${sspT1a}-1-0_gn_205001-209912.nc
+fileURead1bori=/prodigfs/project/input4MIPs/CMIP6/ScenarioMIP/UReading/UReading-CCMI-${sspT1b}-1-0/atmos/mon/vmro3/gn/v20181101/vmro3_input4MIPs_ozone_ScenarioMIP_UReading-CCMI-${sspT1b}-1-0_gn_205001-209912.nc
+fi
+
+fileUReadayear=${ici}o3_URead_${sspT1a}_${year}.nc
+fileUReadbyear=${ici}o3_URead_${sspT1b}_${year}.nc
+
+fileUReadatronq=${ici}o3_URead_${sspT1a}_${year}_tronq.nc
+fileUReadbtronq=${ici}o3_URead_${sspT1b}_${year}_tronq.nc
+
+# Fichier contenant les latitudes UReading exactes
+latURead=${ici}lat_URead.nc
+
+filevsa=${ici}o3_${sspT2}_${year}_vs_${sspT1a}.nc
+filevsb=${ici}o3_${sspT2}_${year}_vs_${sspT1b}.nc
+
+echo "... Selecting year "${year}" in Tier-1 UReading ssp "${sspT1a}"..."
+# Extraction de l'annÃ©e d'intÃ©rÃªt du ssp-Tier1 de UReading
+cdo selyear,${year} ${fileURead1aori} ${fileUReadayear}
+# Pas besoin de la seconde borne Tier1 si on est dans le cas ssp119 :
+if [ ${sspT1a} != ${sspT1b} ]
+then
+echo "... Selecting year "${year}" in Tier-1 UReading ssp "${sspT1b}"..."
+cdo selyear,${year} ${fileURead1bori} ${fileUReadbyear}
+fi
+
+rm -f ${latURead}
+ncks -s 'lat=double(lat)' ${fileUReadayear} ${latURead}
+
+echo "... Chopping up pressure levels, retaining only > 3 Pa..."
+# Troncature des niveaux verticaux
+cdo sellevel,1000,925,850,800,780,750,700,650,600,500,450,400,350,300,285,250,200,170,150,130,115,100,90,80,70,60,50,40,35,30,25,20,15,10,7,5,4,3,2,1.5,1,0.7,0.5,0.4,0.3,0.2,0.15,0.1,0.07,0.05,0.04,0.03 ${fileUReadayear} ${fileUReadatronq}
+# On vient d'enlever les 14 niveaux supÃ©rieurs du fichier UReading. 66-14 = 52 niveaux restent.
+# On fait la mÃªme chose pour la seconde borne Tier1, sauf dans le cas ssp119, oÃ¹ on n'en a pas besoin.
+if [ ${sspT1a} != ${sspT1b} ]
+then
+cdo sellevel,1000,925,850,800,780,750,700,650,600,500,450,400,350,300,285,250,200,170,150,130,115,100,90,80,70,60,50,40,35,30,25,20,15,10,7,5,4,3,2,1.5,1,0.7,0.5,0.4,0.3,0.2,0.15,0.1,0.07,0.05,0.04,0.03 ${fileUReadbyear} ${fileUReadbtronq}
+fi
+
+# Remappage horizontal des fichiers CNRM
+fileCNRMT1arem=${ici}o3_CNRM_${sspT1a}_${year}_remap.nc
+fileCNRMT1brem=${ici}o3_CNRM_${sspT1b}_${year}_remap.nc
+fileCNRMT2rem=${ici}o3_CNRM_${sspT2}_${year}_remap.nc
+
+echo "... Remapping CNRM file for "${sspT2}" and year "${year}"..."
+cdo remapcon,grid_URead_144x96 ${fileCNRMT2} ${fileCNRMT2rem}
+echo "... Remapping CNRM file for "${sspT1a}" and year "${year}"..."
+cdo remapcon,grid_URead_144x96 ${fileCNRMT1a} ${fileCNRMT1arem}
+rm -f ${fileCNRMT1a}
+if [ ${sspT1a} != ${sspT1b} ]
+then
+echo "... Remapping CNRM file for "${sspT1b}" and year "${year}"..."
+cdo remapcon,grid_URead_144x96 ${fileCNRMT1b} ${fileCNRMT1brem}
+fi
+
+# Fix des latitudes, qui ne tombent pas juste
+echo "... Fixing latitudes..."
+ncks -3 -A -v lat ${latURead} ${fileCNRMT2rem}
+ncks -3 -A -v lat ${latURead} ${fileCNRMT1arem}
+if [ ${sspT1a} != ${sspT1b} ]
+then
+ncks -3 -A -v lat ${latURead} ${fileCNRMT1brem}
+fi
+
+# La rÃšgle de trois en soi : l'ozone s'appelle "o3" dans les fichiers CNRM, et "vmro3" dans les fichiers UReading.
+echo "... Performing scaling : "${sspT2}" created by scaling against "${sspT1a}", for year "${year}"..."
+cdo -div -mul -selname,o3 ${fileCNRMT2rem} -selname,vmro3 ${fileUReadatronq} -selname,o3 ${fileCNRMT1arem} ${filevsa}
+if [ ${sspT1a} != ${sspT1b} ]
+then
+echo "... Performing scaling : "${sspT2}" created by scaling against "${sspT1b}", for year "${year}"..."
+cdo -div -mul -selname,o3 ${fileCNRMT2rem} -selname,vmro3 ${fileUReadbtronq} -selname,o3 ${fileCNRMT1brem} ${filevsb}
+fi
+
+# Fin des boucles ssp-Tier2
+done
+
+
+# PHASE 3
+# "CHAPEAUTAGE"
+
+if [ ${year} -le 2049 ]
+then
+fileURead1ori=/prodigfs/project/input4MIPs/CMIP6/ScenarioMIP/UReading/UReading-CCMI-${sspT1a}-1-0/atmos/mon/vmro3/gn/v20181101/vmro3_input4MIPs_ozone_ScenarioMIP_UReading-CCMI-${sspT1a}-1-0_gn_201501-204912.nc
+fileURead2ori=/prodigfs/project/input4MIPs/CMIP6/ScenarioMIP/UReading/UReading-CCMI-${sspT1b}-1-0/atmos/mon/vmro3/gn/v20181101/vmro3_input4MIPs_ozone_ScenarioMIP_UReading-CCMI-${sspT1b}-1-0_gn_201501-204912.nc
+else
+fileURead1ori=/prodigfs/project/input4MIPs/CMIP6/ScenarioMIP/UReading/UReading-CCMI-${sspT1a}-1-0/atmos/mon/vmro3/gn/v20181101/vmro3_input4MIPs_ozone_ScenarioMIP_UReading-CCMI-${sspT1a}-1-0_gn_205001-209912.nc
+fileURead2ori=/prodigfs/project/input4MIPs/CMIP6/ScenarioMIP/UReading/UReading-CCMI-${sspT1b}-1-0/atmos/mon/vmro3/gn/v20181101/vmro3_input4MIPs_ozone_ScenarioMIP_UReading-CCMI-${sspT1b}-1-0_gn_205001-209912.nc
+fi
+
+fileURead1year=${ici}o3_URead_${sspT1a}_${year}.nc
+fileURead2year=${ici}o3_URead_${sspT1b}_${year}.nc
+
+echo "... Preparing hat for year "${year}"..."
+if [ ${sspT1a} != ${sspT1b} ]
+then
+cdo selyear,${year} ${fileURead1ori} ${fileURead1year}
+cdo selyear,${year} ${fileURead2ori} ${fileURead2year}
+else
+cdo selyear,${year} ${fileURead1ori} ${fileURead1year}
+fi
+
+fileURead1chapo=${ici}o3_URead_${sspT1a}_${year}_chapo.nc
+fileURead2chapo=${ici}o3_URead_${sspT1b}_${year}_chapo.nc
+
+if [ ${sspT1a} != ${sspT1b} ]
+then
+cdo sellevel,0.02,0.015,0.01,0.007,0.005,0.004,0.003,0.002,0.0015,0.001,0.0008,0.0005,0.0003,0.0001 ${fileURead1year} ${fileURead1chapo}
+cdo sellevel,0.02,0.015,0.01,0.007,0.005,0.004,0.003,0.002,0.0015,0.001,0.0008,0.0005,0.0003,0.0001 ${fileURead2year} ${fileURead2chapo}
+else
+cdo sellevel,0.02,0.015,0.01,0.007,0.005,0.004,0.003,0.002,0.0015,0.001,0.0008,0.0005,0.0003,0.0001 ${fileURead1year} ${fileURead1chapo}
+fi
+
+fileUReadchapo=${ici}o3_URead_${year}_chapo.nc
+
+# Moyenne des chapeaux des deux ssp
+if [ ${sspT1a} != ${sspT1b} ]
+then
+echo "... Averaging hats for Tier 1 UReading "${sspT1a}" and "${sspT1b}", for "${year}"..."
+cdo mulc,.5 -add -selname,vmro3 ${fileURead1chapo} -selname,vmro3 ${fileURead2chapo} ${fileUReadchapo}
+else
+echo "... Using solely hat from Tier 1 UReading "${sspT1a}", for year "${year}"..."
+cp ${fileURead1chapo} ${fileUReadchapo}
+fi
+
+echo "... Cleaning up..."
+if [ ${sspT1a} != ${sspT1b} ]
+then
+rm -f ${fileURead1year} ${fileURead2year} ${fileURead1chapo} ${fileURead2chapo}
+else
+rm -f ${fileURead1year} ${fileURead1chapo}
+fi
+
+# Ajout du chapeau aux ssp Tier2 crÃ©Ã©s Ã  l'Ã©tape 2
+
+for sspT2 in "ssp119" "ssp434" "ssp460" "ssp534-over"
+do
+
+# choix des ssp-Tier1 associÃ©s
+if [ ${sspT2} = "ssp119" ]
+then
+sspT1a=ssp126
+sspT1b=ssp126
+fi
+if [ ${sspT2} = "ssp434" ]
+then
+sspT1a=ssp126
+sspT1b=ssp245
+fi
+if [ ${sspT2} = "ssp460" ]
+then
+sspT1a=ssp245
+sspT1b=ssp370
+fi
+if [ ${sspT2} = "ssp534-over"Â ]
+then
+sspT1a=ssp126
+sspT1b=ssp245
+fi
+
+filevsa=${ici}o3_${sspT2}_${year}_vs_${sspT1a}.nc
+filevsb=${ici}o3_${sspT2}_${year}_vs_${sspT1b}.nc
+
+fileUReadfinala=${repout}${sspT2}/o3_${sspT2}_${year}_vs_${sspT1a}.nc
+fileUReadfinalb=${repout}${sspT2}/o3_${sspT2}_${year}_vs_${sspT1b}.nc
+
+echo "... Renaming O3 variable and fixing NetCDF-3 inconsistency..."
+ncrename -v o3,vmro3 ${filevsa}
+ncks -3 ${filevsa} ${ici}tmp.nc
+mv ${ici}tmp.nc ${filevsa}
+if [ ${sspT1a} != ${sspT1b} ]
+then
+ncrename -v o3,vmro3 ${filevsb}
+ncks -3 ${filevsb} ${ici}tmp.nc
+mv ${ici}tmp.nc ${filevsb}
+fi
+
+echo "... Topping newly created "${sspT2}" (against "${sspT1a}") for year "${year}" with hat..."
+cdo merge ${filevsa} ${fileUReadchapo} ${fileUReadfinala}
+if [ ${sspT1a} != ${sspT1b} ]
+then
+echo "... Topping newly created "${sspT2}" (against "${sspT1b}") for year "${year}" with hat..."
+cdo merge ${filevsb} ${fileUReadchapo} ${fileUReadfinalb}
+fi
+
+# Je moyenne les deux versions
+echo "... Merging "${sspT2}" against "${sspT1a}" and "${sspT2}" against "${sspT1b}" into a single file..."
+filemerge=${repout}${sspT2}/o3_${sspT2}_${year}.nc
+if [ ${sspT1a} != ${sspT1b} ]
+then
+cdo add -mulc,.5 -selname,vmro3 ${fileUReadfinala} -mulc,.5 -selname,vmro3 ${fileUReadfinalb} ${filemerge}
+else
+cp ${fileUReadfinala} ${filemerge}
+fi
+
+# Petite entourloupe finale pour avoir prÃ©cisÃ©ment les niveaux de pression d'origine (entre-temps, on est passÃ©s en double avec petite erreur d'arrondi...)
+echo "... Fixing z-axis for year "${year}"..."
+cdo setzaxis,URead_axis_en_hPa ${filemerge} ${repout}tmp.nc
+mv ${repout}tmp.nc ${filemerge}
+
+# Fix des latitudes
+ncks -3 -A -v lat ${latURead} ${filemerge}
+
+# Fin de la boucle sur les sspT2
+done
+
+rm -f ${fileUReadchapo}
+
+# MÃ©nage
+rm -f ${ici}o3*${year}*.nc ${ici}ssp*${year}*.nc
+
+# Fin de la grande boucle sur les annÃ©es
+done
