source: trunk/Mars/hawc/quality.C@ 20014

Last change on this file since 20014 was 20003, checked in by tbretz, 5 years ago
New functions, more plots.
File size: 16.2 KB
Line 
1#include <iostream>
2#include "TROOT.h"
3#include "TSystem.h"
4#include "TGraph.h"
5#include "TArrayD.h"
6#include "TCanvas.h"
7#include "TH1.h"
8#include "TAxis.h"
9#include "TText.h"
10#include "TLine.h"
11#include "TStyle.h"
12
13#include "fits.h"
14#include "MTime.h"
15
16#ifndef __CLING__
17#include <algorithm>
18#include <functional>
19#include <vector>
20#endif
21
22#include <numeric>
23
24Bool_t Contains(vector<double> *vec, Double_t t0, Double_t range=0)
25{
26 auto &arr0 = vec[0];
27 auto &arr1 = vec[1];
28 auto &arr2 = vec[2];
29
30 for (size_t i=0; i<arr0.size(); i++)
31 {
32 Double_t start = arr1[i];
33 Double_t stop = arr2[i];
34
35 if (stop>start+305./24/3600)
36 stop = start+305./24/3600;
37
38 if (t0>start-range && t0<stop+range)
39 {
40 //if (fmod(t0, 1)>4./24 && fmod(t0,1)<4.1/24)
41 // cout << t0-start << " " <<t0 << " " << stop-t0 << " " << start-15779 << " " << stop-15779 << " " << (*arr0)[i] << endl;
42 return kTRUE;
43 }
44 }
45
46 return arr0.empty();
47}
48
49Int_t ReadRuns(vector<double> *vec, TString fname)
50{
51 fname += ".FAD_CONTROL_RUNS.fits";
52
53 fits file(fname.Data());
54 if (!file)
55 {
56 cerr << fname << ": " << gSystem->GetError() << endl;
57 return -2;
58 }
59
60 //cout << fname << endl;
61
62 Double_t time;
63 uint32_t qos, stats[2];
64
65 if (!file.SetPtrAddress("Time", &time))
66 return -1;
67 if (!file.SetPtrAddress("QoS", &qos))
68 return -1;
69 if (!file.SetPtrAddress("stats", stats))
70 return -1;
71
72 int32_t run = -1;
73 double start=0;
74
75 while (file.GetNextRow())
76 {
77 switch (qos)
78 {
79 case 1:
80 start = time;
81 run = stats[0];
82 break;
83 case 0:
84 if (uint32_t(run)==stats[1] && stats[0]==stats[1])
85 {
86 vec[0].push_back(run);
87 vec[1].push_back(start);
88 vec[2].push_back(time);
89 //cout << setw(3) << ": " << MTime(start) << " -- " << MTime(time) << endl;
90 }
91 run = -1;
92 break;
93 }
94 }
95
96 return 0;
97}
98
99
100Int_t PlotThresholds(double start, vector<double> *vec, TString fname)
101{
102 fname += ".FTU_CONTROL_DAC.fits";
103
104 fits file(fname.Data());
105 if (!file)
106 {
107 cerr << fname << ": " << gSystem->GetError() << endl;
108 return -2;
109 }
110
111 //cout << fname << endl;
112
113 Double_t time;
114 UShort_t th[4];
115
116 if (!file.SetPtrAddress("Time", &time))
117 return -1;
118
119 if (!file.SetPtrAddress("dac_ch", th))
120 return -1;
121
122 TGraph g;
123 g.SetName("Threshold");
124 while (file.GetNextRow())
125 {
126 if (Contains(vec, time, 10./(24*3600)))
127 g.SetPoint(g.GetN(), time*24*3600, th[0]);
128 }
129
130
131 gROOT->SetSelectedPad(0);
132 gPad->GetCanvas()->cd(1);
133 gPad->SetGrid();
134 gPad->SetTopMargin(0);
135 gPad->SetRightMargin(0.001);
136 gPad->SetLeftMargin(0.04);
137 gPad->SetBottomMargin(0);
138
139 TH1C frame("Threshold_frame", "", 1000, (start+0.45)*24*3600, (start+1.025)*24*3600);
140 frame.SetMinimum(300);
141 frame.SetMaximum(1200);
142 frame.SetStats(kFALSE);
143 frame.GetXaxis()->SetTimeDisplay(true);
144 frame.GetXaxis()->SetTimeFormat("%H:%M %F1995-01-01 00:00:00 GMT");
145 frame.GetXaxis()->SetLabelSize(0.12);
146 frame.GetYaxis()->SetLabelSize(0.1);
147 frame.GetYaxis()->SetTitle("THRESHOLD");
148 frame.GetYaxis()->CenterTitle();
149 frame.GetYaxis()->SetTitleOffset(0.2);
150 frame.GetYaxis()->SetTitleSize(0.1);
151 frame.DrawCopy();
152
153 g.SetMarkerStyle(kFullDotMedium);
154 g.DrawClone("P");
155
156 return 0;
157}
158
159Int_t PlotCurrent(double start, vector<double> *vec, TString fname)
160{
161 fname += ".BIAS_CONTROL_AVERAGE_DATA.fits";
162
163 fits file(fname.Data());
164 if (!file)
165 {
166 cerr << fname << ": " << gSystem->GetError() << endl;
167 return -2;
168 }
169
170 //cout << fname << endl;
171
172 Double_t time;
173 Float_t Iavg[64];
174 Float_t Tavg[64];
175 Float_t Tpsu[2];
176 Float_t Vavg[64];
177 Float_t Uavg[64];
178
179 if (!file.SetPtrAddress("Time", &time))
180 return -1;
181
182 if (!file.SetPtrAddress("Iavg", Iavg))
183 return -1;
184
185 if (!file.SetPtrAddress("Tavg", Tavg))
186 return -1;
187
188 if (!file.SetPtrAddress("Tavg_psu", Tpsu))
189 return -1;
190
191 if (!file.SetPtrAddress("Vavg", Vavg))
192 return -1;
193
194 bool has_uset = true;
195 try
196 {
197 has_uset = file.SetPtrAddress("Uavg", Uavg);
198 }
199 catch (const exception&)
200 {
201 }
202
203 TGraph g1;
204 TGraph g2;
205 TGraph g3;
206 TGraph g4;
207 TGraph g5;
208 g1.SetName("CurrentAvg");
209 g2.SetName("TemperatureAvg");
210 g3.SetName("TempPSUavg");
211 g4.SetName("VoltageAvg");
212 g5.SetName("UsetAvg");
213
214 while (file.GetNextRow())
215 if (Contains(vec, time))
216 {
217 sort(Iavg, Iavg+64);
218 sort(Tavg, Tavg+64);
219 sort(Vavg, Vavg+64);
220 if (has_uset)
221 sort(Uavg, Uavg+64);
222
223 g1.SetPoint(g1.GetN(), time*24*3600, Iavg[31]);
224 g2.SetPoint(g2.GetN(), time*24*3600, Tavg[31]);
225 g3.SetPoint(g3.GetN(), time*24*3600, (Tpsu[0]+Tpsu[1])/2);
226 g4.SetPoint(g4.GetN(), time*24*3600, Vavg[31]);
227 if (has_uset)
228 g5.SetPoint(g5.GetN(), time*24*3600, Uavg[31]);
229 }
230
231
232 gROOT->SetSelectedPad(0);
233 gPad->GetCanvas()->cd(2);
234 gPad->SetGrid();
235 gPad->SetTopMargin(0);
236 gPad->SetRightMargin(0.001);
237 gPad->SetLeftMargin(0.04);
238 gPad->SetBottomMargin(0);
239
240 TH1C frame1("Current_frame", "", 1000, (start+0.45)*24*3600, (start+1.025)*24*3600);
241 frame1.SetMinimum(0);
242 frame1.SetMaximum(1000);
243 frame1.SetStats(kFALSE);
244 frame1.GetXaxis()->SetTimeDisplay(true);
245 frame1.GetXaxis()->SetTimeFormat("%H:%M %F1995-01-01 00:00:00 GMT");
246 frame1.GetXaxis()->SetLabelSize(0.12);
247 frame1.GetYaxis()->SetLabelSize(0.1);
248 frame1.GetYaxis()->SetTitle("CURRENT [#mu A]");
249 frame1.GetYaxis()->CenterTitle();
250 frame1.GetYaxis()->SetTitleOffset(0.2);
251 frame1.GetYaxis()->SetTitleSize(0.1);
252 frame1.DrawCopy();
253
254 g1.SetMarkerStyle(kFullDotMedium);
255 g1.DrawClone("P");
256
257
258
259 gROOT->SetSelectedPad(0);
260 gPad->GetCanvas()->cd(5);
261 gPad->SetGrid();
262 gPad->SetTopMargin(0);
263 gPad->SetRightMargin(0.001);
264 gPad->SetLeftMargin(0.04);
265 gPad->SetBottomMargin(0);
266
267 TH1C frame4("Voltage_frame", "", 1000, (start+0.45)*24*3600, (start+1.025)*24*3600);
268 frame4.SetMinimum(29);
269 frame4.SetMaximum(30);
270 frame4.SetStats(kFALSE);
271 frame4.GetXaxis()->SetTimeDisplay(true);
272 frame4.GetXaxis()->SetTimeFormat("%H:%M %F1995-01-01 00:00:00 GMT");
273 frame4.GetXaxis()->SetLabelSize(0.12);
274 frame4.GetYaxis()->SetLabelSize(0.1);
275 frame4.GetYaxis()->SetTitle("BIAS VOLTAGE [V]");
276 frame4.GetYaxis()->CenterTitle();
277 frame4.GetYaxis()->SetTitleOffset(0.2);
278 frame4.GetYaxis()->SetTitleSize(0.1);
279 frame4.DrawCopy();
280
281 if (has_uset)
282 {
283 g5.SetMarkerColor(kBlue);
284 g5.SetMarkerStyle(kFullDotMedium);
285 g5.DrawClone("P");
286 }
287
288 g4.SetMarkerStyle(kFullDotMedium);
289 g4.DrawClone("P");
290
291
292 gROOT->SetSelectedPad(0);
293 gPad->GetCanvas()->cd(6);
294 gPad->SetGrid();
295 gPad->SetTopMargin(0);
296 gPad->SetRightMargin(0.001);
297 gPad->SetLeftMargin(0.04);
298 gPad->SetBottomMargin(0);
299
300 TH1C frame2("Temperature_frame", "", 1000, (start+0.45)*24*3600, (start+1.025)*24*3600);
301 frame2.SetMinimum(0);
302 frame2.SetMaximum(65);
303 frame2.SetStats(kFALSE);
304 frame2.GetXaxis()->SetTimeDisplay(true);
305 frame2.GetXaxis()->SetTimeFormat("%H:%M %F1995-01-01 00:00:00 GMT");
306 frame2.GetXaxis()->SetLabelSize(0.12);
307 frame2.GetYaxis()->SetLabelSize(0.1);
308 frame2.GetYaxis()->SetTitle("TEMPERATURE [#circ C]");
309 frame2.GetYaxis()->CenterTitle();
310 frame2.GetYaxis()->SetTitleOffset(0.2);
311 frame2.GetYaxis()->SetTitleSize(0.1);
312 frame2.DrawCopy();
313
314 g2.SetMarkerStyle(kFullDotMedium);
315 g2.DrawClone("P");
316
317 g3.SetMarkerColor(kRed);
318 g3.SetMarkerStyle(kFullDotMedium);
319 g3.DrawClone("P");
320
321 return 0;
322}
323
324Int_t PlotDrsTemp(double /*start*/, vector<double> *vec, TString fname)
325{
326 fname += ".FAD_CONTROL_TEMPERATURE.fits";
327
328 fits file(fname.Data());
329 if (!file)
330 {
331 cerr << fname << ": " << gSystem->GetError() << endl;
332 return -2;
333 }
334
335 //cout << fname << endl;
336
337 Double_t time;
338 uint16_t cnt;
339 Float_t temp[160];
340
341 if (!file.SetPtrAddress("Time", &time))
342 return -1;
343
344 if (!file.SetPtrAddress("cnt", &cnt))
345 return -1;
346
347 if (!file.SetPtrAddress("temp", temp))
348 return -1;
349
350 TGraph g1, g2;
351 g1.SetName("TempDRS1");
352 g2.SetName("TempDRS2");
353
354 while (file.GetNextRow())
355 if (Contains(vec, time))
356 {
357 g1.SetPoint(g1.GetN(), time*24*3600, accumulate(temp, temp+4, 0.)/4);
358 g2.SetPoint(g2.GetN(), time*24*3600, accumulate(temp+4, temp+8, 0.)/4);
359 }
360
361
362 gROOT->SetSelectedPad(0);
363 gPad->GetCanvas()->cd(6);
364
365 g1.SetMarkerColor(kBlue);
366 g1.SetMarkerStyle(kFullDotMedium);
367 g1.DrawClone("P");
368
369 g2.SetMarkerColor(kBlue);
370 g2.SetMarkerStyle(kFullDotMedium);
371 g2.DrawClone("P");
372
373 return 0;
374}
375
376Int_t PlotPatchRate(double /*start*/, vector<double> *vec, TString fname)
377{
378 fname += ".FTU_CONTROL_DATA.fits";
379
380 fits file(fname.Data());
381 if (!file)
382 {
383 cerr << fname << ": " << gSystem->GetError() << endl;
384 return -2;
385 }
386
387 //cout << fname << endl;
388
389 Double_t time;
390 uint32_t qos;
391 uint32_t counter_ch[8];
392 float dt_sec;
393
394 if (!file.SetPtrAddress("QoS", &qos))
395 return -1;
396 if (!file.SetPtrAddress("Time", &time))
397 return -1;
398
399 if (!file.SetPtrAddress("counter_ch", counter_ch))
400 return -1;
401 if (!file.SetPtrAddress("dt_sec", &dt_sec))
402 return -1;
403
404 TGraph g1, g2;
405 g1.SetName("MinPatchRate");
406 g2.SetName("MaxPatchRate");
407
408 while (file.GetNextRow())
409 if (Contains(vec, time))
410 {
411 if (qos==0) // Only automatic reports
412 continue;
413
414 sort(counter_ch, counter_ch+8);
415
416 g1.SetPoint(g1.GetN(), time*24*3600, counter_ch[1]/dt_sec);
417 g2.SetPoint(g2.GetN(), time*24*3600, counter_ch[7]/dt_sec);
418 }
419
420
421 gROOT->SetSelectedPad(0);
422 gPad->GetCanvas()->cd(3);
423 /*
424 gPad->SetGrid();
425 gPad->SetTopMargin(0);
426 gPad->SetBottomMargin(0);
427 gPad->SetRightMargin(0.001);
428 gPad->SetLeftMargin(0.04);
429
430 TH1C frame1("Rate_frame", "", 1000, (start+0.45)*24*3600, (start+1.025)*24*3600);
431 frame1.SetMinimum(0);
432 frame1.SetMaximum(90);
433 frame1.SetStats(kFALSE);
434 frame1.GetXaxis()->SetTimeDisplay(true);
435 frame1.GetXaxis()->SetTimeFormat("%H:%M %F1995-01-01 00:00:00 GMT");
436 frame1.GetXaxis()->SetLabelSize(0.12);
437 frame1.GetYaxis()->SetLabelSize(0.1);
438 frame1.GetYaxis()->SetTitle("TRIGGER RATE");
439 frame1.GetYaxis()->CenterTitle();
440 frame1.GetYaxis()->SetTitleOffset(0.2);
441 frame1.GetYaxis()->SetTitleSize(0.1);
442 frame1.DrawCopy();*/
443
444 g1.SetMarkerColor(kGreen);
445 g1.SetMarkerStyle(kFullDotMedium);
446 g1.DrawClone("P");
447
448 g2.SetMarkerColor(kRed);
449 g2.SetMarkerStyle(kFullDotMedium);
450 g2.DrawClone("P");
451
452 return 0;
453}
454
455Int_t PlotCameraRate(double start, vector<double> *vec, TString fname)
456{
457 fname += ".FTM_CONTROL_DATA.fits";
458
459 fits file(fname.Data());
460 if (!file)
461 {
462 cerr << fname << ": " << gSystem->GetError() << endl;
463 return -2;
464 }
465
466 //cout << fname << endl;
467
468 Double_t time;
469 uint32_t trg_counter;
470 uint32_t dead_time, run_time;
471 float temp;
472
473 if (!file.SetPtrAddress("Time", &time))
474 return -1;
475
476 if (!file.SetPtrAddress("trg_counter", &trg_counter))
477 return -1;
478 if (!file.SetPtrAddress("dead_time", &dead_time))
479 return -1;
480 if (!file.SetPtrAddress("run_time", &run_time))
481 return -1;
482 if (!file.SetPtrAddress("temp", &temp))
483 return -1;
484
485 TGraph g1, g2, g3;
486 g1.SetName("TriggerRate");
487 g2.SetName("RelOnTime");
488 g3.SetName("TempFTM");
489
490
491 uint32_t prev = 0;
492
493 while (file.GetNextRow())
494 if (Contains(vec, time))
495 {
496 g1.SetPoint(g1.GetN(), time*24*3600, 1e8*(trg_counter-prev)/(run_time-dead_time));
497 g2.SetPoint(g2.GetN(), time*24*3600, 1.-float(dead_time)/run_time);
498 g3.SetPoint(g2.GetN(), time*24*3600, temp);
499 prev = trg_counter;
500 }
501
502
503 gROOT->SetSelectedPad(0);
504 gPad->GetCanvas()->cd(3);
505 gPad->SetGrid();
506 gPad->SetTopMargin(0);
507 gPad->SetBottomMargin(0);
508 gPad->SetRightMargin(0.001);
509 gPad->SetLeftMargin(0.04);
510
511 TH1C frame1("Rate_frame", "", 1000, (start+0.45)*24*3600, (start+1.025)*24*3600);
512 frame1.SetMinimum(0);
513 frame1.SetMaximum(90);
514 frame1.SetStats(kFALSE);
515 frame1.GetXaxis()->SetTimeDisplay(true);
516 frame1.GetXaxis()->SetTimeFormat("%H:%M %F1995-01-01 00:00:00 GMT");
517 frame1.GetXaxis()->SetLabelSize(0.12);
518 frame1.GetYaxis()->SetLabelSize(0.1);
519 frame1.GetYaxis()->SetTitle("TRIGGER RATE [Hz]");
520 frame1.GetYaxis()->CenterTitle();
521 frame1.GetYaxis()->SetTitleOffset(0.2);
522 frame1.GetYaxis()->SetTitleSize(0.1);
523 frame1.DrawCopy();
524
525 g1.SetMarkerStyle(kFullDotMedium);
526 g1.DrawClone("P");
527
528
529
530 gROOT->SetSelectedPad(0);
531 gPad->GetCanvas()->cd(4);
532 gPad->SetGrid();
533 gPad->SetTopMargin(0);
534 gPad->SetBottomMargin(0);
535 gPad->SetRightMargin(0.001);
536 gPad->SetLeftMargin(0.04);
537
538 TH1C frame2("OnTime_frame", "", 1000, (start+0.45)*24*3600, (start+1.025)*24*3600);
539 frame2.SetMinimum(0);
540 frame2.SetMaximum(1);
541 frame2.SetStats(kFALSE);
542 frame2.GetXaxis()->SetTimeDisplay(true);
543 frame2.GetXaxis()->SetTimeFormat("%H:%M %F1995-01-01 00:00:00 GMT");
544 frame2.GetXaxis()->SetLabelSize(0.12);
545 frame2.GetYaxis()->SetLabelSize(0.1);
546 frame2.GetYaxis()->SetTitle("RELATIVE ON TIME");
547 frame2.GetYaxis()->CenterTitle();
548 frame2.GetYaxis()->SetTitleOffset(0.2);
549 frame2.GetYaxis()->SetTitleSize(0.1);
550 frame2.DrawCopy();
551
552 g2.SetMarkerStyle(kFullDotMedium);
553 g2.DrawClone("P");
554
555 gROOT->SetSelectedPad(0);
556 gPad->GetCanvas()->cd(6);
557
558 g3.SetMarkerColor(kGreen);
559 g3.SetMarkerStyle(kFullDotMedium);
560 g3.DrawClone("P");
561
562 return 0;
563}
564
565int quality(const char *basepath="/data/aux", UInt_t y=0, UInt_t m=0, UInt_t d=0, const char *outpath=0, const bool all=true)
566{
567 // To get correct dates in the histogram you have to add
568 // the MJDREF offset (should be 40587) and 9131.
569
570 if (y==0)
571 {
572 UInt_t nt = MTime(MTime(-1).GetMjd()-1.5).GetNightAsInt();
573 y = nt/10000;
574 m = (nt/100)%100;
575 d = nt%100;
576
577 cout << y << "/" << m << "/" << d << endl;
578 }
579
580 TString fname=Form("%s/%04d/%02d/%02d/%04d%02d%02d", basepath, y, m, d, y, m, d);
581
582 cout << "quality" << endl;
583 cout << "-------" << endl;
584 cout << endl;
585 cout << "Night: " << Form("%04d-%02d-%02d", y, m, d) << endl;
586 cout << endl;
587
588 const double start = MTime(y, m, d).GetMjd()-40587;
589
590 vector<double> runs[3]; // { &run, &beg, &end };
591 if (!all)
592 ReadRuns(runs, fname);
593
594 //check if the sqm was already installed on the telescope
595 TCanvas *c = new TCanvas("quality", Form("Quality %04d/%02d/%02d", y, m, d), 1280, 1120);
596 c->Divide(1, 7, 1e-5, 1e-5);
597
598 PlotThresholds(start, runs, fname);
599 PlotCurrent(start, runs, fname);
600 PlotDrsTemp(start, runs, fname);
601 PlotCameraRate(start, runs, fname);
602 PlotPatchRate(start, runs, fname);
603
604 gROOT->SetSelectedPad(0);
605 gPad->GetCanvas()->cd(7);
606 gPad->SetTopMargin(0);
607 gPad->SetBottomMargin(0.99);
608 gPad->SetRightMargin(0.001);
609 gPad->SetLeftMargin(0.04);
610
611 TH1C frame1("Axis", "", 1000, (start+0.45)*24*3600, (start+1.025)*24*3600);
612 frame1.SetMinimum(0);
613 frame1.SetMaximum(1);
614 frame1.SetStats(kFALSE);
615 frame1.GetXaxis()->SetTimeDisplay(true);
616 frame1.GetXaxis()->SetTimeFormat("%H:%M %F1995-01-01 00:00:00");// 18:00:00 GMT);
617 frame1.GetXaxis()->SetLabelSize(0.2);
618 frame1.GetXaxis()->SetTitleSize(0.2);
619 frame1.GetXaxis()->SetTitle("HAWC Local Time");
620 frame1.GetXaxis()->CenterTitle();
621 frame1.DrawCopy();
622
623 if (outpath)
624 c->SaveAs(Form("%s/quality-%04d%02d%02d.pdf", outpath, y, m, d));
625
626 return 0;
627}
628
629int quality(UInt_t y, UInt_t m, UInt_t d, const char *outpath=0, const bool all=true)
630{
631 return quality("/data/aux", y, m, d, outpath, all);
632}
Note: See TracBrowser for help on using the repository browser.