32 #include <THnSparse.h>
35 #include <TObjString.h>
36 #include <Riostream.h>
40 #include <TGraphAsymmErrors.h>
44 #include <THashList.h>
45 #include <TParameter.h>
62 void PlotMuonEfficiencyVsX(TString var, TString fileNameData, TString fileNameSave, Bool_t saveEdges, Bool_t print, Bool_t draw);
64 TString fileNameData, TString fileNameSave, Bool_t print, Bool_t draw);
66 void PlotMuonEfficiencyVsXY(TString xVar, TString yVar, TString fileNameData, TString fileNameSave, Bool_t draw, Bool_t rap = kFALSE);
68 void PlotMuonEfficiency(TString fileNameData, TString fileNameSave, Bool_t saveEdges, Bool_t print, Bool_t draw);
69 void PlotMuonEfficiencyVsRun(TString runList, TString fileNameData, TString fileNameSave, Bool_t print, Bool_t draw);
76 Bool_t
GetChamberEfficiency(THnSparse &TT, THnSparse &TD, TArrayD &chEff, TArrayD chEffErr[2], Bool_t printError = kFALSE);
77 void GetDEEfficiency(THnSparse &TT, THnSparse &TD, TGraphAsymmErrors &effVsDE);
78 void ComputeStationEfficiency(TArrayD &chEff, TArrayD chEffErr[2], Int_t iSt, Double_t &stEff, Double_t stEffErr[2]);
79 void GetStationEfficiency(TArrayD &chEff, TArrayD chEffErr[2], Int_t iSt, TGraphAsymmErrors *effVsX[3], Int_t ip, Double_t xp);
81 void GetStation45Efficiency(TArrayD &chEff, TArrayD chEffErr[2], TGraphAsymmErrors *effVsX[3], Int_t ip, Double_t xp);
82 void ComputeTrackingEfficiency(Double_t stEff[6], Double_t stEffErr[6][2], Double_t &spectroEff, Double_t spectroEffErr[2]);
84 TGraphAsymmErrors *effVsX[3], Int_t ip, Double_t xp, Bool_t print = kFALSE);
86 TGraphAsymmErrors &effVsX, Int_t ip, Double_t xp);
90 TGraphAsymmErrors*
CreateGraph(
const char* name,
const char*
title,
int value=-1);
91 void BeautifyGraph(TGraphAsymmErrors &g,
const char* xAxisName,
const char* yAxisName);
92 void BeautifyGraphs(TObjArray& array,
const char* xAxisName,
const char* yAxisName);
93 void SetRunLabel(TGraphAsymmErrors &g, Int_t irun,
const TList& runs);
94 void SetRunLabel(TObjArray& array, Int_t irun,
const TList& runs);
99 TString fileNameWeights =
"",
100 TString objNameExtension =
"",
101 TString fileNameData =
"AnalysisResults.root",
102 TString fileNameSave =
"efficiency_new.root")
137 void PlotMuonEfficiencyVsX(TString var, TString fileNameData, TString fileNameSave, Bool_t saveEdges, Bool_t print, Bool_t draw)
141 printf(
"plotting efficiency versus %s...\n", var.Data());
144 if (var ==
"centrality") xDim = 1;
145 else if (var ==
"pt") xDim = 2;
146 else if (var ==
"y") xDim = 3;
147 else if (var ==
"phi") xDim = 4;
148 else if (var ==
"charge") xDim = 5;
150 printf(
"incorrect variable. Choices are centrality, pt, y, phi and charge.\n");
155 TFile *
file =
new TFile(fileNameData.Data(),
"read");
156 if (!file || !file->IsOpen()) {
157 printf(
"cannot open file %s \n",fileNameData.Data());
160 TList *listTT =
static_cast<TList*
>(file->FindObjectAny(Form(
"TotalTracksPerChamber%s",
gObjNameExtension.Data())));
161 TList *listTD =
static_cast<TList*
>(file->FindObjectAny(Form(
"TracksDetectedPerChamber%s",
gObjNameExtension.Data())));
162 THnSparse *TT =
static_cast<THnSparse*
>(listTT->At(10));
163 THnSparse *TD =
static_cast<THnSparse*
>(listTD->At(10));
166 TGraphAsymmErrors *effVsX[3] = {0x0, 0x0, 0x0};
167 TString nameAdd[3] = {
"",
"Low",
"Up"};
168 TString titleAdd[3] = {
"",
" - lower limit",
" - upper limit"};
169 for (Int_t i = 0; i < 1 || (saveEdges && i < 3); ++i) {
170 effVsX[i] =
CreateGraph(Form(
"trackingEffVs%s%s",var.Data(),nameAdd[i].Data()),
171 Form(
"Measured tracking efficiency versus %s%s",var.Data(),titleAdd[i].Data()));
183 for (Int_t ix = 1; ix <= TT->GetAxis(xDim)->GetNbins(); ix++) {
185 if (print) cout << var.Data() <<
" " << TT->GetAxis(xDim)->GetBinLowEdge(ix) <<
"-" << TT->GetAxis(xDim)->GetBinUpEdge(ix) <<
":" << endl;
188 TT->GetAxis(xDim)->SetRange(ix, ix);
189 TD->GetAxis(xDim)->SetRange(ix, ix);
195 TGraphAsymmErrors *dummy[3] = {0x0, 0x0, 0x0};
196 GetTrackingEfficiency(chEff, chEffErr, dummy, effVsX, ix-1, TT->GetAxis(xDim)->GetBinCenter(ix), print);
201 effVsX[0]->SetPoint(ix-1,TT->GetAxis(xDim)->GetBinCenter(ix),1.);
202 effVsX[0]->SetPointError(ix-1,0.,0.,1.,0.);
207 effVsX[1]->SetPoint(ix-1,TT->GetAxis(xDim)->GetBinCenter(ix),0.);
208 effVsX[1]->SetPointError(ix-1,0.,0.,0.,0.);
211 effVsX[2]->SetPoint(ix-1,TT->GetAxis(xDim)->GetBinCenter(ix),1.);
212 effVsX[2]->SetPointError(ix-1,0.,0.,0.,0.);
224 for (Int_t i = 0; i < 1 || (saveEdges && i < 3); ++i)
BeautifyGraph(*effVsX[i], var.Data(),
"efficiency");
227 new TCanvas(Form(
"cTrackingEffVs%s",var.Data()), Form(
"Measured tracking efficiency versus %s",var.Data()),1000,400);
228 effVsX[0]->DrawClone(
"ap");
232 file =
new TFile(fileNameSave.Data(),
"update");
233 for (Int_t i = 0; i < 1 || (saveEdges && i < 3); ++i) effVsX[i]->Write(0x0, TObject::kOverwrite);
237 for (Int_t i = 0; i < 1 || (saveEdges && i < 3); ++i)
delete effVsX[i];
244 TString fileNameData, TString fileNameSave, Bool_t print, Bool_t draw)
248 printf(
"plotting integrated efficiency versus %s...\n", var.Data());
253 printf(
"Cannot compute integrated efficiency without run-by-run weights\n");
258 ifstream inFile(runList.Data());
259 if (!inFile.is_open()) {
260 printf(
"cannot open file %s\n",runList.Data());
265 TGraphAsymmErrors *intEffVsX =
CreateGraph(Form(
"integratedTrackingEffVs%s",var.Data()),
266 Form(
"Integrated tracking efficiency versus %s",var.Data()));
275 while (!inFile.eof()) {
279 currRun.ReadLine(inFile,kTRUE);
280 if(currRun.IsNull() || !currRun.IsDec())
continue;
281 Int_t run = currRun.Atoi();
283 printf(
"run %d: ", run);
286 TString dataFile = Form(
"runs/%d/%s", run, fileNameData.Data());
287 TString outFile = Form(
"runs/%d/%s", run, fileNameSave.Data());
291 TParameter<Double_t> *weight =
static_cast<TParameter<Double_t>*
>(
runWeights->FindObject(currRun.Data()));
293 printf(
"weight not found for run %s\n", currRun.Data());
296 Double_t w = weight->GetVal();
300 TFile *
file =
new TFile(outFile.Data(),
"read");
301 if (!file || !file->IsOpen()) {
302 printf(
"cannot open file\n");
305 TGraphAsymmErrors *effVsX[2];
306 effVsX[0] =
static_cast<TGraphAsymmErrors*
>(file->FindObjectAny(Form(
"trackingEffVs%sLow",var.Data())));
307 effVsX[1] =
static_cast<TGraphAsymmErrors*
>(file->FindObjectAny(Form(
"trackingEffVs%sUp",var.Data())));
308 if (!effVsX[0] || !effVsX[1]) {
309 printf(
"trackingEffVs%sLow(Up) object not found\n", var.Data());
315 n = effVsX[0]->GetN();
316 x.Set(n, effVsX[0]->GetX());
317 for (Int_t i = 0; i < 2; ++i) {
321 }
else if (n != effVsX[0]->GetN()) {
322 printf(
"number of points in graph trackingEffVs%sLow(Up) for run %d is different than from previous runs\n", var.Data(), run);
328 for (Int_t ix = 0; ix < n; ++ix) {
329 Double_t ieffErr[2] = {effVsX[0]->GetErrorYlow(ix), effVsX[1]->GetErrorYhigh(ix)};
330 for (Int_t i = 0; i < 2; ++i) {
331 rec[i][ix] += w*effVsX[i]->GetY()[ix];
332 effErr[i][ix] += w2*ieffErr[i]*ieffErr[i];
346 for (Int_t ix = 0; ix < n; ++ix) {
348 intEffVsX->SetPoint(ix, x[ix], rec[1][ix]/gen);
349 intEffVsX->SetPointError(ix, 0., 0., (rec[1][ix]-rec[0][ix]+TMath::Sqrt(effErr[0][ix]))/gen, TMath::Sqrt(effErr[1][ix])/gen);
355 for (Int_t ix = 0; ix < n; ++ix) {
357 printf(
"impossible to integrate, all weights = 0 or unknown ?!?\n");
359 intEffVsX->SetPoint(ix, x[ix], -1.);
360 intEffVsX->SetPointError(ix, 0., 0., 0., 0.);
370 new TCanvas(Form(
"cIntegratedTrackingEffVs%s",var.Data()), Form(
"Integrated tracking efficiency versus %s",var.Data()),1000,400);
371 intEffVsX->DrawClone(
"ap");
375 TFile *
file =
new TFile(fileNameSave.Data(),
"update");
376 intEffVsX->Write(0x0, TObject::kOverwrite);
386 void PlotMuonEfficiencyVsXY(TString xVar, TString yVar, TString fileNameData, TString fileNameSave, Bool_t draw, Bool_t rap)
390 printf(
"plotting efficiency versus %s/%s...\n", xVar.Data(), yVar.Data());
393 if (xVar ==
"centrality") xDim = 1;
394 else if (xVar ==
"pt") xDim = 2;
395 else if (xVar ==
"y") xDim = 3;
396 else if (xVar ==
"phi") xDim = 4;
397 else if (xVar ==
"charge") xDim = 5;
399 printf(
"incorrect variable. Choices are centrality, pt, y, phi and charge.\n");
403 if (yVar ==
"centrality") yDim = 1;
404 else if (yVar ==
"pt") yDim = 2;
405 else if (yVar ==
"y") yDim = 3;
406 else if (yVar ==
"phi") yDim = 4;
407 else if (yVar ==
"charge") yDim = 5;
409 printf(
"incorrect variable. Choices are centrality, pt, y, phi and charge.\n");
414 TFile *
file =
new TFile(fileNameData.Data(),
"read");
415 if (!file || !file->IsOpen()) {
416 printf(
"cannot open file %s\n",fileNameData.Data());
419 TList *listTT =
static_cast<TList*
>(file->FindObjectAny(Form(
"TotalTracksPerChamber%s",
gObjNameExtension.Data())));
420 TList *listTD =
static_cast<TList*
>(file->FindObjectAny(Form(
"TracksDetectedPerChamber%s",
gObjNameExtension.Data())));
421 THnSparse *TT =
static_cast<THnSparse*
>(listTT->At(10));
422 THnSparse *TD =
static_cast<THnSparse*
>(listTD->At(10));
425 Int_t nxBins = TT->GetAxis(xDim)->GetNbins();
426 Int_t nyBins = TT->GetAxis(yDim)->GetNbins();
427 TH2F *effVsXY =
new TH2F(Form(
"trackingEffVs%s-%s",xVar.Data(),yVar.Data()),
428 Form(
"Measured tracking efficiency versus %s and %s",xVar.Data(),yVar.Data()),
429 nxBins, TT->GetAxis(xDim)->GetBinLowEdge(1), TT->GetAxis(xDim)->GetBinUpEdge(nxBins),
430 nyBins, TT->GetAxis(yDim)->GetBinLowEdge(1), TT->GetAxis(yDim)->GetBinUpEdge(nyBins));
431 effVsXY->SetDirectory(0);
442 for (Int_t ix = 1; ix <= nxBins; ++ix) {
445 TT->GetAxis(xDim)->SetRange(ix, ix);
446 TD->GetAxis(xDim)->SetRange(ix, ix);
448 for (Int_t iy = 1; iy <= nyBins; ++iy) {
451 TT->GetAxis(yDim)->SetRange(iy, iy);
452 TD->GetAxis(yDim)->SetRange(iy, iy);
458 TGraphAsymmErrors *dummy[3] = {0x0, 0x0, 0x0};
462 effVsXY->Fill(TT->GetAxis(xDim)->GetBinCenter(ix),TT->GetAxis(yDim)->GetBinCenter(iy),chEff[0]);
463 effVsXY->SetBinError(ix,iy,TMath::Max(chEffErr[0][0], chEffErr[1][0]));
468 effVsXY->Fill(TT->GetAxis(xDim)->GetBinCenter(ix),TT->GetAxis(yDim)->GetBinCenter(iy),0.);
469 effVsXY->SetBinError(ix,iy,1.);
481 effVsXY->GetXaxis()->SetTitle(xVar.Data());
482 effVsXY->GetXaxis()->CenterTitle(kTRUE);
483 effVsXY->GetXaxis()->SetLabelFont(22);
484 effVsXY->GetXaxis()->SetTitleFont(22);
485 effVsXY->GetYaxis()->SetTitle(yVar.Data());
486 effVsXY->GetYaxis()->CenterTitle(kTRUE);
487 effVsXY->GetYaxis()->SetLabelFont(22);
488 effVsXY->GetYaxis()->SetTitleFont(22);
489 effVsXY->GetZaxis()->SetTitle(
"efficiency");
490 effVsXY->GetZaxis()->SetLabelFont(22);
491 effVsXY->GetZaxis()->SetTitleFont(22);
494 new TCanvas(Form(
"cTrackingEffVs%s-%s",xVar.Data(),yVar.Data()), Form(
"Measured tracking efficiency versus %s and %s",xVar.Data(),yVar.Data()),700,600);
495 effVsXY->DrawClone(
"surf1");
499 file =
new TFile(fileNameSave.Data(),
"update");
500 effVsXY->Write(0x0, TObject::kOverwrite);
503 if (yDim == 3 && rap) {
505 TH2F* effVsXYrap =
new TH2F();
506 TString rapName = Form(
"trackingEffVs%s-%sRapBins", xVar.Data(), yVar.Data());
507 TString rapTitle = Form(
"Measured tracking efficiency versus %s and %s", xVar.Data(), yVar.Data());
508 effVsXYrap->SetTitle(rapTitle.Data());
509 effVsXYrap->SetName(rapName.Data());
511 Double_t xBinEdge[nxBins+1];
512 for (Int_t xbin = 0; xbin <= nxBins; ++xbin)
513 xBinEdge[xbin] = effVsXY->GetXaxis()->GetBinLowEdge(xbin+1);
515 Double_t yBinEdge[nyBins+1];
516 for (Int_t ybin = 0; ybin <= nyBins; ++ybin)
517 yBinEdge[ybin] = 2*TMath::ATan(TMath::Exp((effVsXY->GetYaxis()->GetBinLowEdge(ybin+1))));
519 effVsXYrap->SetBins(nxBins, xBinEdge, nyBins, yBinEdge);
521 for (Int_t xbin = 1; xbin <= nxBins; ++xbin)
522 for (Int_t ybin = 1; ybin <= nyBins; ++ybin)
523 effVsXYrap->SetBinContent(xbin, ybin, effVsXY->GetBinContent(xbin,ybin));
525 effVsXYrap->Write(0x0, TObject::kOverwrite);
540 void PlotMuonEfficiency(TString fileNameData, TString fileNameSave, Bool_t saveEdges, Bool_t print, Bool_t draw)
544 printf(
"plotting efficiency...\n");
547 TFile *
file =
new TFile(fileNameData.Data(),
"read");
548 if (!file || !file->IsOpen()) {
549 printf(
"cannot open file %s \n",fileNameData.Data());
552 TList *listTT =
static_cast<TList*
>(file->FindObjectAny(Form(
"TotalTracksPerChamber%s",
gObjNameExtension.Data())));
553 TList *listTD =
static_cast<TList*
>(file->FindObjectAny(Form(
"TracksDetectedPerChamber%s",
gObjNameExtension.Data())));
554 THnSparse *TT =
static_cast<THnSparse*
>(listTT->At(10));
555 THnSparse *TD =
static_cast<THnSparse*
>(listTD->At(10));
558 TGraphAsymmErrors *effVsCh[3] = {0x0, 0x0, 0x0};
559 TGraphAsymmErrors *effVsSt[3] = {0x0, 0x0, 0x0};
560 TString nameAdd[3] = {
"",
"Low",
"Up"};
561 TString titleAdd[3] = {
"",
" - lower limit",
" - upper limit"};
562 for (Int_t i = 0; i < 1 || (saveEdges && i < 3); ++i) {
564 Form(
"Measured efficiency per chamber (0 = spectro)%s",titleAdd[i].
Data()));
566 Form(
"Measured efficiency per station (6 = st4&5)%s",titleAdd[i].
Data()));
587 effVsCh[0]->SetPoint(0,0.,1.);
588 effVsCh[0]->SetPointError(0,0.,0.,1.,0.);
593 effVsCh[1]->SetPoint(0,0.,0.);
594 effVsCh[1]->SetPointError(0,0.,0.,0.,0.);
597 effVsCh[2]->SetPoint(0,0.,1.);
598 effVsCh[2]->SetPointError(0,0.,0.,0.,0.);
605 for (Int_t iCh = 1; iCh < 11; iCh++) {
606 for (Int_t i = 0; i < 1 || (saveEdges && i < 3); ++i) {
607 effVsCh[i]->SetPoint(iCh,iCh,chEff[iCh]);
608 effVsCh[i]->SetPointError(iCh,0.,0.,chEffErr[0][iCh],chEffErr[1][iCh]);
616 for (Int_t i = 0; i < 1 || (saveEdges && i < 3); ++i) {
617 effVsCh[i]->GetXaxis()->Set(12, -0.5, 10.5);
618 effVsCh[i]->GetXaxis()->SetNdivisions(11);
620 effVsSt[i]->GetXaxis()->Set(7, 0.5, 6.5);
621 effVsSt[i]->GetXaxis()->SetNdivisions(6);
626 TCanvas *c =
new TCanvas(
"cEfficiency",
"Measured tracking efficiency" , 1000, 400);
628 gROOT->SetSelectedPad(c->cd(1));
629 effVsCh[0]->DrawClone(
"ap");
630 gROOT->SetSelectedPad(c->cd(2));
631 effVsSt[0]->DrawClone(
"ap");
635 file =
new TFile(fileNameSave.Data(),
"update");
636 for (Int_t i = 0; i < 1 || (saveEdges && i < 3); ++i) effVsCh[i]->Write(0x0, TObject::kOverwrite);
637 for (Int_t i = 0; i < 1 || (saveEdges && i < 3); ++i) effVsSt[i]->Write(0x0, TObject::kOverwrite);
641 for (Int_t i = 0; i < 1 || (saveEdges && i < 3); ++i) {
654 printf(
"plotting efficiency versus run...\n");
657 ifstream inFile(runList.Data());
658 if (!inFile.is_open()) {
659 printf(
"cannot open file %s\n",runList.Data());
664 TObjArray chamberVsRunGraphs;
665 chamberVsRunGraphs.SetOwner(kTRUE);
666 for ( Int_t iCh = 1; iCh < 11; ++iCh)
667 chamberVsRunGraphs.Add(
CreateGraph(
"effCh%dVsRun",
"Measured efficiency for chamber %d versus run",iCh));
669 TObjArray stationVsRunGraphs[3];
670 TGraphAsymmErrors *trkVsRun[3];
671 TString nameAdd[3] = {
"",
"Low",
"Up"};
672 TString titleAdd[3] = {
"",
" - lower limit",
" - upper limit"};
673 for (Int_t i = 0; i < 3; ++i) {
675 stationVsRunGraphs[i].SetOwner(kTRUE);
676 for ( Int_t iSt = 1; iSt < 6; ++iSt)
677 stationVsRunGraphs[i].Add(
CreateGraph(Form(
"effSt%%dVsRun%s",nameAdd[i].
Data()),
678 Form(
"Measured efficiency for station %%d versus run%s",titleAdd[i].
Data()),iSt));
679 stationVsRunGraphs[i].Add(
CreateGraph(Form(
"effSt4&5VsRun%s",nameAdd[i].
Data()),
680 Form(
"Measured efficiency for station 4&5 versus run%s",titleAdd[i].
Data())));
682 trkVsRun[i] =
CreateGraph(Form(
"trackingEffVsRun%s",nameAdd[i].
Data()),
683 Form(
"Measured tracking efficiency versus run%s", titleAdd[i].
Data()));
692 while (!inFile.eof()) {
696 currRun.ReadLine(inFile,kTRUE);
697 if(currRun.IsNull())
continue;
698 runs.AddLast(
new TObjString(currRun));
701 Int_t run = currRun.Atoi();
703 printf(
"run %d: ", run);
706 TString dataFile = Form(
"runs/%d/%s", run, fileNameData.Data());
707 TString outFile = Form(
"runs/%d/%s", run, fileNameSave.Data());
710 TFile *
file =
new TFile(outFile.Data(),
"read");
711 if (!file || !file->IsOpen()) {
712 printf(
"cannot open file\n");
717 for (Int_t i = 0; i < 3; ++i) {
720 TGraphAsymmErrors *effVsCh;
721 TGraphAsymmErrors *effVsSt;
722 effVsCh =
static_cast<TGraphAsymmErrors*
>(file->FindObjectAny(Form(
"chamberEff%s",nameAdd[i].
Data())));
723 effVsSt =
static_cast<TGraphAsymmErrors*
>(file->FindObjectAny(Form(
"stationEff%s",nameAdd[i].
Data())));
724 if (!effVsCh || !effVsSt) {
725 printf(
"efficiency graph not found\n");
730 if (i == 0)
for ( Int_t iCh = 0; iCh < 10; ++iCh) {
731 TGraphAsymmErrors *g =
static_cast<TGraphAsymmErrors*
>(chamberVsRunGraphs.UncheckedAt(iCh));
732 g->SetPoint(irun,irun,effVsCh->GetY()[iCh+1]);
733 g->SetPointError(irun,0.,0.,effVsCh->GetErrorYlow(iCh+1),effVsCh->GetErrorYhigh(iCh+1));
737 for ( Int_t iSt = 0; iSt < 6; ++iSt) {
738 TGraphAsymmErrors *g =
static_cast<TGraphAsymmErrors*
>(stationVsRunGraphs[i].UncheckedAt(iSt));
739 g->SetPoint(irun,irun,effVsSt->GetY()[iSt]);
740 g->SetPointError(irun,0.,0.,effVsSt->GetErrorYlow(iSt),effVsSt->GetErrorYhigh(iSt));
744 trkVsRun[i]->SetPoint(irun,irun,effVsCh->GetY()[0]);
745 trkVsRun[i]->SetPointError(irun,0.,0.,effVsCh->GetErrorYlow(0),effVsCh->GetErrorYhigh(0));
758 for (Int_t i = 0; i < 3; ++i) {
766 new TCanvas(
"cTrackingEffVsRun",
"Tracking efficiency versus run",1000,400);
767 trkVsRun[0]->DrawClone(
"ap");
771 TFile*
file =
new TFile(fileNameSave.Data(),
"update");
772 chamberVsRunGraphs.Write(
"ChamberEffVsRun", TObject::kOverwrite | TObject::kSingleKey);
773 for (Int_t i = 0; i < 3; ++i)
774 stationVsRunGraphs[i].Write(Form(
"StationEffVsRun%s",nameAdd[i].
Data()), TObject::kOverwrite | TObject::kSingleKey);
775 for (Int_t i = 0; i < 3; ++i) trkVsRun[i]->Write(0x0, TObject::kOverwrite);
779 for (Int_t i = 0; i < 3; ++i)
delete trkVsRun[i];
789 printf(
"plotting integrated efficiency...\n");
794 printf(
"Cannot compute integrated efficiency without run-by-run weights\n");
799 TFile *
file =
new TFile(fileNameSave.Data(),
"update");
800 if (!file || !file->IsOpen()) {
801 printf(
"cannot open file\n");
804 TObjArray *chamberVsRunGraphs =
static_cast<TObjArray*
>(file->FindObjectAny(
"ChamberEffVsRun"));
805 TObjArray *stationVsRunGraphs[2];
806 stationVsRunGraphs[0] =
static_cast<TObjArray*
>(file->FindObjectAny(
"StationEffVsRunLow"));
807 stationVsRunGraphs[1] =
static_cast<TObjArray*
>(file->FindObjectAny(
"StationEffVsRunUp"));
808 TGraphAsymmErrors *trkVsRun[2];
809 trkVsRun[0] =
static_cast<TGraphAsymmErrors*
>(file->FindObjectAny(
"trackingEffVsRunLow"));
810 trkVsRun[1] =
static_cast<TGraphAsymmErrors*
>(file->FindObjectAny(
"trackingEffVsRunUp"));
811 if (!chamberVsRunGraphs || !stationVsRunGraphs[0] || !stationVsRunGraphs[1] || !trkVsRun[0] || !trkVsRun[1]) {
812 printf(
"object not found --> you must first plot the efficiency versus run\n");
817 TGraphAsymmErrors *effVsCh =
CreateGraph(
"integratedChamberEff",
"Integrated efficiency per chamber (0 = spectro)");
818 TGraphAsymmErrors *effVsSt =
CreateGraph(
"integratedStationEff",
"Integrated efficiency per station (6 = st4&5)");
824 for ( Int_t iCh = 0; iCh < 10; ++iCh) {
825 TGraphAsymmErrors *g =
static_cast<TGraphAsymmErrors*
>(chamberVsRunGraphs->UncheckedAt(iCh));
830 for ( Int_t iSt = 0; iSt < 6; ++iSt) {
831 TGraphAsymmErrors *gLow =
static_cast<TGraphAsymmErrors*
>(stationVsRunGraphs[0]->UncheckedAt(iSt));
832 TGraphAsymmErrors *gUp =
static_cast<TGraphAsymmErrors*
>(stationVsRunGraphs[1]->UncheckedAt(iSt));
838 for (Int_t iCh = 1; iCh < 11; ++iCh) {
839 cout <<
"Efficiency chamber " << iCh <<
" : ";
840 cout << effVsCh->GetY()[iCh] <<
" + " << effVsCh->GetErrorYhigh(iCh) <<
" - " << effVsCh->GetErrorYlow(iCh) << endl;
842 for (Int_t iSt = 0; iSt < 6; ++iSt) {
843 if (iSt < 5) cout <<
"Station " << iSt+1 <<
" = ";
844 else cout <<
"Station 45 = ";
845 cout << effVsSt->GetY()[iSt] <<
" + " << effVsSt->GetErrorYhigh(iSt) <<
" - " << effVsSt->GetErrorYlow(iSt) << endl;
847 cout <<
"Total tracking efficiency : ";
848 cout << effVsCh->GetY()[0] <<
" + " << effVsCh->GetErrorYhigh(0) <<
" - " << effVsCh->GetErrorYlow(0) << endl << endl;
852 effVsCh->GetXaxis()->Set(12, -0.5, 10.5);
853 effVsCh->GetXaxis()->SetNdivisions(11);
855 effVsSt->GetXaxis()->Set(7, 0.5, 6.5);
856 effVsSt->GetXaxis()->SetNdivisions(6);
860 TCanvas *c =
new TCanvas(
"cIntegratedEfficiency",
"Integrated tracking efficiency" , 1000, 400);
862 gROOT->SetSelectedPad(c->cd(1));
863 effVsCh->DrawClone(
"ap");
864 gROOT->SetSelectedPad(c->cd(2));
865 effVsSt->DrawClone(
"ap");
869 effVsCh->Write(0x0, TObject::kOverwrite);
870 effVsSt->Write(0x0, TObject::kOverwrite);
885 printf(
"plotting efficiency per DE...\n");
888 TFile *
file =
new TFile(fileNameData.Data(),
"read");
889 if (!file || !file->IsOpen()) {
890 printf(
"cannot open file %s \n",fileNameData.Data());
893 TList *listTT =
static_cast<TList*
>(file->FindObjectAny(Form(
"TotalTracksPerChamber%s",
gObjNameExtension.Data())));
894 TList *listTD =
static_cast<TList*
>(file->FindObjectAny(Form(
"TracksDetectedPerChamber%s",
gObjNameExtension.Data())));
897 TObjArray chamberVsDEGraphs;
898 chamberVsDEGraphs.SetOwner(kTRUE);
899 TGraphAsymmErrors *gCh[10];
901 TObjArray stationVsDEGraphs[3];
902 TString nameAdd[3] = {
"",
"Low",
"Up"};
903 TString titleAdd[3] = {
"",
" - lower limit",
" - upper limit"};
904 for (Int_t i = 0; i < 3; ++i) stationVsDEGraphs[i].SetOwner(kTRUE);
905 TGraphAsymmErrors *gSt[6][3];
906 for (Int_t iSt = 0; iSt < 6; ++iSt)
for (Int_t i = 0; i < 3; ++i) gSt[iSt][i] = 0x0;
909 for (Int_t iCh = 0; iCh < 10; ++iCh) {
912 chamberVsDEGraphs.Add(
CreateGraph(
"effCh%dVsDE",
"Measured efficiency for chamber %d per DE",iCh+1));
913 gCh[iCh] =
static_cast<TGraphAsymmErrors*
>(chamberVsDEGraphs.UncheckedAt(iCh));
916 THnSparse *TT =
static_cast<THnSparse*
>(listTT->At(iCh));
917 THnSparse *TD =
static_cast<THnSparse*
>(listTD->At(iCh));
937 for (Int_t iSt = 0; iSt < 3; ++iSt) {
940 for (Int_t i = 0; i < 1 || (saveEdges && i < 3); ++i) {
941 stationVsDEGraphs[i].Add(
CreateGraph(Form(
"effSt%%dVsDE%s",nameAdd[i].
Data()),
942 Form(
"Measured efficiency for station %%d per DE%s",titleAdd[i].
Data()),iSt+1));
943 gSt[iSt][i] =
static_cast<TGraphAsymmErrors*
>(stationVsDEGraphs[i].UncheckedAt(iSt));
947 Int_t
nDE = gCh[2*iSt]->GetN();
948 for (Int_t iDE = 0; iDE <
nDE; ++iDE) {
951 for (Int_t iCh = 0; iCh < 2; ++iCh) {
952 chEff[2*iSt+iCh+1] = gCh[2*iSt+iCh]->GetY()[iDE];
953 chEffErr[0][2*iSt+iCh+1] = gCh[2*iSt+iCh]->GetErrorYlow(iDE);
954 chEffErr[1][2*iSt+iCh+1] = gCh[2*iSt+iCh]->GetErrorYhigh(iDE);
965 for (Int_t i = 0; i < 1 || (saveEdges && i < 3); ++i) {
966 for (Int_t iSt = 3; iSt < 5; ++iSt) {
967 stationVsDEGraphs[i].Add(
CreateGraph(Form(
"effSt%%dVsDE%s",nameAdd[i].
Data()),
968 Form(
"Measured efficiency for station %%d per DE%s",titleAdd[i].
Data()),iSt+1));
969 gSt[iSt][i] =
static_cast<TGraphAsymmErrors*
>(stationVsDEGraphs[i].UncheckedAt(iSt));
971 stationVsDEGraphs[i].Add(
CreateGraph(Form(
"effSt4&5VsDE%s",nameAdd[i].
Data()),
972 Form(
"Measured efficiency for station 4&5 per DE%s",titleAdd[i].
Data())));
973 gSt[5][i] =
static_cast<TGraphAsymmErrors*
>(stationVsDEGraphs[i].UncheckedAt(5));
977 Int_t
nDE = gCh[6]->GetN();
978 for (Int_t iDE = 0; iDE <
nDE; ++iDE) {
981 for (Int_t iCh = 6; iCh < 10; ++iCh) {
982 chEff[iCh+1] = gCh[iCh]->GetY()[iDE];
983 chEffErr[0][iCh+1] = gCh[iCh]->GetErrorYlow(iDE);
984 chEffErr[1][iCh+1] = gCh[iCh]->GetErrorYhigh(iDE);
988 for (Int_t iSt = 3; iSt < 5; ++iSt)
GetStationEfficiency(chEff, chEffErr, iSt, gSt[iSt], iDE, iDE);
996 for (Int_t iCh = 0; iCh < 10; ++iCh) {
997 Int_t nDE = gCh[iCh]->GetN();
998 gCh[iCh]->GetXaxis()->Set(nDE+1, -0.5, nDE-0.5);
999 gCh[iCh]->GetXaxis()->SetNdivisions(nDE);
1001 BeautifyGraphs(chamberVsDEGraphs,
"Detection Element",
"efficiency");
1002 for (Int_t i = 0; i < 1 || (saveEdges && i < 3); ++i) {
1003 for (Int_t iSt = 0; iSt < 6; ++iSt) {
1004 Int_t nDE = gSt[iSt][i]->GetN();
1005 gSt[iSt][i]->GetXaxis()->Set(nDE+1, -0.5, nDE-0.5);
1006 gSt[iSt][i]->GetXaxis()->SetNdivisions(nDE);
1008 BeautifyGraphs(stationVsDEGraphs[i],
"Detection Element",
"efficiency");
1013 file =
new TFile(fileNameSave.Data(),
"update");
1014 chamberVsDEGraphs.Write(
"ChamberEffVsDE", TObject::kOverwrite | TObject::kSingleKey);
1015 for (Int_t i = 0; i < 1 || (saveEdges && i < 3); ++i)
1016 stationVsDEGraphs[i].Write(Form(
"StationEffVsDE%s",nameAdd[i].
Data()), TObject::kOverwrite | TObject::kSingleKey);
1027 printf(
"plotting efficiency per DE versus run...\n");
1030 ifstream inFile(runList.Data());
1031 if (!inFile.is_open()) {
1032 printf(
"cannot open file %s\n",runList.Data());
1037 TObjArray deVsRunGraphs;
1038 deVsRunGraphs.SetOwner(kTRUE);
1039 TObjArray stDEVsRunGraphs[3];
1040 TString nameAdd[3] = {
"",
"Low",
"Up"};
1041 TString titleAdd[3] = {
"",
" - lower limit",
" - upper limit"};
1042 for (Int_t i = 0; i < 3; ++i) stDEVsRunGraphs[i].SetOwner(kTRUE);
1043 Bool_t createGraph = kTRUE;
1050 while (!inFile.eof()) {
1054 currRun.ReadLine(inFile,kTRUE);
1055 if(currRun.IsNull())
continue;
1056 runs.AddLast(
new TObjString(currRun));
1059 Int_t run = currRun.Atoi();
1061 printf(
"run %d: ", run);
1064 TString dataFile = Form(
"runs/%d/%s", run, fileNameData.Data());
1065 TString outFile = Form(
"runs/%d/%s", run, fileNameSave.Data());
1068 TFile *
file =
new TFile(outFile.Data(),
"read");
1069 if (!file || !file->IsOpen()) {
1070 printf(
"cannot open file\n");
1075 TObjArray *chamberVsDEGraphs =
static_cast<TObjArray*
>(file->FindObjectAny(
"ChamberEffVsDE"));
1076 if (!chamberVsDEGraphs) {
1077 printf(
"efficiency graph not found\n");
1081 Int_t currentDE = 0;
1084 for ( Int_t iCh = 0; iCh < 10; ++iCh) {
1087 TGraphAsymmErrors *g =
static_cast<TGraphAsymmErrors*
>(chamberVsDEGraphs->UncheckedAt(iCh));
1088 Int_t
nDE = g->GetN();
1091 for (Int_t iDE = 0; iDE <
nDE; ++iDE) {
1094 if (createGraph) deVsRunGraphs.Add(
CreateGraph(
"effDE%dVsRun",
"Measured efficiency for DE %d versus run",100*(iCh+1)+iDE));
1095 TGraphAsymmErrors *gDE=
static_cast<TGraphAsymmErrors*
>(deVsRunGraphs.UncheckedAt(currentDE++));
1098 gDE->SetPoint(irun,irun,g->GetY()[iDE]);
1099 gDE->SetPointError(irun,0.,0.,g->GetErrorYlow(iDE),g->GetErrorYhigh(iDE));
1106 for (Int_t i = 0; i < 3; ++i) {
1109 TObjArray *stationVsDEGraphs =
static_cast<TObjArray*
>(file->FindObjectAny(Form(
"StationEffVsDE%s",nameAdd[i].
Data())));
1110 if (!stationVsDEGraphs) {
1111 printf(
"efficiency graph not found\n");
1115 Int_t currentStDE = 0;
1118 for ( Int_t iSt = 0; iSt < 6; ++iSt) {
1121 TGraphAsymmErrors *g =
static_cast<TGraphAsymmErrors*
>(stationVsDEGraphs->UncheckedAt(iSt));
1122 Int_t
nDE = g->GetN();
1125 for (Int_t iDE = 0; iDE <
nDE; ++iDE) {
1129 TString sSt = (iSt<5) ? Form(
"%d",iSt+1) :
"4&5";
1130 stDEVsRunGraphs[i].Add(
CreateGraph(Form(
"effSt%sDE%%dVsRun%s",sSt.Data(),nameAdd[i].Data()),
1131 Form(
"Measured efficiency for DE %%d in station %s versus run%s",sSt.Data(),titleAdd[i].Data()),iDE));
1133 TGraphAsymmErrors *gDE=
static_cast<TGraphAsymmErrors*
>(stDEVsRunGraphs[i].UncheckedAt(currentStDE++));
1136 gDE->SetPoint(irun,irun,g->GetY()[iDE]);
1137 gDE->SetPointError(irun,0.,0.,g->GetErrorYlow(iDE),g->GetErrorYhigh(iDE));
1147 createGraph = kFALSE;
1156 for (Int_t i = 0; i < 3; ++i) {
1162 TFile*
file =
new TFile(fileNameSave.Data(),
"update");
1163 deVsRunGraphs.Write(
"DEEffVsRun", TObject::kOverwrite | TObject::kSingleKey);
1164 for (Int_t i = 0; i < 3; ++i)
1165 stDEVsRunGraphs[i].Write(Form(
"DEEffPerStationVsRun%s",nameAdd[i].
Data()), TObject::kOverwrite | TObject::kSingleKey);
1176 printf(
"plotting integrated efficiency per DE...\n");
1181 printf(
"Cannot compute integrated efficiency without run-by-run weights\n");
1186 TFile *
file =
new TFile(fileNameSave.Data(),
"update");
1187 if (!file || !file->IsOpen()) {
1188 printf(
"cannot open file\n");
1191 TObjArray *deVsRunGraphs =
static_cast<TObjArray*
>(file->FindObjectAny(
"DEEffVsRun"));
1192 TObjArray *stDEVsRunGraphs[2];
1193 stDEVsRunGraphs[0] =
static_cast<TObjArray*
>(file->FindObjectAny(
"DEEffPerStationVsRunLow"));
1194 stDEVsRunGraphs[1] =
static_cast<TObjArray*
>(file->FindObjectAny(
"DEEffPerStationVsRunUp"));
1195 if (!deVsRunGraphs || !stDEVsRunGraphs[0] || !stDEVsRunGraphs[1]) {
1196 printf(
"object not found --> you must first plot the efficiency versus run\n");
1201 TObjArray chamberVsDEGraphs;
1202 chamberVsDEGraphs.SetOwner(kTRUE);
1203 for ( Int_t iCh = 0; iCh < 10; ++iCh)
1204 chamberVsDEGraphs.Add(
CreateGraph(
"integratedEffCh%dVsDE",
"Integrated efficiency for chamber %d per DE",iCh+1));
1206 TObjArray stationVsDEGraphs;
1207 stationVsDEGraphs.SetOwner(kTRUE);
1208 for ( Int_t iSt = 0; iSt < 5; ++iSt)
1209 stationVsDEGraphs.Add(
CreateGraph(
"integratedEffSt%dVsDE",
"Integrated efficiency for station %d per DE",iSt+1));
1210 stationVsDEGraphs.Add(
CreateGraph(
"integratedEffSt4&5VsDE",
"Integrated efficiency for station 4&5 per DE"));
1213 TIter nextDE(deVsRunGraphs);
1214 TGraphAsymmErrors *gDE = 0x0;
1215 while ((gDE = static_cast<TGraphAsymmErrors*>(nextDE()))) {
1219 sscanf(gDE->GetName(),
"effDE%dVsRun", &deId);
1220 Int_t iCh = deId/100-1;
1221 Int_t iDE = deId%100;
1224 TGraphAsymmErrors *g =
static_cast<TGraphAsymmErrors*
>(chamberVsDEGraphs.UncheckedAt(iCh));
1230 Int_t ng = stDEVsRunGraphs[0]->GetEntries();
1231 for (Int_t ig = 0; ig < ng; ++ig) {
1233 TGraphAsymmErrors *gDELow =
static_cast<TGraphAsymmErrors*
>(stDEVsRunGraphs[0]->UncheckedAt(ig));
1234 TGraphAsymmErrors *gDEUp =
static_cast<TGraphAsymmErrors*
>(stDEVsRunGraphs[1]->UncheckedAt(ig));
1238 if (strstr(gDELow->GetName(),
"4&5")) {
1240 sscanf(gDELow->GetName(),
"effSt4&5DE%dVsRunLow", &iDE);
1242 sscanf(gDELow->GetName(),
"effSt%dDE%dVsRunLow", &iSt, &iDE);
1247 TGraphAsymmErrors *g =
static_cast<TGraphAsymmErrors*
>(stationVsDEGraphs.UncheckedAt(iSt));
1253 for ( Int_t iCh = 0; iCh < 10; ++iCh) {
1254 TGraphAsymmErrors *g =
static_cast<TGraphAsymmErrors*
>(chamberVsDEGraphs.UncheckedAt(iCh));
1255 Int_t
nDE = g->GetN();
1256 g->GetXaxis()->Set(nDE+1, -0.5, nDE-0.5);
1257 g->GetXaxis()->SetNdivisions(nDE);
1259 BeautifyGraphs(chamberVsDEGraphs,
"Detection Element",
"efficiency");
1260 for ( Int_t iSt = 0; iSt < 6; ++iSt) {
1261 TGraphAsymmErrors *g =
static_cast<TGraphAsymmErrors*
>(stationVsDEGraphs.UncheckedAt(iSt));
1262 Int_t
nDE = g->GetN();
1263 g->GetXaxis()->Set(nDE+1, -0.5, nDE-0.5);
1264 g->GetXaxis()->SetNdivisions(nDE);
1266 BeautifyGraphs(stationVsDEGraphs,
"Detection Element",
"efficiency");
1269 chamberVsDEGraphs.Write(
"IntegratedChamberEffVsDE", TObject::kOverwrite | TObject::kSingleKey);
1270 stationVsDEGraphs.Write(
"IntegratedStationEffVsDE", TObject::kOverwrite | TObject::kSingleKey);
1283 TH1D *TTdraw = TT.Projection(0,
"e");
1284 TH1D *TDdraw = TD.Projection(0,
"e");
1287 TGraphAsymmErrors *efficiency = (TTdraw->GetEntries() > 0) ?
new TGraphAsymmErrors(TDdraw, TTdraw,Form(
"%se0",
effErrMode)) : 0x0;
1288 Bool_t ok = (efficiency);
1293 Bool_t missingEff = kFALSE;
1295 for (Int_t i = 0; i < 10; i++) {
1297 if (TTdraw->GetBinContent(i+1) > 0) {
1299 chEff[i+1] = efficiency->GetY()[i];
1300 chEffErr[0][i+1] = efficiency->GetErrorYlow(i);
1301 chEffErr[1][i+1] = efficiency->GetErrorYhigh(i);
1306 chEffErr[0][i+1] = 0.;
1307 chEffErr[1][i+1] = 0.;
1315 if (missingEff && printError) cout <<
"efficiency partially unknown" << endl;
1319 for (Int_t i = 0; i < 10; i++) {
1322 chEffErr[0][i+1] = 0.;
1323 chEffErr[1][i+1] = 0.;
1327 if (printError) cout <<
"efficiency unknown" << endl << endl;
1347 TH1D *TTdraw = TT.Projection(0,
"e");
1348 TH1D *TDdraw = TD.Projection(0,
"e");
1351 TGraphAsymmErrors *efficiency = (TTdraw->GetEntries() > 0) ?
new TGraphAsymmErrors(TDdraw, TTdraw,Form(
"%se0",
effErrMode)) : 0x0;
1352 Int_t
nDE = TTdraw->GetNbinsX();
1356 for (Int_t iDE = 0; iDE <
nDE; ++iDE) {
1358 if (TTdraw->GetBinContent(iDE+1) > 0) {
1360 effVsDE.SetPoint(iDE,iDE,efficiency->GetY()[iDE]);
1361 effVsDE.SetPointError(iDE,0,0,efficiency->GetErrorYlow(iDE),efficiency->GetErrorYhigh(iDE));
1365 effVsDE.SetPoint(iDE,iDE,-1);
1366 effVsDE.SetPointError(iDE,0,0,0,0);
1374 for (Int_t iDE = 0; iDE <
nDE; ++iDE) {
1376 effVsDE.SetPoint(iDE,iDE,-1);
1377 effVsDE.SetPointError(iDE,0,0,0,0);
1396 stEff = 1.-(1.-chEff[2*iSt+1])*(1.-chEff[2*iSt+2]);
1398 Double_t d1 = (1. - chEff[2*iSt+2]); d1 *= d1;
1399 Double_t d2 = (1. - chEff[2*iSt+1]); d2 *= d2;
1401 for (Int_t i = 0; i < 2; ++i) {
1402 Double_t s1 = chEffErr[i][2*iSt+1] * chEffErr[i][2*iSt+1];
1403 Double_t s2 = chEffErr[i][2*iSt+2] * chEffErr[i][2*iSt+2];
1404 stEffErr[i] = TMath::Sqrt(d1*s1 + d2*s2 + s1*s2);
1407 stEffErr[0] = TMath::Min(stEff, stEffErr[0]);
1408 stEffErr[1] = TMath::Min(1.-stEff, stEffErr[1]);
1414 void GetStationEfficiency(TArrayD &chEff, TArrayD chEffErr[2], Int_t iSt, TGraphAsymmErrors *effVsX[3], Int_t ip, Double_t xp)
1418 if (chEff[2*iSt+1] >= 0 && chEff[2*iSt+2] >= 0) {
1421 Double_t stEff, stEffErr[2];
1425 for (Int_t i = 0; i < 3; ++i) {
1427 effVsX[i]->SetPoint(ip,xp,stEff);
1428 effVsX[i]->SetPointError(ip,0.,0.,stEffErr[0],stEffErr[1]);
1434 Double_t edge[2] = {0., 1.};
1435 TArrayD chEffEdge[2];
1436 Double_t stEffEdge[2];
1437 Double_t stEffEdgeErr[2][2];
1439 for (Int_t i = 0; i < 2; ++i) {
1442 chEffEdge[i].Set(11);
1443 for (Int_t iCh = 1; iCh < 3; ++iCh) chEffEdge[i][2*iSt+iCh] = (chEff[2*iSt+iCh] < 0) ? edge[i] : chEff[2*iSt+iCh];
1450 effVsX[i+1]->SetPoint(ip,xp,stEffEdge[i]);
1451 effVsX[i+1]->SetPointError(ip,0.,0.,stEffEdgeErr[i][0],stEffEdgeErr[i][1]);
1458 effVsX[0]->SetPoint(ip,xp,stEffEdge[1]);
1459 effVsX[0]->SetPointError(ip,0.,0.,stEffEdge[1]-stEffEdge[0]+stEffEdgeErr[0][0],stEffEdgeErr[1][1]);
1472 st45Eff = chEff[7]*chEff[8]*chEff[9] + chEff[7]*chEff[8]*chEff[10] + chEff[7]*chEff[9]*chEff[10] + chEff[8]*chEff[9]*chEff[10] - 3.*chEff[7]*chEff[8]*chEff[9]*chEff[10];
1474 Double_t d1 = chEff[8]*chEff[9] + chEff[8]*chEff[10] + chEff[9]*chEff[10] - 3.*chEff[8]*chEff[9]*chEff[10]; d1 *= d1;
1475 Double_t d2 = chEff[7]*chEff[9] + chEff[7]*chEff[10] + chEff[9]*chEff[10] - 3.*chEff[7]*chEff[9]*chEff[10]; d2 *= d2;
1476 Double_t d3 = chEff[7]*chEff[8] + chEff[7]*chEff[10] + chEff[8]*chEff[10] - 3.*chEff[7]*chEff[8]*chEff[10]; d3 *= d3;
1477 Double_t d4 = chEff[7]*chEff[8] + chEff[7]*chEff[9] + chEff[8]*chEff[9] - 3.*chEff[7]*chEff[8]*chEff[9]; d4 *= d4;
1478 Double_t d12 = chEff[9] + chEff[10] - 3.*chEff[9]*chEff[10]; d12 *= d12;
1479 Double_t d13 = chEff[8] + chEff[10] - 3.*chEff[8]*chEff[10]; d13 *= d13;
1480 Double_t d14 = chEff[8] + chEff[9] - 3.*chEff[8]*chEff[9]; d14 *= d14;
1481 Double_t d23 = chEff[7] + chEff[10] - 3.*chEff[7]*chEff[10]; d23 *= d23;
1482 Double_t d24 = chEff[7] + chEff[9] - 3.*chEff[7]*chEff[9]; d24 *= d24;
1483 Double_t d34 = chEff[7] + chEff[8] - 3.*chEff[7]*chEff[8]; d34 *= d34;
1484 Double_t d123 = 1. - 3.*chEff[10]; d123 *= d123;
1485 Double_t d124 = 1. - 3.*chEff[9]; d124 *= d124;
1486 Double_t d134 = 1. - 3.*chEff[8]; d134 *= d134;
1487 Double_t d234 = 1. - 3.*chEff[7]; d234 *= d234;
1488 Double_t d1234 = - 3.; d1234 *= d1234;
1490 for (Int_t i = 0; i < 2; ++i) {
1491 Double_t s1 = chEffErr[i][7] * chEffErr[i][7];
1492 Double_t s2 = chEffErr[i][8] * chEffErr[i][8];
1493 Double_t s3 = chEffErr[i][9] * chEffErr[i][9];
1494 Double_t s4 = chEffErr[i][10] * chEffErr[i][10];
1495 st45EffErr[i] = TMath::Sqrt(d1*s1 + d2*s2 + d3*s3 + d4*s4 + d12*s1*s2 + d13*s1*s3 + d14*s1*s4 + d23*s2*s3 + d24*s2*s4 + d34*s3*s4 + d123*s1*s2*s3 + d124*s1*s2*s4 + d134*s1*s3*s4 + d234*s2*s3*s4 + d1234*s1*s2*s3*s4);
1498 st45EffErr[0] = TMath::Min(st45Eff, st45EffErr[0]);
1499 st45EffErr[1] = TMath::Min(1.-st45Eff, st45EffErr[1]);
1509 if (chEff[7] >= 0 && chEff[8] >= 0 && chEff[9] >= 0 && chEff[10] >= 0) {
1512 Double_t stEff, stEffErr[2];
1516 for (Int_t i = 0; i < 3; ++i) {
1518 effVsX[i]->SetPoint(ip,xp,stEff);
1519 effVsX[i]->SetPointError(ip,0.,0.,stEffErr[0],stEffErr[1]);
1525 Double_t edge[2] = {0., 1.};
1526 TArrayD chEffEdge[2];
1527 Double_t stEffEdge[2];
1528 Double_t stEffEdgeErr[2][2];
1530 for (Int_t i = 0; i < 2; ++i) {
1533 chEffEdge[i].Set(11);
1534 for (Int_t iCh = 7; iCh < 11; ++iCh) chEffEdge[i][iCh] = (chEff[iCh] < 0) ? edge[i] : chEff[iCh];
1541 effVsX[i+1]->SetPoint(ip,xp,stEffEdge[i]);
1542 effVsX[i+1]->SetPointError(ip,0.,0.,stEffEdgeErr[i][0],stEffEdgeErr[i][1]);
1549 effVsX[0]->SetPoint(ip,xp,stEffEdge[1]);
1550 effVsX[0]->SetPointError(ip,0.,0.,stEffEdge[1]-stEffEdge[0]+stEffEdgeErr[0][0],stEffEdgeErr[1][1]);
1564 for (Int_t iSt = 0; iSt < 6; iSt++) de[iSt][0] = stEff[iSt]*stEff[iSt];
1568 spectroEff = stEff[0] * stEff[1] * stEff[2] * stEff[3] * stEff[4];
1570 for (Int_t i = 0; i < 2; i++) {
1572 for (Int_t iSt = 0; iSt < 6; iSt++) de[iSt][1] = stEffErr[iSt][i]*stEffErr[iSt][i];
1574 spectroEffErr[i] = 0.;
1575 for (Int_t j = 1; j < 32; j++) {
1576 Double_t sigmaAdd = 1.;
1577 for (Int_t iSt = 0; iSt < 5; iSt++) sigmaAdd *= de[iSt][TESTBIT(j,iSt)];
1578 spectroEffErr[i] += sigmaAdd;
1580 spectroEffErr[i] = TMath::Sqrt(spectroEffErr[i]);
1586 spectroEff = stEff[0] * stEff[1] * stEff[2] * stEff[5];
1588 for (Int_t i = 0; i < 2; i++) {
1590 for (Int_t iSt = 0; iSt < 6; iSt++) de[iSt][1] = stEffErr[iSt][i]*stEffErr[iSt][i];
1592 spectroEffErr[i] = 0.;
1593 for (Int_t j = 1; j < 16; j++) {
1594 Double_t sigmaAdd = de[5][TESTBIT(j,3)];
1595 for (Int_t iSt = 0; iSt < 3; iSt++) sigmaAdd *= de[iSt][TESTBIT(j,iSt)];
1596 spectroEffErr[i] += sigmaAdd;
1598 spectroEffErr[i] = TMath::Sqrt(spectroEffErr[i]);
1604 spectroEffErr[0] = TMath::Min(spectroEff, spectroEffErr[0]);
1605 spectroEffErr[1] = TMath::Min(1.-spectroEff, spectroEffErr[1]);
1612 TGraphAsymmErrors *effVsX[3], Int_t ip, Double_t xp, Bool_t print)
1617 Double_t stEffErr[6][2];
1620 Bool_t allEffKnown = kTRUE;
1621 for (Int_t iCh = 1; iCh < 11 && allEffKnown; ++iCh)
if (chEff[iCh] < 0) allEffKnown = kFALSE;
1630 Double_t spectroEff, spectroEffErr[2];
1632 chEff[0] = spectroEff;
1633 chEffErr[0][0] = spectroEffErr[0];
1634 chEffErr[1][0] = spectroEffErr[1];
1637 for (Int_t i = 0; i < 3; ++i) {
1639 effVsX[i]->SetPoint(ip,xp,chEff[0]);
1640 effVsX[i]->SetPointError(ip,0.,0.,chEffErr[0][0],chEffErr[1][0]);
1642 if (effVsSt[i])
for (Int_t iSt = 0; iSt < 6; ++iSt) {
1643 effVsSt[i]->SetPoint(iSt,iSt+1,stEff[iSt]);
1644 effVsSt[i]->SetPointError(iSt,0.,0.,stEffErr[iSt][0],stEffErr[iSt][1]);
1650 Double_t edge[2] = {0., 1.};
1651 TArrayD chEffEdge[2];
1652 Double_t stEffEdge[2][6];
1653 Double_t stEffEdgeErr[2][6][2];
1654 Double_t spectroEffEdge[2], spectroEffEdgeErr[2][2];
1656 for (Int_t i = 0; i < 2; ++i) {
1659 chEffEdge[i].Set(11);
1660 for (Int_t iCh = 1; iCh < 11; ++iCh) chEffEdge[i][iCh] = (chEff[iCh] < 0) ? edge[i] : chEff[iCh];
1663 for (Int_t iSt = 0; iSt < 5; ++iSt)
ComputeStationEfficiency(chEffEdge[i], chEffErr, iSt, stEffEdge[i][iSt], stEffEdgeErr[i][iSt]);
1671 effVsX[i+1]->SetPoint(ip,xp,spectroEffEdge[i]);
1672 effVsX[i+1]->SetPointError(ip,0.,0.,spectroEffEdgeErr[i][0],spectroEffEdgeErr[i][1]);
1674 if (effVsSt[i+1])
for (Int_t iSt = 0; iSt < 6; ++iSt) {
1675 effVsSt[i+1]->SetPoint(iSt,iSt+1,stEffEdge[i][iSt]);
1676 effVsSt[i+1]->SetPointError(iSt,0.,0.,stEffEdgeErr[i][iSt][0],stEffEdgeErr[i][iSt][1]);
1682 for (Int_t iSt = 0; iSt < 6; ++iSt) {
1683 stEff[iSt] = stEffEdge[1][iSt];
1684 stEffErr[iSt][0] = stEffEdge[1][iSt] - stEffEdge[0][iSt] + stEffEdgeErr[0][iSt][0];
1685 stEffErr[iSt][1] = stEffEdgeErr[1][iSt][1];
1689 chEff[0] = spectroEffEdge[1];
1690 chEffErr[0][0] = spectroEffEdge[1] - spectroEffEdge[0] + spectroEffEdgeErr[0][0];
1691 chEffErr[1][0] = spectroEffEdgeErr[1][1];
1695 effVsX[0]->SetPoint(ip,xp,chEff[0]);
1696 effVsX[0]->SetPointError(ip,0.,0.,chEffErr[0][0],chEffErr[1][0]);
1698 if (effVsSt[0])
for (Int_t iSt = 0; iSt < 6; ++iSt) {
1699 effVsSt[0]->SetPoint(iSt,iSt+1,stEff[iSt]);
1700 effVsSt[0]->SetPointError(iSt,0.,0.,stEffErr[iSt][0],stEffErr[iSt][1]);
1707 for (Int_t iCh = 1; iCh < 11; ++iCh) {
1708 cout <<
"Efficiency chamber " << iCh <<
" : ";
1709 cout << chEff[iCh] <<
" + " << chEffErr[1][iCh] <<
" - " << chEffErr[0][iCh] << endl;
1711 for (Int_t iSt = 0; iSt < 6; ++iSt) {
1712 if (iSt < 5) cout <<
"Station " << iSt+1 <<
" = ";
1713 else cout <<
"Station 45 = ";
1714 cout << stEff[iSt] <<
" + " << stEffErr[iSt][1] <<
" - " << stEffErr[iSt][0] << endl;
1716 cout <<
"Total tracking efficiency : ";
1717 cout << chEff[0] <<
" + " << chEffErr[1][0] <<
" - " << chEffErr[0][0] << endl << endl;
1725 TGraphAsymmErrors &effVsX, Int_t ip, Double_t xp)
1731 effVsX.SetPoint(ip,xp,-1.);
1732 effVsX.SetPointError(ip,0.,0.,0.,0.);
1737 Double_t rec[2] = {0., 0.};
1739 Double_t intEffErr[2] = {0., 0.};
1743 Int_t nRuns = effVsRunLow.GetN();
1744 for (Int_t iRun = 0; iRun < nRuns; ++iRun) {
1747 TString sRun = effVsRunLow.GetXaxis()->GetBinLabel(iRun+1);
1748 TParameter<Double_t> *weight =
static_cast<TParameter<Double_t>*
>(
runWeights->FindObject(sRun.Data()));
1750 printf(
"weight not found for run %s\n", sRun.Data());
1753 Double_t w = weight->GetVal();
1757 Double_t eff[2] = {effVsRunLow.GetY()[iRun], effVsRunUp.GetY()[iRun]};
1758 Double_t effErr[2] = {effVsRunLow.GetErrorYlow(iRun), effVsRunUp.GetErrorYhigh(iRun)};
1759 if (eff[0] < 0. || eff[1] < 0.) {
1760 printf(
"no efficiency measurement --> use 0(1) ± 0 as lower(upper) limit\n");
1769 for (Int_t i = 0; i < 2; ++i) {
1771 intEffErr[i] += w2*effErr[i]*effErr[i];
1777 if (gen > 0. && ok) {
1779 effVsX.SetPoint(ip,xp,rec[1]/gen);
1780 effVsX.SetPointError(ip,0.,0.,(rec[1]-rec[0]+TMath::Sqrt(intEffErr[0]))/gen,TMath::Sqrt(intEffErr[1])/gen);
1784 if (gen <= 0.) printf(
"impossible to integrate, all weights = 0 or unknown ?!?\n");
1785 else printf(
"efficiency never measured --> return -1 ± 0\n");
1787 effVsX.SetPoint(ip,xp,-1.);
1788 effVsX.SetPointError(ip,0.,0.,0.,0.);
1804 ifstream inFile(fileName.Data());
1805 if (!inFile.is_open()) {
1806 printf(
"cannot open file %s\n", fileName.Data());
1814 while (! inFile.eof() ) {
1816 line.ReadLine(inFile,kTRUE);
1817 if(line.IsNull())
continue;
1819 TObjArray *param = line.Tokenize(
" ");
1820 if (param->GetEntries() != 2) {
1821 printf(
"bad input line %s", line.Data());
1825 Int_t run = ((TObjString*)param->UncheckedAt(0))->String().Atoi();
1827 printf(
"invalid run number: %d", run);
1831 Float_t weight = ((TObjString*)param->UncheckedAt(1))->String().Atof();
1833 printf(
"invalid weight: %g", weight);
1837 runWeights->Add(
new TParameter<Double_t>(Form(
"%d",run), weight));
1854 printf(
"Wrong Pt range, ptMax must be higher than ptMin\n");
1858 if (charge < -1 || charge > 1) {
1859 printf(
"Selected charge must be 0, 1 or -1\n");
1867 Int_t lowBin = SparseData.GetAxis(1)->FindBin(
centMin);
1868 Int_t upBin = SparseData.GetAxis(1)->FindBin(
centMax);
1869 SparseData.GetAxis(1)->SetRange(lowBin, upBin);
1872 lowBin = SparseData.GetAxis(2)->FindBin(
ptMin);
1873 if (
ptMax == -1) { upBin = SparseData.GetAxis(2)->GetNbins()+1; }
1874 else { upBin = SparseData.GetAxis(2)->FindBin(
ptMax);}
1875 SparseData.GetAxis(2)->SetRange(lowBin, upBin);
1879 lowBin = SparseData.GetAxis(5)->FindBin(
charge);
1880 SparseData.GetAxis(5)->SetRange(lowBin, lowBin);
1891 TGraphAsymmErrors* g =
new TGraphAsymmErrors();
1895 g->SetName(Form(name,value));
1896 g->SetTitle(Form(title,value));
1912 void BeautifyGraph(TGraphAsymmErrors &g,
const char* xAxisName,
const char* yAxisName)
1918 g.SetMarkerStyle(20);
1919 g.SetMarkerSize(0.7);
1920 g.SetMarkerColor(2);
1922 g.GetXaxis()->SetTitle(xAxisName);
1923 g.GetXaxis()->CenterTitle(kTRUE);
1924 g.GetXaxis()->SetLabelFont(22);
1925 g.GetXaxis()->SetTitleFont(22);
1927 g.GetYaxis()->SetTitle(yAxisName);
1928 g.GetYaxis()->CenterTitle(kTRUE);
1929 g.GetYaxis()->SetLabelFont(22);
1930 g.GetYaxis()->SetTitleFont(22);
1940 for ( Int_t i = 0; i <= array.GetLast(); ++i)
1941 BeautifyGraph(*static_cast<TGraphAsymmErrors*>(array.UncheckedAt(i)), xAxisName, yAxisName);
1951 g.GetXaxis()->Set(irun+1, -0.5, irun+0.5);
1953 TIter nextRun(&runs);
1954 TObjString *srun = 0x0;
1956 while ((srun = static_cast<TObjString*>(nextRun()))) {
1957 g.GetXaxis()->SetBinLabel(iirun, srun->GetName());
1969 for (Int_t i = 0; i <= array.GetLast(); ++i)
1970 SetRunLabel(*static_cast<TGraphAsymmErrors*>(array.UncheckedAt(i)), irun, runs);
void ComputeStationEfficiency(TArrayD &chEff, TArrayD chEffErr[2], Int_t iSt, Double_t &stEff, Double_t stEffErr[2])
void MuonTrackingEfficiency(TString runList="runList.txt", TString fileNameWeights="", TString objNameExtension="", TString fileNameData="AnalysisResults.root", TString fileNameSave="efficiency_new.root")
void GetStationEfficiency(TArrayD &chEff, TArrayD chEffErr[2], Int_t iSt, TGraphAsymmErrors *effVsX[3], Int_t ip, Double_t xp)
void GetTrackingEfficiency(TArrayD &chEff, TArrayD chEffErr[2], TGraphAsymmErrors *effVsSt[3], TGraphAsymmErrors *effVsX[3], Int_t ip, Double_t xp, Bool_t print=kFALSE)
void ComputeTrackingEfficiency(Double_t stEff[6], Double_t stEffErr[6][2], Double_t &spectroEff, Double_t spectroEffErr[2])
void IntegrateMuonEfficiency(TGraphAsymmErrors &effVsRunLow, TGraphAsymmErrors &effVsRunUp, TGraphAsymmErrors &effVsX, Int_t ip, Double_t xp)
void PlotMuonEfficiencyVsRun(TString runList, TString fileNameData, TString fileNameSave, Bool_t print, Bool_t draw)
void PlotMuonEfficiency(TString fileNameData, TString fileNameSave, Bool_t saveEdges, Bool_t print, Bool_t draw)
void PlotIntegratedMuonEfficiency(TString fileNameWeights, TString fileNameSave, Bool_t print, Bool_t draw)
void PlotMuonEfficiencyVsXY(TString xVar, TString yVar, TString fileNameData, TString fileNameSave, Bool_t draw, Bool_t rap=kFALSE)
Bool_t moreTrackCandidates
void LoadRunWeights(TString fileName)
void GetDEEfficiency(THnSparse &TT, THnSparse &TD, TGraphAsymmErrors &effVsDE)
TString gObjNameExtension
void SetRunLabel(TGraphAsymmErrors &g, Int_t irun, const TList &runs)
void PlotIntegratedMuonEfficiencyVsX(TString var, TString runList, TString fileNameWeights, TString fileNameData, TString fileNameSave, Bool_t print, Bool_t draw)
void BeautifyGraphs(TObjArray &array, const char *xAxisName, const char *yAxisName)
TGraphAsymmErrors * CreateGraph(const char *name, const char *title, int value=-1)
void PlotIntegratedMuonEfficiencyPerDE(TString fileNameWeights, TString fileNameSave)
void ComputeStation45Efficiency(TArrayD &chEff, TArrayD chEffErr[2], Double_t &st45Eff, Double_t st45EffErr[2])
Bool_t Data(TH1F *h, Double_t *rangefit, Bool_t writefit, Double_t &sgn, Double_t &errsgn, Double_t &bkg, Double_t &errbkg, Double_t &sgnf, Double_t &errsgnf, Double_t &sigmafit, Int_t &status)
const Char_t * effErrMode
void BeautifyGraph(TGraphAsymmErrors &g, const char *xAxisName, const char *yAxisName)
void PlotMuonEfficiencyVsX(TString var, TString fileNameData, TString fileNameSave, Bool_t saveEdges, Bool_t print, Bool_t draw)
void PlotMuonEfficiencyPerDE(TString fileNameData, TString fileNameSave, Bool_t saveEdges)
void PlotMuonEfficiencyPerDEVsRun(TString runList, TString fileNameData, TString fileNameSave)
void GetStation45Efficiency(TArrayD &chEff, TArrayD chEffErr[2], TGraphAsymmErrors *effVsX[3], Int_t ip, Double_t xp)
Bool_t GetChamberEfficiency(THnSparse &TT, THnSparse &TD, TArrayD &chEff, TArrayD chEffErr[2], Bool_t printError=kFALSE)
void SetCentPtCh(THnSparse &SparseData)