
#################################### FORECAST DEBIAS ###############################################
echo "------------------------------------------------"
echo "Bias Correct NCEP Global Ensemble Forecast & GFS Forecast"
echo "------------------------------------------------"
echo "History: March 2006 - First implementation of this new script"
echo "         March 2007 - Modified for GEFS upgrade (20 members)"
echo "         May   2007 - Add GFS bias correction"
echo "         Aug   2007 - Add hybrid method to combine bias-corrected GFS and bias-corrected GEFS"
echo "AUTHOR: Bo Cui  (wx20cb)"
####################################################################################################

### To submit this job for T00Z, T06Z T12Z and T18Z, four cycles per day
                                                                                                                    
### need pass the values of PDY, CYC, DATA, COMIN, COM, COMOUTBC, COMOUTAN and COMOUTWT

################################################################
# define exec variable, and entry grib utility 
################################################################

export NDATE=$EXECUTIL/ndate

export ENSANOMALY=$USHGEFS/gefs_climate_anomaly.sh
export ENSWEIGHTS=$USHGEFS/gefs_weights.sh

export pgm=gefs_debias
. prep_step

########################################################
### define the days for searching bias estimation backup
########################################################
###
ymdh=${PDY}${cyc}
export PDYm8=`$NDATE -192 $ymdh | cut -c1-8`
export PDYm9=`$NDATE -216 $ymdh | cut -c1-8`
export PDYm10=`$NDATE -240 $ymdh | cut -c1-8`
export PDYm11=`$NDATE -264 $ymdh | cut -c1-8`
export PDYm12=`$NDATE -288 $ymdh | cut -c1-8`
export PDYm13=`$NDATE -312 $ymdh | cut -c1-8`
export PDYm14=`$NDATE -336 $ymdh | cut -c1-8`
export PDYm15=`$NDATE -360 $ymdh | cut -c1-8`
export PDYm16=`$NDATE -384 $ymdh | cut -c1-8`
export PDYm17=`$NDATE -408 $ymdh | cut -c1-8`
export PDYm18=`$NDATE -432 $ymdh | cut -c1-8`

##########################################################################
# bias correct NCEP global ensemble for each forecast time and each member
##########################################################################

#export memberlist="p01 p02 p03 p04 p05 p06 p07 p08 p09 p10 p11 p12 p13 p14 p15 p16 p17 p18 p19 p20"

for nens in $memberlist
do

  hourlist=" 00  06  12  18  24  30  36  42  48  54  60  66  72  78  84  90  96 \
            102 108 114 120 126 132 138 144 150 156 162 168 174 180 186 192 198 \
            204 210 216 222 228 234 240 246 252 258 264 270 276 282 288 294 300 \
            306 312 318 324 330 336 342 348 354 360 366 372 378 384"

  if [ "$nens" = "gfs" ]; then
    hourlist=" 00  06  12  18  24  30  36  42  48  54  60  66  72  78  84  90  96 \
              102 108 114 120 126 132 138 144 150 156 162 168 174 180 "
  fi

for nfhrs in $hourlist
  do

###
#  set the index ( exist of bias estimation ) as default, 0
###

    cstart_ens=0
    cstart_gfs=0
    if_gfs=0

