#!/bin/bash
#
EXPNO=test31
EXPERIMENTSYEAR=2000
EXPERIMENTEYEAR=2001

# Restart:
RESTARTYEAR=2000                  # Set to ${EXPERIMENTSYEAR} for initial run 
RESTARTMONTH=01                   # Month (01, 02 etc.);  Set to 01 for initial run
#
EXPERIMENTEYEAR1=`expr ${EXPERIMENTEYEAR} + 1`

#set -ex
set -e
export MPLARGS="verbose=1;check=0"
export VPP_MBX_SIZE=64000000
export VPP_STATS=10
#export GMON_OUT_PREFIX=gmon.out-
#
# Job file to run echam model on brutus cluster
# =============================================
#
F_RECLUNIT=BYTE ; export F_RECLUNIT
MPIPROGINF=detail ; export MPIPROGINF
F_SYSLEN=600 ; export F_SYSLEN
#
NCPUS=48
NPROCA=8
NPROCB=6
#
NPROMA=16
#
RES=42
LEV=39
#
# Pressure levels for output files <EXPNO>_chem_m_<YYYY><MM>.nc
# =============================================================
rm -f afterburner_levels.nml
cat > afterburner_levels.nml << EOL 
&AFTERBURNER_LEVELS
  LEVELS = 1, 5, 10, 30, 50, 100, 160, 250, 400, 630, 1000,  1585, 2512,  3981, 
  6310, 10000, 15000, 20000, 25000, 35000, 45000, 55000, 65000, 75000, 80000, 
  85000, 90000, 95000, 98000, 100000 /
EOL
#
LSOCOL=.TRUE.     # Switch for SOCOL model
LCHEM=.TRUE.      # Switch for chemistry module
#
ECHAM5DIR=/cluster/home/stenkea/SOCOLv3.1
WORKDIR=/cluster/work/climate/stenkea/echam5
RUNDIR=/cluster/scratch/stenkea/echam5
DAT=$RUNDIR/run/run$EXPNO
DPATH=$DAT/

if [ ! -d $DAT ]
then
  mkdir $DAT
fi

#
################################################################################
#
cd $DPATH               # output and rerun files are written into $DPATH
#
MODEL=$ECHAM5DIR/bin/echam5  # <--- path to the executable
#
INI=$WORKDIR/data/T${RES}
INIAMIP=$WORKDIR/data/T${RES}/amip2
INIHADLEY=$WORKDIR/data/T${RES}/hadley
INISOCOL=$WORKDIR/data/SOCOL/T${RES}
#
################################################################################
#
# specification of file structures
#
# stop execution after the first run time error
F_ERRCNT=0
export F_ERRCNT
#
################################################################################
#
rm -f unit.?? sst* ice* rrtadata *.codes atmout chem_initial
#
ln -s  ${INI}/T${RES}L${LEV}_jan_spec.nc      unit.23
ln -s  ${INI}/T${RES}_jan_surf.nc             unit.24
#
################################################################################
#
#  Initial file for chemical species (MEZON):

ln -s  ${INISOCOL}/chem_initial/T${RES}L${LEV}_restc_197412+tropchem+ods.nc  chem_initial
#
################################################################################
## for climatological sst and ice (LAMIP=F) use:
#ln -s  ${INIAMIP}/T${RES}_amip2sst_clim.nc    unit.20
#ln -s  ${INIAMIP}/T${RES}_amip2sic_clim.nc    unit.96
##
# for AMIP (variable) sst and ice (LAMIP=T) use:

YY=`expr ${RESTARTYEAR} - 1`
while [ ${YY} -le ${EXPERIMENTEYEAR1} ]; do
    ln -s  ${INIHADLEY}/T${RES}_HadISST_sst_${YY}.nc sst${YY}
    ln -s  ${INIHADLEY}/T${RES}_HadISST_sic_${YY}.nc ice${YY}
    YY=`expr ${YY} + 1`
done
#
################################################################################
#
ln -s  ${INI}/T${RES}_O3clim2.nc    unit.21
ln -s  ${INISOCOL}/misc/OH_ACCENT_T${RES}.nc unit.22
ln -s  ${INI}/T${RES}_VLTCLIM.nc    unit.90
ln -s  ${INI}/T${RES}_VGRATCLIM.nc  unit.91
ln -s  ${INI}/T${RES}_TSLCLIM2.nc   unit.92
ln -s  ${INI}/surrta_data           rrtadata
#
################################################################################
#
# Boundary conditions for SOCOL:

