X-Git-Url: https://git.donarmstrong.com/?a=blobdiff_plain;f=bp_bin%2Fremove_primers;h=1f7b7b216ab3fc9f089848b9303771d62d36a77f;hb=48bea5c28b89dc5586d0bddb338ccd6ba23aa1f9;hp=4e4e2e7fe48abdb101f4de2ae0a7e59e1097a1d3;hpb=bc41bf1f0f46fbf4c6468a0bc55bdb093c496084;p=biopieces.git diff --git a/bp_bin/remove_primers b/bp_bin/remove_primers index 4e4e2e7..1f7b7b2 100755 --- a/bp_bin/remove_primers +++ b/bp_bin/remove_primers @@ -32,6 +32,9 @@ require 'pp' require 'maasha/biopieces' require 'maasha/seq' +require 'maasha/seq/backtrack' + +include BackTrack casts = [] casts << {:long=>'forward', :short=>'f', :type=>'string', :mandatory=>true, :default=>nil, :allowed=>nil, :disallowed=>nil} @@ -44,25 +47,39 @@ options = Biopieces.options_parse(ARGV, casts) Biopieces.open(options[:stream_in], options[:stream_out]) do |input, output| input.each do |record| - if record.has_key? :SEQ + if record[:SEQ] + forward = false + reverse = false seq = Seq.new_bp(record) - pos = 0 - if forward = seq.patscan(options[:forward].to_s, pos, options[:mismatches], options[:insertions], options[:deletions]) - record[:FORWARD_POS] = forward.last.pos - record[:FORWARD_LEN] = forward.last.length - pos = forward.last.pos + forward.last.length - seq.subseq!(pos) + seq.patscan(options[:forward].to_s, 0, options[:mismatches], options[:insertions], options[:deletions]) do |match| + record[:FORWARD_POS] = match.pos + record[:FORWARD_LEN] = match.length + pos = match.pos + match.length + len = seq.length - pos + seq.subseq!(pos, len) if len > 0 + forward = true + break end - if reverse = seq.patscan(options[:reverse].to_s, pos, options[:mismatches], options[:insertions], options[:deletions]) - record[:REVERSE_POS] = reverse.first.pos - record[:REVERSE_LEN] = reverse.first.length - pos = reverse.first.pos - seq.subseq!(0, pos) + seq.patscan(options[:reverse].to_s, 0, options[:mismatches], options[:insertions], options[:deletions]) do |match| + record[:REVERSE_POS] = match.pos + record[:REVERSE_LEN] = match.length + pos = 0 + len = match.pos + + if len == 0 + seq.seq = "" + seq.qual = "" if seq.qual + else + seq.subseq!(pos, len) + end + + reverse = true + break end - if pos > 0 + if forward or reverse record.merge!(seq.to_bp) end