code_ruby/lib/maasha/cigar.rb

   1 # Copyright (C) 2007-2011 Martin A. Hansen.
   2
   3 # This program is free software; you can redistribute it and/or
   4 # modify it under the terms of the GNU General Public License
   5 # as published by the Free Software Foundation; either version 2
   6 # of the License, or (at your option) any later version.
   7
   8 # This program is distributed in the hope that it will be useful,
   9 # but WITHOUT ANY WARRANTY; without even the implied warranty of
  10 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  11 # GNU General Public License for more details.
  12
  13 # You should have received a copy of the GNU General Public License
  14 # along with this program; if not, write to the Free Software
  15 # Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
  16
  17 # http://www.gnu.org/copyleft/gpl.html
  18
  19 # >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>><<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
  20
  21 # This software is part of the Biopieces framework (www.biopieces.org).
  22
  23 # >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>><<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
  24
  25 # Methods to handle extended CIGAR format used in the
  26 # SAM format version v1.4-r962 - April 17, 2011
  27
  28 # http://samtools.sourceforge.net/SAM1.pdf
  29
  30 require 'pp'
  31
  32 # Error class for all exceptions to do with CIGAR.
  33 class CigarError < StandardError; end
  34
  35 # Class to parse and write SAM files.
  36 class Cigar
  37   attr_accessor :cigar
  38
  39   def initialize(cigar)
  40     @cigar = cigar
  41
  42     check_cigar
  43   end
  44
  45   def each
  46     cigar.scan(/(\d+)([MIDNSHPX=])/).each do |len, op|
  47       yield [len.to_i, op]
  48     end
  49   end
  50
  51   def length
  52     total = 0
  53
  54     self.each do |len, op|
  55       total += len unless op == 'I'
  56     end
  57
  58     total
  59   end
  60
  61   def matches
  62     total = 0
  63
  64     self.each do |len, op|
  65       total += len if op == 'M'
  66     end
  67
  68     total
  69   end
  70
  71   def insertions
  72     total = 0
  73
  74     self.each do |len, op|
  75       total += len if op == 'I'
  76     end
  77
  78     total
  79   end
  80
  81   def deletions
  82     total = 0
  83
  84     self.each do |len, op|
  85       total += len if op == 'D'
  86     end
  87
  88     total
  89   end
  90
  91   def clip_hard
  92     total = 0
  93
  94     self.each do |len, op|
  95       total += len if op == 'H'
  96     end
  97
  98     total
  99   end
 100
 101   def clip_soft
 102     total = 0
 103
 104     self.each do |len, op|
 105       total += len if op == 'S'
 106     end
 107
 108     total
 109   end
 110
 111   private
 112
 113   def check_cigar
 114     unless  cigar =~ /^(\*|([0-9]+[MIDNSHPX=])+)$/
 115       raise CigarError, "Bad cigar format: #{cigar}"
 116     end
 117
 118     if cigar.gsub(/^[0-9]+H|[0-9]+H$/, "").match('H')
 119       raise CigarError, "Bad cigar with internal H: #{cigar}"
 120     end
 121
 122     if cigar.gsub(/^[0-9]+H|[0-9]+H$/, "").gsub(/^[0-9]+S|[0-9]+S$/, "").match('S')
 123       raise CigarError, "Bad cigar with internal S: #{cigar}"
 124     end
 125   end
 126 end