1 /****
2 DIAMOND protein aligner
3 Copyright (C) 2013-2017 Benjamin Buchfink <buchfink@gmail.com>
4 
5 This program is free software: you can redistribute it and/or modify
6 it under the terms of the GNU General Public License as published by
7 the Free Software Foundation, either version 3 of the License, or
8 (at your option) any later version.
9 
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.
14 
15 You should have received a copy of the GNU General Public License
16 along with this program.  If not, see <http://www.gnu.org/licenses/>.
17 ****/
18 
19 #ifndef INTERVAL_PARTITION_H_
20 #define INTERVAL_PARTITION_H_
21 
22 #include <assert.h>
23 #include <map>
24 #include <limits.h>
25 #include "interval.h"
26 
27 struct IntervalNode
28 {
IntervalNodeIntervalNode29 	IntervalNode():
30 		count(0),
31 		min_score(INT_MAX),
32 		max_score(0)
33 	{}
IntervalNodeIntervalNode34 	IntervalNode(int count, int min_score, int max_score):
35 		count(count),
36 		min_score(min_score),
37 		max_score(max_score)
38 	{
39 	}
addIntervalNode40 	IntervalNode add(int score, int cap) const
41 	{
42 		return IntervalNode(count + 1, count < cap ? std::min(min_score, score) : min_score, std::max(max_score, score));
43 	}
44 	int count, min_score, max_score;
45 };
46 
47 struct IntervalPartition : protected std::map<int, IntervalNode>
48 {
49 
50 	struct MaxScore {};
51 	struct MinScore {};
52 
53 	struct Iterator
54 	{
IteratorIntervalPartition::Iterator55 		Iterator(const_iterator i, const_iterator j, const IntervalPartition &parent) :
56 			i_(i),
57 			j_(j),
58 			parent_(parent)
59 		{
60 		}
goodIntervalPartition::Iterator61 		bool good() const
62 		{
63 			return i_ != parent_.end();
64 		}
65 		Iterator& operator++()
66 		{
67 			i_ = j_;
68 			if(j_ != parent_.end())
69 				++j_;
70 			return *this;
71 		}
72 		std::pair<interval, IntervalNode> operator*() const
73 		{
74 			return std::make_pair(interval(i_->first, j_ == parent_.end() ? INT_MAX : j_->first), i_->second);
75 		}
76 	private:
77 		const_iterator i_, j_;
78 		const IntervalPartition &parent_;
79 	};
80 
IntervalPartitionIntervalPartition81 	IntervalPartition(int cap):
82 		cap(cap)
83 	{
84 		(*this)[0] = IntervalNode();
85 	}
86 
insertIntervalPartition87 	void insert(interval k, int score)
88 	{
89 		iterator i = lower_bound(k.begin_);
90 		if (i == end())
91 			i = std::map<int, IntervalNode>::insert(std::make_pair(k.begin_, IntervalNode())).first;
92 		else if (i->first != k.begin_) {
93 			//assert(i != std::map<int, IntervalNode>::begin());
94 			i--;
95 			i = std::map<int, IntervalNode>::insert(std::make_pair(k.begin_, i->second)).first;
96 		}
97 		IntervalNode last;
98 		while (i != end() && i->first < k.end_) {
99 			last = i->second;
100 			i->second = i->second.add(score, cap);
101 			++i;
102 		}
103 		if (i == end() || i->first != k.end_)
104 			(*this)[k.end_] = last;
105 	}
106 
coveredIntervalPartition107 	int covered(interval k) const
108 	{
109 		Iterator i = begin(k.begin_);
110 		std::pair<interval, IntervalNode> l;
111 		int c = 0;
112 		while (i.good() && (l = *i).first.begin_ < k.end_) {
113 			if (l.second.count >= cap)
114 				c += k.overlap(l.first);
115 			++i;
116 		}
117 		return c;
118 	}
119 
coveredIntervalPartition120 	int covered(interval k, int max_score, const MaxScore&) const
121 	{
122 		Iterator i = begin(k.begin_);
123 		std::pair<interval, IntervalNode> l;
124 		int c = 0;
125 		while (i.good() && (l = *i).first.begin_ < k.end_) {
126 			if (l.second.max_score >= max_score)
127 				c += k.overlap(l.first);
128 			++i;
129 		}
130 		return c;
131 	}
132 
coveredIntervalPartition133 	int covered(interval k, int min_score, const MinScore&) const
134 	{
135 		Iterator i = begin(k.begin_);
136 		std::pair<interval, IntervalNode> l;
137 		int c = 0;
138 		while (i.good() && (l = *i).first.begin_ < k.end_) {
139 			if (l.second.count >= cap && l.second.min_score >= min_score)
140 				c += k.overlap(l.first);
141 			++i;
142 		}
143 		return c;
144 	}
145 
min_scoreIntervalPartition146 	int min_score(interval k) const
147 	{
148 		Iterator i = begin(k.begin_);
149 		std::pair<interval, IntervalNode> l;
150 		int s = INT_MAX;
151 		while (i.good() && (l = *i).first.begin_ < k.end_) {
152 			if (l.second.count < cap)
153 				return 0;
154 			s = std::min(s, l.second.min_score);
155 			++i;
156 		}
157 		return s;
158 	}
159 
max_scoreIntervalPartition160 	int max_score(interval k) const
161 	{
162 		Iterator i = begin(k.begin_);
163 		std::pair<interval, IntervalNode> l;
164 		int s = INT_MAX;
165 		while (i.good() && (l = *i).first.begin_ < k.end_) {
166 			s = std::min(s, l.second.max_score);
167 			++i;
168 		}
169 		assert(s != INT_MAX);
170 		return s;
171 	}
172 
beginIntervalPartition173 	Iterator begin(int p) const
174 	{
175 		const_iterator i = lower_bound(p), j;
176 		if (i == end() || i->first != p) {
177 			//assert(i != std::map<int, IntervalNode>::begin());
178 			j = i;
179 			i--;
180 		}
181 		else {
182 			j = i;
183 			++j;
184 		}
185 		return Iterator(i, j, *this);
186 	}
187 
188 	const int cap;
189 
190 };
191 
192 #endif