source: trunk/Mars/hawc/processing/currents.C@ 20033

Last change on this file since 20033 was 20030, checked in by tbretz, 5 years ago
Do not do any splitting of the aux data in the script files, instead add the commas in the root macros (simplification)
File size: 3.6 KB
Line 
1
2
3void currents(const char *fname, double beg=0, double end=100000)
4{
5 fits file(fname);
6
7 //file.PrintColumns();
8 //file.PrintKeys();
9
10 Double_t time;
11 Double_t jd;
12
13 float Imon[64];
14 float Tmon[64];
15 float Umon[64];
16 float Tpsu[2];
17 file.SetPtrAddress("Time", &time);
18 file.SetPtrAddress("Imon", Imon);
19 file.SetPtrAddress("Tmon", Tmon);
20 file.SetPtrAddress("Vmon", Umon);
21 file.SetPtrAddress("Tpsu", Tpsu);
22
23 UInt_t offset = file.GetUInt("MJDREF");
24 if (beg < 30000)
25 beg+=offset;
26 if (end < 30000)
27 end+=offset;
28
29 double AvgUmed = 0;
30 double AvgTmed = 0;
31 double AvgImed = 0;
32 double AvgUavg = 0;
33 double AvgTavg = 0;
34 double AvgIavg = 0;
35 double AvgUrms = 0;
36 double AvgTrms = 0;
37 double AvgIrms = 0;
38 double AvgUdev = 0;
39 double AvgTdev = 0;
40 double AvgIdev = 0;
41 double AvgPsu = 0;
42/*
43 double RmsUmed = 0;
44 double RmsTmed = 0;
45 double RmsImed = 0;
46 double RmsUavg = 0;
47 double RmsTavg = 0;
48 double RmsIavg = 0;
49 double RmsUrms = 0;
50 double RmsTrms = 0;
51 double RmsIrms = 0;
52 double RmsPsu = 0;
53*/
54 int count = 0;
55
56 double Umed = 0;
57 double Tmed = 0;
58 double Imed = 0;
59 double Udev = 0;
60 double Tdev = 0;
61 double Idev = 0;
62 double Uavg = 0;
63 double Tavg = 0;
64 double Iavg = 0;
65 double Urms = 0;
66 double Trms = 0;
67 double Irms = 0;
68
69 double psu = 0;
70
71 double diff = -1;
72
73 while (file.GetNextRow())
74 {
75 time += offset;
76 jd = time + 2400000.5;
77
78 if (time>end)
79 break;
80
81 Udev = MMath::MedianDev(64, Umon, Umed);
82 Tdev = MMath::MedianDev(64, Tmon, Tmed);
83 Idev = MMath::MedianDev(64, Imon, Imed);
84 Uavg = TMath::Mean(64, Umon);
85 Tavg = TMath::Mean(64, Tmon);
86 Iavg = TMath::Mean(64, Imon);
87 Urms = TMath::RMS(64, Umon);
88 Trms = TMath::RMS(64, Tmon);
89 Irms = TMath::RMS(64, Imon);
90
91 psu = (Tpsu[0]+Tpsu[1])/2;
92
93 if (time<beg)
94 {
95 diff = beg-time;
96 continue;
97 }
98
99 AvgUmed += Umed;
100 AvgTmed += Tmed;
101 AvgImed += Imed;
102 AvgUavg += Uavg;
103 AvgTavg += Tavg;
104 AvgIavg += Iavg;
105 AvgUrms += Urms;
106 AvgTrms += Trms;
107 AvgIrms += Irms;
108 AvgUdev += Udev;
109 AvgTdev += Tdev;
110 AvgIdev += Idev;
111 AvgPsu += psu;
112
113 count ++;
114 }
115
116 if (count==0)
117 {
118 if (diff<5./24/3600)
119 return;
120
121 cout << "result ";
122 cout << Umed << ", ";
123 cout << Uavg << ", ";
124 cout << Udev << ", ";
125 cout << Urms << ", ";
126
127 cout << Imed << ", ";
128 cout << Iavg << ", ";
129 cout << Idev << ", ";
130 cout << Irms << ", ";
131
132 cout << Tmed << ", ";
133 cout << Tavg << ", ";
134 cout << Tdev << ", ";
135 cout << Trms << ", ";
136
137 cout << psu;
138 cout << endl;
139 return;
140 }
141
142 AvgUmed /= count;
143 AvgImed /= count;
144 AvgTmed /= count;
145
146 AvgUavg /= count;
147 AvgIavg /= count;
148 AvgTavg /= count;
149
150 AvgUrms /= count;
151 AvgIrms /= count;
152 AvgTrms /= count;
153
154 AvgUdev /= count;
155 AvgIdev /= count;
156 AvgTdev /= count;
157
158 AvgPsu /= count;
159
160
161 cout << "result";
162 cout << AvgUmed << ", ";
163 cout << AvgUavg << ", ";
164 cout << AvgUdev << ", ";
165 cout << AvgUrms << ", ";
166
167 cout << AvgImed << ", ";
168 cout << AvgIavg << ", ";
169 cout << AvgIdev << ", ";
170 cout << AvgIrms << ", ";
171
172 cout << AvgTmed << ", ";
173 cout << AvgTavg << ", ";
174 cout << AvgTdev << ", ";
175 cout << AvgTrms << ", ";
176
177 cout << AvgPsu;
178
179 cout << endl;
180}
Note: See TracBrowser for help on using the repository browser.