Skip to content
Snippets Groups Projects
Commit 0df75e5e authored by Aaron's avatar Aaron
Browse files

reserve space up front for chromsweep

parent a5b2d831
No related branches found
No related tags found
No related merge requests found
...@@ -232,6 +232,7 @@ void BedIntersect::IntersectBed() { ...@@ -232,6 +232,7 @@ void BedIntersect::IntersectBed() {
ChromSweep sweep = ChromSweep(_bedA, _bedB, _sameStrand, _diffStrand); ChromSweep sweep = ChromSweep(_bedA, _bedB, _sameStrand, _diffStrand);
pair<BED, vector<BED> > hit_set; pair<BED, vector<BED> > hit_set;
hit_set.second.reserve(10000);
while (sweep.Next(hit_set)) { while (sweep.Next(hit_set)) {
processHits(hit_set.first, hit_set.second); processHits(hit_set.first, hit_set.second);
} }
......
...@@ -49,8 +49,8 @@ ChromSweep::ChromSweep(string &bedAFile, string &bedBFile) ...@@ -49,8 +49,8 @@ ChromSweep::ChromSweep(string &bedAFile, string &bedBFile)
_qy_lineNum = 0; _qy_lineNum = 0;
_db_lineNum = 0; _db_lineNum = 0;
_hits.reserve(1000); _hits.reserve(100000);
_cache.reserve(1000); _cache.reserve(100000);
_bedA = new BedFile(bedAFile); _bedA = new BedFile(bedAFile);
_bedB = new BedFile(bedBFile); _bedB = new BedFile(bedBFile);
...@@ -76,18 +76,8 @@ void ChromSweep::ScanCache() { ...@@ -76,18 +76,8 @@ void ChromSweep::ScanCache() {
while (c != _cache.end()) while (c != _cache.end())
{ {
if ((_curr_qy.chrom == c->chrom) && !(after(_curr_qy, *c))) { if ((_curr_qy.chrom == c->chrom) && !(after(_curr_qy, *c))) {
if (overlaps(_curr_qy.start, _curr_qy.end, c->start, c->end) > 0) { if (IsValidHit(_curr_qy, *c)) {
bool strands_are_same = (_curr_qy.strand == c->strand); _hits.push_back(*c);
// test for necessary strandedness
if ( (_sameStrand == false && _diffStrand == false)
||
(_sameStrand == true && strands_are_same == true)
||
(_diffStrand == true && strands_are_same == false)
)
{
_hits.push_back(*c);
}
} }
++c; ++c;
} }
...@@ -99,7 +89,7 @@ void ChromSweep::ScanCache() { ...@@ -99,7 +89,7 @@ void ChromSweep::ScanCache() {
} }
bool ChromSweep::ChromChange() bool ChromSweep::ChromChange()
{ {
// the files are on the same chrom // the files are on the same chrom
if ((_curr_qy.chrom == _curr_db.chrom) || (_db_status == BED_INVALID) || (_qy_status == BED_INVALID)) { if ((_curr_qy.chrom == _curr_db.chrom) || (_db_status == BED_INVALID) || (_qy_status == BED_INVALID)) {
...@@ -135,6 +125,24 @@ bool ChromSweep::ChromChange() ...@@ -135,6 +125,24 @@ bool ChromSweep::ChromChange()
} }
} }
bool ChromSweep::IsValidHit(const BED &query, const BED &db) {
// do we have an overlap in the DB?
if (overlaps(query.start, query.end, db.start, db.end) > 0) {
// Now test for necessary strandedness.
bool strands_are_same = (query.strand == db.strand);
if ( (_sameStrand == false && _diffStrand == false)
||
(_sameStrand == true && strands_are_same == true)
||
(_diffStrand == true && strands_are_same == false)
)
{
return true;
}
}
return false;
}
bool ChromSweep::Next(pair<BED, vector<BED> > &next) { bool ChromSweep::Next(pair<BED, vector<BED> > &next) {
if (!_bedA->Empty()) { if (!_bedA->Empty()) {
...@@ -145,19 +153,8 @@ bool ChromSweep::Next(pair<BED, vector<BED> > &next) { ...@@ -145,19 +153,8 @@ bool ChromSweep::Next(pair<BED, vector<BED> > &next) {
// 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 (!_bedB->Empty() && _curr_qy.chrom == _curr_db.chrom && !(after(_curr_db, _curr_qy))) while (!_bedB->Empty() && _curr_qy.chrom == _curr_db.chrom && !(after(_curr_db, _curr_qy)))
{ {
// do we have an overlap in the DB? if (IsValidHit(_curr_qy, _curr_db)) {
if (overlaps(_curr_qy.start, _curr_qy.end, _curr_db.start, _curr_db.end) > 0) { _hits.push_back(_curr_db);
// Now test for necessary strandedness.
bool strands_are_same = (_curr_qy.strand == _curr_db.strand);
if ( (_sameStrand == false && _diffStrand == false)
||
(_sameStrand == true && strands_are_same == true)
||
(_diffStrand == true && strands_are_same == false)
)
{
_hits.push_back(_curr_db);
}
} }
_cache.push_back(_curr_db); _cache.push_back(_curr_db);
_db_status = _bedB->GetNextBed(_curr_db, _db_lineNum); _db_status = _bedB->GetNextBed(_curr_db, _db_lineNum);
......
...@@ -83,6 +83,7 @@ private: ...@@ -83,6 +83,7 @@ private:
void ScanCache(); void ScanCache();
bool ChromChange(); bool ChromChange();
bool IsValidHit(const BED &query, const BED &db);
}; };
#endif /* CHROMSWEEP_H */ #endif /* CHROMSWEEP_H */
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment