rewrite of FASTQ internals

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

diff --git a/code_ruby/lib/maasha/seq.rb b/code_ruby/lib/maasha/seq.rb
index 4a0af47ed0f5fde9d5ffe87ca19a3caec087f6d6..c8cc9f1e9943959c1b3ae1009cce6ba098ef0308 100644 (file)
--- a/code_ruby/lib/maasha/seq.rb
+++ b/code_ruby/lib/maasha/seq.rb
@@ -73,7 +73,7 @@ class SeqError < StandardError; end
 
 class Seq
   # Quality scores bases
-  SCORE_BASE = 64
+  SCORE_BASE = 33
   SCORE_MIN  = 0
   SCORE_MAX  = 40
 
@@ -550,61 +550,51 @@ class Seq
   def qual_base33?
     self.qual.match(/[!-:]/) ? true : false
   end
-
-  # Method that determines if a quality score string can be
-  # absolutely identified as base 64.
+ 
+  # Method that determines if a quality score string may be base 64.
   def qual_base64?
     self.qual.match(/[K-h]/) ? true : false
   end
 
-  # Method to determine if a quality score is valid.
+  # Method to determine if a quality score is valid accepting only 0-40 range.
   def qual_valid?(encoding)
     raise SeqError, "Missing qual" if self.qual.nil?
 
-    case encoding.downcase
-    when "sanger"     then return true if self.qual.match(/^[!-~]*$/)
-    when "454"        then return true if self.qual.match(/^[@-~]*$/)
-    when "solexa"     then return true if self.qual.match(/^[;-~]*$/)
-    when "illumina13" then return true if self.qual.match(/^[@-~]*$/)
-    when "illumina15" then return true if self.qual.match(/^[@-~]*$/)
-    when "illumina18" then return true if self.qual.match(/^[!-~]*$/)
+    case encoding
+    when :base_33 then return true if self.qual.match(/^[!-I]*$/)
+    when :base_64 then return true if self.qual.match(/^[@-h]*$/)
     else raise SeqError, "unknown quality score encoding: #{encoding}"
     end
 
     false
   end
 
-  # Method to convert quality scores inbetween formats.
-  # Sanger     base 33, range  0-40 
-  # 454        base 64, range  0-40 
-  # Solexa     base 64, range -5-40 
-  # Illumina13 base 64, range  0-40 
-  # Illumina15 base 64, range  0-40 
-  # Illumina18 base 33, range  0-41 
-  def convert_scores!(from, to)
-    unless from == to
-      na_qual = NArray.to_na(self.qual, "byte")
+  # Method to coerce quality scores to be within the 0-40 range.
+  def qual_coerce!(encoding)
+    raise SeqError, "Missing qual" if self.qual.nil?
 
-      case from.downcase
-      when "sanger"     then na_qual -= 33
-      when "454"        then na_qual -= 64
-      when "solexa"     then na_qual -= 64
-      when "illumina13" then na_qual -= 64
-      when "illumina15" then na_qual -= 64
-      when "illumina18" then na_qual -= 33
-      else raise SeqError, "unknown quality score encoding: #{from}"
-      end
+    case encoding
+    when :base_33 then self.qual.tr!("[J-~]", "I")
+    when :base_64 then self.qual.tr!("[i-~]", "h")
+    else raise SeqError, "unknown quality score encoding: #{encoding}"
+    end 
 
-      case to.downcase
-      when "sanger"     then na_qual += 33
-      when "454"        then na_qual += 64
-      when "solexa"     then na_qual += 64
-      when "illumina13" then na_qual += 64
-      when "illumina15" then na_qual += 64
-      when "illumina18" then na_qual += 33
-      else raise SeqError, "unknown quality score encoding: #{to}"
-      end
+    self
+  end
 
+  # Method to convert quality scores.
+  def qual_convert!(from, to)
+    raise SeqError, "unknown quality score encoding: #{from}" unless from == :base_33 or from == :base_64
+    raise SeqError, "unknown quality score encoding: #{to}"   unless to   == :base_33 or to   == :base_64
+
+    if from == :base_33 and to == :base_64
+      na_qual   = NArray.to_na(self.qual, "byte")
+      na_qual  += 64 - 33
+      self.qual = na_qual.to_s
+    elsif from == :base_64 and to == :base_33
+      self.qual.tr!("[;-?]", "@")  # Handle negative Solexa values from -5 to -1 (set these to 0).
+      na_qual   = NArray.to_na(self.qual, "byte")
+      na_qual  -= 64 - 33
       self.qual = na_qual.to_s
     end
 
@@ -617,6 +607,7 @@ class Seq
 
     na_qual = NArray.to_na(self.qual, "byte")
     na_qual -= SCORE_BASE
+
     na_qual.mean
   end