t/primaryqual.t

   1 # -*-Perl-*-
   2 ## Bioperl Test Harness Script for Modules
   3 ## $Id$
   4 #
   5 # modeled after the t/Allele.t test script
   6
   7 use strict;
   8 use vars qw($DEBUG);
   9 $DEBUG = $ENV{'BIOPERLDEBUG'};
  10 my $verbose = -1 unless $DEBUG;
  11 BEGIN {
  12     # to handle systems with no installed Test module
  13     # we include the t dir (where a copy of Test.pm is located)
  14     # as a fallback
  15     eval { require Test; };
  16     if( $@ ) {
  17         use lib 't';
  18     }
  19     use Test;
  20     plan tests => 31;
  21 }
  22
  23 END {
  24     unlink qw(batch_write_qual.qual write_qual.qual);
  25
  26 }
  27 # redirect STDERR to STDOUT
  28 open (STDERR, ">&STDOUT");
  29 use Bio::Root::IO;
  30 use Bio::SeqIO;
  31 use Bio::Seq::SeqWithQuality;
  32 use Bio::Seq::PrimaryQual;
  33
  34 my $string_quals = "10 20 30 40 50 40 30 20 10";
  35 print("Quals are $string_quals\n") if($DEBUG);
  36 my $qualobj = Bio::Seq::PrimaryQual->new( '-qual' => $string_quals,
  37                                           '-id'  => 'QualityFragment-12',
  38                                           '-accession_number' => 'X78121',
  39                                           );
  40 ok($qualobj);
  41 ok($qualobj->display_id, 'QualityFragment-12');
  42 ok($qualobj->accession_number, 'X78121');
  43
  44 my @q2 = split/ /,$string_quals;
  45 $qualobj = Bio::Seq::PrimaryQual->new
  46     ( '-qual'             => \@q2,
  47       '-primary_id'       =>    'chads primary_id',
  48       '-desc'             =>    'chads desc',
  49       '-accession_number' => 'chads accession_number',
  50       '-id'               =>    'chads id'
  51       );
  52
  53 ok($qualobj->primary_id, 'chads primary_id');
  54 my $rqual = $qualobj->qual();
  55 ok(ref($rqual) eq "ARRAY");
  56
  57 my $newqualstring = "50 90 1000 20 12 0 0";
  58
  59 $qualobj->qual($newqualstring);
  60 my $retrieved_quality = $qualobj->qual();
  61 my $retrieved_quality_string = join(' ', @$retrieved_quality);
  62 ok($retrieved_quality_string,$newqualstring);
  63
  64 my @newqualarray = split/ /,$newqualstring;
  65 $qualobj->qual(\@newqualarray);
  66 $retrieved_quality = $qualobj->qual();
  67 $retrieved_quality_string = join(' ',@$retrieved_quality);
  68 ok($retrieved_quality_string,$newqualstring);
  69
  70 eval {
  71     $qualobj->qual("chad");
  72 };
  73 ok($@ =~ /not look healthy/);
  74
  75 eval { $qualobj->qual(""); };
  76 ok(!$@);
  77
  78 eval { $qualobj->qual(" 4"); };
  79 ok(!$@);
  80
  81 ok($qualobj->length(),2 );
  82 $qualobj->qual("10 20 30 40 50 40 30 20 10");
  83 my @subquals = @{$qualobj->subqual(3,6);};
  84 ok(@subquals, 4);
  85      # chad, note to self, evaluate border conditions
  86 ok ("30 20 10" eq join(' ',@{$qualobj->subqual(7,9)}));
  87
  88
  89
  90 my @false_comparator = qw(30 40 70 40);
  91 my @true_comparator = qw(30 40 50 40);
  92 ok(!&compare_arrays(\@subquals,\@true_comparator));
  93
  94 eval { $qualobj->subqual(-1,6); };
  95 ok($@ =~ /EX/ );
  96 eval { $qualobj->subqual(1,6); };
  97 ok(!$@);
  98 eval { $qualobj->subqual(1,9); };
  99 ok(!$@);
 100 eval { $qualobj->subqual(9,1); };
 101 ok($@ =~ /EX/ );
 102
 103
 104 ok($qualobj->display_id() eq "chads id");
 105 $qualobj->display_id("chads new display_id");
 106 ok($qualobj->display_id() eq "chads new display_id");
 107
 108 ok($qualobj->accession_number(), "chads accession_number");
 109 $qualobj->accession_number("chads new accession_number");
 110 ok($qualobj->accession_number(), "chads new accession_number");
 111 ok($qualobj->primary_id(), "chads primary_id");
 112 $qualobj->primary_id("chads new primary_id");
 113 ok($qualobj->primary_id(), "chads new primary_id");
 114
 115 ok($qualobj->desc(), "chads desc");
 116 $qualobj->desc("chads new desc");
 117 ok($qualobj->desc(), "chads new desc");
 118 ok($qualobj->display_id(), "chads new display_id");
 119 $qualobj->display_id("chads new id");
 120 ok($qualobj->display_id(), "chads new id");
 121
 122 my $in_qual  = Bio::SeqIO->new(-file => "<" . Bio::Root::IO->catfile("t","data","qualfile.qual") ,
 123                                '-format' => 'qual',
 124                                '-verbose' => $verbose);
 125 ok($in_qual);
 126 my $pq = $in_qual->next_seq();
 127 ok($pq->qual()->[99], '39'); # spot check boundary
 128 ok($pq->qual()->[100], '39'); # spot check boundary
 129
 130 my $out_qual = Bio::SeqIO->new('-file'    => ">write_qual.qual",
 131                                '-format'  => 'qual',
 132                                '-verbose' => $verbose);
 133 $out_qual->write_seq(-source    =>      $pq);
 134
 135 my $swq545 = Bio::Seq::SeqWithQuality->new (    -seq    =>      "ATA",
 136                                                 -qual   =>      $pq
 137                                         );
 138 $out_qual->write_seq(-source    =>      $swq545);
 139
 140
 141
 142 $in_qual = Bio::SeqIO->new('-file' => Bio::Root::IO->catfile("t","data","qualfile.qual") ,
 143                            '-format' => 'qual',
 144                            '-verbose' => $verbose);
 145
 146 my $out_qual2 = Bio::SeqIO->new('-file'    => ">batch_write_qual.qual",
 147                                 '-format'  => 'qual',
 148                                 '-verbose' => $verbose);
 149
 150 while ( my $batch_qual = $in_qual->next_seq() ) {
 151         $out_qual2->write_seq(-source   =>      $batch_qual);
 152 }
 153
 154 sub display {
 155     if($DEBUG ) {
 156         my @quals;
 157         print("I saw these in qualfile.qual:\n") ;
 158         while ( my $qual = $in_qual->next_seq() ) {
 159             # ::dumpValue($qual);
 160             print($qual->display_id()."\n");
 161             @quals = @{$qual->qual()};
 162             print("(".scalar(@quals).") quality values.\n");
 163         }
 164     }
 165 }
 166
 167 # dumpValue($qualobj);
 168
 169 sub compare_arrays {
 170     my ($a1,$a2) = @_;
 171     return 1 if (scalar(@{$a1}) != scalar(@{$a2}));
 172     my ($v1,$v2,$diff,$curr);
 173     for ($curr=0;$curr<scalar(@{$a1});$curr++){
 174         return 1 if ($a1->[$curr] ne $a2->[$curr]);
 175     }
 176     return 0;
 177 }