#!/bin/ksh -fx
#####################################################################
#  UNIX Script Documentation Block
#                      .
# Script name:         exaqm_premaq.sh
# Script description:  Run premaq emissions processing for CMAQ
#
# Author:      Marina Tsidulko and Pius Lee Org: NP22 Date: 2004-03-31
#
# Abstract: This script runs premaq emissions processing
#
# Script history log:
# 2003-07-03    Marina Tsidulko, initial version
# 2003-07-15    Julia Zhu, modified for production  
# 2004-03-31    Pius Lee, netcdf independence and gfs ozone blending
#####################################################################
set -xa
msg="JOB $job HAS BEGUN"
postmsg "$jlogfile" "$msg"

cd $DATA

## Specify the input pollutant data directory:
if [ ${cycle} = 't00z' -o ${cycle} = 't06z' ]
then
   yyyy=`echo $PDYm1 |cut -c1-4`
   year_month=`echo $PDYm1 |cut -c1-6`
elif [ ${cycle} = 't12z' -o ${cycle} = 't18z' ]
then
   yyyy=`echo $PDY |cut -c1-4`
   year_month=`echo $PDY |cut -c1-6`
fi

export LOGFILE=chkreads.log

rm -f $LOGFILE

#------------------ names for tests ----------------
INPPLTDRY=$INPEMIDIR/$yyyy
INPPLTDR=$INPEMIDIR/$year_month
#---------------------------------------------------
job="PREMAQ CONUS"
INGRID=138
NLAYS=22      # '22' (collapse from 60), '60' (full nmm level)
OUTGRID=CONUS # 'NEUS' (Northeast), 'EASTUS' (3x), or 'CONUS' (5x)
EXPT=invert
DDATE=01May
PROG=premaq
export MOBILE_PREMAQ_VERSION=2005
export BINIO=T      # binary I/O?
MKGRID=T     # make GRID files?
MKEMIS=T     # make emissions file?
#SE MKBC=T       # make time-varying BCs?
MKBC=F       # make time-varying BCs?
PX_RA=F      # If F => use NAM calculated aerodynamic resistance
             #    T => use Pleim calculated Ra using NAM parameters
MBL_T=F      # If T => use Temperature dependent scheme for Mobile Em
             #    F => No temperature dependence 
             #      if F also must comment out line 494 "echo 'MGTS_L'
GAS_CHEM=cb4 # select from   ( cb4, cb5)
INV_BS=2005aa # select from   (2005aa)

typeset -u OUTGRID=$OUTGRID

export fhr=00
export soilhr=06   ##06 for operational run
let "soilhr=soilhr+1"

case $cyc in 
 00) export endfhr=06
     export endday=$PDY
     export soil_input=$COMINm1/aqm.t18z.soil06.bin
     export nmm_00=$COMINm1/aqm.t18z.nmm06.tm00;;
 06) export endfhr=48
     export endday=$PDYp2
     export soil_input=$COMIN/aqm.t00z.soil06.bin
     export nmm_00=$COMIN/aqm.t00z.nmm06.tm00;;
 12) export endfhr=48
     export endday=$PDYp2
     export soil_input=$COMIN/aqm.t06z.soil06.bin
     export nmm_00=$COMIN/aqm.t06z.nmm06.tm00;;
 18) export endfhr=06
     export endday=$PDYp1
     export soil_input=$COMIN/aqm.t12z.soil06.bin
     export nmm_00=$COMIN/aqm.t12z.nmm06.tm00;;
esac
 
let "endcyc=cyc+endfhr"
while [ $endcyc -ge 24 ]
do
  let "endcyc=endcyc-24"
done
typeset -Z2 endcyc
echo $endcyc 

if [ -s $soil_input ]
then
   echo "=====>  The Previous Cycle Ran Successfully"
   export START=WARM
   export INITIAL_RUN=N
else
   export START=COLD
   export INITIAL_RUN=Y
fi

if [ $INITIAL_RUN = "N" ]
then
   echo 'INITIAL_RUN = ' $INITIAL_RUN
   ln -s $soil_input soilinp.$PDY.${OUTGRID}.bin
else
   echo 'INITIAL_RUN = ' $INITIAL_RUN
