#!/bin/bash -vx # Par exemple: # ./make_axe2.sh V6.5AIS8 # ----------------------------------------------------------------- local=`pwd` #source ~/env_Multi_atlas.sh #module purge #module load gcc/9.4.0 #module load netcdf-fortran/4.5.3-parallel-netcdf #module load cdo/2.0.6 #module load nco/5.0.1 #module load ncl/6.6.2 #module load anaconda3-py/2020.11 #module load ferret/7.6.0 source ${local}/env_Multi_atlas_axe2.sh set -vx pr_yr=1 pr_da=1 # ----------------------------------------------------------------- # Based on /home/fabric/LMDZ/UTILS/config.sh # ----------------------------------------------------------------- runstxt=~/${local}/runs.txt while test -n "${1}"; do case $1 in -runstxt) runstxt=$2 ; shift ;; -h) echo "Example: ./make_axe2.sh [-runstxt runs.description.file] V6.5AIS8" ;; *) comp=$1 ;; esac shift done WRK=$MULTIDIR/AXE2/WORK$$ if [ -d $WRK ] ; then WRK=$WRK$$ ; fi mkdir -p $WRK ; cd $WRK COMP_D=$MULTIDIR/$comp if [ ! -d $COMP_D ] ; then echo "$COMP_D not found" && exit ; fi # ----------------------------------------------------------------- # Liste des simulations de la comparasion $comp # ----------------------------------------------------------------- liste_simy="" DEF_FILE=$COMP_D/def.txt for s in `awk ' {print $1} ' $DEF_FILE` ; do liste_sim="$liste_sim $s" ; done for s in `awk ' {print $1"_"$2} ' $DEF_FILE` ; do liste_simy="$liste_simy $s" ; done echo liste_sim: $liste_sim echo liste_simy: $liste_simy # ----------------------------------------------------------------- # ----------------------------------------------------------------- # Analyse de la distribution longitudinale des precipitation # sur ocean et continent # ----------------------------------------------------------------- # ----------------------------------------------------------------- echo OK0 outputdir=$COMP_D/AXE2/PR_YR if [ -d $outputdir/PNG ] ; then pr_yr=0 ; fi if [ $pr_yr = 1 -o 0 = 0 ] ; then echo OK1 outputdir=$COMP_D/AXE2/PR_YR for sub in [ PNG NC TXT ] ; do mkdir -p $outputdir/$sub ; done echo OK2 $liste_sim listey="" for simy in $liste_simy ; do sim=`echo $simy | sed -e 's/_[0-9][0-9][0-9][0-9]_[0-9][0-9][0-9][0-9]//'` years=`echo $simy | sed -e 's/'$sim'_//'` simdir=`grep -w "^$sim " $runstxt | awk ' { print $2 } '` ncks -v precip $simdir/ATM/Analyse/SE/${sim}_SE_${years}_1M_histmth.nc -O grid.nc cdo timavg -remapcon,grid.nc ${MAIN_CMOR}/OBS/pr.nc gpcp.nc cdo timavg -remapcon,grid.nc ${MAIN}/TRMM/prm2000-2009.nc trmm.nc ###################################################### # Sauvegarde d'un masque oceanique pour le diag MJO ###################################################### if [ ! -f $WRK/pourc_oce.nc ] ; then ncks -d time_counter,1 -v pourc_oce "$simdir/ATM/Analyse/SE/${sim}_SE_${years}_1M_histmth.nc" $WRK/pourc_oce.nc ; fi echo ncks -d time_counter,1 -v pourc_oce "$simdir/ATM/Analyse/SE/${sim}_SE_${years}_1M_histmth.nc" $WRK/pourc_oce.nc cat <| tmpyr.jnl use "$simdir/ATM/Analyse/SE/${sim}_SE_${years}_1M_histmth.nc" use gpcp.nc use trmm.nc let con=pourc_ter[d=1,l=1:12@ave]+pourc_lic[d=1,l=1:12@ave] let oce=pourc_oce[d=1,l=1:12@ave]+pourc_sic[d=1,l=1:12@ave] let ocem=if ( oce[l=1] ge 99. ) then 1 let conm=if ( con[l=1] ge 99. ) then 1 let mixm=if ( oce[l=1] gt 1. AND con[l=1] gt 1. ) then 1 let totm=tsol[d=1,l=1]*0.+1. let filtre=var go tmp0.jnl TOT,8,totm frame/file=tot.gif let filtre=if ( oce[l=1] ge 99.9 ) then var*oce/100. go tmp0.jnl OCE,9,ocem frame/file=oce.gif let filtre=if ( con[l=1] ge 99.9 ) then var*con/100. go tmp0.jnl CONT,10,conm frame/file=con.gif let filtre=if ( con[l=1] gt 0.1 AND oce[l=1] gt 0.1 ) then var go tmp0.jnl MIXT,11,mixm frame/file=mix.gif list/x=-10:10/y=12:18 "PRlmdzAMMA",86400*precip[d=1,i=@ave,j=@ave,l=1:12@ave] list/x=-10:10/y=12:18 "PRtrmmAMMA",r[d=3,i=@ave,j=@ave,l=1:12@ave] list/x=-10:10/y=12:18 "PRgpcpAMMA",86400*pr[d=2,i=@ave,j=@ave,l=1:12@ave] eod cat <| tmp0.jnl let mask=\$3 let var=86400*precip[d=1]*mask[i=@ave,j=@ave] plot/line=\$2/vlim=0:10/title="$sim $years pr \$1 (mm/d)" filtre[i=@ave,l=1:12@ave] list/y=-50:50 "PRlmdz\$1",filtre[i=@ave,j=@ave,l=1:12@ave] let var=r[d=3]*mask plot/line=1/o/title="TRMM" filtre[i=@ave,l=1] list/y=-50:50 "PRtrmm\$1",filtre[i=@ave,j=@ave,l=1:12@ave] let var=86400*pr[d=2]*mask plot/line=7/o/title="GPCP" filtre[i=@ave,l=1:12@ave] list/y=-50:50 "PRgpcp\$1",filtre[i=@ave,j=@ave,l=1:12@ave] quit eod \rm gotmpyr.jnl cat<gotmpyr.jnl go tmpyr.jnl quit eof echo ferret $outputdir/TXT/${sim}_${years} ferret $outputdir/TXT/${sim}_${years} for type in con mix tot oce ; do convert -density 144 $type.gif $outputdir/PNG/${sim}_${years}_$type.png listey="$listey ${sim}_${years}_$type" done done if [ -f $COMP_D/entete.html ] ; then cat $COMP_D/entete.html > $outputdir/index.html ; fi ~/Multi_atlas/concat_html.sh $COMP_D/AXE2/PR_YR/PNG "OK" "$listey" 4 >> $outputdir/index.html fi # pr_yr=1 # ----------------------------------------------------------------- # ----------------------------------------------------------------- # Analyse de la variabilité jour à jour des pluies # ----------------------------------------------------------------- # ----------------------------------------------------------------- outputdir=$COMP_D/AXE2/PR_DAY # if [ -d $outputdir ] ; then pr_da = 0 ; fi if [ $pr_da = 1 ] ; then listey="" outputdir=$COMP_D/AXE2/PR_DAY for sub in [ PNG NC TXT ] ; do mkdir -p $outputdir/$sub ; done # PLOT SIMULATION RESULTS first=1 for simy in $liste_simy ; do sim=`echo $simy | sed -e 's/_[0-9][0-9][0-9][0-9]_[0-9][0-9][0-9][0-9]//'` years=`echo $simy | sed -e 's/'$sim'_//'` yi=`echo $years | cut -d_ -f1` yf=`echo $years | cut -d_ -f2` echo FICHIERS POUR LA SIMULATION $sim annees $years simdir=`grep -w "^$sim " $runstxt | awk ' { print $2 } '` run=`basename $simdir ` prd=`ls $simdir/ATM/Analyse/TS_DA/${run}_${yi}*_${yf}*_1D_precip.nc | tail -1` prc=`ls $simdir/ATM/Analyse/TS_DA/${run}_${yi}*_${yf}*_1D_pluc.nc | tail -1` prl=`ls $simdir/ATM/Analyse/TS_DA/${run}_${yi}*_${yf}*_1D_plul.nc | tail -1` ln -s $prd $outputdir/NC/${run}_${yi}_${yf}.nc prm=${MAIN_SE}/ORIG/${run}_SE_${years}_1M_histmth.nc if [ ! -f ${prm} ]; then prm=$simdir/ATM/Analyse/SE/${run}_SE_${years}_1M_histmth.nc fi if [ $first = 1 ] ; then remap=$prm ; first=0 trmm=$outputdir/NC/TRMM.nc if [ ! -f $trmm ] ; then ncks -v precip $prm -O $outputdir/NC/grid.nc ###cdo remapcon,$outputdir/NC/grid.nc /thredds/ipsl/fabric/lmdz/TRMM/pr2000-2009.nc $trmm cdo remapcon,$outputdir/NC/grid.nc ${MAIN}/TRMM/pr2000-2009.nc $trmm fi ts_da="TRMM" fi echo PRD $prd if [ -f "$prd" ] ; then ts_da="$ts_da ${run}_${yi}_${yf}" ; fi done for sim in $ts_da ; do echo Tracer de $sim if [ $sim = TRMM ] ; then var=r else var='86400*precip' fi nom=`echo $sim | sed -e 's/_/ /g'` cat <| tmp.jnl use "$WRK/pourc_oce.nc" use "$outputdir/NC/$sim.nc" DEFINE VIEWPORT/XLIM=0.,1./YLIM= .5,1. V1 DEFINE VIEWPORT/XLIM=0.,1./YLIM= 0.2,0.7 V2 DEFINE VIEWPORT/XLIM=0.,1./YLIM= 0.,0.35 V3 set memory/size=1000 let var=$var let dv=var-var[l=@sbx:30] let dv2=dv*dv let dv2oce=if ( pourc_oce[d=1,l=1] ge 99. ) then dv2 reg/y=-50:50 fill/pal=rain_cmyk/nolab/set/lev=(3,15,3)(Inf) (dv2[l=@ave])^0.5 list "PRd1-30",(dv2[i=@ave,j=@ave,l=@ave])^0.5 ! ppl axlsze,0.14,0.14 ppl fill label/nouser 3.,6.7,0,0,0.24 "$nom" label/nouser 3.,6.2,0,0,0.24 "Std deviation, pr(1d)-pr(30d)" go land frame/file=hf.gif let dv=var[l=@sbx:20]-var[l=@sbx:120] fill/pal=rain_cmyk/nolab/set/lev=(1,5,1)(Inf) (dv2[l=@ave])^0.5 list "PRd20-120",(dv2[i=@ave,j=@ave,l=@ave])^0.5 list/x=50:180/y=-20:20 "PRd20-120MJO",(dv2oce[i=@ave,j=@ave,l=@ave])^0.5 ppl fill go land label/nouser 3.,6.7,0,0,0.24 "$nom" label/nouser 3,6.2,0,0,0.24 "Std deviation, pr(20d)-pr(120d)" frame/file=is.gif let pr1= if ( $var gt 1 ) then 1 else 0 fill/pal=rain_cmyk/nolab/set/lev=(0.1,1,0.15) pr1[l=@ave] list "PRnd1",pr1[i=@ave,j=@ave,l=@ave] ppl fill label/nouser 3.,6.7,0,0,0.24 "$nom" label/nouser 3.,6.2,0,0,0.24 "time frac. wth daily pr. above 1mm/day" go land frame/file=nd.gif plot/line=7/vlim=0:9/title="$nom, PR (mm/day)" ${var}[i=@ave,l=@ave] let pr10= if ( $var gt 10 ) then $var else 0 let pr20= if ( $var gt 20 ) then $var else 0 let pr50= if ( $var gt 50 ) then $var else 0 let pr100= if ( $var gt 100 ) then $var else 0 let pr200= if ( $var gt 200 ) then $var else 0 plot/o/line=8/title="PR>10mm/d" pr10[i=@ave,l=@ave] plot/o/line=9/title="PR>20mm/d" pr20[i=@ave,l=@ave] plot/o/line=10/title="PR>50mm/d" pr50[i=@ave,l=@ave] plot/o/line=11/title="PR>100mm/d" pr100[i=@ave,l=@ave] plot/o/line=12/title="PR>200mm/d" pr200[i=@ave,l=@ave] frame/file=zon.gif let dpr10=$var-pr10 let dpr20=pr10-pr20 let dpr50=pr20-pr50 ! list "PRTOT",${var}[i=@ave,j=@ave,l=@ave] list "PR0-10",dpr10[i=@ave,j=@ave,l=@ave] list "PR10-20",dpr20[i=@ave,j=@ave,l=@ave] list "PR20-50",dpr50[i=@ave,j=@ave,l=@ave] list "PR50",pr50[i=@ave,j=@ave,l=@ave] quit eod \rm gotmp.jnl cat<gotmp.jnl go tmp.jnl quit eof ferret $outputdir/TXT/$sim for type in hf is nd zon ; do convert $type.gif $outputdir/PNG/${sim}_$type.png listey="$listey ${sim}_$type" done done echo $ts_da if [ -f $COMP_D/entete.html ] ; then cat $COMP_D/entete.html > $outputdir/index.html ; fi ~/Multi_atlas/concat_html.sh $COMP_D/AXE2/PR_DAY/PNG "OK" "$listey" 4 >> $outputdir/index.html fi # pr_da=1 ####################################################################### # Bar chart bilan ####################################################################### cd $COMP_D/AXE2 liste_simy="" DEF_FILE=$COMP_D/def.txt for s in `awk ' {print $1"_"$2 } ' $DEF_FILE` ; do liste_simy="$liste_simy $s" ; done echo $liste_simy ref=1 mkdir -p XMGR for sim in $liste_simy ; do if [ $ref = 1 ] ; then clims="lmdz gpcp trmm" else clims=lmdz fi ref=0 for clim in $clims ; do \rm -f $clim ; for type in TOT OCE CONT MIXT AMMA ; do grep $type PR_YR/TXT/$sim | grep $clim | awk ' { print $2 } ' >> $clim ; done ; done for clim in $clims ; do for type in TOT OCE CONT MIXT AMMA ; do grep $type PR_YR/TXT/$sim | grep $clim | awk ' { print $2 } ' ; done ; done \mv -f lmdz XMGR/$sim awk ' { print $2 } ' PR_DAY/TXT/$sim >> XMGR/$sim done mv trmm XMGR/TRMM mv gpcp XMGR/GPCP awk ' { print $2 } ' PR_DAY/TXT/TRMM >> XMGR/TRMM awk ' { print $2 * 0 } ' PR_DAY/TXT/TRMM >> XMGR/GPCP cd XMGR nsims=`echo $liste_simy | wc -w | awk ' { print $1 + 2 } '` size=`echo $nsims | awk ' { print 3.5 / $1 } '` cat <| precip.param g0 type Chart world -0.5, 0, 12.5, 8 yaxis label "PR (mm/day) " yaxis tick major 1 xaxis ticklabel font 0 xaxis ticklabel color 1 xaxis tick place both xaxis tick spec type both xaxis tick spec 13 xaxis tick major 0, 0 xaxis ticklabel 0, "TOT" xaxis tick major 1, 1 xaxis ticklabel 1, "OCE" xaxis tick major 2, 2 xaxis ticklabel 2, "CONT" xaxis tick major 3, 3 xaxis ticklabel 3, "MIXT" xaxis tick major 4, 4 xaxis ticklabel 4, "AMMA" xaxis tick major 5, 5 xaxis ticklabel 5, "SIG1-30" xaxis tick major 6, 6 xaxis ticklabel 6, "SIG20-120" xaxis tick major 7, 7 xaxis ticklabel 7, "MJO" xaxis tick major 8, 8 xaxis ticklabel 8, "FD(PR>1)" xaxis tick major 9, 9 xaxis ticklabel 9, "<10" xaxis tick major 10, 10 xaxis ticklabel 10, "10-20" xaxis tick major 11, 11 xaxis ticklabel 11, "20-50" xaxis tick major 12, 12 xaxis ticklabel 12, ">50" xaxis tick major 13, 13 xaxis ticklabel char size 1.3 xaxis ticklabel angle 90 legend 0.8, 0.85 eod color=( 1 1 1 2 3 4 5 6 8 9 10 11 12 13 14 15 ) pattern=( 10 6 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 ) is=0 while [ $is != $nsims ] ; do cat <> precip.param s$is type bar s$is symbol fill pattern ${pattern[$is]} s$is symbol color ${color[$is]} s$is symbol fill color ${color[$is]} s$is symbol size $size s$is line linestyle 0 eod is=$(( $is + 1 )) done xmgrace TRMM GPCP $liste_simy -param precip.param -legend load -hardcopy -hdevice EPS -printfile tmp.eps epstopdf tmp.eps convert -density 144 tmp.pdf tmp.png