source: trunk/Mars/fact/processing/currents.C@ 15658

Last change on this file since 15658 was 15222, checked in by Daniela Dorner, 12 years ago
added (macro to extract currents from the file CALIBRATED_CURRENTS)
File size: 1.9 KB
Line 
1void currents(const char *fname, double beg=0, double end=100000)
2{
3 fits file(fname);
4
5 //file.PrintColumns();
6 //file.PrintKeys();
7
8 Double_t time;
9 Float_t med, dev;
10 file.SetPtrAddress("Time", &time);
11 file.SetPtrAddress("I_med", &med);
12 file.SetPtrAddress("I_dev", &dev);
13
14 UInt_t offset = file.GetUInt("MJDREF");
15 if (beg < 30000)
16 beg+=offset;
17 if (end < 30000)
18 end+=offset;
19
20 double med_beg = 0;
21 double med_end = 0;
22 double med_avg = 0;
23 double med_rms = 0;
24 double dev_avg = 0;
25 double dev_rms = 0;
26 int cnt = 0;
27 int cnt_beg = 0;
28 int cnt_end = 0;
29
30 double med_last = -1;
31 double dev_last = -1;
32 double diff = -1;
33
34 while (file.GetNextRow())
35 {
36 time += offset;
37
38 if (time>end)
39 break;
40
41 if (time<beg)
42 {
43 med_last = med;
44 dev_last = dev;
45 diff = beg-time;
46 continue;
47 }
48
49 med_avg += med;
50 med_rms += med*med;
51 dev_avg += dev;
52 dev_rms += dev*dev;
53
54 if (time<beg+15./24/3600)
55 {
56 med_beg += med;
57 cnt_beg++;
58 }
59 if (time>end-15./24/3600)
60 {
61 med_end += med;
62 cnt_end++;
63 }
64
65 cnt++;
66 }
67
68 if (cnt==0)
69 {
70 if (diff<5./24/3600)
71 return;
72
73 cout << "result " << med_last << " 0 " << dev_last << " 0 0 0" << endl;
74 return;
75 }
76
77 med_beg = cnt_beg>0 ? med_beg/cnt_beg : -1;
78 med_end = cnt_end>0 ? med_end/cnt_end : -1;
79
80 med_avg /= cnt;
81 med_rms /= cnt;
82 med_rms = sqrt(med_rms-med_avg*med_avg);
83
84 dev_avg /= cnt;
85 dev_rms /= cnt;
86 dev_rms = sqrt(dev_rms-dev_avg*dev_avg);
87
88 cout << "result " << med_avg << " " << med_rms << " " << dev_avg << " " << dev_rms << " " << med_beg << " " << med_end << endl;
89}
Note: See TracBrowser for help on using the repository browser.