fi

  while [ $fhr -le $endfhr ]
   do
    ln -s $COMNAM/aqm.${cycle}.nmm${fhr}.tm00  $DATA/.
    let "fhr=fhr+1"
    typeset -Z2 fhr
   done

   rm -f $DATA/aqm.${cycle}.nmm00.tm00
   cp $COMNAM/aqm.${cycle}.nmm00.tm00 $DATA/nmm00
   ${utilexec}/wgrib.aqm $DATA/nmm00  | awk '{if ($1 > 0 && $1 < 848) print $0}' FS=':' | ${utilexec}/wgrib.aqm -i -grib -append $DATA/nmm00 -o $DATA/aqm.${cycle}.nmm00.tm00
   ${utilexec}/wgrib.aqm -d 848 -grib -append $nmm_00 -o $DATA/aqm.${cycle}.nmm00.tm00
   ${utilexec}/wgrib.aqm $DATA/nmm00  | awk '{if ($1 > 848 && $1 < 909) print $0}' FS=':' | ${utilexec}/wgrib.aqm -i -grib -append $DATA/nmm00 -o $DATA/aqm.${cycle}.nmm00.tm00

SYYYY=`echo $PDY |cut -c1-4`
SMM=`echo $PDY |cut -c5-6`
SDD=`echo $PDY |cut -c7-8`
EYYYY=`echo $endday |cut -c1-4`
EMM=`echo $endday |cut -c5-6`
EDD=`echo $endday |cut -c7-8`

if [[ $OUTGRID = 'NEUS' ]]
then
  CoordName=AQF_LamCon_12      # 16-character maximum
  GridName=AQF_Northeast_US    # 16-character maximum
  BTRIM=0
  X0=1
  Y0=1
  NCOLS=20
  NROWS=20
elif [[ $OUTGRID = 'EASTUS' ]]
then
  CoordName=AQF_RPO            # 16-character maximum
  GridName=AQF_East_US         # 16-character maximum
  BTRIM=-1
  X0=199
  Y0=1
  NCOLS=268
  NROWS=259
  LPRT_COL=0
  LPRT_ROW=0  
elif [[ $OUTGRID = 'CONUS' ]]
then
  CoordName=AQF_RPO            # 16-character maximum
  GridName=AQF_CONUS           # 16-character maximum
  BTRIM=-1
  X0=25
  Y0=1
  NCOLS=442
  NROWS=265
  LPRT_COL=0
  LPRT_ROW=0
  VGTOP_GD=200.0             # P of top (=mb*100)
else
  print Unexpected output grid specified
  print OUTGRID is $OUTGRID
  exit 111
fi

if [ $NLAYS = 22 ]; then
   CTMLAYS="1.000000, 0.991971, 0.975778, 0.959038, 0.941624,  \
            0.923362, 0.903749, 0.882275, 0.858604, 0.830028,  \
            0.791646, 0.744810, 0.663481, 0.572346, 0.472587,  \
            0.363190, 0.245169, 0.123430,  \
            1.000000, 0.820492, 0.650769, 0.464416, 0.242600"
   NETAC=18
elif [ $NLAYS = 29 ]; then
   CTMLAYS="1.000000, 0.991971, 0.983896, 0.975778, 0.967493,  \
            0.959038, 0.941624, 0.923362, 0.903749, 0.882275,  \
            0.858604, 0.830028, 0.791646, 0.744810, 0.691719,  \
            0.634061, 0.572346, 0.506741, 0.437247, 0.363190,  \
            0.245169, 0.123430,  \
            1.000000, 0.879922, 0.762043, 0.650769, 0.552946,  \
            0.464416, 0.362189, 0.242600"
   NETAC=22
elif [ $NLAYS = 56 ]; then
   CTMLAYS="1.000000, 0.991971, 0.983896, 0.975778, 0.967493,  \
            0.959038, 0.950415, 0.941624, 0.932661, 0.923362,  \
            0.913724, 0.903749, 0.893265, 0.882275, 0.870777,  \
            0.858604, 0.845076, 0.830028, 0.811937, 0.791646,  \
            0.768820, 0.744810, 0.719109, 0.691719, 0.663481,  \
            0.634061, 0.603796, 0.572346, 0.540051, 0.506741,  \
            0.472587, 0.437247, 0.400727, 0.363190, 0.324469,  \
            0.285074, 0.245169, 0.204758, 0.164179, 0.123430,  \
            0.082512, 0.041425,  \
            1.000000, 0.939839, 0.879922, 0.820492, 0.762043,  \
            0.704572, 0.650769, 0.599902, 0.552946, 0.508926,  \
            0.464416, 0.417216, 0.362189, 0.302762, 0.242600"
   NETAC=42
