@@ -139,9 +139,9 @@ void DBG::subgraph() {
139
139
140
140
std::vector<std::function<bool ()>> jobs;
141
141
std::array<uint16_t , 2 > mapRange = {0 ,0 };
142
-
142
+
143
143
while (mapRange[1 ] < mapCount) {
144
-
144
+
145
145
mapRange = computeMapRange (mapRange);
146
146
loadMapRange (mapRange);
147
147
@@ -154,30 +154,10 @@ void DBG::subgraph() {
154
154
jobs.clear ();
155
155
156
156
deleteMapRange (mapRange);
157
-
157
+
158
158
}
159
159
lg.verbose (" Merging subgraphs" );
160
160
mergeSubgraphs ();
161
- lg.verbose (" Searching graph" );
162
- if (userInput.travAlgorithm == " best-first" ) {
163
- if (userInput.kmerDepth == -1 )
164
- userInput.kmerDepth = userInput.kmerLen ; // unidirectional search
165
- bestFirst ();
166
- }else if (userInput.travAlgorithm == " traversal" ) {
167
- if (userInput.kmerDepth == -1 )
168
- userInput.kmerDepth = std::ceil ((float )userInput.kmerLen /2 ); // kmer search is in both directions
169
- traversal ();
170
- }else {
171
- fprintf (stderr, " Cannot find input algorithm (%s). Terminating.\n " , userInput.travAlgorithm .c_str ());
172
- exit (EXIT_FAILURE);
173
- }
174
- lg.verbose (" Remove missing edges" );
175
- removeMissingEdges ();
176
- lg.verbose (" Computing summary graph" );
177
- summary (*DBGsubgraph);
178
- lg.verbose (" Generating GFA" );
179
- DBGgraphToGFA ();
180
-
181
161
}
182
162
183
163
void DBG::summary (ParallelMap32color& DBGsubgraph) {
@@ -307,6 +287,21 @@ bool DBG::DBGsubgraphFromSegment(InSegment *inSegment, std::array<uint16_t, 2> m
307
287
return true ;
308
288
}
309
289
290
+ void DBG::searchGraph () {
291
+ if (userInput.travAlgorithm == " best-first" ) {
292
+ if (this ->userInput .kmerDepth == -1 )
293
+ userInput.kmerDepth = userInput.kmerLen ; // unidirectional search
294
+ bestFirst ();
295
+ }else if (userInput.travAlgorithm == " traversal" ) {
296
+ if (userInput.kmerDepth == -1 )
297
+ userInput.kmerDepth = std::ceil ((float )userInput.kmerLen /2 ); // kmer search is in both directions
298
+ traversal ();
299
+ }else {
300
+ fprintf (stderr, " Cannot find input algorithm (%s). Terminating.\n " , userInput.travAlgorithm .c_str ());
301
+ exit (EXIT_FAILURE);
302
+ }
303
+ }
304
+
310
305
void DBG::traversal () {
311
306
312
307
ParallelMap32color candidates, newCandidates;
@@ -426,25 +421,34 @@ ParallelMap32color DBG::traversalPass(ParallelMap32color* subgraph, std::array<u
426
421
void DBG::bestFirst () {
427
422
428
423
ParallelMap32color* candidates = new ParallelMap32color;
429
- uint32_t explored = 0 , total = DBGsubgraph->size ();
424
+ uint64_t explored = 0 , total = DBGsubgraph->size ();
430
425
std::array<uint16_t , 2 > mapRange;
431
426
ParallelMap32color* DBGsubgraphCpy = new ParallelMap32color;
427
+ bool * visited = new bool [DBGsubgraph->size ()]{false };
428
+
432
429
while (explored < total) {
433
430
434
431
mapRange = {0 ,0 };
432
+ uint64_t i = 0 ;
435
433
436
434
while (mapRange[1 ] < mapCount) {
437
435
438
436
mapRange = computeMapRange (mapRange);
439
437
loadMapRange (mapRange);
440
438
for (auto pair : *DBGsubgraph) {
441
- auto results = dijkstra (pair, mapRange);;
442
- explored += results.first ;
443
- if (results.first ) {
444
- candidates->insert (results.second .begin (), results.second .end ());
445
- DBGsubgraphCpy->insert (pair);
446
- // DBGsubgraph->erase(pair.first);
439
+
440
+ if (!visited[i]) {
441
+
442
+ auto results = dijkstra (pair, mapRange);;
443
+ explored += results.first ;
444
+ if (results.first ) {
445
+ candidates->insert (results.second .begin (), results.second .end ());
446
+ DBGsubgraphCpy->insert (pair);
447
+ // DBGsubgraph->erase(pair.first);
448
+ visited[i] = true ;
449
+ }
447
450
}
451
+ ++i;
448
452
// std::cout<<DBGsubgraphCpy.size()<<std::endl;
449
453
}
450
454
deleteMapRange (mapRange);
@@ -454,6 +458,7 @@ void DBG::bestFirst() {
454
458
delete DBGsubgraph;
455
459
DBGsubgraph = DBGsubgraphCpy;
456
460
delete candidates;
461
+ delete[] visited;
457
462
}
458
463
459
464
std::pair<bool ,ParallelMap32color> DBG::dijkstra (std::pair<uint64_t ,DBGkmer32color> source, std::array<uint16_t , 2 > mapRange) {
@@ -473,12 +478,12 @@ std::pair<bool,ParallelMap32color> DBG::dijkstra(std::pair<uint64_t,DBGkmer32col
473
478
int16_t depth = 0 ;
474
479
bool direction = true , isFw;
475
480
476
- while (Q.size () > 0 && depth < userInput.kmerDepth + 1 ) { // The main loop
481
+ while (Q.size () > 0 && depth < userInput.kmerDepth + 1 ) { // the main loop
477
482
explored = false ; // if there are still node in the queue we cannot be done
478
483
ParallelMap *map;
479
484
// ParallelMap32 *map32;
480
485
481
- std::pair<const uint64_t , DBGkmer32>* u = Q.extractMin (); // Remove and return best vertex
486
+ std::pair<const uint64_t , DBGkmer32>* u = Q.extractMin (); // remove and return best vertex
482
487
auto got = prev.find (u->first ); // check direction
483
488
if (got != prev.end ()) {
484
489
direction = got->second .second ;
@@ -622,6 +627,8 @@ void DBG::removeMissingEdges() {
622
627
}
623
628
}
624
629
}
630
+ lg.verbose (" Computing summary graph" );
631
+ summary (*DBGsubgraph);
625
632
}
626
633
627
634
0 commit comments