Skip to content

Commit

Permalink
Implement version 8 of seqStore. Working, but not complete and not we…
Browse files Browse the repository at this point in the history 
…ll tested.
  • Loading branch information
@brianwalenz
brianwalenz committed Jul 5, 2019
1 parent 3c35b83 commit 40763d8
Show file tree
Hide file tree
Showing 95 changed files with 3,043 additions and 3,739 deletions.
19 changes: 9 additions & 10 deletions src/bogart/AS_BAT_ReadInfo.C
View file
Original file line number Diff line number Diff line change
Expand Up @@ -32,11 +32,11 @@ ReadInfo::ReadInfo(const char *seqStorePath,
const char *prefix,
uint32 minReadLen) {

sqStore *seqStore = sqStore::sqStore_open(seqStorePath);
sqStore *seqStore = new sqStore(seqStorePath);

_numBases = 0;
_numReads = seqStore->sqStore_getNumReads();
_numLibraries = seqStore->sqStore_getNumLibraries();
_numReads = seqStore->sqStore_lastReadID();
_numLibraries = seqStore->sqStore_lastLibraryID();

_readStatus = new ReadStatus [_numReads + 1];

Expand All @@ -53,24 +53,23 @@ ReadInfo::ReadInfo(const char *seqStorePath,
uint32 numLoaded = 0;

for (uint32 fi=1; fi<=_numReads; fi++) {
sqRead *read = seqStore->sqStore_getRead(fi);
uint32 iid = read->sqRead_readID();
uint32 len = read->sqRead_sequenceLength();
uint32 len = seqStore->sqStore_getReadLength(fi);

if (len < minReadLen) {
if ((len < minReadLen) ||
(seqStore->sqStore_isIgnoredRead(fi))) {
numSkipped++;
continue;
}

_numBases += len;

_readStatus[iid].readLength = len;
_readStatus[iid].libraryID = read->sqRead_libraryID();
_readStatus[fi].readLength = len;
_readStatus[fi].libraryID = seqStore->sqStore_getLibraryIDForRead(fi);

numLoaded++;
}

seqStore->sqStore_close();
delete seqStore;

if (minReadLen > 0)
writeStatus("ReadInfo()-- Using %d reads, ignoring %u reads less than " F_U32 " bp long.\n",
Expand Down
6 changes: 3 additions & 3 deletions src/bogart/addReadsToUnitigs.C
View file
Original file line number Diff line number Diff line change
Expand Up @@ -386,7 +386,7 @@ main(int argc, char **argv) {
//

fprintf(stderr, "Processing mate pairs, updating seqStore.\n");
seqStore = sqStore::sqStore_open(seqName, false, true); // last arg - TRUE - writable
seqStore = new sqStore(seqName, false, true); // last arg - TRUE - writable

uint32 unpaired = 0;
uint32 multiple = 0;
Expand Down Expand Up @@ -458,7 +458,7 @@ main(int argc, char **argv) {
pairsToDiff++;
}

seqStore->sqStore_close();
delete seqStore;

fprintf(stderr, "Will NOT add %u pairs - one read failed to map.\n", unpaired);
fprintf(stderr, "Will NOT add %u pairs - multiple mappings.\n", multiple);
Expand All @@ -469,7 +469,7 @@ main(int argc, char **argv) {
// Open stores. seqStore cannot be opened for writing, because then we can't loadall.
//

seqStore = sqStore::sqStore_open(seqName, false, false); // last arg - false - not writable
seqStore = new sqStore(seqName, false, false); // last arg - false - not writable
tigStore = new MultiAlignStore(tigName, tigVers, 0, 0, true, true, false); // Write back to the same version

if (loadall) {
Expand Down
4 changes: 2 additions & 2 deletions src/bogart/analyzeBest.C
View file
Original file line number Diff line number Diff line change
Expand Up @@ -86,7 +86,7 @@ main(int argc, char **argv) {

fprintf(stderr, "Loading read to library mapping.\n");

sqStore *seq = sqStore::sqStore_open(seqName, false, false);
sqStore *seq = new sqStore(seqName, false, false);
gkStream *str = new gkStream(seq, 0, 0, GKFRAGMENT_INF);
sqRead fr;

Expand Down Expand Up @@ -193,7 +193,7 @@ main(int argc, char **argv) {
dovePerLib[i], (dovePerLib[i] == 0) ? 0.0 : 100.0 * dovePerLib[i] / tot);
}

seq->sqStore_close();
delete seq;

AS_UTL_closeFile(be, bEdge);
AS_UTL_closeFile(bc, bCont);
Expand Down
54 changes: 24 additions & 30 deletions src/correction/falconsense.C
View file
Original file line number Diff line number Diff line change
Expand Up @@ -83,20 +83,19 @@ loadReadList(char *readListName, uint32 iidMin, uint32 iidMax, set<uint32> &read



sqReadData *
sqRead *
loadReadData(uint32 readID,
sqStore *seqStore,
map<uint32, sqRead *> &reads,
map<uint32, sqReadData *> &datas) {
map<uint32, sqRead *> &reads) {

if (datas.count(readID) == 0) {
datas[readID] = new sqReadData;
seqStore->sqStore_loadReadData(readID, datas[readID]);
if (reads.count(readID) == 0) {
reads[readID] = new sqRead;
seqStore->sqStore_getRead(readID, reads[readID]);
}

assert(datas[readID] != NULL);
assert(reads[readID] != NULL);

return(datas[readID]);
return(reads[readID]);
}


Expand All @@ -106,7 +105,6 @@ generateFalconConsensus(falconConsensus *fc,
tgTig *layout,
sqCache *seqCache,
map<uint32, sqRead *> &reads,
map<uint32, sqReadData *> &datas,
bool trimToAlign,
uint32 minOlapLength) {

Expand Down Expand Up @@ -232,16 +230,12 @@ generateFalconConsensus(falconConsensus *fc,

;

// Clean up. Remvoe all the reads[] and datas[] we've loaded.
// Clean up. Remvoe all the reads[] we've loaded.

for (map<uint32, sqRead *>::iterator it=reads.begin(); it != reads.end(); ++it)
delete it->second;

for (map<uint32, sqReadData *>::iterator it=datas.begin(); it != datas.end(); ++it)
delete it->second;

reads.clear();
datas.clear();

delete fd;
delete [] evidence;
Expand Down Expand Up @@ -426,7 +420,7 @@ main(int argc, char **argv) {

if (seqName) {
fprintf(stderr, "-- Opening seqStore '%s'.\n", seqName);
seqStore = sqStore::sqStore_open(seqName);
seqStore = new sqStore(seqName);
seqCache = new sqCache(seqStore, sqRead_raw);
}

Expand All @@ -436,24 +430,24 @@ main(int argc, char **argv) {
}

if ((seqStore) &&
(seqStore->sqStore_getNumReads() < idMax)) // Limit the range of processing to the
idMax = seqStore->sqStore_getNumReads(); // number of reads in the store.
(seqStore->sqStore_lastReadID() < idMax)) // Limit the range of processing to the
idMax = seqStore->sqStore_lastReadID(); // number of reads in the store.

loadReadList(readListName, idMin, idMax, readList); // Further limit to a set of good reads.

// Open any import or export files.

FILE *exportFile = NULL;
FILE *importFile = NULL;
writeBuffer *exportFile = NULL;
readBuffer *importFile = NULL;

if (exportName) {
fprintf(stderr, "-- Opening export file '%s'.\n", exportName);
exportFile = AS_UTL_openOutputFile(exportName);
exportFile = new writeBuffer(exportName, "w");
}

if (importName) {
fprintf(stderr, "-- Opening import file '%s'.\n", importName);
importFile = AS_UTL_openInputFile(importName);
importFile = new readBuffer(importName);
}

// Open logging and summary files
Expand All @@ -476,7 +470,6 @@ main(int argc, char **argv) {

falconConsensus *fc = new falconConsensus(minOutputCoverage, minOutputLength, minOlapIdentity, minOlapLength, restrictToOverlap);
map<uint32, sqRead *> reads;
map<uint32, sqReadData *> datas;

if (memoryLimit == 0) {
fprintf(stdout, " read read evidence corrected\n");
Expand All @@ -494,12 +487,11 @@ main(int argc, char **argv) {
FILE *importedLayouts = AS_UTL_openOutputFile(importName, '.', "layout", (importName != NULL));
FILE *importedReads = AS_UTL_openOutputFile(importName, '.', "fasta", (importName != NULL));

while (layout->importData(importFile, reads, datas, NULL, NULL) == true) {
while (layout->importData(importFile, reads, NULL, NULL) == true) {
generateFalconConsensus(fc,
layout,
seqCache,
reads,
datas,
trimToAlign,
minOlapLength);

Expand Down Expand Up @@ -552,7 +544,7 @@ main(int argc, char **argv) {
//

else if (memoryLimit > 0) {
uint32 lastID = seqStore->sqStore_getNumReads();
uint32 lastID = seqStore->sqStore_lastReadID();
uint32 *readLens = new uint32 [lastID + 1];
uint32 *readRefs = new uint32 [lastID + 1];

Expand All @@ -565,7 +557,7 @@ main(int argc, char **argv) {
// cacheEntry - storage internal to the cache.

for (uint32 ii=0; ii <= lastID; ii++) {
readLens[ii] = 12 + seqStore->sqStore_getRead(ii)->sqRead_sequenceLength(sqRead_raw) / 4 + 4 + sizeof(sqCacheEntry); // Round up, and 3 extra uint32.
readLens[ii] = 12 + seqStore->sqStore_getReadLength(ii, sqRead_raw) / 4 + 4 + sizeof(sqCacheEntry); // Round up, and 3 extra uint32.
readRefs[ii] = 0;
}

Expand Down Expand Up @@ -596,6 +588,9 @@ main(int argc, char **argv) {
exit(1);
}

//fprintf(stderr, "readsPerBatch %u\n", readsPerBatch);
//fprintf(stderr, "memoryLimit %f GB\n", memoryLimit / 1024.0 / 1024.0 / 1024.0);

fprintf(batFile, "batch bgnID endID nReads memory (base memory %.3f GB)\n", memUsedBase / 1024.0 / 1024.0 / 1024.0);
fprintf(batFile, "----- --------- --------- ------- -------\n");

Expand Down Expand Up @@ -708,7 +703,6 @@ main(int argc, char **argv) {
layout,
seqCache,
reads,
datas,
trimToAlign,
minOlapLength);

Expand All @@ -735,15 +729,15 @@ main(int argc, char **argv) {
AS_UTL_closeFile(seqFile);
AS_UTL_closeFile(batFile);

AS_UTL_closeFile(exportFile);
AS_UTL_closeFile(importFile);
delete exportFile;
delete importFile;

delete fc;
delete corStore;

delete seqCache;

seqStore->sqStore_close();
delete seqStore;

fprintf(stderr, "\n");
fprintf(stderr, "Bye.\n");
Expand Down
4 changes: 2 additions & 2 deletions src/correction/filterCorrectionLayouts.C
View file
Original file line number Diff line number Diff line change
Expand Up @@ -444,12 +444,12 @@ main(int argc, char **argv) {

sqRead_setDefaultVersion(sqRead_raw);

sqStore *seqStore = sqStore::sqStore_open(seqStoreName);
sqStore *seqStore = new sqStore(seqStoreName);
tgStore *corStore = new tgStore(corStoreName, 1);

falconConsensus *fc = new falconConsensus(minOutputCoverage, minOutputLength, 0, 0); // For memory estimtes

uint32 numReads = seqStore->sqStore_getNumReads();
uint32 numReads = seqStore->sqStore_lastReadID();

readStatus *status = new readStatus [numReads + 1];

Expand Down
10 changes: 5 additions & 5 deletions src/correction/filterCorrectionOverlaps.C
View file
Original file line number Diff line number Diff line change
Expand Up @@ -197,7 +197,7 @@ main(int argc, char **argv) {

sqRead_setDefaultVersion(sqRead_raw);

sqStore *seqStore = sqStore::sqStore_open(seqStoreName);
sqStore *seqStore = new sqStore(seqStoreName);

ovStore *ovlStore = new ovStore(ovlStoreName, seqStore);
ovStoreHistogram *ovlHisto = ovlStore->getHistogram();
Expand All @@ -208,7 +208,7 @@ main(int argc, char **argv) {
uint32 ovlMax = 0;
ovOverlap *ovl = NULL;

uint16 *scores = new uint16 [seqStore->sqStore_getNumReads() + 1];
uint16 *scores = new uint16 [seqStore->sqStore_lastReadID() + 1];
uint16 scoreExact = 0;
uint16 scoreEstim = 0;

Expand All @@ -227,7 +227,7 @@ main(int argc, char **argv) {
//fprintf(stdout, "-------- ------ ------\n");
}

for (uint32 id=0; id <= seqStore->sqStore_getNumReads(); id++) {
for (uint32 id=0; id <= seqStore->sqStore_lastReadID(); id++) {
scores[id] = UINT16_MAX;

if (numOlaps[id] == 0) {
Expand Down Expand Up @@ -258,7 +258,7 @@ main(int argc, char **argv) {
}

if (scoreFile)
writeToFile(scores, "scores", seqStore->sqStore_getNumReads() + 1, scoreFile);
writeToFile(scores, "scores", seqStore->sqStore_lastReadID() + 1, scoreFile);

AS_UTL_closeFile(scoreFile, scoreFileName);
AS_UTL_closeFile(logFile, logFileName);
Expand All @@ -270,7 +270,7 @@ main(int argc, char **argv) {
delete ovlHisto;
delete ovlStore;

seqStore->sqStore_close();
delete seqStore;

if (noStats == true)
exit(0);
Expand Down
10 changes: 5 additions & 5 deletions src/correction/generateCorrectionLayouts.C
View file
Original file line number Diff line number Diff line change
Expand Up @@ -62,7 +62,7 @@ loadThresholds(sqStore *seqStore,
char *scoreName,
uint32 expectedCoverage,
FILE *scoFile) {
uint32 numReads = seqStore->sqStore_getNumReads();
uint32 numReads = seqStore->sqStore_lastReadID();
uint16 *olapThresh = new uint16 [numReads + 1];

if (scoreName != NULL)
Expand Down Expand Up @@ -303,11 +303,11 @@ main(int argc, char **argv) {

sqRead_setDefaultVersion(sqRead_raw);

sqStore *seqStore = sqStore::sqStore_open(seqName);
sqStore *seqStore = new sqStore(seqName);
ovStore *ovlStore = new ovStore(ovlName, seqStore);
tgStore *corStore = new tgStore(corName);

uint32 numReads = seqStore->sqStore_getNumReads();
uint32 numReads = seqStore->sqStore_lastReadID();

// Threshold the range of reads to operate on.

Expand Down Expand Up @@ -346,7 +346,7 @@ main(int argc, char **argv) {
tgTig *layout = new tgTig;

layout->_tigID = rr;
layout->_layoutLen = seqStore->sqStore_getRead(rr)->sqRead_sequenceLength(sqRead_raw);
layout->_layoutLen = seqStore->sqStore_getReadLength(rr, sqRead_raw);

generateLayout(layout,
olapThresh,
Expand All @@ -369,7 +369,7 @@ main(int argc, char **argv) {
delete corStore;
delete ovlStore;

seqStore->sqStore_close();
delete seqStore;

fprintf(stderr, "Bye.\n");

Expand Down
Loading

0 comments on commit 40763d8

Please sign in to comment.