elif [ $NLAYS = 60 ]; then
   CTMLAYS="1.000000, 0.991971, 0.983896, 0.975778, 0.967493,  \
            0.959038, 0.950415, 0.941624, 0.932661, 0.923362,  \
            0.913724, 0.903749, 0.893265, 0.882275, 0.870777,  \
            0.858604, 0.845076, 0.830028, 0.811937, 0.791646,  \
            0.768820, 0.744810, 0.719109, 0.691719, 0.663481,  \
            0.634061, 0.603796, 0.572346, 0.540051, 0.506741,  \
            0.472587, 0.437247, 0.400727, 0.363190, 0.324469,  \
            0.285074, 0.245169, 0.204758, 0.164179, 0.123430,  \
            0.082512, 0.041425,  \
            1.000000, 0.939839, 0.879922, 0.820492, 0.762043,  \
            0.704572, 0.650769, 0.599902, 0.552946, 0.508926,  \
            0.464416, 0.417216, 0.362189, 0.302762, 0.242600,  \
            0.181950, 0.121299, 0.060648, 0.000000"
   NETAC=42
else
   print Unexpected ctmlay layers specified
   print NLAYS is $NLAYS
   exit 111
fi

export fhr=00
while [ $fhr -le $endfhr ]
do
fcsthrs[$fhr]=$fhr
InMetFiles[$fhr]=$DATA/meso.AQFNAM$fhr
let "fhr=fhr+1"
typeset -Z2 fhr
done

for fhr in ${fcsthrs[*]}
do
   ln -s -f $DATA/aqm.${cycle}.nmm$fhr.tm00 meso.AQFNAM$fhr
done

AQF_START=${SYYYY}-${SMM}-${SDD}-${cyc}:00:00.0000
AQF_END=${EYYYY}-${EMM}-${EDD}-${endcyc}:00:00.0000

INTVL=60
METVERT=1

#-------------------------------------------------------------------
#  varilables of the namelist "userdefs"
#
#  As in WRF-NMM (60 MDL, 61 INTERFACE)
#     "METLAYH", the first  41 value is "AEta2"
#                the second 19 value is "AEta1"
#
#     "CTMLAYS", the first  41 value is "ETA2"
#                the second 20 value is "Eta1"
#
#  *** CAUTION: "METLAYH" & "CTMLAYS" are BOTTOM UP (surface 1st)
#-------------------------------------------------------------------

METLAYH="0.995986, 0.987934, 0.979837, 0.971636, 0.963266,  \
         0.954727, 0.946020, 0.937142, 0.928012, 0.918543,  \
         0.908736, 0.898507, 0.887770, 0.876526, 0.864690,  \
         0.851840, 0.837552, 0.820983, 0.801791, 0.780233,  \
         0.756815, 0.731959, 0.705414, 0.677600, 0.648771,  \
         0.618929, 0.588071, 0.556199, 0.523396, 0.489664,  \
         0.454917, 0.418987, 0.381959, 0.343830, 0.304772,  \
         0.265121, 0.224964, 0.184469, 0.143805, 0.102971,  \
         0.061969, 0.020713,  \
         0.969919, 0.909880, 0.850207, 0.791267, 0.733308,  \
         0.677671, 0.625335, 0.576424, 0.530936, 0.486671,  \
         0.440816, 0.389702, 0.332475, 0.272681, 0.212275,  \
         0.151624, 0.090974, 0.030324"

NETA2=42          # number of lower sigma Eta2 level
NETAC=18          # number of lower sigma Eta2 level in CTM

IWIND=0           # way of vertical wind correction (top-down or bottom-up)

BTRIM=-1
METVERT=1         # vertical grid info, 1:sigma, 2:nmm sigma-p
                  # ONLY "1" AT THIS MOMENT (it covers both "1" & "2")

