#!/usr/bin/env bash
#zjxu@simm.ac.cn
#Create: 7/14/2015
#Function: ligprep to keep the initial 3D structure. ligprep was composed of applyhtreat, desalter, neutralizer, epik, stereoizer, ring_conf, premin, bmin step by step. To keep the origin 3D structure, only the first 4 steps shoud be retained.
#Update: remove -retain in neutralizer module. 7/14/2015
 
  USAGE()
  {
  echo "-------------------------------------------------------"
  echo "./01_keep3D.sh -l [ligand.mol2]"
  echo ""
  echo "for example:"
  echo "./01_keep3D.sh -l ligand.mol2"
  echo ""
  echo "-------------------------------------------------------"

  }

  program_state=0

  while getopts l: OPT
  do
  case $OPT in
  l) ligand=${OPTARG}
  echo "ligand file is ${ligand}"
  program_state=1
  ;;
  \?) USAGE
  exit
  ;;
  esac
  done

  if [ $program_state -eq 0 ]
  then USAGE; exit 1
  fi

  #--------------------------------------------------------------------------
  # do some tests

  if [ ! -f ${ligand} ]
  then
  echo "can't find file ${ligand}"
  USAGE
  exit
  fi

mol2convert -imol2 ${ligand} -omae ${ligand%.mol2}.mae
if [ ! -f ${ligand%.mol2}.mae ]
  then echo "mol2convert failed for ${ligand}"
  exit
fi

#applyhtreat -WAIT  ${ligand%.mol2}.mae ${ligand%.mol2}_htout.maegz -t 'All-atom with No-Lp' all 1 > ${ligand%.mol2}_ht.log 2>&1
applyhtreat ${ligand%.mol2}.mae ${ligand%.mol2}_htout.maegz 'All-atom with No-Lp' all 1 > ${ligand%.mol2}_ht.log 2>&1
if [ ! -f ${ligand%.mol2}_htout.maegz ]
  then echo "applyhtreat failed for ${ligand}"
  exit
fi

desalter ${ligand%.mol2}_htout.maegz ${ligand%.mol2}_dsout.maegz > ${ligand%.mol2}_ds.log 2>&1
if [ ! -f ${ligand%.mol2}_dsout.maegz ]
  then echo "desalter failed for ${ligand}"
  exit
fi

#neutralizer -retain ${ligand%.mol2}_dsout.maegz ${ligand%.mol2}_nuout.maegz > ${ligand%.mol2}_nu.log 2>&1
neutralizer ${ligand%.mol2}_dsout.maegz ${ligand%.mol2}_nuout.maegz > ${ligand%.mol2}_nu.log 2>&1
if [ ! -f ${ligand%.mol2}_nuout.maegz ]
  then echo "neutralizer failed for ${ligand}"
  exit
fi

#epik:
#  -ma <number>   Structures containing more than <number> atoms will not
#                 be adjusted.
#                 Default: 150
#  -ms <number>   Maximum number of generated structures per
#                 input structure.
#                 Default: 16
#  -nt            Do not tautomerize
#  -tn <number>   maximum number of tautomers.
#                 Default 8.
# -NO_REDIRECT      Do not use job control and print Epik messages to the
#                   screen.

#epik -WAIT -NO_REDIRECT -ph 7.0 -pht 0.0 -tn 8 -ma 200 -imae ${ligand%.mol2}_nuout.maegz -omae ${ligand%.mol2}_epikout.maegz > ${ligand%.mol2}_epik.log 2>&1
epik -WAIT -NO_REDIRECT -ph 7.35 -pht 0.0 -nt -ma 200 -ms 1 -imae ${ligand%.mol2}_nuout.maegz -omae ${ligand%.mol2}_epikout.mae > ${ligand%.mol2}_epik.log 2>&1
if [ ! -f ${ligand%.mol2}_epikout.mae ]
  then echo "epik failed for ${ligand}"
  exit
fi

echo "01_keep3D.sh for ${ligand} succeed."