#!/usr/local/apps/perl5.004/bin/perl
##!D:\Perl5\bin\perl.exe # for use on PC systems

##################################################################
#
#  match.pl - a script that uses a genetic sequence database file 
#       to create a file containing a list of positives and
#	negatives.
#
#  Author:  David Martinez
#       San Diego Supercomputer Center
#            (damartin@sdsc.edu)        August 1997
#
# The file contains:
# original database file name 
# Number of true positives
# Number of true negatives
#
##################################################################

########################## GLOBALS ###############################
# list of possible input files:
#
#   bmotra2x.seq   gbphg.seq@     gbvrt.seq@     lepqbrep.seq   ysctrf.seq
#   gbmam.seq@     gbrna.seq@     qbeta.seq
#
# name of the input file
#$infile = "../seq/bmotra2x.seq";
#$infile = "../seq/gbphg.seq";
$infile = "../seq/gbrna.seq";
#$infile = "../seq/gbvrt.seq";
#$infile = "../seq/lepqbrep.seq";
#$infile = "../seq/qbeta.seq";
#$infile = "../seq/ysctrf.seq";
##################################################################

# for error checking
if ( !(-e $infile) ) {
   print "sequence database file $infile doesn't exist!\n";
   print "Bailing out..\n";
   exit(1);
}
print "using sequence database file: $infile\n";

open(INFILE, "<$infile") || die "Unable to open filename \"$!\".\n";

# need query filename
#
#
#$queryfile="demo.qry";
$queryfile="RNaseP.qry";
#$queryfile="Xlaevis5srna.qry";
print "Using $queryfile queryfile\n";
##################################################################

@temp1= split('\/',$infile);
@temp1= split('\.',$temp1[2]);
@temp2= split('\.',$queryfile);

# name of the match output file
#$outfilename = "match/gbrna_Xlaevis5srna.mat"; # database file _ qry file
#$outfilename = "match/gbrna_demo.mat";
#$outfilename = "match/gbrna_RNaseP.mat";

$outfilename = "match/".$temp1[0]."_".$temp2[0].".mat";
print "Creating $outfilename matchfile\n";

$posval="1";
$negval=0;

$numlocus=0;
$n=0;
$m=0;

# target string information
#-----------------------------------------------------------
#	Used for Xlaevis5sXX.qry files
if ($queryfile eq "Xlaevis5srna.qry") {
   $target0[0]="5S";
   $target0[1]="ribosomal";
   $target0[2]="RNA";
   $tar_str0=$target0[0].$target0[1].$target0[2];

   $target1[0]="5S";
   $target1[1]="ribosomal";
   $target1[2]="RNA";
   $tar_str1=$tar_str0;
}
#
#-----------------------------------------------------------
#	Used for demo.qry file
if ($queryfile eq "demo.qry") {
   $target0[0]="trna";
   $tar_str0=$target0[0];

   $target1[0]="transfer";
   $target1[1]="rna";
   $tar_str1=$target1[0].$target1[1];
}
#
#-----------------------------------------------------------
#	Used for RNaseP.qry file
if ($queryfile eq "RNaseP.qry") {
   $target0[0]="rnase";
   $target0[1]="p";
   $tar_str0=$target0[0].$target0[1];

   $target1[0]="RNase";
   $target1[1]="P";
   $tar_str1=$target1[0].$target1[1];
}
#-----------------------------------------------------------


$keyword="";
$new=0;

########################### MAIN #################################

#-------------------------------------------------------------
# start creating output file
#-------------------------------------------------------------
if ( !open(OUTFILE,"> $outfilename") ) {
   print "Error opening an output file $outfilename: $!\n";
   print "Bailing out..\n";
   exit(1);
}


## Start with sequence database file
  while ($line = <INFILE>) {
    $line =~ m/^\s*(\S+)/s;
    if ($1 eq "LOCUS") {
	$numlocus++;
	if ($new) {
	    if ((lc $keyword eq lc $tar_str0)|(lc $keyword eq lc $tar_str1)) {
		$result=$posval;
		$n++;  
	    }
	    else {
		$result=$negval;
		$m++;
	    }

	    $correctval{$key}=$result;
	    print OUTFILE "key  $key  result  $result\n";
	} ## end if ($new)

	$new=1;
	$kw=0; # reset keywords flag
	$keyword="";
        ##  We're at the beginning of a sequence
        ##  Parse sequence name
	@fields = split(' ',$line);

	## set key to LOCUS name
	$key = $fields[1];  ## set key to LOCUS name

    } ## end if ($1 eq LOCUS)

    elsif ( $kw==0 ) {
        ##  Parse sequence keywords
	@fields = split(' ',$line);
	# Get number of elements
	$len_fields=$#fields+1;
	## reconstruct keyword string
	## cycle thru all values of $i based upon number 
	## of elements in the $fields array

	for ($i=1;$i<$len_fields && $kw==0;$i++) {
	   $keyword="";

	   ## remove non-alpha/num from string
	   $fields[$i]=~s/[^\w]//g;
## check for all target strings within sequence database

	   $index=0;
	   $junk=$fields[$i];
	   $start=$i;
	   $seq=$start;

	   while ((lc $junk eq lc $target0[$index]) && $seq<$len_fields) { 
		$index++;
		$keyword.=$junk;
		$seq++;
		$fields[$seq]=~s/[^\w]//g;
		$junk=$fields[$seq];
	   } # end while

           ##
           ## set flag once we have found entire target string
           ##
           if (lc $keyword eq lc $tar_str0) {
                ## keyword is found
                ## Go to next Locus
                $kw=1;
                last;
           } # end if

	   $seq=$start;
           $junk=$fields[$i];
	   $index=0;

           while ((lc $junk eq lc $target1[$index]) && $seq<$len_fields) {
                $index++;
                $keyword.=$junk;
                $seq++;
                $fields[$seq]=~s/[^\w]//g;
                $junk=$fields[$seq];
           } # end while

	   ##
	   ## set flag once we have found entire target string
	   ##
	   if (lc $keyword eq lc $tar_str1) {
		## keyword is found
		## Go to next Locus
		$kw=1;
		last;
	   } # end if

	} # end for
  } # end elsif
  } # end while

##
##  Do for last LOCUS found
    if ((lc $keyword eq lc $tar_str0)|(lc $keyword eq lc $tar_str1)) {
	$result=$posval;
	$n++;  
    }
    else {
	$result=$negval;
	$m++;
    }
    $correctval{$key}=$result;
    print OUTFILE "key  $key  result  $result\n\n\n";

print OUTFILE "LOCUS $numlocus values within $infile\n";
print OUTFILE "TRUE POSITIVES $n sequences with $tar_str0, $tar_str1\n";
print OUTFILE "TRUE NEGATIVES $m sequences with $tar_str0, $tar_str1\n";

close OUTFILE;

print "DONE\n";