fixed SAM bug

git-svn-id: http://biopieces.googlecode.com/svn/trunk@756 74ccb610-7750-0410-82ae-013aeee3265d

diff --git a/code_perl/Maasha/SAM.pm b/code_perl/Maasha/SAM.pm
index e794d75f69b066b8c94ca4419505166efb5315f2..be1002aa4d026258e094bd357a4077dfaa6dd6be 100644 (file)
--- a/code_perl/Maasha/SAM.pm
+++ b/code_perl/Maasha/SAM.pm
@@ -143,6 +143,11 @@ sub sam2biopiece
     $record->{ 'READ_FAIL' }      = $record->{ 'FLAG' } & 0x0200;
     $record->{ 'READ_BAD' }       = $record->{ 'FLAG' } & 0x0400;
 
+    if ( $record->{ 'READ_PAIRED' } ) {
+        $record->{ 'BLOCK_COUNT' } = 2;
+    } else {
+        $record->{ 'BLOCK_COUNT' } = 1;
+    }
 
     # extra fields - (hate hate hate SAM!!!)
     
@@ -175,9 +180,10 @@ sub sam2align
 
     # Returns a string.
 
-    my ( $offset, $len, $op, $i, $nt, $nt2, @align, $aln_str );
+    my ( $offset, $len, $op, $i, $nt, $nt2, @nts, $insertions, @align, $aln_str );
 
-    $offset = 0;
+    $offset     = 0;
+    $insertions = 0;
 
     while ( length( $cigar ) > 0 )
     {
@@ -206,31 +212,41 @@ sub sam2align
 
     while ( $nm > 0 )
     {
-        if ( $md =~ s/^(\d+)// )
+        if ( $md =~ s/^(\d+)// )   # match
         {
             $offset += $1;
         }
-        elsif ( $md =~ s/^(\w)// )
+        elsif( $md =~ s/^\^([ATCGN]+)// )  # insertions
         {
-            $nt = $1;
+            @nts = split //, $1;
 
-            $nt2 = substr $seq, $offset, 1;
+            foreach $nt ( @nts )
+            {
+                $nt2 = substr $seq, $offset, 1;
 
-            push @align, [ $offset, "$nt>$nt2" ];
+                push @align, [ $offset, "$nt2>-" ];
+        
+                $offset++;
 
-            $offset++;
+                $insertions++;
 
-            $nm--;
+                $nm--;
+            }
         }
-        elsif( $md =~ s/^\^// )
+        elsif ( $md =~ s/^([ATCGN]+)// ) # mismatch
         {
-            $nt = substr $seq, $offset, 1;
+            @nts = split //, $1;
 
-            push @align, [ $offset, "$nt>-" ];
-        
-            $offset++;
+            foreach $nt ( @nts )
+            {
+                $nt2 = substr $seq, $offset - $insertions, 1;
+
+                push @align, [ $offset, "$nt>$nt2" ];
 
-            $nm--;
+                $offset++;
+
+                $nm--;
+            }
         }
         else
         {
@@ -238,7 +254,6 @@ sub sam2align
         }
     }
 
-
     $aln_str = "";
 
     @align = sort { $a->[ 0 ] <=> $b->[ 0 ] } @align;