· 9 years ago · Oct 31, 2016, 06:22 PM
1/*
2 * parse
3 * Date: Oct-31-2016
4 * Author : Gabriel Renaud gabriel.reno [at sign here ] gmail.com
5 *
6 */
7
8#include <iostream>
9#include <fstream>
10#include <gzstream.h>
11
12#include "utils.h"
13
14using namespace std;
15
16template <typename T1, typename T2>
17struct less_second {
18 typedef pair<T1, T2> type;
19 bool operator ()(type const& a, type const& b) const {
20 return a.second > b.second;
21 }
22};
23
24string return3rdstring(const string & toparse){
25 int field=0;
26 bool inWhite=false;
27 string temp="";
28
29 for(unsigned int i=0;i<toparse.size();i++){
30 if(inWhite){
31 if(toparse[i] == ' '){
32
33 }else{
34 temp = toparse[i];
35 inWhite=false;
36 }
37 }else{//not in white
38 if(toparse[i] == ' '){
39 field++;
40 if(field == 4){
41 return temp;
42 }
43 inWhite=true;
44 }else{
45 temp+=toparse[i];
46 }
47 }
48 }
49
50 cerr<<"ERROR"<<endl;
51 exit(1);
52}
53
54
55int main (int argc, char *argv[]) {
56
57 string line;
58 igzstream myFile;
59 string filename = string(argv[1]);
60 double logcutoff = destringify<double>(string(argv[2]));
61 unsigned int lengthCutoff = destringify<unsigned int >(string(argv[3]));
62
63 myFile.open(filename.c_str(), ios::in);
64 bool blockPassed=false;
65 bool blockFirst =true ;
66 map<string,int> chr2count;
67 map<string,unsigned int> chr2AlLength;
68
69 string lastChr="";
70 string lastScaff="";
71
72 if (myFile.good()){
73 while ( getline (myFile,line)){
74 //cout<<line<<endl;
75 if(line[0] == '#'){
76 if(line.size()>=7){
77 //cout<<">"<<line.substr(0,7)<<"<"<<endl;
78 if(line.substr(0,7) == "# batch"){
79 //return 1;
80
81 //dump
82 if(!chr2count.empty()){
83
84 // for (map<string,int>::iterator it=chr2count.begin(); it!=chr2count.end(); ++it)
85 // cout << lastScaff<<"\t"<<it->first << "\t" << it->second << endl;
86
87 // vector<pair<string, int> > chr2countCpy(chr2count.begin(), chr2count.end());
88 // sort(chr2countCpy.begin(), chr2countCpy.end(), less_second<string, int>());
89 // for(unsigned int i=0;i<chr2countCpy.size();i++){
90 // cout<<chr2countCpy[i].first<<"\t"<<chr2countCpy[i].second<<endl;
91 // }
92 vector<pair<string, unsigned int> > chr2AlLengthCpy(chr2AlLength.begin(), chr2AlLength.end());
93 sort(chr2AlLengthCpy.begin(), chr2AlLengthCpy.end(), less_second<string, unsigned int>());
94 for(unsigned int i=0;i<chr2AlLengthCpy.size();i++){
95 cout<<lastScaff<<"\t"<<chr2AlLengthCpy[i].first<<"\t"<<chr2AlLengthCpy[i].second<<"\t"<<chr2count[chr2AlLengthCpy[i].first]<<endl;
96 }
97
98
99 }
100 //re-init
101 chr2count.clear();
102 chr2AlLength.clear();
103
104 }
105 }
106 }
107
108 if(line[0] == 'a'){
109 blockPassed=false;
110 blockFirst =true;
111
112 vector<string> tokens = allTokens(line,' ');
113 string eval = tokens[3];
114
115 if(eval[0] != 'E' ||
116 eval[1] != '=' ){
117 cerr << "Unable to open file "<<filename<<endl;
118 return 1;
119 }
120
121 eval = eval.substr(2);
122 //cout<<eval<<endl;
123
124 string lengthScore = tokens[1];
125
126
127
128 double eVal = destringify<double>(eval);
129 if(eVal==0){
130 blockPassed=true;
131 }else{
132 double logeVal = -1.0*log(eVal);
133 //cout<<eVal<<"\t"<<log(eVal)<<endl;
134 if(logeVal > logcutoff){
135 blockPassed=true;
136 }
137 }
138 //cout<<endl<<line<<" "<<blockPassed<<endl;
139
140 }
141
142 if(line[0] == 's'){
143 if( blockFirst ){
144 //cout<<blockPassed<<"schr "<<line<<endl;
145 if(blockPassed){
146 vector<string> tokens = allTokens(line,' ');
147 lastChr = tokens[1];
148
149 string lengthAl = return3rdstring(line);
150 unsigned int lengthAlUI= destringify<unsigned int>(lengthAl);
151
152
153 if(lengthAlUI>=lengthCutoff){
154
155 chr2count[lastChr]++;
156 chr2AlLength[lastChr] += lengthAlUI;
157 //cout<<"l=\t"<<lengthAl<<"\t"<<lastChr<<"\t"<<chr2AlLength[lastChr]<<endl;
158 }
159 }
160
161 blockFirst=false;
162 }else{//scaf
163 //cout<<"sscaf "<<line<<endl;
164 vector<string> tokens = allTokens(line,' ');
165 lastScaff=tokens[1];
166 //cout<<line<<endl;
167 blockFirst=true;
168 }
169 }
170
171 }
172 myFile.close();
173 }else{
174 cerr << "Unable to open file "<<filename<<endl;
175 return 1;
176 }
177
178 if(!chr2count.empty()){
179
180 // for (map<string,int>::iterator it=chr2count.begin(); it!=chr2count.end(); ++it)
181 // cout << lastScaff<<"\t"<<it->first << "\t" << it->second << endl;
182
183 // vector<pair<string, int> > chr2countCpy(chr2count.begin(), chr2count.end());
184 // sort(chr2countCpy.begin(), chr2countCpy.end(), less_second<string, int>());
185 // for(unsigned int i=0;i<chr2countCpy.size();i++){
186 // cout<<chr2countCpy[i].first<<"\t"<<chr2countCpy[i].second<<endl;
187 // }
188 vector<pair<string, unsigned int> > chr2AlLengthCpy(chr2AlLength.begin(), chr2AlLength.end());
189 sort(chr2AlLengthCpy.begin(), chr2AlLengthCpy.end(), less_second<string, unsigned int>());
190 for(unsigned int i=0;i<chr2AlLengthCpy.size();i++){
191 cout<<lastScaff<<"\t"<<chr2AlLengthCpy[i].first<<"\t"<<chr2AlLengthCpy[i].second<<"\t"<<chr2count[chr2AlLengthCpy[i].first]<<endl;
192 }
193
194
195 }
196
197
198 return 0;
199}