18 #include "SamStatistics.h"
22 SamStatistics::SamStatistics()
27 SamStatistics::~SamStatistics()
32 void SamStatistics::reset()
35 myMappedReadCount = 0;
36 myPairedReadCount = 0;
37 myProperPairedReadCount = 0;
39 myMappedReadBases = 0;
41 myQCFailureReadCount = 0;
45 bool SamStatistics::updateStatistics(
SamRecord& samRecord)
54 uint16_t flag = samRecord.
getFlag();
57 if(SamFlag::isMapped(flag))
60 myMappedReadBases += readLen;
62 if(SamFlag::isPaired(flag))
65 if(SamFlag::isProperPair(flag))
67 ++myProperPairedReadCount;
70 if(SamFlag::isDuplicate(flag))
74 if(SamFlag::isQCFailure(flag))
76 ++myQCFailureReadCount;
80 myBaseCount += readLen;
86 void SamStatistics::print()
88 double DIVIDE_UNITS = 1000000;
89 std::string units =
"(e6)";
91 std::cerr << std::fixed << std::setprecision(2);
94 if(myReadCount < DIVIDE_UNITS)
101 std::cerr <<
"TotalReads" << units <<
"\t"
102 << myReadCount/DIVIDE_UNITS << std::endl;
103 std::cerr <<
"MappedReads" << units <<
"\t"
104 << myMappedReadCount/DIVIDE_UNITS << std::endl;
105 std::cerr <<
"PairedReads" << units <<
"\t"
106 << myPairedReadCount/DIVIDE_UNITS << std::endl;
107 std::cerr <<
"ProperPair" << units <<
"\t"
108 << myProperPairedReadCount/DIVIDE_UNITS << std::endl;
109 std::cerr <<
"DuplicateReads" << units <<
"\t"
110 << myDupReadCount/DIVIDE_UNITS << std::endl;
111 std::cerr <<
"QCFailureReads" << units <<
"\t"
112 << myQCFailureReadCount/DIVIDE_UNITS << std::endl;
113 std::cerr << std::endl;
116 std::cerr <<
"MappingRate(%)\t"
117 << 100 * myMappedReadCount/(double)myReadCount << std::endl;
118 std::cerr <<
"PairedReads(%)\t"
119 << 100 * myPairedReadCount/(double)myReadCount << std::endl;
120 std::cerr <<
"ProperPair(%)\t"
121 << 100 * myProperPairedReadCount/(double)myReadCount << std::endl;
122 std::cerr <<
"DupRate(%)\t"
123 << 100 * myDupReadCount/(double)myReadCount << std::endl;
124 std::cerr <<
"QCFailRate(%)\t"
125 << 100 * myQCFailureReadCount/(double)myReadCount << std::endl;
126 std::cerr << std::endl;
129 std::cerr <<
"TotalBases" << units <<
"\t"
130 << myBaseCount/DIVIDE_UNITS << std::endl;
131 std::cerr <<
"BasesInMappedReads" << units <<
"\t"
132 << myMappedReadBases/DIVIDE_UNITS << std::endl;
Class providing an easy to use interface to get/set/operate on the fields in a SAM/BAM record.
uint16_t getFlag()
Get the flag (FLAG).
int32_t getReadLength()
Get the length of the read.