last bunch of changes before jump to new version
[qmc.git] / analyze.cpp
blob6bf5f47aa35493a3af68d7c2237403681a541456
1 /*
2 * Copyright (c) 2009 Mauro Iazzi
4 * Permission is hereby granted, free of charge, to any person
5 * obtaining a copy of this software and associated documentation
6 * files (the "Software"), to deal in the Software without
7 * restriction, including without limitation the rights to use,
8 * copy, modify, merge, publish, distribute, sublicense, and/or sell
9 * copies of the Software, and to permit persons to whom the
10 * Software is furnished to do so, subject to the following
11 * conditions:
13 * The above copyright notice and this permission notice shall be
14 * included in all copies or substantial portions of the Software.
16 * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
17 * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES
18 * OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
19 * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT
20 * HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY,
21 * WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
22 * FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR
23 * OTHER DEALINGS IN THE SOFTWARE.
26 #include <string>
27 #include <vector>
28 #include <iostream>
29 #include <fstream>
30 #include <sstream>
31 #include <cmath>
32 #include <cstdlib>
34 using namespace std;
36 std::pair<double, double> compute_stats (const std::vector<double> &d) {
37 double avg = 0.0, var = 0.0;
38 const int N = d.size();
39 for (int i=0;i<N;i++) {
40 avg += d[i];
41 var += d[i]*d[i];
43 avg /= N;
44 var /= (N);
45 var -= avg*avg;
46 var = sqrt(var);
47 return std::pair<double, double>(avg, var);
50 std::string analyze_data_file (std::string filename, const int skip, const int len, int w = -1) {
51 std::vector<double> sums;
52 std::string ret = "";
54 int n = -skip;
55 double norm = 0.0;
56 std::string line_string = "";
57 std::string state = "";
58 std::ifstream f(filename.c_str(), std::fstream::in);
59 while (!f.eof()) {
60 double energy = 0.0;
61 double weight = 0.0;
62 getline(f, line_string);
63 std::stringstream line(line_string);
64 if (line_string[0]=='#') {
65 if (line_string[1]=='!') {
66 if (line_string.substr(2, 8)=="prepend ") {
67 ret.insert(0, " ");
68 ret.insert(0, line_string.substr(10, line_string.size()));
69 } else {
72 } else {
73 n++;
74 if (n>=0) {
75 line >> state >> energy;
76 weight = 1.0;
77 for (int j=0;j<w;j++) {
78 line >> weight;
80 norm += weight;
81 if (n%len==0) {
82 if (n>0) sums[sums.size()-1] /= norm;
83 sums.push_back(0.0);
84 norm = 0.0;
86 sums[n/len] += energy*(w?weight:1.0);
90 f.close();
91 sums.pop_back(); // last one may be incomplete.
92 // time to do some thing with the sums (now they are averages)
93 std::pair<double, double> s = compute_stats(sums);
94 std::stringstream rs;
95 rs << s.first << " " << s.second;
96 ret.append(rs.str());
97 return ret;
100 #define NEXTOPT(a, n) ( (a[i][2]=='\0')?(a[++i]):&(a[i][2]) )
102 int main (int argc, char **argv) {
103 int skip = 1000;
104 int len = 1000;
105 int w_col = -1;
106 vector<string> files;
107 for (int i=1;i<argc;i++) {
108 if (argv[i][0]=='-') {
109 switch (argv[i][1]) {
110 case 'w':
111 w_col = atoi(NEXTOPT(argv, i));
112 break;
113 case 's':
114 skip = atoi(NEXTOPT(argv, i));
115 break;
116 case 'l':
117 len = atoi(NEXTOPT(argv, i));
118 break;
120 } else {
121 files.push_back(argv[i]);
124 for (vector<string>::iterator iter = files.begin();iter != files.end();iter++) {
125 cout << analyze_data_file(*iter, skip, len, w_col) << endl;
127 return 0;