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
33 #ifndef SEQAN_HEADER_GRAPH_ALIGN_TCOFFEE_PROGRESSIVE_H
34 #define SEQAN_HEADER_GRAPH_ALIGN_TCOFFEE_PROGRESSIVE_H
35
36 namespace SEQAN_NAMESPACE_MAIN
37 {
38
39 //////////////////////////////////////////////////////////////////////////////
40 // Progressive Alignment
41 //////////////////////////////////////////////////////////////////////////////
42
43 //////////////////////////////////////////////////////////////////////////////
44
45 template<typename TStringSet, typename TCargo, typename TSpec, typename TPosition, typename TSequence>
46 inline void
_buildLeafString(Graph<Alignment<TStringSet,TCargo,TSpec>> const & g,TPosition const pos,TSequence & alignSeq)47 _buildLeafString(Graph<Alignment<TStringSet, TCargo, TSpec> > const& g,
48 TPosition const pos,
49 TSequence& alignSeq)
50 {
51 SEQAN_CHECKPOINT
52 typedef Graph<Alignment<TStringSet, TCargo, TSpec> > TGraph;
53 typedef typename Size<TGraph>::Type TSize;
54 typedef typename Id<TGraph>::Type TId;
55 typedef typename VertexDescriptor<TGraph>::Type TVertexDescriptor;
56 typedef typename Value<TSequence>::Type TVertexString;
57
58 //TVertexDescriptor nilVertex = getNil<TVertexDescriptor>();
59 TStringSet& str = stringSet(g);
60 TId seqId = positionToId(str, pos);
61 TSize lenRoot = length(str[pos]);
62 TSize i = 0;
63 while(i<lenRoot) {
64 TVertexDescriptor nextVertex = findVertex(const_cast<TGraph&>(g), seqId, i);
65 //SEQAN_ASSERT(nextVertex != nilVertex);
66 //if (nextVertex == nilVertex) {
67 // std::cout << "Warning: Nil Vertex" << std::endl;
68 // TSize j = i + 1;
69 // while ((j < lenRoot) && (findVertex(const_cast<TGraph&>(g), seqId, j) == nilVertex)) ++j;
70 // nextVertex = addVertex(const_cast<TGraph&>(g), seqId, i, j-i);
71 //}
72 TVertexString vs;
73 appendValue(vs, nextVertex);
74 appendValue(alignSeq, vs, Generous());
75 i += fragmentLength(g, nextVertex);
76 }
77 }
78
79 //////////////////////////////////////////////////////////////////////////////
80
81 template<typename TStringSet, typename TCargo, typename TSpec, typename TSegmentString, typename TOutGraph>
82 inline void
_createAlignmentGraph(Graph<Alignment<TStringSet,TCargo,TSpec>> const & g,TSegmentString & alignSeq,TOutGraph & gOut)83 _createAlignmentGraph(Graph<Alignment<TStringSet, TCargo, TSpec> > const& g,
84 TSegmentString& alignSeq,
85 TOutGraph& gOut)
86 {
87 SEQAN_CHECKPOINT
88 typedef Graph<Alignment<TStringSet, TCargo, TSpec> > TGraph;
89 typedef typename Size<TGraph>::Type TSize;
90 typedef typename VertexDescriptor<TGraph>::Type TVertexDescriptor;
91 typedef String<TVertexDescriptor> TVertexString;
92
93 // Create the alignment graph
94 TSize alignSeqLen = length(alignSeq);
95 for(TSize i = 0; i<alignSeqLen;++i) {
96 TVertexString& alignSeq_i = alignSeq[i];
97 TSize len_i = length(alignSeq_i);
98 for(TSize j=0; j<len_i; ++j) {
99 TVertexDescriptor v = alignSeq_i[j];
100 SEQAN_ASSERT(fragmentBegin(g,v) < length(getValueById(stringSet(g), sequenceId(g,v))));
101 SEQAN_ASSERT(fragmentLength(g,v) > 0);
102 SEQAN_ASSERT(fragmentBegin(g,v) + fragmentLength(g,v) <= length(getValueById(stringSet(g), sequenceId(g,v))));
103 TVertexDescriptor l = addVertex(gOut, sequenceId(g, v), fragmentBegin(g,v), fragmentLength(g,v));
104 //std::cout << l << label(gOut, l) << ',';
105 TSize count = 1;
106 for(TSize k = j; k>0; --k) {
107 //SEQAN_ASSERT(fragmentLength(gOut,l) == fragmentLength(gOut,l - count));
108 addEdge(gOut, (TVertexDescriptor) (l - count), (TVertexDescriptor) l);
109 ++count;
110 }
111 }
112 //std::cout << std::endl;
113 }
114 }
115
116 //////////////////////////////////////////////////////////////////////////////
117
118 /*!
119 * @fn progressiveAlignment
120 * @headerfile <seqan/graph_msa.h>
121 * @brief Perform a progressive multiple sequence alignment (MSA).
122 *
123 * @signature void progressiveAlignment(inputGraph, guideTree, outputGraph);
124 *
125 * @param[in] inputGraph A @link AlignmentGraph @endlink with multiple sequence information.
126 * @param[in] guideTree A @link Tree @endlink to use as the guide tree.
127 * @param[out] outputGraph An @link AlignmentGraph @endlink for the final MSA.
128 */
129
130 template<typename TStringSet, typename TCargo, typename TSpec, typename TGuideTree, typename TOutGraph>
131 inline void
progressiveAlignment(Graph<Alignment<TStringSet,TCargo,TSpec>> & g,TGuideTree & tree,TOutGraph & gOut)132 progressiveAlignment(Graph<Alignment<TStringSet, TCargo, TSpec> >& g,
133 TGuideTree& tree,
134 TOutGraph& gOut)
135 {
136 SEQAN_CHECKPOINT
137 typedef Graph<Alignment<TStringSet, TCargo, TSpec> > TGraph;
138 typedef typename Size<TGraph>::Type TSize;
139 typedef typename VertexDescriptor<TGuideTree>::Type TVertexDescriptor;
140 typedef typename Iterator<TGuideTree, BfsIterator>::Type TBfsIterator;
141 typedef typename Iterator<TGuideTree, AdjacencyIterator>::Type TAdjacencyIterator;
142 typedef String<TVertexDescriptor> TVertexString;
143 typedef String<TVertexString> TSegmentString;
144
145 // Initialization
146 TVertexDescriptor rootVertex = getRoot(tree);
147 TSize nVertices = numVertices(tree);
148
149 // Vertices in reversed bfs order
150 TVertexString vertices;
151 resize(vertices, nVertices);
152
153 // All Strings of Strings of vertices for each node of the guide tree
154 String<TSegmentString> segString;
155 resize(segString, nVertices);
156
157 // Walk through the tree in bfs order
158 typedef typename Iterator<TVertexString, Standard>::Type TVertexIter;
159 TVertexIter itVert = begin(vertices, Standard());
160 TVertexIter itVertEnd = end(vertices, Standard());
161 --itVertEnd;
162 TBfsIterator bfsIt(tree, rootVertex);
163 for(;!atEnd(bfsIt);goNext(bfsIt), --itVertEnd)
164 *itVertEnd = *bfsIt;
165
166 // Progressive alignment
167 itVert = begin(vertices, Standard());
168 itVertEnd = end(vertices, Standard());
169 for(;itVert != itVertEnd; ++itVert) {
170 if(isLeaf(tree, *itVert)) _buildLeafString(g, *itVert, segString[*itVert]);
171 else {
172 // Align the two children (Binary tree)
173 TAdjacencyIterator adjIt(tree, *itVert);
174 TVertexDescriptor child1 = *adjIt; goNext(adjIt);
175 heaviestCommonSubsequence(g, segString[child1], segString[*adjIt], segString[*itVert]);
176 clear(segString[child1]);
177 clear(segString[*adjIt]);
178 }
179 }
180
181 // Create the alignment graph
182 _createAlignmentGraph(g, segString[rootVertex], gOut);
183 }
184
185
186
187 //////////////////////////////////////////////////////////////////////////////
188 // Progressive Matching, just for testing purposes
189 //////////////////////////////////////////////////////////////////////////////
190
191
192 //////////////////////////////////////////////////////////////////////////////
193
194 template<typename TStringSet, typename TCargo, typename TSpec, typename TSegmentString, typename TEdgeMap, typename TOutGraph>
195 inline void
_createMatchingGraph(Graph<Alignment<TStringSet,TCargo,TSpec>> const & g,TSegmentString & alignSeq,TEdgeMap & edgeMap,TOutGraph & gOut,TEdgeMap & edgeMapOut)196 _createMatchingGraph(Graph<Alignment<TStringSet, TCargo, TSpec> > const& g,
197 TSegmentString& alignSeq,
198 TEdgeMap& edgeMap,
199 TOutGraph& gOut,
200 TEdgeMap& edgeMapOut)
201 {
202 typedef Graph<Alignment<TStringSet, TCargo, TSpec> > TGraph;
203 typedef typename Size<TGraph>::Type TSize;
204 typedef typename VertexDescriptor<TGraph>::Type TVertexDescriptor;
205 typedef typename EdgeDescriptor<TGraph>::Type TEdgeDescriptor;
206 typedef String<TVertexDescriptor> TVertexString;
207
208 // Initialization
209 TVertexDescriptor nilVertex = getNil<TVertexDescriptor>();
210
211 // Create the matching graph
212 clear(edgeMapOut);
213 TSize alignSeqLen = length(alignSeq);
214 for(TSize i = 0; i<alignSeqLen;++i) {
215 TVertexString& alignSeq_i = alignSeq[i];
216 TSize len_i = length(alignSeq_i);
217 for(TSize j=0; j<len_i - 1; ++j) {
218 for(TSize k=j+1; k<len_i; ++k) {
219 TVertexDescriptor v1 = getValue(alignSeq_i, j);
220 TVertexDescriptor v2 = getValue(alignSeq_i, k);
221 if ((v1 == nilVertex) || (v2 == nilVertex)) continue;
222
223 TVertexDescriptor v1New = findVertex(gOut, sequenceId(g, v1), fragmentBegin(g,v1));
224 if (v1New == nilVertex) v1New = addVertex(gOut, sequenceId(g, v1), fragmentBegin(g,v1), fragmentLength(g,v1));
225 TVertexDescriptor v2New = findVertex(gOut, sequenceId(g, v2), fragmentBegin(g,v2));
226 if (v2New == nilVertex) v2New = addVertex(gOut, sequenceId(g, v2), fragmentBegin(g,v2), fragmentLength(g,v2));
227
228 TEdgeDescriptor e = findEdge(g, v1, v2);
229 addEdge(gOut, v1New, v2New, cargo(e));
230 appendValue(edgeMapOut, property(edgeMap, e));
231 }
232 }
233 }
234 }
235
236 //////////////////////////////////////////////////////////////////////////////
237
238 template<typename TStringSet, typename TCargo, typename TSpec, typename TString, typename TOutString>
239 inline TCargo
heaviestMatching(Graph<Alignment<TStringSet,TCargo,TSpec>> const & g,TString const & str1,TString const & str2,TOutString & align)240 heaviestMatching(Graph<Alignment<TStringSet, TCargo, TSpec> > const& g,
241 TString const& str1,
242 TString const& str2,
243 TOutString& align)
244 {
245 typedef Graph<Alignment<TStringSet, TCargo, TSpec> > TGraph;
246 typedef __int64 TLargeSize;
247 typedef typename Size<TStringSet>::Type TSize;
248 typedef typename VertexDescriptor<TGraph>::Type TVertexDescriptor;
249 typedef typename Iterator<TGraph, OutEdgeIterator>::Type TOutEdgeIterator;
250 typedef typename Value<TString>::Type TVertexSet;
251
252 TSize m = length(str1); // How many sets of vertex descriptors in seq1
253 TSize n = length(str2); // How many sets of vertex descriptors in seq1
254
255 // Size of the sequences
256 // Note for profile alignments every member of the sequence is a String!!! of vertex descriptors
257 TCargo divider = (TCargo) length(str1[0]) * (TCargo) length(str2[0]);
258 TVertexDescriptor nilVertex = getNil<TVertexDescriptor>();
259
260 // Fill the vertex to position map for str1
261 // Remember for each vertex descriptor the position in the sequence
262 typedef std::map<TVertexDescriptor, TSize> TVertexToPosMap;
263 typedef typename TVertexToPosMap::const_iterator TVertexToPosMapIter;
264 TVertexToPosMap map;
265 typedef typename Iterator<TString const>::Type TStringIterConst;
266 typedef typename Iterator<TVertexSet const>::Type TVertexSetIterConst;
267 TStringIterConst itStrEnd1 = end(str1);
268 TSize pos = 0;
269 for(TStringIterConst itStr1 = begin(str1);itStr1 != itStrEnd1;++itStr1, ++pos) {
270 TVertexSetIterConst itVEnd = end(getValue(itStr1));
271 for(TVertexSetIterConst itV = begin(getValue(itStr1));itV != itVEnd;++itV) {
272 if (*itV != nilVertex) map.insert(std::make_pair(*itV, pos));
273 }
274 }
275
276 // Given a full graph, what positions are occupied?
277 typedef std::set<TLargeSize> TOccupiedPositions;
278 TOccupiedPositions occupiedPositions;
279 TStringIterConst itStrEnd2 = end(str2);
280 TSize posItStr2 = 0;
281 for(TStringIterConst itStr2 = begin(str2);itStr2 != itStrEnd2;++itStr2, ++posItStr2) {
282 TVertexSetIterConst itVEnd = end(getValue(itStr2));
283 for(TVertexSetIterConst itV = begin(getValue(itStr2));itV != itVEnd;++itV) {
284 if (*itV != nilVertex) {
285 TOutEdgeIterator itOut(g, *itV);
286 for(;!atEnd(itOut); ++itOut) {
287 // Target vertex must be in the map
288 TVertexToPosMapIter pPos = map.find(targetVertex(itOut));
289 if (pPos != map.end()) occupiedPositions.insert( (TLargeSize) (pPos->second) * (TLargeSize) n + (TLargeSize) (n - posItStr2 - 1) );
290 }
291 }
292 }
293 }
294 // Map the occupied positions to slots
295 typedef std::map<TLargeSize, TSize> TPositionToSlotMap;
296 TPositionToSlotMap posToSlotMap;
297 TSize counter = 0;
298 for(typename TOccupiedPositions::const_iterator setIt = occupiedPositions.begin();setIt != occupiedPositions.end(); ++setIt, ++counter) {
299 posToSlotMap.insert(std::make_pair(*setIt, counter));
300 }
301 occupiedPositions.clear();
302
303 // Walk through str2 and fill in the weights of the actual edges
304 typedef String<TCargo> TWeights;
305 typedef typename Iterator<TWeights>::Type TWeightsIter;
306 TWeights weights;
307 resize(weights, posToSlotMap.size(),0);
308 posItStr2 = 0;
309 for(TStringIterConst itStr2 = begin(str2);itStr2 != itStrEnd2;++itStr2, ++posItStr2) {
310 TVertexSetIterConst itVEnd = end(getValue(itStr2));
311 for(TVertexSetIterConst itV = begin(getValue(itStr2));itV != itVEnd;++itV) {
312 if (*itV != nilVertex) {
313 TOutEdgeIterator itOut(g, *itV);
314 for(;!atEnd(itOut); ++itOut) {
315 // Target vertex must be in the map
316 TVertexToPosMapIter pPos = map.find(targetVertex(itOut));
317 if (pPos != map.end()) weights[posToSlotMap[ (TLargeSize) (pPos->second) * (TLargeSize) n + (TLargeSize) (n - posItStr2 - 1) ]] += (TCargo) cargo(*itOut);
318 }
319 }
320 }
321 }
322 map.clear();
323 // Average weights
324 TWeightsIter itWeights = begin(weights);
325 TWeightsIter itWeightsEnd = begin(weights);
326 for(;itWeights != itWeightsEnd; ++itWeights) *itWeights /= divider;
327
328 // Create the match graph
329 Graph<Undirected<void> > matchGraph;
330 typedef typename VertexDescriptor<Graph<Undirected<void> > >::Type TVD;
331 String<bool> vertexMap;
332 resize(vertexMap, m + n, false);
333 for(TSize i = 0; i < m; ++i) addVertex(matchGraph);
334 for(TSize i = 0; i < n; ++i) value(vertexMap, addVertex(matchGraph)) = true;
335 for(typename TPositionToSlotMap::const_iterator mapIt = posToSlotMap.begin();mapIt != posToSlotMap.end(); ++mapIt) {
336 TLargeSize pos = mapIt->first;
337 TSize i = (TSize) (pos / (TLargeSize) n); // Get the index in str1
338 TSize j = n - 1 - (TSize) (pos % (TLargeSize) n); // Get the index in str2
339 addEdge(matchGraph, i, m+j);
340 }
341
342 // Compute the best matching
343 typedef String<Pair<TVD, TVD> > TEdges;
344 TEdges edges;
345
346 //std::fstream strm;
347 //strm.open("Z:\\my_graph.dot", std::ios_base::out | std::ios_base::trunc);
348 //writeRecords(strm,matchGraph,DotDrawing());
349 //strm.close();
350
351 TCargo val = weightedBipartiteMatching(matchGraph, vertexMap, weights, edges);
352
353 // Retrieve the aligned segments
354 TSize seqsInStr1 = length(str1[0]);
355 TSize seqsInStr2 = length(str2[0]);
356 typedef typename Value<TString>::Type TVertexSet;
357 //typedef typename Iterator<TString const, Rooted>::Type TStringIter;
358 //typedef typename Iterator<TString, Rooted>::Type TSIter;
359 typedef typename Iterator<TVertexSet const, Rooted>::Type TVertexSetIter;
360 //typedef typename Iterator<TVertexSet, Rooted>::Type TIter;
361 clear(align);
362 // Retrieve all matches
363 String<bool> matchedVertices;
364 resize(matchedVertices, n + m, false);
365 typedef typename Iterator<TEdges>::Type TEdgesIter;
366 TEdgesIter itEdges = begin(edges);
367 TEdgesIter itEdgesEnd = end(edges);
368 for(;itEdges != itEdgesEnd; ++itEdges) {
369 TSize i = (value(itEdges)).i1;
370 TSize j = (value(itEdges)).i2 - m;
371 value(matchedVertices, i) = true;
372 value(matchedVertices, m+j) = true;
373 TVertexSet tmp;
374 resize(tmp, (seqsInStr1 + seqsInStr2), nilVertex);
375 TVertexSetIter itVEnd = end(value(str1,i));
376 TSize count = 0;
377 for(TVertexSetIter itV = begin(value(str1,i));itV != itVEnd;++itV) {
378 tmp[count] = *itV;
379 ++count;
380 }
381 TVertexSetIter itVEnd2 = end(value(str2,j));
382 for(TVertexSetIter itV2 = begin(value(str2,j));itV2 != itVEnd2;++itV2) {
383 tmp[count] = *itV2;
384 ++count;
385 }
386 appendValue(align, tmp);
387 }
388 typedef typename Iterator<Graph<Undirected<void> >, VertexIterator>::Type TVertexIterator;
389 TVertexIterator itVertex(matchGraph);
390 for(;!atEnd(itVertex);++itVertex) {
391 if (getProperty(matchedVertices, *itVertex) == true) continue;
392 if (*itVertex < m) {
393 TSize i = *itVertex;
394 TVertexSet tmp;
395 resize(tmp, (seqsInStr1 + seqsInStr2), nilVertex);
396 TVertexSetIter itVEnd = end(value(str1, i));
397 TSize count = 0;
398 for(TVertexSetIter itV = begin(value(str1, i));itV != itVEnd;++itV) {
399 tmp[count] = *itV;
400 ++count;
401 }
402 appendValue(align, tmp);
403 } else {
404 TSize j = *itVertex - m;
405 TVertexSet tmp;
406 resize(tmp, (seqsInStr1 + seqsInStr2), nilVertex);
407 TVertexSetIter itVEnd = end(value(str2, j));
408 TSize count = 0;
409 for(TVertexSetIter itV = begin(value(str2, j));itV != itVEnd;++itV) {
410 tmp[count] = *itV;
411 ++count;
412 }
413 appendValue(align, tmp);
414 }
415 }
416
417 return val;
418 }
419
420 //////////////////////////////////////////////////////////////////////////////
421
422 template<typename TStringSet, typename TCargo, typename TSpec, typename TGuideTree, typename TVertexDescriptor, typename TSequence>
423 inline void
_recursiveProgressiveMatching(Graph<Alignment<TStringSet,TCargo,TSpec>> & g,TGuideTree & tree,TVertexDescriptor const root,TSequence & alignSeq)424 _recursiveProgressiveMatching(Graph<Alignment<TStringSet, TCargo, TSpec> >& g,
425 TGuideTree& tree,
426 TVertexDescriptor const root,
427 TSequence& alignSeq)
428 {
429 //typedef Graph<Alignment<TStringSet, TCargo, TSpec> > TGraph;
430 //typedef typename Size<TGraph>::Type TSize;
431 //typedef typename Id<TGraph>::Type TId;
432 typedef typename Iterator<TGuideTree, AdjacencyIterator>::Type TAdjacencyIterator;
433 //typedef typename Iterator<TGuideTree, DfsPreorder>::Type TDfsPreorderIterator;
434
435 if(isLeaf(tree, root)) {
436 _buildLeafString(g, root, alignSeq);
437 } else {
438 // Align the two children (Binary tree)
439 typedef String<String<TVertexDescriptor> > TSegmentString;
440 TSegmentString seq1;
441 TSegmentString seq2;
442 TAdjacencyIterator adjIt(tree, root);
443 _recursiveProgressiveMatching(g,tree, *adjIt, seq1);
444 goNext(adjIt);
445 _recursiveProgressiveMatching(g,tree, *adjIt, seq2);
446
447 heaviestMatching(g, seq1, seq2, alignSeq);
448
449 //// Debug Code
450 //for(TSize i = 0; i<length(alignSeq);++i) {
451 // std::cout << '(';
452 // for(TSize j=0; j<length(alignSeq[i]);++j) {
453 // std::cout << getValue(alignSeq[i], j) << ',';
454 // }
455 // std::cout << ')' << ',';
456 //}
457 //std::cout << std::endl;
458 }
459 }
460
461 //////////////////////////////////////////////////////////////////////////////
462
463 template<typename TStringSet, typename TCargo, typename TSpec, typename TGuideTree, typename TEdgeMap, typename TOutGraph>
464 inline void
progressiveMatching(Graph<Alignment<TStringSet,TCargo,TSpec>> & g,TGuideTree & tree,TEdgeMap & edgeMap,TOutGraph & gOut,TEdgeMap & edgeMapOut)465 progressiveMatching(Graph<Alignment<TStringSet, TCargo, TSpec> >& g,
466 TGuideTree& tree,
467 TEdgeMap& edgeMap,
468 TOutGraph& gOut,
469 TEdgeMap& edgeMapOut)
470 {
471 typedef Graph<Alignment<TStringSet, TCargo, TSpec> > TGraph;
472 //typedef typename Size<TGraph>::Type TSize;
473 typedef typename VertexDescriptor<TGraph>::Type TVertexDescriptor;
474 typedef String<TVertexDescriptor> TVertexString;
475 typedef String<TVertexString> TSegmentString;
476
477 // Perform progressive alignment
478 TSegmentString alignSeq;
479 _recursiveProgressiveMatching(g,tree,getRoot(tree),alignSeq);
480
481 // Create the alignment graph
482 _createMatchingGraph(g, alignSeq, edgeMap, gOut, edgeMapOut);
483 }
484
485
486 }// namespace SEQAN_NAMESPACE_MAIN
487
488 #endif //#ifndef SEQAN_HEADER_...
489