Skip to content
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
85 changes: 85 additions & 0 deletions src/fastqchunkindex.cpp
Original file line number Diff line number Diff line change
@@ -0,0 +1,85 @@
#include "fastqchunkindex.h"
#include <fcntl.h>
#include <unistd.h>
#include <sys/stat.h>
#include <cstring>

// Scan buffer size: 4MB
static const size_t SCAN_BUF_SIZE = 1 << 22;

bool FastqChunkIndex::build(const std::string& filename, int packSize) {
mFd = open(filename.c_str(), O_RDONLY);
if (mFd < 0)
return false;

struct stat st;
if (fstat(mFd, &st) != 0 || !S_ISREG(st.st_mode)) {
::close(mFd);
mFd = -1;
return false;
}
mFileSize = st.st_size;
if (mFileSize == 0) {
::close(mFd);
mFd = -1;
return false;
}

#ifdef __linux__
posix_fadvise(mFd, 0, mFileSize, POSIX_FADV_SEQUENTIAL);
#endif

// Every FASTQ record is 4 lines, so one pack = packSize * 4 newlines
const int newlinesPerPack = packSize * 4;

mOffsets.clear();
mOffsets.push_back(0); // first pack starts at byte 0

char* buf = new char[SCAN_BUF_SIZE];
size_t fileOffset = 0;
int nlCount = 0;

while (fileOffset < mFileSize) {
size_t toRead = SCAN_BUF_SIZE;
if (fileOffset + toRead > mFileSize)
toRead = mFileSize - fileOffset;

ssize_t n = pread(mFd, buf, toRead, fileOffset);
if (n <= 0)
break;

// Scan for newlines using memchr for speed
const char* p = buf;
const char* end = buf + n;
while (p < end) {
const char* nl = (const char*)memchr(p, '\n', end - p);
if (!nl)
break;
nlCount++;
if (nlCount == newlinesPerPack) {
// The next pack starts right after this newline
size_t offset = fileOffset + (nl - buf) + 1;
if (offset < mFileSize)
mOffsets.push_back(offset);
nlCount = 0;
}
p = nl + 1;
}

fileOffset += n;
}

// Sentinel: the end of the last pack is fileSize
if (mOffsets.back() != mFileSize)
mOffsets.push_back(mFileSize);

delete[] buf;
return true;
}

FastqChunkIndex::~FastqChunkIndex() {
if (mFd >= 0) {
::close(mFd);
mFd = -1;
}
}
41 changes: 41 additions & 0 deletions src/fastqchunkindex.h
Original file line number Diff line number Diff line change
@@ -0,0 +1,41 @@
#ifndef FASTQ_CHUNK_INDEX_H
#define FASTQ_CHUNK_INDEX_H

#include <vector>
#include <string>
#include <cstddef>

// Scans an uncompressed FASTQ file and builds a pack offset index.
// Each entry marks the byte offset where a pack of PACK_SIZE reads begins.
// The last pack may contain fewer than PACK_SIZE reads.
class FastqChunkIndex {
public:
// Scan the file, recording a byte offset every (4 * packSize) newlines.
// Returns false if the file cannot be opened or is not a regular file.
bool build(const std::string& filename, int packSize);

// Number of packs
size_t packCount() const { return mOffsets.size() > 0 ? mOffsets.size() - 1 : 0; }

// Byte offset of pack i
size_t packStart(size_t i) const { return mOffsets[i]; }

// Byte length of pack i
size_t packLength(size_t i) const { return mOffsets[i + 1] - mOffsets[i]; }

// File descriptor (caller must not close it; owned by this object)
int fd() const { return mFd; }

// Total file size
size_t fileSize() const { return mFileSize; }

// Clean up
~FastqChunkIndex();

private:
std::vector<size_t> mOffsets; // packCount+1 entries: [0]=0, [last]=fileSize
int mFd = -1;
size_t mFileSize = 0;
};

#endif
117 changes: 117 additions & 0 deletions src/fastqchunkparser.cpp
Original file line number Diff line number Diff line change
@@ -0,0 +1,117 @@
#include "fastqchunkparser.h"
#include "common.h"
#include <cstring>
#include <vector>
#include <iostream>