###
#  GEFS bias estimation entry, if_gfs: 1 = for GFS high resolution forecast
###

    ibias_ens=geavg.t${cyc}z.pgrba_mef${nfhrs}

    if [ "$nens" = "gfs" ]; then
      ibias_ens=gegfs.t${cyc}z.pgrba_mef${nfhrs}
      if_gfs=1
    fi

    if [ -s $COM/gefs.$PDY/${cyc}/pgrba/$ibias_ens ]; then
      cp $COM/gefs.$PDY/${cyc}/pgrba/$ibias_ens $ibias_ens
    elif [ -s $COM/gefs.$PDYm1/${cyc}/pgrba/$ibias_ens ]; then
      cp $COM/gefs.$PDYm1/${cyc}/pgrba/$ibias_ens $ibias_ens
    elif [ -s $COM/gefs.$PDYm2/${cyc}/pgrba/$ibias_ens ]; then
      cp $COM/gefs.$PDYm2/${cyc}/pgrba/$ibias_ens $ibias_ens
    elif [ -s $COM/gefs.$PDYm3/${cyc}/pgrba/$ibias_ens ]; then
      cp $COM/gefs.$PDYm3/${cyc}/pgrba/$ibias_ens $ibias_ens
    elif [ -s $COM/gefs.$PDYm4/${cyc}/pgrba/$ibias_ens ]; then
      cp $COM/gefs.$PDYm4/${cyc}/pgrba/$ibias_ens $ibias_ens
    elif [ -s $COM/gefs.$PDYm5/${cyc}/pgrba/$ibias_ens ]; then
      cp $COM/gefs.$PDYm5/${cyc}/pgrba/$ibias_ens $ibias_ens
    elif [ -s $COM/gefs.$PDYm6/${cyc}/pgrba/$ibias_ens ]; then
      cp $COM/gefs.$PDYm6/${cyc}/pgrba/$ibias_ens $ibias_ens
    elif [ -s $COM/gefs.$PDYm7/${cyc}/pgrba/$ibias_ens ]; then
      cp $COM/gefs.$PDYm7/${cyc}/pgrba/$ibias_ens $ibias_ens
    elif [ -s $COM/gefs.$PDYm8/${cyc}/pgrba/$ibias_ens ]; then
      cp $COM/gefs.$PDYm8/${cyc}/pgrba/$ibias_ens $ibias_ens
    elif [ -s $COM/gefs.$PDYm9/${cyc}/pgrba/$ibias_ens ]; then
      cp $COM/gefs.$PDYm9/${cyc}/pgrba/$ibias_ens $ibias_ens
    elif [ -s $COM/gefs.$PDYm10/${cyc}/pgrba/$ibias_ens ]; then
      cp $COM/gefs.$PDYm10/${cyc}/pgrba/$ibias_ens $ibias_ens
    elif [ -s $COM/gefs.$PDYm11/${cyc}/pgrba/$ibias_ens ]; then
      cp $COM/gefs.$PDYm11/${cyc}/pgrba/$ibias_ens $ibias_ens
    elif [ -s $COM/gefs.$PDYm12/${cyc}/pgrba/$ibias_ens ]; then
      cp $COM/gefs.$PDYm12/${cyc}/pgrba/$ibias_ens $ibias_ens
    elif [ -s $COM/gefs.$PDYm13/${cyc}/pgrba/$ibias_ens ]; then
      cp $COM/gefs.$PDYm13/${cyc}/pgrba/$ibias_ens $ibias_ens
    elif [ -s $COM/gefs.$PDYm14/${cyc}/pgrba/$ibias_ens ]; then
      cp $COM/gefs.$PDYm14/${cyc}/pgrba/$ibias_ens $ibias_ens
    elif [ -s $COM/gefs.$PDYm15/${cyc}/pgrba/$ibias_ens ]; then
      cp $COM/gefs.$PDYm15/${cyc}/pgrba/$ibias_ens $ibias_ens
    elif [ -s $COM/gefs.$PDYm16/${cyc}/pgrba/$ibias_ens ]; then
      cp $COM/gefs.$PDYm16/${cyc}/pgrba/$ibias_ens $ibias_ens
    else
      echo " There is no Bias Estimation at " ${nfhrs} 
      cstart_ens=1
    fi

