]> git.donarmstrong.com Git - biopieces.git/blobdiff - code_perl/Maasha/Align.pm
adding bzip2 support in ruby
[biopieces.git] / code_perl / Maasha / Align.pm
index 986a9b7a9bd1246619cd14d547ee754860685ea3..262ee48255634929b0e7383224952b3800fc88aa 100644 (file)
@@ -28,6 +28,7 @@ package Maasha::Align;
 # >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>><<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
 
 
+use warnings;
 use strict;
 use Data::Dumper;
 use IPC::Open2;
@@ -38,8 +39,8 @@ use Maasha::Seq;
 use vars qw ( @ISA @EXPORT );
 
 use constant {
-    HEAD  => 0,
-    SEQ   => 1,
+    SEQ_NAME => 0,
+    SEQ      => 1,
 };
 
 @ISA = qw( Exporter );
@@ -157,7 +158,7 @@ sub align_print_multi
     $consensus = align_consensus( $entries ) if not $no_cons;
 
     unshift @{ $entries }, $ruler if not $no_ruler;
-    push    @{ $entries }, $consensus;
+    push    @{ $entries }, $consensus if not $no_cons;
 
     align_print( $entries, $fh, $wrap );
 }
@@ -180,7 +181,7 @@ sub align_print
 
     $max = 0;
 
-    map { $max = length $_->[ HEAD ] if length $_->[ HEAD ] > $max } @{ $entries };
+    map { $max = length $_->[ SEQ_NAME ] if length $_->[ SEQ_NAME ] > $max } @{ $entries };
 
     $blocks = align_wrap( $entries, $wrap );
 
@@ -188,12 +189,12 @@ sub align_print
     {
         foreach $entry ( @{ $block } )
         {
-            $entry->[ HEAD ] =~ s/stats|ruler|consensus//;                                                                          
+            $entry->[ SEQ_NAME ] =~ s/stats|ruler|consensus//;                                                                          
 
             if ( $fh ) {
-                print $fh $entry->[ HEAD ], " " x ( $max + 3 - length $entry->[ HEAD ] ), $entry->[ SEQ ], "\n";                                 
+                print $fh $entry->[ SEQ_NAME ], " " x ( $max + 3 - length $entry->[ SEQ_NAME ] ), $entry->[ SEQ ], "\n";                                 
             } else {
-                print $entry->[ HEAD ], " " x ( $max + 3 - length $entry->[ HEAD ] ), $entry->[ SEQ ], "\n";                                 
+                print $entry->[ SEQ_NAME ], " " x ( $max + 3 - length $entry->[ SEQ_NAME ] ), $entry->[ SEQ ], "\n";                                 
             }
         }
     }
@@ -225,7 +226,7 @@ sub align_wrap
 
         for ( $c = 0; $c < @{ $entries }; $c++ )
         {
-            if ( $entries->[ $c ]->[ HEAD ] eq "ruler" )
+            if ( $entries->[ $c ]->[ SEQ_NAME ] eq "ruler" )
             {
                 $ruler = substr $entries->[ $c ]->[ SEQ ], $i, $wrap;
 
@@ -241,7 +242,7 @@ sub align_wrap
             }
             else
             {
-                push @lines, [ $entries->[ $c ]->[ HEAD ], substr $entries->[ $c ]->[ SEQ ], $i, $wrap ];
+                push @lines, [ $entries->[ $c ]->[ SEQ_NAME ], substr $entries->[ $c ]->[ SEQ ], $i, $wrap ];
             }
         }
 
@@ -271,18 +272,18 @@ sub align_pins
 
     $type ||= Maasha::Seq::seq_guess_type( $entry1->[ SEQ ] );
 
