/[PAMELA software]/PamCut/Collections/EffRigCollection/EffRigCollection.cpp
ViewVC logotype

Diff of /PamCut/Collections/EffRigCollection/EffRigCollection.cpp

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 1.2 by pam-fi, Thu Sep 24 18:17:26 2009 UTC revision 1.4 by pam-fi, Thu Oct 1 10:33:49 2009 UTC
# Line 9  Line 9 
9    
10  #include "EffRigCollection.h"  #include "EffRigCollection.h"
11    
12    #include "TGraphAsymmErrors.h"
13    
14  extern "C" {  extern "C" {
15  bool efficiency_(Int_t*, Int_t*, Double_t*, Double_t*, Double_t*);  bool efficiency_(Int_t*, Int_t*, Double_t*, Double_t*, Double_t*);
16  }  }
17    
18  EffRigCollection::EffRigCollection(const char *collectionName, TString outFileBase, TString rigBinsFile, bool absRig) :  EffRigCollection::EffRigCollection(const char *collectionName, TString outFileBase, TString rigBinsFile, int errMethod) :
19    EffCollection(collectionName, outFileBase), _bins(0), _selVector(0), _detVector(0), _outUp(0), _outDown(0) {    EffCollection(collectionName, outFileBase, errMethod), _bins(0), _selVector(0), _detVector(0), _outUp(0), _outDown(0) {
20    
21    ifstream rigBinListFile;    ifstream rigBinListFile;
22    rigBinListFile.open(rigBinsFile);    rigBinListFile.open(rigBinsFile);
# Line 39  int EffRigCollection::ApplyCut(PamLevel2 Line 41  int EffRigCollection::ApplyCut(PamLevel2
41      // Check if the event is inside the rigidity range      // Check if the event is inside the rigidity range
42      // NOTE: at this point a TrkPhSinCut should be already performed,      // NOTE: at this point a TrkPhSinCut should be already performed,
43      //       since we are going to retrieve rigidity.      //       since we are going to retrieve rigidity.
44      float rig;  
45      if (_absRig) {      float rig = event->GetTrack(0)->GetTrkTrack()->GetRigidity();
46        rig = event->GetTrack(0)->GetTrkTrack()->GetRigidity();  
     }  
     else  
       rig = 1. / event->GetTrack(0)->GetTrkTrack()->GetDeflection();  
47      if (rig >= _bins[0]) {      if (rig >= _bins[0]) {
48        int i = 1;        int i = 1;
49        while (rig >= _bins[i] && i < (int) _bins.size()) {        while (rig >= _bins[i] && i < (int) _bins.size()) {
# Line 91  void EffRigCollection::Finalize() { Line 90  void EffRigCollection::Finalize() {
90    for (unsigned int i = 0; i < _selVector.size(); i++) {    for (unsigned int i = 0; i < _selVector.size(); i++) {
91      sel[i] = (Int_t) _selVector[i];      sel[i] = (Int_t) _selVector[i];
92      det[i] = (Int_t) _detVector[i];      det[i] = (Int_t) _detVector[i];
     efficiency_(&(sel[i]), &(det[i]), &(eff[i]), &(errLow[i]), &(errHigh[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 == EFFERR_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          for (int j = 0; j < det[i]; j++)
110            pass.Fill((binning[i + 1] + binning[i]) / 2.);
111          for (int j = 0; j < sel[i]; j++)
112            total.Fill((binning[i + 1] + binning[i]) / 2.);
113        }
114    
115        errGraph.BayesDivide(&pass, &total);
116        Double_t graphX;
117        double currBin;
118        int currPoint = 0;
119    
120        for (unsigned int i = 0; i < _selVector.size(); i++) {
121          errGraph.GetPoint(currPoint, graphX, eff[i]);
122          currBin = (binning[i + 1] + binning[i]) / 2.;
123          if (currBin == graphX) {
124            if (_selVector[i] < 8) {
125              eff[i] = 1.1;
126              errLow[i] = errHigh[i] = 0.;
127              errGraph.SetPoint(currPoint, graphX, 1.1);
128              float halfBin = (binning[i + 1] - binning[i]) / 2.;
129              errGraph.SetPointError(currPoint, halfBin, halfBin, 0., 0.);
130            }
131            currPoint++;
132          }
133        }
134    
135      }
136      if (_errMethod == EFFERR_SERGIO) {
137        for (unsigned int i = 0; i < _selVector.size(); i++) {
138          efficiency_(&(sel[i]), &(det[i]), &(eff[i]), &(errLow[i]), &(errHigh[i]));
139    
140          float R = (_bins[i] + _bins[i + 1]) / 2.;
141          float halfBin = (_bins[i + 1] - _bins[i]) / 2.;
142          errGraph.SetPoint(i, R, eff[i]);
143          errGraph.SetPointError(i, halfBin, halfBin, errLow[i], errHigh[i]);
144    
145        }
146    
147      }
148    
149    // Write the output files    // Write the output files
150    if (_outFileBase != "") {    if (_outFileBase != "") {
151      ofstream outTextFile((_outFileBase + "-eff-rig.txt").Data(), ios_base::out);      ofstream outTextFile((_outFileBase + "-" + GetName() + "-rig.txt").Data(), ios_base::out);
152      streamsize newPrec = 4;      streamsize newPrec = 4;
153      outTextFile.precision(newPrec);      outTextFile.precision(newPrec);
154      outTextFile.setf(ios::fixed, ios::floatfield);      outTextFile.setf(ios::fixed, ios::floatfield);
155      for (unsigned int i = 0; i < _selVector.size(); i++) {      for (unsigned int i = 0; i < _selVector.size(); i++) {
156        outTextFile << setw(10) << _detVector[i] << setw(10) << _selVector[i];        outTextFile << setw(10) << _bins[i] << setw(10) << _bins[i + 1] << setw(10) << _detVector[i] << setw(10)
157              << _selVector[i];
158        if (_selVector[i] != 0)        if (_selVector[i] != 0)
159          outTextFile << setw(10) << eff[i] << setw(10) << errLow[i] << setw(10) << errHigh[i] << "\n";          outTextFile << setw(10) << eff[i] << setw(10) << errLow[i] << setw(10) << errHigh[i] << "\n";
160        else        else
161          outTextFile << setw(10) << 0. << setw(10) << 0. << setw(10) << 0. << endl;          outTextFile << setw(10) << 0. << setw(10) << 0. << setw(10) << 0. << endl;
162    
163      }      }
164    
165      outTextFile.close();      outTextFile.close();
166    
167        TFile outRootFile((_outFileBase + "-" + GetName() + "-rig.root"), "RECREATE");
168        outRootFile.cd();
169        errGraph.Write();
170        outRootFile.Close();
171    
172    }    }
173    
174  }  }
Legend:
Removed from v.1.2  
changed lines
  Added in v.1.4
  ViewVC Help
Powered by ViewVC 1.1.23