###
#  GFS bias estimation entry
###

    ibias_gfs=gegfs.t${cyc}z.pgrba_mef${nfhrs}

    if [ -s $COM/gefs.$PDY/${cyc}/pgrba/$ibias_gfs ]; then
      cp $COM/gefs.$PDY/${cyc}/pgrba/$ibias_gfs $ibias_gfs
    elif [ -s $COM/gefs.$PDYm1/${cyc}/pgrba/$ibias_gfs ]; then
      cp $COM/gefs.$PDYm1/${cyc}/pgrba/$ibias_gfs $ibias_gfs
    elif [ -s $COM/gefs.$PDYm2/${cyc}/pgrba/$ibias_gfs ]; then
      cp $COM/gefs.$PDYm2/${cyc}/pgrba/$ibias_gfs $ibias_gfs
    elif [ -s $COM/gefs.$PDYm3/${cyc}/pgrba/$ibias_gfs ]; then
      cp $COM/gefs.$PDYm3/${cyc}/pgrba/$ibias_gfs $ibias_gfs
    elif [ -s $COM/gefs.$PDYm4/${cyc}/pgrba/$ibias_gfs ]; then
      cp $COM/gefs.$PDYm4/${cyc}/pgrba/$ibias_gfs $ibias_gfs
    elif [ -s $COM/gefs.$PDYm5/${cyc}/pgrba/$ibias_gfs ]; then
      cp $COM/gefs.$PDYm5/${cyc}/pgrba/$ibias_gfs $ibias_gfs
    elif [ -s $COM/gefs.$PDYm6/${cyc}/pgrba/$ibias_gfs ]; then
      cp $COM/gefs.$PDYm6/${cyc}/pgrba/$ibias_gfs $ibias_gfs
    elif [ -s $COM/gefs.$PDYm7/${cyc}/pgrba/$ibias_gfs ]; then
      cp $COM/gefs.$PDYm7/${cyc}/pgrba/$ibias_gfs $ibias_gfs
    elif [ -s $COM/gefs.$PDYm8/${cyc}/pgrba/$ibias_gfs ]; then
      cp $COM/gefs.$PDYm8/${cyc}/pgrba/$ibias_gfs $ibias_gfs
    elif [ -s $COM/gefs.$PDYm9/${cyc}/pgrba/$ibias_gfs ]; then
      cp $COM/gefs.$PDYm9/${cyc}/pgrba/$ibias_gfs $ibias_gfs
    elif [ -s $COM/gefs.$PDYm10/${cyc}/pgrba/$ibias_gfs ]; then
      cp $COM/gefs.$PDYm10/${cyc}/pgrba/$ibias_gfs $ibias_gfs
    elif [ -s $COM/gefs.$PDYm11/${cyc}/pgrba/$ibias_gfs ]; then
      cp $COM/gefs.$PDYm11/${cyc}/pgrba/$ibias_gfs $ibias_gfs
    elif [ -s $COM/gefs.$PDYm12/${cyc}/pgrba/$ibias_gfs ]; then
      cp $COM/gefs.$PDYm12/${cyc}/pgrba/$ibias_gfs $ibias_gfs
    elif [ -s $COM/gefs.$PDYm13/${cyc}/pgrba/$ibias_gfs ]; then
      cp $COM/gefs.$PDYm13/${cyc}/pgrba/$ibias_gfs $ibias_gfs
    elif [ -s $COM/gefs.$PDYm14/${cyc}/pgrba/$ibias_gfs ]; then
      cp $COM/gefs.$PDYm14/${cyc}/pgrba/$ibias_gfs $ibias_gfs
    elif [ -s $COM/gefs.$PDYm15/${cyc}/pgrba/$ibias_gfs ]; then
      cp $COM/gefs.$PDYm15/${cyc}/pgrba/$ibias_gfs $ibias_gfs
    elif [ -s $COM/gefs.$PDYm16/${cyc}/pgrba/$ibias_gfs ]; then
      cp $COM/gefs.$PDYm16/${cyc}/pgrba/$ibias_gfs $ibias_gfs
    else
      echo " There is no GFS Bias Estimation at " ${nfhrs} 
      cstart_gfs=1
    fi

    echo "&message"  >input.$nfhrs.$nens
    echo " icstart_ens=${cstart_ens}," >> input.$nfhrs.$nens
    echo " icstart_gfs=${cstart_gfs}," >> input.$nfhrs.$nens
    echo " if_gfs=${if_gfs}," >> input.$nfhrs.$nens
    echo " nfhr=${nfhrs}," >> input.$nfhrs.$nens
    echo "/" >>input.$nfhrs.$nens

