CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
List of all members | Public Member Functions | Private Member Functions | Private Attributes
ErrorsPropagationAnalyzer Class Reference

#include <MuonAnalysis/MomentumScaleCalibration/plugins/ErrorsPropagationAnalyzer.cc>

Inheritance diagram for ErrorsPropagationAnalyzer:
edm::EDAnalyzer edm::EDConsumerBase

Public Member Functions

 ErrorsPropagationAnalyzer (const edm::ParameterSet &)
 
 ~ErrorsPropagationAnalyzer ()
 
- Public Member Functions inherited from edm::EDAnalyzer
void callWhenNewProductsRegistered (std::function< void(BranchDescription const &)> const &func)
 
 EDAnalyzer ()
 
ModuleDescription const & moduleDescription () const
 
std::string workerType () const
 
virtual ~EDAnalyzer ()
 
- Public Member Functions inherited from edm::EDConsumerBase
 EDConsumerBase ()
 
ProductHolderIndexAndSkipBit indexFrom (EDGetToken, BranchType, TypeID const &) const
 
void itemsMayGet (BranchType, std::vector< ProductHolderIndexAndSkipBit > &) const
 
void itemsToGet (BranchType, std::vector< ProductHolderIndexAndSkipBit > &) const
 
std::vector
< ProductHolderIndexAndSkipBit >
const & 
itemsToGetFromEvent () const
 
void labelsForToken (EDGetToken iToken, Labels &oLabels) const
 
bool registeredToConsume (ProductHolderIndex, bool, BranchType) const
 
bool registeredToConsumeMany (TypeID const &, BranchType) const
 
void updateLookup (BranchType iBranchType, ProductHolderIndexHelper const &)
 
virtual ~EDConsumerBase ()
 

Private Member Functions

virtual void analyze (const edm::Event &, const edm::EventSetup &)
 
void drawHistograms (const TProfile *histo, const TProfile *histoPlusErr, const TProfile *histoMinusErr, const TString &type, const TString &yLabel)
 
virtual void endJob ()
 
void fillHistograms ()
 
void fillValueError ()
 
double massResolution (const lorentzVector &mu1, const lorentzVector &mu2, const std::vector< double > &parval, const double &sigmaPt1, const double &sigmaPt2)
 Modified method to take into account the error. More...
 
double massResolution (const lorentzVector &mu1, const lorentzVector &mu2, double *parval, const double &sigmaPt1, const double &sigmaPt2)
 

Private Attributes

bool debug_
 
std::vector< int > errorFactors_
 
std::vector< double > errors_
 
int etaBins_
 
double etaMax_
 
double etaMaxCut_
 
double etaMin_
 
double etaMinCut_
 
uint32_t maxEvents_
 
TString outputFileName_
 
std::vector< double > parameters_
 
int ptBins_
 
double ptMax_
 
double ptMaxCut_
 
double ptMin_
 
double ptMinCut_
 
int resolFitType_
 
TProfile * sigmaMassOverMassVsEta_
 
TProfile * sigmaMassOverMassVsEtaMinusErr_
 
TProfile * sigmaMassOverMassVsEtaPlusErr_
 
TProfile * sigmaMassOverMassVsPt_
 
TProfile * sigmaMassOverMassVsPtMinusErr_
 
TProfile * sigmaMassOverMassVsPtPlusErr_
 
TProfile * sigmaMassVsEta_
 
TProfile * sigmaMassVsEtaMinusErr_
 
TProfile * sigmaMassVsEtaPlusErr_
 
TProfile * sigmaMassVsPt_
 
TProfile * sigmaMassVsPtMinusErr_
 
TProfile * sigmaMassVsPtPlusErr_
 
TProfile * sigmaPtVsEta_
 
TProfile * sigmaPtVsEtaDiff_
 
TProfile * sigmaPtVsEtaMinusErr_
 
TProfile * sigmaPtVsEtaPlusErr_
 
TProfile * sigmaPtVsPt_
 
TProfile * sigmaPtVsPtDiff_
 
TProfile * sigmaPtVsPtMinusErr_
 
TProfile * sigmaPtVsPtPlusErr_
 
TString treeFileName_
 
std::vector< double > valueMinusError_
 
std::vector< double > valuePlusError_
 

Additional Inherited Members

- Public Types inherited from edm::EDAnalyzer
typedef EDAnalyzer ModuleType
 
- Static Public Member Functions inherited from edm::EDAnalyzer
static const std::string & baseType ()
 
static void fillDescriptions (ConfigurationDescriptions &descriptions)
 
static void prevalidate (ConfigurationDescriptions &)
 
- Protected Member Functions inherited from edm::EDConsumerBase
template<typename ProductType , BranchType B = InEvent>
EDGetTokenT< ProductType > consumes (edm::InputTag const &tag)
 
EDGetToken consumes (const TypeToGet &id, edm::InputTag const &tag)
 
template<BranchType B>
EDGetToken consumes (TypeToGet const &id, edm::InputTag const &tag)
 
ConsumesCollector consumesCollector ()
 Use a ConsumesCollector to gather consumes information from helper functions. More...
 
template<typename ProductType , BranchType B = InEvent>
void consumesMany ()
 
void consumesMany (const TypeToGet &id)
 
template<BranchType B>
void consumesMany (const TypeToGet &id)
 
template<typename ProductType , BranchType B = InEvent>
EDGetTokenT< ProductType > mayConsume (edm::InputTag const &tag)
 
EDGetToken mayConsume (const TypeToGet &id, edm::InputTag const &tag)
 
template<BranchType B>
EDGetToken mayConsume (const TypeToGet &id, edm::InputTag const &tag)
 

Detailed Description

Description: <one line="" class="" summary>="">

Implementation: <Notes on="" implementation>="">

Definition at line 49 of file ErrorsPropagationAnalyzer.h.

Constructor & Destructor Documentation

ErrorsPropagationAnalyzer::ErrorsPropagationAnalyzer ( const edm::ParameterSet iConfig)
explicit

Definition at line 6 of file ErrorsPropagationAnalyzer.cc.