# Greenhouse gases and ODSs:
############################
rm -f co2_y ch4_y n2o_y odsmem_y odscll_y odscls_y odsbr_y

ln -s  ${INISOCOL}/../GHG_ODS/co2_socol_1950-2100.txt         co2_y
ln -s  ${INISOCOL}/../GHG_ODS/ch4_socol_1950-2100.txt         ch4_y
ln -s  ${INISOCOL}/../GHG_ODS/n2o_socol_1950-2100.txt         n2o_y
ln -s  ${INISOCOL}/../GHG_ODS/odsmem_socol_1950-2100.txt      odsmem_y
ln -s  ${INISOCOL}/../GHG_ODS/odscls_socol_1950-2100.txt      odscls_y
ln -s  ${INISOCOL}/../GHG_ODS/odscll_socol_1950-2100.txt      odscll_y
ln -s  ${INISOCOL}/../GHG_ODS/odsbr_socol_1950-2100.txt       odsbr_y

# Solar irradiation, Schumann-Runge bands/Lyman alpha line heating 
# parameterisation, lookup-tables for photolysis rates:
#######################################################
rm -f sun_irrad_y sun_par_y photolysis???? photolysis_mean photolysis_delta_e_corr 

ln -s ${INISOCOL}/../CCMI_solar_data/irrad_echam5/sun_irrad_C1_1950_2011.txt       sun_irrad_y 
ln -s ${INISOCOL}/../CCMI_solar_data/extra_heating/sun_par_C1_1950-2011.txt        sun_par_y

YY=`expr ${RESTARTYEAR} - 1`
while [ ${YY} -le ${EXPERIMENTEYEAR1} ]; do
  ln -s ${INISOCOL}/../CCMI_solar_data/photolysis/phot_C1_${YY}.nc photolysis${YY} 
  YY=`expr ${YY} + 1`
done

ln -s ${INISOCOL}/../CCMI_solar_data/photolysis/phot_C1_mean_1977-1998.nc photolysis_mean # mean 1979-1998 (two solar cycles)

# Delta-E photolysis corrections:
ln -s ${INISOCOL}/../sun/L${LEV}_photolysis_delta_e_corr_socol.nc photolysis_delta_e_corr

#CCMI
# Ionization data
#######################################################
rm -f geomag_data_1600_2100 ????_ionrate-nitrate.txt monthly_mean_ap gcr_data

ln -s ${INISOCOL}/../ionization/geomag_data_1600_2100.txt geomag_data_1600_2100
ln -s ${INISOCOL}/../ionization/gcr_data_1600_2100_constant.nc gcr_data
ln -s ${INISOCOL}/../CCMI_solar_data/ap_index/ap_annual_REF-C2.txt  monthly_mean_ap

YY=`expr ${RESTARTYEAR} - 1`
while [ ${YY} -le ${EXPERIMENTEYEAR1} ]; do
    ln -s ${INISOCOL}/../CCMI_solar_data/ionrate/${YY}_ionrate_spe_REF_C2.txt  ${YY}_ionrate-nitrate.txt
    YY=`expr ${YY} + 1`
done

# Stratospheric aerosols:
#########################
rm -f strataer*

YY=`expr ${RESTARTYEAR} - 1`
while [ ${YY} -le ${EXPERIMENTEYEAR1} ]; do
    ln -s  ${INISOCOL}/strat_aerosols/CMIP6_pw_v4/T${RES}L${LEV}_CMIP6_3l_v4.0.0_pw_${YY}.nc strataer${YY}
    YY=`expr ${YY} + 1`
done

ln -s  ${INISOCOL}/strat_aerosols/CMIP6_pw_v4/T${RES}L${LEV}_CMIP6_3l_v4.0.0_pw_2000.nc strataerbg

# Troposperic aerosol climatology:
##################################
rm -f tropoaer*

YY=`expr ${RESTARTYEAR} - 1`
 while [ ${YY} -le ${EXPERIMENTEYEAR1} ]; do
     ln -s  ${INISOCOL}/tropo_aerosols/T${RES}L${LEV}_tropoaero_socol_clim.nc tropoaer${YY} #option for interannually varying tropospheric aerosols
     YY=`expr ${YY} + 1`
 done
 
