updated to v0.7.5

src/ChangeLog

@@ -1,3 +1,9 @@
+[0.7.5]
+- fixed a compiler error for Clang
+- compiles on OS X
+- fixed a bug that would not produce the correct sort order for HSPs with e-value 0
+- option --run-len defaults to 1 for sequences shorther than 30 letters
+
 [0.7.4]
 - removed OpenMP dependency
 - using explicit Makefile instead of autoconf

src/Makefile

@@ -1,11 +1,21 @@
-WARN=-Wall
+WARN=-Wall -Wno-uninitialized
 C=gcc
 CXX=g++
 CFLAGS=-O3 -DNDEBUG
 CXXFLAGS=-O3 -DNDEBUG -Iboost/include $(WARN)
-LINKFLAGS=-Lboost/lib
-LIBS=-Wl,-Bstatic -lboost_thread -lboost_system -lboost_timer -lboost_chrono -lboost_iostreams -lboost_program_options -Wl,-Bdynamic -lpthread -lz -lrt
-PROG=diamond
+
+ifneq ($(wildcard boost/.),)
+	LIBS=boost/lib/libboost_thread.a boost/lib/libboost_system.a boost/lib/libboost_chrono.a boost/lib/libboost_timer.a boost/lib/libboost_iostreams.a boost/lib/libboost_program_options.a
+else
+	LIBS=-lboost_thread -lboost_system -lboost_timer -lboost_chrono -lboost_iostreams -lboost_program_options
+endif
+
+LIBS+=-lpthread -lz
+
+ifneq (,$(findstring Darwin,$(shell uname)))
+else
+LIBS+=-lrt
+endif
 
 GCC_VER_GTE42 := $(shell echo `gcc -dumpversion | cut -f1-2 -d.` \>= 4.2 | bc )
 ifeq ($(GCC_VER_GTE42),1)
@@ -36,6 +46,6 @@ search/align.h	search/align_range.h  search/align_ungapped.h  search/collision.h
 	$(C) $(CFLAGS) -c -o $@ $<
 
 diamond : $(OBJS)
-	$(CXX) $(LINKFLAGS) -o ../bin/$@ $(OBJS) $(LIBS)
+	$(CXX) -o ../bin/$@ $(OBJS) $(LIBS)
 
-all : $(PROG)
+all : diamond

src/align/align_sequence.h

@@ -22,10 +22,8 @@ Author: Benjamin Buchfink
 #define ALIGN_SEQUENCE_H_
 
 #include <vector>
-#include <boost/thread/tss.hpp>
 #include "../dp/floating_sw.h"
 
-using boost::thread_specific_ptr;
 using std::vector;
 
 template<typename _val, typename _locr, typename _locl>
