#!/bin/tcsh -f
#
#
# Sort, scale and merge the data from a mosflm mtz file
#
#
if (${0:h} == ${0}) then
set cmd = ""
else
set cmd = "${0:h}/"
endif
#
#
source ${cmd}makesymbols
#
#
source symbols
#
#
# Define the usage specific to integ
#
#
set input = '[[mtzfile] mtzfile ...]'
#
#
set help =  (title            \
             outputroot       \
             highres          \
             lowres           \
             merge            \
             addscale         \
             nresidue         \
             anomalous        \
             exclude          \
             sdadd)
#
#
source ${cmd}cmdargs
#
#
# Interpret the unmatched arguments
#
#
if (${#argv} >= 1) then 
set merge = (${argv})
endif
#
#
# Save the symbols
#
#
source ${cmd}savesymbols
#
#
# Do any symbols have defaults?
#
#
source ${cmd}defaults
#
#
# Interrupt action
#
#
onintr interrupt
#
#
# Check for required input
#
#
if (${?merge} == 0) then

echo " "
echo " ERROR: mtz files to merge not specified, program exit. "
echo " "
echo " Usage:" ${0} " mtzfile [mtzfile ...] "
echo " "

kill -9 $$

endif
#
#
# Make the batch numbers unique
#
#
set num = 1
while (${num} <= ${#merge})
if (${num} == 1) then
if (${#merge} == 1) then
if (${merge[1]} !~ [/\\]*) then
set filewd = ${cwd}/
else
set filewd = 
endif
if (-l ${merge[1]}) then
cp ${filewd}${merge[1]} ${filewd}${merge[1]}.unique
else
ln -sf ${filewd}${merge[1]} ${filewd}${merge[1]}.unique
endif
else
cp ${merge[${num}]} ${merge[${num}]}.unique
endif
else
@ add = (${num} - 1) * 1000
rebatch HKLIN ${merge[${num}]} HKLOUT ${merge[${num}]}.unique << eof
BATCH ADD ${add}
eof
endif
@ num = ${num} + 1
end
#
#
# Make the sortmtz command file
#
#
echo "H K L M/ISYM BATCH I SIGI" > ${tmpfile}.sortmtz
set num = 1
while (${num} <= ${#merge})
echo "${merge[${num}]}.unique" >> ${tmpfile}.sortmtz
@ num = ${num} + 1
end
echo "" >> ${tmpfile}.sortmtz
#
#
# Make the scala command file
#
#
set noscales      = 0
set nointensities = 0
set norun         = 0
set noreject      = 0
if (${?addscale}) then
set num = 1
while (${num} <= ${#addscale})
if (${noscales} == 0) then
@ noscales = ("${addscale[${num}]}" =~ [Ss][Cc][Aa][Ll][Ee]*)
endif
if (${nointensities} == 0) then
@ nointensities = ("${addscale[${num}]}" =~ [Ii][Nn][Tt][Ee][Nn]*)
endif
if (${norun} == 0) then
@ norun    = ("${addscale[${num}]}" =~ [Rr][Uu][Nn]*)
endif
if (${noreject} == 0) then
@ noreject = ("${addscale[${num}]}" =~ [Rr][Ee][Jj][Ee]*)
endif
@ num = ${num} + 1
end
endif
#
#
echo "" > ${tmpfile}.scala
if (${norun} == 0) then
set num = 1
while (${num} <= ${#merge})
@ addstart =  (${num} - 1) * 1000
@ addend   = ((${num} - 1) * 1000) + 999
echo -n "run ${num} include batch ${addstart} to ${addend}" >> ${tmpfile}.scala
set ex = 1
set excluded = 0
if (${?exclude}) then
while (${ex} <= ${#exclude})
if ((${addstart} <= ${exclude[${ex}]}) && (${addend} >= ${exclude[${ex}]})) then
if (${excluded} == 0) then
echo -n " exclude" >> ${tmpfile}.scala
set excluded = 1
endif
echo -n " ${exclude[${ex}]}" >> ${tmpfile}.scala
endif
@ ex = ${ex} + 1
end
endif
echo "" >> ${tmpfile}.scala
@ num = ${num} + 1
end
endif
#
#
if (${noscales} == 0) then
echo "scales brotation spacing 20 batch" >> ${tmpfile}.scala
endif
#
#
if (${nointensities} == 0) then
echo "intensities partials" >> ${tmpfile}.scala
endif
#
#
if (${?addscale}) then
set num = 1
while (${num} <= ${#addscale})
echo "${addscale[${num}]}" >> ${tmpfile}.scala
@ num = ${num} + 1
end
endif
#
#
if (${?highres}) then
echo "resolution   HIGH ${highres}" >> ${tmpfile}.scala
endif
if (${?lowres}) then
echo "resolution   LOW ${lowres}" >> ${tmpfile}.scala
endif
#
#
echo "cycles 	   40" >> ${tmpfile}.scala
#
#
if (${?sdadd}) then
echo "sdcorrection 1.5 ${sdadd}" >> ${tmpfile}.scala
else
echo "sdcorrection 1.5 0.02" >> ${tmpfile}.scala
endif
#
#
if (${noreject} == 0) then
echo "reject       merge 4" >> ${tmpfile}.scala
endif
#
#
if (${?anomalous}) then
echo "anomalous    ${anomalous}" >> ${tmpfile}.scala
else
echo "anomalous    off" >> ${tmpfile}.scala
endif
#
#
echo "bins         20" >> ${tmpfile}.scala
#echo "name project jcsg crystal tm0000 dataset wave"  >> ${tmpfile}.scala
#
# make a second scala file to run "OUTPUT POLISH UNMERGED ORIGINAL", for jcsg usage
cp ${tmpfile}.scala ${tmpfile}.scala_outHKL
echo "OUTPUT POLISH UNMERGED ORIGINAL" >> ${tmpfile}.scala_outHKL
#
# Make the truncate command file
#
#
echo "" > ${tmpfile}.truncate
#
#
if (${?title}) then
echo "title    ${title}" >> ${tmpfile}.truncate
endif
#
#
if (${?nresidue}) then
echo "nresidue ${nresidue}" >> ${tmpfile}.truncate
endif
#
#
# Sort the data
#
#
sortmtz HKLOUT ${tmpfile}.scratch1 < ${tmpfile}.sortmtz
#
#
# Delete the .unique files
#
#
set num = 1
while (${num} <= ${#merge})
rm -f "${merge[${num}]}.unique"
@ num = ${num} + 1
end
#
#
# Get the current version number
#
#
set logfiles = `ls ${outputroot}.scale_log*`
set version = 0
set num = 1
while (${num} <= ${#logfiles})
echo -n "set tmp = " > ${tmpfile}
echo "(${logfiles[${num}]:e})" | sed "s/[a-z_]/ /g" >> ${tmpfile}
source ${tmpfile}
rm -f ${tmpfile}
if (${tmp} > ${version}) then
@ version = ${tmp}
endif
@ num = ${num} + 1
end
@ version = ${version} + 1
#
#
# Create the log file
#
#
echo "" > ${outputroot}.scale_log${version}
#
#
# Update the link to the log file
#
#
if (-e ${outputroot}.scale_log) then
if (-l ${outputroot}.scale_log) then
rm -f ${outputroot}.scale_log
ln -s ${outputroot}.scale_log${version} ${outputroot}.scale_log
endif
else
ln -s ${outputroot}.scale_log${version} ${outputroot}.scale_log
endif

echo "" > ${outputroot}.truncate_log${version}
if (-e ${outputroot}.truncate_log) then
if (-l ${outputroot}.truncate_log) then
rm -f ${outputroot}.truncate_log
ln -s ${outputroot}.truncate_log${version} ${outputroot}.truncate_log
endif
else
ln -s ${outputroot}.truncate_log${version} ${outputroot}.truncate_log
endif

#
#
# Follow the progress of the scaling
#
#
tail -f ${outputroot}.scale_log${version} | grep -E "Mean and maximum shift" &
#
#
# Calculate and apply the scale factors
#
#
cp ${tmpfile}.scala scala.inp
scala HKLIN  ${tmpfile}.scratch1 \
      HKLOUT ${tmpfile}.scratch2 >> ${outputroot}.scale_log${version} < ${tmpfile}.scala
#
#
if (-e ANOMPLOT) then
if (${?anomalous}) then
if (${anomalous} !~ [Oo][Ff][Ff]) then
mv ANOMPLOT ${outputroot}.anomplot
endif
endif
endif
if (-e NORMPLOT) then
mv NORMPLOT ${outputroot}.normplot
endif
if (-e ROGUES) then
mv ROGUES ${outputroot}.rogues
endif
if (-e SCALES) then
mv SCALES ${outputroot}.scales
endif
#
# run it again to generate HKL formatted file, for jcsg
rm -rf SCALEPACK
scala HKLIN  ${tmpfile}.scratch1 \
      HKLOUT ${outputroot}.scala_umerged.mtz >/dev/null  < ${tmpfile}.scala_outHKL
#cleaning up
rm -rf ANOMPLOT NORMPLOT ROGUES SCALES
mv SCALEPACK ${outputroot}.scala_umerged.hkl
#
# Create summary file (only on an alpha)
#
#
if (${HOSTTYPE} == "alpha") then

echo "" > ${outputroot}.scale_sum
echo "  Scale factors analysed by Batch for each run" >> ${outputroot}.scale_sum
echo "  ============================================" >> ${outputroot}.scale_sum
echo "" >> ${outputroot}.scale_sum
grep -E -p"========" \
               "Scales v rotation range" ${outputroot}.scale_log${version} \
                                        >> ${outputroot}.scale_sum
#
#
echo "" >> ${outputroot}.scale_sum
grep -E -p"[I=][N=][T=][E=][N=][S=][I=]" \
               "Rmerge against resolution" ${outputroot}.scale_log${version} \
                                        >> ${outputroot}.scale_sum
#
#
echo "  Completeness and multiplicity, including reflections measured only once" \
     >> ${outputroot}.scale_sum
echo "  =======================================================================" \
     >> ${outputroot}.scale_sum
echo "" >> ${outputroot}.scale_sum
grep -E -p"=======" "TABLE: +Completeness" ${outputroot}.scale_log${version} \
                                        >> ${outputroot}.scale_sum

else

if (-e ${outputroot}.scale_sum) then

rm -f ${outputroot}.scale_sum

endif
endif
#
#
# Convert the averaged intensities into mean amplitudes
#
#

# added by QP to save scala mtz file and to separate scale/truncate logfile
cp ${tmpfile}.scratch2  ${outputroot}_scala.mtz

#
truncate HKLIN  ${tmpfile}.scratch2 \
         HKLOUT ${outputroot}_F.mtz > ${outputroot}.truncate_log${version} < ${tmpfile}.truncate
#
#
if (${HOSTTYPE} == "alpha") then

echo "  Wilson Plot" >> ${outputroot}.scale_sum
echo "  ===========" >> ${outputroot}.scale_sum
echo "" >> ${outputroot}.scale_sum
grep -E -p"[*i][*n][*c][*l][*u][*d][*e]" \
               "TABLE: Wilson Plot" ${outputroot}.scale_log${version} \
                                        >> ${outputroot}.scale_sum

endif
#
#
interrupt:
#
kill %+
#
#
# Delete all the scratch files
#
#
rm -f ${tmpfile}.sortmtz
rm -f ${tmpfile}.scala
rm -f ${tmpfile}.truncate
rm -f ${tmpfile}.scratch[12]