~ubuntu-branches/ubuntu/trusty/abyss/trusty-proposed

« back to all changes in this revision

Viewing changes to MergePaths/MergePaths.cpp

• browse files at revision 7
• compare with another revision
• download diff
• download tarball
• view history from revision 7

Show diffs side-by-side

added

removed

MergePaths/MergePaths.cpp

1

1

#include "config.h"

2

2

#include "Common/Options.h"

3

 

#include "ContigGraph.h"

4

 

#include "ContigGraphAlgorithms.h"

5

3

#include "ContigID.h"

6

4

#include "ContigLength.h"

7

5

#include "ContigPath.h"

8

 

#include "DirectedGraph.h"

9

 

#include "DotIO.h"

10

 

#include "GraphAlgorithms.h"

11

 

#include "GraphUtil.h"

 

6

#include "Functional.h" // for mem_var

12

7

#include "IOUtil.h"

13

8

#include "Iterator.h"

14

9

#include "Uncompress.h"

 

10

#include "Graph/Assemble.h"

 

11

#include "Graph/ContigGraph.h"

 

12

#include "Graph/DirectedGraph.h"

 

13

#include "Graph/DotIO.h"

 

14

#include "Graph/GraphAlgorithms.h"

 

15

#include "Graph/GraphUtil.h"

15

16

#include <boost/tuple/tuple.hpp>

16

17

#include <algorithm>

17

18

#include <cassert>

53

54

"  -k, --kmer=KMER_SIZE  k-mer size\n"

54

55

"  -s, --seed-length=L   minimum length of a seed contig [0]\n"

55

56

"  -o, --out=FILE        write result to FILE\n"

56

 

"      --greedy          use the greedy algorithm [default]\n"

57

 

"      --no-greedy       do not use the greedy algorithm\n"

 

57

"      --no-greedy       use the non-greedy algorithm [default]\n"

 

58

"      --greedy          use the greedy algorithm\n"

58

59

"  -g, --graph=FILE      write the path overlap graph to FILE\n"

59

60

"  -j, --threads=N       use N parallel threads [1]\n"

60

61

"  -v, --verbose         display verbose output\n"

72

73

        static unsigned seedLen;

73

74

 

74

75

        /** Use a greedy algorithm. */

75

 

        static int greedy = true;

 

76

        static int greedy;

76

77

 

77

78

        /** Write the path overlap graph to this file. */

78

79

        static string graphPath;

233

234

        assert(it != paths.end());

234

235

        ContigPath path = it->second;

235

236

        if (u.sense())

236

 

                path.reverseComplement();

 

237

                reverseComplement(path.begin(), path.end());

237

238

        return path;

238

239

}

239

240

 

256

257

 

257

258

                ContigPath path2 = path2It->second;

258

259

                if (seed2.sense())

259

 

                        path2.reverseComplement();

 

260

                        reverseComplement(path2.begin(), path2.end());

260

261

                addOverlapEdge(gout, seed2, seed1, path1, seed2, path2);

261

262

        }

262

263

}

279

280

 

280

281

                ContigPath path2 = path2It->second;

281

282

                if (pivot.sense())

282

 

                        path2.reverseComplement();

 

283

                        reverseComplement(path2.begin(), path2.end());

283

284

                ContigPath consensus = align(path, path2, pivot);

284

285

                if (consensus.empty()) {

285

286

                        invalid.push_back(pivot);

311

312

        assert(path1It != paths.end());

312

313

        ContigPath path(path1It->second);

313

314

        if (seedPath.front().sense())

314

 

                path.reverseComplement();

 

315

                reverseComplement(path.begin(), path.end());

315

316

        if (opt::verbose > 1)

316

317

#pragma omp critical(cout)

317

318

                cout << "\n* " << seedPath << '\n'

324

325

                assert(path2It != paths.end());

325

326

                ContigPath path2 = path2It->second;

326

327

                if (seed2.sense())

327

 

                        path2.reverseComplement();

 

328

                        reverseComplement(path2.begin(), path2.end());

328

329

 

329

330

                ContigNode pivot

330

331

                        = find(path.begin(), path.end(), seed2) != path.end()

455

456

                        continue;

456

457

                ContigPath path2 = path2It->second;

457

458

                if (pivot.sense())

458

 

                        path2.reverseComplement();

 

459

                        reverseComplement(path2.begin(), path2.end());

459

460

                ContigPath consensus = align(path, path2, pivot);

460

461

                if (consensus.empty())

461

462

                        continue;

672

673

                PathGraph& pathGraph, ContigPathMap& paths)

673

674

{

674

675

        ContigPaths seedPaths;

675

 

        assemble(pathGraph, back_inserter(seedPaths));

 

676

        assembleDFS(pathGraph, back_inserter(seedPaths));

676

677

        ContigPaths mergedPaths = mergeSeedPaths(paths, seedPaths);

677

678

        if (opt::verbose > 1)

678

679

                cout << '\n';

682

683

                        it1 != seedPaths.end(); ++it1) {

683

684

                const ContigPath& path(mergedPaths[it1 - seedPaths.begin()]);

684

685

                ContigPath pathrc(path);

685

 

                pathrc.reverseComplement();

 

686

                reverseComplement(pathrc.begin(), pathrc.end());

686

687

                for (ContigPath::const_iterator it2 = it1->begin();

687

688

                                it2 != it1->end(); ++it2) {

688

689

                        ContigNode seed(*it2);

707

708

/** Read a set of paths from the specified file. */

708

709

static ContigPathMap readPaths(const string& filePath)

709

710

{

 

711

        if (opt::verbose > 0)

 

712

                cerr << "Reading `" << filePath << "'..." << endl;

710

713

        ifstream in(filePath.c_str());

711

714

        assert_good(in, filePath);

712

715

 

830

833

                omp_set_num_threads(opt::threads);

831

834

#endif

832

835

 

 

836

        if (opt::verbose > 0)

 

837

                cerr << "Reading `" << argv[optind] << "'..." << endl;

833

838

        g_contigLengths = readContigLengths(argv[optind++]);

834

839

        ContigPathMap originalPathMap = readPaths(argv[optind++]);

835

840

 

1252

1257

                const ContigPath& path1, const ContigPath& path2,

1253

1258

                ContigNode pivot, dir_type& orientation)

1254

1259

{

1255

 

        if (path1 == path2 && &path1 != &path2) {

1256

 

                // These two paths are identical. Ignore the trivial alignment

1257

 

                // when aligning a path to itself.

 

1260

        if (&path1 == &path2) {

 

1261

                // Ignore the trivial alignment when aligning a path to

 

1262

                // itself.

 

1263

        } else if (path1 == path2) {

 

1264

                // These two paths are identical.

1258

1265

                orientation = DIR_B;

1259

1266

                return path1;

 

1267

        } else {

 

1268

                ContigPath::const_iterator it

 

1269

                        = search(path1.begin(), path1.end(),

 

1270

                                path2.begin(), path2.end());

 

1271

                if (it != path1.end()) {

 

1272

                        // path2 is subsumed in path1.

 

1273

                        // Determine the orientation of the edge.

 

1274

                        orientation

 

1275

                                = it == path1.begin() ? DIR_R

 

1276

                                : it + path2.size() == path1.end() ? DIR_F

 

1277

                                : DIR_B;

 

1278

                        return path1;

 

1279

                }

1260

1280

        }

1261

1281

 

1262

1282

        // Find a suitable pivot.

• Older »

Loggerhead is a web-based interface for Breezy
Version: 2.0.1