-
Notifications
You must be signed in to change notification settings - Fork 2
/
Copy pathDstatCounter.h
99 lines (77 loc) · 2.17 KB
/
DstatCounter.h
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
/*
* DstatCounter
* Date: Oct-15-2012
* Author : Gabriel Renaud [email protected]
*
*/
#ifndef DstatCounter_h
#define DstatCounter_h
#include <iostream>
#include <ostream>
#include <vector>
#include <math.h>
#include <sys/time.h>
#include <stdlib.h>
#include "libgab.h"
//#include <iostream>
//#include <fstream>
using namespace std;
class DstatCounter{
private:
public:
//counter of mutations
unsigned int counterAncAnc;
unsigned int counterAncDer;
unsigned int counterDerAnc;
unsigned int counterDerDer;
DstatCounter(); //constructor
// vector<double> performBoostraps() const;
void reinitializedCounters();
string headerForCount() const;
double computeDST() const;
DstatCounter & operator+=(const DstatCounter & other);
DstatCounter & operator-=(const DstatCounter & other);
friend ostream& operator<<(ostream& os, const DstatCounter & ct){
/* vector<double> boostraps = ct.performBoostraps(); */
double dst= (double(ct.counterAncDer)-double(ct.counterDerAnc))/(double(ct.counterAncDer)+double(ct.counterDerAnc));
os<<ct.counterAncAnc<<"\t"
<<ct.counterAncDer<<"\t"
<<ct.counterDerAnc<<"\t"
<<ct.counterDerDer<<"\t"
//(ADDA - DADA )/ (ADDA+DADA )
<<dst<<"\t" ;
/* pair<double,double> res = computeMeanSTDDEV(boostraps); */
/* os<< ((res.first-dst)/res.second); */
return os;
}
friend istream &operator>>(istream &is , DstatCounter &ct){
//cerr<<"DstatCounter in"<<endl;
double dst;
//istringstream dst;
is>>
ct.counterAncAnc>>
ct.counterAncDer>>
ct.counterDerAnc>>
ct.counterDerDer;
//get dst but will be discarded
string dstStr;
char c;
is.get(c);
if(c != '\t'){
cerr<<"A single dstat should have 5 fields AA\tAD\tDA\tDD\tdst"<<endl;
}
while (is.get(c)){ // loop getting single characters
if(c == '\t') break;
dstStr+=c;
//cerr<<"#"<<c<<"#"<<endl;
}
if(dstStr == "-nan"){
dst = -1.0*numeric_limits<double>::quiet_NaN();
}else{
dst = destringify<double>(dstStr);
}
//cerr<<ct.counterAncAnc<<"\t"<<ct.counterAncDer<<"\t"<<ct.counterDerAnc<<"\t"<<ct.counterDerDer<<"\t"<<dst<<endl;;
return is;
}
};
#endif