set -A InMM ${InMetFiles[*]}
FILE_GD=$DATA/GRIDDESC
#--------------------------------------------------------------------
# Create namelist with user definitions.
#--------------------------------------------------------------------
if [ -f $DATA/files_list ]
then
  rm $DATA/files_list
fi

export fhr=00
while [ $fhr -le $endfhr ]
do
echo "\"${InMetFiles[$fhr]}\"," >> $DATA/files_list
let "fhr=fhr+1"
typeset -Z2 fhr
done

Files_List=`cat $DATA/files_list`

MACHTYPE=`uname`
if [[ $MACHTYPE = "AIX" ]]
then
  Marker="/"
else
  Marker="&END"
fi
cat > $DATA/namelist.${PROG} << EOF05

 &FILENAMES
  file_gd    = "$FILE_GD"
  files_in  = $Files_List
  makegrid   = .${MKGRID}.
  makeemis   = .${MKEMIS}.
  makebc     = .${MKBC}.
 $Marker

 &USERDEFS
  aqf_start  = "$AQF_START"
  aqf_end    = "$AQF_END"
  intvl      =  $INTVL
  iwind      =  $IWIND
  vgtop_gd   =  $VGTOP_GD
  coordnam   = "$CoordName"
  grdnam     = "$GridName"
  metlayh    =  $METLAYH
  ctmlays    =  $CTMLAYS
  btrim      =  $BTRIM
  metvert    =  $METVERT
  neta2      =  $NETA2
  netac      =  $NETAC
  lprt_col   =  $LPRT_COL
  lprt_row   =  $LPRT_ROW
 $Marker

 &WINDOWDEFS
  x0         =  $X0
  y0         =  $Y0
  ncolsin    =  $NCOLS
  nrowsin    =  $NROWS
 $Marker

EOF05

rm -f fort.*
ln -s -f $FILE_GD                   fort.4
ln -s -f ./namelist.${PROG}         fort.8

integer NUMFIL=0
for fil in $InMetFiles
do
  ((NN = NUMFIL + 10))
  ln -s -f $fil fort.$NN
  ((NUMFIL = NUMFIL + 1))
done

#---------------------------------------------------------------------
# Set output file names and other miscellaneous environment variables.
#---------------------------------------------------------------------

export IOAPI_CHECK_HEADERS=F
export EXECUTION_ID=PREMAQ
export NOTCDF_FATAL=F
export BNDY_GASC_1=BIN:$FIXaqm/aqm_bc_constanteast_${OUTGRID}.bin
if [[ $MKGRID = "T" ]] then
  export GRID_BDY_2D_G1=BIN:$DATA/GRIDBDY2D
  export GRID_BDY_3D_G1=BIN:$DATA/GRIDBDY3D
  export GRID_CRO_2D_G1=BIN:$DATA/GRIDCRO2D
  export GRID_CRO_3D_G1=BIN:$DATA/GRIDCRO3D
  export GRID_DOT_2D_G1=BIN:$DATA/GRIDDOT2D
fi
export MET_BDY_2D_G1=BIN:$DATA/METBDY2D
export MET_BDY_3D_G1=BIN:$DATA/METBDY3D
export MET_CRO_2D_G1=BIN:$DATA/METCRO2D
export MET_CRO_3D_G1=BIN:$DATA/METCRO3D
export MET_DOT_3D_G1=BIN:$DATA/METDOT3D
export CHEM_BDY_3D_G1=BIN:$DATA/CHEMBDY3D
export CHEM_CRO_3D_G1=BIN:$DATA/CHEMCRO3D

#---------------------------------------------------------------------
# Set input file names.
#---------------------------------------------------------------------
export LANDA_G1=BIN:$FIXaqm/aqm_LANDA_${OUTGRID}.bin
if [[ $MKBC = "T" ]]
then
export BNDY_GASC_1=BIN:$FIXaqm/aqm_bc_constanteast_${OUTGRID}.bin
fi

#---------------------------------------------------------------------
# Set BEIS311 environment variables and file names.
#---------------------------------------------------------------------

