#!/bin/bash
# script for starting xleap.sh
# Stefan Henrich 29.03.2005
#
version="combine2go 13.09.2005"
# 05.09.2005	minimization of multiple ligands in separated files (pdb) 
#				together with a receptor model (pdb);
#				getting netcharge of SD file
# 13.09.2005	minimization of multiple docking solutions in one sdf file 
#				together with receptor model (pdb); command line options
#
# give pdb as argument
#
# required programs:
# grepsdtag.linux
# xleap.sh
# select_sdf
#
# local_start.sh | tee xleap.log
###########################################
COMBINEHOME="/home/henricsn/combine2go/"
sdf=$COMBINEHOME"sdf/PDB_uPA_300805_3D.sdf"              #library SD file
template=$COMBINEHOME"parameter/protmin.dielc4.template" #directory for minimize template
prgdir=$COMBINEHOME"programs/"                           #directory for programs
outdir="./"
#sdf=$COMBINEHOME"sdf/final-residue.sdf"                  #directory for SD file
##################################################################
##################################################################
# DO NOT CHANGE ANYTHING BELOW THIS LINE
# UNLESS YOU KNOW WHAT YOU ARE DOING
##################################################################
##################################################################
until [ -z "$1" ]; do
	case $1 in
		--help|-h )
			echo options:
			echo "-h       help"
			echo "-vers    version ($version)"
			echo "-sdf     library SD file ($sdf)"
			echo "-proj    project path ($COMBINEHOME)"
			echo "-out     output directory ($outdir)"
			echo "-inrec   input directory of receptor ($indir_receptor)"
			echo "-inlig   input directory of ligand ($indir_ligand)"
			echo "-model   receptor model (pdb file)"
			echo "-ligbase ligand base name (eg. *.prot.ligand.align.pdb)"
			echo "-multi   SD file with multiple GOLD docking solutions"
			echo "-list    list of ligand pdb files (eg. \"A01.pdb A02.pdb A03.pdb\")" 
			echo "-min     minimization protocol ($template)"
			exit 0
		;;
		--version|-vers )
			echo $version
			exit 0
		;;
		--librarysdf|-sdf )
			#library SD file
			shift
			sdf=$1
		;;
		--projectfolder|-proj )
			shift
			COMBINEHOME=$1
		;;
		--output|-out )
			shift
			outdir=$1
		;;
		--indirreceptor|-inrec )
			shift
			indir_receptor=$1
		;;
		--indirligand|-inlig )
			shift
			indir_ligand=$1
		;;
		--model|-model )
			shift
			model=$1
		;;
		--ligbase|-ligbase )
			shift
			ligbase=$1
		;;
		--multisolutions|-multi )
			shift
			multi=$1
		;;
		--list|-list )
			shift
			list=$1
		;;
		--minimize|-min )
			shift
			template=$1
		;;
		* ) 
			echo no command is given
		;;
	esac 

	shift

	if [ "$#" = "0" ]; then
		break  
	fi
done 

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

indir_receptor=${indir_receptor:=$COMBINEHOME"data/urokinase/pdb/model/"}
indir_ligand=${indir_ligand:=$COMBINEHOME"data/urokinase/pdb/ligand/"}
model=${model:=$indir_receptor"model5.pdb"}                        #receptor model
##################################################################
###### search for ligands with $ligbase
if [ "$ligbase" != "" ]; then
	getmodifpdb=$(find $indir_ligand -name "*$ligbase")
fi #ligbase

###### multiple ligands in SD file
if [ "$multi" != "" ]; then
	maxsdf=$(grep '$$$$' $multi|wc|awk '{ print $1 } ')
	echo number of ligands: $maxsdf
	counter=1
	until [ "$maxsdf" -le "$counter" ]; do
		echo select_sdf -start $counter -end $counter \< $multi \> temp.sdf
		select_sdf -start $counter -end $counter < $multi > temp.sdf
		### calculate netcharge
		charge_numb=0;netcharge=0;charge=0
		chargeline=$(grep "M  CHG" temp.sdf)
		charged_atoms=`echo $chargeline| awk '{ print $3 } '`
		until [ "$charged_atoms" -le "$charge_numb" ]; do
			column=$( echo "5 + $charge_numb * 2" |bc -l)
			#variable column of bash will be used in awk
			charge=`echo $chargeline| awk -v column="$column" '{ total = $column } END {print $column}'`
			netcharge=$( echo "$netcharge + $charge" |bc -l)
			let charge_numb+=1
		done
		residue=$(extract_prop_sdf -property_name residue < temp.sdf)
		echo residue $residue netcharge $netcharge
		if [ "$residue" = "undef" ]; then
			id=$(extract_prop_sdf -property_name ID < temp.sdf)
			echo $id
			if [ $id = "undef"  ]; then
				residue="xxx"
			else
				id="000"$id
				#take the last three numbers of id
				residue=${res:$(echo "${#id} -3"|bc -l):3}
			fi
			id=""
		fi #ID

		### start xleap.sh
		echo $prgdir"xleap.sh"	-gold -d $outdir -m $template -notemp -bres \
								-prot $model -pdb $residue -alt $counter -n $residue  \
								-sdf temp.sdf -nc $netcharge

		$prgdir"xleap.sh"		-gold -d $outdir -m $template -notemp -bres \
								-prot $model -pdb $residue -alt $counter -n $residue  \
								-sdf temp.sdf -nc $netcharge
		

		let counter+=1

	done #maxsdf

	###
	#rm -f sol_*.temp
##################################################################
else

	FILES="$getmodifpdb"
	
	###### list of ligand pdb files
	if [ "$list" != "" ]; then
		FILES=$list
	fi #list
	
	### here can you give a list of files:
	#FILES="$indir_ligand/A04.pdb \
	#"
	######

	for file in $FILES
	do
		echo
		echo ---------------------------------------------------------------------------
		echo
		
		###grep residue name of ligand pdb file
		line=$(grep "^HETATM" $file|uniq -w6)
		residue=${line:17:3}
		
		###create sdf and calculate netcharge
		$prgdir"grepsdtag.linux" -t residue $residue < $sdf > temp.sdf
		charge_numb=0;netcharge=0;charge=0
		chargeline=$(grep "M  CHG" temp.sdf)
		charged_atoms=`echo $chargeline| awk '{ print $3 } '`
		until [ "$charged_atoms" -le "$charge_numb" ]; do
			column=$( echo "5 + $charge_numb * 2" |bc -l)
			#variable column of bash will be used in awk
			charge=`echo $chargeline| awk -v column="$column" '{ total = $column } END {print $column}'`
			netcharge=$( echo "$netcharge + $charge" |bc -l)
			let charge_numb+=1
		done
		
		###name for output files
		counter=0
		while [ -e $residue"."$counter".min.pdb" ]
		do 
			let counter+=1
		done
		echo  $residue netcharge $netcharge counter $counter
		
		
		###start xleap.sh
		echo $prgdir"xleap.sh" -d $outdir -pdb $residue -n $residue -prot $model -lig $file \
		-alt $counter -nc $netcharge -sdf temp.sdf -m $template -w -notemp -bres -noedit
		$prgdir"xleap.sh" -d $outdir -pdb $residue -n $residue -prot $model -lig $file \
		-alt $counter -nc $netcharge -sdf temp.sdf -m $template -w -notemp -bres -noedit
		
		#echo > $residue"."$counter".min.pdb"
	done
		rm -f temp.sdf
fi #multi

	exit

#grep -A6 "loading mimimized molecule in pymol:" sxleap010805_2.log > rmsd3.log