References gather_cfg::cout, errorFactors_, errors_, etaBins_, etaMax_, etaMin_, cmsRelvalreport::exit, fillValueError(), edm::ParameterSet::getParameter(), parameters_, ptBins_, ptMax_, ptMin_, sigmaMassOverMassVsEta_, sigmaMassOverMassVsEtaMinusErr_, sigmaMassOverMassVsEtaPlusErr_, sigmaMassOverMassVsPt_, sigmaMassOverMassVsPtMinusErr_, sigmaMassOverMassVsPtPlusErr_, sigmaMassVsEta_, sigmaMassVsEtaMinusErr_, sigmaMassVsEtaPlusErr_, sigmaMassVsPt_, sigmaMassVsPtMinusErr_, sigmaMassVsPtPlusErr_, sigmaPtVsEta_, sigmaPtVsEtaDiff_, sigmaPtVsEtaMinusErr_, sigmaPtVsEtaPlusErr_, sigmaPtVsPt_, sigmaPtVsPtDiff_, sigmaPtVsPtMinusErr_, and sigmaPtVsPtPlusErr_.

6  :
7  treeFileName_( iConfig.getParameter<std::string>("InputFileName") ),
8  resolFitType_( iConfig.getParameter<int>("ResolFitType") ),
9  maxEvents_( iConfig.getParameter<int>("MaxEvents") ),
10  outputFileName_( iConfig.getParameter<std::string>("OutputFileName") ),
11  ptBins_( iConfig.getParameter<int>("PtBins") ),
12  ptMin_( iConfig.getParameter<double>("PtMin") ),
13  ptMax_( iConfig.getParameter<double>("PtMax") ),
14  etaBins_( iConfig.getParameter<int>("EtaBins") ),
15  etaMin_( iConfig.getParameter<double>("EtaMin") ),
16  etaMax_( iConfig.getParameter<double>("EtaMax") ),
17  debug_( iConfig.getParameter<bool>("Debug") ),
18  ptMinCut_( iConfig.getUntrackedParameter<double>("PtMinCut", 0.) ),
19  ptMaxCut_( iConfig.getUntrackedParameter<double>("PtMaxCut", 999999.) ),
20  etaMinCut_( iConfig.getUntrackedParameter<double>("EtaMinCut", 0.) ),
21  etaMaxCut_( iConfig.getUntrackedParameter<double>("EtaMaxCut", 100.) )
22 {
23  parameters_ = iConfig.getParameter<std::vector<double> >("Parameters");
24  errors_ = iConfig.getParameter<std::vector<double> >("Errors");
25  errorFactors_ = iConfig.getParameter<std::vector<int> >("ErrorFactors");
26 
27  if( (parameters_.size() != errors_.size()) || (parameters_.size() != errorFactors_.size()) ) {
28  std::cout << "Error: parameters and errors have different number of values" << std::endl;
29  exit(1);
30  }
31 
33 
34  sigmaPtVsPt_ = new TProfile("sigmaPtVsPtProfile", "sigmaPtVsPt", ptBins_, ptMin_, ptMax_);
35  sigmaPtVsPtPlusErr_ = new TProfile("sigmaPtVsPtPlusErrProfile", "sigmaPtVsPtPlusErr", ptBins_, ptMin_, ptMax_);
36  sigmaPtVsPtMinusErr_ = new TProfile("sigmaPtVsPtMinusErrProfile", "sigmaPtVsPtMinusErr", ptBins_, ptMin_, ptMax_);
37 
38  sigmaPtVsEta_ = new TProfile("sigmaPtVsEtaProfile", "sigmaPtVsEta", etaBins_, etaMin_, etaMax_);
39  sigmaPtVsEtaPlusErr_ = new TProfile("sigmaPtVsEtaPlusErrProfile", "sigmaPtVsEtaPlusErr", etaBins_, etaMin_, etaMax_);
40  sigmaPtVsEtaMinusErr_ = new TProfile("sigmaPtVsEtaMinusErrProfile", "sigmaPtVsEtaMinusErr", etaBins_, etaMin_, etaMax_);
41 
42  sigmaMassVsPt_ = new TProfile("sigmaMassVsPtProfile", "sigmaMassVsPt", ptBins_, ptMin_, ptMax_);
43  sigmaMassVsPtPlusErr_ = new TProfile("sigmaMassVsPtPlusErrProfile", "sigmaMassVsPtPlusErr", ptBins_, ptMin_, ptMax_);
44  sigmaMassVsPtMinusErr_ = new TProfile("sigmaMassVsPtMinusErrProfile", "sigmaMassVsPtMinusErr", ptBins_, ptMin_, ptMax_);
45 
46  sigmaMassVsEta_ = new TProfile("sigmaMassVsEtaProfile", "sigmaMassVsEta", etaBins_, etaMin_, etaMax_);
47  sigmaMassVsEtaPlusErr_ = new TProfile("sigmaMassVsEtaPlusErrProfile", "sigmaMassVsEtaPlusErr", etaBins_, etaMin_, etaMax_);
48  sigmaMassVsEtaMinusErr_ = new TProfile("sigmaMassVsEtaMinusErrProfile", "sigmaMassVsEtaMinusErr", etaBins_, etaMin_, etaMax_);
49 
50  sigmaMassOverMassVsPt_ = new TProfile("sigmaMassOverMassVsPtProfile", "sigmaMassOverMassVsPt", ptBins_, ptMin_, ptMax_);
51  sigmaMassOverMassVsPtPlusErr_ = new TProfile("sigmaMassOverMassVsPtPlusErrProfile", "sigmaMassOverMassVsPtPlusErr", ptBins_, ptMin_, ptMax_);
52  sigmaMassOverMassVsPtMinusErr_ = new TProfile("sigmaMassOverMassVsPtMinusErrProfile", "sigmaMassOverMassVsPtMinusErr", ptBins_, ptMin_, ptMax_);
53 
54  sigmaMassOverMassVsEta_ = new TProfile("sigmaMassOverMassVsEtaProfile", "sigmaMassOverMassVsEta", etaBins_, etaMin_, etaMax_);
55  sigmaMassOverMassVsEtaPlusErr_ = new TProfile("sigmaMassOverMassVsEtaPlusErrProfile", "sigmaMassOverMassVsEtaPlusErr", etaBins_, etaMin_, etaMax_);
56  sigmaMassOverMassVsEtaMinusErr_ = new TProfile("sigmaMassOverMassVsEtaMinusErrProfile", "sigmaMassOverMassVsEtaMinusErr", etaBins_, etaMin_, etaMax_);
57 
58  sigmaPtVsPtDiff_ = new TProfile("sigmaPtVsPtDiffProfile", "sigmaPtVsPtDiff", ptBins_, ptMin_, ptMax_);
59  sigmaPtVsEtaDiff_ = new TProfile("sigmaPtVsEtaDiffProfile", "sigmaPtVsEtaDiff", etaBins_, etaMin_, etaMax_);
60 }
T getParameter(std::string const &) const
T getUntrackedParameter(std::string const &, T const &) const
tuple cout
Definition: gather_cfg.py:121
ErrorsPropagationAnalyzer::~ErrorsPropagationAnalyzer ( )

Definition at line 77 of file ErrorsPropagationAnalyzer.cc.

References drawHistograms(), fillHistograms(), download_sqlite_cfg::outputFile, outputFileName_, sigmaMassOverMassVsEta_, sigmaMassOverMassVsEtaMinusErr_, sigmaMassOverMassVsEtaPlusErr_, sigmaMassOverMassVsPt_, sigmaMassOverMassVsPtMinusErr_, sigmaMassOverMassVsPtPlusErr_, sigmaMassVsEta_, sigmaMassVsEtaMinusErr_, sigmaMassVsEtaPlusErr_, sigmaMassVsPt_, sigmaMassVsPtMinusErr_, sigmaMassVsPtPlusErr_, sigmaPtVsEta_, sigmaPtVsEtaDiff_, sigmaPtVsEtaMinusErr_, sigmaPtVsEtaPlusErr_, sigmaPtVsPt_, sigmaPtVsPtDiff_, sigmaPtVsPtMinusErr_, and sigmaPtVsPtPlusErr_.

78 {
79  gROOT->SetStyle("Plain");
80 
82 
83  TFile * outputFile = new TFile(outputFileName_, "RECREATE");
84 
85  drawHistograms(sigmaPtVsEta_, sigmaPtVsEtaPlusErr_, sigmaPtVsEtaMinusErr_, "sigmaPtVsEta", "#sigma(Pt)/Pt");
86  drawHistograms(sigmaPtVsPt_, sigmaPtVsPtPlusErr_, sigmaPtVsPtMinusErr_, "sigmaPtVsPt", "#sigma(Pt)/Pt");
87 
90 
93 
94  sigmaPtVsPtDiff_->Write();
95  sigmaPtVsEtaDiff_->Write();
96 
97  outputFile->Write();
98  outputFile->Close();
99 }
void drawHistograms(const TProfile *histo, const TProfile *histoPlusErr, const TProfile *histoMinusErr, const TString &type, const TString &yLabel)

Member Function Documentation

void ErrorsPropagationAnalyzer::analyze ( const edm::Event iEvent,
const edm::EventSetup iSetup 
)
privatevirtual

Implements edm::EDAnalyzer.

Definition at line 211 of file ErrorsPropagationAnalyzer.cc.

212 {
213 }
void ErrorsPropagationAnalyzer::drawHistograms ( const TProfile *  histo,
const TProfile *  histoPlusErr,
const TProfile *  histoMinusErr,
const TString &  type,
const TString &  yLabel 
)
private

Definition at line 101 of file ErrorsPropagationAnalyzer.cc.

References svgfig::canvas(), i, sigmaMassOverMassVsEta_, sigmaMassOverMassVsEtaMinusErr_, sigmaMassOverMassVsEtaPlusErr_, sigmaMassVsEta_, sigmaMassVsEtaMinusErr_, sigmaMassVsEtaPlusErr_, sigmaPtVsEta_, sigmaPtVsEtaMinusErr_, sigmaPtVsEtaPlusErr_, indexGen::title, and makeHLTPrescaleTable::values.

Referenced by ~ErrorsPropagationAnalyzer().

