| 1 |
/* |
| 2 |
* EffRigCollection.cpp |
| 3 |
* |
| 4 |
* Created on: 10/ago/2009 |
| 5 |
* Author: Nicola Mori |
| 6 |
*/ |
| 7 |
|
| 8 |
/*! @file EffRigCollection.cpp The EffRigCollection class implementation file. */ |
| 9 |
|
| 10 |
#include "EffRigCollection.h" |
| 11 |
|
| 12 |
#include "TGraphAsymmErrors.h" |
| 13 |
|
| 14 |
extern "C" { |
| 15 |
bool efficiency_(Int_t*, Int_t*, Double_t*, Double_t*, Double_t*); |
| 16 |
} |
| 17 |
|
| 18 |
EffRigCollection::EffRigCollection(const char *collectionName, TString outFileBase, TString rigBinsFile, int errMethod) : |
| 19 |
EffCollection(collectionName, outFileBase, errMethod), _bins(0), _selVector(0), _detVector(0), _outUp(0), _outDown(0) { |
| 20 |
|
| 21 |
ifstream rigBinListFile; |
| 22 |
rigBinListFile.open(rigBinsFile); |
| 23 |
|
| 24 |
TString auxString; |
| 25 |
while (!rigBinListFile.eof()) { |
| 26 |
rigBinListFile >> auxString; |
| 27 |
if (auxString != "") { |
| 28 |
_bins.push_back(auxString.Atof()); |
| 29 |
} |
| 30 |
} |
| 31 |
rigBinListFile.close(); |
| 32 |
|
| 33 |
_selVector.resize(_bins.size() - 1, 0); |
| 34 |
_detVector.resize(_bins.size() - 1, 0); |
| 35 |
} |
| 36 |
|
| 37 |
int EffRigCollection::ApplyCut(PamLevel2 *event) { |
| 38 |
|
| 39 |
_nEv++; |
| 40 |
if (_selCollection.ApplyCut(event) == CUTOK) { |
| 41 |
// Check if the event is inside the rigidity range |
| 42 |
// NOTE: at this point a TrkPhSinCut should be already performed, |
| 43 |
// since we are going to retrieve rigidity. |
| 44 |
|
| 45 |
float rig = event->GetTrack(0)->GetTrkTrack()->GetRigidity(); |
| 46 |
|
| 47 |
if (rig >= _bins[0]) { |
| 48 |
int i = 1; |
| 49 |
while (rig >= _bins[i] && i < (int) _bins.size()) { |
| 50 |
i++; |
| 51 |
} |
| 52 |
i--; |
| 53 |
|
| 54 |
if (i < (int) (_selVector.size())) { |
| 55 |
_selVector[i]++; |
| 56 |
_sel++; |
| 57 |
if (_detCollection.ApplyCut(event) == CUTOK) { |
| 58 |
_detVector[i]++; |
| 59 |
_det++; |
| 60 |
_nGood++; |
| 61 |
return CUTOK; |
| 62 |
} |
| 63 |
} |
| 64 |
else |
| 65 |
_outUp++; |
| 66 |
|
| 67 |
} |
| 68 |
else { |
| 69 |
_outDown++; |
| 70 |
return 0; |
| 71 |
} |
| 72 |
} |
| 73 |
|
| 74 |
return 0; |
| 75 |
} |
| 76 |
|
| 77 |
void EffRigCollection::Finalize() { |
| 78 |
|
| 79 |
// Print the report |
| 80 |
EffCollection::Finalize(); |
| 81 |
cout << " Events below the minimum rigidity: " << _outDown << "\n"; |
| 82 |
cout << " Events above the maximum rigidity: " << _outUp << "\n"; |
| 83 |
|
| 84 |
// Compute the error |
| 85 |
Int_t sel[_selVector.size()]; |
| 86 |
Int_t det[_detVector.size()]; |
| 87 |
Double_t eff[_selVector.size()]; |
| 88 |
Double_t errLow[_selVector.size()]; |
| 89 |
Double_t errHigh[_selVector.size()]; |
| 90 |
for (unsigned int i = 0; i < _selVector.size(); i++) { |
| 91 |
sel[i] = (Int_t) _selVector[i]; |
| 92 |
det[i] = (Int_t) _detVector[i]; |
| 93 |
} |
| 94 |
|
| 95 |
TGraphAsymmErrors errGraph; |
| 96 |
errGraph.SetName(GetName()); |
| 97 |
errGraph.SetMarkerColor(kRed); |
| 98 |
errGraph.SetMarkerStyle(7); |
| 99 |
errGraph.GetYaxis()->SetRangeUser(0, 1.2); |
| 100 |
|
| 101 |
if (_errMethod == EFFRIG_ROOT) { |
| 102 |
double binning[_bins.size()]; |
| 103 |
for (unsigned int i = 0; i < _bins.size(); i++) |
| 104 |
binning[i] = _bins[i]; |
| 105 |
|
| 106 |
TH1I pass("pass", "pass", _bins.size() - 1, binning); |
| 107 |
TH1I total("total", "total", _bins.size() - 1, binning); |
| 108 |
for (unsigned int i = 0; i < _selVector.size(); i++) { |
| 109 |
cout << (binning[i + 1] + binning[i]) / 2. << " " << det[i] << " " << sel[i] << endl; |
| 110 |
for (int j = 0; j < det[i]; j++) |
| 111 |
pass.Fill((binning[i + 1] + binning[i]) / 2.); |
| 112 |
for (int j = 0; j < sel[i]; j++) |
| 113 |
total.Fill((binning[i + 1] + binning[i]) / 2.); |
| 114 |
} |
| 115 |
|
| 116 |
errGraph.BayesDivide(&pass, &total); |
| 117 |
Double_t graphX; |
| 118 |
double currBin; |
| 119 |
int currPoint = 0; |
| 120 |
//cout << errGraph.GetN() << " " << pass.GetNbinsX() << " " << total.GetNbinsX() << endl; |
| 121 |
for (unsigned int i = 0; i < _selVector.size(); i++) { |
| 122 |
errGraph.GetPoint(currPoint, graphX, eff[i]); |
| 123 |
currBin = (binning[i + 1] + binning[i]) / 2.; |
| 124 |
if (currBin == graphX) { |
| 125 |
if (_selVector[i] < 8) { |
| 126 |
eff[i] = 1.1; |
| 127 |
errLow[i] = errHigh[i] = 0.; |
| 128 |
errGraph.SetPoint(currPoint, graphX, 1.1); |
| 129 |
float halfBin = (binning[i + 1] - binning[i]) / 2.; |
| 130 |
errGraph.SetPointError(currPoint, halfBin, halfBin, 0., 0.); |
| 131 |
} |
| 132 |
currPoint++; |
| 133 |
} |
| 134 |
} |
| 135 |
|
| 136 |
} |
| 137 |
if (_errMethod == EFFRIG_SERGIO) { |
| 138 |
for (unsigned int i = 0; i < _selVector.size(); i++) { |
| 139 |
efficiency_(&(sel[i]), &(det[i]), &(eff[i]), &(errLow[i]), &(errHigh[i])); |
| 140 |
} |
| 141 |
|
| 142 |
} |
| 143 |
|
| 144 |
// Write the output files |
| 145 |
if (_outFileBase != "") { |
| 146 |
ofstream outTextFile((_outFileBase + "-eff-rig.txt").Data(), ios_base::out); |
| 147 |
streamsize newPrec = 4; |
| 148 |
outTextFile.precision(newPrec); |
| 149 |
outTextFile.setf(ios::fixed, ios::floatfield); |
| 150 |
for (unsigned int i = 0; i < _selVector.size(); i++) { |
| 151 |
outTextFile << setw(10) << _detVector[i] << setw(10) << _selVector[i]; |
| 152 |
if (_selVector[i] != 0) |
| 153 |
outTextFile << setw(10) << eff[i] << setw(10) << errLow[i] << setw(10) << errHigh[i] << "\n"; |
| 154 |
else |
| 155 |
outTextFile << setw(10) << 0. << setw(10) << 0. << setw(10) << 0. << endl; |
| 156 |
|
| 157 |
} |
| 158 |
|
| 159 |
outTextFile.close(); |
| 160 |
if (_errMethod == EFFRIG_ROOT) { |
| 161 |
TFile outRootFile((_outFileBase + "-eff-rig.root"), "RECREATE"); |
| 162 |
outRootFile.cd(); |
| 163 |
errGraph.Write(); |
| 164 |
outRootFile.Close(); |
| 165 |
} |
| 166 |
} |
| 167 |
|
| 168 |
} |