export SOILOUT_TIME=${soilhr}
export PROMPTFLAG=F
export BIOSW_YN=N
export SUMMER_YN=Y
export B3GRD=BIN:$FIXaqm/aqm_static_b3grd_${OUTGRID}.bin
export GSPRO=$FIXaqm/aqm_static_gspro
export B3GTS_L=BIN:$DATA/aqm_b3gts_l.${PDY}.${OUTGRID}.bin
export SOILINP=BIN:$DATA/soilinp.${PDY}.${OUTGRID}.bin
export SOILOUT=BIN:$DATA/soilout.${PDY}.${OUTGRID}.bin

#---------------------------------------------------------------------
# Set laypoint environment variables and file names.
#
# Note: Program aborts before reaching here for non-defined layers
#---------------------------------------------------------------------

export SMK_SOURCE='P'

if [ $NLAYS = 22 ]; then
   export SMK_EMLAYS=15
   export REP_LAYER_MAX=15
elif [ $NLAYS = 29 ]; then
   export SMK_EMLAYS=19
   export REP_LAYER_MAX=19
elif [ $NLAYS = 56 ]; then
   export SMK_EMLAYS=34
   export REP_LAYER_MAX=34
elif [ $NLAYS = 60 ]; then
   export SMK_EMLAYS=34
   export REP_LAYER_MAX=34
fi

export SMK_SPECELEV_YN=N
export HOUR_PLUMEDATA_YN=N
export VELOC_RECALC=N
export IOAPI_ISPH=19
export INSPH_CODE=19
export OUTZONE=0
export UNIFORM_TPROF_YN=N
export DAY_SPECIFIC_YN=N
export HOUR_SPECIFIC_YN=N
export SMK_AVEDAY_YN=N
export RENORM_TRPOF=T
export REPORT_DEFAULTS=N
export ZONE4WM=T
export OUTPUT='.'


export PNTS=$INPPLTDRY/pnts.map.txt
export INVNAME1=pnts
export INVNAME2=psrc
export INVTABLE=$FIXaqm/aqm_static_invtable.criteria.nospec.txt

export PLAY=BIN:$DATA/play.${PDY}.1.aqf12_nei1aqf99.bin
rm -f $DATA/play.${PDY}.1.aqf12_nei1aqf99.bin

#---------------------------------------------------------------------
# Set smkmerge environment variables and file names.
#---------------------------------------------------------------------

export MRG_SOURCE='P'
export MRG_CTLMAT_MULT=' '
export MRG_CTLMAT_ADD=' '
export MRG_CTLMAT_REAC=' '
export MRG_TEMPORAL_YN=Y
export MRG_SPCMAT_YN=Y
export MRG_METCHK_YN=Y
export MRG_MARKETPEN_YN=N
export MRG_GRDOUT_YN=Y
export MRG_REPSTA_YN=N
export MRG_REPCNY_YN=N
export MRG_REPINV_YN=Y
export MRG_REPSPC_YN=Y
export MRG_REPCTL_YN=Y
export MRG_LAYERS_YN=Y
export SMK_PING_METHOD=0
export SMK_ASCIIELEV_YN=N
export EXPLICIT_PLUMES_YN=N
export SMK_O3SEASON_YN=Y
export SMK_DEFAULT_TZONE=5
export MRG_BYDAY=' '
export PTREF=$FIXaqm/aqm_static_ptref.m3.us+can.txt          
export PTPRO=$FIXaqm/aqm_static_ptpro.m3.default.us+can.txt  
export HOLIDAYS=$FIXaqm/aqm_static_holidays.txt


export PGMAT=BIN:$INPPLTDRY/aqm_static_pgmat_${OUTGRID}_${yyyy}.bin
export PSMAT=BIN:$INPPLTDRY/aqm_static_psmat_cb4_${yyyy}.bin
export COSTCY=$FIXaqm/costcy_for_2002_19dec2007_v1.txt

export ESCEN='nei1aqf99'

export MRG_GRDOUT_UNIT='mole/s'
export MRG_TOTOUT_UNIT='mole/s'
export PGTS3D_L=BIN:$DATA/pgts3d_l.${PDY}.1.aqf12.nei1aqf99.bin