104 {
105  TH1D * sigmaPtVsEtaTH1D = new TH1D(type, type, histo->GetNbinsX(),
106  histo->GetXaxis()->GetXmin(), histo->GetXaxis()->GetXmax());
107  TH1D * sigmaPtVsEtaPlusErrTH1D = new TH1D(type+"PlusErr", type+"PlusErr", histo->GetNbinsX(),
108  histo->GetXaxis()->GetXmin(), histo->GetXaxis()->GetXmax());
109  TH1D * sigmaPtVsEtaMinusErrTH1D = new TH1D(type+"MinusErr", type+"MinusErr", histo->GetNbinsX(),
110  histo->GetXaxis()->GetXmin(), histo->GetXaxis()->GetXmax());
111 
112  TH1D * sigmaMassVsEtaTH1D = new TH1D(type, type, histo->GetNbinsX(),
113  histo->GetXaxis()->GetXmin(), histo->GetXaxis()->GetXmax());
114  TH1D * sigmaMassVsEtaPlusErrTH1D = new TH1D(type+"PlusErr", type+"PlusErr", histo->GetNbinsX(),
115  histo->GetXaxis()->GetXmin(), histo->GetXaxis()->GetXmax());
116  TH1D * sigmaMassVsEtaMinusErrTH1D = new TH1D(type+"MinusErr", type+"MinusErr", histo->GetNbinsX(),
117  histo->GetXaxis()->GetXmin(), histo->GetXaxis()->GetXmax());
118 
119  TH1D * sigmaMassOverMassVsEtaTH1D = new TH1D(type, type, histo->GetNbinsX(),
120  histo->GetXaxis()->GetXmin(), histo->GetXaxis()->GetXmax());
121  TH1D * sigmaMassOverMassVsEtaPlusErrTH1D = new TH1D(type+"PlusErr", type+"PlusErr", histo->GetNbinsX(),
122  histo->GetXaxis()->GetXmin(), histo->GetXaxis()->GetXmax());
123  TH1D * sigmaMassOverMassVsEtaMinusErrTH1D = new TH1D(type+"MinusErr", type+"MinusErr", histo->GetNbinsX(),
124  histo->GetXaxis()->GetXmin(), histo->GetXaxis()->GetXmax());
125 
126  TCanvas * canvas = new TCanvas("canvas"+type, "canvas"+type, 1000, 800);
127  for( int iBin = 1; iBin <= histo->GetNbinsX(); ++iBin ) {
128  sigmaPtVsEtaTH1D->SetBinContent(iBin, histo->GetBinContent(iBin));
129  sigmaPtVsEtaPlusErrTH1D->SetBinContent(iBin, histoPlusErr->GetBinContent(iBin));
130  sigmaPtVsEtaMinusErrTH1D->SetBinContent(iBin, histoMinusErr->GetBinContent(iBin));
131  }
132  int numBins = sigmaPtVsEtaTH1D->GetNbinsX();
133  // Draw TGraph with asymmetric errors
134  Double_t * values = sigmaPtVsEtaTH1D->GetArray();
135  Double_t * valuesPlus = sigmaPtVsEtaPlusErrTH1D->GetArray();
136  Double_t * valuesMinus = sigmaPtVsEtaMinusErrTH1D->GetArray();
137  double * posErrors = new double[numBins];
138  double * negErrors = new double[numBins];
139 
140  TGraphAsymmErrors * graphAsymmErrors = new TGraphAsymmErrors(sigmaPtVsEtaTH1D);
141  TGraph * graph = new TGraph(sigmaPtVsEtaTH1D);
142 
143  for( int i=1; i<=numBins; ++i ) {
144  // std::cout << "filling " << i << std::endl;
145  posErrors[i-1] = valuesPlus[i] - values[i];
146  if( valuesMinus[i] < 0 ) negErrors[i-1] = values[i];
147  else negErrors[i-1] = values[i] - valuesMinus[i];
148 
149  graphAsymmErrors->SetPointEYlow(i-1, negErrors[i-1]);
150  graphAsymmErrors->SetPointEYhigh(i-1, posErrors[i-1]);
151  }
152 
153  canvas->Draw();
154 
155  graphAsymmErrors->SetTitle("");
156  graphAsymmErrors->SetFillColor(kGray);
157  graphAsymmErrors->Draw("A2");
158  TString title(type);
159  if( type == "Eta" ) title = "#eta";
160  graphAsymmErrors->GetXaxis()->SetTitle(title);
161  graphAsymmErrors->GetYaxis()->SetTitle("yLabel");
162  graph->Draw();
163 
164  // graph->SetLineColor(kGray);
165  // graph->SetMarkerColor(kBlack);
166  // graph->Draw("AP");
167 
168 // if( debug_ ) {
169 // sigmaPtVsEtaPlusErrTH1D->SetLineColor(kRed);
170 // sigmaPtVsEtaMinusErrTH1D->SetLineColor(kBlue);
171 // }
172 // else {
173 // sigmaPtVsEtaPlusErrTH1D->SetFillColor(kGray);
174 // sigmaPtVsEtaPlusErrTH1D->SetLineColor(kWhite);
175 // sigmaPtVsEtaMinusErrTH1D->SetFillColor(kWhite);
176 // sigmaPtVsEtaMinusErrTH1D->SetLineColor(kWhite);
177 // }
178 // sigmaPtVsEtaPlusErrTH1D->Draw();
179 // sigmaPtVsEtaTH1D->Draw("SAME");
180 // sigmaPtVsEtaMinusErrTH1D->Draw("SAME");
181 
182  sigmaPtVsEtaPlusErrTH1D->Write();
183  sigmaPtVsEtaTH1D->Write();
184  sigmaPtVsEtaMinusErrTH1D->Write();
185 
186  sigmaPtVsEtaPlusErr_->Write();
187  sigmaPtVsEta_->Write();
188  sigmaPtVsEtaMinusErr_->Write();
189 
190  // Mass
191  sigmaMassVsEtaPlusErrTH1D->Write();
192  sigmaMassVsEtaTH1D->Write();
193  sigmaMassVsEtaMinusErrTH1D->Write();
194 
195  sigmaMassOverMassVsEtaPlusErrTH1D->Write();
196  sigmaMassOverMassVsEtaTH1D->Write();
197  sigmaMassOverMassVsEtaMinusErrTH1D->Write();
198 
199  sigmaMassVsEtaPlusErr_->Write();
200  sigmaMassVsEta_->Write();
201  sigmaMassVsEtaMinusErr_->Write();
202 
204  sigmaMassOverMassVsEta_->Write();
206 
207  canvas->Write();
208 }
type
Definition: HCALResponse.h:21
int i
Definition: DBlmapReader.cc:9
def canvas
Definition: svgfig.py:481
Definition: adjgraph.h:12
virtual void ErrorsPropagationAnalyzer::endJob ( void  )
inlineprivatevirtual

Reimplemented from edm::EDAnalyzer.

Definition at line 61 of file ErrorsPropagationAnalyzer.h.

61 {};
void ErrorsPropagationAnalyzer::fillHistograms ( void  )
private

Definition at line 286 of file ErrorsPropagationAnalyzer.cc.

References gather_cfg::cout, MuScleFitUtils::debug, debug_, MuScleFitUtils::debugMassResol_, errors_, etaMaxCut_, etaMinCut_, i, massResolution(), maxEvents_, parameters_, ptMaxCut_, ptMinCut_, RootTreeHandler::readTree(), MuScleFitUtils::resfind, resolFitType_, MuScleFitUtils::resolutionFunction, resolutionFunctionService(), resolutionFunctionVecService(), sigmaMassOverMassVsEta_, sigmaMassOverMassVsEtaMinusErr_, sigmaMassOverMassVsEtaPlusErr_, sigmaMassOverMassVsPt_, sigmaMassOverMassVsPtMinusErr_, sigmaMassOverMassVsPtPlusErr_, sigmaMassVsEta_, sigmaMassVsEtaMinusErr_, sigmaMassVsEtaPlusErr_, sigmaMassVsPt_, sigmaMassVsPtMinusErr_, sigmaMassVsPtPlusErr_, resolutionFunctionBase< T >::sigmaPt(), resolutionFunctionBase< T >::sigmaPtError(), sigmaPtVsEta_, sigmaPtVsEtaDiff_, sigmaPtVsEtaMinusErr_, sigmaPtVsEtaPlusErr_, sigmaPtVsPt_, sigmaPtVsPtDiff_, sigmaPtVsPtMinusErr_, sigmaPtVsPtPlusErr_, SigmaPtDiff::squaredDiff(), and treeFileName_.

Referenced by ~ErrorsPropagationAnalyzer().