using namespace std;

// Find next newline in buf[pos..len), return position of char after \n.
// Handles \r\n. Returns len if no newline found.
static inline size_t nextLine(const char* buf, size_t pos, size_t len) {
const char* p = (const char*)memchr(buf + pos, '\n', len - pos);
if (!p)
return len;
return (p - buf) + 1;
}

// Extract a line from buf[start..end) stripping trailing \r\n
static inline void extractLine(const char* buf, size_t start, size_t lineEnd, string* out) {
size_t end = lineEnd;
// strip trailing \n and \r
while (end > start && (buf[end - 1] == '\n' || buf[end - 1] == '\r'))
end--;
out->assign(buf + start, end - start);
}

ReadPack* FastqChunkParser::parse(const char* buf, size_t len, int tid, ReadPool* pool, bool phred64) {
// Pre-count records: count newlines / 4
int estimatedReads = 0;
{
const char* p = buf;
const char* end = buf + len;
while (p < end) {
p = (const char*)memchr(p, '\n', end - p);
if (!p) break;
estimatedReads++;
p++;
}
estimatedReads /= 4;
}
if (estimatedReads <= 0)
estimatedReads = 1;

Read** data = new Read*[estimatedReads];
int count = 0;
size_t pos = 0;

while (pos < len) {
// Skip empty lines
while (pos < len && (buf[pos] == '\n' || buf[pos] == '\r'))
pos++;
if (pos >= len)
break;

// Line 1: name (starts with @)
size_t nameStart = pos;
size_t nameEnd = nextLine(buf, pos, len);
if (nameStart >= len || buf[nameStart] != '@')
break;
pos = nameEnd;

// Line 2: sequence
if (pos >= len) break;
size_t seqStart = pos;
size_t seqEnd = nextLine(buf, pos, len);
pos = seqEnd;

// Line 3: strand (+)
if (pos >= len) break;
size_t strandStart = pos;
size_t strandEnd = nextLine(buf, pos, len);
pos = strandEnd;

// Line 4: quality
if (pos >= len && strandEnd >= len) break;
size_t qualStart = pos;
size_t qualEnd = nextLine(buf, pos, len);
pos = qualEnd;

// Build Read object
Read* r = nullptr;
if (pool)
r = pool->getOne();

if (r) {
extractLine(buf, nameStart, nameEnd, r->mName);
extractLine(buf, seqStart, seqEnd, r->mSeq);
extractLine(buf, strandStart, strandEnd, r->mStrand);
extractLine(buf, qualStart, qualEnd, r->mQuality);
} else {
string* name = new string();
string* seq = new string();
string* strand = new string();
string* quality = new string();
extractLine(buf, nameStart, nameEnd, name);
extractLine(buf, seqStart, seqEnd, seq);
extractLine(buf, strandStart, strandEnd, strand);
extractLine(buf, qualStart, qualEnd, quality);
r = new Read(name, seq, strand, quality, phred64);
}

if (count >= estimatedReads) {
// Shouldn't happen, but be safe
Read** newData = new Read*[estimatedReads * 2];
memcpy(newData, data, sizeof(Read*) * count);
delete[] data;
data = newData;
estimatedReads *= 2;
}
data[count++] = r;
}

ReadPack* pack = new ReadPack;
pack->data = data;
pack->count = count;
return pack;
}
17 changes: 17 additions & 0 deletions src/fastqchunkparser.h
Original file line number Diff line number Diff line change
@@ -0,0 +1,17 @@
#ifndef FASTQ_CHUNK_PARSER_H
#define FASTQ_CHUNK_PARSER_H

#include "read.h"
#include "readpool.h"

// Parse a raw byte buffer (from pread) into a ReadPack.
// The buffer must start at a record boundary and contain complete records.
class FastqChunkParser {
public:
// Parse buf[0..len) into Read* objects, return a ReadPack.
// If pool is non-null, try to reuse Read objects from the pool.
// tid is the thread id for the pool.
static ReadPack* parse(const char* buf, size_t len, int tid, ReadPool* pool, bool phred64);
};

#endif
Loading
Loading