Commit a5c240d8 authored by Paul Igor Costea's avatar Paul Igor Costea
Browse files

Refactored SNP calling/determination code for clarity.

parent b101615c
......@@ -587,26 +587,31 @@ std::map<std::string, char> CodonMap(codons, codons + sizeof codons / sizeof cod
std::string oldCodon = "";
bool write = false;
for (int i=0; i<snps.length(); ++i) {
//Skip same base
if (snps[i] == base) continue;
std::ostringstream internal;
char check[3];
sprintf(check,"%c%c",snps[i],toupper(snps[i]));
snpCount = getSum(map,check);
bool writeThis = false;
std::string* sToWrite = NULL;
if (snpCount >= 2*CALLING_THRESHOLD) {//Could be that this is a sample specific SNP.
if (isPersonalSNP(map,nrSamples,check)) {
writeThis = true;
sToWrite = &indiv;
}
}
if ((snpCount >= CALLING_THRESHOLD) && (snpCount >= cov*CALLING_MIN_FRACTION)) {
write = true;
writeThis = true;
sToWrite = &s;
}
//Skip same base
if (snps[i] == base) continue;
std::ostringstream internal;
char check[3];
sprintf(check,"%c%c",snps[i],toupper(snps[i]));
snpCount = getSum(map,check);//How many times do we see this position as a variant in the entire population?
bool writeThis = false;
std::string* sToWrite = NULL;
if ((snpCount >= CALLING_THRESHOLD) && (snpCount >= cov*CALLING_MIN_FRACTION)) { // This is a "common" variant
write = true;
writeThis = true;
sToWrite = &s;
} else { //May be a individual variant
for (int i=1; i<=nrSamples; ++i) {
int s = getSum(map,snp,i);
if (s >= CALLING_THRESHOLD) {//We observe the variant "CALLING_THRESHOLD" times in one sample.
writeThis = true;
sToWrite = &indiv;
break;
}
}
}
if (writeThis) {
if ((hasGenes) && (isInGene)) {
//TODO: consider multiple genes!
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment