CMS 3D CMS Logo
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
Functions
ResponseFitter.cc File Reference
#include <iostream>
#include <string.h>
#include <fstream>
#include <cmath>
#include <TFile.h>
#include <TH1F.h>
#include <TMath.h>
#include <TKey.h>
#include <TList.h>
#include "Utilities.h"

Go to the source code of this file.

Functions

int main (int argc, char **argv)
 

Function Documentation

int main ( int  argc,
char **  argv 
)

Definition at line 14 of file ResponseFitter.cc.

References EnergyCorrector::c, alignmentValidation::c1, CalculateCorrection(), CalculateResponse(), CommandLine::check(), gather_cfg::cout, alignCSCRings::e, GetMEAN(), GetMPV(), CommandLine::getValue(), CommandLine::getVector(), h, HistoExists(), EcalCondDB::inf, j, relval_steps::key, mergeVDriftHistosByStation::name, NETA, GetRecoTauVFromDQM_MC_cff::next, CommandLine::parse(), CommandLine::print(), alignCSCRings::r, and AlCaHLTBitMon_QueryRunRegistry::string.

15 {
16  CommandLine c1;
17  c1.parse(argc,argv);
18 
19  std::string HistoFilename = c1.getValue<string>("HistoFilename");
20  std::string FitterFilename = c1.getValue<string>("FitterFilename");
21  std::string L3ResponseTxtFilename = c1.getValue<string>("L3ResponseTxtFilename");
22  std::string L3CorrectionTxtFilename = c1.getValue<string>("L3CorrectionTxtFilename");
23  std::string L3OutputROOTFilename = c1.getValue<string>("L3OutputROOTFilename");
24  std::string L2CorrectionTxtFilename = c1.getValue<string>("L2CorrectionTxtFilename");
25  std::string L2OutputROOTFilename = c1.getValue<string>("L2OutputROOTFilename");
26  bool UseRatioForResponse = c1.getValue<bool>("UseRatioForResponse");
27  std::vector<double> pt_vec = c1.getVector<double>("RefPtBoundaries");
28  std::vector<double> eta_vec = c1.getVector<double>("EtaBoundaries");
29  if (!c1.check()) return 0;
30  c1.print();
32  const int MAX_NETA = 83;
33  const int MAX_NREFPTBINS = 30;
34  int NRefPtBins = pt_vec.size()-1;
35  int NETA = eta_vec.size()-1;
36  char name[100];
37  int j,etabin;
38  double e,mR,eR,sR,seR,mRBarrel,eRBarrel,sRBarrel,seRBarrel,r,c;
39  double mCaloPt,eCaloPt,sCaloPt,mRefPt,eRefPt,sRefPt;
40  double mRefPtEtaBin,eRefPtEtaBin,sRefPtEtaBin,mCaloPtEtaBin,eCaloPtEtaBin,sCaloPtEtaBin;
41  double EtaBoundaries[MAX_NETA],RefPtBoundaries[MAX_NREFPTBINS];
42  std::vector<std::string> HistoNamesList;
43  for(j=0;j<=NRefPtBins;j++)
44  RefPtBoundaries[j] = pt_vec[j];
45  for(j=0;j<=NETA;j++)
46  EtaBoundaries[j] = eta_vec[j];
47  TFile *inf;//Input file containing the response histograms.
48  TFile *outf;//Output file containing the fitter results.
49  TH1F *BarrelResponse;//Histogram with the barrel response in RefPt bins.
50  TH1F *BarrelCorrection;//Histogram with the barrel correction in RefPt bins.
51  TH1F *MeanRefPt_Barrel;//Histogram with the barrel average RefPt in RefPt bins.
52  TH1F *MeanCaloPt_Barrel;//Histogram with the barrel average CaloPt in RefPt bins.
53  TH1F *MeanRefPt_EtaBin[MAX_NETA];//Histograms with the average RefPt in Eta & RefPt bins.
54  TH1F *MeanCaloPt_EtaBin[MAX_NETA];//Histograms with the average CaloPt in Eta & RefPt bins.
55  TH1F *ResponseVsEta_RefPt[MAX_NREFPTBINS];//Histograms with the average response vs Eta in RefPt bins.
56  TH1F *CorrectionVsEta_RefPt[MAX_NREFPTBINS];//Histograms with the average correction vs Eta in RefPt bins.
57  TH1F *h;//Auxilary histogram;
58  TKey *key;
60  inf = new TFile(HistoFilename.c_str(),"r");
61  if (inf->IsZombie()) return(0);
62  TIter next(inf->GetListOfKeys());
63  while ((key = (TKey*)next()))
64  HistoNamesList.push_back(key->GetName());
65  outf = new TFile(FitterFilename.c_str(),"RECREATE");
66  TDirectory *dir_Response = (TDirectory*)outf->mkdir("FittedHistograms");//Directory in output file to store the fitted histograms.
67  BarrelResponse = new TH1F("Response","Response",NRefPtBins,RefPtBoundaries);
68  BarrelCorrection = new TH1F("Correction","Correction",NRefPtBins,RefPtBoundaries);
69  MeanRefPt_Barrel = new TH1F("MeanRefPt","MeanRefPt",NRefPtBins,RefPtBoundaries);
70  MeanCaloPt_Barrel = new TH1F("MeanCaloPt","MeanCaloPt",NRefPtBins,RefPtBoundaries);
71  if (NETA>1)//multiple eta bins: used for L2+L3 correction calculation
72  {
73  std::cout<<"************* Fitting Response Histograms in multiple Eta bins. ************"<<std::endl;
74  for(etabin=0;etabin<NETA;etabin++)
75  {
76  sprintf(name,"MeanRefPt_Eta%d",etabin);
77  MeanRefPt_EtaBin[etabin] = new TH1F(name,name,NRefPtBins,RefPtBoundaries);
78  sprintf(name,"MeanCaloPt_Eta%d",etabin);
79  MeanCaloPt_EtaBin[etabin] = new TH1F(name,name,NRefPtBins,RefPtBoundaries);
80  }
81  for (j=0; j<NRefPtBins; j++)//loop over RefPt bins
82  {
83  std::cout<<"RefJetPt Bin: ["<<RefPtBoundaries[j]<<","<<RefPtBoundaries[j+1]<<"] GeV"<<std::endl;
84  sprintf(name,"ptRef_RefPt%d_Barrel",j);
85  if (!HistoExists(HistoNamesList,name)) return(0);
86  h = (TH1F*)inf->Get(name);
87  GetMEAN(h,mRefPt,eRefPt,sRefPt);
88  sprintf(name,"ptCalo_RefPt%d_Barrel",j);
89  if (!HistoExists(HistoNamesList,name)) return(0);
90  h = (TH1F*)inf->Get(name);
91  GetMEAN(h,mCaloPt,eCaloPt,sCaloPt);
92  sprintf(name,"Response_RefPt%d_Barrel",j);
93  if (!HistoExists(HistoNamesList,name)) return(0);
94  h = (TH1F*)inf->Get(name);
95  GetMPV(name,h,dir_Response,mRBarrel,eRBarrel,sRBarrel,seRBarrel);
97  MeanRefPt_Barrel->SetBinContent(j+1,mRefPt);
98  MeanRefPt_Barrel->SetBinError(j+1,eRefPt);
100  MeanCaloPt_Barrel->SetBinContent(j+1,mCaloPt);
101  MeanCaloPt_Barrel->SetBinError(j+1,eCaloPt);
103  CalculateResponse(UseRatioForResponse,mRefPt,eRefPt,mRBarrel,eRBarrel,r,e);
104  BarrelResponse->SetBinContent(j+1,r);
105  BarrelResponse->SetBinError(j+1,e);
107  CalculateCorrection(UseRatioForResponse,mRefPt,eRefPt,mRBarrel,eRBarrel,c,e);
108  BarrelCorrection->SetBinContent(j+1,c);
109  BarrelCorrection->SetBinError(j+1,e);
111  sprintf(name,"Response_vs_Eta_RefPt%d",j);
112  ResponseVsEta_RefPt[j] = new TH1F(name,name,NETA,EtaBoundaries);
113  sprintf(name,"Correction_vs_Eta_RefPt%d",j);
114  CorrectionVsEta_RefPt[j] = new TH1F(name,name,NETA,EtaBoundaries);
115  for(etabin=0;etabin<NETA;etabin++)//loop over eta bins
116  {
118  sprintf(name,"Response_RefPt%d_Eta%d",j,etabin);
119  if (!HistoExists(HistoNamesList,name)) return(0);
120  h = (TH1F*)inf->Get(name);
121  GetMPV(name,h,dir_Response,mR,eR,sR,seR);
122  sprintf(name,"ptRef_RefPt%d_Eta%d",j,etabin);
123  if (!HistoExists(HistoNamesList,name)) return(0);
124  h = (TH1F*)inf->Get(name);
125  GetMEAN(h,mRefPtEtaBin,eRefPtEtaBin,sRefPtEtaBin);
126  sprintf(name,"ptCalo_RefPt%d_Eta%d",j,etabin);
127  if (!HistoExists(HistoNamesList,name)) return(0);
128  h = (TH1F*)inf->Get(name);
129  GetMEAN(h,mCaloPtEtaBin,eCaloPtEtaBin,sCaloPtEtaBin);
131  MeanRefPt_EtaBin[etabin]->SetBinContent(j+1,mRefPtEtaBin);
132  MeanRefPt_EtaBin[etabin]->SetBinError(j+1,eRefPtEtaBin);
134  MeanCaloPt_EtaBin[etabin]->SetBinContent(j+1,mCaloPtEtaBin);
135  MeanCaloPt_EtaBin[etabin]->SetBinError(j+1,eCaloPtEtaBin);
137  CalculateResponse(UseRatioForResponse,mRefPtEtaBin,eRefPtEtaBin,mR,eR,r,e);
138  ResponseVsEta_RefPt[j]->SetBinContent(etabin+1,r);
139  ResponseVsEta_RefPt[j]->SetBinError(etabin+1,e);
141  CalculateCorrection(UseRatioForResponse,mRefPtEtaBin,eRefPtEtaBin,mR,eR,c,e);
142  CorrectionVsEta_RefPt[j]->SetBinContent(etabin+1,c);
143  CorrectionVsEta_RefPt[j]->SetBinError(etabin+1,e);
144  }//end of EtaBin loop
145  }// end of Pt loop
146  }
147  else//single eta bin: used for L3 correction calculation
148  {
149  std::cout<<"************* Fitting Response Histograms in single eta bin. ************"<<std::endl;
150  for (j=0; j<NRefPtBins; j++)//loop over Pt bins
151  {
152  std::cout<<"RefJetPt Bin: ["<<RefPtBoundaries[j]<<","<<RefPtBoundaries[j+1]<<"] GeV"<<std::endl;
153  sprintf(name,"ptRef_RefPt%d",j);
154  if (!HistoExists(HistoNamesList,name)) return(0);
155  h = (TH1F*)inf->Get(name);
156  GetMEAN(h,mRefPt,eRefPt,sRefPt);
157  sprintf(name,"ptCalo_RefPt%d",j);
158  if (!HistoExists(HistoNamesList,name)) return(0);
159  h = (TH1F*)inf->Get(name);
160  GetMEAN(h,mCaloPt,eCaloPt,sCaloPt);
161  sprintf(name,"Response_RefPt%d",j);
162  if (!HistoExists(HistoNamesList,name)) return(0);
163  h = (TH1F*)inf->Get(name);
164  GetMPV(name,h,dir_Response,mRBarrel,eRBarrel,sRBarrel,seRBarrel);
166  MeanRefPt_Barrel->SetBinContent(j+1,mRefPt);
167  MeanRefPt_Barrel->SetBinError(j+1,eRefPt);
169  MeanCaloPt_Barrel->SetBinContent(j+1,mCaloPt);
170  MeanCaloPt_Barrel->SetBinError(j+1,eCaloPt);
172  CalculateResponse(UseRatioForResponse,mRefPt,eRefPt,mRBarrel,eRBarrel,r,e);
173  BarrelResponse->SetBinContent(j+1,r);
174  BarrelResponse->SetBinError(j+1,e);
176  CalculateCorrection(UseRatioForResponse,mRefPt,eRefPt,mRBarrel,eRBarrel,c,e);
177  BarrelCorrection->SetBinContent(j+1,c);
178  BarrelCorrection->SetBinError(j+1,e);
179  }// end of Pt loop
180  }
182  outf->cd();
183  MeanRefPt_Barrel->Write();
184  MeanCaloPt_Barrel->Write();
185  BarrelResponse->Write();
186  BarrelCorrection->Write();
187  if (NETA>1)
188  {
189  for(etabin=0;etabin<NETA;etabin++)
190  {
191  MeanRefPt_EtaBin[etabin]->Write();
192  MeanCaloPt_EtaBin[etabin]->Write();
193  }
194  for(j=0;j<NRefPtBins;j++)
195  {
196  ResponseVsEta_RefPt[j]->Write();
197  CorrectionVsEta_RefPt[j]->Write();
198  }
199  }
200  outf->Close();
201 }
CalculateResponse
void CalculateResponse(bool UseRatioForResponse, double x, double ex, double y, double ey, double &r, double &e)
Definition: Utilities.h:452
NETA
const int NETA
Definition: CaloTowerBoundries.h:4
CommandLine::check
bool check()
Definition: Utilities.h:240
CommandLine::getVector
std::vector< T > getVector(const std::string &name)
Definition: Utilities.h:144
AlCaHLTBitMon_QueryRunRegistry.string
string string
Definition: AlCaHLTBitMon_QueryRunRegistry.py:255
EcalCondDB.inf
string inf
Definition: EcalCondDB.py:94
j
int j
Definition: DBlmapReader.cc:9
CommandLine::print
void print()
Definition: Utilities.h:262
h
The Signals That Services Can Subscribe To This is based on ActivityRegistry h
Helper function to determine trigger accepts.
Definition: Activities.doc:4
CommandLine::getValue
T getValue(const std::string &name)
Definition: Utilities.h:74
relval_steps.key
string key
FastSim: produces sample of signal events, overlayed with premixed minbias events.
Definition: relval_steps.py:736
GetMPV
void GetMPV(char name[100], TH1F *histo, TDirectory *Dir, double &peak, double &error, double &sigma, double &err_sigma)
Definition: Utilities.h:368
dir2webdir.argc
tuple argc
Definition: dir2webdir.py:38
GetMEAN
void GetMEAN(TH1F *histo, double &peak, double &error, double &sigma)
Definition: Utilities.h:435
gather_cfg.cout
tuple cout
Definition: gather_cfg.py:121
GetRecoTauVFromDQM_MC_cff.next
tuple next
Definition: GetRecoTauVFromDQM_MC_cff.py:30
CalculateCorrection
void CalculateCorrection(bool UseRatioForResponse, double x, double ex, double y, double ey, double &c, double &e)
Definition: Utilities.h:474
CommandLine::parse
bool parse(int argc, char **argv)
Definition: Utilities.h:199
HistoExists
bool HistoExists(std::vector< std::string > LIST, std::string hname)
Definition: Utilities.h:496
Generated for CMSSW Reference Manual by doxygen 1.8.5