NVBIO
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
alignments.h
Go to the documentation of this file.
1 /*
2  * nvbio
3  * Copyright (c) 2011-2014, NVIDIA CORPORATION. All rights reserved.
4  *
5  * Redistribution and use in source and binary forms, with or without
6  * modification, are permitted provided that the following conditions are met:
7  * * Redistributions of source code must retain the above copyright
8  * notice, this list of conditions and the following disclaimer.
9  * * Redistributions in binary form must reproduce the above copyright
10  * notice, this list of conditions and the following disclaimer in the
11  * documentation and/or other materials provided with the distribution.
12  * * Neither the name of the NVIDIA CORPORATION nor the
13  * names of its contributors may be used to endorse or promote products
14  * derived from this software without specific prior written permission.
15  *
16  * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND
17  * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
18  * WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
19  * DISCLAIMED. IN NO EVENT SHALL NVIDIA CORPORATION BE LIABLE FOR ANY
20  * DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
21  * (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
22  * LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND
23  * ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
24  * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
25  * SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
26  */
27 
28 //#include <nvBowtie/bowtie2/cuda/defs.h>
29 #include <nvbio/basic/cuda/arch.h>
30 #include <nvbio/basic/cached_iterator.h>
31 #include <nvbio/basic/packedstream.h>
32 #include <nvbio/basic/pod.h>
33 #include <algorithm>
34 
35 #pragma once
36 
37 namespace nvbio {
38 namespace io {
39 
42 
46 enum MDS_OP {
47  MDS_MATCH = 0,
48  MDS_MISMATCH = 1,
49  MDS_INSERTION = 2,
50  MDS_DELETION = 3,
51  MDS_INVALID = 4
52 };
53 
57 struct Cigar
58 {
61  enum Operation {
62  SUBSTITUTION = 0,
63  INSERTION = 1,
64  DELETION = 2,
65  SOFT_CLIPPING = 3,
66  };
67 
68  NVBIO_FORCEINLINE NVBIO_HOST_DEVICE
69  Cigar() {}
70 
71  NVBIO_FORCEINLINE NVBIO_HOST_DEVICE
72  Cigar(const uint8 type, const uint16 len) : m_type( type ), m_len( len ) {}
73 
74  uint16 m_type:2, m_len:14;
75 };
76 
80 struct Alignment
81 {
82  NVBIO_FORCEINLINE NVBIO_HOST_DEVICE
83  static uint32 max_ed() { return 255u; }
84 
85  NVBIO_FORCEINLINE NVBIO_HOST_DEVICE
86  static int32 max_score() { return (1 << 17) - 1; }
87 
88  NVBIO_FORCEINLINE NVBIO_HOST_DEVICE
89  static int32 min_score() { return -((1 << 17) - 1); }
90 
91  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE Alignment() {}
92  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE Alignment(
93  const uint32 pos,
94  const uint32 ed,
95  const int32 score,
96  const uint32 rc,
97  const uint32 mate = 0,
98  const bool paired = false,
99  const bool discordant = false)
100  {
101  //NVBIO_CUDA_ASSERT( ed < 1024u );
102  //NVBIO_CUDA_ASSERT( score >= min_score() && score <= max_score() );
103  m_align = pos;
104  m_ed = ed;
105  m_score = score < 0 ? uint32( -score ) : uint32( score );
106  m_score_sgn = score < 0 ? 1u : 0u;
107  m_rc = rc;
108  m_mate = mate;
109  m_paired = paired ? 1u : 0u;
110  m_discordant = discordant ? 1u : 0u;
111  }
112 
113  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE int32 score() const { return m_score_sgn ? -int32( m_score ) : int32( m_score ); }
114  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE bool is_aligned() const { return m_align != uint32(-1); }
115  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE uint32 alignment() const { return m_align; }
116  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE bool is_rc() const { return m_rc; }
117  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE uint32 ed() const { return m_ed; }
118  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE uint32 sink() const { return m_ed; }
119  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE uint32 mate() const { return m_mate; }
120  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE bool is_paired() const { return (m_paired == true) && is_aligned(); }
121  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE bool is_unpaired() const { return (m_paired == false) && is_aligned(); }
122  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE bool is_concordant() const { return (m_paired == true) && (m_discordant == false); }
123  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE bool is_discordant() const { return (m_paired == true) && (m_discordant == true); }
124 
125  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
126  static Alignment invalid() { return Alignment( uint32(-1), max_ed(), max_score(), 0u, 0u, false ); }
127 
128  uint32 m_score_sgn : 1, m_score:17, m_ed:10, m_rc:1, m_mate:1, m_paired:1, m_discordant:1;
129  uint32 m_align;
130 };
131 
135 struct AlignmentCompare
136 {
137  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
138  bool operator() (const Alignment f, const Alignment s) const { return (f.m_score > s.m_score); }
139 };
140 
144 struct BestAlignments
145 {
146  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
147  BestAlignments() {}
148 
149  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
150  BestAlignments(const Alignment& a1, const Alignment& a2) : m_a1( a1 ), m_a2( a2 ) {}
151 
152  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE bool is_aligned() const { return m_a1.is_aligned(); }
153  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE bool has_second() const { return m_a2.is_aligned(); }
154  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE int32 best_score() const { return m_a1.score(); }
155  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE uint32 best_ed() const { return m_a1.ed(); }
156  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE uint32 best_alignment_pos() const { return m_a1.alignment(); }
157  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE int32 second_score() const { return m_a2.score(); }
158  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE uint32 second_ed() const { return m_a2.ed(); }
159  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE uint32 second_alignment_pos() const { return m_a2.alignment(); }
160 
161  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE const Alignment& best() const { return m_a1; }
162  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE const Alignment& second_best() const { return m_a2; }
163 
164  template <uint32 I>
165  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
166  const Alignment& alignment() const { return I == 0 ? m_a1 : m_a2; }
167 
168  template <uint32 I>
169  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
170  Alignment& alignment() { return I == 0 ? m_a1 : m_a2; }
171 
172  template <uint32 I>
173  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
174  uint32 alignment_pos() const { return I == 0 ? best_alignment_pos() : second_alignment_pos(); }
175 
176  template <uint32 I>
177  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
178  int32 score() const { return I == 0 ? best_score() : second_score(); }
179 
180  template <uint32 I>
181  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
182  uint32 ed() const { return I == 0? best_ed() : second_ed(); }
183 
184  template <uint32 I>
185  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
186  bool is_aligned() const { return I == 0 ? m_a1.is_aligned() : m_a2.is_aligned(); }
187 
188  template <uint32 I>
189  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
190  bool is_rc() const { return I == 0 ? m_a1.is_rc() : m_a2.is_rc(); }
191 
192  Alignment m_a1;
193  Alignment m_a2;
194 };
195 
199 struct PairedAlignments
200 {
201  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
202  PairedAlignments() {}
203 
204  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
205  PairedAlignments(const Alignment& a, const Alignment& o) : m_a( a ), m_o( o ) {}
206 
207  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE bool is_aligned() const { return m_a.is_aligned() && m_o.is_aligned(); }
208  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE bool is_paired() const { return m_a.is_paired(); }
209  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE bool is_concordant() const { return m_a.is_concordant(); }
210  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE bool is_discordant() const { return m_a.is_discordant(); }
211  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE int32 score() const { return m_a.score() + m_o.score(); }
212  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE int32 ed() const { return m_a.ed() + m_o.ed(); }
213 
214  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
215  Alignment& mate(const uint32 m) { return m == m_a.mate() ? m_a : m_o; }
216 
217  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
218  Alignment mate(const uint32 m) const { return m == m_a.mate() ? m_a : m_o; }
219 
220  Alignment m_a;
221  Alignment m_o;
222 };
223 
227 struct BestPairedAlignments
228 {
229  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
230  BestPairedAlignments() {}
231 
232  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
233  BestPairedAlignments(const BestAlignments& a, const BestAlignments& o) : m_a1( a.m_a1 ), m_a2( a.m_a2 ), m_o1( o.m_a1 ), m_o2( o.m_a2 ) {}
234 
235  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
236  BestPairedAlignments(const BestAlignments& a) : m_a1( a.m_a1 ), m_a2( a.m_a2 ), m_o1( Alignment::invalid() ), m_o2( Alignment::invalid() ) {}
237 
238  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE bool is_paired() const { return m_a1.is_paired(); }
239  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE bool is_aligned() const { return m_a1.is_aligned() /*&& m_o1.is_aligned()*/; }
240  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE bool has_second_paired() const { return m_a2.is_paired(); }
241  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE bool has_second_unpaired() const { return m_a2.is_unpaired(); }
242  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE bool has_second() const { return is_paired() ? has_second_paired() : m_a2.is_aligned(); }
243  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE int32 best_score() const { return m_a1.score() + (is_paired() ? m_o1.score() : 0); }
244  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE uint32 best_ed() const { return m_a1.ed() + (is_paired() ? m_o1.ed() : 0u); }
245  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE uint32 best_alignment_pos() const { return m_a1.alignment(); }
246  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE int32 second_score() const { return m_a2.score() + (has_second_paired() ? m_o2.score() : 0); }
247  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE uint32 second_ed() const { return m_a2.ed() + (has_second_paired() ? m_o2.ed() : 0u); }
248  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE uint32 second_alignment_pos() const { return m_a2.alignment(); }
249 
250  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE BestAlignments best_anchor() const { return BestAlignments( m_a1, m_a2 ); }
251  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE BestAlignments best_opposite() const { return BestAlignments( m_o1, m_o2 ); }
252 
253  template <uint32 I>
254  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
255  PairedAlignments pair() const { return I == 0 ? PairedAlignments( m_a1, m_o1 ) : PairedAlignments( m_a2, m_o2 ); }
256 
257  template <uint32 I>
258  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
259  const Alignment& alignment() const { return I == 0 ? m_a1 : m_a2; }
260 
261  template <uint32 I>
262  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
263  Alignment& alignment() { return I == 0 ? m_a1 : m_a2; }
264 
265  template <uint32 I>
266  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
267  const Alignment& opposite_alignment() const { return I == 0 ? m_o1 : m_o2; }
268 
269  template <uint32 I>
270  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
271  Alignment& opposite_alignment() { return I == 0 ? m_o1 : m_o2; }
272 
273  template <uint32 I>
274  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
275  uint32 alignment_pos() const { return I == 0 ? best_alignment_pos() : second_alignment_pos(); }
276 
277  template <uint32 I>
278  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
279  uint32 opposite_alignment_pos() const { return I == 0 ? m_o1.alignment() : m_o2.alignment(); }
280 
281  template <uint32 I>
282  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
283  int32 score() const { return I == 0 ? best_score() : second_score(); }
284 
285  template <uint32 I>
286  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
287  uint32 ed() const { return I == 0? best_ed() : second_ed(); }
288 
289  template <uint32 I>
290  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
291  bool is_aligned() const { return I == 0 ? m_a1.is_aligned() : m_a2.is_aligned(); }
292 
293  template <uint32 I>
294  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
295  bool is_rc() const { return I == 0 ? m_a1.is_rc() : m_a2.is_rc(); }
296 
297  template <uint32 I>
298  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
299  bool is_opposite_rc() const { return I == 0 ? m_o1.is_rc() : m_o2.is_rc(); }
300 
301  template <uint32 I>
302  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
303  uint32 anchor_mate() const { return I == 0 ? m_a1.mate() : m_a2.mate(); }
304 
305  template <uint32 I>
306  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
307  uint32 opposite_mate() const { return I == 0 ? m_o1.mate() : m_o2.mate(); }
308 
309  template <uint32 I>
310  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
311  int32 anchor_score() const { return I == 0 ? m_a1.score() : m_a2.score(); }
312 
313  template <uint32 I>
314  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
315  uint32 anchor_ed() const { return I == 0? m_a1.ed() : m_a2.ed(); }
316 
317  template <uint32 I>
318  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
319  int32 opposite_score() const { return I == 0 ? m_o1.score() : m_o2.score(); }
320 
321  template <uint32 I>
322  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
323  uint32 opposite_ed() const { return I == 0 ? m_o1.ed() : m_o2.ed(); }
324 
325  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
326  BestAlignments anchor() const { return BestAlignments( m_a1, m_a2 ); }
327 
328  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
329  BestAlignments opposite() const { return BestAlignments( m_o1, m_o2 ); }
330 
331  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
332  BestAlignments mate(const uint32 m) const
333  {
334  return BestAlignments(
335  m == m_a1.mate() ? m_a1 : m_o1,
336  m == m_a2.mate() ? m_a2 : m_o2 );
337  }
338 
339  Alignment m_a1;
340  Alignment m_a2;
341  Alignment m_o1;
342  Alignment m_o2;
343 };
344 
345 struct has_second
346 {
347  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
348  bool operator() (const BestAlignments& op) { return op.has_second(); }
349 
350  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
351  bool operator() (const BestPairedAlignments& op) { return op.has_second(); }
352 };
353 
354 struct has_second_paired
355 {
356  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
357  bool operator() (const BestAlignments& op) { return op.has_second() && op.alignment<1>().is_paired(); }
358 
359  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
360  bool operator() (const BestPairedAlignments& op) { return op.has_second_paired(); }
361 };
362 
363 struct has_second_unpaired
364 {
365  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
366  bool operator() (const BestAlignments& op) { return op.has_second() && op.alignment<1>().is_unpaired() ; }
367 
368  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
369  bool operator() (const BestPairedAlignments& op) { return op.has_second_unpaired(); }
370 };
371 
372 struct is_paired
373 {
374  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
375  bool operator() (const BestAlignments& op) { return op.alignment<0>().is_paired(); }
376 
377  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
378  bool operator() (const Alignment& op) { return op.is_paired(); }
379 };
380 struct is_unpaired
381 {
382  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
383  bool operator() (const BestAlignments& op) { return op.alignment<0>().is_unpaired(); }
384 
385  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
386  bool operator() (const Alignment& op) { return op.is_unpaired(); }
387 };
388 struct is_concordant
389 {
390  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
391  bool operator() (const BestAlignments& op) { return op.alignment<0>().is_concordant(); }
392 
393  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
394  bool operator() (const Alignment& op) { return op.is_concordant(); }
395 };
396 struct is_discordant
397 {
398  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
399  bool operator() (const BestAlignments& op) { return op.alignment<0>().is_discordant(); }
400 
401  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
402  bool operator() (const Alignment& op) { return op.is_discordant(); }
403 };
404 struct is_not_concordant // note: !concordant = discordant | unpaired - i.e. this is different from is_discordant
405 {
406  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
407  bool operator() (const BestAlignments& op) { return op.alignment<0>().is_concordant() == false; }
408 
409  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
410  bool operator() (const Alignment& op) { return op.is_concordant() == false; }
411 };
412 struct is_aligned
413 {
414  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
415  bool operator() (const BestAlignments& op) { return op.alignment<0>().is_aligned(); }
416 
417  NVBIO_HOST_DEVICE NVBIO_FORCEINLINE
418  bool operator() (const Alignment& op) { return op.is_aligned(); }
419 };
420 
421 // check whether two alignments are distinct
422 //
423 NVBIO_FORCEINLINE NVBIO_HOST_DEVICE bool distinct_alignments(
424  const uint32 pos1,
425  const bool rc1,
426  const uint32 pos2,
427  const bool rc2,
428  const uint32 dist);
429 
430 // check whether two paired-end alignments are distinct
431 //
432 NVBIO_FORCEINLINE NVBIO_HOST_DEVICE bool distinct_alignments(
433  const uint32 apos1,
434  const uint32 opos1,
435  const bool arc1,
436  const bool orc1,
437  const uint32 apos2,
438  const uint32 opos2,
439  const bool arc2,
440  const bool orc2,
441  const uint32 dist);
442 
443 // check whether two alignments are distinct
444 //
445 // NOTE: this uses Alignment::sink() which is only valid during alignment,
446 // as sink() is unioned with ed().
447 //
448 NVBIO_FORCEINLINE NVBIO_HOST_DEVICE bool distinct_alignments(
449  const Alignment& p1,
450  const Alignment& p2,
451  const uint32 dist = 1);
452 
453 // check whether two paired-end alignments are distinct
454 //
455 // NOTE: this uses Alignment::sink() which is only valid during alignment,
456 // as sink() is unioned with ed().
457 //
458 NVBIO_FORCEINLINE NVBIO_HOST_DEVICE bool distinct_alignments(
459  const PairedAlignments& p1,
460  const PairedAlignments& p2);
461 
462 // check whether two paired-end alignments are distinct
463 //
464 // NOTE: this uses Alignment::sink() which is only valid during alignment,
465 // as sink() is unioned with ed().
466 //
467 NVBIO_FORCEINLINE NVBIO_HOST_DEVICE bool distinct_alignments(
468  const PairedAlignments& p1,
469  const PairedAlignments& p2,
470  const uint32 dist);
471 
473 
474 } // namespace io
475 
476 // pod_type specializations
477 template <> struct pod_type<io::Alignment> { typedef uint2 type; };
478 template <> struct pod_type<io::BestAlignments> { typedef uint4 type; };
479 
480 } // namespace nvbio
481 
482 #include <nvbio/io/alignments_inl.h>
Generated on Wed Feb 25 2015 08:32:59 for NVBIO by   doxygen 1.8.4