NVBIO
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
score_all_inl.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 
32 
33 namespace nvbio {
34 namespace bowtie2 {
35 namespace cuda {
36 
37 namespace detail {
38 
41 
44 
47 
53 template <typename AlignerType, typename PipelineType>
54 struct AllScoreStream : public AlignmentStreamBase<SCORE_STREAM,AlignerType,PipelineType>
55 {
57  typedef typename base_type::context_type context_type;
59 
71  const uint32 _band_len,
72  const PipelineType _pipeline,
73  const AlignerType _aligner,
74  const ParamsPOD _params,
75  const uint32 _buffer_offset,
76  const uint32 _buffer_size,
77  uint32* _counter) :
78  base_type( _pipeline, _aligner, _params ), m_band_len( _band_len ),
79  m_buffer_offset( _buffer_offset ),
80  m_buffer_size( _buffer_size ),
81  m_counter( _counter ) {}
82 
86  uint32 max_pattern_length() const { return base_type::m_pipeline.reads.max_sequence_len(); }
87 
91  uint32 max_text_length() const { return base_type::m_pipeline.reads.max_sequence_len() + m_band_len; }
92 
96  uint32 size() const { return base_type::m_pipeline.hits_queue_size; }
97 
102  const uint32 i,
103  context_type* context) const
104  {
105  context->idx = base_type::m_pipeline.idx_queue[i];
106 
107  // initialize the sink
108  context->sink = aln::BestSink<int32>();
109 
110  // fetch the hit to process
111  HitReference<HitQueuesDeviceView> hit = base_type::m_pipeline.scoring_queues.hits[ context->idx ];
112 
113  // setup the read info
114  context->mate = 0;
115  context->read_rc = hit.seed.rc;
116  context->read_id = hit.read_id;
117  context->read_range = base_type::m_pipeline.reads.get_range( context->read_id );
118 
119  // setup the genome range
120  const uint32 g_pos = hit.loc;
121 
122  const uint32 read_len = context->read_range.y - context->read_range.x;
123  context->genome_begin = g_pos > m_band_len/2 ? g_pos - m_band_len/2 : 0u;
124  context->genome_end = nvbio::min( context->genome_begin + m_band_len + read_len, base_type::m_pipeline.genome_length );
125 
126  // setup the minimum score
127  context->min_score = base_type::m_pipeline.scoring_scheme.min_score( read_len );
128  return true;
129  }
130 
134  void output(
135  const uint32 i,
136  const context_type* context) const
137  {
138  const aln::BestSink<int32> sink = context->sink;
139 
140  // append all valid alignments to a list
141  if (sink.score >= context->min_score)
142  {
143  HitReference<HitQueuesDeviceView> hit = base_type::m_pipeline.scoring_queues.hits[ context->idx ];
144 
145  #if defined(NVBIO_DEVICE_COMPILATION)
146  const uint32 slot = atomicAdd( m_counter, 1u );
147  #else
148  const uint32 slot = (*m_counter)++;
149  #endif
150  base_type::m_pipeline.buffer_read_info[ (m_buffer_offset + slot) % m_buffer_size ] = hit.read_id; //(hit.read_id << 1) | hit.seed.rc;
151  base_type::m_pipeline.buffer_alignments[ (m_buffer_offset + slot) % m_buffer_size ] = io::Alignment( hit.loc, 0u, sink.score, hit.seed.rc );
152  // TODO: rewrite hit.loc
153  // hit.loc = context->genome_begin;
154  NVBIO_CUDA_DEBUG_PRINT_IF( base_type::m_params.debug.show_score( context->read_id, (sink.score >= context->min_score) ), "score: %d (rc[%u], pos[%u], [qid %u])]\n", sink.score, context->read_rc, context->genome_begin, i );
155  }
156  }
157 
162 };
163 
167 template <typename aligner_type, typename pipeline_type>
169  const uint32 band_len,
170  const pipeline_type& pipeline,
171  const aligner_type aligner,
172  const ParamsPOD params,
173  const uint32 buffer_offset,
174  const uint32 buffer_size,
175  uint32* counter)
176 {
177  const uint32 static_band_len =
178  (band_len < 4) ? 3u :
179  (band_len < 8) ? 7u :
180  (band_len < 16) ? 15u :
181  31u;
182 
183  typedef AllScoreStream<aligner_type,pipeline_type> stream_type;
184 
185  stream_type stream(
186  static_band_len,
187  pipeline,
188  aligner,
189  params,
190  buffer_offset,
191  buffer_size,
192  counter );
193 
194  //typedef aln::DeviceThreadScheduler scheduler_type;
195  typedef aln::DeviceStagedThreadScheduler scheduler_type;
196 
197  if (band_len < 4)
198  {
200 
201  batch.enact( stream, pipeline.dp_buffer_size, pipeline.dp_buffer );
202  }
203  else if (band_len < 8)
204  {
206 
207  batch.enact( stream, pipeline.dp_buffer_size, pipeline.dp_buffer );
208  }
209  else if (band_len < 16)
210  {
212 
213  batch.enact( stream, pipeline.dp_buffer_size, pipeline.dp_buffer );
214  }
215  else
216  {
218 
219  batch.enact( stream, pipeline.dp_buffer_size, pipeline.dp_buffer );
220  }
221 }
222 
226 
227 } // namespace detail
228 
229 //
230 // execute a batch of single-ended banded-alignment score calculations, best mapping
231 //
232 // \param band_len alignment band length
233 // \param pipeline all mapping pipeline
234 // \param params alignment params
235 // \param buffer_offset ring buffer offset
236 // \param buffer_size ring buffer size
237 // \return number of valid alignments
238 template <typename scheme_type>
240  const uint32 band_len,
241  const AllMappingPipelineState<scheme_type>& pipeline,
242  const ParamsPOD params,
243  const uint32 ring_offset,
244  const uint32 ring_size)
245 {
246  thrust::device_vector<uint32> counter(1,0);
247 
248  if (params.alignment_type == LocalAlignment)
249  {
251  band_len,
252  pipeline,
253  pipeline.scoring_scheme.local_aligner(),
254  params,
255  ring_offset,
256  ring_size,
257  nvbio::device_view(counter) );
258  }
259  else
260  {
262  band_len,
263  pipeline,
264  pipeline.scoring_scheme.end_to_end_aligner(),
265  params,
266  ring_offset,
267  ring_size,
268  nvbio::device_view(counter) );
269  }
270  //cudaDeviceSynchronize();
271  //nvbio::cuda::check_error("score kernel");
272  return counter[0];
273 }
274 
275 } // namespace cuda
276 } // namespace bowtie2
277 } // namespace nvbio