3 # Copyright (C) 2007-2011 Martin A. Hansen.
5 # This program is free software; you can redistribute it and/or
6 # modify it under the terms of the GNU General Public License
7 # as published by the Free Software Foundation; either version 2
8 # of the License, or (at your option) any later version.
10 # This program is distributed in the hope that it will be useful,
11 # but WITHOUT ANY WARRANTY; without even the implied warranty of
12 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
13 # GNU General Public License for more details.
15 # You should have received a copy of the GNU General Public License
16 # along with this program; if not, write to the Free Software
17 # Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
19 # http://www.gnu.org/copyleft/gpl.html
21 # >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>><<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
23 # This program is part of the Biopieces framework (www.biopieces.org).
25 # >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> DESCRIPTION <<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
27 # Usearch sequences in the stream against a specified database.
29 # >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>><<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
31 require 'maasha/biopieces'
32 require 'maasha/fasta'
37 def initialize(infile, outfile, options)
44 # Method to execute database search.
46 @command << "usearch --query #{@infile} --db #{@options[:database]} --userout #{@outfile}"
47 @command << "--userfields target+tloz+thiz+query+bits+strand"
48 @command << "--id #{@options[:identity]}" if @options.has_key? :identity
49 @command << "--evalue #{@options[:e_val]}" if @options.has_key? :e_val
54 # Method to parse a Useach .uc file and for each line of data
55 # yield a Biopiece record.
59 File.open(@outfile, mode="r") do |ios|
60 ios.gets # skip comment line
61 ios.each_line do |line|
62 fields = line.chomp.split("\t")
64 record[:REC_TYPE] = "USEARCH"
65 record[:S_ID] = fields[0]
66 record[:S_BEG] = fields[1].to_i
67 record[:S_END] = fields[2].to_i
68 record[:Q_ID] = fields[3]
69 record[:SCORE] = fields[4].to_f
70 record[:STRAND] = fields[5]
81 # Method to execute a command using a system() call.
82 # The command is composed of bits from the @command variable.
84 @command.unshift "nice -n 19"
86 @command << "> /dev/null 2>&1" unless @options[:verbose]
87 command = @command.join(" ")
89 raise "Command failed: #{command}" unless $?.success?
96 casts << {:long=>'database', :short=>'d', :type=>'file!', :mandatory=>false, :default=>nil, :allowed=>nil, :disallowed=>nil}
97 casts << {:long=>'identity', :short=>'i', :type=>'float', :mandatory=>false, :default=>nil, :allowed=>nil, :disallowed=>nil}
98 casts << {:long=>'e_val', :short=>'e', :type=>'float', :mandatory=>false, :default=>nil, :allowed=>nil, :disallowed=>nil}
100 options = Biopieces.options_parse(ARGV, casts)
102 raise ArgumentError, "--identity or --e_val must be specified" unless options[:identity] or options[:e_val]
104 tmpdir = Biopieces.mktmpdir
105 infile = File.join(tmpdir, "in.fna")
106 outfile = File.join(tmpdir, "out.uc")
108 Biopieces.open(options[:stream_in], options[:stream_out]) do |input, output|
109 Fasta.open(infile, mode="w") do |fasta_io|
110 input.each_record do |record|
113 if record.has_key? :SEQ_NAME and record.has_key? :SEQ
114 fasta_io.puts Seq.new_bp(record).to_fasta
119 uc = Usearch.new(infile, outfile, options)
129 # >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>><<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<