Commit 1b4ecc01 authored by Aaron's avatar Aaron
Browse files

fix corner case in chrom sweep. addresses GC issue #141

parent 6a853e2f
......@@ -246,14 +246,17 @@ void BedIntersect::IntersectBed() {
// treat the BED as a single "block"
if (_obeySplits == false)
FindOverlaps(a, hits);
// split the BED12 into blocks and look for overlaps in each discrete block
// split the BED12 into blocks and look for
// overlaps in each discrete block
else {
// find the hits that overlap with the full span of the blocked BED
// find the hits that overlap with the
// full span of the blocked BED
_bedB->allHits(a.chrom, a.start, a.end, a.strand,
hits, _sameStrand, _diffStrand,
0.0, false);
// break a into discrete blocks, as we need to
// measure overlap with the individual blocks, not the full span.
// measure overlap with the individual blocks,
// not the full span.
bedVector a_blocks;
GetBedBlocks(a, a_blocks);
// find the overlaps between the block in A and B
......@@ -275,12 +278,14 @@ void BedIntersect::IntersectBed() {
pair<BED, vector<BED> > hit_set;
hit_set.second.reserve(10000);
while (sweep.Next(hit_set)) {
if (_obeySplits == false)
processHits(hit_set.first, hit_set.second);
if (_obeySplits == false) {
processHits(hit_set.first, hit_set.second);
}
else {
bedVector a_blocks;
GetBedBlocks(hit_set.first, a_blocks);
FindBlockedOverlaps(hit_set.first, a_blocks, hit_set.second, false);
FindBlockedOverlaps(hit_set.first, a_blocks,
hit_set.second, false);
}
}
}
......@@ -351,7 +356,8 @@ void BedIntersect::IntersectBam(string bamFile) {
hits, _sameStrand, _diffStrand,
_overlapFraction, _reciprocal);
// find the overlaps between the block in A and B
overlapsFound = FindBlockedOverlaps(bed, bed_blocks, hits, _bamOutput);
overlapsFound =
FindBlockedOverlaps(bed, bed_blocks, hits, _bamOutput);
}
else if ((_bamOutput == false) && (_obeySplits == false))
{
......
......@@ -72,6 +72,7 @@ void ChromSweep::ScanCache() {
vector<BED>::iterator c = _cache.begin();
while (c != _cache.end())
{
if ((_curr_qy.chrom == c->chrom) && !(after(_curr_qy, *c))) {
if (IsValidHit(_curr_qy, *c)) {
_hits.push_back(*c);
......@@ -91,9 +92,11 @@ bool ChromSweep::ChromChange()
if (_curr_qy.chrom == _curr_db.chrom) {
return false;
}
// the query is ahead of the database. fast-forward the database to catch-up.
else if (_curr_qy.chrom > _curr_db.chrom) {
while (_db->GetNextBed(_curr_db, true) && _curr_db.chrom < _curr_qy.chrom)
// the query is ahead of the database.
// fast-forward the database to catch-up.
else if ((_curr_qy.chrom > _curr_db.chrom) && (!_db->Empty())) {
while (_db->GetNextBed(_curr_db, true) &&
_curr_db.chrom < _curr_qy.chrom)
{
}
_cache.clear();
......@@ -158,7 +161,8 @@ bool ChromSweep::Next(pair<BED, vector<BED> > &next) {
if (ChromChange() == false) {
// scan the database cache for hits
ScanCache();
// advance the db until we are ahead of the query. update hits and cache as necessary
// advance the db until we are ahead of the query.
// update hits and cache as necessary
while (!_db->Empty() &&
_curr_qy.chrom == _curr_db.chrom &&
!(after(_curr_db, _curr_qy)))
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment