NVBIO
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
batched_stream.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 #pragma once
29 
30 #include <nvbio/basic/types.h>
32 #include <nvbio/alignment/utils.h>
35 
36 namespace nvbio {
37 namespace aln {
38 
41 
53 template <typename stream_type>
55 {
56  typedef typename stream_type::aligner_type aligner_type;
60 
63  ScoreStream(const stream_type _stream, uint8* _columns, uint8* _checkpoints, const uint32 _stride) :
64  stream(_stream),
65  columns( (cell_type*)_columns ),
66  checkpoints( (cell_type*)_checkpoints ),
67  stride( _stride )
68  {}
69 
73  uint32 size() const { return stream.size(); }
74 
79  {
80  return checkpoint_type( checkpoints + i, stride );
81  }
82 
86  column_type column(const uint32 i) const
87  {
88  return column_type( columns + i, stride );
89  }
90 
93  template <typename score_unit_type>
95  void get(const uint32 job_id, score_unit_type* unit, const uint2 queue_slot) const
96  {
97  // setup the given work unit
98  unit->setup( job_id, queue_slot.x, stream );
99  }
100 
101  stream_type stream;
105 };
106 
125 template <typename stream_type, typename derived_type>
127 {
128  static const uint32 WINDOW_SIZE = 32;
129 
130  typedef typename stream_type::context_type context_type;
131  typedef typename stream_type::strings_type strings_type;
132  typedef typename stream_type::aligner_type aligner_type;
133 
134  // constructor
136  void setup(const uint32 _job_id, const uint32 _queue_slot, const stream_type& stream)
137  {
138  job_id = _job_id;
139  queue_slot = _queue_slot;
140  window_begin = 0u;
141 
142  // load the alignment context
143  valid = stream.init_context( job_id, &context );
144 
145  // load the strings to be aligned
146  if (valid)
147  {
148  const uint32 len = equal<typename aligner_type::algorithm_tag,PatternBlockingTag>() ?
149  stream.pattern_length( job_id, &context ) :
150  stream.text_length( job_id, &context );
151  stream.load_strings( job_id, 0, len, &context, &strings );
152  }
153  }
154 
155  // run method
157  bool run(const ScoreStream<stream_type>& score_stream)
158  {
159  if (valid)
160  {
161  // compute the end of the current DP matrix window
162  const uint32 len = equal<typename aligner_type::algorithm_tag,PatternBlockingTag>() ?
163  strings.pattern.length() :
164  strings.text.length();
165  const uint32 window_end = nvbio::min( window_begin + WINDOW_SIZE, len );
166 
167  valid = (static_cast<derived_type*>(this))->execute(
168  score_stream,
169  queue_slot,
170  strings.pattern,
171  strings.quals,
172  strings.text,
173  context.min_score,
174  window_begin,
175  window_end,
176  context.sink );
177 
178  if (window_end >= len)
179  valid = false;
180 
181  // update the pattern window
182  window_begin = window_end;
183  }
184 
185  // check whether we are done processing the pattern
186  if (valid == false)
187  {
188  score_stream.stream.output( job_id, &context );
189  return false;
190  }
191  return true;
192  }
193 
198  volatile uint32 window_begin;
199  volatile bool valid;
200 };
201 
205 template <typename stream_type>
208  stream_type, // pass the stream_type
209  StagedScoreUnit<stream_type> > // specify self as the derived_type
210 {
212 
213  // execute method
214  template <
215  typename pattern_string,
216  typename qual_string,
217  typename text_string,
218  typename sink_type>
220  bool execute(
221  const ScoreStream<stream_type>& score_stream,
222  const uint32 queue_slot,
223  const pattern_string pattern,
224  const qual_string quals,
225  const text_string text,
226  const int32 min_score,
227  const uint32 window_begin,
228  const uint32 window_end,
229  sink_type& sink)
230  {
231  // fetch this job's temporary storage
232  column_type column = score_stream.column( queue_slot );
233 
234  // score the current DP matrix window
235  return alignment_score(
236  score_stream.stream.aligner(),
237  pattern,
238  quals,
239  text,
240  min_score,
241  window_begin,
242  window_end,
243  sink,
244  column );
245  }
246 };
247 
251 template <uint32 BAND_LEN, typename stream_type>
254  stream_type, // pass the stream_type
255  BandedScoreUnit<BAND_LEN,stream_type> > // specify self as the derived_type
256 {
260 
261  // execute method
262  template <
263  typename pattern_string,
264  typename qual_string,
265  typename text_string,
266  typename sink_type>
268  bool execute(
269  const ScoreStream<stream_type>& score_stream,
270  const uint32 queue_slot,
271  const pattern_string pattern,
272  const qual_string quals,
273  const text_string text,
274  const int32 min_score,
275  const uint32 window_begin,
276  const uint32 window_end,
277  sink_type& sink)
278  {
279  // fetch this job's temporary storage
280  column_type column = score_stream.column( queue_slot );
281 
282  // score the current DP matrix window
283  return banded_alignment_score<BAND_LEN>(
284  score_stream.stream.aligner(),
285  pattern,
286  quals,
287  text,
288  min_score,
289  window_begin,
290  window_end,
291  sink,
292  column );
293  }
294 };
295 
297 
299 
301 
302 } // namespace alignment
303 } // namespace nvbio