ln -s  ${INISOCOL}/tropo_aerosols/T${RES}L${LEV}_tropoaero_socol_clim.nc tropoaer

# CO and NOx emissions:
#######################
rm -f co_nox_emiss_surf_aircr???? nox_lightning lightn_fac

YY=`expr ${RESTARTYEAR} - 1`
while [ ${YY} -le ${EXPERIMENTEYEAR1} ]; do
    ln -s  ${INISOCOL}/co_nox/T${RES}L${LEV}_co_nox_emiss_surf_aircr_socol_${YY}_ccmi.nc co_nox_emiss_surf_aircr${YY}
    YY=`expr ${YY} + 1`
done
ln -s  ${INISOCOL}/co_nox/T${RES}L${LEV}_nox_lightning_socol.nc nox_lightning
#CCMI
ln -s  ${INISOCOL}/co_nox/T${RES}_fac_lightning_socol.nc lightn_fac #eth_af_update
#  Input files for tropospheric chemistry:
##########################################
rm -f  c5h8_emissions* ch2o_emissions* ch3cooh_emissions* hcooh_emissions photolysis_rates_messy 
 
  YY=`expr ${RESTARTYEAR} - 1`
  while [ ${YY} -le ${EXPERIMENTEYEAR1} ]; do
    ln -s ${INISOCOL}/tropchem/ACCMIP_emissions_Isoprene_T${RES}_2000.nc c5h8_emissions${YY} 
    ln -s ${INISOCOL}/tropchem/ACCMIP_emissions_Formaldehyde_T${RES}_2000.nc ch2o_emissions${YY} 
    ln -s ${INISOCOL}/tropchem/ACCMIP_emissions_CH3COOH_T${RES}_2000.nc ch3cooh_emissions${YY} 
    YY=`expr ${YY} + 1`
   done
 
#
ln -s ${INISOCOL}/tropchem/ACCMIP_emissions_Isoprene_T${RES}_2000.nc c5h8_emissions 
ln -s ${INISOCOL}/tropchem/ACCMIP_emissions_Formaldehyde_T${RES}_2000.nc ch2o_emissions 
ln -s ${INISOCOL}/tropchem/ACCMIP_emissions_CH3COOH_T${RES}_2000.nc ch3cooh_emissions
#
ln -s ${INISOCOL}/tropchem/IPCC_emissions_HCOOH_2000_T${RES}.nc hcooh_emissions 
#
ln -s ${INISOCOL}/tropchem/Photolysis_tropchem_T${RES}.nc photolysis_rates_messy
#
 
# NMVOC emissions:
########################
rm -f nmvoc_emiss_biogenic???? nmvoc_emiss_bb???? nmvoc_emiss_anthrop????
rm -f nmvoc_emiss_biogenic nmvoc_emiss_bb nmvoc_emiss_anthrop

YY=`expr ${RESTARTYEAR} - 1`
while [ ${YY} -le ${EXPERIMENTEYEAR1} ]; do
    ln -s ${INISOCOL}/nmhc/ACCMIP_emissions_biogenic_NMHC_T${RES}_2000.nc  nmvoc_emiss_biogenic${YY} 
    ln -s ${INISOCOL}/nmhc/ACCMIP_emissions_anthropogenic_NMHC_T${RES}_2000.nc  nmvoc_emiss_anthrop${YY} 
    ln -s ${INISOCOL}/nmhc/ACCMIP_emissions_biomassburning_NMHC_T${RES}_2000.nc   nmvoc_emiss_bb${YY} 
    YY=`expr ${YY} + 1`
done

    ln -s ${INISOCOL}/nmhc/ACCMIP_emissions_biogenic_NMHC_T${RES}_2000.nc nmvoc_emiss_biogenic 
    ln -s ${INISOCOL}/nmhc/ACCMIP_emissions_anthropogenic_NMHC_T${RES}_2000.nc nmvoc_emiss_anthrop 
    ln -s ${INISOCOL}/nmhc/ACCMIP_emissions_biomassburning_NMHC_T${RES}_2000.nc nmvoc_emiss_bb 
# QBO:
######
rm -f qbo????

