worked on sam - added 1 to phred score check in seq.rb

[biopieces.git] / code_ruby / lib / maasha / sam.rb

diff --git a/code_ruby/lib/maasha/sam.rb b/code_ruby/lib/maasha/sam.rb
index 79aa5ad26353212c243be40b7804b58ddeef1893..272d8fc977f859141f2f9b3100f24013099a4fae 100644 (file)
--- a/code_ruby/lib/maasha/sam.rb
+++ b/code_ruby/lib/maasha/sam.rb
@@ -29,6 +29,7 @@
 require 'pp'
 require 'maasha/filesys'
 require 'maasha/seq'
+require 'maasha/cigar'
 
 REGEX_HEADER  = Regexp.new(/^@[A-Za-z][A-Za-z](\t[A-Za-z][A-Za-z0-9]:[ -~]+)+$/)
 REGEX_COMMENT = Regexp.new(/^@CO\t.*/)
@@ -88,13 +89,10 @@ class Sam < Filesys
   # http://code.google.com/p/biopieces/wiki/KissFormat
   def self.align_descriptors(sam)
     offset     = 0
-    insertions = 0
     align      = []
 
     # Insertions
-    sam[:CIGAR].scan(/([0-9]+)([MIDNSHPX=])/).each do |len, op|
-      len = len.to_i
-
+    sam[:CIGAR].each do |len, op|
       if op == 'I'
         (0 ... len).each_with_index do |i|
           nt = sam[:SEQ].seq[offset + i]
@@ -346,7 +344,6 @@ class Sam < Filesys
     check_rname(rname)
     check_pos(pos)
     check_mapq(mapq)
-    check_cigar(cigar, seq)
     check_rnext(rnext)
     check_pnext(pnext)
     check_tlen(tlen)
@@ -359,7 +356,7 @@ class Sam < Filesys
     entry[:RNAME] = rname
     entry[:POS]   = pos
     entry[:MAPQ]  = mapq
-    entry[:CIGAR] = cigar
+    entry[:CIGAR] = Cigar.new(cigar)
     entry[:RNEXT] = rnext
     entry[:PNEXT] = pnext
     entry[:TLEN]  = tlen
@@ -425,44 +422,6 @@ class Sam < Filesys
     raise SamError, "Bad mapq: #{mapq}" unless (0 .. 2**8 - 1).include? mapq
   end
 
-  # Method to check cigar string.
-  def check_cigar(cigar, seq)
-    raise SamError, "Bad cigar: #{cigar}" unless cigar =~ /^(\*|([0-9]+[MIDNSHPX=])+)$/
-
-    unless cigar == '*'
-      check_cigar_hard_clip(cigar)
-      check_cigar_soft_clip(cigar)
-      check_cigar_seq_len(cigar, seq) unless seq == '*'
-    end
-  end
-
-  # Method to check cigar hard clipping only at ends.
-  def check_cigar_hard_clip(cigar)
-    if cigar.gsub(/^[0-9]+H|[0-9]+H$/, "").match('H')
-      raise SamError, "Bad cigar with internal H: #{cigar}"
-    end
-  end
-
-  # Method to check cigar soft clipping only at ends or H.
-  def check_cigar_soft_clip(cigar)
-    if cigar.gsub(/^[0-9]+H|[0-9]+H$/, "").gsub(/^[0-9]+S|[0-9]+S$/, "").match('S')
-      raise SamError, "Bad cigar with internal S: #{cigar}"
-    end
-  end
-
-  # Method to check cigar length matches sequence length.
-  def check_cigar_seq_len(cigar, seq)
-    cigar_len = 0
-
-    cigar.scan(/([0-9]+)([MIDNSHPX=])/).each do |len, op|
-      cigar_len += len.to_i if op =~ /[MISX=]/
-    end
-
-    if cigar_len != seq.length
-      raise SamError, "cigar and sequence length mismatch: #{cigar_len} != #{seq.length}"
-    end
-  end
-
   # Method to check if rnext, when not '*' or '='
   # and @SQ header lines are present, is located
   # in the header hash.