source: branches/Corsika7500Compatibility/fact/plots/quality.C@ 19366

Last change on this file since 19366 was 18260, checked in by mknoetig, 10 years ago
quality.C had a small bug that droppend the x-label fro the axis when no sqm data was present ( before 20140823)
File size: 25.5 KB
Line 
1#include <algorithm>
2#include <functional>
3
4Bool_t Contains(TArrayD **vec, Double_t t0, Double_t range=0)
5{
6 TArrayD *arr0 = vec[0];
7 TArrayD *arr1 = vec[1];
8 TArrayD *arr2 = vec[2];
9
10 for (int i=0; i<arr0->GetSize(); i++)
11 {
12 Double_t start = (*arr1)[i];
13 Double_t stop = (*arr2)[i];
14
15 if (stop>start+305./24/3600)
16 stop = start+305./24/3600;
17
18 if (t0>start-range && t0<stop+range)
19 //{
20 // if (fmod(t0, 1)>4./24 && fmod(t0,1)<4.1/24)
21 // cout << t0-start << " " <<t0 << " " << stop-t0 << " " << start-15779 << " " << stop-15779 << " " << (*arr0)[i] << endl;
22 return kTRUE;
23 //}
24 }
25
26 return arr0->GetSize()==0;
27}
28
29Int_t PlotThresholds(TArrayD **vec, TString fname)
30{
31 fname += ".RATE_CONTROL_THRESHOLD.fits";
32
33 fits file(fname.Data());
34 if (!file)
35 {
36 cerr << fname << ": " << gSystem->GetError() << endl;
37 return -2;
38 }
39
40 //cout << fname << endl;
41
42 Double_t time;
43 UShort_t th;
44
45 if (!file.SetPtrAddress("Time", &time))
46 return -1;
47
48 if (!file.SetPtrAddress("threshold", &th))
49 return -1;
50
51 TGraph g;
52 g.SetName("Threshold");
53
54 while (file.GetNextRow())
55 {
56 if (Contains(vec, time, 10./(24*3600)))
57 g.SetPoint(g.GetN(), time*24*3600, th);
58 }
59
60 g.SetMinimum(281);
61 g.SetMarkerStyle(kFullDotMedium);
62 g.GetXaxis()->SetTimeDisplay(true);
63 g.GetXaxis()->SetTimeFormat("%H:%M %F1995-01-01 00:00:00 GMT");
64 g.GetXaxis()->SetLabelSize(0.12);
65 g.GetYaxis()->SetLabelSize(0.1);
66 g.GetYaxis()->SetTitle("THRESHOLD");
67 g.GetYaxis()->SetTitleOffset(0.2);
68 g.GetYaxis()->SetTitleSize(0.1);
69 g.DrawClone("AP");
70
71 return 0;
72}
73
74#include <vector>
75#include <pair>
76
77vector<pair<double, Nova::EquPosn>> vecp;
78
79Nova::EquPosn FindPointing(Double_t time)
80{
81 for (int i=0; i<vecp.size(); i++)
82 if (time<vecp[i].first)
83 {
84 if (i==0)
85 return Nova::EquPosn();
86 else
87 return vecp[i-1].second;
88 }
89
90 return vecp[vecp.size()-1].second;
91}
92
93Float_t Prediction(Double_t time)
94{
95 Double_t jd = time + 40587 + 2400000.5;
96
97 // Sun properties
98 Nova::EquPosn sun = Nova::GetSolarEquCoords(jd);
99 Nova::ZdAzPosn hrzs = Nova::GetHrzFromEqu(sun, jd);
100
101 // Get source position
102 Nova::EquPosn pos = FindPointing(time);
103
104
105 // Moon properties
106 Nova::EquPosn moon = Nova::GetLunarEquCoords(jd, 0.01);
107 Nova::HrzPosn hrzm = Nova::GetHrzFromEqu(moon, jd);
108 double disk = Nova::GetLunarDisk(jd);
109
110 // Derived moon properties
111 double angle = Nova::GetAngularSeparation(moon, pos);
112 double edist = Nova::GetLunarEarthDist(jd)/384400;
113
114 // Current prediction
115 double sin_malt = hrzm.alt<0 ? 0 : sin(hrzm.alt*TMath::DegToRad());
116 double cos_mdist = cos(angle*TMath::DegToRad());
117 double sin_szd = sin(hrzs.zd*TMath::DegToRad());
118
119 double c0 = pow(disk, 2.63);
120 double c1 = pow(sin_malt, 0.60);
121 double c2 = pow(edist, -2.00);
122 double c3 = exp(0.67*cos_mdist*cos_mdist*cos_mdist*cos_mdist);
123 double c4 = exp(-97.8+105.8*sin_szd*sin_szd);
124
125 double cur = 6.2 + 95.7*c0*c1*c2*c3 + c4;
126
127 return cur;
128}
129
130Int_t ReadSources(TString fname)
131{
132 fname += ".DRIVE_CONTROL_SOURCE_POSITION.fits";
133
134 fits file(fname.Data());
135 if (!file)
136 {
137 cerr << fname << ": " << gSystem->GetError() << endl;
138 return -2;
139 }
140
141 Double_t time, ra, dec;
142 if (!file.SetPtrAddress("Time", &time))
143 return -1;
144 if (!file.SetPtrAddress("Ra_cmd", &ra))
145 return -1;
146 if (!file.SetPtrAddress("Dec_cmd", &dec))
147 return -1;
148
149 while (file.GetNextRow())
150 {
151 Nova::EquPosn p;
152 p.ra = ra*15;
153 p.dec = dec;
154
155 vecp.push_back(make_pair(time, p));
156 }
157
158 return 0;
159}
160
161Int_t PlotCurrent(TArrayD **vec, TString fname)
162{
163 Int_t rc = ReadSources(fname);
164 if (rc<0)
165 return rc;
166
167 fname += ".FEEDBACK_CALIBRATED_CURRENTS.fits";
168
169 fits file(fname.Data());
170 if (!file)
171 {
172 cerr << fname << ": " << gSystem->GetError() << endl;
173 return -2;
174 }
175
176 //cout << fname << endl;
177
178 Double_t time;
179 Float_t Imed;
180 Float_t Idev;
181 Float_t I[416];
182
183 if (!file.SetPtrAddress("Time", &time))
184 return -1;
185
186 if (!file.SetPtrAddress("I_med", &Imed))
187 return -1;
188
189 if (!file.SetPtrAddress("I_dev", &Idev))
190 return -1;
191
192 if (!file.SetPtrAddress("I", I))
193 return -1;
194
195 TGraph g1;
196 TGraph g2;
197 TGraph g3;
198 TGraph g4;
199 TGraph g5;
200 g1.SetName("CurrentsMed");
201 g2.SetName("Prediction");
202 g3.SetName("CurrentsDev");
203 g4.SetName("CurrentsMax-4");
204 g5.SetName("CurrentsMax");
205
206 while (file.GetNextRow())
207 if (Contains(vec, time))
208 {
209 // crazy pixels
210 I[66] = 0;
211 I[191] = 0;
212 I[193] = 0;
213
214 sort(I, I+320);
215
216 g1.SetPoint(g1.GetN(), time*24*3600, Imed);
217 g2.SetPoint(g2.GetN(), time*24*3600, Prediction(time));
218 g3.SetPoint(g3.GetN(), time*24*3600, Imed+Idev);
219 g4.SetPoint(g4.GetN(), time*24*3600, I[315]);
220 g5.SetPoint(g5.GetN(), time*24*3600, I[319]);
221 }
222
223 g1.SetMinimum(0);
224 g1.SetMaximum(149);
225
226 g1.SetMarkerStyle(kFullDotMedium);
227 g1.GetXaxis()->SetTimeDisplay(true);
228 g1.GetXaxis()->SetTimeFormat("%H:%M %F1995-01-01 00:00:00 GMT");
229 g1.GetXaxis()->SetLabelSize(0.12);
230 g1.GetYaxis()->SetLabelSize(0.1);
231 g1.GetYaxis()->SetTitle("CURRENT");
232 g1.GetYaxis()->SetTitleOffset(0.2);
233 g1.GetYaxis()->SetTitleSize(0.1);
234 g1.DrawClone("AP");
235
236 g5.SetMarkerColor(kGray);
237 g5.DrawClone("P");
238
239 g4.SetMarkerColor(kGray+1);
240 g4.DrawClone("P");
241
242 g3.SetMarkerColor(kGray+2);
243 g3.DrawClone("P");
244
245 g2.SetMarkerColor(kBlue);
246 g2.SetMarkerStyle(kFullDotMedium);
247 g2.DrawClone("P");
248
249 g1.DrawClone("P");
250
251 return 0;
252}
253
254Int_t PlotRate(TArrayD **vec, TString fname)
255{
256 fname += ".FTM_CONTROL_TRIGGER_RATES.fits";
257
258 fits file(fname.Data());
259 if (!file)
260 {
261 cerr << fname << ": " << gSystem->GetError() << endl;
262 return -2;
263 }
264
265 //cout << fname << endl;
266
267 Double_t time;
268 Float_t rate;
269 Float_t ontime, elapsed;
270
271 if (!file.SetPtrAddress("Time", &time))
272 return -1;
273
274 if (!file.SetPtrAddress("TriggerRate", &rate))
275 return -1;
276 if (!file.SetPtrAddress("OnTime", &ontime))
277 return -1;
278 if (!file.SetPtrAddress("ElapsedTime", &elapsed))
279 return -1;
280
281 TGraph g1, g2;
282 g1.SetName("TriggerRate");
283 g2.SetName("RelOnTime");
284
285 while (file.GetNextRow())
286 if (Contains(vec, time))
287 {
288 if (rate>=0)
289 {
290 g1.SetPoint(g1.GetN(), time*24*3600, rate);
291 g2.SetPoint(g2.GetN(), time*24*3600, ontime/elapsed);
292 }
293 }
294
295 g1.SetMinimum(0);
296 g1.SetMaximum(269);
297 g1.SetMarkerStyle(kFullDotMedium);
298 g1.GetXaxis()->SetTimeDisplay(true);
299 g1.GetXaxis()->SetTimeFormat("%H:%M %F1995-01-01 00:00:00 GMT");
300 g1.GetXaxis()->SetLabelSize(0.12);
301 g1.GetYaxis()->SetLabelSize(0.1);
302 g1.GetYaxis()->SetTitle("TRIGGER RATE");
303 g1.GetYaxis()->SetTitleOffset(0.2);
304 g1.GetYaxis()->SetTitleSize(0.1);
305 g1.DrawClone("AP");
306
307 gROOT->SetSelectedPad(0);
308 gPad->GetCanvas()->cd(4);
309
310 gPad->SetGrid();
311 gPad->SetTopMargin(0);
312 gPad->SetBottomMargin(0);
313 gPad->SetRightMargin(0.001);
314 gPad->SetLeftMargin(0.04);
315
316 g2.SetMinimum(0);
317 g2.SetMaximum(1);
318 g2.SetMarkerStyle(kFullDotMedium);
319 g2.GetXaxis()->SetTimeDisplay(true);
320 g2.GetXaxis()->SetTimeFormat("%H:%M %F1995-01-01 00:00:00 GMT");
321 g2.GetXaxis()->SetLabelSize(0.12);
322 g2.GetYaxis()->SetLabelSize(0.1);
323 g2.GetYaxis()->SetTitle("RELATIVE ON TIME");
324 g2.GetYaxis()->SetTitleOffset(0.2);
325 g2.GetYaxis()->SetTitleSize(0.1);
326 g2.DrawClone("AP");
327
328 return 0;
329}
330
331void PlotRateQC(UInt_t night, MSQLServer &serv)
332{
333 TString query =
334 "LEFT JOIN AnalysisResultsRunLP USING(fNight, fRunID) "
335 "WHERE fRunTypeKey=1 AND NOT ISNULL (AnalysisResultsRunLP.fNumEvtsAfterCleaning) AND fNight=";
336 query += night;
337
338 TTree *t = serv.GetTree("RunInfo", query);
339 if (!t)
340 return;
341
342 int save = gErrorIgnoreLevel;
343 gErrorIgnoreLevel = kFatal;
344
345 gROOT->SetSelectedPad(0);
346 gPad->GetCanvas()->cd(3);
347
348 t->Draw("AnalysisResultsRunLP.fNumEvtsAfterCleaning/AnalysisResultsRunLP.fOnTimeAfterCuts:(RunInfo.fRunStart+RunInfo.fRunStop)/2+9131*24*3600", "", "same");
349 TGraph *g = (TGraph*)gPad->GetPrimitive("Graph");
350 if (g)
351 {
352 g->SetName("CleaningRate");
353 g->SetMarkerColor(kRed);
354 g->SetMarkerStyle(29);//kFullDotMedium);
355 }
356
357 t->Draw("AnalysisResultsRunLP.fNumEvtsAfterQualCuts/AnalysisResultsRunLP.fOnTimeAfterCuts:(RunInfo.fRunStart+RunInfo.fRunStop)/2+9131*24*3600", "", "same");
358 g = (TGraph*)gPad->GetPrimitive("Graph");
359 if (g)
360 {
361 g->SetName("RateAfterQC");
362 g->SetMarkerColor(kBlue);
363 g->SetMarkerStyle(29);//kFullDotMedium);
364 }
365
366 gErrorIgnoreLevel = save;
367}
368
369Int_t PlotSqm(TArrayD **vec, TString fname)
370{
371 fname += ".SQM_CONTROL_DATA.fits";
372
373 fits file(fname.Data());
374 if (!file)
375 {
376 cerr << fname << ": " << gSystem->GetError() << endl;
377 return -2;
378 }
379
380 //cout << fname << endl;
381
382 Double_t time;
383 Float_t mag; // magnitude
384
385 if (!file.SetPtrAddress("Time", &time))
386 return -1;
387 if (!file.SetPtrAddress("Mag", &mag))
388 return -1;
389
390 //const int marker_style = kFullDotMedium;
391 const int marker_style = kDot;
392
393 TGraph g1;
394 g1.SetName("SQM");
395
396
397 bool found_first_time = false;
398 while (file.GetNextRow())
399 if (Contains(vec, time))
400 {
401 g1.SetPoint(g1.GetN(), time*24*3600, mag);
402 }
403
404
405 g1.SetMinimum(19.0);
406 g1.SetMaximum(21.5);
407 g1.SetMarkerColor(kBlack);
408 g1.SetMarkerStyle(marker_style);
409 g1.GetXaxis()->SetTimeDisplay(true);
410 g1.GetXaxis()->SetTimeFormat("%H:%M %F1995-01-01 00:00:00 GMT");
411 g1.GetXaxis()->SetLabelSize(0.12);
412 g1.GetYaxis()->SetLabelSize(0.1);
413 g1.GetYaxis()->SetTitle("SQM [mag]");
414 g1.GetYaxis()->SetTitleOffset(0.2);
415 g1.GetYaxis()->SetTitleSize(0.1);
416 g1.DrawClone("AP");
417
418 return 0;
419}
420
421Int_t PlotSqmLinear(TArrayD **vec, TString fname)
422{
423 fname += ".SQM_CONTROL_DATA.fits";
424
425 fits file(fname.Data());
426 if (!file)
427 {
428 cerr << fname << ": " << gSystem->GetError() << endl;
429 return -2;
430 }
431
432 //cout << fname << endl;
433
434 Double_t time;
435 Float_t mag; // magnitude
436
437 if (!file.SetPtrAddress("Time", &time))
438 return -1;
439 if (!file.SetPtrAddress("Mag", &mag))
440 return -1;
441
442 //const int marker_style = kFullDotMedium;
443 const int marker_style = kDot;
444
445 TGraph g1;
446 g1.SetName("SQM");
447
448
449 bool found_first_time = false;
450 while (file.GetNextRow())
451 if (Contains(vec, time))
452 {
453 Double_t mag_double = 4.4 * pow(10, 20) * pow( 10, mag*(-0.4));
454 g1.SetPoint(g1.GetN(), time*24*3600, mag_double);
455 }
456
457
458 g1.SetMinimum(1.e12);
459 g1.SetMaximum(5.e12);
460 g1.SetMarkerColor(kBlack);
461 g1.SetMarkerStyle(marker_style);
462 g1.GetXaxis()->SetTimeDisplay(true);
463 g1.GetXaxis()->SetTimeFormat("%H:%M %F1995-01-01 00:00:00 GMT");
464 g1.GetXaxis()->SetLabelSize(0.12);
465 g1.GetYaxis()->SetLabelSize(0.1);
466 g1.GetYaxis()->SetTitle("SQM lin [phot./(s sr m^2)]");
467 g1.GetYaxis()->SetTitleOffset(0.2);
468 g1.GetYaxis()->SetTitleSize(0.1);
469 g1.DrawClone("AP");
470
471 return 0;
472}
473
474Int_t PlotHumidity(TArrayD **vec, TString fname)
475{
476 fname += ".FSC_CONTROL_HUMIDITY.fits";
477
478 fits file(fname.Data());
479 if (!file)
480 {
481 cerr << fname << ": " << gSystem->GetError() << endl;
482 return -2;
483 }
484
485 //cout << fname << endl;
486
487 Double_t time;
488 Float_t H[4];
489
490 if (!file.SetPtrAddress("Time", &time))
491 return -1;
492 if (!file.SetPtrAddress("H", H))
493 return -1;
494
495 //const int marker_style = kFullDotMedium;
496 const int marker_style = kDot;
497
498 TGraph g1;
499 TGraph g2;
500 TGraph g3;
501 TGraph g4;
502 //TGraph g5;
503 g1.SetName("H0");
504 g2.SetName("H1");
505 g3.SetName("H2");
506 g4.SetName("H3");
507 //g5.SetName("PFmini");
508
509
510 Double_t first_time, last_time;
511 bool found_first_time = false;
512 while (file.GetNextRow())
513 if (Contains(vec, time))
514 {
515 if (!found_first_time){
516 first_time = time*24*3600;
517 found_first_time = true;
518 }
519 g1.SetPoint(g1.GetN(), time*24*3600, H[0]);
520 g2.SetPoint(g2.GetN(), time*24*3600, H[1]);
521 g3.SetPoint(g3.GetN(), time*24*3600, H[2]);
522 g4.SetPoint(g4.GetN(), time*24*3600, H[3]);
523 //g5.SetPoint(g5.GetN(), time*24*3600, I[319]);
524 last_time = time*24*3600;
525 }
526
527
528 g1.SetMinimum(10);
529 g1.SetMaximum(80);
530 g1.SetMarkerColor(kAzure);
531 g1.SetMarkerStyle(marker_style);
532 g1.GetXaxis()->SetTimeDisplay(true);
533 g1.GetXaxis()->SetTimeFormat("%H:%M %F1995-01-01 00:00:00 GMT");
534 g1.GetXaxis()->SetLabelSize(0.12);
535 g1.GetYaxis()->SetLabelSize(0.1);
536 g1.GetYaxis()->SetTitle("HUMITIDY");
537 g1.GetYaxis()->SetTitleOffset(0.2);
538 g1.GetYaxis()->SetTitleSize(0.1);
539 g1.DrawClone("AP");
540
541 g2.SetMarkerColor(kAzure+1);
542 g2.SetMarkerStyle(marker_style);
543 g2.DrawClone("P");
544
545 g3.SetMarkerColor(kAzure+3);
546 g3.SetMarkerStyle(marker_style);
547 g3.DrawClone("P");
548
549 g4.SetMarkerColor(kAzure+6);
550 g4.SetMarkerStyle(marker_style);
551 g4.DrawClone("P");
552
553 //g5.SetMarkerColor(kAzure+1);
554 //g5.SetMarkerStyle(kFullDotMedium);
555 //g5.DrawClone("P");
556
557 g1.DrawClone("P");
558
559 TLine l1(first_time-600, 40, last_time+600, 40);
560 l1.SetLineColor(kOrange);
561 l1.DrawClone();
562 TText t1(first_time-600, 41, "Please, note in logbook");
563 t1.SetTextSize(0.1);
564 t1.DrawClone();
565
566
567 TLine l2(first_time-600, 55, last_time+600, 55);
568 l2.SetLineColor(kRed);
569 l2.DrawClone();
570 TText t2(first_time-600, 56, "Please, report to fact-online");
571 t2.SetTextSize(0.1);
572 t2.DrawClone();
573
574
575 return 0;
576}
577
578Int_t PlotHumidity2(TArrayD **vec, TString fname)
579{
580 fname += ".PFMINI_CONTROL_DATA.fits";
581
582 fits file(fname.Data());
583 if (!file)
584 {
585 cerr << fname << ": " << gSystem->GetError() << endl;
586 return -2;
587 }
588
589 //cout << fname << endl;
590
591 Double_t time;
592 Float_t H;
593
594 if (!file.SetPtrAddress("Time", &time))
595 return -1;
596 if (!file.SetPtrAddress("Humidity", &H))
597 return -1;
598
599 const int marker_style = kFullDotMedium;
600 //const int marker_style = kDot;
601
602 TGraph g1;
603 g1.SetName("PFmini");
604
605
606 while (file.GetNextRow())
607 if (Contains(vec, time))
608 {
609 g1.SetPoint(g1.GetN(), time*24*3600, H);
610 }
611
612 g1.SetMarkerStyle(marker_style);
613 g1.SetMarkerColor(kGreen);
614 g1.DrawClone("P");
615 return 0;
616}
617
618Int_t PlotPointing(TArrayD **vec, TString fname)
619{
620 fname += ".DRIVE_CONTROL_POINTING_POSITION.fits";
621
622 fits file(fname.Data());
623 if (!file)
624 {
625 cerr << fname << ": " << gSystem->GetError() << endl;
626 return -2;
627 }
628
629 //cout << fname << endl;
630
631 Double_t time;
632 Double_t zd;
633
634 if (!file.SetPtrAddress("Time", &time))
635 return -1;
636 if (!file.SetPtrAddress("Zd", &zd))
637 return -1;
638
639 TGraph g;
640 g.SetName("Zd");
641
642 while (file.GetNextRow())
643 if (Contains(vec, time))
644 g.SetPoint(g.GetN(), time*24*3600, 90-zd);
645
646 g.SetMinimum(1);
647 g.SetMaximum(90);
648 g.SetMarkerStyle(kFullDotMedium);
649 g.GetXaxis()->SetTimeDisplay(true);
650 g.GetXaxis()->SetTimeFormat("%H:%M %F1995-01-01 00:00:00 GMT");
651 g.GetXaxis()->SetLabelSize(0.12);
652 g.GetYaxis()->SetLabelSize(0.1);
653 g.GetYaxis()->SetTitle("ELEVATION");
654 g.GetYaxis()->SetTitleOffset(0.2);
655 g.GetYaxis()->SetTitleSize(0.1);
656 g.DrawClone("AP");
657
658 return 0;
659}
660
661Int_t PlotTemperature1(TArrayD **vec, TString fname)
662{
663 fname += ".TEMPERATURE_DATA.fits";
664
665 fits file(fname.Data());
666 if (!file)
667 {
668 cerr << fname << ": " << gSystem->GetError() << endl;
669 return -2;
670 }
671
672 //cout << fname << endl;
673
674 Double_t time;
675 Float_t temp;
676
677 if (!file.SetPtrAddress("Time", &time))
678 return -1;
679 if (!file.SetPtrAddress("T", &temp))
680 return -1;
681
682 TGraph g;
683 g.SetName("ContainerTemp");
684
685 while (file.GetNextRow())
686 if (Contains(vec, time))
687 g.SetPoint(g.GetN(), time*24*3600, temp);
688
689 g.SetMinimum(-5);
690 g.SetMaximum(49);
691 g.SetMarkerStyle(kFullDotMedium);
692 g.SetMarkerColor(kRed);
693 g.GetXaxis()->SetTimeDisplay(true);
694 g.GetXaxis()->SetTimeFormat("%H:%M %F1995-01-01 00:00:00 GMT");
695 g.GetXaxis()->SetLabelSize(0.1);
696 g.GetYaxis()->SetLabelSize(0.1);
697 g.GetYaxis()->SetTitle("TEMPERATURE");
698 g.GetYaxis()->SetTitleOffset(0.2);
699 g.GetYaxis()->SetTitleSize(0.1);
700 g.DrawClone("AP");
701
702 return 0;
703}
704
705Int_t PlotTemperature2(TArrayD **vec, TString fname)
706{
707 fname += ".MAGIC_WEATHER_DATA.fits";
708
709 fits file(fname.Data());
710 if (!file)
711 {
712 cerr << fname << ": " << gSystem->GetError() << endl;
713 return -2;
714 }
715
716 //cout << fname << endl;
717
718 Double_t time;
719 Float_t temp;
720
721 if (!file.SetPtrAddress("Time", &time))
722 return -1;
723 if (!file.SetPtrAddress("T", &temp))
724 return -1;
725
726 TGraph g;
727 g.SetName("OutsideTemp");
728
729 while (file.GetNextRow())
730 if (Contains(vec, time))
731 g.SetPoint(g.GetN(), time*24*3600, temp);
732
733 g.SetMarkerStyle(kFullDotMedium);
734 g.SetMarkerColor(kBlue);
735 g.DrawClone("P");
736
737 return 0;
738}
739
740Int_t PlotTemperature3(TArrayD **vec, TString fname)
741{
742 fname += ".FSC_CONTROL_TEMPERATURE.fits";
743
744 fits file(fname.Data());
745 if (!file)
746 {
747 cerr << fname << ": " << gSystem->GetError() << endl;
748 return -2;
749 }
750
751 //cout << fname << endl;
752
753 Double_t time;
754 Float_t temp[31];
755
756 if (!file.SetPtrAddress("Time", &time))
757 return -1;
758 if (!file.SetPtrAddress("T_sens", temp))
759 return -1;
760
761 TGraph g, g1, g2;
762 g.SetName("SensorTempAvg");
763 g1.SetName("SensorTempMin");
764 g2.SetName("SensorTempMax");
765
766 while (file.GetNextRow())
767 if (Contains(vec, time))
768 {
769 float min = 100;
770 float max = -100;
771 double avg = 0;
772 int num = 0;
773 for (int i=0; i<31; i++)
774 if (temp[i]!=0)
775 {
776 avg += temp[i];
777 num++;
778
779 min = TMath::Min(min, temp[i]);
780 max = TMath::Max(max, temp[i]);
781 }
782
783 g.SetPoint(g.GetN(), time*24*3600, avg/num);
784 g1.SetPoint(g1.GetN(), time*24*3600, min);
785 g2.SetPoint(g2.GetN(), time*24*3600, max);
786 }
787
788 g.SetMarkerStyle(kFullDotMedium);
789 g.DrawClone("P");
790
791 /*
792 g1.SetLineWidth(1);
793 g1.DrawClone("L");
794
795 g2.SetLineWidth(1);
796 g2.DrawClone("L");
797 */
798 return 0;
799}
800
801Int_t PlotTemperature4(TArrayD **vec, TString fname)
802{
803 fname += ".FAD_CONTROL_TEMPERATURE.fits";
804
805 fits file(fname.Data());
806 if (!file)
807 {
808 cerr << fname << ": " << gSystem->GetError() << endl;
809 return -2;
810 }
811
812 //cout << fname << endl;
813
814 Double_t time;
815 Float_t temp[160];
816
817 if (!file.SetPtrAddress("Time", &time))
818 return -1;
819// if (!file.SetPtrAddress("Data1", temp) &&
820// !file.SetPtrAddress("temp", temp))
821 if (!file.SetPtrAddress("temp", temp))
822 return -1;
823
824 Int_t num = file.GetN("temp")==0 ? file.GetN("Data1") : file.GetN("temp");
825 Int_t beg = num==82 ? 2 : 0;
826
827 TGraphErrors g1;
828 TGraph g2,g3;
829
830 g1.SetName("FadTempAvg");
831 g2.SetName("FadTempMin");
832 g3.SetName("FadTempMax");
833
834 while (file.GetNextRow())
835 if (Contains(vec, time))
836 {
837 double avg = 0;
838 double rms = 0;
839 float min = 100;
840 float max = -100;
841 for (int i=beg; i<num; i++)
842 {
843 avg += temp[i];
844 rms += temp[i]*temp[i];
845
846 min = TMath::Min(min, temp[i]);
847 max = TMath::Max(max, temp[i]);
848 }
849
850 avg /= num-beg;
851 rms /= num-beg;
852
853 g1.SetPoint(g1.GetN(), time*24*3600, avg);
854 g1.SetPointError(g1.GetN()-1, 0, sqrt(rms-avg*avg));
855
856 g2.SetPoint(g2.GetN(), time*24*3600, min);
857 g3.SetPoint(g3.GetN(), time*24*3600, max);
858 }
859
860 g1.SetLineColor(kGreen);
861 g1.DrawClone("[]");
862
863 g2.SetLineColor(kGreen);
864 g2.SetLineWidth(1);
865 g2.DrawClone("L");
866
867 g3.SetLineColor(kGreen);
868 g3.SetLineWidth(1);
869 g3.DrawClone("L");
870
871 return 0;
872}
873
874int quality(UInt_t y=0, UInt_t m=0, UInt_t d=0, const char *outpath="quality")
875{
876 // To get correct dates in the histogram you have to add
877 // the MJDREF offset (should be 40587) and 9131.
878
879 if (y==0)
880 {
881 UInt_t nt = MTime(MTime(-1).GetMjd()-1.5).GetNightAsInt();
882 y = nt/10000;
883 m = (nt/100)%100;
884 d = nt%100;
885
886 cout << y << "/" << m << "/" << d << endl;
887 }
888
889 TString fname=Form("/fact/aux/%04d/%02d/%02d/%04d%02d%02d", y, m, d, y, m, d);
890
891 UInt_t night = MTime(y, m, d, 0).GetNightAsInt();
892
893 MSQLMagic serv("sql.rc");
894 Bool_t con = serv.IsConnected();
895
896 cout << "quality" << endl;
897 cout << "-------" << endl;
898 cout << endl;
899 if (con)
900 {
901 cout << "Connected to " << serv.GetName() << endl;
902 cout << endl;
903 }
904 cout << "Night: " << Form("%04d-%02d-%02d", y, m, d) << endl;
905 cout << endl;
906
907 TArrayD run, beg, end;
908
909 TArrayD *runs[3] = { &run, &beg, &end };
910
911 if (con)
912 {
913 TString query;
914 query += "SELECT fRunID, fRunStart, fRunStop FROM RunInfo";
915 query += " WHERE fNight=";
916 query += night;
917 query += " AND fRunTypeKey=1 ORDER BY fRunStart";
918
919 TSQLResult *res = serv.Query(query);
920 if (!res)
921 return 1;
922
923 run.Set(res->GetRowCount());
924 beg.Set(res->GetRowCount());
925 end.Set(res->GetRowCount());
926
927 Int_t n = 0;
928
929 TSQLRow *row = 0;
930 while ((row=res->Next()))
931 {
932 run[n] = atoi((*row)[0]);
933 beg[n] = MTime((*row)[1]).GetMjd()-40587;
934 end[n] = MTime((*row)[2]).GetMjd()-40587;
935 n++;
936 delete row;
937 }
938
939 delete res;
940
941 if (n==0)
942 cout << "WARNING - No data runs in db, displaying all data." << endl;
943 else
944 cout << "Num: " << n << "\n" << endl;
945 }
946
947 //check if the sqm was already installed on the telescope
948 TCanvas *c = NULL;
949 TString datestring = Form("%04d%02d%02d", y, m, d);
950 if( datestring.Atoi() > 20140723 ) {
951 TCanvas *c = new TCanvas("quality", Form("Quality %04d/%02d/%02d", y, m, d), 1280, 1280);
952 c->Divide(1, 8, 1e-5, 1e-5);
953 }else{
954 TCanvas *c = new TCanvas("quality", Form("Quality %04d/%02d/%02d", y, m, d), 1280, 1120);
955 c->Divide(1, 7, 1e-5, 1e-5);
956 }
957
958 gROOT->SetSelectedPad(0);
959
960 c->cd(1);
961 gPad->SetGrid();
962 gPad->SetTopMargin(0);
963 gPad->SetRightMargin(0.001);
964 gPad->SetLeftMargin(0.04);
965 gPad->SetBottomMargin(0);
966 cout << PlotThresholds(runs, fname) << endl;
967
968 gROOT->SetSelectedPad(0);
969 c->cd(2);
970 gPad->SetGrid();
971 gPad->SetTopMargin(0);
972 gPad->SetRightMargin(0.001);
973 gPad->SetLeftMargin(0.04);
974 gPad->SetBottomMargin(0);
975 cout << PlotCurrent(runs, fname) << endl;
976
977 gROOT->SetSelectedPad(0);
978 c->cd(3);
979 gPad->SetGrid();
980 gPad->SetTopMargin(0);
981 gPad->SetBottomMargin(0);
982 gPad->SetRightMargin(0.001);
983 gPad->SetLeftMargin(0.04);
984 cout << PlotRate(runs, fname) << endl;
985 cout << PlotRateQC(night, serv) << endl;
986
987 gROOT->SetSelectedPad(0);
988 c->cd(5);
989 gPad->SetGrid();
990 gPad->SetTopMargin(0);
991 gPad->SetBottomMargin(0);
992 gPad->SetRightMargin(0.001);
993 gPad->SetLeftMargin(0.04);
994 cout << PlotPointing(runs, fname) << endl;
995
996 gROOT->SetSelectedPad(0);
997 c->cd(6);
998 gPad->SetGrid();
999 gPad->SetTopMargin(0);
1000 gPad->SetRightMargin(0.001);
1001 gPad->SetLeftMargin(0.04);
1002 gPad->SetBottomMargin(0);
1003 cout << PlotTemperature1(runs, fname) << endl;
1004 cout << PlotTemperature2(runs, fname) << endl;
1005 cout << PlotTemperature3(runs, fname) << endl;
1006 cout << PlotTemperature4(runs, fname) << endl;
1007
1008 gROOT->SetSelectedPad(0);
1009 c->cd(7);
1010 gPad->SetGrid();
1011 gPad->SetTopMargin(0);
1012 gPad->SetBottomMargin(0);
1013 gPad->SetRightMargin(0.001);
1014 gPad->SetLeftMargin(0.04);
1015 cout << PlotHumidity(runs, fname) << endl;
1016 cout << PlotHumidity2(runs, fname) << endl;
1017
1018 //check if the sqm was already installed
1019 if( datestring.Atoi() > 20140723 ) {
1020 gROOT->SetSelectedPad(0);
1021 c->cd(7);
1022 gPad->SetGrid();
1023 gPad->SetTopMargin(0);
1024 gPad->SetBottomMargin(0);
1025 gPad->SetRightMargin(0.001);
1026 gPad->SetLeftMargin(0.04);
1027 cout << PlotHumidity(runs, fname) << endl;
1028 cout << PlotHumidity2(runs, fname) << endl;
1029
1030 gROOT->SetSelectedPad(0);
1031 c->cd(8);
1032 gPad->SetGrid();
1033 gPad->SetTopMargin(0);
1034 gPad->SetRightMargin(0.001);
1035 gPad->SetLeftMargin(0.04);
1036 cout << PlotSqm(runs, fname) << endl;
1037 }else{
1038
1039 gROOT->SetSelectedPad(0);
1040 c->cd(7);
1041 gPad->SetGrid();
1042 gPad->SetTopMargin(0);
1043 gPad->SetRightMargin(0.001);
1044 gPad->SetLeftMargin(0.04);
1045 cout << PlotHumidity(runs, fname) << endl;
1046 cout << PlotHumidity2(runs, fname) << endl;
1047 }
1048
1049 c->SaveAs(Form("%s/%04d%02d%02d.png", outpath, y, m, d));
1050
1051 return 0;
1052}
Note: See TracBrowser for help on using the repository browser.