source: lmdz_wrf/trunk/tools/obs-sim_Comparison.bash @ 1991

Last change on this file since 1991 was 1991, checked in by lfita, 6 years ago

Improvements and more Header oriented names of figures
  • Property svn:executable set to *
File size: 40.8 KB
Line 
1#!/bin/bash
2# Shell script to compare observations with simulations using PyNCplot
3#
4#   observations: different files from create_OBSnetcdf.py and/or UWyoming_snd_nc.py
5#   simulations: multiple outputs from different models and runs
6#
7
8# Name of the file with the configuration
9if test $1 = '-h'; then
10  echo "************************************************"
11  echo "***            Script to compare             ***"
12  echo "*** sounding and single-station observations ***"
13  echo "***         with multiple simulations        ***"
14  echo "************************************************"
15  echo "obs-sim_Comparison [ConfFile](configuration file)"
16else
17  rootsh=`pwd`
18
19  configfname=$1
20
21####### ###### ##### #### ### ## #
22
23function uploadvars() {
24# Function to upload variables to the system from an ASCII file as:
25#   [varname] = [value]
26  fileval=$1
27  errormsg='ERROR -- error -- ERROR -- error'
28
29  if test ! -f ${fileval}; then
30    echo ${errormsg}
31    echo "  "${fname}": file '"${fileval}"' does not exist!!"
32    exit
33  fi
34
35  Nlin=`wc -l ${fileval} | awk '{print $1}'`
36
37  ilin=1
38  while test ${ilin} -le ${Nlin}; do
39    line=`head -n ${ilin} ${fileval} | tail -n 1`
40    varname=`echo ${line} | tr '=' ' ' | awk '{print $1}'`
41    value=`echo ${line} | tr '=' ' ' | awk '{print $2}'`
42    Lvarname=`expr length ${varname}'0'`
43
44    if test ${Lvarname} -gt 1 && test ! ${varname:0:1} = '#'; then
45      export ${varname}=${value}
46    fi
47    ilin=`expr ${ilin} + 1`
48  done 
49}
50
51uploadvars ${configfname}
52echo "END upload -- end UPLOAD -- END upload -- end UPLOAD"
53
54# files from scratch
55if test; then ${fscratch} = 'true'; fscratch=true
56else fscratch=false; fi
57
58# figures from scratch
59if test; then ${gscratch} = 'true'; gscratch=true
60else gscratch=false; fi
61
62function num_hex(){
63# Fucntion to transform a number to a hexagonal basis
64#   num: number to transform
65  fname='num_hex'
66  num=$1
67
68  if test ${num} -gt 255; then
69    echo ${errormsg}
70    echo "  "${fname}": value "${num}" too large !!"
71    exit
72  fi
73  if test ${num} -gt 16; then
74    xviA=`expr ${num} / 16`
75    numA=`expr ${xviA} '*' 16`
76    xviB=`expr ${num} - ${numA}`
77  else
78    xviA=0
79    xviB=${num}
80  fi
81
82  if test ${xviA} -gt 9; then
83    case ${xviA} in
84      10)
85        hexA='A'
86      ;;
87      11)
88        hexA='B'
89      ;;
90      12)
91        hexA='C'
92      ;;
93      13)
94        hexA='D'
95      ;;
96      14)
97        hexA='E'
98      ;;
99      15)
100        hexA='F'
101      ;;
102    esac
103  else
104    hexA=${xviA}
105  fi
106
107  if test ${xviB} -gt 9; then
108    case ${xviB} in
109      10)
110        hexB='A'
111      ;;
112      11)
113        hexB='B'
114      ;;
115      12)
116        hexB='C'
117      ;;
118      13)
119        hexB='D'
120      ;;
121      14)
122        hexB='E'
123      ;;
124      15)
125        hexB='F'
126      ;;
127    esac
128  else
129    hexB=${xviB}
130  fi
131
132  echo ${hexA}${hexB}
133}
134
135function WindRose_plot() {
136# Function to plot a Wind rose
137#   it: time-step to use from the file
138#   sndstid: id of the sounding station
139#   sndstn: name of the sounding station to appear in plot ('!' for spaces)
140#   timeS: String with the format of the actual time ('!' for spaces)
141#   ofigdir: folder for the output of the figure
142#   timefS: date of the plot to be used for the figure name
143#   filen: file to use
144
145    valuesfig=$1
146
147    it=`echo ${valuesfig} | tr '#' ' ' | awk '{print $1}'`
148    sndstid=`echo ${valuesfig} | tr '#' ' ' | awk '{print $2}'`
149    sndstn=`echo ${valuesfig} | tr '#' ' ' | awk '{print $3}'`
150    timeS=`echo ${valuesfig} | tr '#' ' ' | awk '{print $4}'`
151    ofigdir=`echo ${valuesfig} | tr '#' ' ' | awk '{print $5}'`
152    timefS=`echo ${valuesfig} | tr '#' ' ' | awk '{print $6}'`
153    filen=`echo ${valuesfig} | tr '#' ' ' | awk '{print $7}'`
154
155    values='pres|-1;time|'${it}':linepoint;multicol;pres;auto;auto;rainbow;auto:'
156    values=${values}${sndstn}'!sounding!WindRose!on!'${timeS}'!local!time:'${kfig}
157    values=${values}':cardinals:False:WindRose:True'
158    ofign=${ofigdir}'/WindRose_'${sndstid}'_'${timefS}'.'${kfig}
159    if ${gscratch}; then rm ${ofign} >& /dev/null; fi
160    if test ! -f ${ofign}; then
161      echo "  Wind Rose on: "${timeS}
162      python ${pyHOME}/drawing.py -o draw_WindRose -S ${values} -v ua,va -f ${filen} >& /dev/null
163      if test $? -ne 0; then
164        echo ${errmsg}
165        echo "  python failed!!"
166        echo python ${pyHOME}/drawing.py -o draw_WindRose -S ${values} -v ua,va -f ${filen}
167        exit
168      fi
169      convert -trim WindRose.png ${ofign}
170      echo "* "${ofign} >> ${ofilefigs}
171      echo python ${pyHOME}/drawing.py -o draw_WindRose -S "'"${values}"'" -v ua,va -f "'"${filen}"'" >> ${ofilefigs}
172      echo " " >> ${ofilefigs}
173    fi
174}
175
176function SkewT_logP_plot() {
177# Function to plot a SkewT_logP plot
178#   it: time step to use from file
179#   sndstid: id of the sounding station
180#   sndstn: name of the sounding station to appear in plot ('!' for spaces)
181#   timeS: String with the format of the actual time ('!' for spaces)
182#   ofigdir: folder for the output of the figure
183#   timefS: date of the plot to be used for the figure name
184#   filen: file to use
185
186    valuesfig=$1
187
188    it=`echo ${valuesfig} | tr '#' ' ' | awk '{print $1}'`
189    sndstid=`echo ${valuesfig} | tr '#' ' ' | awk '{print $2}'`
190    sndstn=`echo ${valuesfig} | tr '#' ' ' | awk '{print $3}'`
191    timeS=`echo ${valuesfig} | tr '#' ' ' | awk '{print $4}'`
192    ofigdir=`echo ${valuesfig} | tr '#' ' ' | awk '{print $5}'`
193    timefS=`echo ${valuesfig} | tr '#' ' ' | awk '{print $6}'`
194    filen=`echo ${valuesfig} | tr '#' ' ' | awk '{print $7}'`
195 
196    values='time|'${it}',pres|-1:auto:auto:'${sndstn}'!sounding!on!'${timeS}'!local!'
197    values=${values}'time:'${kfig}':True'
198    ofign=${ofigdir}'/SkewT-logP_obs_ta-tda_'${sndstid}'_'${timefS}'.'${kfig}
199    if ${gscratch}; then rm ${ofign} >& /dev/null; fi
200    if test ! -f ${ofign}; then
201      echo "  Sounding on: "$(echo ${timeS} | tr '!' ' ')
202      python ${pyHOME}/drawing.py -o draw_SkewT -S ${values} -v ta,tda,pres -f ${filen} >& /dev/null
203      if test $? -ne 0; then
204        echo ${errmsg}
205        echo "  python failed!!"
206        echo python ${pyHOME}/drawing.py -o draw_SkewT -S ${values} -v ta,tda,pres -f ${filen}
207        exit
208      fi
209      convert -trim SkewT.png ${ofign}
210      echo "* "${ofign} >> ${ofilefigs}
211      echo python ${pyHOME}/drawing.py -o draw_SkewT -S "'"${values}"'" -v ta,tda,pres -f "'"${filen}"'" >> ${ofilefigs}
212      echo " " >> ${ofilefigs}
213    fi
214}
215
216function multi_soundings_plot() {
217# Function to plot a multi_sounding figure
218#   sndstid: id of the sounding station
219#   sndstn: name of the sounding station to appear in plot ('!' for spaces)
220#   tatda_evol_labs: labels for the figure
221#   tatda_evol_cols: colors of the lines for the figure
222#   timeS3: date in the figure
223#   ofigdir: folder for the output figure
224
225  valuesfig=$1
226
227  sndstid=`echo ${valuesfig} | tr '@' ' ' | awk '{print $1}'`
228  sndstn=`echo ${valuesfig} | tr '@' ' ' | awk '{print $2}'`
229  tatda_evol_labs=`echo ${valuesfig} | tr '@' ' ' | awk '{print $3}'`
230  tatda_evol_cols=`echo ${valuesfig} | tr '@' ' ' | awk '{print $4}'`
231  tatda_evol_files=`echo ${valuesfig} | tr '@' ' ' | awk '{print $5}'`
232  timeS3=`echo ${valuesfig} | tr '@' ' ' | awk '{print $6}'`
233  ofigdir=`echo ${valuesfig} | tr '@' ' ' | awk '{print $7}'`
234
235  tatda_evol_values='auto:auto:multilines!'${tatda_evol_labs}'!'${tatda_evol_cols}'!'
236  tatda_evol_values=${tatda_evol_values}'-!,!2:0,auto:'${ststn}'!sounding!evolution!on!'
237  tatda_evol_values=${tatda_evol_values}${timeS3}'!local!time:'${kfig}':True'
238
239  ofign=${ofigdir}'/SkewT-logP_obs_evol_'${stid}'.'${kfig}
240  if ${gscratch}; then rm ${ofign} >& /dev/null; fi
241  if test ! -f ${ofign}; then
242    echo "  Sounding on: "$(echo ${timeS} | tr '!' ' ')
243    python ${pyHOME}/drawing.py -o draw_multi_SkewT -S ${tatda_evol_values} -f ${tatda_evol_files} >& /dev/null
244    if test $? -ne 0; then
245      echo ${errmsg}
246      echo "  python failed!!"
247      python ${pyHOME}/drawing.py -o draw_multi_SkewT -S ${tatda_evol_values} -f ${tatda_evol_files}
248      exit
249    fi
250    convert -trim multi_SkewT.png ${ofign}
251    echo "* "${ofign} >> ${ofilefigs}
252    echo python ${pyHOME}/drawing.py -o draw_multi_SkewT -S ${tatda_evol_values} -f ${tatda_evol_files} >> ${ofilefigs}
253    echo " " >> ${ofilefigs}
254  fi
255}
256
257function multi_SkewT_logP_plot() {
258# Function to plot a multi_SkewT_logP plot with multiple lines
259#   sndstid: id of the sounding station
260#   sndstn: name of the sounding station to appear in plot ('!' for spaces)
261#   complabs: labels for the plot
262#   compcols: colors for the plot
263#   complines: type of lines for the plot
264#   compfiles: files to use for the plot
265#   dateS: String with the format of the actual time ('!' for spaces)
266#   ofigdir: folder for the output of the figure
267#   datefS: date of the plot to be used for the figure name
268
269    valuesfig=$1
270    sndstid=`echo ${valuesfig} | tr '@' ' ' | awk '{print $1}'`
271    sndstn=`echo ${valuesfig} | tr '@' ' ' | awk '{print $2}'`
272    complabs=`echo ${valuesfig} | tr '@' ' ' | awk '{print $3}'`
273    compcols=`echo ${valuesfig} | tr '@' ' ' | awk '{print $4}'`
274    complines=`echo ${valuesfig} | tr '@' ' ' | awk '{print $5}'`
275    compfiles=`echo ${valuesfig} | tr '@' ' ' | awk '{print $6}'`
276    dateS=`echo ${valuesfig} | tr '@' ' ' | awk '{print $7}'`
277    ofigdir=`echo ${valuesfig} | tr '@' ' ' | awk '{print $8}'`
278    datefS=`echo ${valuesfig} | tr '@' ' ' | awk '{print $9}'`
279
280    values='auto:auto:multilines!'${complabs}'!'${compcols}'!'
281    values=${values}${complines}'!,!2:0,auto:'${ststn}'!sounding!obs!,!sim!'
282    values=${values}'comparison!on!'${dateS}'!UTC:png:True'
283
284    ofign=${ofigdir}'/SkewT-logP_obs-sim_step_'${sndstid}'_'${datefS}'.'${kfig}
285    if ${gscratch}; then rm ${ofign} >& /dev/null; fi
286    if test ! -f ${ofign}; then
287      echo "  Sounding on: "$(echo ${dateS} | tr '!' ' ')
288      python ${pyHOME}/drawing.py -o draw_multi_SkewT -S ${values} -f ${compfiles} >& /dev/null
289      if test $? -ne 0; then
290        echo ${errmsg}
291        echo "  python failed!!"
292        echo python ${pyHOME}/drawing.py -o draw_multi_SkewT -S ${values} -f ${compfiles}
293        exit
294      fi
295      convert -trim multi_SkewT.png ${ofign}
296      echo "* "${ofign} >> ${ofilefigs}
297      echo python ${pyHOME}/drawing.py -o draw_multi_SkewT -S ${values} -f ${compfiles} >> ${ofilefigs}
298      echo " " >> ${ofilefigs}
299    fi
300}
301
302function Homoeller_SkewT_map_plot() {
303# Function to plot a Hovmoeller SkewT_map
304#   sndstid: id of the sounding station
305#   sndstn: name of the sounding station to appear in plot ('!' for spaces)
306#   sndv: variable to plot
307#   ofigdir: folder for the output of the figure
308#   datefS: date of the plot to be used for the figure name
309#   timevals: values to use for the time-axis in the plot
310#   simfilen: name of the simulation file at the given point
311#   sndfilen: name of the sounding file
312#   expl: label of the experiment
313#   expn: name of the experiment
314#   cbar: color bar
315
316  valuesfig=$1
317  sndstid=`echo ${valuesfig} | tr '#' ' ' | awk '{print $1}'`
318  sndstn=`echo ${valuesfig} | tr '#' ' ' | awk '{print $2}'`
319  sndv=`echo ${valuesfig} | tr '#' ' ' | awk '{print $3}'`
320  ofigdir=`echo ${valuesfig} | tr '#' ' ' | awk '{print $4}'`
321  datefS=`echo ${valuesfig} | tr '#' ' ' | awk '{print $5}'`
322  timevals=`echo ${valuesfig} | tr '#' ' ' | awk '{print $6}'`
323  simfilen=`echo ${valuesfig} | tr '#' ' ' | awk '{print $7}'`
324  sndfilen=`echo ${valuesfig} | tr '#' ' ' | awk '{print $8}'`
325  expl=`echo ${valuesfig} | tr '#' ' ' | awk '{print $9}'`
326  expn=`echo ${valuesfig} | tr '#' ' ' | awk '{print $10}'`
327  cbar=`echo ${valuesfig} | tr '#' ' ' | awk '{print $11}'`
328
329  cfiles=${simfilen}';'${sndv}';time;pres;time|-1,pres|-1@'
330  cfiles=${cfiles}${sndfilen}';'${sndv}';time;pres;time|-1,pres|-1'
331  cvalues=${sndv}';y;auto|'${timevals}';Vfix,auto,50.,auto'
332  cvalues=${cvalues}';'${cbar}',auto,auto;Srange,Srange;auto;obs!,!'${expn}'!'
333  cvalues=${cvalues}${sndstn}'!sounding;'${kfig}';flip@y;None;yes'
334
335  ofign=${ofigdir}'/SkewT-logP_obs-sim_evol_'${expl}'_'${sndstid}'_'${sndv}'.'${kfig}
336  if ${gscratch}; then rm ${ofign} >& /dev/null; fi
337  if test ! -f ${ofign}; then
338    python ${pyHOME}/drawing.py -o draw_2D_shad_contdisc_time -f ${cfiles} -S ${cvalues}
339    if test $? -ne 0; then
340      echo ${errmsg}
341      echo "  python failed!!"
342      echo python ${pyHOME}/drawing.py -o draw_2D_shad_contdisc_time -f ${cfiles} -S ${cvalues}
343      exit
344    fi
345    convert -trim 2Dshad_obs-sim_comparison_time.png ${ofign}
346    echo "* "${ofign} >> ${ofilefigs}
347    echo python ${pyHOME}/drawing.py -o draw_2D_shad_contdisc_time -f ${cfiles} -S ${cvalues} >> ${ofilefigs}
348    echo " " >> ${ofilefigs}
349  fi
350}
351
352function multiple_time_series_plot(){
353# Function to plot multiple time-series
354#   sti: id of the surface station
355#   stn: name of the surface station to appear in plot ('!' for spaces)
356#   sfcv: name of variable to plot
357#   ofigdir: folder for the output of the figure
358#   timevals: values of the time-axis to use
359#   files: files to use
360#   labs: labels of the legend
361#   lcol: color of lines to use
362#   nsfc: minimum value to plot
363#   xsfc: maximum value to plot
364
365    valuesfig=$1
366    stid=`echo ${valuesfig} | tr '@' ' ' | awk '{print $1}'`
367    stn=`echo ${valuesfig} | tr '@' ' ' | awk '{print $2}'`
368    sfcv=`echo ${valuesfig} | tr '@' ' ' | awk '{print $3}'`
369    ofigdir=`echo ${valuesfig} | tr '@' ' ' | awk '{print $4}'`
370    timevals=`echo ${valuesfig} | tr '@' ' ' | awk '{print $5}'`
371    files=`echo ${valuesfig} | tr '@' ' ' | awk '{print $6}'`
372    labs=`echo ${valuesfig} | tr '@' ' ' | awk '{print $7}'`
373    lcol=`echo ${valuesfig} | tr '@' ' ' | awk '{print $8}'`
374    nsfc=`echo ${valuesfig} | tr '@' ' ' | awk '{print $9}'`
375    xsfc=`echo ${valuesfig} | tr '@' ' ' | awk '{print $10}'`
376
377    timen=`echo ${timevals} | tr '|' ' ' | awk '{print $3}'`
378    fmtT=`echo ${timevals} | tr '|' ' ' | awk '{print $1"|"$2}'`
379
380    values='ststimes,time;y;'${timen}';auto;'${labs}';'${sfcv}';'${sfcv}
381    values=${values}'|evolution|at|'$(echo ${stn} | tr '!' '|')';'${nsfc}','${xsfc}';'
382    values=${values}'time|'${fmtT}';0|12;'
383    values=${values}${kfig}';-;'${lcol}';.;2.;2.;all;-1;True'
384
385    ofign=${ofigdir}'/obs-sim_evol_ts_'${sfcv}_${sti}'.'${kfig}
386    if ${gscratch}; then rm ${ofign} >& /dev/null; fi
387    if test ! -f ${ofign}; then
388      python ${pyHOME}/drawing.py -o draw_lines_time -f ${files} -S ${values} -v ${sfcv}
389      if test $? -ne 0; then
390        echo ${errmsg}
391        echo "  python failed!!"
392        python ${pyHOME}/drawing.py -o draw_lines_time -f ${files} -S ${values} -v ${sfcv}
393        exit
394      fi
395      convert -trim lines_time.${kfig} ${ofign}
396      echo "* "${ofign} >> ${ofilefigs}
397      echo python ${pyHOME}/drawing.py -o draw_lines_time -f ${files} -S ${values} -v ${sfcv} >> ${ofilefigs}
398      echo " "${ofign} >> ${ofilefigs}
399    fi
400}
401
402function shad_contdisc_map_plot(){
403# Function to plot a map of a continuos and a discontinuos variable
404#   expl: label of the experiment
405#   expn: name of the experiment in the title of the plot
406#   sfcv: name of variable to plot
407#   ofigdir: folder for the output of the figure
408#   obsfiles: observtional file to use
409#   simfile: simulation file to use
410#   it: time-step from simulation
411#   oit: time-step from observations
412#   timeS: format of time to appear in the title
413#   timefS: format of time to be used for the file name
414#   nsfc: minimum value to plot
415#   xsfc: maximum value to plot
416#   mapv: map values
417#   mapcover: cover of the map
418
419    valuesfig=$1
420
421    expl=`echo ${valuesfig} | tr '#' ' ' | awk '{print $1}'`
422    expn=`echo ${valuesfig} | tr '#' ' ' | awk '{print $2}'`
423    sfcv=`echo ${valuesfig} | tr '#' ' ' | awk '{print $3}'`
424    ofigdir=`echo ${valuesfig} | tr '#' ' ' | awk '{print $4}'`
425    obsfile=`echo ${valuesfig} | tr '#' ' ' | awk '{print $5}'`
426    simfile=`echo ${valuesfig} | tr '#' ' ' | awk '{print $6}'`
427    it=`echo ${valuesfig} | tr '#' ' ' | awk '{print $7}'`
428    oit=`echo ${valuesfig} | tr '#' ' ' | awk '{print $8}'`
429    timeS=`echo ${valuesfig} | tr '#' ' ' | awk '{print $9}'`
430    timefS=`echo ${valuesfig} | tr '#' ' ' | awk '{print $10}'`
431    nsfc=`echo ${valuesfig} | tr '#' ' ' | awk '{print $11}'`
432    xsfc=`echo ${valuesfig} | tr '#' ' ' | awk '{print $12}'`
433    mapv=`echo ${valuesfig} | tr '#' ' ' | awk '{print $13}'`
434    mapcover=`echo ${valuesfig} | tr '#' ' ' | awk '{print $14}'`
435
436    CFvarvals=`python $pyHOME/generic.py -o variables_values -S ${sfcv}`
437    cbar=`echo ${CFvarvals} | tr ':' ' ' | awk '{print $7}'`
438
439    cfiles=${simfile}';'${sfcv}';XLONG;XLAT;Time|'${it}',time|'${it}','
440    cfiles=${cfiles}'west_east|-1,south_north|-1@'${obsfile}';'${sfcv}
441    cfiles=${cfiles}';stslon;stslat;time|'${oit}',lon|-1,lat|-1'
442    cvalues=${sfcv}':west_east,south_north:auto:'${cbar}',auto,auto:'${nsfc}','
443    cvalues=${cvalues}${xsfc}':''auto:obs!,!'${expn}'!'${sfcv}'!on!'${timeS}'!UTC:'
444    cvalues=${cvalues}${kfig}':None:'${mapv}':'${mapcover}':yes'
445
446    ofign=${ofigdir}'/map_obs-sim_map_'${timefS}'_'${sfcv}'_'${expl}'.'${kfig}
447    if ${gscratch}; then rm ${ofign} >& /dev/null; fi
448    if test ! -f ${ofign}; then
449      python ${pyHOME}/drawing.py -o draw_2D_shad_contdisc -f ${cfiles} -S ${cvalues}
450      if test $? -ne 0; then
451        echo ${errmsg}
452        echo "  python failed!!"
453        echo python ${pyHOME}/drawing.py -o draw_2D_shad_contdisc -f ${cfiles} -S ${cvalues}
454      fi
455      convert -trim 2Dshad_obs-sim_comparison.png ${ofign} 
456      echo "* "${ofign} >> ${ofilefigs}
457      echo python ${pyHOME}/drawing.py -o draw_2D_shad_contdisc -f ${cfiles} -S ${cvalues} >> ${ofilefigs}
458      echo " " >> ${ofilefigs}
459    fi
460}
461
462#######    #######
463## MAIN
464    #######
465
466mkdir -p ${odir}
467
468ofilefigs=${wdir}/allins_figures.inf
469if ${fscratch}; then 
470  echo "*** Instruction from all figures" > ${ofilefigs}
471fi
472
473###
474## Preparing simulations
475###
476
477# sounding diags
478diagns=`echo ${wrfsnddiags} | tr ':' ' '`
479snddiagd='time@time,bottom_top@P,south_north@XLAT,west_east@XLONG'
480diagvals=''
481idiag=1
482for diagn in ${diagns}; do
483  diagv=`cat $pyHOME/diagnostics.inf | grep ${diagn}',' | tr ',' ' ' |               \
484    awk '{print $2"|"$3}'`
485  if test ${idiag} -eq 1; then
486    diagvals=${diagv}
487  else
488    diagvals=${diagvals}','${diagv}
489  fi
490  idiag=`expr ${idiag} + 1`
491done
492
493if test ${sfcvars} != 'None'; then
494# surface diags
495diagns=`echo ${wrfsfcdiags} | tr ':' ' '`
496sfcdiagd='time@time,south_north@XLAT,west_east@XLONG'
497diagvals=''
498idiag=1
499iadv=1
500for diagn in ${diagns}; do
501  diagv=`cat $pyHOME/diagnostics.inf | grep WRF | grep ${diagn}',' | tr ',' ' ' |    \
502    awk '{print $2"|"$3}'`
503  Ldiagv=`expr length $(echo ${diagv} | tr ' ' '@')0`
504  if test ${Ldiagv} -lt 2; then
505    echo ${errmsg}
506    echo "  diagnostic '"${diagn}"' (for WRF) does not exist !!"
507    if test ${iadv} -eq 1; then advs=${diagn}; iadv=`expr ${iadv} + 1`
508    else advs=${adv}','${diagn}; fi
509  else
510    if test ${idiag} -eq 1; then
511      diagvals=${diagv}
512      idiag=`expr ${idiag} + 1`
513    else
514      diagvals=${diagvals}','${diagv}
515    fi
516  fi
517done
518fi
519
520exps=`echo ${sims} | tr ':' ' '`
521
522iexp=1
523for exp in ${exps}; do
524  expf=`echo ${exp} | tr ',' ' ' | awk '{print $1}'`
525  expl=`echo ${exp} | tr ',' ' ' | awk '{print $2}'`
526  expn=`echo ${exp} | tr ',' ' ' | awk '{print $3}'`
527  expc=`echo ${exp} | tr ',' ' ' | awk '{print $4}'`
528
529  # Joining files
530  simjoinselvars=${odir}'/simout_vars_'${expl}'.nc'
531  if ${fscratch}; then rm ${simjoinselvars}; fi
532  if test ! -f ${simjoinselvars}; then
533    values=${simdir}'/'${expf}','${simtimedimn}','${simtimevarn}
534    HMT=${simHMT}
535    python ${pyHOME}/nc_var.py -o netcdf_fold_concatenation_HMT -S ${values}         \
536     -f ${HMT} -v ${simvars}
537    if test $? -ne 0; then
538      echo ${errmsg}
539      echo "  python failed!!"
540      echo python ${pyHOME}/nc_var.py -o netcdf_fold_concatenation_HMT -S ${values}  \
541        -f ${HMT} -v ${simvars}
542      exit
543    fi
544    mv netcdf_fold_concatenated_HMT.nc ${simjoinselvars}
545  fi
546
547  # Computing sounding diagnostics
548  simdiagsf=${odir}'/simout_snddiags_'${expl}'.nc'
549  if ${fscratch}; then rm ${simdiagsf}; fi
550  if test ! -f ${simdiagsf}; then
551    python ${pyHOME}/diagnostics.py -f ${simjoinselvars} -d ${snddiagd} -v ${diagvals}
552    if test $? -ne 0; then
553      echo ${errmsg}
554      echo "  python failed!!"
555      echo python ${pyHOME}/diagnostics.py -f ${simjoinselvars} -d ${snddiagd} -v ${diagvals}
556      exit
557    fi
558    mv diagnostics.nc ${simdiagsf}
559    #To deg
560    python ${pyHOME}/nc_var.py -o valmod -S subc,273.15 -f ${simdiagsf} -v ta
561    #To Pa
562    python ${pyHOME}/nc_var.py -o valmod -S divc,100. -f ${simdiagsf} -v pres
563  fi
564
565  if test ${sfcvars} != 'None'; then
566  # Computing surface diagnostics
567  simsfcdiagf=${odir}'/simout_sfcdiags_'${expl}'.nc'
568  if ${fscratch}; then rm ${simsfcdiagf}; fi
569  if test ! -f ${simsfcdiagf}; then
570    python ${pyHOME}/diagnostics.py -f ${simjoinselvars} -d ${sfcdiagd} -v ${diagvals}
571    if test $? -ne 0; then
572      echo ${errmsg}
573      echo "  python failed!!"
574      echo python ${pyHOME}/diagnostics.py -f ${simjoinselvars} -d ${sfcdiagd} -v ${diagvals}
575      exit
576    fi
577    mv diagnostics.nc ${simsfcdiagf}
578    #To rad
579    python ${pyHOME}/nc_var.py -o valmod -S mulc,0.0174532925199 -f ${simsfcdiagf} -v wds
580    python ${pyHOME}/nc_var.py -o varaddattr -S 'units|rad' -f ${simsfcdiagf} -v wds
581  fi
582
583  #Adding non-diagnostic variables
584  if test ${iadv} -ne 1; then
585    varsadd=`echo ${advs} | tr ',' ' '`
586    for vadd in ${varsadd}; do
587      python ${pyHOME}/nc_var.py -o fvaradd -S ${simjoinselvars},${vadd} -f ${simsfcdiagf}
588      if test $? -ne 0; then
589        echo ${errmsg}
590        echo "  python failed!!"
591        echo python ${pyHOME}/nc_var.py -o fvaradd -S ${simjoinselvars},${vadd} -f ${simsfcdiagf}
592        exit
593        rm ${simsfcdiagf}
594      fi
595      # CF varname
596      CFvarn=`python $pyHOME/generic.py -o variables_values -S ${vadd} | tr ':' ' ' |  \
597        awk '{print $1}'`
598      if test $? -ne 0; then
599        echo ${errmsg}
600        echo "  python failed!!"
601        echo python $pyHOME/generic.py -o variables_values -S ${vadd}
602        exit
603        rm ${simsfcdiagf}
604      fi
605      python $pyHOME/nc_var.py -o chvarname -S ${CFvarn} -f ${simsfcdiagf} -v ${vadd}
606      if test $? -ne 0; then
607        echo ${errmsg}
608        echo "  python failed!!"
609        echo python ${pyHOME}/nc_var.py -o chvarname -S ${CFvarn} -f ${simsfcdiagf} -v ${vadd}
610        exit
611        rm ${simsfcdiagf}
612      fi
613    done
614  fi
615  fi
616
617  iexp=`expr ${iexp} + 1`
618
619# End experiments
620done
621
622mkdir -p ${ofigdir}
623
624###
625## Soundings figures
626###
627
628# Files in sounding folder (from UWyoming_netcdf.py)
629sndfiles=`ls -1 ${snddir}/${sndHf}*`
630
631echo "Soundings files:"${sndfiles}
632
633isnd=1
634for sndf in ${sndfiles}; do
635  stid=`python $pyHOME/nc_var.py -o grattr -f ${sndf} -S Station_number`
636
637  # Getting wind components
638  ouavasndf=${odir}'/sounding_uava_'${stid}'.nc'
639  if ${fscratch}; then rm ${ouavasndf}; fi
640  if test ! -f ${ouavasndf}; then
641    python ${pyHOME}/diagnostics.py -f ${sndf} -d 'time@time,pres@pres' -v 'uavaFROMwswd|ws@wd'
642    if test $? -ne 0; then
643      echo ${errmsg}
644      echo "  python failed!!"
645      echo python ${pyHOME}/diagnostics.py -f ${sndf} -d 'time@time,pres@pres' -v 'uavaFROMwswd|ws@wd'
646      exit
647    fi 
648    mv diagnostics.nc ${ouavasndf}
649  fi
650  if test ${isnd} -eq 1; then sndstids=${stid}'@'${sndf}
651  else sndstids=${sndstids}':'${stid}'@'${sndf}; fi
652  isnd=`expr ${isnd} + 1`
653done
654
655sndids=`echo ${sndstids} | tr ':' ' '`
656
657for stidf in ${sndids}; do
658  stid=`echo ${stidf} | tr '@' ' ' | awk '{print $1}'`
659  sndorigfn=`echo ${stidf} | tr '@' ' ' | awk '{print $2}'`
660  echo ${stid}" ..."
661  it=0
662  tatda_evol_values=''
663  tatda_evol_files=''
664  sndtimes=''
665  snddates=''
666  sndFdates=''
667  snddimt=`python ${pyHOME}/nc_var.py -o itime -S CFtime -f ${sndorigfn} -v time | \
668     grep 'dimt:' | awk '{print $2-1}'`
669
670  # Time characteristics
671  tunits=`python ${pyHOME}/nc_var.py -o ivattrs -v time -f ${sndorigfn} | \
672    grep -v allvattrs | grep units | awk '{print $3}'`
673  Tu=`echo ${tunits} | tr '!' ' ' | awk '{print $1}'`
674  Dref0=`echo ${tunits} | tr '!' ' ' | awk '{print $3}'`
675  Dref=${Dref0:0:4}${Dref0:5:2}${Dref0:8:2}
676
677  # Getting values for evolution plot
678  col=`echo ${icolta} | tr ':' ' ' | awk '{print $1}'`
679  col2=`echo ${icolta} | tr ':' ' ' | awk '{print $2}'`
680  while test ${it} -le ${snddimt}; do
681    timev=`python $pyHOME/nc_var.py -o varout -S time:${it} -f ${sndorigfn} -v time | \
682      grep NC | awk '{printf ("%d", $2)}'`
683    timeS=`date +${figsndTfmt} -d"${Dref} ${timev} ${Tu}"`
684    timeS2=`date +%d"$^{"%H"}$" -d"${Dref} ${timev} ${Tu}"`
685    timeS3=`date +%Y"-"%b -d"${Dref} ${timev} ${Tu}"`
686    timeS4=`date +%Y -d"${Dref} ${timev} ${Tu}"`
687    echo "  "${it}" "${timev}" "${Dref}" '"${timev}"' "${Tu}" "$(echo ${timeS} | tr '!' ' ')
688    timefS=`date +%Y%m%d%H%M%S -d"${Dref} ${timev} ${Tu}"`
689
690    # for now...
691    sndstn=${stid}
692
693    WindRose_plot ${it}'#'${stid}'#'${sndstn}'#'${timeS}'#'${ofigdir}'#'${timefS}'#'${ouavasndf}
694
695    SkewT_logP_plot ${it}'#'${stid}'#'${sndstn}'#'${timeS}'#'${ofigdir}'#'${timefS}'#'${sndorigfn}
696
697    colH=`num_hex ${col}`
698    colH2=`num_hex ${col2}`
699    if test ${it} -eq 0; then
700      sndtimes=${timev}
701      snddates=${timeS}
702      sndFdates=${timefS}
703      tatda_evol_files=${sndorigfn}':pres|-1,time|'${it}':ta,pres;'
704      tatda_evol_files=${tatda_evol_files}${sndorigfn}':pres|-1,time|'${it}':tda,pres'
705      tatda_evol_labs='ta_'${timeS2}',tda'
706      tatda_evol_cols='#'${colH}${colH2}${colH2}',#'${colH2}${colH2}${colH}
707    else
708      sndtimes=${sndtimes}','${timev}
709      snddates=${snddates}','${timeS}
710      sndFdates=${sndFdates}','${timefS}
711      tatda_evol_files=${tatda_evol_files}';'${sndorigfn}':pres|-1,time|'${it}':ta,pres'
712      tatda_evol_files=${tatda_evol_files}';'${sndorigfn}':pres|-1,time|'${it}':tda,pres'
713      if test $(expr ${it} % 2) -eq 0; then 
714        tatda_evol_labs=${tatda_evol_labs}','${timeS2}',None'
715      else 
716        tatda_evol_labs=${tatda_evol_labs}',None,None'
717      fi
718      tatda_evol_cols=${tatda_evol_cols}',#'${colH}${colH2}${colH2}',#'${colH2}${colH2}${colH}
719    fi
720    col=`expr ${col} - ${ddcol}`
721    col2=`expr ${col2} - ${ddcol}`
722
723    #exit
724    it=`expr ${it} + 1`
725  done
726
727  vals=${stid}'@'${sndstn}'@'${tatda_evol_labs}'@'${tatda_evol_cols}'@'${tatda_evol_files}'@'${timeS3}'@'${ofigdir}
728  multi_soundings_plot ${vals}
729
730# End of soundings
731done
732
733# snd-sims comparison
734
735for stidf in ${sndids}; do
736  stid=`echo ${stidf} | tr '@' ' ' | awk '{print $1}'`
737  sndorigfn=`echo ${stidf} | tr '@' ' ' | awk '{print $2}'`
738  echo "snd: "${stid}
739
740  sndlon=`ncdump -h ${sndorigfn} | grep Station_longitude | awk '{print $3}'`
741  sndlat=`ncdump -h ${sndorigfn} | grep Station_latitude | awk '{print $3}'`
742  sndts=`echo ${sndtimes} | tr ',' ' '`
743
744  # for now...
745  sndstn=${stid}
746
747  # Time characteristics
748  tunits=`python ${pyHOME}/nc_var.py -o ivattrs -v time -f ${sndorigfn} |            \
749    grep -v allvattrs | grep units | awk '{print $3}'`
750  Tu=`echo ${tunits} | tr '!' ' ' | awk '{print $1}'`
751  Dref0=`echo ${tunits} | tr '!' ' ' | awk '{print $3}'`
752  Dref=${Dref0:0:4}${Dref0:5:2}${Dref0:8:2}
753
754  obsdimt=`python ${pyHOME}/nc_var.py -o itime -S CFtime, -v time -f ${sndorigfn} |  \
755    grep dimt | awk '{print $2-1}'`
756
757  # Plotting
758  it=0
759  while test ${it} -le ${obsdimt}; do
760
761    otimev=`python $pyHOME/nc_var.py -o varout -S time:${it} -f ${sndorigfn} -v time |\
762      grep NC | awk '{printf ("%d", $2)}'`
763    dateS=`date +${figsndTfmt} -d"${Dref} ${otimev} ${Tu}"`
764    timeS2=`date +%d"$^{"%H"}$" -d"${Dref} ${otimev} ${Tu}"`
765    timeS3=`date +%Y"-"%b -d"${Dref} ${otimev} ${Tu}"`
766    timeS4=`date +%Y -d"${Dref} ${otimev} ${Tu}"`
767    datefS=`date +%Y%m%d%H%M%S -d"${Dref} ${otimev} ${Tu}"`
768
769    iexp=1
770    for exp in ${exps}; do
771      expf=`echo ${exp} | tr ',' ' ' | awk '{print $1}'`
772      expl=`echo ${exp} | tr ',' ' ' | awk '{print $2}'`
773      expn=`echo ${exp} | tr ',' ' ' | awk '{print $3}'`
774      expc=`echo ${exp} | tr ',' ' ' | awk '{print $4}'`
775
776      simjoinselvars=${odir}'/simout_vars_'${expl}'.nc'
777      simdiagsf=${odir}'/simout_snddiags_'${expl}'.nc'
778
779      gridsndv=`python ${pyHOME}/nc_var.py -o get_point -f ${simjoinselvars}         \
780        -S 'XLONG:XLAT:Time|0' -v ${sndlon}','${sndlat} | tr ' ' '!'`
781      xsnd=`echo ${gridsndv} | tr '!' ' ' | awk '{print $1}' | tr ',' ' ' |          \
782        awk '{print $1}'`
783      ysnd=`echo ${gridsndv} | tr '!' ' ' | awk '{print $1}' | tr ',' ' ' |          \
784        awk '{print $2}'`
785      dsnd=`echo ${gridsndv} | tr '!' ' ' | awk '{print $2}'`
786      if test ${it} -eq 0 && test ${iexp} -eq 1; then
787        echo "    Equivalent simulated sounding grid point: "${gridsnd}" distance: "\
788          ${dsnd}
789      fi
790
791      # Getting values at the sounding point
792      simsndptf=${odir}'/simout_vars_sndpt_'${expl}'.nc'
793      if ${fscratch}; then rm ${simsndptf} >& /dev/null; fi
794      if test ! -f ${simsndptf}; then
795        values='time,0,-1,1@bottom_top,0,-1,1@south_north,'${ysnd}','${ysnd}',1@'
796        values=${values}'west_east,'${xsnd}','${xsnd}',1'
797        python ${pyHOME}/nc_var.py -o DataSetSection_multidims -f ${simdiagsf}       \
798          -S ${values} -v all
799        if test $? -ne 0; then
800          echo ${errmsg}
801          echo "  python failed!!"
802          echo python ${pyHOME}/nc_var.py -o DataSetSection_multidims -f ${simdiagsf}\
803            -S ${values} -v all
804          exit
805        fi
806        ofilen=`ls -rt1 ${odir} | tail -n 1`
807        mv ${odir}/${ofilen} ${simsndptf}
808      fi
809
810      # Time information from simulation
811      simdimt=`python ${pyHOME}/nc_var.py -o itime -S CFtime, -v time                \
812        -f ${simdiagsf} | grep dimt | awk '{print $2}'`
813      simtunits=`python ${pyHOME}/nc_var.py -o ivattrs -v time -f ${simdiagsf} |     \
814        grep -v allvattrs | grep units | awk '{print $3}'`
815      simTu=`echo ${tunits} | tr '!' ' ' | awk '{print $1}'`
816      simDref0=`echo ${tunits} | tr '!' ' ' | awk '{print $3}'`
817      simDref=${simDref0:0:4}${simDref0:5:2}${simDref0:8:2}
818
819      # Time step from simulated file
820      timestepv=`python ${pyHOME}/nc_var.py -o get_time -f ${simsndptf}               \
821        -S ${otimev}';'${tunits} -v time | tr ' ' '@'`
822      timestep=`echo ${timestepv} | tr '@' ' ' | awk '{print $1}'`
823      if test $? -ne 0; then
824        echo ${errmsg}
825        echo "  python failed!!"
826        echo python ${pyHOME}/nc_var.py -o get_time -f ${simsndptf}                  \
827          -S ${otimev}';'${tunits} -v time
828        exit
829      fi
830
831      if test ${iexp} -eq 1; then
832        echo "     "${it}": it="${timestep}" "$(echo ${dateS} | tr '!' ' ')" "${datefS}
833        tacomplabs='$ta^{obs}$,$ta^{'${expn}'}$'
834        tdacomplabs='$tda^{obs}$,$tda^{'${expn}'}$'
835        tacompcols='#FF6464,'${expc}
836        tdacompcols='#6464FF,'${expc}
837        tacompltyp='-,-'
838        tdacompltyp='-,-.'
839        tacompfiles=${sndorigfn}':pres|-1,time|'${it}':ta,pres;'
840        tacompfiles=${tacompfiles}${simsndptf}':bottom_top|-1,time|'${timestep}':'
841        tacompfiles=${tacompfiles}'ta,pres'
842        tdacompfiles=${sndorigfn}':pres|-1,time|'${it}':tda,pres;'
843        tdacompfiles=${tdacompfiles}${simsndptf}':bottom_top|-1,time|'${timestep}':'
844        tdacompfiles=${tdacompfiles}'tda,pres'
845      else
846        tacomplabs=${tacomplabs}',$ta^{'${expn}'}$'
847        tdacomplabs=${tdacomplabs}',$tda^{'${expn}'}$'
848        tacompcols=${tacompcols}','${expc}
849        tdacompcols=${tdacompcols}','${expc}
850        tacompltyp=${tacompltyp}',-'
851        tdacompltyp=${tdacompltyp}',-.'
852        tacompfiles=${tacompfiles}';'${simsndptf}':bottom_top|-1,time|'${timestep}':'
853        tacompfiles=${tacompfiles}'ta,pres'
854        tdacompfiles=${tdacompfiles}';'${simsndptf}':bottom_top|-1,time|'${timestep}
855        tdacompfiles=${tdacompfiles}':tda,pres'
856      fi
857
858      # Hovmoeller maps
859      sndvs=`echo ${sndvars} | tr ':' ' '`
860      for sndv in ${sndvs}; do
861
862        CFvarvals=`python $pyHOME/generic.py -o variables_values -S ${sndv}`
863        cbar=`echo ${CFvarvals} | tr ':' ' ' | awk '{print $7}'`
864
865        fivals=${stid}'#'${sndstn}'#'${sndv}'#'${ofigdir}'#'${datefS}'#'${fmtTts}
866        fivals=${fivals}'#'${simsndptf}'#'${sndorigfn}'#'${expl}'#'${expn}'#'${cbar}
867
868        Homoeller_SkewT_map_plot ${fivals}
869
870      # end sounding vars
871      done
872
873      iexp=`expr ${iexp} + 1`
874    # End experiments
875    done
876
877    complabs=${tacomplabs}','${tdacomplabs}
878    compcols=${tacompcols}','${tdacompcols}
879    compltyp=${tacompltyp}','${tdacompltyp}
880    compfiles=${tacompfiles}';'${tdacompfiles}
881
882    fivals=${stid}'@'${sndstn}'@'${complabs}'@'${compcols}'@'${compltyp}'@'
883    fivals=${fivals}${compfiles}'@'${dateS}'@'${ofigdir}'@'${datefS}'@'${ofigdir}
884
885    multi_SkewT_logP_plot ${fivals}
886
887    it=`expr ${it} + 1`
888  # End sounding time-steps
889  done
890
891# End of soundings stations
892done
893
894###
895## SFC figures
896###
897
898# Joining all observations
899obssfcfile=${odir}'/all_single-stations.nc'
900if test ! -f ${wdir}/${sfcobsdir}/${obssfcfile}; then
901  python ${pyHOME}/nc_var.py -o join_singlestation_obsfiles -S                       \
902    ${sfcobsdir}':OBSnetcdf' -v all
903  if test $? -ne 0; then
904    echo ${errmsg}
905    echo "  python failed!!"
906    echo  python ${pyHOME}/nc_var.py -o join_singlestation_obsfiles \
907      -S ${sfcobsdir}':OBSnetcdf' -v all
908  fi
909  mv joined_singlestations.nc ${obssfcfile}
910fi
911
912# getting min,max values
913if test ${sfcvars} != 'None'; then
914fcs=`echo ${sfcvars} | tr ':' ' '`
915
916iv=1
917for sfcpv in ${sfcs}; do
918  sfcv=`echo ${sfcpv} | tr '@' ' ' | awk '{print $1}'`
919  nsfc=`echo ${sfcpv} | tr '@' ' ' | awk '{print $2}'`
920  xsfc=`echo ${sfcpv} | tr '@' ' ' | awk '{print $3}'`
921  if test ${nsfc} = 'FROMobs'; then
922    fstats=`python $pyHOME/nc_var.py -o field_stats -S full,1.e+20,None              \
923      -f ${obssfcfile} -v ${sfcv} | grep MAT | grep ${sfcv}`
924    nsfc=`echo ${fstats} | awk '{print $3}'`
925  elif test ${nsfc} = 'FROMsims'; then
926    iexp=1
927    for exp in ${exps}; do
928      expl=`echo ${exp} | tr ',' ' ' | awk '{print $2}'`
929      simsfcdiagf=${odir}'/simout_sfcdiags_'${expl}'.nc'
930   
931      fstats=`python $pyHOME/nc_var.py -o field_stats -S full,1.e+20,None            \
932        -f ${simsfcdiagf} -v ${sfcv} | grep MAT | grep ${sfcv}`
933      if test ${iexp} -eq 1; then
934        nsfc=`echo ${fstats} | awk '{print $3}'`
935      else
936        nsfc2=`echo ${fstats} | awk '{print $3}'`
937        nsfcf=`echo ${nsfc}" "${nsfc2} | awk '{if ($1 < $2) {print $1;} else {print $2;}}'`
938      fi
939    done
940  elif test ${nsfc} = 'FROMobssims'; then
941    fstats=`python $pyHOME/nc_var.py -o field_stats -S full,1.e+20,None              \
942      -f ${obssfcfile} -v ${sfcv} | grep MAT | grep ${sfcv}`
943    nsfc=`echo ${fstats} | awk '{print $3}'`
944    for exp in ${exps}; do
945      expl=`echo ${exp} | tr ',' ' ' | awk '{print $2}'`
946      simsfcdiagf=${odir}'/simout_sfcdiags_'${expl}'.nc'
947   
948      fstats=`python $pyHOME/nc_var.py -o field_stats -S full,1.e+20,None            \
949        -f ${simsfcdiagf} -v ${sfcv} | grep MAT | grep ${sfcv}`
950      nsfc2=`echo ${fstats} | awk '{print $3}'`
951      nsfc=`echo ${nsfc}" "${nsfc2} | awk '{if ($1 < $2) {print $1;} else {print $2;}}'`
952    done
953  fi
954  if test ${xsfc} = 'FROMobs'; then
955    fstats=`python $pyHOME/nc_var.py -o field_stats -S full,1.e+20,None              \
956      -f ${obssfcfile} -v ${sfcv} | grep MAT | grep ${sfcv}`
957    xsfc=`echo ${fstats} | awk '{print $4}'`
958  elif test ${xsfc} = 'FROMsims'; then
959    iexp=1
960    for exp in ${exps}; do
961      expl=`echo ${exp} | tr ',' ' ' | awk '{print $2}'`
962      simsfcdiagf=${odir}'/simout_sfcdiags_'${expl}'.nc'
963   
964      fstats=`python $pyHOME/nc_var.py -o field_stats -S full,1.e+20,None            \
965        -f ${simsfcdiagf} -v ${sfcv} | grep MAT | grep ${sfcv}`
966      if test ${iexp} -eq 1; then
967        xsfc=`echo ${fstats} | awk '{print $4}'`
968      else
969        xsfc2=`echo ${fstats} | awk '{print $4}'`
970        xsfcf=`echo ${xsfc}" "${xsfc2} | awk '{if ($1 > $2) {print $1;} else {print $2;}}'`
971      fi
972    done
973  elif test ${nsfc} = 'FROMobssims'; then
974    fstats=`python $pyHOME/nc_var.py -o field_stats -S full,1.e+20,None              \
975      -f ${obssfcfile} -v ${sfcv} | grep MAT | grep ${sfcv}`
976    xsfc=`echo ${fstats} | awk '{print $4}'`
977    iexp=1
978    for exp in ${exps}; do
979      expl=`echo ${exp} | tr ',' ' ' | awk '{print $2}'`
980      simsfcdiagf=${odir}'/simout_sfcdiags_'${expl}'.nc'
981   
982      fstats=`python $pyHOME/nc_var.py -o field_stats -S full,1.e+20,None            \
983        -f ${simsfcdiagf} -v ${sfcv} | grep MAT | grep ${sfcv}`
984      xsfc2=`echo ${fstats} | awk '{print $4}'`
985      xsfcf=`echo ${xsfc}" "${xsfc2} | awk '{if ($1 > $2) {print $1;} else {print $2;}}'`
986    done
987  fi
988  if test $iv -eq 1; then
989    newsfcvars=${sfcv}'@'${nsfc}'@'${xsfc}
990  else
991    newsfcvars=${newsfcvars}':'${sfcv}'@'${nsfc}'@'${xsfc}
992  fi
993  iv=`expr ${iv} + 1`
994done
995sfcs=`echo ${newsfcvars} | tr ':' ' '`
996echo "new sfc variables: "${sfcs}
997
998# Time-series
999
1000# stations
1001# FROM: https://stackoverflow.com/questions/2961635/using-awk-to-print-all-columns-from-the-nth-to-the-last
1002stsids=`python ${pyHOME}/nc_var.py -o varout -v id -f \
1003  ${obssfcfile} -S 'Nstations:-1|time:0' | awk '{$1=""; printf("%d:",substr($0,2))}'`
1004stsnames=`python ${pyHOME}/nc_var.py -o varout -v stsname -f \
1005  ${obssfcfile} -S Nstations:-1 | awk '{$1=""; print substr($0,2)}' | tr '\n' ':' | \
1006  tr ' ' '!'`
1007stslons=`python ${pyHOME}/nc_var.py -o varout -v stslon -f \
1008  ${obssfcfile} -S Nstations:-1 | awk '{$1=""; print substr($0,2)}' | tr '\n' ':'`
1009stslats=`python ${pyHOME}/nc_var.py -o varout -v stslat -f \
1010  ${obssfcfile} -S Nstations:-1 | awk '{$1=""; print substr($0,2)}' | tr '\n' ':'`
1011
1012Nsts=`echo ${stsnames} | tr ':' ' ' | wc -w | awk '{print $1}'`
1013ist=1
1014while test ${ist} -le ${Nsts}; do
1015  sti=`echo ${stsids} | tr ':' '\n' | head -n ${ist} | tail -n 1`
1016  stn=`echo ${stsnames} | tr ':' '\n' | head -n ${ist} | tail -n 1`
1017  stl=`echo ${stslons} | tr ':' '\n' | head -n ${ist} | tail -n 1`
1018  stL=`echo ${stslats} | tr ':' '\n' | head -n ${ist} | tail -n 1`
1019  ist_1=`expr ${ist} - 1`
1020
1021  sfccompfiles=${obssfcfile}'%time|-1;Nstations|'${ist_1}
1022  sfccomplabs='obs'
1023  sfccompcols='#000000'
1024
1025  iexp=1
1026  for exp in ${exps}; do
1027    expl=`echo ${exp} | tr ',' ' ' | awk '{print $2}'`
1028    expn=`echo ${exp} | tr ',' ' ' | awk '{print $3}'`
1029    expc=`echo ${exp} | tr ',' ' ' | awk '{print $4}'`
1030
1031    simsfcdiagf=${odir}'/simout_sfcdiags_'${expl}'.nc'
1032
1033    gridsfcv=`python ${pyHOME}/nc_var.py -o get_point -f ${simsfcdiagf}              \
1034      -S 'XLONG:XLAT:Time|0' -v ${stl}','${stL} | tr ' ' '!'`
1035    xgrid=`echo ${gridsfcv} | tr '!' ' ' | awk '{print $1}' | tr ',' ' ' |           \
1036      awk '{print $1}'`
1037    ygrid=`echo ${gridsfcv} | tr '!' ' ' | awk '{print $1}' | tr ',' ' ' |           \
1038      awk '{print $2}'`
1039    distgrid=`echo ${gridsfcv} | tr '!' ' ' | awk '{print $2}'`
1040
1041    sfcstats=${sti}'@'${stn}'@'${stl}'@'${stL}'@'${gridsfc}
1042
1043    sfccompfiles=${sfccompfiles}','${simsfcdiagf}'%time|-1;south_north|'${ygrid}';west_east|'
1044    sfccompfiles=${sfccompfiles}${xgrid}
1045    sfccomplabs=${sfccomplabs}','${expn}
1046    sfccompcols=${sfccompcols}','${expc}
1047  done
1048  # end experiments
1049
1050  for sfcpv in ${sfcs}; do
1051    sfcv=`echo ${sfcpv} | tr '@' ' ' | awk '{print $1}'`
1052    nsfc=`echo ${sfcpv} | tr '@' ' ' | awk '{print $2}'`
1053    xsfc=`echo ${sfcpv} | tr '@' ' ' | awk '{print $3}'`
1054
1055    fivals=${sti}'@'${stn}'@'${sfcv}'@'${ofigdir}'@'${fmtTts}'@'${sfccompfiles}'@'
1056    fivals=${fivals}${sfccomplabs}'@'${sfccompcols}'@'${nsfc}'@'${xsfc}
1057
1058    multiple_time_series_plot ${fivals}
1059
1060  done
1061  # end observations
1062
1063  ist=`expr $ist + 1`
1064done
1065# end sfc statinos
1066
1067# Cont-disc comparison
1068iexp=1
1069for exp in ${exps}; do
1070  expl=`echo ${exp} | tr ',' ' ' | awk '{print $2}'`
1071  expn=`echo ${exp} | tr ',' ' ' | awk '{print $3}'`
1072  expc=`echo ${exp} | tr ',' ' ' | awk '{print $4}'`
1073
1074  simsfcdiagf=${odir}'/simout_sfcdiags_'${expl}'.nc'
1075
1076  # Time information from simulation
1077  simdimt=`python ${pyHOME}/nc_var.py -o itime -S CFtime, -v time -f ${simsfcdiagf}  \
1078    | grep dimt | awk '{print $2}'`
1079  simtunits=`python ${pyHOME}/nc_var.py -o ivattrs -v time -f ${simsfcdiagf} |       \
1080    grep -v allvattrs | grep units | awk '{print $3}'`
1081  simTu=`echo ${tunits} | tr '!' ' ' | awk '{print $1}'`
1082  simDref0=`echo ${tunits} | tr '!' ' ' | awk '{print $3}'`
1083  simDref=${simDref0:0:4}${simDref0:5:2}${simDref0:8:2}
1084
1085  it=0
1086  while test ${it} -le ${simdimt}; do
1087    simtstep=`python ${pyHOME}/nc_var.py -o varout -f ${simsfcdiagf} -S              \
1088      'time:'${it} -v time | grep NC | awk '{printf ("%d",$2)}'`
1089    if test $? -ne 0; then
1090      echo ${errmsg}
1091      echo "  python failed!!"
1092      echo python ${pyHOME}/nc_var.py -o varout -f ${simsfcdiagf} -S 'time:'${it} -v time
1093      exit
1094    fi
1095
1096    timeS=`date +${figsndTfmt} -d"${simDref} ${simtstep} ${simTu}"`
1097    echo ${it}" "${simtstep}" "${simDref}" "${simTu}" "$(echo ${timeS} | tr '!' ' ')
1098    timefS=`date +%Y%m%d%H%M%S -d"${simDref} ${simtstep} ${simTu}"`
1099
1100    # Looking for obs time
1101    oitv=`python ${pyHOME}/nc_var.py -o get_time -f ${obssfcfile}                    \
1102      -S ${simtstep}';'${tunits} -v ststimes | grep -v differ | tr ' ' '@'`
1103    oit=`echo ${oitv} | tr '@' ' ' | awk '{print $1}'`
1104
1105    for sfcnxv in ${sfcs}; do
1106      sfcv=`echo ${sfcnxv} | tr '@' ' ' | awk '{print $1}'`
1107      nsfc=`echo ${sfcnxv} | tr '@' ' ' | awk '{print $2}'`
1108      xsfc=`echo ${sfcnxv} | tr '@' ' ' | awk '{print $3}'`
1109
1110      fivals=${expl}'#'${expn}'#'${sfcv}'#'${ofigdir}'#'${obssfcfile}'#'${simsfcdiagf}
1111      fivals=${fivals}'#'${it}'#'${oit}'#'${timeS}'#'${timefS}'#'${nsfc}'#'${xsfc}'#'
1112      fivals=${fivals}${mapv}'#'${mapcover}
1113
1114      shad_contdisc_map_plot ${fivals}
1115
1116    # end of variables
1117    done
1118
1119  it=`expr ${it} + 1`
1120  # end of time-steps
1121  done
1122
1123# end of simulations
1124done
1125fi
1126
1127fi
Note: See TracBrowser for help on using the repository browser.