#include "NonCall.h" NonCall NonCalls::aggregateAll(void) { NonCall aggregate; bool first = true; for (NonCalls::const_iterator nc = this->begin(); nc != this->end(); ++nc) { for (map >::const_iterator p = nc->second.begin(); p != nc->second.end(); ++p) { for (map::const_iterator s = p->second.begin(); s != p->second.end(); ++s) { const NonCall& nonCall = s->second; aggregate.refCount += nonCall.refCount; aggregate.altCount += nonCall.altCount; aggregate.reflnQ += nonCall.reflnQ; aggregate.altlnQ += nonCall.altlnQ; aggregate.nCount += 1; if (first) { aggregate.minDepth = nonCall.refCount + nonCall.altCount; first = false; } else { aggregate.minDepth = min(aggregate.minDepth, nonCall.refCount + nonCall.altCount); } } } } return aggregate; } void NonCalls::aggregatePerSample(map& perSample) { set seen; for (NonCalls::const_iterator nc = this->begin(); nc != this->end(); ++nc) { for (map >::const_iterator p = nc->second.begin(); p != nc->second.end(); ++p) { for (map::const_iterator s = p->second.begin(); s != p->second.end(); ++s) { const string& name = s->first; const NonCall& nonCall = s->second; NonCall& aggregate = perSample[name]; aggregate.refCount += nonCall.refCount; aggregate.altCount += nonCall.altCount; aggregate.reflnQ += nonCall.reflnQ; aggregate.altlnQ += nonCall.altlnQ; aggregate.nCount += 1; if (!seen.count(name)) { aggregate.minDepth = nonCall.refCount + nonCall.altCount; seen.insert(name); } else { aggregate.minDepth = min(aggregate.minDepth, nonCall.refCount + nonCall.altCount); } } } } } void NonCalls::record(const string& seqName, long pos, const Samples& samples) { map& site = (*this)[seqName][pos]; for (Samples::const_iterator s = samples.begin(); s != samples.end(); ++s) { // tally ref and non-ref alleles const string& name = s->first; const Sample& sample = s->second; NonCall& noncall = site[name]; for (Sample::const_iterator a = sample.begin(); a != sample.end(); ++a) { const vector& alleles = a->second; for (vector::const_iterator o = alleles.begin(); o != alleles.end(); ++o) { Allele& allele = **o; if (allele.isReference()) { ++noncall.refCount; noncall.reflnQ += allele.lnquality; } else { ++noncall.altCount; noncall.altlnQ += allele.lnquality; } } } } } pair NonCalls::firstPos(void) { const string& startChrom = begin()->first; long startPos = begin()->second.begin()->first; return make_pair(startChrom, startPos); } pair NonCalls::lastPos(void) { const string& endChrom = rbegin()->first; long endPos = rbegin()->second.rbegin()->first; return make_pair(endChrom, endPos); }