diff --git a/SCRIPT/mk_mht.bash b/SCRIPT/mk_mht.bash index 5dd7602..62faa56 100755 --- a/SCRIPT/mk_mht.bash +++ b/SCRIPT/mk_mht.bash @@ -39,4 +39,4 @@ else fi # extract only 26.5 -ncks -d y,793,793 $FILEOUT nemo_${RUN_NAME}o_${FREQ}_${TAG}_mht_265.nc +ncks -O -d y,793,793 $FILEOUT nemo_${RUN_NAME}o_${FREQ}_${TAG}_mht_265.nc diff --git a/SCRIPT/plot_time_series.py b/SCRIPT/plot_time_series.py index 94e9688..747cdfb 100644 --- a/SCRIPT/plot_time_series.py +++ b/SCRIPT/plot_time_series.py @@ -5,9 +5,11 @@ import matplotlib.colors as colors import sys import re +import datetime as dt import argparse import pandas as pd import matplotlib.dates as mdates +import matplotlib.ticker as ticker # def class runid class run(object): @@ -28,45 +30,49 @@ def load_time_series(self, cfile, cvar): for kf,cf in enumerate(cfile): try: ncid = nc.Dataset(cf) - except Exception as e: - print 'issue in trying to load file : '+cf - print e - sys.exit(42) - ncvtime = ncid.variables[ctime] - if 'units' in ncvtime.ncattrs(): - cunits = ncvtime.units - else: - cunits = "seconds since 1900-01-01 00:00:00" - - # define calendar - if 'calendar' in ncvtime.ncattrs(): - ccalendar = ncvtime.calendar - else: - ccalendar = "noleap" - time = nc.num2date(ncid.variables[ctime][:].squeeze(), cunits, ccalendar) + ncvtime = ncid.variables[ctime] + if 'units' in ncvtime.ncattrs(): + cunits = ncvtime.units + else: + cunits = "seconds since 1900-01-01 00:00:00" - # convert to proper datetime object - if isinstance(time,list): - ntime = time.shape[0] - else: - ntime = 1 - - timeidx=[None]*ntime - for itime in range(0, ntime): + # define calendar + if 'calendar' in ncvtime.ncattrs(): + ccalendar = ncvtime.calendar + else: + ccalendar = "noleap" + time = nc.num2date(ncid.variables[ctime][:].squeeze(), cunits, ccalendar) + + # convert to proper datetime object if isinstance(time,list): - timeidx[itime] = np.datetime64(time[itime],'D') + ntime = time.shape[0] else: - timeidx[itime] = np.datetime64(time,'D') + ntime = 1 - # build series - cnam=get_varname(cf,cvar) - df[kf] = pd.Series(ncid.variables[cnam][:].squeeze()*self.sf, index = timeidx, name = self.name) + timeidx=[None]*ntime + for itime in range(0, ntime): + if isinstance(time,list): + timeidx[itime] = np.datetime64(time[itime],'D') + else: + timeidx[itime] = np.datetime64(time,'D') + + # build series + cnam=get_varname(cf,cvar) + df[kf] = pd.Series(ncid.variables[cnam][:].squeeze()*self.sf, index = timeidx, name = self.name) + + except Exception as e: + print 'issue in trying to load file : '+cf + print e + sys.exit(42) + # build dataframe self.ts = pd.DataFrame(pd.concat(df)).sort_index() self.mean = self.ts[self.name].mean() self.std = self.ts[self.name].std() + self.min = self.ts[self.name].min() + self.max = self.ts[self.name].max() def __str__(self): return 'runid = {}, name = {}, line = {}, color = {}'.format(self.runid, self.name, self.line, self.color) @@ -108,9 +114,10 @@ def find_key(char, fid): def load_argument(): # deals with argument parser = argparse.ArgumentParser() - parser.add_argument("-runid", metavar='runid list' , help="used to look information in runid.db" , type=str, nargs='+' , required=True) + parser.add_argument("-runid", metavar='runid list' , help="used to look information in runid.db" , type=str, nargs='+' , required=True ) parser.add_argument("-f" , metavar='file list' , help="file list to plot (default is runid_var.nc)" , type=str, nargs='+' , required=False) - parser.add_argument("-var" , metavar='var list' , help="variable to look for in the netcdf file ./runid_var.nc", type=str, nargs='+' , required=True) + parser.add_argument("-var" , metavar='var list' , help="variable to look for in the netcdf file ./runid_var.nc", type=str, nargs='+' , required=True ) + parser.add_argument("-varf" , metavar='var list' , help="variable to look for in the netcdf file ./runid_var.nc", type=str, nargs='+' , required=False) parser.add_argument("-title", metavar='title' , help="subplot title (associated with var)" , type=str, nargs='+' , required=False) parser.add_argument("-dir" , metavar='directory of input file' , help="directory of input file" , type=str, nargs=1 , required=False, default=['./']) parser.add_argument("-sf" , metavar='scale factor', help="scale factor" , type=float, nargs=1 , required=False, default=[1]) @@ -225,6 +232,11 @@ def main(): plt.figure(figsize=np.array([210, 210]) / 25.4) # need to deal with multivar + mintime=dt.date.max + maxtime=dt.date.min + ymin=-sys.float_info.max + ymax=sys.float_info.max + for ivar, cvar in enumerate(args.var): ax[ivar] = plt.subplot(nvar, 1, ivar+1) # load obs @@ -246,6 +258,16 @@ def main(): if len(cfile)==0: print 'no file found with this pattern '+args.dir[0]+'/'+runid+'/'+fglob sys.exit(42) + elif args.varf: + # in case only one file pattern given + if len(args.varf) == 1 : + fglob = args.varf[0] + else : + fglob = args.varf[ivar] + cfile = glob.glob(args.dir[0]+'/'+runid+'/'+fglob) + if len(cfile)==0: + print 'no file found with this pattern '+args.dir[0]+'/'+runid+'/'+fglob + sys.exit(42) else: cfile = glob.glob(args.dir[0]+'/'+runid+'_'+cvar+'.nc') if len(cfile)==0: @@ -254,14 +276,32 @@ def main(): run_lst[irun].load_time_series(cfile, cvar) ts_lst[irun] = run_lst[irun].ts - lg = ts_lst[irun].plot(ax=ax[ivar], legend=False, style=run_lst[irun].line,color=run_lst[irun].color,label=run_lst[irun].name, linewidth=2, rot=0) + lg = ts_lst[irun].plot(ax=ax[ivar], legend=False, style=run_lst[irun].line,color=run_lst[irun].color,label=run_lst[irun].name, x_compat=True, linewidth=2, rot=0) + # + # limit of time axis + mintime=min([mintime,ts_lst[irun].index[0].to_pydatetime().date()]) + maxtime=max([maxtime,ts_lst[irun].index[-1].to_pydatetime().date()]) # set title if (args.title): ax[ivar].set_title(args.title[ivar],fontsize=20) - - # rm xaxis label - ax[ivar].xaxis.set_major_locator(mdates.YearLocator()) + + # set x axis + nlabel=5 + ndays=(maxtime-mintime).days + nyear=ndays/365 + if nyear < 10: + nyt=1 + elif 10<=nyear<50: + nyt=5 + elif 50<=nyear<100: + nyt=10 + else: + nyt=100 + nmt=ts_lst[irun].index[-1].to_pydatetime().date().month + ndt=ts_lst[irun].index[-1].to_pydatetime().date().day + + ax[ivar].xaxis.set_major_locator(mdates.YearLocator(nyt,month=nmt,day=1)) ax[ivar].tick_params(axis='both', labelsize=16) if (ivar != nvar-1): ax[ivar].set_xticklabels([]) @@ -270,9 +310,9 @@ def main(): for lt in ax[ivar].get_xticklabels(): lt.set_ha('center') -# ax[ivar].set_xticklabels(ha='center') rmin[ivar],rmax[ivar]=get_ybnd(run_lst,obs_min[ivar],obs_max[ivar]) + print rmin[ivar], rmax[ivar] ax[ivar].set_ylim([rmin[ivar],rmax[ivar]]) ax[ivar].grid() diff --git a/param.bash b/param.bash index 90c95ef..8e5153d 100644 --- a/param.bash +++ b/param.bash @@ -12,4 +12,33 @@ EXEPATH=${HOME}/VALSO/ # SCRIPT location SCRPATH=${HOME}/VALSO/SCRIPT/ +#RUNVALSO=0 +#RUNVALGLO=0 +#RUNALL=1 +# +#if [[ $RUNVALSO == 1 ]]; then +# runACC=1 #acc ts +# runMLD=1 #mld ts +# runBSF=1 #gyre ts +# runBOT=1 #bottom TS ts +#elif [[ $RUNVALGLO == 1 ]]; then +# runAMOC=1 +# runAMHT=1 +# runSIE=1 +# runSST=1 +#elif [[ $RUNALL == 1 ]]; then +# runACC=1 #acc ts +# runMLD=1 #mld ts +# runBSF=1 #gyre ts +# runBOT=1 #bottom TS ts +# runAMOC=1 +# runAMHT=1 +# runSIE=1 +# runSST=1 +#else +# echo 'need to define what you want in param.bash; exit 42' +# exit 42 +#fi + + module load gcc/8.1.0 mpi/mpich/3.2.1/gnu/8.1.0 hdf5/1.8.20/gnu/8.1.0 netcdf/4.6.1/gnu/8.1.0 diff --git a/run_all.bash b/run_all.bash index 206d530..d3bbb4c 100755 --- a/run_all.bash +++ b/run_all.bash @@ -28,6 +28,7 @@ for RUNID in `echo $RUNIDS`; do TAG=${YEAR}1201-$((YEAR+1))1201 TAG09=${YEAR}0901-${YEAR}1001 TAG02=${YEAR}0201-${YEAR}0301 + TAG02=${YEAR}0301-${YEAR}0401 # get data mooVyid=$(sbatch --job-name=moo_${YEAR}_V --output=${JOBOUT_PATH}/moo_${YEAR}_V ${SCRPATH}/get_data.bash $CONFIG $RUNID 1y $TAG grid-V | awk '{print $4}') # for mk_trp and mk_psi @@ -35,6 +36,7 @@ for RUNID in `echo $RUNIDS`; do mooTyid=$(sbatch --job-name=moo_${YEAR}_T --output=${JOBOUT_PATH}/moo_${YEAR}_T ${SCRPATH}/get_data.bash $CONFIG $RUNID 1y $TAG grid-T | awk '{print $4}') # for mk_bot.bash mooT09mid=$(sbatch --job-name=moo_${YEAR}_T --output=${JOBOUT_PATH}/moo_${YEAR}_T09 ${SCRPATH}/get_data.bash $CONFIG $RUNID 1m $TAG09 grid-T | awk '{print $4}') # for mk_mxl.bash mooT02mid=$(sbatch --job-name=moo_${YEAR}_T --output=${JOBOUT_PATH}/moo_${YEAR}_T02 ${SCRPATH}/get_data.bash $CONFIG $RUNID 1m $TAG02 grid-T | awk '{print $4}') # for mk_mxl.bash + mooT03mid=$(sbatch --job-name=moo_${YEAR}_T --output=${JOBOUT_PATH}/moo_${YEAR}_T03 ${SCRPATH}/get_data.bash $CONFIG $RUNID 1m $TAG03 grid-T | awk '{print $4}') # for mk_mxl.bash # run cdftools # scheduler option @@ -54,15 +56,15 @@ for RUNID in `echo $RUNIDS`; do # mld # VALSO - sbatchrunopt="--dependency=afterany:$mooTmid --job-name=SO_mxl_${TAG}_${RUNID} --output=${JOBOUT_PATH}/mxl_${TAG}.out" + sbatchrunopt="--dependency=afterany:$mooT09mid --job-name=SO_mxl_${TAG}_${RUNID} --output=${JOBOUT_PATH}/mxl_${TAG}.out" sbatch ${sbatchschopt} ${sbatchrunopt} ${SCRPATH}/mk_mxl.bash $CONFIG $RUNID $TAG09 1m > /dev/null 2>&1 & njob=$((njob+1)) # botT/S # VALSO - #sbatchrunopt="--dependency=afterany:$mooTyid --job-name=SO_bot_${TAG}_${RUNID} --output=${JOBOUT_PATH}/bot_${TAG}.out" - #sbatch ${sbatchschopt} ${sbatchrunopt} ${SCRPATH}/mk_bot.bash $CONFIG $RUNID $TAG 1y > /dev/null 2>&1 & - #njob=$((njob+1)) + sbatchrunopt="--dependency=afterany:$mooTyid --job-name=SO_bot_${TAG}_${RUNID} --output=${JOBOUT_PATH}/bot_${TAG}.out" + sbatch ${sbatchschopt} ${sbatchrunopt} ${SCRPATH}/mk_bot.bash $CONFIG $RUNID $TAG 1y > /dev/null 2>&1 & + njob=$((njob+1)) # moc # VALGLO @@ -94,6 +96,11 @@ for RUNID in `echo $RUNIDS`; do sbatch ${sbatchschopt} ${sbatchrunopt} ${SCRPATH}/mk_sie.bash $CONFIG $RUNID $TAG02 1m > /dev/null 2>&1 & njob=$((njob+1)) + # VALGLO/VALSO + sbatchrunopt="--dependency=afterany:$mooT03mid --job-name=GLO_sie_${TAG}_${RUNID} --output=${JOBOUT_PATH}/sie_${TAG03}.out" + sbatch ${sbatchschopt} ${sbatchrunopt} ${SCRPATH}/mk_sie.bash $CONFIG $RUNID $TAG03 1m > /dev/null 2>&1 & + njob=$((njob+1)) + # sie # VALGLO/VALSO sbatchrunopt="--dependency=afterany:$mooT09mid --job-name=GLO_sie_${TAG}_${RUNID} --output=${JOBOUT_PATH}/sie_${TAG09}.out" diff --git a/run_plot_VALGLO.bash b/run_plot_VALGLO.bash index 710d58b..d9f477a 100755 --- a/run_plot_VALGLO.bash +++ b/run_plot_VALGLO.bash @@ -1,5 +1,7 @@ #!/bin/bash +set -x + if [ $# -eq 0 ] ; then echo 'need a [KEYWORD] (will be inserted inside the figure title and output name) and a list of id [RUNIDS RUNID ...] (definition of line style need to be done in RUNID.db)'; exit; fi module load scitools/production-os41-1 @@ -48,12 +50,12 @@ convert ${KEY}_fig06.png -crop 1240x1040+0+0 tmp06.png # ARC SIE echo 'plot 02 SIE time series' -python2.7 SCRIPT/plot_time_series.py -noshow -runid $RUNIDS -f '*sie*0201*.nc' -var NExnsidc SExnsidc -sf 0.001 -title "SIE Arctic m02 [1e6 km2] : ${KEY}" "SIE Antarctic m02 [1e6 km2] : ${KEY}" -dir ${DATPATH} -o ${KEY}_fig07 -obs OBS/ARC_sie02_obs.txt OBS/ANT_sie02_obs.txt +python2.7 SCRIPT/plot_time_series.py -noshow -runid $RUNIDS -varf '*sie*0301*.nc' '*sie*0901*.nc' -var NExnsidc NExnsidc -sf 0.001 -title "SIE Arctic m03 [1e6 km2] : ${KEY}" "SIE Arctic m09 [1e6 km2] : ${KEY}" -dir ${DATPATH} -o ${KEY}_fig07 -obs OBS/ARC_sie03_obs.txt OBS/ARC_sie09_obs.txt if [[ $? -ne 0 ]]; then exit 42; fi convert ${KEY}_fig07.png -crop 1240x1040+0+0 tmp07.png echo 'plot 09 SIE time series' -python2.7 SCRIPT/plot_time_series.py -noshow -runid $RUNIDS -f '*sie*0901*.nc' -var NExnsidc SExnsidc -sf 0.001 -title "SIE Arctic m09 [1e6 km2] : ${KEY}" "SIE Antarctic m09 [1e6 km2] : ${KEY}" -dir ${DATPATH} -o ${KEY}_fig08 -obs OBS/ARC_sie09_obs.txt OBS/ANT_sie09_obs.txt +python2.7 SCRIPT/plot_time_series.py -noshow -runid $RUNIDS -varf '*sie*0201*.nc' '*sie*0901*.nc' -var SExnsidc SExnsidc -sf 0.001 -title "SIE Antarctic m02 [1e6 km2] : ${KEY}" "SIE Antarctic m09 [1e6 km2] : ${KEY}" -dir ${DATPATH} -o ${KEY}_fig08 -obs OBS/ANT_sie02_obs.txt OBS/ANT_sie09_obs.txt if [[ $? -ne 0 ]]; then exit 42; fi convert ${KEY}_fig08.png -crop 1240x1040+0+0 tmp08.png