Forked from anonymous/get_count_with_sufficient_coverage.cc
Last active
November 4, 2015 15:47
-
-
Save gatoravi/454abf932c0b262058ab to your computer and use it in GitHub Desktop.
Parse 'samtools depth' and obtain lines with sufficient coverage. Assumes two samples in the output. Same as, awk '$3 >= cutoff && $4 >= cutoff'
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
#include <stdint.h> | |
#include <cstdlib> | |
#include <iostream> | |
#include <fstream> | |
#include <sstream> | |
using namespace std; | |
void usage() { | |
cerr << "Usage: ./get_count_with_sufficient_coverage samtools_depth_file cutoff[INT]"; | |
cerr << "Assumes 2 samples."; | |
cerr << endl; | |
return; | |
} | |
int count_number_satisfying_lines(string depth_file, int cutoff) { | |
ifstream fin; | |
fin.open(depth_file.c_str()); | |
string line; | |
string chr; | |
uint32_t pos, s1_readcount, s2_readcount; | |
int lines_satisfying = 0; | |
while(getline(fin, line)) { | |
stringstream ss(line); | |
ss >> chr >> pos >> s1_readcount >> s2_readcount; | |
if(s1_readcount >= cutoff && s2_readcount >= cutoff) { | |
lines_satisfying++; | |
} | |
} | |
return lines_satisfying; | |
} | |
int main(int argc, char* argv[]) { | |
if (argc < 3) { | |
usage(); | |
return 1; | |
} | |
string samtools_depth_file = string(argv[1]); | |
int cutoff = atoi(argv[2]); | |
cerr << "\nCutoff " << cutoff; | |
cerr << "\nSamtools file: " << samtools_depth_file; | |
cerr << "\nNumber of lines satisfying cutoff: "; | |
cout << count_number_satisfying_lines(samtools_depth_file, cutoff); | |
cout << endl; | |
return 0; | |
} |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment