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 DPScout for the banded chain alignment algorithm.
35 // The ScoutState is used to keep track of the initialization values for
36 // the next matrix.
37 // ==========================================================================
38 
39 #ifndef INCLUDE_SEQAN_SEEDS_BANDED_CHAIN_ALIGNMENT_SCOUT_H_
40 #define INCLUDE_SEQAN_SEEDS_BANDED_CHAIN_ALIGNMENT_SCOUT_H_
41 
42 namespace seqan {
43 
44 // ============================================================================
45 // Forwards
46 // ============================================================================
47 
48 // ============================================================================
49 // Tags, Classes, Enums
50 // ============================================================================
51 
52 // ----------------------------------------------------------------------------
53 // Tag BandedChainAlignmentScout
54 // ----------------------------------------------------------------------------
55 
56 struct BandedChainAlignmentScout_;
57 typedef Tag<BandedChainAlignmentScout_> BandedChainAlignmentScout;
58 
59 // ----------------------------------------------------------------------------
60 // Class BandedChainAlignmentScoutState
61 // ----------------------------------------------------------------------------
62 
63 // Used to determine the scout state of the BandedChainAlignmentScout.
64 template <typename TSpec>
65 struct BandedChainAlignmentScoutState;
66 
67 // ----------------------------------------------------------------------------
68 // Class DPScoutState_
69 // ----------------------------------------------------------------------------
70 
71 // Stores the state of the algorithm.
72 // It keeps track of the initialization values for horziontal and vertical
73 // direction of the current matrix and the next matrix that intersects with
74 // the current active matrix.
75 template <typename TDPCell>
76 class DPScoutState_<BandedChainAlignmentScoutState<TDPCell> >
77 {
78 public:
79     typedef Triple<unsigned, unsigned, TDPCell> TInitCell;
80     typedef std::set<TInitCell> TInitializationCellSet;
81 
82     unsigned int _horizontalNextGridOrigin;
83     unsigned int _verticalNextGridOrigin;
84 
85     String<TDPCell> _horizontalInitCurrentMatrix;
86     String<TDPCell> _verticalInitCurrentMatrix;
87     String<TDPCell> _horizontalInitNextMatrix;
88     String<TDPCell> _verticalInitNextMatrix;
89     TInitializationCellSet _nextInitializationCells;
90 
DPScoutState_()91     DPScoutState_() : _horizontalNextGridOrigin(0u), _verticalNextGridOrigin(0u)
92     {}
93 
94 };
95 
96 // ----------------------------------------------------------------------------
97 // Class DPScout_
98 // ----------------------------------------------------------------------------
99 
100 template<typename TDPCell>
101 class DPScout_<TDPCell, BandedChainAlignmentScout>
102 {
103 public:
104     typedef typename ScoutStateSpecForScout_<DPScout_>::Type TDPScoutStateSpec;
105     typedef DPScoutState_<TDPScoutStateSpec> TScoutState;
106     typedef String<unsigned int> TMaxHostPositionString;
107     typedef typename Value<TDPCell>::Type TScoreValue;
108 
109 //    TScoreValue _maxScore;      // the maximal score detected
110     TDPCell _maxScore;
111     TScoutState * _dpScoutStatePtr;
112     TMaxHostPositionString _maxHostPositions;    // the array containing all positions that have a maximal score
113 
DPScout_()114     DPScout_() : _maxScore(),
115                 _dpScoutStatePtr(0),
116                 _maxHostPositions()
117     {}
118 
DPScout_(TScoutState & scoutState)119     DPScout_(TScoutState & scoutState) :
120                 _maxScore(),
121                 _dpScoutStatePtr(&scoutState),
122                 _maxHostPositions()
123     {}
124 };
125 
126 // ============================================================================
127 // Metafunctions
128 // ============================================================================
129 
130 // ----------------------------------------------------------------------------
131 // Metafunction ScoutSpecForAlignmentAlgorithm_
132 // ----------------------------------------------------------------------------
133 
134 template <typename TSpec, typename TDPMatrixLocation>
135 struct ScoutSpecForAlignmentAlgorithm_<BandedChainAlignment_<TSpec, TDPMatrixLocation> >
136 {
137     typedef BandedChainAlignmentScout Type;
138 };
139 
140 // ----------------------------------------------------------------------------
141 // Metafunction ScoutStateSpecForScout_
142 // ----------------------------------------------------------------------------
143 
144 template <typename TDPCell>
145 struct ScoutStateSpecForScout_<DPScout_<TDPCell, BandedChainAlignmentScout> >
146 {
147     typedef BandedChainAlignmentScoutState<TDPCell> Type;
148 };
149 
150 // ============================================================================
151 // Functions
152 // ============================================================================
153 
154 // ----------------------------------------------------------------------------
155 // Function _rinitScout()
156 // ----------------------------------------------------------------------------
157 
158 template <typename TDPCell>
159 inline void
160 _reinitScout(DPScout_<TDPCell, BandedChainAlignmentScout> & dpScout)
161 {
162     //typedef typename Value<TDPCell>::Type TScoreValue;
163 
164     dpScout._maxScore = TDPCell();
165     resize(dpScout._maxHostPositions, 1, 0);
166 }
167 
168 // ----------------------------------------------------------------------------
169 // Function _reinitScoutState()
170 // ----------------------------------------------------------------------------
171 
172 template <typename TDPCell, typename TPosH, typename TPosV, typename TSizeCurrInit, typename TSizeNextInit>
173 inline void
174 _reinitScoutState(DPScoutState_<BandedChainAlignmentScoutState<TDPCell> > & scoutState,
175                   TPosH const & originNextMatrixH,
176                   TPosV const & originNextMatrixV,
177                   TSizeCurrInit const & sizeCurrMatrixInitH,
178                   TSizeCurrInit const & sizeCurrMatrixInitV,
179                   TSizeNextInit const & sizeNextMatrixInitH,
180                   TSizeNextInit const & sizeNextMatrixInitV)
181 {
182     typedef DPScoutState_<BandedChainAlignmentScoutState<TDPCell> > TDPScoutState;
183     typedef typename TDPScoutState::TInitializationCellSet TInitCellSet;
184     typedef typename TInitCellSet::iterator TInitCellSetIterator;
185 
186     scoutState._horizontalNextGridOrigin = originNextMatrixH;
187     scoutState._verticalNextGridOrigin = originNextMatrixV;
188 
189     // initialize the current initialization values of the tracker
190     arrayFill(begin(scoutState._horizontalInitCurrentMatrix), end(scoutState._horizontalInitCurrentMatrix), TDPCell());
191     arrayFill(begin(scoutState._verticalInitCurrentMatrix), end(scoutState._verticalInitCurrentMatrix), TDPCell());
192     arrayFill(begin(scoutState._horizontalInitNextMatrix), end(scoutState._horizontalInitNextMatrix), TDPCell());
193     arrayFill(begin(scoutState._verticalInitNextMatrix), end(scoutState._verticalInitNextMatrix), TDPCell());
194 
195     // check if the value needs to be resized ... (can be longer but not smaller)
196     if ((TSizeCurrInit) length(scoutState._horizontalInitCurrentMatrix) < sizeCurrMatrixInitH)
197         resize(scoutState._horizontalInitCurrentMatrix, sizeCurrMatrixInitH, TDPCell());
198     if ((TSizeCurrInit) length(scoutState._verticalInitCurrentMatrix) < sizeCurrMatrixInitV)
199         resize(scoutState._verticalInitCurrentMatrix, sizeCurrMatrixInitV, TDPCell());
200 
201     if ((TSizeNextInit) length(scoutState._horizontalInitNextMatrix) < sizeNextMatrixInitH)
202         resize(scoutState._horizontalInitNextMatrix, sizeNextMatrixInitH, TDPCell());
203     if ((TSizeNextInit) length(scoutState._verticalInitNextMatrix) < sizeNextMatrixInitV)
204         resize(scoutState._verticalInitNextMatrix, sizeNextMatrixInitV, TDPCell());
205 
206     // Parsing the set to get the values.
207     for (TInitCellSetIterator it = scoutState._nextInitializationCells.begin();  // We need to plant the initialization values here. At the moment we only put the old values here.
208         it != scoutState._nextInitializationCells.end(); ++it)
209     {
210 //        std::cerr << "TInitCell == " << it->i1 << " " << it->i2 << " " << _scoreOfCell(it->i3) << std::endl;
211         if (it->i1 == 0)
212             scoutState._verticalInitCurrentMatrix[it->i2] = it->i3;
213         if (it->i2 == 0)
214             scoutState._horizontalInitCurrentMatrix[it->i1] = it->i3;
215 
216     }
217 }
218 
219 // ----------------------------------------------------------------------------
220 // Function _setScoutState()
221 // ----------------------------------------------------------------------------
222 
223 // Sets the state of the scout.
224 template <typename TDPCell, typename TStateSpec>
225 inline void
226 _setScoutState(DPScout_<TDPCell, BandedChainAlignmentScout> & dpScout,
227                DPScoutState_<BandedChainAlignmentScoutState<TStateSpec> > & state)
228 {
229     dpScout._dpScoutStatePtr = & state;
230 }
231 
232 // ----------------------------------------------------------------------------
233 // Function _scoutBestScore()
234 // ----------------------------------------------------------------------------
235 
236 
237 template <typename TDPCell, typename TTraceMatrixNavigator>
238 inline void
239 _scoutBestScore(DPScout_<TDPCell, BandedChainAlignmentScout> & dpScout,
240                 TDPCell const & dpCell,
241                 TTraceMatrixNavigator const & navigator,
242                 bool isLastColumn,
243                 bool isLastRow,
244                 bool trackNextInitColumn,
245                 bool trackNextInitRow)
246 {
247     // Store value for vertical initialization of next grid.
248     if(trackNextInitColumn)
249         dpScout._dpScoutStatePtr->_verticalInitNextMatrix[coordinate(navigator, +DPMatrixDimension_::VERTICAL) -
250                                                        dpScout._dpScoutStatePtr->_verticalNextGridOrigin] = dpCell;
251     // Store value for horizontal initialization of next grid.
252     if (trackNextInitRow)
253         dpScout._dpScoutStatePtr->_horizontalInitNextMatrix[coordinate(navigator, +DPMatrixDimension_::HORIZONTAL) -
254                                                          dpScout._dpScoutStatePtr->_horizontalNextGridOrigin] = dpCell;
255 
256     // Now we have to track the optimal score from the last column or row.
257     if (isLastColumn || isLastRow)
258     {
259         if (_scoreOfCell(dpCell) >= _scoreOfCell(dpScout._maxScore))
260         {
261             if (_scoreOfCell(dpCell) == _scoreOfCell(dpScout._maxScore))
262                 appendValue(dpScout._maxHostPositions, position(navigator));
263             else
264             {
265                 resize(dpScout._maxHostPositions, 1);
266                 dpScout._maxHostPositions[0] = position(navigator);
267                 dpScout._maxScore = dpCell;
268             }
269         }
270     }
271 }
272 
273 template <typename TDPCell, typename TTraceMatrixNavigator>
274 inline void
275 _scoutBestScore(DPScout_<TDPCell, BandedChainAlignmentScout> &,
276                 TDPCell const &,
277                 TTraceMatrixNavigator const &,
278                 bool isLastColumn = false,
279                 bool isLastRow = false)
280 {
281     (void) isLastColumn;
282     (void) isLastRow;
283     //no-op
284 }
285 
286 // Delegate.
287 template <typename TDPCell, typename TTraceMatrixNavigator, typename TIsLastColumn, typename TIsLastRow>
288 inline void
289 _scoutBestScore(DPScout_<TDPCell, BandedChainAlignmentScout> & dpScout,
290                 TDPCell const & activeCell,
291                 TTraceMatrixNavigator const & navigator,
292                 TIsLastColumn const & /**/,
293                 TIsLastRow const & /**/)
294 {
295     _scoutBestScore(dpScout, activeCell, navigator, TIsLastColumn::VALUE, TIsLastRow::VALUE);
296 }
297 
298 // ----------------------------------------------------------------------------
299 // Function maxScore()
300 // ----------------------------------------------------------------------------
301 
302 template <typename TDPCell>
303 inline typename Value<TDPCell>::Type &
304 maxScore(DPScout_<TDPCell, BandedChainAlignmentScout> & scout)
305 {
306     return _scoreOfCell(scout._maxScore);
307 }
308 
309 template <typename TDPCell>
310 inline typename Value<TDPCell>::Type const &
311 maxScore(DPScout_<TDPCell, BandedChainAlignmentScout> const & scout)
312 {
313     return _scoreOfCell(scout._maxScore);
314 }
315 
316 // ----------------------------------------------------------------------------
317 // Function maxHostPositions()
318 // ----------------------------------------------------------------------------
319 
320 template <typename TDPCell>
321 inline typename DPScout_<TDPCell, BandedChainAlignmentScout>::TMaxHostPositionString const &
322 maxHostPositions(DPScout_<TDPCell, BandedChainAlignmentScout> const & scout)
323 {
324     return scout._maxHostPositions;
325 }
326 
327 template <typename TDPCell>
328 inline typename DPScout_<TDPCell, BandedChainAlignmentScout>::TMaxHostPositionString &
329 maxHostPositions(DPScout_<TDPCell, BandedChainAlignmentScout> & scout)
330 {
331     return scout._maxHostPositions;
332 }
333 
334 // ----------------------------------------------------------------------------
335 // Function maxHostPosition()
336 // ----------------------------------------------------------------------------
337 
338 template <typename TDPCell>
339 inline unsigned int
340 maxHostPosition(DPScout_<TDPCell, BandedChainAlignmentScout> const & scout)
341 {
342     return scout._maxHostPositions[0];
343 }
344 
345 // ----------------------------------------------------------------------------
346 // Function nextMatrixBeginH()
347 // ----------------------------------------------------------------------------
348 
349 template <typename TDPCell>
350 inline unsigned int
351 _nextMatrixBeginH(DPScout_<TDPCell, BandedChainAlignmentScout> const & scout)
352 {
353     return scout._posH;
354 }
355 
356 // ----------------------------------------------------------------------------
357 // Function nextMatrixBeginV()
358 // ----------------------------------------------------------------------------
359 
360 template <typename TDPCell>
361 inline unsigned int
362 _nextMatrixBeginV(DPScout_<TDPCell, BandedChainAlignmentScout> const & scout)
363 {
364     return scout._posV;
365 }
366 
367 }  // namespace seqan
368 
369 #endif  // #ifndef INCLUDE_SEQAN_SEEDS_BANDED_CHAIN_ALIGNMENT_SCOUT_H_
370