rm -f $DATA/pgts3d_l.${PDY}.1.aqf12.nei1aqf99.bin
#---------------------------------------------------------------------
# Set mobile source processing variables.
#---------------------------------------------------------------------
export COEFS_EVAPEXH=BIN:$INPPLTDRY/aqm_static_coef_evapexh.cb4.${OUTGRID}.${yyyy}.bin
export MGTS_L=BIN:$DATA/aqm_mgts_l.${PDY}.${OUTGRID}.bin
export LOOKUP_PM=BIN:$INPPLTDRY/aqm_lookup_table.${OUTGRID}.${yyyy}.bin
if [ ${cycle} = 't00z' -o ${cycle} = 't06z' ]
then
   export AGTS_L=BIN:$INPPLTDR/aqm_pollutant_${GAS_CHEM}_${PDYm1}_agts_${INV_BS}.${OUTGRID}.bin
   export CGTS_L=BIN:$INPPLTDR/aqm_pollutant_${GAS_CHEM}_${PDYm1}_mgts_${INV_BS}.${OUTGRID}.bin
elif [ ${cycle} = 't12z' -o ${cycle} = 't18z' ]
then
   export AGTS_L=BIN:$INPPLTDR/aqm_pollutant_${GAS_CHEM}_${PDY}_agts_${INV_BS}.${OUTGRID}.bin
   export CGTS_L=BIN:$INPPLTDR/aqm_pollutant_${GAS_CHEM}_${PDY}_mgts_${INV_BS}.${OUTGRID}.bin
fi

#---------------------------------------------------------------------
# Set mrggid processing variables.
#---------------------------------------------------------------------

export FILELIST=$DATA/filelist.aqf
if [ -f $FILELIST ]
then
  rm -f $FILELIST
fi
if [ $MBL_T = "F" ]; then
  echo 'AGTS_L'   >> $FILELIST
  echo 'B3GTS_L'  >> $FILELIST
  echo 'PGTS3D_L' >> $FILELIST
  echo 'CGTS_L'  >> $FILELIST # No need MGTS_L since no Temp dependence
else
  echo 'AGTS_L'   >> $FILELIST
  echo 'B3GTS_L'  >> $FILELIST
  echo 'PGTS3D_L' >> $FILELIST
  echo 'MGTS_L'  >> $FILELIST   
  echo 'CGTS_L'  >> $FILELIST
fi

export OUTFILE=BIN:$DATA/emission.${PDY}.m6.bin
#---------------------------------------------------------------------
# Execute PREMAQ.
#---------------------------------------------------------------------
export pgm=aqm_premaq
. prep_step

startmsg
$EXECaqm/aqm_premaq_5xwrf >> $pgmout 2>errfile
export err=$?;err_chk


let "soilhr=soilhr-1"
typeset -Z2 soilhr
if [ "$SENDCOM" = 'YES' ]
then
     cp $DATA/aqm.${cycle}.nmm00.tm00 $COMOUT/aqm.${cycle}.nmm00.tm00
     cp $DATA/GRIDDESC  $COMOUT/aqm_griddesc05
     cp $DATA/GRIDCRO2D $COMOUT/aqm.${cycle}.grdcro2d.bin
     cp $DATA/GRIDDOT2D $COMOUT/aqm.${cycle}.grddot2d.bin
     cp $DATA/METBDY3D $COMOUT/aqm.${cycle}.metbdy3d.bin
     cp $DATA/METCRO2D $COMOUT/aqm.${cycle}.metcro2d.bin
     cp $DATA/METCRO3D $COMOUT/aqm.${cycle}.metcro3d.bin
     cp $DATA/METDOT3D $COMOUT/aqm.${cycle}.metdot3d.bin
     cp $DATA/CHEMBDY3D $COMOUT/aqm.${cycle}.chembdy3d.bin
     cp $DATA/CHEMCRO3D $COMOUT/aqm.${cycle}.chemcro3d.bin
     cp /com/output/${envir}/today/aqm_*${cyc}.* $COMOUT/.
     cp $DATA/emission.${PDY}.m6.bin            $COMOUT/aqm.${cycle}.emission.bin
     cp $DATA/soilout.${PDY}.${OUTGRID}.bin     $COMOUT/aqm.${cycle}.soil${soilhr}.bin
fi

echo EXITING $0

########################################################

msg='ENDED NORMALLY.'
#postmsg "$jlogfile" "$msg"

################## END OF SCRIPT #######################