#!/bin/csh
# This script is written by K.Hirata
# Please inform me if you find bugs in this script
# 2006.11.30 K.Hirata e-mail:hirata@spring8.or.jp
#
# Version up on 20080511 by K.Hirata
# Version up on 20080523 by K.Hirata
# -> ARPwARP model building cycle was added.
# -> Thanks to Dr. Makino
#
# What you must do
# Modify variants in the header of this script
# This script extracts cell parameters from 'merged' .sca file.
#

# Making scalepack script
# Parameters to be modified 
###################################################################
set SYMM='P43212'
set SCA_FILENAME=`ls *sca | sed -e 's/@//g'`
set PROJECT='iom'
set PHASE_MIN_RESOL=25
set PHASE_MAX_RESOL=1.8
set NUM_OF_TRY=100
set NUM_OF_FIND=3
set ANOM_ATOM=Br
set SOLVENT=0.42
set NUM_OF_DM=100
set NRESIDUES=129
set SEQ_FILE='/isilon/MOVE/BLBench/3.Lysozyme/lys.pir'
# End of parameters to be modified 
###################################################################

## for merged .sca file
set tmp=`head -3 $SCA_FILENAME|tail -1`
set CELL=`echo $tmp | awk '{print $1,$2,$3,$4,$5,$6}'`
echo $CELL

### SHELX part ####
shelxc $PROJECT << eof >shelxc.log
SAD $SCA_FILENAME
FIND $NUM_OF_FIND
SFAC $ANOM_ATOM
NTRY $NUM_OF_TRY
cell $CELL
SPAG $SYMM
SHELL $PHASE_MIN_RESOL $PHASE_MAX_RESOL
eof

echo "SHELXC finished"

shelxd ${PROJECT}_fa > shelxd.log

echo "SHELXD finished"
SHELXE:
shelxe ${PROJECT} ${PROJECT}_fa -s${SOLVENT} -m${NUM_OF_DM} > shelxe_o.log
shelxe ${PROJECT} ${PROJECT}_fa -s${SOLVENT} -m${NUM_OF_DM} -i > shelxe_i.log

#### PLOT shelxe results
grep cycle shelxe_i.log | grep wt | awk '{print $6}' | sed -e 's/,//g' > plot.dat
echo "" >> plot.dat
echo "" >> plot.dat
grep cycle shelxe_o.log | grep wt | awk '{print $6}' | sed -e 's/,//g' >> plot.dat

gnuplot << eof > gnuplot.log
set terminal postscript enhanced color eps
set output "shelx_plot.eps"
set title "SHELXE i/o"
set xlabel "cycle"
set ylabel "Contrast"
plot \
 "plot.dat" i 0 u 1 w lp lt 1 lw 3.0 title "Inverse",\
 "plot.dat" i 1 u 1 w lp lt 2 lw 3.0 title "Original"
eof

SELECTION:
echo GNUPLOT finished for SHELXE
## Select correct solution
set which=`choose_shelxe_chiral`
set filename=$PROJECT$which

if (-e ${filename}.phs) then
echo " The file ${filename} is available"
else
echo " The file ${filename} is not available (did you use the extension ?)"
echo ""
exit
endif

# Run f2mtz
echo " Running f2mtz"
f2mtz hklin ${filename}.phs hklout ${filename}.mtz <<END >${filename}_f2mtz.log
cell $CELL
symmetry $SYMM
labout H K L FP FOM PHIB SIGFP
CTYPOUT H H H F W P Q
END

grep "Number of Reflections " ${filename}_f2mtz.log > info
grep "Space Group" ${filename}_f2mtz.log >> info
echo -n " * Resolution Range: " >> info
grep " A )" ${filename}_f2mtz.log >> info
# Run cad
cad hklin1 ${filename}.mtz hklout ${filename}_cad.mtz > ${filename}_cad.log <<eof-cad
LABIN FILE 1 E1=FP E2=FOM E3=PHIB E4=SIGFP
LABOUT E1=FP E2=FOM E3=PHIB E4=SIGFP
eof-cad

# Free-R flag
uniqueify -p 0.05 ${filename}_cad.mtz ${filename}_free.mtz > ${filename}_freer.log

# Model building with ARPwARP
MODEL_BUILDING:
set filename='iom'
##
## Auto model building
##
$warpbin/auto_warp.sh \
datafile `pwd`/${filename}_free.mtz \
nresidues $NRESIDUES \
workdir `pwd` \
seqin $SEQ_FILE \
fp FP \
sigfp SIGFP \
phibest PHIB \
fom FOM \
cycles 50  > arp_start.log

END:
