BmnRoot
Loading...
Searching...
No Matches
BmnTriggerEfficiencyRun8.cxx
Go to the documentation of this file.
1#include "BmnTriggerEfficiencyRun8.h"
2
3#include "BmnTrigDigit.h"
4#include "BmnTrigWaveDigit.h"
5#include "FairLogger.h"
6
7#include <TClonesArray.h>
8#include <stdlib.h>
9
10BmnTriggerEfficiencyRun8::BmnTriggerEfficiencyRun8()
11 : isValid(kTRUE)
12 , fNBd(-1)
13 , fEffFile1(nullptr)
14 , fEffFile2(nullptr)
15 , fEffFile3(nullptr)
16 , fHBd(nullptr)
17 , fHFd(nullptr)
18 , fHCctSyst(nullptr)
19 , fHFrf(nullptr)
20 , fDstFile(nullptr)
21 , fDstTree(nullptr)
22 , fBTrigInfo(nullptr)
23 , fBPrimaryVertex(nullptr)
24 , fRunId(-1.)
25 , fBC2IntCut(-1.)
26 , fFDTrigWindowLeftEdge(-1)
27 , fFDPeakLimit(-1)
28{
29 InitEfficiencyInput();
30}
31
32BmnTriggerEfficiencyRun8::BmnTriggerEfficiencyRun8(TString dstFileName)
33 : isValid(kTRUE)
34 , fNBd(-1)
35 , fEffFile1(nullptr)
36 , fEffFile2(nullptr)
37 , fEffFile3(nullptr)
38 , fHBd(nullptr)
39 , fHFd(nullptr)
40 , fHCctSyst(nullptr)
41 , fHFrf(nullptr)
42 , fDstFile(nullptr)
43 , fDstTree(nullptr)
44 , fBTrigInfo(nullptr)
45 , fBPrimaryVertex(nullptr)
46 , fRunId(-1.)
47 , fBC2IntCut(-1.)
48 , fFDTrigWindowLeftEdge(-1)
49 , fFDPeakLimit(-1)
50{
51 InitEfficiencyInput();
52
53 fDstFile = TFile::Open(dstFileName, "READONLY");
54 // fDstFile = new TFile(dstFileName, "READ");
55 if (!fDstFile || !fDstFile->IsOpen() || fDstFile->IsZombie()) {
56 LOGF(error, "Cannot open file %s", dstFileName.Data());
57 isValid = kFALSE;
58 return;
59 }
60 fDstTree = (TTree*)fDstFile->Get("bmndata");
61 if (!fDstTree) {
62 LOGF(error, "Cannot find tree bmndata in the file");
63 isValid = kFALSE;
64 return;
65 }
66
67 Int_t ret = fDstTree->SetBranchAddress("BmnTrigInfo.", &fBTrigInfo);
68 if ((ret < 0) || (ret > 2)) {
69 LOGF(error, "Cannot match branch BmnTrigInfo.");
70 isValid = kFALSE;
71 return;
72 }
73 ret = fDstTree->SetBranchAddress("PrimaryVertex.", &fBPrimaryVertex);
74 if ((ret < 0) || (ret > 2)) {
75 LOGF(error, "Cannot match branch PrimaryVertex.");
76 isValid = kFALSE;
77 return;
78 }
79
80 // Get cuts for this run id
81 Ssiz_t startCh = dstFileName.Index("mpd_run_Top_");
82 if (startCh == -1) {
83 LOGF(error, "Dst file name must contain \"mpd_run_Top_<RunId>\"");
84 isValid = kFALSE;
85 return;
86 }
87 startCh = dstFileName.Index("mpd_run_Top_") + TString("mpd_run_Top_").Sizeof() - 1;
88 fRunId = (TString(dstFileName(startCh, 4))).Atoi();
89 LOGF(info, "runId = %d", fRunId);
90
91 SetBC2IntCutByRun(fRunId);
92 SetFDTrigWindowLeftEdgeByRun(fRunId);
93 SetFDPeakLimitByRun(fRunId);
94
95 LOGF(info, "File %s successfully opened", dstFileName.Data());
96}
97
98BmnTriggerEfficiencyRun8::BmnTriggerEfficiencyRun8(Int_t runId, BmnTrigInfoDst* trigInfo, CbmVertex* primaryVertex)
99 : isValid(kTRUE)
100 , fNBd(-1)
101 , fEffFile1(nullptr)
102 , fEffFile2(nullptr)
103 , fEffFile3(nullptr)
104 , fHBd(nullptr)
105 , fHFd(nullptr)
106 , fHCctSyst(nullptr)
107 , fHFrf(nullptr)
108 , fDstFile(nullptr)
109 , fDstTree(nullptr)
110 , fBTrigInfo(trigInfo)
111 , fBPrimaryVertex(primaryVertex)
112 , fRunId(runId)
113 , fBC2IntCut(-1.)
114 , fFDTrigWindowLeftEdge(-1)
115 , fFDPeakLimit(-1)
116{
117 InitEfficiencyInput();
118
119 LOGF(info, "runId = %d", fRunId);
120
121 SetBC2IntCutByRun(fRunId);
122 SetFDTrigWindowLeftEdgeByRun(fRunId);
123 SetFDPeakLimitByRun(fRunId);
124}
125
126Double_t BmnTriggerEfficiencyRun8::GetBDEfficiency(Int_t runId, Int_t ntrPV)
127{
128 LOGF(debug, "BmnTriggerEfficiencyRun8::GetBDEfficiency(): Start");
129 return GetEfficiencyOrError(fHBd, runId, ntrPV, kTRUE);
130}
131
132Double_t BmnTriggerEfficiencyRun8::GetFDEfficiency(Int_t runId, Int_t ntrPV)
133{
134 LOGF(debug, "BmnTriggerEfficiencyRun8::GetFDEfficiency(): Start");
135 return GetEfficiencyOrError(fHFd, runId, ntrPV, kTRUE);
136}
137
138Double_t BmnTriggerEfficiencyRun8::GetTriggerEfficiency(Int_t runId, Int_t ntrPV)
139{
140 LOGF(debug, "BmnTriggerEfficiencyRun8::GetTriggerEfficiency(): Start");
141
142 Double_t bdEff = GetEfficiencyOrError(fHBd, runId, ntrPV, kTRUE);
143 if (bdEff < 0)
144 return -1.;
145
146 return bdEff * GetEfficiencyOrError(fHFd, runId, ntrPV, kTRUE);
147}
148
149Bool_t BmnTriggerEfficiencyRun8::GetTriggerEfficiencyWithErrors(Int_t runId,
150 Int_t ntrPV,
151 Double_t& eff,
152 Double_t& statErr,
153 Double_t& systErr,
154 Double_t& errSystErr)
155{
156 statErr = systErr = errSystErr = -1.;
157
158 eff = GetTriggerEfficiency(runId, ntrPV);
159 if (eff < 0.)
160 return kFALSE;
161 else if (eff == 0.) {
162 statErr = systErr = errSystErr = 0.;
163 return kTRUE;
164 }
165
166 Double_t relStatErrBD =
167 GetEfficiencyOrError(fHBd, runId, ntrPV, kFALSE) / GetEfficiencyOrError(fHBd, runId, ntrPV, kTRUE);
168 Double_t relStatErrFD =
169 GetEfficiencyOrError(fHFd, runId, ntrPV, kFALSE) / GetEfficiencyOrError(fHFd, runId, ntrPV, kTRUE);
170 statErr = eff * sqrt(relStatErrBD * relStatErrBD + relStatErrFD * relStatErrFD);
171
172 systErr = abs(GetEfficiencyOrError(fHCctSyst, runId, ntrPV, kTRUE));
173
174 errSystErr = GetEfficiencyOrError(fHCctSyst, runId, ntrPV, kFALSE);
175
176 return kTRUE;
177}
178
179Bool_t BmnTriggerEfficiencyRun8::GetTriggerEfficiencyWithErrors(Double_t& eff,
180 Double_t& statErr,
181 Double_t& systErr,
182 Double_t& errSystErr)
183{
184 if (!fBPrimaryVertex) {
185 LOGF(debug, "Primary Vertex branch should be defined");
186 return kFALSE;
187 }
188 Int_t pvntr = fBPrimaryVertex->GetNTracks();
189
190 return GetTriggerEfficiencyWithErrors(fRunId, pvntr, eff, statErr, systErr, errSystErr);
191}
192
193Bool_t BmnTriggerEfficiencyRun8::GetFluxRejectionFactorWithError(Int_t runId, Double_t& factor, Double_t& err)
194{
195 factor = err = -1.;
196
197 if (!IsValid()) {
198 LOGF(error, "This BmnTriggerEfficiencyRun8 object is not valid");
199 return kFALSE;
200 }
201
202 Int_t bin = fHFrf->FindBin(runId);
203 LOGF(debug, "BmnTriggerEfficiencyRun8::GetFluxRejectionFactorWithError(): bin = %d", bin);
204 Int_t binx, biny, binz;
205 fHFrf->GetBinXYZ(bin, binx, biny, binz);
206 LOGF(debug, "BmnTriggerEfficiencyRun8::GetFluxRejectionFactorWithError(): binx = %d, biny = %d, binz = %d", binx,
207 biny, binz);
208
209 if ((binx < 1) || (binx > fHFrf->GetNbinsX())) {
210 LOGF(error, "Run id is out of analysis range");
211 return kFALSE;
212 }
213
214 factor = fHFrf->GetBinContent(bin);
215 if (factor <= 0.)
216 return kFALSE;
217
218 err = fHFrf->GetBinError(bin);
219
220 return kTRUE;
221}
222
223Bool_t BmnTriggerEfficiencyRun8::GetFluxRejectionFactorWithError(Double_t& factor, Double_t& err)
224{
225 return GetFluxRejectionFactorWithError(fRunId, factor, err);
226}
227
228Bool_t BmnTriggerEfficiencyRun8::IsEventAnalysable()
229{
230 Double_t pvx, pvy, pvz;
231 Int_t pvntr;
232 pvx = fBPrimaryVertex->GetX();
233 pvy = fBPrimaryVertex->GetY();
234 pvz = fBPrimaryVertex->GetZ();
235 pvntr = fBPrimaryVertex->GetNTracks();
236
237 LOGF(debug, "(pvx,pvy,pvz,pvntr) = (%f,%f,%f,%d)", pvx, pvy, pvz, pvntr);
238
239 TVector3 vtxPos;
240 vtxPos.SetXYZ(0.3179, 0., -0.0169);
241 if ((abs(pvz - vtxPos.Z()) > 1.) || (pow(pvx - vtxPos.X(), 2) + pow(pvy - vtxPos.Y(), 2) > pow(1.5, 2))
242 || (pvntr < 2))
243 {
244 LOGF(debug, "Event primary vertex does not pass analysis conditions");
245 return kFALSE;
246 }
247
248 Double_t pv_chi2tondf = fBPrimaryVertex->GetChi2() / fBPrimaryVertex->GetNDF();
249 if ((pv_chi2tondf < 0.1) || (pv_chi2tondf > 10)) {
250 LOGF(debug, "Event primary vertex does not pass analysis conditions");
251 return kFALSE;
252 }
253
254 Int_t bc2inttw40 = GetBC2Int(258, 267);
255 if ((bc2inttw40 > fBC2IntCut) || (fBC2IntCut < 0)) {
256 LOGF(debug, "This event does not pass BC2 analysis conditions for pile up");
257 return kFALSE;
258 }
259
260 auto TrigPattern = fBTrigInfo->GetInputSignalsAR();
261 auto bitCCT2 = (TrigPattern >> 7) & 1;
262 if ((bitCCT2 != 1) || !IsBDTriggered() || !IsFDTriggered()) {
263 LOGF(debug, "This event does not pass trigger conditions");
264 return kFALSE;
265 }
266
267 return kTRUE;
268}
269
270Bool_t BmnTriggerEfficiencyRun8::IsEventAnalysable(Int_t evId)
271{
272 if ((evId >= fDstTree->GetEntries()) || (evId < 0)) {
273 LOGF(debug, "Event id is out of range for this file");
274 return kFALSE;
275 }
276
277 fDstTree->GetEntry(evId);
278
279 return IsEventAnalysable();
280}
281
282Double_t BmnTriggerEfficiencyRun8::GetEfficiencyOrError(TH2D* hist, Int_t runId, Int_t ntrPV, Bool_t isEff)
283{
284 if (ntrPV < 2) {
285 LOGF(error, "Number of tracks in PV is out of analysis range");
286 return -1.;
287 }
288 if (!IsValid()) {
289 LOGF(error, "This BmnTriggerEfficiencyRun8 object is not valid");
290 return -1.;
291 }
292
293 Int_t bin = hist->FindBin(runId, ntrPV);
294 LOGF(debug, "BmnTriggerEfficiencyRun8::GetEfficiencyOrError(): bin = %d", bin);
295 Int_t binx, biny, binz;
296 hist->GetBinXYZ(bin, binx, biny, binz);
297 LOGF(debug, "BmnTriggerEfficiencyRun8::GetEfficiencyOrError(): binx = %d, biny = %d, binz = %d", binx, biny, binz);
298
299 if ((binx < 1) || (binx > hist->GetNbinsX())) {
300 LOGF(error, "Run id is out of analysis range");
301 return -1.;
302 }
303
304 if (biny > hist->GetNbinsY()) {
305 if (isEff)
306 return hist->GetBinContent(binx, hist->GetNbinsY());
307 else
308 return hist->GetBinError(binx, hist->GetNbinsY());
309 }
310
311 Double_t ret;
312
313 if (isEff)
314 ret = hist->GetBinContent(bin);
315 else
316 ret = hist->GetBinError(bin);
317
318 return ret;
319}
320
321void BmnTriggerEfficiencyRun8::InitEfficiencyInput()
322{
323 TString fileName1 = TString(getenv("VMCWORKDIR")) + TString("/parameters/triggerefficiency/EffBdFd_Run8_v2.root");
324 fEffFile1 = new TFile(fileName1, "READ");
325 if (!fEffFile1->IsOpen() || fEffFile1->IsZombie()) {
326 LOGF(error, "Cannot open file %s", fileName1.Data());
327 isValid = kFALSE;
328 return;
329 }
330
331 fHBd = (TH2D*)fEffFile1->Get("hbd");
332 if (!fHBd) {
333 LOGF(error, "Cannot find histogram hbd in the file");
334 isValid = kFALSE;
335 return;
336 }
337 fHFd = (TH2D*)fEffFile1->Get("hfd");
338 if (!fHFd) {
339 LOGF(error, "Cannot find histogram hfd in the file");
340 isValid = kFALSE;
341 return;
342 }
343
344 LOGF(info, "File %s successfully opened", fileName1.Data());
345
346 TString fileName2 = TString(getenv("VMCWORKDIR")) + TString("/parameters/triggerefficiency/EffCctSyst_Run8.root");
347 fEffFile2 = new TFile(fileName2, "READ");
348 if (!fEffFile2->IsOpen() || fEffFile2->IsZombie()) {
349 LOGF(error, "Cannot open file %s", fileName2.Data());
350 isValid = kFALSE;
351 return;
352 }
353
354 fHCctSyst = (TH2D*)fEffFile2->Get("herr");
355 if (!fHCctSyst) {
356 LOGF(error, "Cannot find histogram herr in the file");
357 isValid = kFALSE;
358 return;
359 }
360
361 LOGF(info, "File %s successfully opened", fileName2.Data());
362
363 TString fileName3 =
364 TString(getenv("VMCWORKDIR")) + TString("/parameters/triggerefficiency/FluxRejectionFactor_Run8.root");
365 fEffFile3 = new TFile(fileName3, "READ");
366 if (!fEffFile3->IsOpen() || fEffFile3->IsZombie()) {
367 LOGF(error, "Cannot open file %s", fileName3.Data());
368 isValid = kFALSE;
369 return;
370 }
371
372 fHFrf = (TH1D*)fEffFile3->Get("hfrc_bt");
373 if (!fHFrf) {
374 LOGF(error, "Cannot find histogram hfrc_bt in the file");
375 isValid = kFALSE;
376 return;
377 }
378
379 LOGF(info, "File %s successfully opened", fileName3.Data());
380
381 return;
382}
383
384void BmnTriggerEfficiencyRun8::SetBC2IntCutByRun(Int_t runId)
385{
386 if (runId < 7840)
387 fBC2IntCut = 9000.;
388 else if (runId < 7984)
389 fBC2IntCut = 11000.;
390 else if (runId < 8005)
391 fBC2IntCut = 14000.;
392 else if (runId < 8011)
393 fBC2IntCut = 13000.;
394 else if (runId < 8030)
395 fBC2IntCut = 12000.;
396 else if (runId < 8060)
397 fBC2IntCut = 10000.;
398 else
399 fBC2IntCut = 11000.;
400
401 return;
402}
403
404void BmnTriggerEfficiencyRun8::SetFDTrigWindowLeftEdgeByRun(Int_t runId)
405{
406 fFDTrigWindowLeftEdge = 279 - FDTrigWindow;
407 if (runId < 7069)
408 fFDTrigWindowLeftEdge += 2;
409 else if (runId < 7114)
410 fFDTrigWindowLeftEdge += 3;
411 else if (runId < 7237)
412 fFDTrigWindowLeftEdge += 2;
413 else if (runId < 7560)
414 fFDTrigWindowLeftEdge += 0;
415 else if (runId < 7594)
416 fFDTrigWindowLeftEdge += 3;
417
418 return;
419}
420
421void BmnTriggerEfficiencyRun8::SetFDPeakLimitByRun(Int_t runId)
422{
423 fFDPeakLimit = 4250;
424 if (runId < 7069)
425 fFDPeakLimit = 3500;
426 else if (runId < 7114)
427 fFDPeakLimit = 4250;
428 else if (runId < 7148)
429 fFDPeakLimit = 3500;
430 else if (runId < 7237)
431 fFDPeakLimit = 2500;
432 else if (runId < 7258)
433 fFDPeakLimit = 3500;
434 else if (runId < 7429)
435 fFDPeakLimit = 4000;
436 else if (runId < 7486)
437 fFDPeakLimit = 3750;
438 else if (runId < 7560)
439 fFDPeakLimit = 4000;
440 else if (runId < 7594)
441 fFDPeakLimit = 3500;
442
443 return;
444}
445
446Int_t BmnTriggerEfficiencyRun8::GetBC2Int(Int_t leftEl, Int_t rightEl)
447{
448 Int_t bcIntTW = 0;
449
450 // Search for a number of BC1 signals between leftEdge and rightEdge
451 TClonesArray* bc2Digits = fBTrigInfo->GetBC2Digits();
452 if (bc2Digits->GetEntriesFast() != 1) {
453 LOGF(error, "Unexpected behaviour. bc2Digits->GetEntriesFast() = %d", bc2Digits->GetEntriesFast());
454 return -1;
455 }
456 BmnTrigWaveDigit* bc2Digit = (BmnTrigWaveDigit*)bc2Digits->UncheckedAt(0);
457 for (Int_t i = leftEl; i <= rightEl; ++i) {
458 if (bc2Digit->GetShortValue()[i] > 100)
459 bcIntTW += bc2Digit->GetShortValue()[i];
460 }
461
462 return bcIntTW;
463}
464
465Bool_t BmnTriggerEfficiencyRun8::IsBDTriggered()
466{
467 Bool_t bdTrig = kFALSE;
468
469 // Search for a number of BD signals inside the trigger window (+-50 ns)
470 TClonesArray* bdDigits = fBTrigInfo->GetBDDigits();
471 Int_t nBd = 0;
472 // Int_t aBd[40] = {0};
473 for (Int_t iBd = 0; iBd < bdDigits->GetEntriesFast(); iBd++) {
474 BmnTrigDigit* bdDigit = (BmnTrigDigit*)bdDigits->UncheckedAt(iBd);
475 if ((bdDigit->GetMod() < 0) || (bdDigit->GetMod() > 39)) {
476 LOGF(error, "Unexpected BD module number %d", bdDigit->GetMod());
477 return kFALSE;
478 }
479 // 7988 - mean=1950., delta=50.
480 // 7842 - mean=1930., delta=30.
481 if (abs(bdDigit->GetTime() - 1930.) < 30.) {
482 // aBd[bdDigit->GetMod()]++;
483 nBd++;
484 }
485 }
486 if (nBd >= 2)
487 bdTrig = kTRUE;
488
489 fNBd = nBd;
490
491 return bdTrig;
492}
493
494Bool_t BmnTriggerEfficiencyRun8::IsFDTriggered()
495{
496 Bool_t fdTrig = kFALSE;
497
498 // Search for a number of FD signals inside the trigger window (+-35 ns)
499 TClonesArray* fdDigits = fBTrigInfo->GetFDDigits();
500 // Double_t fdamp = 0.;
501 Int_t start = fFDTrigWindowLeftEdge;
502 Int_t stop = start + 2 * FDTrigWindow;
503 Int_t peak = -100000;
504 // fdamp = (Double_t)trInfo->GetFDAmp(); // Could be out of trig window
505 BmnTrigWaveDigit* fdDigit0 = (BmnTrigWaveDigit*)fdDigits->UncheckedAt(0);
506 for (Int_t i = start; i <= stop; ++i)
507 if (fdDigit0->GetShortValue()[i] > peak)
508 peak = fdDigit0->GetShortValue()[i];
509
510 // cout << "peak=" << peak << endl;
511 // 7988 (same for 7842) - 4500, 4250 - half the threshold height
512 if (peak < fFDPeakLimit)
513 fdTrig = kTRUE;
514
515 return fdTrig;
516}
sqrt
friend F32vec4 sqrt(const F32vec4 &a)
Definition P4_F32vec4.h:34
i
int i
Definition P4_F32vec4.h:22
BmnTrigDigit::GetTime
Double_t GetTime() const
Definition BmnTrigDigit.h:17
BmnTrigDigit::GetMod
Short_t GetMod() const
Definition BmnTrigDigit.h:15
BmnTrigInfoDst::GetInputSignalsAR
UInt_t GetInputSignalsAR()
Definition BmnTrigInfoDst.h:75
BmnTrigInfoDst::GetBDDigits
TClonesArray * GetBDDigits()
Definition BmnTrigInfoDst.h:49
BmnTrigInfoDst::GetBC2Digits
TClonesArray * GetBC2Digits()
Definition BmnTrigInfoDst.h:48
BmnTrigInfoDst::GetFDDigits
TClonesArray * GetFDDigits()
Definition BmnTrigInfoDst.h:50
BmnTrigWaveDigit::GetShortValue
Short_t * GetShortValue() const
Definition BmnTrigWaveDigit.h:42
BmnTriggerEfficiencyRun8::GetFDEfficiency
Double_t GetFDEfficiency(Int_t runId, Int_t ntrPV)
Definition BmnTriggerEfficiencyRun8.cxx:132
BmnTriggerEfficiencyRun8::GetTriggerEfficiency
Double_t GetTriggerEfficiency(Int_t runId, Int_t ntrPV)
Definition BmnTriggerEfficiencyRun8.cxx:138
BmnTriggerEfficiencyRun8::GetFluxRejectionFactorWithError
Bool_t GetFluxRejectionFactorWithError(Int_t runId, Double_t &factor, Double_t &err)
Definition BmnTriggerEfficiencyRun8.cxx:193
BmnTriggerEfficiencyRun8::GetTriggerEfficiencyWithErrors
Bool_t GetTriggerEfficiencyWithErrors(Int_t runId, Int_t ntrPV, Double_t &eff, Double_t &statErr, Double_t &systErr, Double_t &errSystErr)
Definition BmnTriggerEfficiencyRun8.cxx:149
BmnTriggerEfficiencyRun8::GetBDEfficiency
Double_t GetBDEfficiency(Int_t runId, Int_t ntrPV)
Definition BmnTriggerEfficiencyRun8.cxx:126
CbmVertex::GetZ
Double_t GetZ() const
Definition CbmVertex.h:60
CbmVertex::GetX
Double_t GetX() const
Definition CbmVertex.h:58
CbmVertex::GetNTracks
Int_t GetNTracks() const
Definition CbmVertex.h:63
CbmVertex::GetChi2
Double_t GetChi2() const
Definition CbmVertex.h:61
CbmVertex::GetY
Double_t GetY() const
Definition CbmVertex.h:59
CbmVertex::GetNDF
Int_t GetNDF() const
Definition CbmVertex.h:62