287 {
288  std::cout << "Reading muon pairs from Root Tree in " << treeFileName_ << std::endl;
289  RootTreeHandler rootTreeHandler;
290 
291  typedef std::vector<std::pair<lorentzVector,lorentzVector> > MuonPairVector;
292  MuonPairVector savedPair;
293  std::vector<std::pair<int, int> > evtRun;
294  rootTreeHandler.readTree(maxEvents_, treeFileName_, &savedPair, 0, &evtRun);
295  // rootTreeHandler.readTree(maxEvents, inputRootTreeFileName_, &savedPair, &(MuScleFitUtils::genPair));
296 
301  MuScleFitUtils::resfind = std::vector<int>(6, 0);
302 
303  SigmaPt sigmaPt( parameters_, errors_ );
304  SigmaPtDiff sigmaPtDiff;
305 
306  // Loop on all the pairs
307  unsigned int i = 0;
308  MuonPairVector::iterator it = savedPair.begin();
309  std::cout << "Starting loop on " << savedPair.size() << " muons" << std::endl;
310  for( ; it != savedPair.end(); ++it, ++i ) {
311  double pt1 = it->first.pt();
312  double eta1 = it->first.eta();
313  double pt2 = it->second.pt();
314  double eta2 = it->second.eta();
315 
316  if( debug_ ) {
317  std::cout << "pt1 = " << pt1 << ", eta1 = " << eta1 << ", pt2 = " << pt2 << ", eta2 = " << eta2 << std::endl;
318  }
319  // double fabsEta1 = fabs(eta1);
320  // double fabsEta2 = fabs(eta2);
321 
322  if( pt1 == 0 && pt2 == 0 && eta1 == 0 && eta2 == 0 ) continue;
323 
324 
325  // double sigmaPt1 = sigmaPt( eta1 );
326  // double sigmaPt2 = sigmaPt( eta2 );
327 
328  // double sigmaPtPlusErr1 = sigmaPt1 + sigmaPt.sigma(eta1);
329  // double sigmaPtPlusErr2 = sigmaPt2 + sigmaPt.sigma(eta2);
330  // double sigmaPtMinusErr1 = sigmaPt1 - sigmaPt.sigma(eta1);
331  // double sigmaPtMinusErr2 = sigmaPt2 - sigmaPt.sigma(eta2);
332 
333  double sigmaPt1 = resolutionFunctionForVec->sigmaPt( pt1,eta1,parameters_ );
334  double sigmaPt2 = resolutionFunctionForVec->sigmaPt( pt2,eta2,parameters_ );
335  double sigmaPtPlusErr1 = sigmaPt1 + resolutionFunctionForVec->sigmaPtError( pt1,eta1,parameters_,errors_);
336  double sigmaPtPlusErr2 = sigmaPt2 + resolutionFunctionForVec->sigmaPtError( pt2,eta2,parameters_,errors_ );
337  double sigmaPtMinusErr1 = sigmaPt1 - resolutionFunctionForVec->sigmaPtError( pt1,eta1,parameters_,errors_ );
338  double sigmaPtMinusErr2 = sigmaPt2 - resolutionFunctionForVec->sigmaPtError( pt2,eta2,parameters_,errors_ );
339 
340  // double sigmaMass = MuScleFitUtils::massResolution( it->first, it->second, parameters_ );
341  // double sigmaMassPlusErr = MuScleFitUtils::massResolution( it->first, it->second, valuePlusError_ );
342  // double sigmaMassMinusErr = MuScleFitUtils::massResolution( it->first, it->second, valueMinusError_ );
343  double sigmaMass = massResolution( it->first, it->second, parameters_, sigmaPt1, sigmaPt2 );
344  double sigmaMassPlusErr = massResolution( it->first, it->second, parameters_, sigmaPtPlusErr1, sigmaPtPlusErr2 );
345  double sigmaMassMinusErr = massResolution( it->first, it->second, parameters_, sigmaPtMinusErr1, sigmaPtMinusErr2 );
346 
347  double mass = (it->first + it->second).mass();
348 
349  if( debug_ ) {
350  std::cout << "sigmaPt1 = " << sigmaPt1 << " + " << sigmaPtPlusErr1 << " - " << sigmaPtMinusErr1 << std::endl;
351  std::cout << "sigmaPt2 = " << sigmaPt2 << " + " << sigmaPtPlusErr2 << " - " << sigmaPtMinusErr2 << std::endl;
352  std::cout << "sigmaMass = " << sigmaMass << " + " << sigmaMassPlusErr << " - " << sigmaMassMinusErr << std::endl;
353  }
354 
355  // Protections from nans
356  if( pt1 != pt1 ) continue;
357  if( pt2 != pt2 ) continue;
358  if( eta1 != eta1 ) continue;
359  if( eta2 != eta2 ) continue;
360  if( sigmaPt1 != sigmaPt1 ) continue;
361  if( sigmaPt2 != sigmaPt2 ) continue;
362  if( sigmaPtPlusErr1 != sigmaPtPlusErr1 ) continue;
363  if( sigmaPtPlusErr2 != sigmaPtPlusErr2 ) continue;
364  if( sigmaPtMinusErr1 != sigmaPtMinusErr1 ) continue;
365  if( sigmaPtMinusErr2 != sigmaPtMinusErr2 ) continue;
366  if( sigmaMass != sigmaMass ) continue;
367  if( sigmaMassPlusErr != sigmaMassPlusErr ) continue;
368  if( sigmaMassMinusErr != sigmaMassMinusErr ) continue;
369  if( mass == 0 ) continue;
370 
371  std::cout << "Muon pair number " << i << std::endl;
372 
373  // std::cout << "sigmaPt1 = " << sigmaPt1 << ", sigmaPt2 = " << sigmaPt2 << std::endl;
374  // std::cout << "sigmaPtPlusErr1 = " << sigmaPtPlusErr1 << ", sigmaPtMinusErr2 = " << sigmaPtMinusErr2 << std::endl;
375  // std::cout << "sigmaPtMinusErr1 = " << sigmaPtPlusErr1 << ", sigmaPtMinusErr2 = " << sigmaPtMinusErr2 << std::endl;
376 
377 
378  if( (pt1 >= ptMinCut_ && pt1 <= ptMaxCut_) && (fabs(eta1) >= etaMinCut_ && fabs(eta1) <= etaMaxCut_) ) {
379  sigmaPtVsPt_->Fill(pt1, sigmaPt1);
380  sigmaPtVsPtPlusErr_->Fill(pt1, sigmaPtPlusErr1);
381  sigmaPtVsPtMinusErr_->Fill(pt1, sigmaPtMinusErr1);
382  sigmaPtVsPtDiff_->Fill(pt1, sigmaPtDiff.squaredDiff(eta1));
383  sigmaMassVsPt_->Fill(pt1, sigmaMass*mass);
384  sigmaMassVsPtPlusErr_->Fill(pt1, sigmaMassPlusErr*mass);
385  sigmaMassVsPtMinusErr_->Fill(pt1, sigmaMassMinusErr*mass);
386  sigmaMassOverMassVsPt_->Fill(pt1, sigmaMass);
387  sigmaMassOverMassVsPtPlusErr_->Fill(pt1, sigmaMassPlusErr);
388  sigmaMassOverMassVsPtMinusErr_->Fill(pt1, sigmaMassMinusErr);
389 
390  sigmaPtVsEta_->Fill(eta1, sigmaPt1);
391  sigmaPtVsEtaPlusErr_->Fill(eta1, sigmaPtPlusErr1);
392  sigmaPtVsEtaMinusErr_->Fill(eta1, sigmaPtMinusErr1);
393  sigmaPtVsEtaDiff_->Fill(eta1, sigmaPtDiff.squaredDiff(eta1));
394  sigmaMassVsEta_->Fill(eta1, sigmaMass*mass);
395  sigmaMassVsEtaPlusErr_->Fill(eta1, sigmaMassPlusErr*mass);
396  sigmaMassVsEtaMinusErr_->Fill(eta1, sigmaMassMinusErr*mass);
397  sigmaMassOverMassVsEta_->Fill(eta1, sigmaMass);
398  sigmaMassOverMassVsEtaPlusErr_->Fill(eta1, sigmaMassPlusErr);
399  sigmaMassOverMassVsEtaMinusErr_->Fill(eta1, sigmaMassMinusErr);
400  }
401  if( (pt2 >= ptMinCut_ && pt2 <= ptMaxCut_) && (fabs(eta2) >= etaMinCut_ && fabs(eta2) <= etaMaxCut_) ) {
402  sigmaPtVsPt_->Fill(pt2, sigmaPt2);
403  sigmaPtVsPtPlusErr_->Fill(pt2, sigmaPtPlusErr2);
404  sigmaPtVsPtMinusErr_->Fill(pt2, sigmaPtMinusErr2);
405  sigmaPtVsPtDiff_->Fill(pt2, sigmaPtDiff.squaredDiff(eta2));
406  sigmaMassVsPt_->Fill(pt2, sigmaMass*mass);
407  sigmaMassVsPtPlusErr_->Fill(pt2, sigmaMassPlusErr*mass);
408  sigmaMassVsPtMinusErr_->Fill(pt2, sigmaMassMinusErr*mass);
409  sigmaMassOverMassVsPt_->Fill(pt2, sigmaMass);
410  sigmaMassOverMassVsPtPlusErr_->Fill(pt2, sigmaMassPlusErr);
411  sigmaMassOverMassVsPtMinusErr_->Fill(pt2, sigmaMassMinusErr);
412 
413  sigmaPtVsEta_->Fill(eta2, sigmaPt2);
414  sigmaPtVsEtaPlusErr_->Fill(eta2, sigmaPtPlusErr2);
415  sigmaPtVsEtaMinusErr_->Fill(eta2, sigmaPtMinusErr2);
416  sigmaPtVsEtaDiff_->Fill(eta2, sigmaPtDiff.squaredDiff(eta2));
417  sigmaMassVsEta_->Fill(eta2, sigmaMass*mass);
418  sigmaMassVsEtaPlusErr_->Fill(eta2, sigmaMassPlusErr*mass);
419  sigmaMassVsEtaMinusErr_->Fill(eta2, sigmaMassMinusErr*mass);
420  sigmaMassOverMassVsEta_->Fill(eta2, sigmaMass);
421  sigmaMassOverMassVsEtaPlusErr_->Fill(eta2, sigmaMassPlusErr);
422  sigmaMassOverMassVsEtaMinusErr_->Fill(eta2, sigmaMassMinusErr);
423  }
424  }
425 }
void readTree(const int maxEvents, const TString &fileName, MuonPairVector *savedPair, const int muonType, std::vector< std::pair< int, int > > *evtRun, MuonPairVector *genPair=0)
int i
Definition: DBlmapReader.cc:9
Returns ( sigmaPt/Pt(data)^2 - sigmaPt/Pt(MC)^2 )
Definition: SigmaPtDiff.h:68
static bool debugMassResol_
static int debug
virtual double sigmaPt(const double &pt, const double &eta, const T &parval)=0
resolutionFunctionBase< std::vector< double > > * resolutionFunctionVecService(const int identifier)
Service to build the resolution functor corresponding to the passed identifier when receiving a std::...
Definition: Functions.cc:156
std::vector< std::pair< lorentzVector, lorentzVector > > MuonPairVector
double squaredDiff(const double &eta)
Definition: SigmaPtDiff.h:136
double massResolution(const lorentzVector &mu1, const lorentzVector &mu2, const std::vector< double > &parval, const double &sigmaPt1, const double &sigmaPt2)
Modified method to take into account the error.
resolutionFunctionBase< double * > * resolutionFunctionService(const int identifier)
Service to build the resolution functor corresponding to the passed identifier.
Definition: Functions.cc:116
tuple cout
Definition: gather_cfg.py:121
static std::vector< int > resfind
virtual double sigmaPtError(const double &pt, const double &eta, const T &parval, const T &parError)
Definition: Functions.h:3025
static resolutionFunctionBase< double * > * resolutionFunction
void ErrorsPropagationAnalyzer::fillValueError ( )
private

