#!usr/bin/perl
######################################################################
#
#Ian Odell
#Read Hep C fasta
#Remove whole sequence if it has ambiguous base(s) because
#'bad' sequencing.
#
######################################################################
use strict;
use warnings;
use Bio::SeqIO;
use Bio::Seq;

my $inputfilename = "type1sorted.fasta";
my $outputfilename = ">type1_ambigsRemoved.fasta";

my $in = Bio::SeqIO->new(-file => $inputfilename);

#read fasta and put all sequences into an array.
my @seq_array;
while (my $seq = $in->next_seq()) {
  push @seq_array, $seq;
}

#for each sequence, check if there is an ambiguous base
#if reaches end of seq, add to good array.
my @seqsNOambigs;
my @ambigs = ('U', 'R', 'Y', 'K', 'M', 'S', 'W', 'B', 'D', 'H', 'V', 'N');
my $ambigsize = @ambigs;
my $is_ambig;
my $count;
my $sequence;
my $ambigs = "";
print "@ambigs\n";

foreach my $seq_object (@seq_array) {
  $sequence = $seq_object->seq;
  $is_ambig = 0;
  $count = 0;
  do {
    if ($sequence =~ $ambigs[$count]) { #if the seq does not contain this ambig base
      $is_ambig = 1;		#yes is ambiguous
    }
    $count++;
  }
    until ($is_ambig == 1 or $count == $ambigsize);
  if ($is_ambig == 0) {		#save sequence if it doesn't have any amb bases
    push @seqsNOambigs, $seq_object;
  } else {
    $ambigs .= $seq_object->id . "\n";
  }
}

#write new alignments to their own files
my $out = Bio::SeqIO->new(-file => $outputfilename);
foreach my $seq (@seqsNOambigs) {
  $out->write_seq($seq);
}

print "$ambigs";