-    $blosum = blosum_read() if $type =~ /protein/;
+    $blosum = blosum_read() if $type =~ "PROTEIN";
 
     for ( $i = 0; $i < length $entry1->[ SEQ ]; $i++ )
     {
-        $char1 = substr $entry1->[ SEQ ], $i, 1;
-        $char2 = substr $entry2->[ SEQ ], $i, 1;
+        $char1 = uc substr $entry1->[ SEQ ], $i, 1;
+        $char2 = uc substr $entry2->[ SEQ ], $i, 1;
 
         if ( $blosum and $char1 eq $char2 ) {
             $pins .= $char1;
         } elsif ( $char1 eq $char2 ) {
             $pins .= "|";
-        } elsif ( $blosum and $blosum->{ $char1 }->{ $char2 } > 0 ) {
+        } elsif ( $blosum and defined $blosum->{ $char1 }->{ $char2 } and $blosum->{ $char1 }->{ $char2 } > 0 ) {
             $pins .= "+";
         } else {
             $pins .= " ";
@@ -402,7 +403,7 @@ sub align_consensus
         }
     }
 
-    $entry->[ HEAD ] = "Consensus: $min_sim%";
+    $entry->[ SEQ_NAME ] = "Consensus: $min_sim%";
 
     return wantarray ? @{ $entry } : $entry;
 }
@@ -462,8 +463,8 @@ sub align_tile
     # using pairwise alignments. The result is returned as a list of
     # aligned FASTA entries.
     
-    my ( $ref_entry,   # reference entry as [ HEAD, SEQ ] tuple
-         $q_entries,   # list of [ HEAD, SEQ ] tuples
+    my ( $ref_entry,   # reference entry as [ SEQ_NAME, SEQ ] tuple
+         $q_entries,   # list of [ SEQ_NAME, SEQ ] tuples
          $args,        # argument hash
        ) = @_;
 
@@ -479,16 +480,16 @@ sub align_tile
 
         $type = Maasha::Seq::seq_guess_type( $seq1 );
 
-        if ( $type eq "rna" ) {
+        if ( $type eq "RNA" ) {
             $seq2 = Maasha::Seq::rna_revcomp( $seq1 );
-        } elsif ( $type eq "dna" ) {
+        } elsif ( $type eq "DNA" ) {
             $seq2 = Maasha::Seq::dna_revcomp( $seq1 );
         } else {
             Maasha::Common::error( qq(Bad sequence type->$type) );
         }
 
-        $align1 = Maasha::Align::align_muscle( [ $ref_entry, [ $entry->[ HEAD ] . "_+", $seq1 ] ], "-quiet -maxiters 1" );
-        $align2 = Maasha::Align::align_muscle( [ $ref_entry, [ $entry->[ HEAD ] . "_-", $seq2 ] ], "-quiet -maxiters 1" );
+        $align1 = Maasha::Align::align_muscle( [ $ref_entry, [ $entry->[ SEQ_NAME ] . "_+", $seq1 ] ], "-quiet -maxiters 1" );
+        $align2 = Maasha::Align::align_muscle( [ $ref_entry, [ $entry->[ SEQ_NAME ] . "_-", $seq2 ] ], "-quiet -maxiters 1" );
 
         if ( $args->{ "supress_indels" } )
         {
@@ -509,7 +510,7 @@ sub align_tile
 
             $align1->[ 1 ]->[ SEQ ] =~ s/-{$gaps}$// if $gaps;
             
-            $entry->[ HEAD ] = "$align1->[ 1 ]->[ HEAD ]_$sim1";
+            $entry->[ SEQ_NAME ] = "$align1->[ 1 ]->[ SEQ_NAME ]_$sim1";
             $entry->[ SEQ ]  = $align1->[ 1 ]->[ SEQ ];
 
             push @entries, $entry;
@@ -520,7 +521,7 @@ sub align_tile
 
             $align2->[ 1 ]->[ SEQ ] =~ s/-{$gaps}$// if $gaps;
 
-            $entry->[ HEAD ] = "$align2->[ 1 ]->[ HEAD ]_$sim2";
+            $entry->[ SEQ_NAME ] = "$align2->[ 1 ]->[ SEQ_NAME ]_$sim2";
             $entry->[ SEQ ]  = $align2->[ 1 ]->[ SEQ ];
 
             push @entries, $entry;