#!/usr/bin/perl 
#-w
##############################################################################
#
# Author:Jack Min, on Jan. 8, 2003
############################################ Body of the code 

if ($#ARGV != 2)
{
	print "Usage: seqRetrieve seqFile seqID_file retrieved_seq\n\n";
	exit;
}

print "In progress ...\n";

  $seqFile = '';
  $IDfile  = '';
  $retrievedSeqFile = '';

$seqFile = $ARGV[0];
$IDfile = $ARGV[1];
$retrievedSeqFile = $ARGV[2];
unless( open(IN, $IDfile) ) {
    print STDERR "Cannot open file IDfile\n\n";
    exit;
    }
  @idList = <IN>;
  close IN;

my %idHash = (); # Initalize hash
foreach my $item(@idList)
{
	chomp ($item); # Our input list had an extra \n so couldn't directly input to hash
	@tmp = split(/\t/,$item);
	$id = $tmp[0];
	#$id = $item; 
	$idHash{$id} = 1; # Set key to 1, we need search only values
}

# Step one get new gene
  @seq = ();
# Step two: get sequence data
   unless(open(IN, $seqFile) ) {
   print STDERR "Cannot open file singletsfile\n\n";
    exit;
    }
  #@seq = <IN>;
  #close IN;
  %seqH = (); #seq harsh
  $i = 0; #as a count
  $sequence = '';
  $id = '';
#foreach $line (@seq)
#open (OUT, ">>$retrievedSeqFile");

while (<IN>)
{
	$line = $_;  
	chomp ($line);
	$line =~ s/\*//;
	# discard blank line
    if ($line =~ /^\s*$/)
    {
       	next;                       # discard comment line
    }
    elsif($line =~ /^>/)
    {
          $i++; print $i."\n";
	  if ($i > 1)
	  {
		if (exists $idHash{$id})
		{
			print $i."\n";
			open (OUT, ">>$retrievedSeqFile");
		    print OUT (">".$id."\n".$sequence."\n");
			close OUT;
		}			
	  }		
	@tmp = split(/\|/, $line); #split with one or more space, to get the seqID
	#$id = substr($tmp[0], 1);
	$id = $tmp[1];
	#	$id = substr ($line, 1);
	  	$sequence = '';
    }
       else
       {
         $sequence .= $line;
       }
}
#$seqH{$id} = $sequence;
if (exists $idHash{$id})
{
	open (OUT, ">>$retrievedSeqFile");
	print OUT (">".$id."\n".$sequence."\n");
	close OUT;
}	
close (IN);
#close OUT;

#open (OUT, ">$retrievedSeqFile");
#print OUT (@retrievedSeq);
#close OUT;

print "The job is done!!\n\n";
exit;