Definition at line 62 of file ErrorsPropagationAnalyzer.cc.

References errorFactors_, errors_, i, parameters_, valueMinusError_, and valuePlusError_.

Referenced by ErrorsPropagationAnalyzer().

63 {
64  valuePlusError_.resize(parameters_.size());
65  valueMinusError_.resize(parameters_.size());
66 
67  std::vector<double>::const_iterator parIt = parameters_.begin();
68  std::vector<double>::const_iterator errIt = errors_.begin();
69  std::vector<int>::const_iterator errFactorIt = errorFactors_.begin();
70  int i=0;
71  for( ; parIt != parameters_.end(); ++parIt, ++errIt, ++errFactorIt, ++i ) {
72  valuePlusError_[i] = *parIt + (*errIt)*(*errFactorIt);
73  valueMinusError_[i] = *parIt - (*errIt)*(*errFactorIt);
74  }
75 }
int i
Definition: DBlmapReader.cc:9
std::vector< double > valueMinusError_
std::vector< double > valuePlusError_
double ErrorsPropagationAnalyzer::massResolution ( const lorentzVector mu1,
const lorentzVector mu2,
const std::vector< double > &  parval,
const double &  sigmaPt1,
const double &  sigmaPt2 
)
private

Modified method to take into account the error.

Definition at line 215 of file ErrorsPropagationAnalyzer.cc.

