source: branches/Corsika7500Compatibility/fact/plots/plotratescan.C@ 18454

Last change on this file since 18454 was 18205, checked in by Daniela Dorner, 9 years ago
added information and example ratescan to the plot
File size: 10.5 KB
Line 
1/*void Fit(TGraph &g, double rate_start=3e5, double rate_kink=80, double rate_end=1)
2{
3 double begin = -1;
4 double end = -1;
5
6 double kinkx = -1;
7 double kinky = 0;
8
9 TGraph gl;
10 for (int j=0; j<g.GetN(); j++)
11 {
12 gl.SetPoint(gl.GetN(), g.GetX()[j], log10(g.GetY()[j]));
13 if (g.GetY()[j]<rate_start && begin<0)
14 begin = g.GetX()[j];
15 if (g.GetY()[j]<rate_end && end<0)
16 {
17 end = g.GetX()[j];
18 //if (end<550)
19 // end = 750;
20 }
21 if (g.GetY()[j]<rate_kink && kinkx<0)
22 {
23 kinkx = g.GetX()[j];
24 kinky = g.GetY()[j];
25 }
26 }
27
28 if (end==-1)
29 end = 600;//g.GetX()[g.GetN()-1];
30
31 TF1 func2("f2", "log10("
32 " [1]*exp([2]*(x-[0]))"
33 "+ [1]*exp([3]*(x-[0]))"
34 //"+ [4]*exp([5]*(x-700))"
35 ")", begin, end);
36
37 Double_t v00 = g.Eval(kinkx-150, 0, "S");
38 Double_t v01 = g.Eval(kinkx- 50, 0, "S");
39
40 Double_t v10 = g.Eval(kinkx+100, 0, "S");
41 Double_t v11 = g.Eval(kinkx+200, 0, "S");
42
43 //Double_t v20 = g.Eval(700-50, 0, "S");
44 //Double_t v21 = g.Eval(700, 0, "S");
45 //Double_t v22 = g.Eval(700+50, 0, "S");
46
47 func2.SetParameter(0, kinkx);
48 func2.SetParameter(1, kinky);
49 func2.SetParameter(2, log(v01/v00)/100);
50 func2.SetParameter(3, log(v11/v10)/100);
51
52 func2.FixParameter(4, 0);
53 func2.FixParameter(5, 0);
54 func2.SetLineWidth(1);
55 func2.SetLineColor(kBlue);
56
57 gl.Fit(&func2, "N0", "", begin, end);
58
59 TF1 func("f1",
60 "[0]*exp([1]*(x-[6])) + "
61 "[2]*exp([3]*(x-[6])) + "
62 "[4]*exp([5]*(x-700))", begin, end);
63 func.FixParameter(0, func2.GetParameter(1));
64 func.FixParameter(1, func2.GetParameter(2));
65 func.FixParameter(2, func2.GetParameter(1));
66 func.FixParameter(3, func2.GetParameter(3));
67 func.FixParameter(4, func2.GetParameter(4));
68 func.FixParameter(5, func2.GetParameter(5));
69 func.FixParameter(6, func2.GetParameter(0));
70 func.SetLineWidth(2);
71 func.SetLineColor(kBlue);
72 func.DrawClone("same");
73 func.SetLineWidth(1);
74
75 cout << begin << " -- " << end << " : " << kinkx << endl;
76
77 func.SetLineStyle(kDashed);
78 func.SetRange(0, 1000);
79
80 func.FixParameter(0, 0);
81 func.FixParameter(2, 0);
82 func.FixParameter(4, func2.GetParameter(4));
83 func.DrawClone("same");
84
85 func.FixParameter(0, func2.GetParameter(1));
86 func.FixParameter(2, 0);
87 func.FixParameter(4, 0);
88 func.DrawClone("same");
89
90 func.FixParameter(0, 0);
91 func.FixParameter(2, func2.GetParameter(1));
92 func.FixParameter(4, 0);
93 func.DrawClone("same");
94
95
96 TLatex txt;
97 txt.SetTextColor(kBlue);
98 txt.SetTextAlign(33);
99 txt.DrawLatex(func2.GetParameter(0), func2.GetParameter(1),
100 Form("%.1fHz @ %.0f", func2.GetParameter(1), func2.GetParameter(0)));
101
102 // [1]*exp([2]*(x-[0])) / exp(1) = [1]*exp([3]*(x-[0]))
103 // exp([2]*(x-[0])) / exp(1) = exp([3]*(x-[0]))
104 // exp([2]*(x-[0])) / exp([3]*(x-[0])) = exp(1)
105 // exp([2]*(x-[0]) - [3]*(x-[0])) = exp(1)
106 // [2]*(x-[0]) - [3]*(x-[0]) = 1
107 //
108 // x = 1/([2]-[3]) + [0]
109
110
111 Double_t th = 1. / (func2.GetParameter(3)-func2.GetParameter(2)) + func2.GetParameter(0);
112 Double_t rate = pow(10, func2.Eval(th));
113 Double_t rate_end = pow(10, func2.Eval(end)-0.5);
114
115 txt.SetTextAlign(11);
116 txt.DrawLatex(th+10, rate, Form("%.1fHz @ %.0f", rate, th));
117
118 txt.SetTextAlign(12);
119 txt.DrawLatex(180, 3e6, Form("%.1f", 1./func2.GetParameter(2)));
120 txt.DrawLatex(end, rate_end, Form("%.0f", 1./func2.GetParameter(3)));
121
122 TLine line;
123 line.SetLineStyle(kDashed);
124 line.SetLineColor(kBlue);
125 line.DrawLine(th, 1e-2, th, 1e8);
126}*/
127
128Bool_t plot(MSQLMagic &serv, Long64_t search_id, TGraph &g, TGraph *gb)
129{
130 TString query;
131 query += "SELECT fTimeBegin, fTimeEnd, fVoltageIsOn, fOvervoltage, fCurrentMedMean, fNight, fAzMin, fAzMax, fZdMin, fZdMax ";
132 query += "FROM Ratescan WHERE fRatescanID=";
133 query += search_id;
134
135 TSQLResult *res = serv.Query(query);
136 if (!res)
137 return kFALSE;
138
139 if (res->GetRowCount()!=1)
140 {
141 cout << "ERROR - Row count is not 1." << endl;
142 delete res;
143 return kTRUE;
144 }
145
146 TSQLRow *row = res->Next();
147
148 const char *time_beg = (*row)[0];
149 const char *time_end = (*row)[1];
150 const char *night = (*row)[5];
151 const char *az_beg = (*row)[6];
152 const char *az_end = (*row)[7];
153 const char *zd_beg = (*row)[8];
154 const char *zd_end = (*row)[9];
155
156 int voltage_on = (*row)[2] ? atoi((*row)[2]) : -1;
157 float overvoltage = (*row)[3] ? atof((*row)[3]) : -100;
158 float current = (*row)[4] ? atof((*row)[4]) : -1;
159
160 delete row;
161
162 gROOT->SetSelectedPad(0);
163 TCanvas *c = new TCanvas(time_beg+11, time_beg);
164
165 c->SetGrid();
166 c->SetLogy();
167 c->SetTopMargin(0.01);
168 c->SetRightMargin(0.005);
169
170 TPaveText leg(600, 1e5, 1050, 8e8, "br");
171 leg.SetBorderSize(1);
172 leg.SetFillColor(kWhite);
173 leg.SetTextAlign(12);
174 leg.AddText(Form("Ratescan %d ", search_id));
175 //leg.AddText("");
176 leg.AddText(Form("Begin %s", time_beg));
177 leg.AddText(Form("End %s", time_end));
178 leg.AddText(Form("Az %s#circ to %s#circ", az_beg, az_end));
179 leg.AddText(Form("Zd %s#circ to %s#circ", zd_beg, zd_end));
180 //leg.AddText("");
181 if (voltage_on==0)
182 leg.AddText("Voltage off");
183 else
184 {
185 if (current>=0)
186 leg.AddText(Form("Current <I_{med}> = %.1f #muA", current));
187 if (overvoltage>-70)
188 leg.AddText(Form("Voltage #DeltaU = %+.2f V", overvoltage));
189 }
190
191 TH1D h0("frame", "", 1050, 0, 1050);
192 h0.SetDirectory(0);
193 h0.SetStats(kFALSE);
194 h0.SetXTitle("Threshold [dac counts]");
195 h0.SetYTitle("Trigger rate [Hz]");
196 h0.GetXaxis()->SetLabelSize(0.05);
197 h0.GetXaxis()->SetTitleSize(0.05);
198 h0.GetYaxis()->SetLabelSize(0.05);
199 h0.GetYaxis()->SetTitleSize(0.05);
200 h0.SetMaximum(8e8);
201 h0.SetMinimum(0.61);
202 h0.DrawCopy();
203
204 leg.DrawClone();
205
206 for (int i=0; i<40; i++)
207 {
208 gb[i].SetNameTitle(Form("Board%02d", i), "Board trigger rate (time 40)");
209 gb[i].DrawClone("P");
210 }
211
212 g.SetMarkerStyle(kFullDotMedium);
213 g.SetMarkerColor(kBlue);
214 g.SetLineColor(kBlue);
215 g.SetNameTitle("Camera", "Camera trigger rate");
216 g.DrawClone("PL");
217
218 TGraph gr("good_ratescan_edit.txt", "%lg %lg");
219 gr.SetLineColor(12);
220 gr.DrawClone("L");
221
222 c->Write();
223
224 TString name;
225 name += night;
226 name += "-ratescan ";
227 name += time_beg;
228
229 //c->SaveAs(name+".pdf");
230 //c->SaveAs(name+".eps");
231 //c->SaveAs("/loc_data/analysis/"+name+".png");
232 c->SaveAs(Form("/loc_data/analysis/ratescans/%04d/%02d/%02d/%06d_%d.png", atoi(night)/10000, (atoi(night)/100)%100, atoi(night)%100, atoi(night), search_id));
233
234 delete c;
235}
236
237Int_t plotid(MSQLMagic &serv, fits &file, const char *_search_id, Int_t row)
238{
239 Long64_t search_id = atol(_search_id);
240
241 cout << " " << search_id << endl;
242
243 UInt_t offset = file.GetUInt("MJDREF");
244
245 bool old = file.HasColumn("Data0");
246
247 Double_t *ptime = file.SetPtrAddress("Time");
248 UInt_t *pid = file.SetPtrAddress(old ? "Data0" : "Id");
249 Float_t *rates = file.SetPtrAddress(old ? "Data5" : "BoardRate");
250 UInt_t *pth = file.SetPtrAddress(old ? "Data1" : "Threshold");
251 Float_t *ptrig = file.SetPtrAddress(old ? "Data4" : "TriggerRate");
252 Float_t *pontime = file.SetPtrAddress(old ? "Data3" : "RelOnTime");
253
254 if (!ptime || !pid || !rates || !pth || !ptrig || !pontime)
255 return -1;
256
257 TGraph g;
258 TGraph gb[40];
259
260 while (file.GetRow(row++))
261 {
262 if (pid[0]<search_id)
263 continue;
264
265 if (pid[0]!=search_id)
266 break;
267
268 g.SetPoint(g.GetN(), pth[0], ptrig[0]/pontime[0]);
269 for (int i=0; i<40; i++)
270 gb[i].SetPoint(gb[i].GetN(), pth[0], rates[i]*40);
271 }
272
273 if (g.GetN()>0)
274 plot(serv, search_id, g, gb);
275
276 return row;
277
278
279}
280
281Bool_t plotratescan(MSQLMagic &serv, const char *_night)
282{
283 Long64_t night = atol(_night);
284
285 TString fname = Form("/fact/aux/%04d/%02d/%02d/%06d.RATE_SCAN_DATA.fits",
286 night/10000, (night/100)%100, night%100, night);
287
288 cout << " " << fname << endl;
289
290 fits file(fname.Data());
291 if (!file)
292 {
293 cout << "ERROR - Cannot access " << fname << ": " << gSystem->GetError() << endl;
294 return kFALSE;
295 }
296
297 TString query;
298 query += "SELECT fRatescanID FROM Ratescan WHERE fNight=";
299 query += night;
300 query += " ORDER BY fRatescanID";
301 /*
302 query += "SELECT fRatescanID, fTimeBeg, fTimeEnd, fVoltageIsOn, ";
303 query += " fOvervoltage, fCurrentMedMean, fRateBegin, fRateEnd, ";
304 query += " fThresholdBegin, fThresholdEnd, fNumPoints ";
305 query += "WHERE fNight=";
306 query += night;
307 */
308 TSQLResult *res = serv.Query(query);
309 if (!res)
310 return kFALSE;
311
312 if (res->GetRowCount()==0)
313 {
314 delete res;
315 return kTRUE;
316 }
317
318 TString oname = Form("%06d-ratescan.root", night);
319 //cout << " " << oname << '\n' << endl;
320 TFile rootfile(oname.Data(), "recreate");
321
322 TSQLRow *row = 0;
323
324 Int_t cnt = 0;
325 while ((row=res->Next()) && cnt>=0)
326 {
327 const char *id = (*row)[0];
328 cnt = plotid(serv, file, id, cnt);
329 delete row;
330 }
331
332 delete res;
333
334 cout << endl;
335
336 return cnt>=0;
337}
338
339Bool_t plotratescan(const char *night)
340{
341 MSQLMagic serv("sql.rc");
342 if (!serv.IsConnected())
343 {
344 cout << "ERROR - Connection to database failed." << endl;
345 return 0;
346 }
347
348 cout << "plotratescan" << endl;
349 cout << "------------" << endl;
350 cout << endl;
351 cout << "Connected to " << serv.GetName() << endl;
352 cout << "Night: " << night << endl;
353 cout << endl;
354
355 return plotratescan(serv, night);
356}
357
358Bool_t plotratescan()
359{
360 MSQLMagic serv("sql.rc");
361 if (!serv.IsConnected())
362 {
363 cout << "ERROR - Connection to database failed." << endl;
364 return 0;
365 }
366
367 cout << "plotratescan" << endl;
368 cout << "------------" << endl;
369 cout << endl;
370 cout << "Connected to " << serv.GetName() << endl;
371 cout << endl;
372
373 TSQLResult *res = serv.Query("SELECT fNight FROM Ratescan GROUP BY fNight ORDER BY fNight");
374 if (!res)
375 return kFALSE;
376
377 TSQLRow *row = 0;
378 while ((row=res->Next()))
379 {
380 const char *night = (*row)[0];
381 plotratescan(serv, night);
382 delete row;
383 }
384
385 delete res;
386
387 return kTRUE;
388}
Note: See TracBrowser for help on using the repository browser.