###
#  GEFS bias corrected ensemble forecasting output
###

    ofile=ge${nens}.t${cyc}z.pgrba_bcf${nfhrs}

###
#  check GEFS, GFS and GEFS control forecast file
###
   
    ifile_ens=$COMIN/ge${nens}.t${cyc}z.pgrbaf${nfhrs} 
    ifile_ctl=$COMIN/gec00.t${cyc}z.pgrbaf${nfhrs} 
    ifile_gfs=$COMIN/gegfs.t${cyc}z.pgrbaf${nfhrs} 

    icnt=0
    while [ icnt -le 12 ]; do
      if [ -s $ifile_ens -a -s $ifile_ctl ]; then

        ln -sf $ibias_ens fort.11
        ln -sf $ifile_ens fort.12
        ln -sf $ibias_gfs fort.13
        ln -sf $ifile_gfs fort.14
        ln -sf $ifile_ctl fort.15
        ln -sf $ofile     fort.51

        startmsg
        $EXECGEFS/$pgm   <input.$nfhrs.$nens     > $pgmout.$nfhrs.$nens 2> errfile
        export err=$?;err_chk
 	rm fort.*
        export MEMLIST=$nens
        $ENSANOMALY $PDY$cyc $nfhrs
        if [ "$nens" != "gfs" ]; then
          $ENSWEIGHTS $PDY$cyc $nfhrs
        fi
        icnt=13

#### sendcom  bias corrected forecast
#### and convert to grib2

        if [ "$SENDCOM" = "YES" ]; then
          if [ -s ge${nens}.t${cyc}z.pgrba_bcf${nfhrs} ]; then  
            mv ge${nens}.t${cyc}z.pgrba_bcf${nfhrs} $COMOUTBC/
            if [ ${nfhrs} = 240 ]; then
              $SMSBIN/setev pgrba_bcf240_ge${nens}
            fi
          $EXECUTIL/cnvgrib -g12 -p40 $COMOUTBC/ge${nens}.t${cyc}z.pgrba_bcf${nfhrs} $COMOUTBC_GB2/ge${nens}.t${cyc}z.pgrb2a_bcf${nfhrs}
	  fi
          if [ -s ge${nens}.t${cyc}z.pgrba_anf$nfhrs ]; then
            mv ge${nens}.t${cyc}z.pgrba_anf$nfhrs $COMOUTAN/
          $EXECUTIL/cnvgrib -g12 -p40 $COMOUTAN/ge${nens}.t${cyc}z.pgrba_anf$nfhrs $COMOUTAN_GB2/ge${nens}.t${cyc}z.pgrb2a_anf$nfhrs
          fi
          if [ -s ge${nens}.t${cyc}z.pgrba_wtf$nfhrs ]; then
            mv ge${nens}.t${cyc}z.pgrba_wtf$nfhrs $COMOUTWT/
          $EXECUTIL/cnvgrib -g12 -p40 $COMOUTWT/ge${nens}.t${cyc}z.pgrba_wtf$nfhrs $COMOUTWT_GB2/ge${nens}.t${cyc}z.pgrb2a_wtf$nfhrs
          fi
        fi

      else

        sleep 10
        icnt=`expr $icnt + 1`
	echo $icnt

      fi
    done

  done
done

###
#  final check up and make up
###

