#!/bin/sh set -x cd $DATA set +x ######################################################################### echo "----------------------------------------------------------------" echo "hysplit_conc.sh - Script To Run Hysplit Concentrations and plot." echo "----------------------------------------------------------------" # REVISED 19 May 2005 - hourly output cump2 for volcanoes # 17 Oct 2005 - number particles to MAPTEXT.CFG; MAXPAR # 19 Jun 2006 - 2nd conc grid for RSMC # 01 Dec 2006 - zoom and relcen veriables # 31 Jan 2008 - Google Earth # 01 Feb 2008 - set end sample date/time to zeros # source attribution (backward) # 28 Mar 2008 - lower MAXPAR; change filename cdump1 to cdump # 01 Apr 2008 - change issued to created # 09 Apr 2008 - force units Bq in MAPTEXT for RSMC # 21 Oct 2009 - GDAS-GFS labeling ######################################################################### set -x ################################## # typenum=1 RSMC # typenum=2 Volcano # typenum=3 Homeland Security # typenum=4 Back-tracking (Source Attribution) ################################## ############################################################# # Set Output Name for Dispersion/Concentration Simulation ############################################################# #NOTE FOR HS run: # There will be 2 postscript files created (concplot1 and concplot2) # concplot1 will contain the fine grid concentration plots, while # concplot2 will contain the coarse grid concentration plots # If the GIS option (-a) is turned on (1) for concplot then output # GIS text files will also be created. Each concentration level will # have the following format.... # # GIS_level1-level2_(fine or coarse)_(count).txt # # These files will be read by Arcview and converted to shapefiles #NOTE FOR volcanic ash run: # There will be 9 postscript files created # volcplot06.r0 (ash reduction - none) 6-hour output # volcplot06.r1 (ash reduction - small) # volcplot06.r2 (ash reduction - medium) # volcplot06.r3 (ash reduction - large) # volcplot12.r0 12-hour output for wafs # volcplot12.r1 # volcplot12.r2 # volcplot12.r3 # volcplotae.r0 (alert-ended, 12-hour output for wafs, .rN is required) #NOTE FOR rsmc run: # There will be 1 postscript file created (concplot.ps) ############################################################# # run duration above which should run 3-D particle model # --automatically set when set run duration, and checked at end of setup_hysplit_script # --longrun and neglongrun must be identical in setup_hysplit_script # -- and hysplit_conc.sh # --test is: if [ ${nhrs} -ge ${longrun} -o ${nhrs} -le ${neglongrun} ] longrun=84 neglongrun=`echo "${longrun}*-1" | bc` ############################################################ # Set Up Control File for each type of run ############################################################ # common grid specs ngrd=1 if [ $dhr2 -ne 0 ]; then ngrd=2; fi # grid #1 gctr[1]="$cnlat $cnlon" # center gspc[1]="$dlat $dlon" # spacing gspn[1]="$splat $splon" # span gdir[1]="$DATA/" # directory # expand grid for long-range runs if [ ${nhrs} -ge ${longrun} -o ${nhrs} -le ${neglongrun} ]; then gspn[1]="181. 360." # span fi # grid #2 gctr[2]="$cnlat2 $cnlon2" # center gspc[2]="$dlat2 $dlon2" # spacing gspn[2]="$splat2 $splon2" # span gdir[2]="$DATA/" # directory gfil[2]="cdump2" # filename case $typenum in 1) # RSMC - particle with wet and dry deposition and half-life # 1 concentration grid npol=1 pol[1]=$ident qrate[1]=$qrat # grid specs gfil[1]="cdump" # filename gnlv[1]=2 # number of levels glvl[1]="0 $height" # levels year=`echo $csyr | cut -c3-` # RSMC protocol - samples start at 00z or 12z if [ $dshr -ge 12 ] then gsst[1]="$year $csmo $csda 12 00" else gsst[1]="$year $csmo $csda 00 00" fi gsnd[1]="00 00 00 00 00" # sample stop gint[1]="00 $dhr 00" # sampling interval if [ $ngrd -eq 2 ]; then # number of hours of output (from present time) dfcst=18 # round rlse down to previous multiple of $dhr2 as start time (3 for RSMC) # release start time YR=$dsyr MO=$dsmo DA=$dsda HR=$dshr # sample start time = release start time rounded down if [ ! `expr $HR % $dhr2` -eq 0 ]; then HR=`expr $HR - $HR % $dhr2` if [ $HR -lt 10 ]; then HR="0"$HR; fi fi # present time PYR=`date -u +%Y` PMO=`date -u +%m` PDA=`date -u +%d` PHR=`date -u +%H` # sample end time = (present + 18 hours) rounded up if [ `expr $PHR % $dhr2` -eq 0 ]; then inc=$dfcst NS=`expr $dfcst / $dhr2` # number of samples else inc=`expr $PHR % $dhr2` inc=`expr $dfcst + $dhr2 - $inc` NS=`expr $dfcst / $dhr2 + 1` fi END=`${utilexec}/ndate $inc ${PYR}${PMO}${PDA}${PHR}` EYR=`echo $END | cut -c1-4` EMO=`echo $END | cut -c5-6` EDA=`echo $END | cut -c7-8` EHR=`echo $END | cut -c9-10` # time periods (samples) N1 to N2 will be plotted NHR=`${utilexec}/nhour ${EYR}${EMO}${EDA}${EHR} ${YR}${MO}${DA}${HR}` N2=`expr $NHR / $dhr2` N1=`expr $N2 - $NS + 1` ######################################################################## gnlv[2]=2 # number of levels glvl[2]="0 $height2" # levels gsst[2]="$YR $MO $DA $HR 00" # sample start gsnd[2]="$EYR $EMO $EDA $EHR 00" # sample stop gint[2]="00 $dhr2 00" # sampling interval fi # depositing specs pspc[1]="1.0 1.0 1.0" # particle specs vdep[1]="$dryvl 0.0 0.0 0.0 0.0" # deposition velocity etc. wet[1]="0.0 $wetin $wetlo" # wet removal half[1]=$rhalf # half-life res[1]="0.0" # resuspension ;; 2) # Volcano - 4 particle sizes, with gravitational settling # 2 concentration grids; 1 6-hourly, 2nd 1-hourly # sample start time and run time depend on eruption time # with respect to meteo initialization time npol=4 pol[1]=p006; qrate[1]=0.008 # 0.3 - 1 microns pol[2]=p020; qrate[2]=0.068 # 1 - 3 microns pol[3]=p060; qrate[3]=0.25 # 3 - 10 microns pol[4]=p200; qrate[4]=0.67 # 10 - 30 microns # total emission is 1.0 for i in 1 2 3 4 do qrate[${i}]=`echo "scale=5\n ${qrate[${i}]}/${qhrs}" | bc` done # grid specs for conc grid #1 gfil[1]="cdump" # filename gnlv[1]=3 # number of levels glvl[1]="6096 10668 16764" # levels (needed for volcplot) gsst[1]="00 00 00 00 00" # sample start specified below gsnd[1]="00 00 00 00 00" # sample stop gint[1]="-1 6 00" # sampling interval - hourly concentration every 6-h # start time is met fcst initialization time #YR=`echo ${PDY} | cut -c1-4` #MO=`echo ${PDY} | cut -c5-6` #DA=`echo ${PDY} | cut -c7-8` #HR=${cyc} # grid specs for conc grid #2 (for SAB NAWIPS) # start time is the time of the eruption, if up to 6 hours ago, # or, is 6 hours ago, if the eruption was more than 6 hours ago # present time YR=`date -u +%Y` MO=`date -u +%m` DA=`date -u +%d` HR=`date -u +%H` # 6 hours before present back=6 inc=`expr -1 \* $back` SIXBP=`${utilexec}/ndate $inc ${YR}${MO}${DA}${HR}` syr=`echo $SIXBP | cut -c1-4` smo=`echo $SIXBP | cut -c5-6` sda=`echo $SIXBP | cut -c7-8` shr=`echo $SIXBP | cut -c9-10` # eruption start time # $dsyr # $dsmo # $dsda # $dshr # desired number of forecast hours for hourly output dfcst=18 # extra hours to have if needed extra=3 # dt = difference in hours between present and eruption start # dt>0 ==> eruption before present dt=`${utilexec}/nhour ${YR}${MO}${DA}${HR} $dsyr$dsmo$dsda$dshr` let dt=${dt} # to get rid of leading "0" # dhour = duration of cdump2 sampling if [ ${dt} -ge ${back} ]; then # eruption starts at least X hours before present(BP), start hourlies 6h BP # grid #2 sample start time YR=$syr MO=$smo DA=$sda HR=$shr # grid #2 sample end time let dhour=${dfcst}+${extra}+${back} elif [ ${dt} -ge 0 ]; then # eruption started at the present or up to 6h BP, start hourlies at eruption # grid #2 sample start time YR=$dsyr MO=$dsmo DA=$dsda HR=$dshr # grid #2 sample end time let dhour=${dfcst}+${extra}+${dt} else # don't do hourly plots for eruption starting after present (not real erupt) # ngrd=1 # OR set start time to eruption time and output hourlies for 18 hours only # grid #2 sample start time YR=$dsyr MO=$dsmo DA=$dsda HR=$dshr # grid #2 sample end time dhour=${dfcst} fi # grid #2 sample end time echo $dhour >> DHOUR SMPET=`${utilexec}/ndate $dhour ${YR}${MO}${DA}${HR}` EYR=`echo $SMPET | cut -c1-4` EMO=`echo $SMPET | cut -c5-6` EDA=`echo $SMPET | cut -c7-8` EHR=`echo $SMPET | cut -c9-10` gnlv[2]=3 # number of levels glvl[2]="6096 10668 16764" # levels (want same levels in NAWIPS) gsst[2]="$YR $MO $DA $HR 00" # sample start gsnd[2]="$EYR $EMO $EDA $EHR 00" # sample stop gint[2]="-1 $dhr2 00" # sampling interval ########################################### # sample start time and run time - volcano ########################################### # difference between eruption time (Te) and met fcst init time (Tm) in hours # DT was calculated in hysplit_trajconc.sh # positive for eruption after Tm # negative before ## don't want output plots during analysis period # if Te is at or after 3 hours before Tm if test $DT -ge -3 then # sample start time = eruption time syr=`echo $dsyr | cut -c3-4` smo=$dsmo sda=$dsda shr=$dshr else # if Te is at or after 18-h before Tm and before 3-h before Tm if test $DT -ge -18 -a $DT -lt -3 then # sample start time = Te + 12 inc=12 SMPST=`${utilexec}/ndate $inc $dsyr$dsmo$dsda$dshr` syr=`echo $SMPST | cut -c3-4` smo=`echo $SMPST | cut -c5-6` sda=`echo $SMPST | cut -c7-8` shr=`echo $SMPST | cut -c9-10` elif [ $DT -lt -174 ]; then # could handle earlier eruptions by increasing values for N # as long as have meteorology to support it msg='ERROR Sample start time not configured for these runs' errexit "$msg" else for N in 3 4 5 6 7 8 9 10 11 12 13 14 15 do # for next upgrade - consider simplifying this by just # setting syr/smo/sda/shr to the current time # N=3: -30 --> -18 hours # N=4: -42 --> -30 hours # N=5: -54 --> -42 hours # N=6: -66 --> -54 hours # N=7: -78 --> -66 hours # N=8: -90 --> -78 hours # N=9: -102 --> -90 hours # N=10:-114 -->-102 hours # N=11:-126 -->-114 hours # N=12:-138 -->-126 hours # N=13:-150 -->-138 hours # N=14:-162 -->-150 hours # N=15:-174 -->-162 hours lo=`echo "-12*$N+6" | bc` hi=`echo "-12*($N-1)+6" | bc` if test $DT -ge $lo -a $DT -lt $hi then inc=`echo "12*($N-1)" | bc` SMPST=`${utilexec}/ndate $inc $dsyr$dsmo$dsda$dshr` syr=`echo $SMPST | cut -c3-4` smo=`echo $SMPST | cut -c5-6` sda=`echo $SMPST | cut -c7-8` shr=`echo $SMPST | cut -c9-10` fi done fi fi gsst[1]="$syr $smo $sda $shr 00" # sample start # depositing specs pspc[1]="0.6 2.5 1.0" # particle specs (size in microns, dens, shape) vdep[1]="0.0 0.0 0.0 0.0 0.0" # deposition velocity etc. wet[1]="0.0 0.0 0.0" # wet removal half[1]="0.0" # half-life res[1]="0.0" # resuspension pspc[2]="2.0 2.5 1.0" vdep[2]=${vdep[1]} wet[2]=${wet[1]} half[2]=${half[1]} res[2]=${res[1]} pspc[3]="6.0 2.5 1.0" vdep[3]=${vdep[1]} wet[3]=${wet[1]} half[3]=${half[1]} res[3]=${res[1]} pspc[4]="20.0 2.5 1.0" vdep[4]=${vdep[1]} wet[4]=${wet[1]} half[4]=${half[1]} res[4]=${res[1]} ;; 3) # Homeland Security - gaseous pollutant, non-depositing # 2 concentration grids: 1st 12-hours, after 12-h npol=1 pol[1]=$ident qrate[1]=$qrat # ------------------------------------ # fine grid output 1st $delhr hours # ------------------------------------ # grid specs gfil[1]="cdump" # filename gnlv[1]=2 # number of levels glvl[1]="10 $height" # levels year=`echo $csyr | cut -c3-` gsst[1]="$year $csmo $csda $cshr 00" # sample start gsnd[1]="00 00 00 $delhr 00" # sample end gint[1]="00 $dhr 00" # sampling interval # ----------------------------------------- # coarse grid output beyond $delhr hours # ----------------------------------------- if [ $ngrd -eq 2 ]; then gnlv[2]=2 # number of levels glvl[2]="10 $height2" # levels gsst[2]="00 00 00 $delhr 00" # sample start gsnd[2]="00 00 00 00 00" # sample end gint[2]="00 $dhr2 00" # sampling interval fi # depositing specs if [ ${wetin} -eq 0 -a ${wetlo} -eq 0 ] then pspc[1]="0.0 0.0 0.0" # particle specs 0's --> gas else pspc[1]="1.0 1.0 1.0" # particle specs 1's --> particle fi vdep[1]="$dryvl 0.0 0.0 0.0 0.0" # deposition velocity etc. wet[1]="0.0 $wetin $wetlo" # wet removal half[1]="$rhalf" # half-life res[1]="0.0" # resuspension ;; 4) # Back-tracking (Source Attribution) # no deposition, 1 concentration grid, gas (not particle) npol=1 pol[1]=$ident qrate[1]=$qrat # grid specs gfil[1]="cdump" # filename gnlv[1]=1 # number of levels glvl[1]="$height" # levels year=`echo $csyr | cut -c3-` gsst[1]="$year $csmo $csda $cshr 00" gsnd[1]="00 00 00 00 00" # sample stop gint[1]="00 $dhr 00" # sampling interval # number of hours of output (from present time) dfcst=18 # round rlse down to previous multiple of $dhr2 as start time (3 for RSMC) # release start time YR=$dsyr MO=$dsmo DA=$dsda HR=$dshr # sample start time = release start time rounded down #if [ ! `expr $HR % $dhr2` -eq 0 ]; then # HR=`expr $HR - $HR % $dhr2` # if [ $HR -lt 10 ]; then HR="0"$HR; fi #fi # present time PYR=`date -u +%Y` PMO=`date -u +%m` PDA=`date -u +%d` PHR=`date -u +%H` # sample end time = (present + 18 hours) rounded up #if [ `expr $PHR % $dhr2` -eq 0 ]; then # inc=$dfcst # NS=`expr $dfcst / $dhr2` # number of samples #else # inc=`expr $PHR % $dhr2` # inc=`expr $dfcst + $dhr2 - $inc` # NS=`expr $dfcst / $dhr2 + 1` #fi # END=`${utilexec}/ndate $inc ${PYR}${PMO}${PDA}${PHR}` # EYR=`echo $END | cut -c1-4` # EMO=`echo $END | cut -c5-6` # EDA=`echo $END | cut -c7-8` # EHR=`echo $END | cut -c9-10` # force to zeros for now - need to calc based on request/avail met EYR=00 EMO=00 EDA=00 EHR=00 # time periods (samples) N1 to N2 will be plotted #NHR=`${utilexec}/nhour ${EYR}${EMO}${EDA}${EHR} ${YR}${MO}${DA}${HR}` #N2=`expr $NHR / $dhr2` #N1=`expr $N2 - $NS + 1` ######################################################################## gnlv[2]=2 # number of levels glvl[2]="0 $height2" # levels gsst[2]="$YR $MO $DA $HR 00" # sample start gsnd[2]="$EYR $EMO $EDA $EHR 00" # sample stop gint[2]="00 $dhr2 00" # sampling interval # depositing specs pspc[1]="0.0 0.0 0.0" # particle specs vdep[1]="$dryvl 0.0 0.0 0.0 0.0" # deposition velocity etc. wet[1]="0.0 $wetin $wetlo" # wet removal half[1]=$rhalf # half-life res[1]="0.0" # resuspension ;; esac ######################## # write CONTROL file (this is what hysplit reads) ######################## echo "$ibyr $ibmo $ibda $ibhr" >CONTROL.conc # calculation start echo "2 ">>CONTROL.conc # number of starting locations echo "$olat $olon ${olvl1%.*}">>CONTROL.conc # loc #1 lat lon height echo "$olat $olon ${olvl2%.*}">>CONTROL.conc # loc #2 lat lon height echo "$nhrs ">>CONTROL.conc # total run time echo "0 ">>CONTROL.conc # vertical motion echo "$zdata ">>CONTROL.conc # top of model domain echo "$nfile" >>CONTROL.conc # number of input meteo files # meteo files may appear in any order # hysplit will choose appropriate files in order # see MESSAGE file # meteo files must exist ifile=0 if [ ! "$met1" = "undefined" ] then let ifile=$ifile+1 echo "$DATA/ ">>CONTROL.conc # file #1 directory echo "${met1} ">>CONTROL.conc # file #1 filename fi if [ ! "$met2" = "undefined" ] then let ifile=$ifile+1 echo "$DATA/ ">>CONTROL.conc # file #2 directory echo "${met2} ">>CONTROL.conc # file #2 filename fi if [ ! "$met3" = "undefined" ] then let ifile=$ifile+1 echo "$DATA/ ">>CONTROL.conc # file #3 directory echo "${met3} ">>CONTROL.conc # file #3 filename fi if [ ! "$met4" = "undefined" ] then let ifile=$ifile+1 echo "$DATA/ ">>CONTROL.conc # file #4 directory echo "${met4} ">>CONTROL.conc # file #4 filename fi if [ ! "$met5" = "undefined" ] then let ifile=$ifile+1 echo "$DATA/ ">>CONTROL.conc # file #5 directory echo "${met5} ">>CONTROL.conc # file #5 filename fi if [ ! "$met6" = "undefined" ] then let ifile=$ifile+1 echo "$DATA/ ">>CONTROL.conc # file #5 directory echo "${met6} ">>CONTROL.conc # file #5 filename fi if [ ! "$met7" = "undefined" ] then let ifile=$ifile+1 echo "$DATA/ ">>CONTROL.conc # file #5 directory echo "${met7} ">>CONTROL.conc # file #5 filename fi if [ ! "$met8" = "undefined" ] then let ifile=$ifile+1 echo "$DATA/ ">>CONTROL.conc # file #5 directory echo "${met8} ">>CONTROL.conc # file #5 filename fi if [ ! "$met9" = "undefined" ] then let ifile=$ifile+1 echo "$DATA/ ">>CONTROL.conc # file #5 directory echo "${met9} ">>CONTROL.conc # file #5 filename fi if [ $ifile -ne $nfile ] then msg='ERROR - Unexpected number of meteorology files written to CONTROL.conc' postmsg "$jlogfile" "$msg" export err=1; err_exit fi echo "$npol ">>CONTROL.conc # number of pollutants i=1 while [ $i -le $npol ]; do echo "${pol[$i]} ">>CONTROL.conc # pollutant ID echo "${qrate[$i]} ">>CONTROL.conc # emission rate (per hour) echo "$qhrs ">>CONTROL.conc # hours of emission year=`echo $dsyr | cut -c3-` echo "$year $dsmo $dsda $dshr $dsmn">>CONTROL.conc # release start let i=i+1 done echo "$ngrd ">>CONTROL.conc # number of simultaneous conc grids i=1 while [ $i -le $ngrd ]; do echo "${gctr[$i]} ">>CONTROL.conc # grid i center lat lon echo "${gspc[$i]} ">>CONTROL.conc # grid i grid spacing lat lon echo "${gspn[$i]} ">>CONTROL.conc # grid i grid span lat lon echo "${gdir[$i]} ">>CONTROL.conc # grid i directory echo "${gfil[$i]} ">>CONTROL.conc # grid i filename echo "${gnlv[$i]} ">>CONTROL.conc # grid i number of levels echo "${glvl[$i]} ">>CONTROL.conc # grid i levels echo "${gsst[$i]} ">>CONTROL.conc # grid i sampling start echo "${gsnd[$i]} ">>CONTROL.conc # grid i sampling end echo "${gint[$i]} ">>CONTROL.conc # grid i sampling interval let i=i+1 done echo "$npol ">>CONTROL.conc # number of pollutants depositing i=1 while [ $i -le $npol ]; do echo "${pspc[$i]} ">>CONTROL.conc # particle specs echo "${vdep[$i]} ">>CONTROL.conc # deposition velocity etc echo "${wet[$i]} ">>CONTROL.conc # wet removal echo "${half[$i]} ">>CONTROL.conc # half-life echo "${res[$i]} ">>CONTROL.conc # resuspension let i=i+1 done ############################# # now that control file is written, set gfile[2] to cdump for volcplot ############################# case $typenum in 2) gfil[2]=${gfil[1]};; esac ############################################################################### # Set up namelist parameters for each type of Dispersion/Conc Simulation # only those different from defaults must be given ############################################################################### rm -f SETUP.conc echo "&SETUP" > SETUP.conc echo "INITD = $initd," >> SETUP.conc # initial particle distribution echo "KHMAX = $khmax," >> SETUP.conc # max particle age echo "NUMPAR = $numpar," >> SETUP.conc # number of particles echo "MAXPAR = 50000," >> SETUP.conc # max number of particles echo "ISOT = $isot," >> SETUP.conc # isotropic turbulence echo "NDUMP = $ndump," >> SETUP.conc # write PARDUMP echo "KMSL = $kmsl" >> SETUP.conc # agl/msl if [ ${nhrs} -ge ${longrun} -o ${nhrs} -le ${neglongrun} ]; then echo "DELT = 20" >> SETUP.conc # time step (minutes) fi echo "/" >> SETUP.conc # exit # if just want to see CONTROL.conc and SETUP.conc # but hysplit_post.sh still runs ############################################################ # Run the Concentration Simulation ############################################################ rm -f ${gfil[1]} rm -f ${gfil[2]} pgm=hysplit_hymodelc export pgm;. prep_step startmsg $EXEChysplit/${pgm} conc >> $pgmout 2>errfile_a err=$?; export err;err_chk ############################################################ # verify hymodelc output file exists and filesize greater than just headers ############################################################ if [ ! -f ${gfil[1]} ] then msg='WARNING hysplit_hymodelc output file not created' errexit "$msg" fi echo "`wc -c ${gfil[1]}`" | read fsize tmp if [ $fsize -lt 1000 ] then msg='WARNING concentration graphic not available - problem with hysplit_hymodelc run' errexit "$msg" fi ################################################################## # read timestep from MESSAGE file (1 to 60 minutes) ################################################################## # this is written to MESSAGE in hymodelc.F as: # WRITE(30,*)' NOTICE main: Initial time step (min) ',INT(DT) ################################################################## STEP=`grep "Initial time step (min)" MESSAGE.conc | tail -c4` # convert emission duration in hours to minutes qhrsm=`echo "$qhrs*60" | bc` ################################################################## # if release duration chosen in setup_hysplit is less than the # initial timestep, then duration was actually the timestep instead # need to force the duration for the output plot labels # by re-making MAPTEXT.CFG ################################################################## if [ $STEP -gt $qhrsm ] then # convert from minutes to hours STEP=`echo "scale=4\n $STEP/60" | bc` qhrs=$STEP # qhrs input to volcplot for labeling w=`echo $qhrs | cut -d. -f1` # whole number portion dot=`expr index $qhrs .` # decimal portion if [ dot -gt 0 ]; then # truncate any unnecessary zeros .2500 ==> .25 frac=`echo $qhrs | cut -d. -f2` flen=`expr length $frac` d1=0 d2=0 d3=0 d4=0 if [ $flen -ge 1 ]; then d1=`echo $frac | cut -c1,1` if [ $flen -ge 2 ]; then d2=`echo $frac | cut -c2,2` if [ $flen -ge 3 ]; then d3=`echo $frac | cut -c3,3` if [ $flen -ge 4 ]; then d4=`echo $frac | cut -c4,4` fi fi fi fi if [ $d4 -eq 0 ]; then qhrs=${w}.${d1}${d2}${d3} if [ $d3 -eq 0 ]; then qhrs=${w}.${d1}${d2} if [ $d2 -eq 0 ]; then qhrs=${w}.${d1} if [ $d1 -eq 0 ]; then qhrs=${w}. fi fi fi fi fi msg='NOTICE duration is set to initial time step' postmsg "$jlogfile" "$msg" fi ############################################################################# # Create File to Label Graphics Maps - RSMC conc, Homeland conc # note that concplot uses lines 3,5,8-15 as long as they don't start "Traj" # lines not used for concplot are used in the RSMC coversheet ############################################################################# rm -f hyout.ps case $typenum in 1|3|4) rm -f MAPTEXT.conc echo $title | read tmp tmp1 case $typenum in 1) case $tmp in Test) label="EXERCISE EXERCISE EXERCISE";; Requested) label="REQUESTED SERVICES";; IAEA) label="IAEA NOTIFIED EMERGENCY";; Special) label=$title;; esac ;; 3) case $tmp in Test) label="TEST TEST TEST";; Real) label="As Requested ";; esac ;; 4) case $tmp in Test) label="EXERCISE EXERCISE EXERCISE";; Requested) label="REQUESTED SERVICES";; esac ;; esac # total=`echo "${qrat}*${qhrs}" | bc` # export total=`printf %3.2f $total` rm -f MAPTEXT.conc case $typenum in 1|4) echo "RSMC Washington - NOAA ARL / NCEP" >MAPTEXT.conc;; 3) echo "NOAA NWS NCEP NCO " >MAPTEXT.conc;; esac echo "-----------------------------------------------------------" >>MAPTEXT.conc case $typenum in 1|4) echo "Created: $ihr${imn}UTC $ida/$imo/$iyr (day/month/year) RSMC Washington - NOAA ARL / NCEP" >>MAPTEXT.conc;; 3) echo "Created:$ihr${imn}UTC $imo/$ida/$iyr NOAA NWS NCEP NCO" >>MAPTEXT.conc;; esac echo "-----------------------------------------------------------" >>MAPTEXT.conc echo "Source:${site} lat:${olat} lon:${olon} hgt:${olvl1} to ${olvl2} m" >>MAPTEXT.conc echo "Location:${site} lat:${olat} lon:${olon} " >>MAPTEXT.conc echo "Release Start (YYYY MM DD HH MM):${dsyr} ${dsmo} ${dsda} ${dshr} ${dsmn} " >>MAPTEXT.conc case $typenum in 1) echo "Release ID:${ident} Rate: ${qrat} Bq/hr Duration: ${qhrs} hr Particles: ${numpar}" >>MAPTEXT.conc;; *) echo "Release ID:${ident} Rate: ${qrat}/hr Duration: ${qhrs} hr Particles: ${numpar}" >>MAPTEXT.conc;; esac echo "Traj created $iday $imonth $ida ${ihr}:${imn} UTC ${iyr} " >> MAPTEXT.conc echo "Distribution: Uniform between ${olvl1} and ${olvl2} m AGL" >>MAPTEXT.conc case $typenum in 1|3|4) echo "Dry Deposition Rate:${dryvl} m/s Wet Removal (below/in-cloud):${wetlo}/${wetin}" >>MAPTEXT.conc;; esac echo "Meteorology: $DATlabel " >>MAPTEXT.conc echo "Note: Contour values may change from chart to chart" >>MAPTEXT.conc case $typenum in 1) if [ ${icmt} -eq 1 ]; then echo "Note: RESULTS BASED ON DEFAULT INITIAL VALUES" >>MAPTEXT.conc else echo " " >> MAPTEXT.conc fi ;; 3|4) echo " " >>MAPTEXT.conc;; esac echo "Response: $label " >>MAPTEXT.conc echo "-----------------------------------------------------------" >>MAPTEXT.conc ;; esac ############################################################ # Output graphics ################################################################## # Setup plot program filename, options and post-processing specs # GEMPAK in hysplit_post.sh ################################################################## case $typenum in 1) pgm=hysplit_concplot # output filename concplot.ps PGM=hysplit_concplot map[1]=arlmap map[2]=arlmap gis=0 # GIS? 0=no 1=yes suffix[1]=conc suffix[2]=fine options[1]="-e1 -r3 -uBq -z${zoom1}" # see below #options[2]="-e1 -r3 -uBq -z${zoom2} -n${N1}:${N2}" # current time --> +18 h options[2]="-e1 -r3 -uBq -z${zoom2} " # rlse time --> current time +18 h cp MAPTEXT.conc MAPTEXT.${suffix[2]} ;; 2) pgm=hysplit_volcplot PGM=hysplit_volcplot for imap in 1 2 3 4 5 6 7 8 9; do map[$imap]=arlmap gfil[$imap]=${gfil[1]} # force this to plot all reductions done gis=0 # 6-h plots are for internal use # 12-h plots are for dissemination (wafs) N=0 for tint in 06 12; do numb= kolor=1 # color zoom=${zoom1} # was 80 but now keep it same as +12-h plots if [ $tint -eq 12 ] then numb=-n-2 kolor=0 # B&W (want better resolution in g3-format graphic) zoom=${zoom1} # was 90 fi for red in 0 1 2 3; do let N=$N+1 options[$N]="-1${alert} -2${atyp} -3${red} -4${icmt} -5${site} -6${qhrs} -z${zoom} -k${kolor} -ovolcplot${tint}.r${red} $numb" done done # alert-ended run for dissemination (force alert=1 atyp=2 reduction=0) tint=ae alert=0 atyp=2 red=0 options[9]="-1${alert} -2${atyp} -3${red} -4${icmt} -5${site} -6${qhrs} -z80 -ovolcplot${tint}.r${red} -n-2" ;; 3) pgm=hysplit_concplot PGM=hysplit_concplot map[1]=countymap # fine grid map[2]=arlmap # coarse grid gis=1 suffix[1]=fine suffix[2]=coarse options[1]=" -d2 -s1 -z${zoom1} -c51" if [ ${relcen} -eq 1 ]; then options[1]="${options[1]} -g0" fi options[2]=" -d2 -s1 -z${zoom2} -c1" cp MAPTEXT.conc MAPTEXT.${suffix[1]} cp MAPTEXT.conc MAPTEXT.${suffix[2]} ;; 4) pgm=hysplit_concplot # output filename concplot.ps PGM=hysplit_concplot map[1]=arlmap gis=0 # GIS? 0=no 1=yes suffix[1]=conc options[1]="-u -z${zoom1}" # see below, -u turns off "mass" label #cp MAPTEXT.conc MAPTEXT.${suffix[2]} ;; esac case $typenum in 1) lgrd=$ngrd;; # only set for ngrd=1 2) lgrd=9;; # to do all cases 3) lgrd=$ngrd;; # only set for ngrd=2 4) lgrd=$ngrd;; # only set for ngrd=1 esac # test for GDAS-GFS run to modify volcplot met init time labeling rm -f GDAS-GFS.date # wc -l means count number of lines with a match if [ `grep GDAS MESSAGE.conc | wc -l` -ne 0 ] then if [ `grep GFS MESSAGE.conc | wc -l` -ne 0 ] then FTPDONE=GFS.ftpdone if [ -s ${FTPDONE} ] then YR=`cat ${FTPDONE} | cut -c1-2` if [ ${YR} -lt 10 ] then YR=`cat ${FTPDONE} | cut -c2-2` YR="0"${YR} fi MN=`cat ${FTPDONE} | cut -c3-4` if [ ${MN} -lt 10 ] then MN=`cat ${FTPDONE} | cut -c4-4` MN="0"${MN} fi DY=`cat ${FTPDONE} | cut -c5-6` if [ ${DA} -lt 10 ] then DY=`cat ${FTPDONE} | cut -c6-6` DY="0"${DY} fi HR=`cat ${FTPDONE} | cut -c7-8` if [ ${HR} -lt 10 ] then HR=`cat ${FTPDONE} | cut -c8-8` HR="0"${HR} fi echo "GFS ${YR} ${MN} ${DY} ${HR}" > GDAS-GFS.date else echo "WARNING: GFS.ftpdone file not found, plots will not contain GFS labels" fi fi fi grd=1 while [ $grd -le $lgrd ] do # copy appropriate map boundary file cp $FIXhysplit/hysplit_${map[$grd]} arlmap # do each grid or each volcanic ash reduction level (not both) # -p... is ignored by volcplot # concplot options: # -a = write GIS text file (one file per time period) # -c = contours: (0)-dynamic, 1-fixed, Val-10^(val) # -d = concentration average between levels # -e = exposure # -f = one frame per output file # -o = output filename # -p = output file identifier # -n = number of time periods: (0)-all, numb, min:max, -incr # -r = removal (3=total) # -s = number of species # -u = units # -z = zoom (0-100) # dummy in DHOUR for RSMC and hls (input file to hysplit_concrop for volc ash) case $typenum in 1|3|4) echo "20" > DHOUR ;; esac TMP=${gfil[$grd]} export pgm=hysplit_concrop;. prep_step startmsg case $typenum in 1|3|4) $EXEChysplit/hysplit_concrop -i${gfil[$grd]} -occrop -f0;; 2) $EXEChysplit/hysplit_concrop -i${gfil[$grd]} -occrop -f1;; esac err=$?; export err; err_chk gfil[$grd]=ccrop export pgm=${PGM};. prep_step startmsg $EXEChysplit/${pgm} -i${gfil[$grd]} -p${suffix[$grd]} ${options[$grd]} >> $pgmout 2>errfile_${grd} err=$?; export err; err_chk gfil[$grd]=$TMP let grd=$grd+1 done if [ -s concplot.conc ] then cp concplot.conc concplot.ps fi # add this copy concplot.fine 4-30-07 if [ -s concplot.fine ] then cp concplot.fine concplot.fine.ps fi case $typenum in 3) if [ -s concplot.fine ] then cp concplot.fine concplot_loop.fine fi if [ -s concplot.coarse ] then cp concplot.coarse concplot_loop.coarse fi esac ####################### # Google Earth output ####################### case $typenum in 1|3|4) grd=1 while [ $grd -le $lgrd ] do case $typenum in 1) options[1]="-e1 -r3 -uBq -z10 -a3" options[2]="-e1 -r3 -uBq -z10 -a3";; 3) options[1]="-a3 -d2 -s1 -z10 -c50" options[2]="-a3 -d2 -s1 -z10";; 4) options[1]="-u -z10 -a3" options[2]="-u -z10 -a3";; esac # copy appropriate map boundary file # - don't need map for GIS but get default arlmap so no WARNING msg cp $FIXhysplit/hysplit_arlmap arlmap export pgm=${PGM};. prep_step startmsg $EXEChysplit/${pgm} -i${gfil[$grd]} -p${suffix[$grd]} ${options[$grd]} >> $pgmout 2>errfile_${grd}_a3 err=$?; export err; err_chk # don't need postscript file rm concplot.${suffix[$grd]} if [ -f GELABEL_${suffix[$grd]}.txt ];then export pgm=hysplit_gelabel;. prep_step startmsg $EXEChysplit/${pgm} -p${suffix[$grd]} >> $pgmout 2>errfile_gelabel err=$?; export err; err_chk if [ -f GELABEL_01_${suffix[$grd]}.ps ];then let frame=1 fnum=`printf %2.2d $frame` while [ -f GELABEL_${fnum}_${suffix[$grd]}.ps ]; do ${GSPATH}/gs -dTextAlphaBits=4 -dGraphicsAlphaBits=4 -dNOPAUSE -dSAFER -sDEVICE=pnmraw -q -sOutputFile=- GELABEL_${fnum}_${suffix[$grd]}.ps -c quit | convert +adjoin - -crop 0x0 GIF:GELABEL_${fnum}_${suffix[$grd]}.gif # rm GELABEL_${fnum}_${suffix[$grd]}.ps let frame=$frame+1 fnum=`printf %2.2d $frame` done cp ${FIXhysplit}/hysplit_logocon.gif logocon.gif /usr/bin/zip HYSPLIT_${suffix[$grd]}.kmz HYSPLIT_${suffix[$grd]}.kml logocon.gif GELABEL*${suffix[$grd]}.gif #rm -f HYSPLIT_${suffix[$grd]}.kml #rm -f GELABEL*${suffix[$grd]}.gif fi fi let grd=$grd+1 done ;; esac ########################################################### # For hls, loops made above, ie plots with same concentration contours # now make individual plots with contours set separately for each plot ########################################################### case $typenum in 3) # run barebones concplot to count number of maps to make individual frames grd=1 while [ $grd -le $lgrd ] do optionstmp[$grd]="-d2 -f1 -s1" export pgm=${PGM};. prep_step startmsg $EXEChysplit/${pgm} -i${gfil[$grd]} -p${suffix[$grd]} ${optionstmp[$grd]} >> $pgmout 2>errfile_${grd} err=$?; export err; err_chk NMAPS=0 ls concplot00??.${suffix[$grd]} | while read tmp do let NMAPS=${NMAPS}+1 done rm concplot00??.${suffix[$grd]} # set concplot options options[1]="-d2 -s1 -z${zoom1} -c50" options[2]="-d2 -s1 -z${zoom2}" if [ ${relcen} -eq 1 ]; then options[1]="${options[1]} -g0" fi # copy appropriate map boundary file cp $FIXhysplit/hysplit_${map[$grd]} arlmap # loop through maps, plotting each separately N=1 while [ $N -le $NMAPS ] do ####################### # GIS output ####################### export pgm=${PGM};. prep_step startmsg $EXEChysplit/${pgm} -i${gfil[$grd]} -p${suffix[$grd]} -n${N}:${N} -a${gis} -d2 -s1 -z10 >> $pgmout 2>errfile_${grd} err=$?; export err; err_chk # remove postscript file rm -f concplot.${suffix[$grd]} ####################### # postscript output ####################### # dummy in DHOUR for RSMC and hls (input file to hysplit_concrop for volc ash) case $typenum in 1|3|4) echo "20" > DHOUR ;; esac TMP=${gfil[$grd]} export pgm=hysplit_concrop;. prep_step startmsg $EXEChysplit/hysplit_concrop -i${gfil[$grd]} -occrop -f0 err=$?; export err; err_chk gfil[$grd]=ccrop export pgm=${PGM};. prep_step startmsg $EXEChysplit/${pgm} -i${gfil[$grd]} -p${suffix[$grd]} -n${N}:${N} ${options[$grd]} >> $pgmout 2>errfile_${grd} err=$?; export err; err_chk gfil[$grd]=$TMP if [ $N -lt 10 ]; then n="000${N}" elif [ $N -lt 100 ]; then n="00${N}" elif [ $N -lt 1000 ]; then n="0${N}" fi mv concplot.${suffix[$grd]} concplot${n}.${suffix[$grd]} let N=$N+1 done let grd=$grd+1 done ;; esac ########################################################### # output to /com ########################################################### if [ "$SENDCOM" = "YES" ] then cp CONTROL.conc $COMOUT/hysplit.${cycle}.${site}.${TYP}${seq}.CONTROL.conc cp SETUP.conc $COMOUT/hysplit.${cycle}.${site}.${TYP}${seq}.SETUP.conc cp MESSAGE.conc $COMOUT/hysplit.${cycle}.${site}.${TYP}${seq}.MESSAGE.conc cp hysplit.ini $COMOUT/hysplit.${cycle}.${site}.${TYP}${seq}.hysplit.ini case $typenum in 2) ls volcplot* | while read tmp do cp $tmp $COMOUT/hysplit.${cycle}.${site}.${TYP}${seq}.$tmp done ;; esac if [ -s DHOUR ] then cp DHOUR $COMOUT/hysplit.${cycle}.${site}.${TYP}${seq}.DHOUR fi if [ -s cdump ] then cp cdump $COMOUT/hysplit.${cycle}.${site}.${TYP}${seq}.cdump fi if [ -s cdump2 ] then cp cdump2 $COMOUT/hysplit.${cycle}.${site}.${TYP}${seq}.cdump2 fi if [ -s MAPTEXT.conc ] then cp MAPTEXT.conc $COMOUT/hysplit.${cycle}.${site}.${TYP}${seq}.MAPTEXT.conc fi if [ -s concplot.ps ] then cp concplot.ps $COMOUT/hysplit.${cycle}.${site}.${TYP}${seq}.concplot.ps fi case $typenum in 1|3|4) ls concplot* | while read tmp do cp $tmp $COMOUT/hysplit.${cycle}.${site}.${TYP}${seq}.$tmp done ;; esac case $typenum in 3) ls GIS* | while read tmp do cp $tmp $COMOUT/hysplit.${cycle}.${site}.${TYP}${seq}.$tmp done ;; esac case $typenum in 1|3|4) ls *.kmz | while read tmp do tmpfile=`echo $tmp | cut -f2 -d_` cp $tmp $COMOUT/hysplit.${cycle}.${site}.${TYP}${seq}.$tmpfile done ;; esac fi exit