#!/usr/bin/perl
# 
# This perl script is responsible for finding a certain
# type of connectivity pattern in your xpk files. It will call the
# peakcon program to do this, the reason for this script is to
# allow for additional filtering of patterns which peakcon is not
# capable of doing.
#

#
# Specify the type of connectivity pattern
# Pattern type is residue-residue
#
$TYPE = 'dxx';
$PROG_NAME = 'peakcon';

#
# arguments to this script
#
$BIN_DIR = $ARGV[0];
$INPUT_DIR = $ARGV[1];
$OUTPUT_DIR = $ARGV[2];
$OUTPUT_PREF = $ARGV[3];
$REF = $ARGV[4];
$HNCACB = $ARGV[5];
$CBCACO = $ARGV[6];

#
# we will use the peakcon program with the following input/output 
# RULES tells me what the pattern looks like
#
$PROG = $BIN_DIR.'/'.$PROG_NAME;
$TMP1 = $OUTPUT_DIR.'/'.'tmp.1';
$TMP2 = $OUTPUT_DIR.'/'.'tmp.2';
$TMP3 = $OUTPUT_DIR.'/'.'tmp.3';
$TMP4 = $OUTPUT_DIR.'/'.'tmp.4';
$TMPX = $OUTPUT_DIR.'/'.'tmp.x';
$RULES = $INPUT_DIR.'/'.'rules.'.$TYPE;
$RESULTS = $OUTPUT_DIR.'/'.$OUTPUT_PREF.'.'.$TYPE;
system("rm -f $TMP1 $TMP2 $TMP3 $TMP4 $TMPX");

#
# call peakcon
#
system("$PROG $RULES.00 $REF $HNCACB $CBCACO $TMP1");
system("$PROG $RULES.01 $REF $HNCACB $CBCACO $TMP2");
system("$PROG $RULES.10 $REF $HNCACB $CBCACO $TMP3");
system("$PROG $RULES.11 $REF $HNCACB $CBCACO $TMP4");
system("touch $TMP1 $TMP2 $TMP3 $TMP4");

#
# sort the output (or do whatever else we want to the output) 
# sorting on these chosen fields will most likely place the correct 
# pattern first (when more that one choice available)
#
$C_REF0 = 1;
$C_REF1 = 2;
$C_PROB = 4;
$C_FILTER = 14;
$C_SCA1 = 15;
$C_SCB1 = 16;

#
# Sorting based on fit as opposed to intensity.
#system("cat $TMP1 $TMP2 $TMP3 $TMP4 | sort -r -n -k $C_REF0 -k $C_REF1 -k $C_FILTER -k $C_PROB >$TMPX");

#
# Sorting based on intensity as opposed to frequency
system("cat $TMP1 $TMP2 $TMP3 $TMP4 | sort -r -n -k $C_REF0 -k $C_REF1 -k $C_FILTER -k $C_SCA1 -k $C_SCB1 >$TMPX");


#
# print the first peak of all ref1/ref2 pairs.
#
open(INPUT, "$TMPX");
open(OUTPUT, ">$RESULTS");
printf OUTPUT "#\n";
printf OUTPUT "#i-1   i  Rule Prob ------- Peaklist Ids -------";
printf OUTPUT " Type  Ca1-Int Ca2-Int";
printf OUTPUT "   Ca0   Cb0     N0    H0   Ca1   Cb1     N1    H1\n";
printf OUTPUT "#\n";

$h1 = -1;
$h2 = -1;
while () {
	chomp; @fields = split;
	if ((@fields[0] != $h1) || (@fields[1] != $h2)) {
		$h1 = @fields[0];
		$h2 = @fields[1];
		$rtype = @fields[13];
	}

	if ((@fields[0] == $h1) && (@fields[1] == $h2) && (@fields[13] == $rtype)) {
		printf OUTPUT " %3d %3d", @fields[0], @fields[1];
		printf OUTPUT "  %4s %4.2f", @fields[2], @fields[3];
		printf OUTPUT "  %3d %3d %3d %3d %3d %3d %3d %3d %3d", @fields[4], @fields[5], @fields[6], @fields[7], @fields[8], @fields[9], @fields[10], @fields[11], @fields[12];
		printf OUTPUT " %4d  %6.4g %6.4g", @fields[13], @fields[14], @fields[15];
		printf OUTPUT " %5.2f %5.2f %6.2f %5.2f %5.2f %5.2f %6.2f %5.2f", @fields[16], @fields[17], @fields[18], @fields[19], @fields[20], @fields[21], @fields[22], @fields[23];
		printf OUTPUT "\n";
	}
}

#
# clean up
#
close(INPUT); close(OUTPUT);
system("rm -f $TMP1 $TMP2 $TMP3 $TMP4 $TMPX");