1 // ==========================================================================
2 // SeqAn - The Library for Sequence Analysis
3 // ==========================================================================
4 // Copyright (c) 2006-2015, Knut Reinert, FU Berlin
5 // All rights reserved.
6 //
7 // Redistribution and use in source and binary forms, with or without
8 // modification, are permitted provided that the following conditions are met:
9 //
10 // * Redistributions of source code must retain the above copyright
11 // notice, this list of conditions and the following disclaimer.
12 // * Redistributions in binary form must reproduce the above copyright
13 // notice, this list of conditions and the following disclaimer in the
14 // documentation and/or other materials provided with the distribution.
15 // * Neither the name of Knut Reinert or the FU Berlin nor the names of
16 // its contributors may be used to endorse or promote products derived
17 // from this software without specific prior written permission.
18 //
19 // THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
20 // AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
21 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
22 // ARE DISCLAIMED. IN NO EVENT SHALL KNUT REINERT OR THE FU BERLIN BE LIABLE
23 // FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
24 // DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
25 // SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
26 // CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
27 // LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
28 // OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH
29 // DAMAGE.
30 //
31 // ==========================================================================
32 // Author: Rene Rahn <rene.rahn@fu-berlin.de>
33 // ==========================================================================
34 // Implements the affine gap cost functions.
35 // ==========================================================================
36
37 #ifndef SEQAN_INCLUDE_SEQAN_ALIGN_DP_FORMULA_AFFINE_H_
38 #define SEQAN_INCLUDE_SEQAN_ALIGN_DP_FORMULA_AFFINE_H_
39
40 namespace seqan {
41
42 // ============================================================================
43 // Forwards
44 // ============================================================================
45
46 // ============================================================================
47 // Tags, Classes, Enums
48 // ============================================================================
49
50 // ============================================================================
51 // Metafunctions
52 // ============================================================================
53
54 // ============================================================================
55 // Functions
56 // ============================================================================
57
58 // ----------------------------------------------------------------------------
59 // Function _retrieveTraceAffine [RecursionDirectionAll]
60 // ----------------------------------------------------------------------------
61
62 template <typename TScoreValue>
63 inline TraceBitMap_::TTraceValue
_retrieveTraceAffine(TScoreValue const & globalMax,TScoreValue const & diagScore,TScoreValue const & horiScore,TScoreValue const & horiOpenScore,TScoreValue const & vertiScore,TScoreValue const & vertiOpenScore,RecursionDirectionAll const &)64 _retrieveTraceAffine(TScoreValue const & globalMax,
65 TScoreValue const & diagScore,
66 TScoreValue const & horiScore,
67 TScoreValue const & horiOpenScore,
68 TScoreValue const & vertiScore,
69 TScoreValue const & vertiOpenScore,
70 RecursionDirectionAll const &)
71 {
72 typename TraceBitMap_::TTraceValue traceValue(TraceBitMap_::NONE);
73 _conditionalOrOnEquality(traceValue, globalMax, diagScore, TraceBitMap_::DIAGONAL);
74 _conditionalOrOnEquality(traceValue, globalMax, horiScore, TraceBitMap_::HORIZONTAL);
75 _conditionalOrOnEquality(traceValue, globalMax, horiOpenScore, TraceBitMap_::MAX_FROM_HORIZONTAL_MATRIX);
76 _conditionalOrOnEquality(traceValue, globalMax, vertiScore, TraceBitMap_::VERTICAL);
77 _conditionalOrOnEquality(traceValue, globalMax, vertiOpenScore, TraceBitMap_::MAX_FROM_VERTICAL_MATRIX);
78 return traceValue;
79 }
80
81 // ----------------------------------------------------------------------------
82 // Function _retrieveTraceAffine [RecursionDirectionUpperDiagonal]
83 // ----------------------------------------------------------------------------
84
85 template <typename TScoreValue>
86 inline TraceBitMap_::TTraceValue
_retrieveTraceAffine(TScoreValue const & globalMax,TScoreValue const & diagScore,TScoreValue const & horiScore,TScoreValue const & horiOpenScore,TScoreValue const &,TScoreValue const &,RecursionDirectionUpperDiagonal const &)87 _retrieveTraceAffine(TScoreValue const & globalMax,
88 TScoreValue const & diagScore,
89 TScoreValue const & horiScore,
90 TScoreValue const & horiOpenScore,
91 TScoreValue const &,
92 TScoreValue const &,
93 RecursionDirectionUpperDiagonal const &)
94 {
95 typename TraceBitMap_::TTraceValue traceValue(TraceBitMap_::NONE);
96 _conditionalOrOnEquality(traceValue, globalMax, diagScore, TraceBitMap_::DIAGONAL);
97 _conditionalOrOnEquality(traceValue, globalMax, horiScore, TraceBitMap_::HORIZONTAL);
98 _conditionalOrOnEquality(traceValue, globalMax, horiOpenScore, TraceBitMap_::MAX_FROM_HORIZONTAL_MATRIX);
99 return traceValue;
100 }
101
102 // ----------------------------------------------------------------------------
103 // Function _retrieveTraceAffine [RecursionDirectionLowerDiagonal]
104 // ----------------------------------------------------------------------------
105
106 template <typename TScoreValue>
107 inline TraceBitMap_::TTraceValue
_retrieveTraceAffine(TScoreValue const & globalMax,TScoreValue const & diagScore,TScoreValue const &,TScoreValue const &,TScoreValue const & vertiScore,TScoreValue const & vertiOpenScore,RecursionDirectionLowerDiagonal const &)108 _retrieveTraceAffine(TScoreValue const & globalMax,
109 TScoreValue const & diagScore,
110 TScoreValue const &,
111 TScoreValue const &,
112 TScoreValue const & vertiScore,
113 TScoreValue const & vertiOpenScore,
114 RecursionDirectionLowerDiagonal const &)
115 {
116 typename TraceBitMap_::TTraceValue traceValue(TraceBitMap_::NONE);
117 _conditionalOrOnEquality(traceValue, globalMax, diagScore, TraceBitMap_::DIAGONAL);
118 _conditionalOrOnEquality(traceValue, globalMax, vertiScore, TraceBitMap_::VERTICAL);
119 _conditionalOrOnEquality(traceValue, globalMax, vertiOpenScore, TraceBitMap_::MAX_FROM_VERTICAL_MATRIX);
120 return traceValue;
121 }
122
123 // ----------------------------------------------------------------------------
124 // Function _retrieveTraceAffine [RecursionDirectionHorizontal]
125 // ----------------------------------------------------------------------------
126
127 template <typename TScoreValue>
128 inline TraceBitMap_::TTraceValue
_retrieveTraceAffine(TScoreValue const & globalMax,TScoreValue const &,TScoreValue const & horiScore,TScoreValue const & horiOpenScore,TScoreValue const &,TScoreValue const &,RecursionDirectionHorizontal const &)129 _retrieveTraceAffine(TScoreValue const & globalMax,
130 TScoreValue const &,
131 TScoreValue const & horiScore,
132 TScoreValue const & horiOpenScore,
133 TScoreValue const &,
134 TScoreValue const &,
135 RecursionDirectionHorizontal const &)
136 {
137 typename TraceBitMap_::TTraceValue traceValue(TraceBitMap_::NONE);
138 _conditionalOrOnEquality(traceValue, globalMax, horiScore, TraceBitMap_::HORIZONTAL);
139 _conditionalOrOnEquality(traceValue, globalMax, horiOpenScore, TraceBitMap_::MAX_FROM_HORIZONTAL_MATRIX);
140 return traceValue;
141 }
142
143 // ----------------------------------------------------------------------------
144 // Function _retrieveTraceAffine [RecursionDirectionVertical]
145 // ----------------------------------------------------------------------------
146
147 template <typename TScoreValue>
148 inline TraceBitMap_::TTraceValue
_retrieveTraceAffine(TScoreValue const & globalMax,TScoreValue const &,TScoreValue const &,TScoreValue const &,TScoreValue const & vertiScore,TScoreValue const & vertiOpenScore,RecursionDirectionVertical const &)149 _retrieveTraceAffine(TScoreValue const & globalMax,
150 TScoreValue const &,
151 TScoreValue const &,
152 TScoreValue const &,
153 TScoreValue const & vertiScore,
154 TScoreValue const & vertiOpenScore,
155 RecursionDirectionVertical const &)
156 {
157 typename TraceBitMap_::TTraceValue traceValue(TraceBitMap_::NONE);
158 _conditionalOrOnEquality(traceValue, globalMax, vertiScore, TraceBitMap_::VERTICAL);
159 _conditionalOrOnEquality(traceValue, globalMax, vertiOpenScore, TraceBitMap_::MAX_FROM_VERTICAL_MATRIX);
160 return traceValue;
161 }
162
163 // ----------------------------------------------------------------------------
164 // Function _internalComputeScore [RecursionDirectionDiagonal, AffineGaps]
165 // ----------------------------------------------------------------------------
166
167 template <typename TScoreValue, typename TAffineGaps, typename TTraceValueL, typename TTraceValueGap>
168 inline typename TraceBitMap_::TTraceValue
_internalComputeScore(DPCell_<TScoreValue,TAffineGaps> & activeCell,TScoreValue const & rightCompare,TTraceValueL,TTraceValueGap,TracebackOff const &,RecursionDirectionDiagonal const &)169 _internalComputeScore(DPCell_<TScoreValue, TAffineGaps> & activeCell,
170 TScoreValue const & rightCompare,
171 TTraceValueL,
172 TTraceValueGap,
173 TracebackOff const &,
174 RecursionDirectionDiagonal const &)
175 {
176 if(activeCell._score < rightCompare)
177 activeCell._score = rightCompare;
178 return TraceBitMap_::NONE;
179 }
180
181 template <typename TScoreValue, typename TAffineGaps, typename TTraceValueL, typename TTraceValueGap, typename TGapsPlacement>
182 inline typename TraceBitMap_::TTraceValue
_internalComputeScore(DPCell_<TScoreValue,TAffineGaps> & activeCell,TScoreValue const & rightCompare,TTraceValueL leftTrace,TTraceValueGap gapTrace,TracebackOn<TracebackConfig_<SingleTrace,TGapsPlacement>> const &,RecursionDirectionDiagonal const &)183 _internalComputeScore(DPCell_<TScoreValue, TAffineGaps> & activeCell,
184 TScoreValue const & rightCompare,
185 TTraceValueL leftTrace,
186 TTraceValueGap gapTrace,
187 TracebackOn<TracebackConfig_<SingleTrace, TGapsPlacement> > const &,
188 RecursionDirectionDiagonal const &)
189 {
190 if(activeCell._score <= rightCompare)
191 {
192 activeCell._score = rightCompare;
193 return TraceBitMap_::DIAGONAL | leftTrace;
194 }
195 return leftTrace | gapTrace;
196 }
197
198 template <typename TScoreValue, typename TAffineGaps, typename TTraceValueL, typename TTraceValueGap, typename TGapsPlacement>
199 inline typename TraceBitMap_::TTraceValue
_internalComputeScore(DPCell_<TScoreValue,TAffineGaps> & activeCell,TScoreValue const & rightCompare,TTraceValueL leftTrace,TTraceValueGap gapTrace,TracebackOn<TracebackConfig_<CompleteTrace,TGapsPlacement>> const &,RecursionDirectionDiagonal const &)200 _internalComputeScore(DPCell_<TScoreValue, TAffineGaps> & activeCell,
201 TScoreValue const & rightCompare,
202 TTraceValueL leftTrace,
203 TTraceValueGap gapTrace,
204 TracebackOn<TracebackConfig_<CompleteTrace, TGapsPlacement> > const &,
205 RecursionDirectionDiagonal const &)
206 {
207 if(activeCell._score < rightCompare)
208 {
209 activeCell._score = rightCompare; // Maximal score comes from diagonal.
210 return TraceBitMap_::DIAGONAL | leftTrace; // Return trace for Diagonal.
211 }
212 if (activeCell._score == rightCompare) // Maximal score comes from previous computed directions and diagonal.
213 return leftTrace | TraceBitMap_::DIAGONAL | gapTrace; // Return all directions inclusively the flag indicating max from gap.
214 return leftTrace | gapTrace; // Maximum comes from gap. Return gap value inclusively the flag indicating max from gap.
215 }
216
217 // ----------------------------------------------------------------------------
218 // Function _internalComputeScore [RecursionDirectionHorizontal, AffineGaps]
219 // ----------------------------------------------------------------------------
220
221 template <typename TScoreValue, typename TTraceValueL, typename TTraceValueR>
222 inline typename TraceBitMap_::TTraceValue
_internalComputeScore(DPCell_<TScoreValue,AffineGaps> & activeCell,TScoreValue const & rightCompare,TTraceValueL,TTraceValueR,TracebackOff const &,RecursionDirectionHorizontal const &)223 _internalComputeScore(DPCell_<TScoreValue, AffineGaps> & activeCell,
224 TScoreValue const & rightCompare,
225 TTraceValueL,
226 TTraceValueR,
227 TracebackOff const &,
228 RecursionDirectionHorizontal const &)
229 {
230 if(activeCell._horizontalScore < rightCompare)
231 activeCell._score = activeCell._horizontalScore = rightCompare;
232 else
233 activeCell._score = activeCell._horizontalScore;
234 return TraceBitMap_::NONE;
235 }
236
237 template <typename TScoreValue, typename TTraceValueL, typename TTraceValueR, typename TGapsPlacement>
238 inline typename TraceBitMap_::TTraceValue
_internalComputeScore(DPCell_<TScoreValue,AffineGaps> & activeCell,TScoreValue const & rightCompare,TTraceValueL leftTrace,TTraceValueR rightTrace,TracebackOn<TracebackConfig_<SingleTrace,TGapsPlacement>> const &,RecursionDirectionHorizontal const &)239 _internalComputeScore(DPCell_<TScoreValue, AffineGaps> & activeCell,
240 TScoreValue const & rightCompare,
241 TTraceValueL leftTrace,
242 TTraceValueR rightTrace,
243 TracebackOn<TracebackConfig_<SingleTrace, TGapsPlacement> > const &,
244 RecursionDirectionHorizontal const &)
245 {
246 if(activeCell._horizontalScore < rightCompare)
247 {
248 activeCell._score = activeCell._horizontalScore = rightCompare;
249 return rightTrace;
250 }
251 activeCell._score = activeCell._horizontalScore;
252 return leftTrace;
253 }
254
255 template <typename TScoreValue, typename TTraceValueL, typename TTraceValueR, typename TGapsPlacement>
256 inline typename TraceBitMap_::TTraceValue
_internalComputeScore(DPCell_<TScoreValue,AffineGaps> & activeCell,TScoreValue const & rightCompare,TTraceValueL leftTrace,TTraceValueR rightTrace,TracebackOn<TracebackConfig_<CompleteTrace,TGapsPlacement>> const &,RecursionDirectionHorizontal const &)257 _internalComputeScore(DPCell_<TScoreValue, AffineGaps> & activeCell,
258 TScoreValue const & rightCompare,
259 TTraceValueL leftTrace,
260 TTraceValueR rightTrace,
261 TracebackOn<TracebackConfig_<CompleteTrace, TGapsPlacement> > const &,
262 RecursionDirectionHorizontal const &)
263 {
264 if(activeCell._horizontalScore < rightCompare)
265 {
266 activeCell._score = activeCell._horizontalScore = rightCompare;
267 return rightTrace;
268 }
269 activeCell._score = activeCell._horizontalScore;
270 if (activeCell._horizontalScore == rightCompare)
271 return leftTrace | rightTrace;
272 return leftTrace;
273 }
274
275 // ----------------------------------------------------------------------------
276 // Function _internalComputeScore [RecursionDirectionVertical, AffineGaps]
277 // ----------------------------------------------------------------------------
278
279 template <typename TScoreValue, typename TTraceValueL, typename TTraceValueR>
280 inline typename TraceBitMap_::TTraceValue
_internalComputeScore(DPCell_<TScoreValue,AffineGaps> & activeCell,TScoreValue const & rightCompare,TTraceValueL,TTraceValueR,TracebackOff const &,RecursionDirectionVertical const &)281 _internalComputeScore(DPCell_<TScoreValue, AffineGaps> & activeCell,
282 TScoreValue const & rightCompare,
283 TTraceValueL,
284 TTraceValueR,
285 TracebackOff const &,
286 RecursionDirectionVertical const &)
287 {
288 if(activeCell._verticalScore < rightCompare)
289 activeCell._score = activeCell._verticalScore = rightCompare;
290 else
291 activeCell._score = activeCell._verticalScore;
292 return TraceBitMap_::NONE;
293 }
294
295 template <typename TScoreValue, typename TTraceValueL, typename TTraceValueR, typename TGapsPlacement>
296 inline typename TraceBitMap_::TTraceValue
_internalComputeScore(DPCell_<TScoreValue,AffineGaps> & activeCell,TScoreValue const & rightCompare,TTraceValueL leftTrace,TTraceValueR rightTrace,TracebackOn<TracebackConfig_<SingleTrace,TGapsPlacement>> const &,RecursionDirectionVertical const &)297 _internalComputeScore(DPCell_<TScoreValue, AffineGaps> & activeCell,
298 TScoreValue const & rightCompare,
299 TTraceValueL leftTrace,
300 TTraceValueR rightTrace,
301 TracebackOn<TracebackConfig_<SingleTrace, TGapsPlacement> > const &,
302 RecursionDirectionVertical const &)
303 {
304 if(activeCell._verticalScore < rightCompare)
305 {
306 activeCell._score = activeCell._verticalScore = rightCompare;
307 return rightTrace;
308 }
309 activeCell._score = activeCell._verticalScore;
310 return leftTrace;
311 }
312
313 template <typename TScoreValue, typename TTraceValueL, typename TTraceValueR, typename TGapsPlacement>
314 inline typename TraceBitMap_::TTraceValue
_internalComputeScore(DPCell_<TScoreValue,AffineGaps> & activeCell,TScoreValue const & rightCompare,TTraceValueL leftTrace,TTraceValueR rightTrace,TracebackOn<TracebackConfig_<CompleteTrace,TGapsPlacement>> const &,RecursionDirectionVertical const &)315 _internalComputeScore(DPCell_<TScoreValue, AffineGaps> & activeCell,
316 TScoreValue const & rightCompare,
317 TTraceValueL leftTrace,
318 TTraceValueR rightTrace,
319 TracebackOn<TracebackConfig_<CompleteTrace, TGapsPlacement> > const &,
320 RecursionDirectionVertical const &)
321 {
322 if(activeCell._verticalScore < rightCompare)
323 {
324 activeCell._score = activeCell._verticalScore = rightCompare;
325 return rightTrace;
326 }
327 activeCell._score = activeCell._verticalScore;
328 if (activeCell._verticalScore == rightCompare)
329 return leftTrace | rightTrace;
330 return leftTrace;
331 }
332
333 // ----------------------------------------------------------------------------
334 // Function _internalComputeScore [Vertical vs Horizontal, AffineGaps]
335 // ----------------------------------------------------------------------------
336
337 template <typename TScoreValue>
338 inline typename TraceBitMap_::TTraceValue
_internalComputeScore(DPCell_<TScoreValue,AffineGaps> & activeCell,TracebackOff const &)339 _internalComputeScore(DPCell_<TScoreValue, AffineGaps> & activeCell,
340 TracebackOff const &)
341 {
342 if(activeCell._score < activeCell._horizontalScore)
343 activeCell._score = activeCell._horizontalScore;
344 return TraceBitMap_::NONE;
345 }
346
347 template <typename TScoreValue, typename TGapsPlacement>
348 inline typename TraceBitMap_::TTraceValue
_internalComputeScore(DPCell_<TScoreValue,AffineGaps> & activeCell,TracebackOn<TracebackConfig_<SingleTrace,TGapsPlacement>> const &)349 _internalComputeScore(DPCell_<TScoreValue, AffineGaps> & activeCell,
350 TracebackOn<TracebackConfig_<SingleTrace, TGapsPlacement> > const &)
351 {
352 if(activeCell._score < activeCell._horizontalScore)
353 {
354 activeCell._score = activeCell._horizontalScore;
355 return TraceBitMap_::MAX_FROM_HORIZONTAL_MATRIX;
356 }
357 return TraceBitMap_::MAX_FROM_VERTICAL_MATRIX;
358 }
359
360 template <typename TScoreValue, typename TGapsPlacement>
361 inline typename TraceBitMap_::TTraceValue
_internalComputeScore(DPCell_<TScoreValue,AffineGaps> & activeCell,TracebackOn<TracebackConfig_<CompleteTrace,TGapsPlacement>> const &)362 _internalComputeScore(DPCell_<TScoreValue, AffineGaps> & activeCell,
363 TracebackOn<TracebackConfig_<CompleteTrace, TGapsPlacement> > const &)
364 {
365 if(activeCell._score < activeCell._horizontalScore)
366 {
367 activeCell._score = activeCell._horizontalScore;
368 return TraceBitMap_::MAX_FROM_HORIZONTAL_MATRIX;
369 }
370 if (activeCell._score == activeCell._horizontalScore)
371 return TraceBitMap_::MAX_FROM_VERTICAL_MATRIX | TraceBitMap_::MAX_FROM_HORIZONTAL_MATRIX;
372 return TraceBitMap_::MAX_FROM_VERTICAL_MATRIX;
373 }
374
375 // ----------------------------------------------------------------------------
376 // Function _doComputeScore [RecursionAllDirection, AffineGaps]
377 // ----------------------------------------------------------------------------
378
379 template <typename TScoreValue, typename TSequenceHValue, typename TSequenceVValue, typename TScoringScheme,
380 typename TAlgorithm, typename TTracebackConfig>
381 inline typename TraceBitMap_::TTraceValue
_doComputeScore(DPCell_<TScoreValue,AffineGaps> & activeCell,DPCell_<TScoreValue,AffineGaps> const & previousDiagonal,DPCell_<TScoreValue,AffineGaps> const & previousHorizontal,DPCell_<TScoreValue,AffineGaps> const & previousVertical,TSequenceHValue const & seqHVal,TSequenceVValue const & seqVVal,TScoringScheme const & scoringScheme,RecursionDirectionAll const &,DPProfile_<TAlgorithm,AffineGaps,TTracebackConfig> const &)382 _doComputeScore(DPCell_<TScoreValue, AffineGaps> & activeCell,
383 DPCell_<TScoreValue, AffineGaps> const & previousDiagonal,
384 DPCell_<TScoreValue, AffineGaps> const & previousHorizontal,
385 DPCell_<TScoreValue, AffineGaps> const & previousVertical,
386 TSequenceHValue const & seqHVal,
387 TSequenceVValue const & seqVVal,
388 TScoringScheme const & scoringScheme,
389 RecursionDirectionAll const &,
390 DPProfile_<TAlgorithm, AffineGaps, TTracebackConfig> const &)
391 {
392 typedef typename TraceBitMap_::TTraceValue TTraceValue;
393
394 // Now we have to find a smart version to solve this problem. Which is not as easy I would think.
395
396 activeCell._horizontalScore = _horizontalScoreOfCell(previousHorizontal) +
397 scoreGapExtendHorizontal(scoringScheme, seqHVal, seqVVal);
398 TScoreValue tmpScore = _scoreOfCell(previousHorizontal) + scoreGapOpenHorizontal(scoringScheme, seqHVal, seqVVal);
399
400 TTraceValue tvGap = _internalComputeScore(activeCell, tmpScore, TraceBitMap_::HORIZONTAL, TraceBitMap_::HORIZONTAL_OPEN, TTracebackConfig(), RecursionDirectionHorizontal());
401
402 // Now we can decide for the optimal score in horizontal score or not?
403 activeCell._verticalScore = _verticalScoreOfCell(previousVertical) + scoreGapExtendVertical(scoringScheme, seqHVal, seqVVal);
404 tmpScore = _scoreOfCell(previousVertical) + scoreGapOpenVertical(scoringScheme, seqHVal, seqVVal);
405 tvGap |= _internalComputeScore(activeCell, tmpScore, TraceBitMap_::VERTICAL, TraceBitMap_::VERTICAL_OPEN, TTracebackConfig(), RecursionDirectionVertical());
406
407 // Finds the maximum between the vertical and the horizontal matrix. Stores the flag for coming from a potential direction.
408 TTraceValue tvMax = _internalComputeScore(activeCell, TTracebackConfig()); // Stores from where the maximal score comes.
409 tmpScore = _scoreOfCell(previousDiagonal) + score(scoringScheme, seqHVal, seqVVal);
410 return _internalComputeScore(activeCell, tmpScore, tvGap, tvMax, TTracebackConfig(), RecursionDirectionDiagonal());
411 }
412
413 // ----------------------------------------------------------------------------
414 // Function _doComputeScore [RecursionUpperDiagonalDirection, AffineGaps]
415 // ----------------------------------------------------------------------------
416
417 template <typename TScoreValue, typename TSequenceHValue, typename TSequenceVValue, typename TScoringScheme,
418 typename TAlgorithm, typename TTracebackConfig>
419 inline typename TraceBitMap_::TTraceValue
_doComputeScore(DPCell_<TScoreValue,AffineGaps> & activeCell,DPCell_<TScoreValue,AffineGaps> const & previousDiagonal,DPCell_<TScoreValue,AffineGaps> const & previousHorizontal,DPCell_<TScoreValue,AffineGaps> const &,TSequenceHValue const & seqHVal,TSequenceVValue const & seqVVal,TScoringScheme const & scoringScheme,RecursionDirectionUpperDiagonal const &,DPProfile_<TAlgorithm,AffineGaps,TTracebackConfig> const &)420 _doComputeScore(DPCell_<TScoreValue, AffineGaps> & activeCell,
421 DPCell_<TScoreValue, AffineGaps> const & previousDiagonal,
422 DPCell_<TScoreValue, AffineGaps> const & previousHorizontal,
423 DPCell_<TScoreValue, AffineGaps> const & /*previousVertical*/,
424 TSequenceHValue const & seqHVal,
425 TSequenceVValue const & seqVVal,
426 TScoringScheme const & scoringScheme,
427 RecursionDirectionUpperDiagonal const &,
428 DPProfile_<TAlgorithm, AffineGaps, TTracebackConfig> const &)
429 {
430 typedef typename TraceBitMap_::TTraceValue TTraceValue;
431 activeCell._horizontalScore = _horizontalScoreOfCell(previousHorizontal)
432 + scoreGapExtendHorizontal(scoringScheme, seqHVal, seqVVal);
433
434 activeCell._verticalScore = DPCellDefaultInfinity<DPCell_<TScoreValue, AffineGaps> >::VALUE;
435 TScoreValue tmpScore = _scoreOfCell(previousHorizontal) + scoreGapOpenHorizontal(scoringScheme, seqHVal, seqVVal);
436 TTraceValue tv = _internalComputeScore(activeCell, tmpScore, TraceBitMap_::HORIZONTAL, TraceBitMap_::HORIZONTAL_OPEN, TTracebackConfig(), RecursionDirectionHorizontal());
437 tmpScore = _scoreOfCell(previousDiagonal) + score(scoringScheme, seqHVal, seqVVal);
438 return _internalComputeScore(activeCell, tmpScore, tv, TraceBitMap_::MAX_FROM_HORIZONTAL_MATRIX, TTracebackConfig(), RecursionDirectionDiagonal());
439 }
440
441 // ----------------------------------------------------------------------------
442 // Function _doComputeScore [RecursionDirectionLowerDiagonal, AffineGaps]
443 // ----------------------------------------------------------------------------
444
445 template <typename TScoreValue, typename TSequenceHValue, typename TSequenceVValue, typename TScoringScheme,
446 typename TAlgorithm, typename TTracebackConfig>
447 inline typename TraceBitMap_::TTraceValue
_doComputeScore(DPCell_<TScoreValue,AffineGaps> & activeCell,DPCell_<TScoreValue,AffineGaps> const & previousDiagonal,DPCell_<TScoreValue,AffineGaps> const &,DPCell_<TScoreValue,AffineGaps> const & previousVertical,TSequenceHValue const & seqHVal,TSequenceVValue const & seqVVal,TScoringScheme const & scoringScheme,RecursionDirectionLowerDiagonal const &,DPProfile_<TAlgorithm,AffineGaps,TTracebackConfig> const &)448 _doComputeScore(DPCell_<TScoreValue, AffineGaps> & activeCell,
449 DPCell_<TScoreValue, AffineGaps> const & previousDiagonal,
450 DPCell_<TScoreValue, AffineGaps> const & /*previousHorizontal*/,
451 DPCell_<TScoreValue, AffineGaps> const & previousVertical,
452 TSequenceHValue const & seqHVal,
453 TSequenceVValue const & seqVVal,
454 TScoringScheme const & scoringScheme,
455 RecursionDirectionLowerDiagonal const &,
456 DPProfile_<TAlgorithm, AffineGaps, TTracebackConfig> const &)
457 {
458 typedef typename TraceBitMap_::TTraceValue TTraceValue;
459
460 activeCell._verticalScore = _verticalScoreOfCell(previousVertical) +
461 scoreGapExtendVertical(scoringScheme, seqHVal, seqVVal);
462 TScoreValue tmpScore = _scoreOfCell(previousVertical) + scoreGapOpenVertical(scoringScheme, seqHVal, seqVVal);
463
464 activeCell._horizontalScore = DPCellDefaultInfinity<DPCell_<TScoreValue, AffineGaps> >::VALUE;
465 // This computes the difference between the vertical extend and vertical open.
466 TTraceValue tv = _internalComputeScore(activeCell, tmpScore, TraceBitMap_::VERTICAL, TraceBitMap_::VERTICAL_OPEN, TTracebackConfig(), RecursionDirectionVertical());
467
468 // Up to here, activeCell stores the highest value of vertical or vertical open.
469 tmpScore = _scoreOfCell(previousDiagonal) + score(scoringScheme, seqHVal, seqVVal);
470 return _internalComputeScore(activeCell, tmpScore, tv, TraceBitMap_::MAX_FROM_VERTICAL_MATRIX, TTracebackConfig(), RecursionDirectionDiagonal()); // Now we have this problem. How do we determine if the max comes from the vertical distance.
471 }
472
473 // ----------------------------------------------------------------------------
474 // Function _doComputeScore [RecursionHorizontalDirection]
475 // ----------------------------------------------------------------------------
476
477 template <typename TScoreValue, typename TSequenceHValue, typename TSequenceVValue, typename TScoringScheme,
478 typename TAlgorithm, typename TTracebackConfig>
479 inline typename TraceBitMap_::TTraceValue
_doComputeScore(DPCell_<TScoreValue,AffineGaps> & activeCell,DPCell_<TScoreValue,AffineGaps> const &,DPCell_<TScoreValue,AffineGaps> const & previousHorizontal,DPCell_<TScoreValue,AffineGaps> const &,TSequenceHValue const & seqHVal,TSequenceVValue const & seqVVal,TScoringScheme const & scoringScheme,RecursionDirectionHorizontal const &,DPProfile_<TAlgorithm,AffineGaps,TTracebackConfig> const &)480 _doComputeScore(DPCell_<TScoreValue, AffineGaps> & activeCell,
481 DPCell_<TScoreValue, AffineGaps> const & /*previousDiagonal*/,
482 DPCell_<TScoreValue, AffineGaps> const & previousHorizontal,
483 DPCell_<TScoreValue, AffineGaps> const & /*previousVertical*/,
484 TSequenceHValue const & seqHVal,
485 TSequenceVValue const & seqVVal,
486 TScoringScheme const & scoringScheme,
487 RecursionDirectionHorizontal const &,
488 DPProfile_<TAlgorithm, AffineGaps, TTracebackConfig> const &)
489 {
490 //typedef typename TraceBitMap_::TTraceValue TTraceValue;
491 TScoreValue tmpGapOpenHorizontal = _scoreOfCell(previousHorizontal) + scoreGapOpenHorizontal(scoringScheme, seqHVal, seqVVal);
492 activeCell._horizontalScore = _horizontalScoreOfCell(previousHorizontal) + scoreGapExtendHorizontal(scoringScheme, seqHVal, seqVVal);
493
494 activeCell._verticalScore = DPCellDefaultInfinity<DPCell_<TScoreValue, AffineGaps> >::VALUE;
495 return _internalComputeScore(activeCell, tmpGapOpenHorizontal, TraceBitMap_::HORIZONTAL, TraceBitMap_::HORIZONTAL_OPEN, TTracebackConfig(), RecursionDirectionHorizontal()) | TraceBitMap_::MAX_FROM_HORIZONTAL_MATRIX;
496 }
497
498 // ----------------------------------------------------------------------------
499 // Function _doComputeScore [RecursionVerticalDirection]
500 // ----------------------------------------------------------------------------
501
502 template <typename TScoreValue, typename TSequenceHValue, typename TSequenceVValue, typename TScoringScheme,
503 typename TAlgorithm, typename TTracebackConfig>
504 inline typename TraceBitMap_::TTraceValue
_doComputeScore(DPCell_<TScoreValue,AffineGaps> & activeCell,DPCell_<TScoreValue,AffineGaps> const &,DPCell_<TScoreValue,AffineGaps> const &,DPCell_<TScoreValue,AffineGaps> const & previousVertical,TSequenceHValue const & seqHVal,TSequenceVValue const & seqVVal,TScoringScheme const & scoringScheme,RecursionDirectionVertical const &,DPProfile_<TAlgorithm,AffineGaps,TTracebackConfig> const &)505 _doComputeScore(DPCell_<TScoreValue, AffineGaps> & activeCell,
506 DPCell_<TScoreValue, AffineGaps> const & /*previousDiagonal*/,
507 DPCell_<TScoreValue, AffineGaps> const & /*previousHorizontal*/,
508 DPCell_<TScoreValue, AffineGaps> const & previousVertical,
509 TSequenceHValue const & seqHVal,
510 TSequenceVValue const & seqVVal,
511 TScoringScheme const & scoringScheme,
512 RecursionDirectionVertical const &,
513 DPProfile_<TAlgorithm, AffineGaps, TTracebackConfig> const &)
514 {
515 //typedef typename TraceBitMap_::TTraceValue TTraceValue;
516 TScoreValue tmpGapOpenVertical = _scoreOfCell(previousVertical) + scoreGapOpenVertical(scoringScheme, seqHVal, seqVVal);
517 activeCell._verticalScore = _verticalScoreOfCell(previousVertical) + scoreGapExtendVertical(scoringScheme, seqHVal, seqVVal);
518
519 // Here we distinguish between vertical and vertical open.
520 activeCell._horizontalScore = DPCellDefaultInfinity<DPCell_<TScoreValue, AffineGaps> >::VALUE;
521 return _internalComputeScore(activeCell, tmpGapOpenVertical, TraceBitMap_::VERTICAL, TraceBitMap_::VERTICAL_OPEN, TTracebackConfig(), RecursionDirectionVertical()) | TraceBitMap_::MAX_FROM_VERTICAL_MATRIX;
522 }
523
524 } // namespace seqan
525
526 #endif // #ifndef SEQAN_INCLUDE_SEQAN_ALIGN_DP_FORMULA_AFFINE_H_
527