@@ -41,8 +39,8 @@ void align_sequence(vector<Segment<_val> > &matches,
 	std::sort(begin, end, hit<_locr,_locl>::cmp_normalized_subject);
 	const unsigned q_num (begin->query_);
 	const sequence<const _val> query (query_seqs<_val>::get()[q_num]);
-	unsigned frame = q_num % query_contexts();
-	unsigned query_len = query.length();
+	const unsigned frame = q_num % query_contexts();
+	const unsigned query_len = query.length();
 	padding[frame] = program_options::read_padding<_val>(query_len);
 
 	const Sequence_set<_val> *ref = ref_seqs<_val>::data_;
@@ -61,9 +59,9 @@ void align_sequence(vector<Segment<_val> > &matches,
 				Traceback ());
 		const int score = local.back().score_;
 		std::pair<size_t,size_t> l = ref_seqs<_val>::data_->local_position(i->subject_);
-		matches.push_back(Segment<_val> (score, frame, score_matrix::get().evalue(score, db_letters, query_len), &local.back(), l.first));
-		//matches.back().top_evalue_ = matches.back().evalue_;
+		matches.push_back(Segment<_val> (score, frame, &local.back(), l.first));
 		anchored_transform(local.back(), l.second, i->seed_offset_);
+		stat.inc(Statistics::ALIGNED_QLEN, local.back().query_len_);
 
 		//local.back().print(query, ref_seqs<_val>::get()[l.first]);
 

src/align/link_segments.h

@@ -55,13 +55,13 @@ double link_segments(Segment<_val> &h1, Segment<_val> &h2)
 template<typename _val>
 void link_segments(const typename vector<Segment<_val> >::iterator &begin, const typename vector<Segment<_val> >::iterator &end)
 {
-	double min_ev = begin->evalue_;
+	int max_score = begin->score_;
 	/*for(typename vector<match<_val> >::iterator i=begin; i<end; ++i)
 		if(i->top_evalue_ == -1)
 			for(typename vector<match<_val> >::iterator j=i+1; j<end; ++j)
 				min_ev = std::min(min_ev, link_segments(*i, *j));*/
 	for(typename vector<Segment<_val> >::iterator i=begin; i<end; ++i)
-		i->top_evalue_ = min_ev;
+		i->top_score_ = max_score;
 }
 
 template<typename _val>

src/basic/const.h

@@ -25,7 +25,7 @@ struct Const
 {
 
 	enum {
-		build_version = 53,
+		build_version = 54,
 		build_compatibility = 52,
 		daa_version = 0,
 		seedp_bits = 10,
@@ -46,7 +46,7 @@ struct Const
 
 };
 
-const char* Const::version_string = "0.7.4";
+const char* Const::version_string = "0.7.5";
 const char* Const::program_name = "diamond";
 const char* Const::id_delimiters = " \a\b\f\n\r\t\v";
 

src/basic/match.h

@@ -194,18 +194,16 @@ struct local_match
 template<typename _val>
 struct Segment
 {
-	Segment(unsigned score,
+	Segment(int score,
 			unsigned frame,
-			double evalue = std::numeric_limits<double>::max(),
 			local_match<_val> *traceback = 0,
 			unsigned subject_id = std::numeric_limits<unsigned>::max()):
 		score_ (score),
 		frame_ (frame),
 		traceback_ (traceback),
 		subject_id_ (subject_id),
-		evalue_ (evalue),
 		next_ (0),
-		top_evalue_ (-1)
+		top_score_ (0)
 	{ }
 	Strand strand() const
 	{ return frame_ < 3 ? FORWARD : REVERSE; }
@@ -214,23 +212,22 @@ struct Segment
 	interval subject_range() const
 	{ return traceback_->subject_range(); }
 	bool operator<(const Segment &rhs) const
-	{ return top_evalue_ < rhs.top_evalue_
-			|| (top_evalue_ == rhs.top_evalue_
-			&& (subject_id_ < rhs.subject_id_ || (subject_id_ == rhs.subject_id_ && (evalue_ < rhs.evalue_ || (evalue_ == rhs.evalue_ && score_ > rhs.score_))))); }
+	{ return top_score_ > rhs.top_score_
+			|| (top_score_ == rhs.top_score_
+			&& (subject_id_ < rhs.subject_id_ || (subject_id_ == rhs.subject_id_ && (score_ > rhs.score_ || (score_ == rhs.score_ && traceback_->score_ > rhs.traceback_->score_))))); }
 	static bool comp_subject(const Segment& lhs, const Segment &rhs)
 	{ return lhs.subject_id_ < rhs.subject_id_ || (lhs.subject_id_ == rhs.subject_id_ && lhs.score_ > rhs.score_); }
 	struct Subject
 	{
 		unsigned operator()(const Segment& x) const
 		{ return x.subject_id_; }
 	};
-	unsigned				score_;
+	int						score_;
 	unsigned				frame_;
 	local_match<_val>	   *traceback_;
 	unsigned				subject_id_;
-	double					evalue_;
 	Segment				   *next_;
-	double					top_evalue_;
+	int						top_score_;
 };
 
 #endif /* MATCH_H_ */

src/basic/options.h

@@ -111,7 +111,9 @@ namespace program_options
 	inline unsigned get_run_len(unsigned length)
 	{
 		if(run_len == 0) {
-			if(length < 100)
+			if(length < 30)
+				return 1;
+			else if(length < 100)
 				return 20;
 			else
 				return 40;

src/basic/statistics.h

@@ -27,7 +27,7 @@ struct Statistics
 {
 
 	enum value { SEED_HITS, TENTATIVE_MATCHES0, TENTATIVE_MATCHES1, TENTATIVE_MATCHES2, TENTATIVE_MATCHES3, MATCHES, ALIGNED, GAPPED, DUPLICATES,
-		GAPPED_HITS, QUERY_SEEDS, QUERY_SEEDS_HIT, REF_SEEDS, REF_SEEDS_HIT, QUERY_SIZE, REF_SIZE, OUT_HITS, OUT_MATCHES, COLLISION_LOOKUPS, QCOV, BIAS_ERRORS, SCORE_TOTAL, COUNT };
+		GAPPED_HITS, QUERY_SEEDS, QUERY_SEEDS_HIT, REF_SEEDS, REF_SEEDS_HIT, QUERY_SIZE, REF_SIZE, OUT_HITS, OUT_MATCHES, COLLISION_LOOKUPS, QCOV, BIAS_ERRORS, SCORE_TOTAL, ALIGNED_QLEN, COUNT };
 
 	Statistics()
 	{ memset(data_, 0, sizeof(data_)); }
@@ -59,6 +59,7 @@ struct Statistics
 		log_stream << "Net hits = " << data_[OUT_HITS] << endl;
 		log_stream << "Matches = " << data_[OUT_MATCHES] << endl;
 		log_stream << "Total score = " << data_[SCORE_TOTAL] << endl;
+		log_stream << "Aligned query len = " << data_[ALIGNED_QLEN] << endl;
 		log_stream << "Gapped matches = " << data_[GAPPED] << endl;
 		verbose_stream << "Final matches = " << data_[MATCHES] << endl;
 		verbose_stream << "Queries aligned = " << data_[ALIGNED] << endl;

src/output/daa_record.h

@@ -70,9 +70,7 @@ struct DAA_query_record
 	DAA_query_record(const DAA_file& file, const Binary_buffer &buf):
 		file_ (file),
 		it_ (init(buf))
-	{
-		//cout << sequence<Nucleotide>(source_seq.data(), source_seq.size());
-	}
+	{ }
 
 	Match_iterator begin() const
 	{ return Match_iterator (*this, it_); }
@@ -107,7 +105,8 @@ struct DAA_query_record
 	const Binary_buffer::Iterator it_;
 
 	friend struct DAA_match_record<_val>;
-	friend Binary_buffer::Iterator& operator>>(Binary_buffer::Iterator &it, DAA_match_record<_val> &r);
+
+	template<typename _val2> friend Binary_buffer::Iterator& operator>>(Binary_buffer::Iterator &it, DAA_match_record<_val2> &r);
 
 };
 
@@ -148,32 +147,6 @@ struct DAA_match_record
 
 private:
 
-	friend Binary_buffer::Iterator& operator>>(Binary_buffer::Iterator &it, DAA_match_record &r)
-	{
-		uint32_t subject_id;
-		it >> subject_id;
-		uint8_t flag;
-		it >> flag;
-		it.read_packed(flag & 3, r.score);
-		it.read_packed((flag>>2)&3, r.query_begin);
-		it.read_packed((flag>>4)&3, r.subject_begin);
-		r.transcript.read(it);
-		r.subject_name = r.parent_.file_.ref_name(subject_id);
-		r.total_subject_len = r.parent_.file_.ref_len(subject_id);
-		if(r.parent_.file_.mode() == blastx) {
-			r.frame = (flag&(1<<6)) == 0 ? r.query_begin % 3 : 3+(r.parent_.source_seq.size() - 1 - r.query_begin)%3;
-			r.translated_query_begin = query_translated_begin<_val>(r.query_begin, r.frame, r.parent_.source_seq.size(), true);
-		} else if (r.parent_.file_.mode() == blastp) {
-			r.frame = 0;
-			r.translated_query_begin = r.query_begin;
-		} else {
-			r.frame = (flag&(1<<6)) == 0 ? 0 : 1;
-			r.translated_query_begin = query_translated_begin<_val>(r.query_begin, r.frame, r.parent_.source_seq.size(), false);
-		}
-		r.parse();
-		return it;
-	}
-
 	void parse()
 	{
 		translated_query_len = 0;
@@ -214,6 +187,35 @@ struct DAA_match_record
 
 	const DAA_query_record<_val> &parent_;
 
+	template<typename _val2> friend Binary_buffer::Iterator& operator>>(Binary_buffer::Iterator &it, DAA_match_record<_val2> &r);
+
 };
 
+template<typename _val>
+Binary_buffer::Iterator& operator>>(Binary_buffer::Iterator &it, DAA_match_record<_val> &r)
+{
+	uint32_t subject_id;
+	it >> subject_id;
+	uint8_t flag;
+	it >> flag;
+	it.read_packed(flag & 3, r.score);
+	it.read_packed((flag>>2)&3, r.query_begin);
+	it.read_packed((flag>>4)&3, r.subject_begin);
+	r.transcript.read(it);
+	r.subject_name = r.parent_.file_.ref_name(subject_id);
+	r.total_subject_len = r.parent_.file_.ref_len(subject_id);
+	if(r.parent_.file_.mode() == blastx) {
+		r.frame = (flag&(1<<6)) == 0 ? r.query_begin % 3 : 3+(r.parent_.source_seq.size() - 1 - r.query_begin)%3;
+		r.translated_query_begin = query_translated_begin<_val>(r.query_begin, r.frame, r.parent_.source_seq.size(), true);
+	} else if (r.parent_.file_.mode() == blastp) {
+		r.frame = 0;
+		r.translated_query_begin = r.query_begin;
+	} else {
+		r.frame = (flag&(1<<6)) == 0 ? 0 : 1;
+		r.translated_query_begin = query_translated_begin<_val>(r.query_begin, r.frame, r.parent_.source_seq.size(), false);
+	}
+	r.parse();
+	return it;
+}
+
 #endif /* DAA_RECORD_H_ */

src/output/output_format.h

@@ -40,6 +40,9 @@ template<typename _val>
 struct Blast_tab_format : public Output_format<_val>
 {
 
+	Blast_tab_format()
+	{ }
+
 	virtual void print_match(const DAA_match_record<_val> &r, Text_buffer &out) const
 	{
 		out << r.query_name() << '\t'
@@ -79,6 +82,9 @@ template<typename _val>
 struct Sam_format : public Output_format<_val>
 {
 
+	Sam_format()
+	{ }
+
 	virtual void print_match(const DAA_match_record<_val> &r, Text_buffer &out) const
 	{
 		out << r.query_name() << '\t'

src/run/master_thread.h

@@ -278,6 +278,9 @@ void master_thread()
 	cpu_timer timer2, timer_mapping;
 	timer_mapping.stop();
 
+	if(!check_dir(program_options::tmpdir))
+		throw std::runtime_error("Temporary directory " + program_options::tmpdir + " does not exist or is not a directory. Please use option -t to specify a different directory.");
+
 	task_timer timer ("Opening the database", 1);
 	Database_file db_file;
 	timer.finish();

src/util/binary_file.h

@@ -101,7 +101,7 @@ struct Input_stream : public io::filtering_istream
 	size_t read(_t *ptr, size_t count)
 	{
 		ssize_t n;
-		if((n = io::read(*this, reinterpret_cast<char*>(ptr), sizeof(_t) * count)) != count * sizeof(_t)) {
+		if((n = io::read(*this, reinterpret_cast<char*>(ptr), sizeof(_t) * count)) != (ssize_t)(count * sizeof(_t))) {
 			if(n == EOF)
 				return 0;
 			else if(n >= 0 && this->get() == EOF && (n % sizeof(_t) == 0))
@@ -213,7 +213,7 @@ struct Output_stream : public io::filtering_ostream
 	void write(const _t *ptr, size_t count)
 	{
 		ssize_t n;
-		if((n = io::write(*this, reinterpret_cast<const char*>(ptr), sizeof(_t) * count)) != sizeof(_t) * count)
+		if((n = io::write(*this, reinterpret_cast<const char*>(ptr), sizeof(_t) * count)) != (ssize_t)(sizeof(_t) * count))
 			throw File_write_exception (file_name_.c_str(), sizeof(_t) * count, n);
 	}
 

src/util/hash_table.h

@@ -24,9 +24,7 @@ Author: Benjamin Buchfink
 template<typename T, T value> struct value_compare
 {
 	bool operator()(T x) const
-	{
-		return x == value;
-	}
+	{ return x == value; }
 };
 
 template<typename _K, typename _V, typename _E, typename _H> class hash_table
@@ -43,14 +41,10 @@ template<typename _K, typename _V, typename _E, typename _H> class hash_table
 	hash_table(size_t size):
 		table (new entry[size]),
 		size_ (size)
-	{
-		memset(table, 0, size_ * sizeof(entry));
-	}
+	{ memset(table, 0, size_ * sizeof(entry)); }
 
 	~hash_table()
-	{
-		delete[] table;
-	}
+	{ delete[] table; }
 
 	entry* operator[](_K key) const
 	{
@@ -84,19 +78,6 @@ template<typename _K, typename _V, typename _E, typename _H> class hash_table
 
 private:
 
-	hash_table(const hash_table &):
-		table (NULL),
-		size_ (0)
-	{
-		assert(false);
-	}
-
-	hash_table& operator=(const hash_table &)
-	{
-		assert(false);
-		return *this;
-	}
-
 	entry* get_entry(_K key) const
 	{
 		entry *p = &table[_H()(key) % size_];

src/util/merge_sort.h

@@ -31,7 +31,7 @@ void merge_sort(_it begin, _it end, unsigned n_threads, unsigned level = 0)
 	if(diff <= 1)
 		return;
 
-	if(1 << level >= n_threads) {
+	if(1u << level >= n_threads) {
 		std::sort(begin, end);
 		return;
 	}