:
#chfp_test            version: 20080612
set -euxv

#=======================================================================
#
#     Script chfp_test
#     ----------------
#
#     Purpose     : calculate monthly mean nino 3.4 values from CHFP atmos
#     -------       data and write them into a simple data base.
#
#     Method      : Will access data from the specified data source,
#     ------        either an OPenDAP server or a downloaded netCDF file.
#                   Either process all available data (up to a specified
#                   number of ensemble members), or only new data.
#
#     Usage       : chfp_test
#     -----
#                   Developed to work under ksh. May need modification to
#                   work with other shells. Make sure you set "Local settings"
#                   below to appropriate values before use. NCO software is
#                   available from http://nco.sourceforge.net/ .
#
#     Author      : Tim Stockdale, ECMWF
#     ------
#
#=======================================================================
set  +v

opendap=true
datasource=http://ensembles.ecmwf.int/thredds/dodsC/ensembles/stream1/atmospheric/monthly
#opendap=false
#datasource=CHFPtest.nc

var=ts          # surface temperature
maxens=3        # maximum number of ensemble members to retrieve (or 999 for all)
maxtimes=38     # maximum number of output times to process      (or 9999 for all)
overwrite=true  # false => only download data not already in database


#============ Local settings - please specify your own =================
use nco
WORK=$SCRATCH/chfp_test      # where to run the script
DBROOT=/tmp/net/dbtest       # where to create the nino database
#=======================================================================

mkdir -p $WORK
cd $WORK

# First find out which multi-model ensemble members are available on server

if [[ $opendap = true ]] ; then
  sourcecmd="-p $datasource $var"
else
  sourcecmd="$datasource"
fi

ncks -h -a -v institution,experiment_id,realization,reftime,leadtime \
 -O -o multimodel.nc $sourcecmd

ncks -H -a -v realization \
 -s '%d\n' multimodel.nc > real.list

ncks -H -a -v institution \
 -s '%c' multimodel.nc > i.list
fold -w 15 i.list > inst.list

ncks -H -a -v experiment_id \
 -s '%c' multimodel.nc > e.list
fold -w 4 e.list > expt.list

integer nmmens
nmmens=`wc -l < inst.list`
echo $nmmens
(( nmmens > maxens )) && nmmens=$maxens

ncks -H -a -v reftime \
 -s '%f\n' multimodel.nc > date.list

ncks -H -a -v leadtime \
 -s '%d\n' multimodel.nc > lead.list

integer ntimes
ntimes=`wc -l < date.list`
(( ntimes = ntimes -1 ))
echo $ntimes
(( ntimes > maxtimes )) && ntimes=$maxtimes

# Build list of start dates, to check against existing data

integer nexist
datenum=1
rm -f idate.list
touch idate.list
while [[ $datenum -le $ntimes ]] ; do
  date1=`sed "${datenum}!d" date.list`
  nexist=`wc -l < idate.list`
  nind=1          # general case
  nind=$nexist    # if not mixed up
  match=false
  while [[ $nind -le $nexist ]] ; do
    dexist=`sed "${nind}!d" idate.list `
    if [[ $date1 -eq $dexist ]] ; then
      match=true
    fi
    ((nind=nind+1))
  done
  if [[ $match = false ]] ; then
    echo "$date1" > latest
    cat latest >> idate.list
  fi
  ((datenum=datenum+1))
done

integer nidate
nidate=`wc -l < idate.list`

#----------------------------------------
# Loop over multi-model ensemble members
#----------------------------------------

mynum=1
while [[ $mynum -le $nmmens ]] ; do

# Set the ECMWF "nino database" variables, or abort if data not recognized

  ((odnum=mynum-1))
  head -$mynum inst.list > i1.list
  inst=`tail -1 i1.list`
  fdb_origin=unknown
  [[ $inst = "ECMWF          " ]] && fdb_origin=ecmf
  [[ $inst = "UK Met Office  " ]] && fdb_origin=egrr
  [[ $inst = "Meteo-France   " ]] && fdb_origin=lfpw
  [[ $inst = "IfM-Kiel       " ]] && fdb_origin=ifmk
  if [[ $fdb_origin = unknown ]] ; then
    echo "Institution not recognized"
    exit 1
  fi
  integer number
  head -$mynum real.list > r1.list
  number=`tail -1 r1.list`
  head -$mynum expt.list > e1.list
  expver=`tail -1 e1.list`
  typeset -Z3 number

  mkdir -p $DBROOT/nino/$fdb_origin/$expver/$number || true

#----------------------------------------
# Loop over forecast start dates
#----------------------------------------

  datenum=1
  while [[ $datenum -le nidate ]] ; do

    integer daycount
    daycount=`sed "${datenum}!d" idate.list `

# newdate is an ecmwf utility that we use to calculate the date
# in the form YYYYMMDD. For portability, we replace it with a simple but
# approximate calculation which assumes the initial date is always
# the first of the month (as it is with ENSEMBLES daa).
#   date2=`newdate -D 19500101 +$daycount `
    integer yyyy
    integer mm
    typeset -Z2 mm
    (( yyyy = daycount/365 ))
    (( mm = (daycount-365*yyyy)/29 + 1 ))
    (( yyyy = yyyy + 1950 ))
    date2=${yyyy}${mm}01

    match=false
    if [[ $overwrite = false ]] ; then
      [[ -s $DBROOT/nino/$fdb_origin/$expver/$number/${date2}_${var}amona ]] && match=true
    fi

    if [[ $match = false ]] ; then

# Find start and end "times" of requested data in netCDF file
      ts1=`grep -n $daycount date.list | head -1`
      ts2=`grep -n $daycount date.list | tail -1`
      time1=${ts1%:*}
      time2=${ts2%:*}
      ((time1=time1-1))
      ((time2=time2-1))
      datesubset=" -d time,$time1,$time2"

      ncks -h -a -v $var -d ensemble,$odnum $datesubset \
        -O -o nino_indata $sourcecmd

      ncwa -O -a latitude,longitude -d latitude,-5.0,5.0 -d longitude,190.0,240.0 nino_indata nino_outdata
      ncks -C -H -s '%f\n' -v $var nino_outdata > nino_txt

      mv nino_txt $DBROOT/nino/$fdb_origin/$expver/$number/${date2}_${var}amona

    fi

    (( datenum=datenum+1))
  done   # end of loop over dates

  (( mynum = mynum + 1 ))
done     # end of loop over multi-model ensemble members