From: martinahansen Date: Mon, 9 Nov 2009 12:26:57 +0000 (+0000) Subject: fixed SAM bug X-Git-Url: https://git.donarmstrong.com/?a=commitdiff_plain;h=2729afb13b11b2840f2be0b431a68c435aecb9fe;p=biopieces.git 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 e794d75..be1002a 100644 --- 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;