YY=`expr ${RESTARTYEAR} - 1`
while [ ${YY} -le ${EXPERIMENTEYEAR1} ]; do
    ln -s  ${INISOCOL}/../QBO/QBO_socol_${YY}.nc qbo${YY}
    YY=`expr ${YY} + 1`
done

#  Input files for Drydep
rm -f soilpH surf_para
ln -s  ${INISOCOL}/misc/soilpH_T${RES}.nc  soilpH 
ln -s  ${INISOCOL}/misc/surf_para_T${RES}.nc  surf_para 
#
#
# Wetland parameterisation:
rm -f unit.40 unit.41 unit.43 unit.44

ln -s ${INISOCOL}/misc/global_wetland_frac_T${RES}.nc unit.40
ln -s ${INISOCOL}/misc/HR_T${RES}.nc unit.41
#ln -s tsurf.nc unit.43
ln -s ${INISOCOL}/misc/30008_1995_tsurf_annualmean.nc unit.43

# CMDL CH4 climatology
ln -s ${INISOCOL}/ch4/CH4_clim_T${RES}.nc unit.44


# Prepare rerun files (if necessary):

if [ ${RESTARTYEAR} -gt ${EXPERIMENTSYEAR} -o ${RESTARTMONTH} -gt 01 ]; then
    RERUN=.TRUE.     # Rerun switch; .false. for initial run, .true. for rerun
else
    RERUN=.FALSE.
fi
#
rm -f fortran_error_messages
#
################################################################################
#
# Set namelist and call model:  

if [ ! -s rerun_${EXPNO}_echam ]; then
    #YR=${RESTARTYEAR}
    YR=${EXPERIMENTSYEAR}
    echo "NO RERUN FILE"
else
    #YR=`ncdump -h rerun_${EXPNO}_echam | grep -i vdate | cut -c 12-15`
    YR=$1
fi

echo ${YR}

    for MM in 01 02 03 04 05 06 07 08 09 10 11 12; do

	if [ "${YR}" -gt ${RESTARTYEAR} -o "${MM}" -gt ${RESTARTMONTH} ]; then
	    RERUN=.TRUE.
	fi

        if [ "${YR}" -eq ${RESTARTYEAR} -a "${MM}" -eq ${RESTARTMONTH} ]; then
          CO2=1.00 #For running ensembles 
        else
          CO2=1.00
        fi
		
#  namelist control variables and output control for grid variables
#  spectral variables are written out by default except liq. water
#  for production runs set LABORT=.FALSE.
#
	    cat > namelist.echam << EOF
&RUNCTL
  LSOCOL=${LSOCOL}
  LRESUME=$RERUN,
  OUT_DATAPATH = "$DPATH"
  OUT_EXPNAME  = "$EXPNO"
  OUT_FILETYPE = 2
  TRAC_FILETYPE = 2
  DT_START  = ${EXPERIMENTSYEAR},01,01,12,0,0
  DT_STOP   = ${EXPERIMENTSYEAR},01,31,12,0,0
  NO_CYCLES = 1
  PUTDATA   = 12, 'hours', 'first', 0
  PUTRERUN  = 1, 'months', 'last', 0
  DELTA_TIME = 900
  LAMIP=.TRUE.
  LMIDATM=.TRUE.
  LABORT=.FALSE.
  NPROCA=${NPROCA}
  NPROCB=${NPROCB}
  NPROMA=${NPROMA}
/
&SOCOLCTL
  LCHEM=${LCHEM}
  CYEAR=${CYEAR:-${EXPERIMENTSYEAR}}
  CO2FAC=${CO2}
  lqbonudg=.true.
  nlev_mbc=1
  lsrb_lya_heating = .TRUE.
/
&DYNCTL
  VCHECK=235.
  SPDRAG=0.926E-4
/
&RADCTL
  IAERO = 10
  LGADSRH=.TRUE.
/
EOF

  # Call model:
    date
    mpirun $MODEL

    RERUN=.true.
    
   done

if [ ${YR} -eq ${EXPERIMENTEYEAR} ]; then
    exit
fi

ls -lt

cd $ECHAM5DIR/run

YR=`expr ${YR} + 1`

echo ${YR}

bsub -n $NCPUS -J ${EXPNO} -W 11:59 "./run_socol_test ${YR}"

exit