for nens in $memberlist
do

  hourlist=" 00  06  12  18  24  30  36  42  48  54  60  66  72  78  84  90  96 \
            102 108 114 120 126 132 138 144 150 156 162 168 174 180 186 192 198 \
            204 210 216 222 228 234 240 246 252 258 264 270 276 282 288 294 300 \
            306 312 318 324 330 336 342 348 354 360 366 372 378 384"

  if [ "$nens" = "gfs" ]; then
    hourlist="00  06  12  18  24  30  36  42  48  54  60  66  72  78  84  90  96 \
              102 108 114 120 126 132 138 144 150 156 162 168 174 180 "
  fi

  for nfhrs in $hourlist
  do
    if [ ! -s $COMOUTBC/ge${nens}.t${cyc}z.pgrba_bcf${nfhrs} ]; then  

      rm fort.*
      ibias_ens=geavg.t${cyc}z.pgrba_mef${nfhrs}
      if [ "$nens" = "gfs" ]; then
      ibias_ens=gegfs.t${cyc}z.pgrba_mef${nfhrs}
      fi
      ibias_gfs=gegfs.t${cyc}z.pgrba_mef${nfhrs}
      ifile_ens=$COMIN/ge${nens}.t${cyc}z.pgrbaf${nfhrs} 
      ifile_ctl=$COMIN/gec00.t${cyc}z.pgrbaf${nfhrs} 
      ifile_gfs=$COMIN/gegfs.t${cyc}z.pgrbaf${nfhrs} 
      ofile=ge${nens}.t${cyc}z.pgrba_bcf${nfhrs}

      ln -sf $ibias_ens fort.11
      ln -sf $ifile_ens fort.12
      ln -sf $ibias_gfs fort.13
      ln -sf $ifile_gfs fort.14
      ln -sf $ifile_ctl fort.15
      ln -sf $ofile     fort.51

      startmsg
      $EXECGEFS/$pgm   <input.$nfhrs.$nens    > $pgmout.$nfhrs.$nens 2> errfile
      export err=$?;err_chk
      if [ "$SENDCOM" = "YES" ]; then
        if [ -s $ofile ]; then  
          cp $ofile $COMOUTBC/
        fi
      fi
    fi
    if [ ! -s $COMOUTAN/ge${nens}.t${cyc}z.pgrba_anf$nfhrs ]; then
      rm fort.*
      export MEMLIST=$nens
      $ENSANOMALY $PDY$cyc $nfhrs
      if [ "$SENDCOM" = "YES" ]; then
        if [ -s ge${nens}.t${cyc}z.pgrba_anf$nfhrs ]; then
          mv ge${nens}.t${cyc}z.pgrba_anf$nfhrs $COMOUTAN/
          $EXECUTIL/cnvgrib -g12 -p40 $COMOUTAN/ge${nens}.t${cyc}z.pgrba_anf$nfhrs $COMOUTAN_GB2/ge${nens}.t${cyc}z.pgrb2a_anf$nfhrs
        fi
      fi
    fi
    if [ ! -s $COMOUTWT/ge${nens}.t${cyc}z.pgrba_wtf$nfhrs ]; then
      rm fort.*
      export MEMLIST=$nens
      if [ "$nens" != "gfs" ]; then
        $ENSWEIGHTS $PDY$cyc $nfhrs
      fi
      if [ "$SENDCOM" = "YES" ]; then
        if [ -s ge${nens}.t${cyc}z.pgrba_wtf$nfhrs ]; then
          mv ge${nens}.t${cyc}z.pgrba_wtf$nfhrs $COMOUTWT/
          $EXECUTIL/cnvgrib -g12 -p40 $COMOUTWT/ge${nens}.t${cyc}z.pgrba_wtf$nfhrs $COMOUTWT_GB2/ge${nens}.t${cyc}z.pgrb2a_wtf$nfhrs
        fi
      fi
    fi
  done
###
done

msg="HAS COMPLETED NORMALLY!"
postmsg "$jlogfile" "$msg"