18 #include "glfHandler.h" 19 #include "BaseQualityHelper.h" 21 char glfHandler::translateBase[16] = {0, 1, 2, 0, 3, 0, 0, 0, 4, 0, 0, 0, 0, 0, 0, 0};
22 char glfHandler::backTranslateBase[5] = { 15, 1, 2, 4, 8 };
23 unsigned char glfHandler::nullLogLikelihoods[10] = {0, 0, 0, 0, 0, 0, 0, 0, 0, 0};
24 double glfHandler::nullLikelihoods[10] = {1., 1., 1., 1., 1., 1., 1., 1., 1., 1.};
26 glfHandler::glfHandler()
31 maxPosition = position = endOfSection = 0;
34 glfHandler::~glfHandler()
41 bool glfHandler::Open(
const String & filename)
44 handle =
ifopen(filename,
"rb");
57 return handle != NULL;
60 void glfHandler::OpenStub()
67 maxPosition = 1999999999;
68 position = maxPosition + 1;
71 bool glfHandler::Create(
const String & filename)
85 return handle != NULL;
88 bool glfHandler::isOpen()
90 return handle != NULL;
93 bool glfHandler::ReadHeader()
103 if (
ifread(handle, magicNumber, 4) != 4)
105 errorMsg =
"unexpected end of file";
109 if (magicNumber[0] !=
'G' || magicNumber[1] !=
'L' || magicNumber[2] !=
'F')
111 errorMsg =
"invalid format";
115 if (magicNumber[3] != 3)
117 errorMsg =
"unsupported version";
121 unsigned int headerLength = 0;
123 if (
ifread(handle, &headerLength, 4) != 4)
125 errorMsg =
"unexpected end of file";
129 if (headerLength > 1024 * 1024)
131 errorMsg =
"header too large -- bailing";
135 header.SetLength(headerLength + 1);
136 header[headerLength] = 0;
138 if (headerLength &&
ifread(handle, header.LockBuffer(headerLength + 1), headerLength) != headerLength)
140 errorMsg =
"unexpected end of file";
147 void glfHandler::Close()
153 void glfHandler::Rewind()
166 bool glfHandler::NextSection()
172 maxPosition = 1999999999;
173 position = maxPosition + 1;
177 while (!endOfSection && !
ifeof(handle))
180 endOfSection =
false;
186 if (
ifread(handle, &labelLength,
sizeof(
int)) ==
sizeof(
int))
188 ifread(handle, label.LockBuffer(labelLength+1), labelLength *
sizeof(char));
189 label.UnlockBuffer();
192 ifread(handle, &maxPosition,
sizeof(
int));
194 return ((maxPosition > 0) && !
ifeof(handle));
200 bool glfHandler::NextBaseEntry()
206 result = NextEntry();
208 while (result && data.recordType == 2);
214 bool glfHandler::NextEntry()
220 if (endOfSection || (
ifread(handle, &data, 1) != 1))
224 position = maxPosition + 1;
230 if (position > maxPosition)
233 switch (data.recordType)
237 position = maxPosition + 1;
240 if (
ifread(handle,((
char *) &data) + 1,
sizeof(data) - 1) ==
sizeof(data) - 1)
242 data.refBase = translateBase[data.refBase];
244 for (
int i = 0; i < 10; i++)
245 likelihoods[i] = bQualityConvertor.toDouble(data.lk[i]);
247 position = position + data.offset;
253 position = maxPosition + 1;
256 while (
ifread(handle, ((
char *) &data) + 1,
sizeof(data) - 4) ==
sizeof(data) - 4)
258 data.refBase = translateBase[data.refBase];
260 for (
int i = 0; i < 3; i++)
261 likelihoods[i] = bQualityConvertor.toDouble(data.indel.lk[i]);
263 position = position + data.offset;
265 indelSequence[0].SetLength(abs(data.indel.length[0]) + 1);
266 indelSequence[0][abs(data.indel.length[0])] = 0;
267 if (
ifread(handle, indelSequence[0].LockBuffer(), abs(data.indel.length[0])) != (
unsigned int) abs(data.indel.length[0]))
270 indelSequence[1].SetLength(abs(data.indel.length[1]) + 1);
271 indelSequence[1][abs(data.indel.length[1])] = 0;
272 if (
ifread(handle, indelSequence[1].LockBuffer(), abs(data.indel.length[1])) != (
unsigned int) abs(data.indel.length[1]))
280 position = maxPosition + 1;
290 recordType = rhs.recordType;
294 for (
int i = 0; i < 10; i++)
303 const double * glfHandler::GetLikelihoods(
int pos)
308 return nullLikelihoods;
311 const unsigned char * glfHandler::GetLogLikelihoods(
int pos)
316 return nullLogLikelihoods;
319 char glfHandler::GetReference(
int pos,
char defaultBase)
327 int glfHandler::GetDepth(
int pos)
335 int glfHandler::GetMapQuality(
int pos)
338 return data.mapQuality;
343 void glfHandler::WriteHeader(
const String & headerText)
345 char magicNumber[4] = {
'G',
'L',
'F', 3};
347 ifwrite(handle, magicNumber, 4);
349 unsigned int headerLength = headerText.Length();
351 ifwrite(handle, &headerLength, 4);
352 ifwrite(handle, (
void *)(
const char *) headerText, headerLength);
355 void glfHandler::BeginSection(
const String & sectionLabel,
int sectionLength)
357 int labelLength = sectionLabel.Length() + 1;
359 ifwrite(handle, &labelLength,
sizeof(
int));
360 ifwrite(handle, (
void *)(
const char *) sectionLabel, labelLength);
361 ifwrite(handle, §ionLength,
sizeof(
int));
363 label = sectionLabel;
364 maxPosition = sectionLength;
367 void glfHandler::EndSection()
371 ifwrite(handle, &marker,
sizeof(
char));
374 void glfHandler::WriteEntry(
int outputPosition)
376 data.offset = outputPosition - position;
377 position = outputPosition;
379 switch (data.recordType)
385 data.refBase = backTranslateBase[data.refBase];
386 ifwrite(handle, &data,
sizeof(data));
387 data.refBase = translateBase[data.refBase];
390 data.refBase = backTranslateBase[data.refBase];
391 ifwrite(handle, &data,
sizeof(data) - 3);
392 data.refBase = translateBase[data.refBase];
394 ifwrite(handle, (
void *)(
const char *) indelSequence[0], abs(data.indel.length[0]));
395 ifwrite(handle, (
void *)(
const char *) indelSequence[1], abs(data.indel.length[1]));
unsigned char refBase
"XACMGRSVTWYHKDBN"[ref_base] gives the reference base
unsigned char mapQuality
root mean squared maximum mapping quality for overlapping reads
unsigned char lk[10]
log10 likelihood ratio * 10 for genotypes AA, AC, AG, AT, CC, CG, CT, GG, GT, TT
unsigned int offset
offset of this record from the previous one, in bases
unsigned depth
log10 minimum likelihood * 10 and the number of mapped reads