References AlCaHLTBitMon_ParallelJobs::p.

Referenced by fillHistograms().

220 {
221  double * p = new double[(int)(parval.size())];
222  std::vector<double>::const_iterator it = parval.begin();
223  int id = 0;
224  for ( ; it!=parval.end(); ++it, ++id) {
225  p[id] = *it;
226  }
227  double massRes = massResolution (mu1, mu2, p, sigmaPt1, sigmaPt2);
228  delete[] p;
229  return massRes;
230 }
double massResolution(const lorentzVector &mu1, const lorentzVector &mu2, const std::vector< double > &parval, const double &sigmaPt1, const double &sigmaPt2)
Modified method to take into account the error.
double ErrorsPropagationAnalyzer::massResolution ( const lorentzVector mu1,
const lorentzVector mu2,
double *  parval,
const double &  sigmaPt1,
const double &  sigmaPt2 
)
private

Definition at line 232 of file ErrorsPropagationAnalyzer.cc.

References funct::cos(), resolutionFunctionBase< T >::covPt1Pt2(), create_public_lumi_plots::exp, MuScleFitUtils::mMu2, funct::pow(), MuScleFitUtils::resolutionFunction, resolutionFunctionBase< T >::sigmaCotgTh(), resolutionFunctionBase< T >::sigmaPhi(), funct::sin(), and mathSSE::sqrt().

237 {
238  double mass = (mu1+mu2).mass();
239  double pt1 = mu1.Pt();
240  double phi1 = mu1.Phi();
241  double eta1 = mu1.Eta();
242  double theta1 = 2*atan(exp(-eta1));
243  double pt2 = mu2.Pt();
244  double phi2 = mu2.Phi();
245  double eta2 = mu2.Eta();
246  double theta2 = 2*atan(exp(-eta2));
247 
248  // ATTENTION: need to compute 1/tan(theta) as cos(theta)/sin(theta) because the latter diverges for theta=pi/2
249  // -----------------------------------------------------------------------------------------------------------
250  double dmdpt1 = (pt1/std::pow(sin(theta1),2)*sqrt((std::pow(pt2/sin(theta2),2)+MuScleFitUtils::mMu2)/(std::pow(pt1/sin(theta1),2)+MuScleFitUtils::mMu2))-
251  pt2*(cos(phi1-phi2)+cos(theta1)*cos(theta2)/(sin(theta1)*sin(theta2))))/mass;
252  double dmdpt2 = (pt2/std::pow(sin(theta2),2)*sqrt((std::pow(pt1/sin(theta1),2)+MuScleFitUtils::mMu2)/(std::pow(pt2/sin(theta2),2)+MuScleFitUtils::mMu2))-
253  pt1*(cos(phi2-phi1)+cos(theta2)*cos(theta1)/(sin(theta2)*sin(theta1))))/mass;
254  double dmdphi1 = pt1*pt2/mass*sin(phi1-phi2);
255  double dmdphi2 = pt2*pt1/mass*sin(phi2-phi1);
256  double dmdcotgth1 = (pt1*pt1*cos(theta1)/sin(theta1)*
257  sqrt((std::pow(pt2/sin(theta2),2)+MuScleFitUtils::mMu2)/(std::pow(pt1/sin(theta1),2)+MuScleFitUtils::mMu2)) -
258  pt1*pt2*cos(theta2)/sin(theta2))/mass;
259  double dmdcotgth2 = (pt2*pt2*cos(theta2)/sin(theta2)*
260  sqrt((std::pow(pt1/sin(theta1),2)+MuScleFitUtils::mMu2)/(std::pow(pt2/sin(theta2),2)+MuScleFitUtils::mMu2)) -
261  pt2*pt1*cos(theta1)/sin(theta1))/mass;
262 
263  // Resolution parameters:
264  // ----------------------
265  // double sigma_pt1 = MuScleFitUtils::resolutionFunction->sigmaPt( pt1,eta1,parval );
266  // double sigma_pt2 = MuScleFitUtils::resolutionFunction->sigmaPt( pt2,eta2,parval );
267  double sigma_pt1 = sigmaPt1;
268  double sigma_pt2 = sigmaPt2;
269  double sigma_phi1 = MuScleFitUtils::resolutionFunction->sigmaPhi( pt1,eta1,parval );
270  double sigma_phi2 = MuScleFitUtils::resolutionFunction->sigmaPhi( pt2,eta2,parval );
271  double sigma_cotgth1 = MuScleFitUtils::resolutionFunction->sigmaCotgTh( pt1,eta1,parval );
272  double sigma_cotgth2 = MuScleFitUtils::resolutionFunction->sigmaCotgTh( pt2,eta2,parval );
273  double cov_pt1pt2 = MuScleFitUtils::resolutionFunction->covPt1Pt2( pt1, eta1, pt2, eta2, parval );
274 
275  // Sigma_Pt is defined as a relative sigmaPt/Pt for this reason we need to
276  // multiply it by pt.
277  double mass_res = sqrt(std::pow(dmdpt1*sigma_pt1*pt1,2)+std::pow(dmdpt2*sigma_pt2*pt2,2)+
278  std::pow(dmdphi1*sigma_phi1,2)+std::pow(dmdphi2*sigma_phi2,2)+
279  std::pow(dmdcotgth1*sigma_cotgth1,2)+std::pow(dmdcotgth2*sigma_cotgth2,2)+
280  2*dmdpt1*dmdpt2*cov_pt1pt2);
281 
282  return mass_res;
283 }
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
T sqrt(T t)
Definition: SSEVec.h:48
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
virtual double covPt1Pt2(const double &pt1, const double &eta1, const double &pt2, const double &eta2, const T &parval)
Definition: Functions.h:3031
virtual double sigmaCotgTh(const double &pt, const double &eta, const T &parval)=0
static const double mMu2
virtual double sigmaPhi(const double &pt, const double &eta, const T &parval)=0
Power< A, B >::type pow(const A &a, const B &b)
Definition: Power.h:40
static resolutionFunctionBase< double * > * resolutionFunction

Member Data Documentation

bool ErrorsPropagationAnalyzer::debug_
private

Definition at line 84 of file ErrorsPropagationAnalyzer.h.

Referenced by fillHistograms().

std::vector<int> ErrorsPropagationAnalyzer::errorFactors_
private

Definition at line 90 of file ErrorsPropagationAnalyzer.h.

Referenced by ErrorsPropagationAnalyzer(), and fillValueError().

std::vector<double> ErrorsPropagationAnalyzer::errors_
private
int ErrorsPropagationAnalyzer::etaBins_
private

Definition at line 81 of file ErrorsPropagationAnalyzer.h.

Referenced by ErrorsPropagationAnalyzer().

double ErrorsPropagationAnalyzer::etaMax_
private

Definition at line 83 of file ErrorsPropagationAnalyzer.h.

Referenced by ErrorsPropagationAnalyzer().

double ErrorsPropagationAnalyzer::etaMaxCut_
private

Definition at line 86 of file ErrorsPropagationAnalyzer.h.

Referenced by fillHistograms().

double ErrorsPropagationAnalyzer::etaMin_
private

Definition at line 82 of file ErrorsPropagationAnalyzer.h.

Referenced by ErrorsPropagationAnalyzer().

double ErrorsPropagationAnalyzer::etaMinCut_
private

Definition at line 86 of file ErrorsPropagationAnalyzer.h.

Referenced by fillHistograms().

uint32_t ErrorsPropagationAnalyzer::maxEvents_
private

Definition at line 76 of file ErrorsPropagationAnalyzer.h.

Referenced by fillHistograms().

TString ErrorsPropagationAnalyzer::outputFileName_
private

Definition at line 77 of file ErrorsPropagationAnalyzer.h.

Referenced by ~ErrorsPropagationAnalyzer().

std::vector<double> ErrorsPropagationAnalyzer::parameters_
private
int ErrorsPropagationAnalyzer::ptBins_
private

Definition at line 78 of file ErrorsPropagationAnalyzer.h.

Referenced by ErrorsPropagationAnalyzer().

double ErrorsPropagationAnalyzer::ptMax_
private

Definition at line 80 of file ErrorsPropagationAnalyzer.h.

Referenced by ErrorsPropagationAnalyzer().

double ErrorsPropagationAnalyzer::ptMaxCut_
private

Definition at line 86 of file ErrorsPropagationAnalyzer.h.

Referenced by fillHistograms().

double ErrorsPropagationAnalyzer::ptMin_
private

Definition at line 79 of file ErrorsPropagationAnalyzer.h.

Referenced by ErrorsPropagationAnalyzer().

double ErrorsPropagationAnalyzer::ptMinCut_
private

Definition at line 86 of file ErrorsPropagationAnalyzer.h.

Referenced by fillHistograms().

int ErrorsPropagationAnalyzer::resolFitType_
private

Definition at line 75 of file ErrorsPropagationAnalyzer.h.

Referenced by fillHistograms().

TProfile* ErrorsPropagationAnalyzer::sigmaMassOverMassVsEta_
private
TProfile* ErrorsPropagationAnalyzer::sigmaMassOverMassVsEtaMinusErr_
private
TProfile* ErrorsPropagationAnalyzer::sigmaMassOverMassVsEtaPlusErr_
private
TProfile* ErrorsPropagationAnalyzer::sigmaMassOverMassVsPt_
private
TProfile* ErrorsPropagationAnalyzer::sigmaMassOverMassVsPtMinusErr_
private
TProfile* ErrorsPropagationAnalyzer::sigmaMassOverMassVsPtPlusErr_
private
TProfile* ErrorsPropagationAnalyzer::sigmaMassVsEta_
private
TProfile* ErrorsPropagationAnalyzer::sigmaMassVsEtaMinusErr_
private
TProfile* ErrorsPropagationAnalyzer::sigmaMassVsEtaPlusErr_
private
TProfile* ErrorsPropagationAnalyzer::sigmaMassVsPt_
private
TProfile* ErrorsPropagationAnalyzer::sigmaMassVsPtMinusErr_
private
TProfile* ErrorsPropagationAnalyzer::sigmaMassVsPtPlusErr_
private
TProfile* ErrorsPropagationAnalyzer::sigmaPtVsEta_
private
TProfile* ErrorsPropagationAnalyzer::sigmaPtVsEtaDiff_
private
TProfile* ErrorsPropagationAnalyzer::sigmaPtVsEtaMinusErr_
private
TProfile* ErrorsPropagationAnalyzer::sigmaPtVsEtaPlusErr_
private
TProfile* ErrorsPropagationAnalyzer::sigmaPtVsPt_
private
TProfile* ErrorsPropagationAnalyzer::sigmaPtVsPtDiff_
private
TProfile* ErrorsPropagationAnalyzer::sigmaPtVsPtMinusErr_
private
TProfile* ErrorsPropagationAnalyzer::sigmaPtVsPtPlusErr_
private
TString ErrorsPropagationAnalyzer::treeFileName_
private

Definition at line 74 of file ErrorsPropagationAnalyzer.h.

Referenced by fillHistograms().

std::vector<double> ErrorsPropagationAnalyzer::valueMinusError_
private

Definition at line 93 of file ErrorsPropagationAnalyzer.h.

Referenced by fillValueError().

std::vector<double> ErrorsPropagationAnalyzer::valuePlusError_
private

Definition at line 92 of file ErrorsPropagationAnalyzer.h.

Referenced by fillValueError().