#include <XANADOO/domain/interface/XANAEsdEvent.h>#include <XANADOO/XANAMcInfo/interface/XANAGeneratorEvent.h>#include <XANADOO/XANAMcInfo/interface/XANAGeantEvent.h>#include <XANADOO/XANAMcInfo/interface/XANAGeneratorParticle.h>#include <XANADOO/XANAElectronCandidate/interface/XANAElectronCandidate.h>#include <XANADOO/XANAElectronCandidate/interface/XANAElectronTrack.h>#include <XANADOO/XANAElectronCandidate/interface/XANAElectronGSFTrack.h>#include <XANADOO/XANAVertex/interface/XANAVertex.h>#include <XANADOO/XANATracks/interface/XANATrack.h>#include <XANADOO/XANAnalysisTools/interface/XANAmbiguityResolve.h>#include <XANADOO/XANAnalysisZZ/interface/XANASelectHZZ.h>#include <XANADOO/XANAnalysisZZ/interface/XANASelectMcHZZ.h>#include <XANADOO/XANAnalysisTools/interface/XANAEvtIsolation.h>#include <XANADOO/XANAnalysisTools/interface/XANAEvtRegionIsolation.h>#include <XANADOO/XANAnalysisTools/interface/XANAIsolation.h>#include <XANADOO/XANAnalysisTools/interface/XANAnalysisTrigger.h>#include <XANADOO/XANATriggerInfo/interface/XANATriggerInfo.h>#include <XANADOO/XANAnalysisTools/interface/XANAElectronRecGen.h>#include <XANADOO/XANAnalysisTools/interface/XANAEventSelect.h>#include <XANADOO/XANAnalysisTools/interface/XANAElectronClassification.h>#include <XANADOO/XANAnalysisTools/interface/XANAElectronSelect.h>#include <XANADOO/XANAnalysisTools/interface/XANAElectronCorrection.h>#include <XANADOO/XANAnalysisTools/interface/XANAElectronMomentum.h>#include "CLHEP/Geometry/Point3D.h"#include "TTree.h"#include "TFile.h"#include "TBranch.h"#include "TObjArray.h"#include "TH1.h"#include "TH2.h"#include "TProfile.h"#include "TStopwatch.h"#include "TChain.h"#include <TBits.h>#include <string>Include dependency graph for hzz_4e_analysis.cpp:
Go to the source code of this file.
Functions | |
| void | initSettings (TChain *T, XANAEsdEvent **eventaddress, XANAGeneratorEvent **generatorEventaddress, XANAGeantEvent **geantEventaddress) |
| void | clean (XANAEsdEvent *event, XANAGeneratorEvent *generatorEvent, XANAGeantEvent *geantEvent) |
| int | main (int argc, char **argv) |
|
||||||||||||||||
|
Definition at line 81 of file hzz_4e_analysis.cpp. References XANAGeantEvent::clear(), XANAGeneratorEvent::clear(), and XANAEsdEvent::clear(). Referenced by main().
|
|
||||||||||||||||||||
|
Definition at line 46 of file hzz_4e_analysis.cpp. Referenced by main().
00047 {
00048 // example of how to suppress reading of branches
00049 // has to be put before setting the branch addresses!
00050
00051 // To improve the speed set status for all branches to 0 and then include those needed
00052
00053 // T->SetBranchStatus("Event.*",0);
00054 // T->SetBranchStatus("GeneratorEvent.*",0);
00055 T->SetBranchStatus("GeantEvent.*",0);
00056
00057 //T->SetBranchStatus("Event.caloRecHits_",0);
00058 T->SetBranchStatus("Event.allTrackHits_",0);
00059 // T->SetBranchStatus("Event.emClusters_",0);
00060 T->SetBranchStatus("Event.preshClusters_",0);
00061 T->SetBranchStatus("Event.hadClusters_",0);
00062 T->SetBranchStatus("Event.vertices_",0);
00063 T->SetBranchStatus("Event.muonCandidates_",0);
00064 T->SetBranchStatus("Event.muonTracks_",0);
00065 T->SetBranchStatus("Event.jets_",0);
00066 T->SetBranchStatus("Event.met_",0);
00067
00068 /*
00069 T->SetBranchStatus("Event.electronCandidates_",1);
00070 T->SetBranchStatus("Event.numberOfElectronCandidates_",1);
00071 T->SetBranchStatus("Event.tracks_",1);
00072 T->SetBranchStatus("Event.triggerInfo_",1);
00073 */
00074
00075 // Set the branch addresses: first reconstructed event and then MC information (generator and GEANT event)
00076 T->SetBranchAddress("Event.",eventaddress);
00077 T->SetBranchAddress("GeneratorEvent.",generatorEventaddress);
00078 T->SetBranchAddress("GeantEvent.",geantEventaddress);
00079 }
|
|
||||||||||||
|
Definition at line 88 of file hzz_4e_analysis.cpp. References clean(), XANAGeantEvent::clear(), XANAGeneratorEvent::clear(), XANAEsdEvent::clear(), electrons, XANAElectronSelect::ElSelect_dphi(), XANATrack::getAlgoName(), XANAEsdEvent::getElectronCandidates(), XANAElectronClassification::getElectronClass(), XANASelectMcHZZ::getElectronFromRealZ(), XANASelectHZZ::getElectronFromRealZ(), XANASelectMcHZZ::getElectronFromRealZAftPh(), XANASelectMcHZZ::getElectronFromVirtualZ(), XANASelectHZZ::getElectronFromVirtualZ(), XANASelectMcHZZ::getElectronFromVirtualZAftPh(), XANAElectronRecGen::getElectronsRecGen(), XANAElectronCandidate::getElectronTrack(), XANASuperCluster::getEnergy(), XANASuperCluster::getEnergyScaleFactor(), XANAEventSelect::getGoodElCands(), XANAElectronGSFTrack::getGSFMomentumAtVertexMode(), XANASuperCluster::getHadronicOverEm(), XANASelectMcHZZ::getHiggs(), XANASelectHZZ::getHiggs(), XANAnalysisTrigger::getHLTDecision(), XANAnalysisTrigger::getHLTElecDecision(), XANAEvtIsolation::getIsolationPtSum(), XANAEvtIsolation::getIsolationVector(), XANAnalysisTrigger::getL1Decision(), XANAGeneratorParticle::getMomentum(), XANAElectronCandidate::getMomentumAtVertex(), XANATrack::getMomentumAtVertex(), XANAEsdEvent::getNumberOfElectronCandidates(), XANASuperCluster::getPosition(), XANASelectMcHZZ::getPositronFromRealZ(), XANASelectHZZ::getPositronFromRealZ(), XANASelectMcHZZ::getPositronFromRealZAftPh(), XANASelectMcHZZ::getPositronFromVirtualZ(), XANASelectHZZ::getPositronFromVirtualZ(), XANASelectMcHZZ::getPositronFromVirtualZAftPh(), XANASelectMcHZZ::getRealZ(), XANASelectHZZ::getRealZ(), XANAmbiguityResolve::getResolvedElcands(), XANAElectronSelect::getSelectedElectrons(), XANAElectronCandidate::getSuperCluster(), XANAElectronCandidate::getSuperClusterEnergy(), XANAElectronCandidate::getSuperClusterPosition(), XANAEsdEvent::getSuperClusters(), XANAEsdEvent::getTracks(), XANAEsdEvent::getTriggerInfo(), XANASelectMcHZZ::getVirtualZ(), XANASelectHZZ::getVirtualZ(), initSettings(), XANASelectMcHZZ::isFourElectrons(), XANASelectMcHZZ::isFourTaus(), XANASelectHZZ::isHiggsMassCut(), isolationVector, XANASelectHZZ::isPtCut(), XANASelectHZZ::isRealZMassCut(), XANASelectMcHZZ::isTwoElectronsTwoTaus(), XANASelectHZZ::isValid(), XANASelectHZZ::isVirtualZMassCut(), XANAEventSelect::NumberOfElCandSelect(), XANAmbiguityResolve::ResolveByEoverP(), XANAElectronCandidate::setESuperClusterOverP(), XANAElectronCandidate::setHadronicOverEm(), XANAElectronCandidate::setMomentumAtVertex(), XANAElectronCandidate::setSuperCluster(), XANAElectronCandidate::setSuperClusterEnergy(), superClusters, and tracks.
00089 {
00090
00091 if(!argv[1]) {
00092 cout<<">>> No argument, please give an argument if you want me to do something <<<"<<endl;
00093 cout<<">>> list of possible arguments for the moment: hzz120 hzz150 and tt4e <<<"<<endl;
00094 return 0;
00095 }
00096 cout<< ">>> Hello, welcome in hzz_4e_analysis program <<<"<<endl;
00097
00098 TStopwatch stopWatch;
00099 stopWatch.Start();
00100
00101 // Steering flags
00102 bool momentumCorrection = true;
00103 bool federicoEscaleCorrection = true;
00104 bool egammaEscaleCorrection = true;
00105
00106
00107 // output file definition
00108 char name[80];
00109 if ( federicoEscaleCorrection && momentumCorrection )
00110 sprintf(name,"HZZ_corr_momentum_%s.root",argv[1]);
00111 else if ( federicoEscaleCorrection && !momentumCorrection )
00112 sprintf(name,"HZZ_corr_no_momentum_%s.root",argv[1]);
00113 else if ( !federicoEscaleCorrection && !momentumCorrection )
00114 sprintf(name,"HZZ_no_corr_no_momentum_%s.root",argv[1]);
00115 else if ( !federicoEscaleCorrection && !egammaEscaleCorrection && !momentumCorrection )
00116 sprintf(name,"HZZ_no_corr_at_all_%s.root",argv[1]);
00117
00118
00119 TFile f(name,"RECREATE");
00120 f.cd();
00121
00122 // Histogram booking
00123
00124 // general info
00125 TH1F hElectronNumber("hElectronNumber", "Number of ElectronCandidates", 15, 0., 15.);
00126 TH1F hElectronNumberAfterAmb("hElectronNumberAfterAmb", "Number of ElectronCandidates after ambiguity resolving", 15, 0., 15.);
00127 TH1F hElectronNumberBeforeAfter("hElectronNumberBeforeAfter", "Number of ElectronCandidates before - after ambiguity resolving", 15, 0., 15.);
00128
00129 // isolation
00130 TH1F hIsoNumberTracks("hIsoNumberTracks", "Number of tracks inside cone", 40, 0., 40.);
00131 TH1F hIsoPtTracks("hIsoPtTracks", "Pt summ of all tracks inside cone", 100, 0., 100.);
00132 TH1F hIsoNumberECinsideCone("hIsoNumberECinsideCone", "Number of e candidates inside cone", 10, 0., 10.);
00133 TH1F hIsoPtMax("hIsoPtMax", "Max pT of tracks inside cone", 100, 0., 100.);
00134
00135 TH1F hIsoRegionNumberTracks("hIsoRegionNumberTracks", "Number of tracks inside region", 40, 0., 40.);
00136 TH1F hIsoRegionPtTracks("hIsoRegionPtTracks", "Pt summ of all tracks inside region", 100, 0., 100.);
00137 TH1F hIsoRegionPtMax("hIsoRegionPtMax", "Max pT of tracks inside region", 100, 0., 100.);
00138
00139 // Mc information
00140 TH1F hElectronRecGenPt("hElectronRecGenPt","ElectronCandidate: Rec/Gen pT", 500, 0., 5.);
00141 TH1F hElectronRecGenDr("hElectronRecGenDr","ElectronCandidate: Rec/Gen Delta R", 500, 0., 5.);
00142
00143 // classification
00144 TH1F hElectronClasses("hElectronClasses","Electron Candidates Classification",23, -2., 21.);
00145 TH1F hElectronClass0("hElectronClass0","Electron Candidates Class: 0, Barrel-Golden",100, 0., 1.5);
00146 TH1F hElectronClass1("hElectronClass1","Electron Candidates Class: 1, Barrel-Big Brem",100, 0., 1.5);
00147 TH1F hElectronClass2("hElectronClass2","Electron Candidates Class: 2, Barrel-Narrow",100, 0., 1.5);
00148 TH1F hElectronClass3("hElectronClass3","Electron Candidates Class: 3, Barrel-Showering",100, 0., 1.5);
00149 TH1F hElectronClass4("hElectronClass4","Electron Candidates Class: 4, Barrel-Cracks",100, 0., 1.5);
00150 TH1F hElectronClass10("hElectronClass10","Electron Candidates Class: 0, Endcap-Golden",100, 0., 1.5);
00151 TH1F hElectronClass11("hElectronClass11","Electron Candidates Class: 1, Endcap-Big Brem",100, 0., 1.5);
00152 TH1F hElectronClass12("hElectronClass12","Electron Candidates Class: 2, Endcap-Narrow",100, 0., 1.5);
00153 TH1F hElectronClass13("hElectronClass13","Electron Candidates Class: 3, Endcap-Showering",100, 0., 1.5);
00154
00155 TProfile hElectronClassProfile0("hElectronClassProfile0","Electron Candidates ClassProfile: 0, Barrel-Golden",100, 0., 2.7);
00156 TProfile hElectronClassProfile1("hElectronClassProfile1","Electron Candidates ClassProfile: 1, Barrel-Big Brem",100, 0., 2.7);
00157 TProfile hElectronClassProfile2("hElectronClassProfile2","Electron Candidates ClassProfile: 2, Barrel-Narrow",100, 0., 2.7);
00158 TProfile hElectronClassProfile3("hElectronClassProfile3","Electron Candidates ClassProfile: 3, Barrel-Showering",100, 0., 2.7);
00159 TProfile hElectronClassProfile4("hElectronClassProfile4","Electron Candidates ClassProfile: 4, Barrel-Cracks",100, 0., 2.7);
00160 TProfile hElectronClassProfile10("hElectronClassProfile10","Electron Candidates ClassProfile: 0, Endcap-Golden",100, 0., 2.7);
00161 TProfile hElectronClassProfile11("hElectronClassProfile11","Electron Candidates ClassProfile: 1, Endcap-Big Brem",100, 0., 2.7);
00162 TProfile hElectronClassProfile12("hElectronClassProfile12","Electron Candidates ClassProfile: 2, Endcap-Narrow",100, 0., 2.7);
00163 TProfile hElectronClassProfile13("hElectronClassProfile13","Electron Candidates ClassProfile: 3, Endcap-Showering",100, 0., 2.7);
00164
00165 TH1F hElectronEta("hElectronEta","Electron Candidates Eta",100, 0., 2.7);
00166 TH1F hElectronEtaClass0("hElectronEtaClass0","Electron Candidates Eta Class: 0, Barrel-Golden",100, 0., 2.7);
00167 TH1F hElectronEtaClass1("hElectronEtaClass1","Electron Candidates Eta Class: 1, Barrel-Big Brem",100, 0., 2.7);
00168 TH1F hElectronEtaClass2("hElectronEtaClass2","Electron Candidates Eta Class: 2, Barrel-Narrow",100, 0., 2.7);
00169 TH1F hElectronEtaClass3("hElectronEtaClass3","Electron Candidates Eta Class: 3, Barrel-Showering",100, 0., 2.7);
00170 TH1F hElectronEtaClass4("hElectronEtaClass4","Electron Candidates Eta Class: 4, Barrel-Cracks",100, 0., 2.7);
00171 TH1F hElectronEtaClass10("hElectronEtaClass10","Electron Candidates Eta Class: 0, Endcap-Golden",100, 0., 2.7);
00172 TH1F hElectronEtaClass11("hElectronEtaClass11","Electron Candidates Eta Class: 1, Endcap-Big Brem",100, 0., 2.7);
00173 TH1F hElectronEtaClass12("hElectronEtaClass12","Electron Candidates Eta Class: 2, Endcap-Narrow",100, 0., 2.7);
00174 TH1F hElectronEtaClass13("hElectronEtaClass13","Electron Candidates Eta Class: 3, Endcap-Showering",100, 0., 2.7);
00175
00176 // after corrections
00177 TH1F hElectronClassCorr0("hElectronClassCorr0","Electron Candidates Class After Corr: 0, Barrel-Golden",100, 0., 1.5);
00178 TH1F hElectronClassCorr1("hElectronClassCorr1","Electron Candidates Class After Corr: 1, Barrel-Big Brem",100, 0., 1.5);
00179 TH1F hElectronClassCorr2("hElectronClassCorr2","Electron Candidates Class After Corr: 2, Barrel-Narrow",100, 0., 1.5);
00180 TH1F hElectronClassCorr3("hElectronClassCorr3","Electron Candidates Class After Corr: 3, Barrel-Showering",100, 0., 1.5);
00181 TH1F hElectronClassCorr4("hElectronClassCorr4","Electron Candidates Class After Corr: 4, Barrel-Cracks",100, 0., 1.5);
00182 TH1F hElectronClassCorr10("hElectronClassCorr10","Electron Candidates Class After Corr: 0, Endcap-Golden",100, 0., 1.5);
00183 TH1F hElectronClassCorr11("hElectronClassCorr11","Electron Candidates Class After Corr: 1, Endcap-Big Brem",100, 0., 1.5);
00184 TH1F hElectronClassCorr12("hElectronClassCorr12","Electron Candidates Class After Corr: 2, Endcap-Narrow",100, 0., 1.5);
00185 TH1F hElectronClassCorr13("hElectronClassCorr13","Electron Candidates Class After Corr: 3, Endcap-Showering",100, 0., 1.5);
00186
00187 TProfile hElectronClassProfileCorr0("hElectronClassProfileCorr0","Electron Candidates ClassProfileCorr: 0, Barrel-Golden",100, 0., 2.7);
00188 TProfile hElectronClassProfileCorr1("hElectronClassProfileCorr1","Electron Candidates ClassProfileCorr: 1, Barrel-Big Brem",100, 0., 2.7);
00189 TProfile hElectronClassProfileCorr2("hElectronClassProfileCorr2","Electron Candidates ClassProfileCorr: 2, Barrel-Narrow",100, 0., 2.7);
00190 TProfile hElectronClassProfileCorr3("hElectronClassProfileCorr3","Electron Candidates ClassProfileCorr: 3, Barrel-Showering",100, 0., 2.7);
00191 TProfile hElectronClassProfileCorr4("hElectronClassProfileCorr4","Electron Candidates ClassProfileCorr: 4, Barrel-Cracks",100, 0., 2.7);
00192 TProfile hElectronClassProfileCorr10("hElectronClassProfileCorr10","Electron Candidates ClassProfileCorr: 0, Endcap-Golden",100, 0., 2.7);
00193 TProfile hElectronClassProfileCorr11("hElectronClassProfileCorr11","Electron Candidates ClassProfileCorr: 1, Endcap-Big Brem",100, 0., 2.7);
00194 TProfile hElectronClassProfileCorr12("hElectronClassProfileCorr12","Electron Candidates ClassProfileCorr: 2, Endcap-Narrow",100, 0., 2.7);
00195 TProfile hElectronClassProfileCorr13("hElectronClassProfileCorr13","Electron Candidates ClassProfileCorr: 3, Endcap-Showering",100, 0., 2.7);
00196
00197
00198 // Momentum calculations
00199 TH1F hElectronMomentumEOverP("hElectronMomentumEOverP","E/p",200,0.,2.);
00200 TH1F hElectronMomentumEnergy("hElectronMomentumEnergy","E/Etrue",200,0.,2.);
00201 TH1F hElectronMomentumTrack("hElectronMomentumTrack","p/Etrue",200,0.,2.);
00202 TH1F hElectronMomentumFinal("hElectronMomentumFinal","Efinal/Etrue",200,0.,2.);
00203
00204 TH2F hElectronMomentumEnergy2("hElectronMomentumEnergy2","E/Etrue vs E/p",200,0.,2.,200,0.,2.);
00205 TH2F hElectronMomentumTrack2("hElectronMomentumTrack2","p/Etrue vs E/p",200,0.,2.,200,0.,2.5);
00206 TH2F hElectronMomentumFinal2("hElectronMomentumFinal2","Efinal/Etrue vs E/p",200,0.,2.,200,0.,2.);
00207
00208 // Higgs analysis
00209 TH1F hRealZMass("hRealZMass","Higgs: Real Z mass", 150, 0., 150.);
00210 TH1F hVirtualZMass("hVirtualZMass","Higgs: Virtual Z mass", 150, 0., 150.);
00211 TH1F hHiggsMass("hHiggsMass","Higgs: Higgs mass", 180, 0., 180.);
00212 TH1F hElectronPt("hElectronPt","XANASelectHZZ: Electron pT", 100, 0., 100.);
00213 TH1F hElectronRecGenPtHiggs("hElectronRecGenPtHiggs","ElectronCandidate for Higgs: Rec/Gen pT", 500, 0., 5.);
00214 TH1F hElectronRecGenDrHiggs("hElectronRecGenDrHiggs","ElectronCandidate for Higgs: Rec/Gen Delta R", 500, 0., 5.);
00215
00216 TH1F hRealZMassCuts("hRealZMassCuts","Higgs: Real Z mass, after cuts", 150, 0., 150.);
00217 TH1F hVirtualZMassCuts("hVirtualZMassCuts","Higgs: Virtual Z mass, after cuts", 150, 0., 150.);
00218 TH1F hHiggsMassCuts("hHiggsMassCuts","Higgs: Higgs mass, after cuts", 180, 0., 180.);
00219 TH1F hElectronPtCuts("hElectronPtCuts","XANASelectHZZ: Electron pT, after cuts", 100, 0., 100.);
00220 TH1F hElectronRecGenPtHiggsCuts("hElectronRecGenPtHiggsCuts","ElectronCandidate for Higgs: Rec/Gen pT, after cuts", 500, 0., 5.);
00221 TH1F hElectronRecGenDrHiggsCuts("hElectronRecGenDrHiggsCuts","ElectronCandidate for Higgs: Rec/Gen Delta R, after cuts", 500, 0., 5.);
00222
00223 TH1F hElectronRecGenPtHiggsOrig("hElectronRecGenPtHiggsOrig","ElectronCandidate for Higgs: Rec/Gen (Original) pT", 500, 0., 5.);
00224 TH1F hElectronRecGenDrHiggsOrig("hElectronRecGenDrHiggsOrig","ElectronCandidate for Higgs: Rec/Gen (Original) Delta R", 500, 0., 5.);
00225 TH1F hElectronRecGenPtHiggsOrigAftPh("hElectronRecGenPtHiggsOrigAftPh","ElectronCandidate for Higgs: Rec/Gen (Original) pT, After Photos", 500, 0., 5.);
00226 TH1F hElectronRecGenDrHiggsOrigAftPh("hElectronRecGenDrHiggsOrigAftPh","ElectronCandidate for Higgs: Rec/Gen (Original) Delta R, After Photos", 500, 0., 5.);
00227
00228 TH1F hElectronPhotosAftBefPt("hElectronPhotosAftBefPt","Photos: After/Before pT", 200, 0., 2.);
00229 TH1F hElectronPhotosAftBefE("hElectronPhotosAftBefE","Photos: After/Before E", 200, 0., 2.);
00230
00231 TH1F hRealZMassRecoOrig("hRealZMassRecoOrig","Higgs: Real Z mass (Reconstructed/Original)", 200, 0., 2.);
00232 TH1F hVirtualZMassRecoOrig("hVirtualZMassRecoOrig","Higgs: Virtual Z mass (Reconstructed/Original)", 200, 0., 2.);
00233 TH1F hHiggsMassRecoOrig("hHiggsMassRecoOrig","Higgs: Higgs mass (Reconstructed/Original)", 200, 0., 2.);
00234
00235 // Electrons from Higgs
00236
00237 TH1F hElectronClassesHiggs("hElectronClassesHiggs","Electron Candidates from Higgs - Classification",23, -2., 21.);
00238
00239 TH1F hElectronClass0Higgs("hElectronClass0Higgs","Electron Candidates from Higgs Class: 0, Barrel-Golden",100, 0., 1.5);
00240 TH1F hElectronClass1Higgs("hElectronClass1Higgs","Electron Candidates from Higgs Class: 1, Barrel-Big Brem",100, 0., 1.5);
00241 TH1F hElectronClass2Higgs("hElectronClass2Higgs","Electron Candidates from Higgs Class: 2, Barrel-Narrow",100, 0., 1.5);
00242 TH1F hElectronClass3Higgs("hElectronClass3Higgs","Electron Candidates from Higgs Class: 3, Barrel-Showering",100, 0., 1.5);
00243 TH1F hElectronClass4Higgs("hElectronClass4Higgs","Electron Candidates from Higgs Class: 4, Barrel-Cracks",100, 0., 1.5);
00244 TH1F hElectronClass10Higgs("hElectronClass10Higgs","Electron Candidates from Higgs Class: 0, Endcap-Golden",100, 0., 1.5);
00245 TH1F hElectronClass11Higgs("hElectronClass11Higgs","Electron Candidates from Higgs Class: 1, Endcap-Big Brem",100, 0., 1.5);
00246 TH1F hElectronClass12Higgs("hElectronClass12Higgs","Electron Candidates from Higgs Class: 2, Endcap-Narrow",100, 0., 1.5);
00247 TH1F hElectronClass13Higgs("hElectronClass13Higgs","Electron Candidates from Higgs Class: 3, Endcap-Showering",100, 0., 1.5);
00248
00249 TProfile hElectronClassProfile0Higgs("hElectronClassProfile0Higgs","Electron Candidates from Higgs ClassProfile: 0, Barrel-Golden",100, 0., 2.7);
00250 TProfile hElectronClassProfile1Higgs("hElectronClassProfile1Higgs","Electron Candidates from Higgs ClassProfile: 1, Barrel-Big Brem",100, 0., 2.7);
00251 TProfile hElectronClassProfile2Higgs("hElectronClassProfile2Higgs","Electron Candidates from Higgs ClassProfile: 2, Barrel-Narrow",100, 0., 2.7);
00252 TProfile hElectronClassProfile3Higgs("hElectronClassProfile3Higgs","Electron Candidates from Higgs ClassProfile: 3, Barrel-Showering",100, 0., 2.7);
00253 TProfile hElectronClassProfile4Higgs("hElectronClassProfile4Higgs","Electron Candidates from Higgs ClassProfile: 4, Barrel-Cracks",100, 0., 2.7);
00254 TProfile hElectronClassProfile10Higgs("hElectronClassProfile10Higgs","Electron Candidates from Higgs ClassProfile: 0, Endcap-Golden",100, 0., 2.7);
00255 TProfile hElectronClassProfile11Higgs("hElectronClassProfile11Higgs","Electron Candidates from Higgs ClassProfile: 1, Endcap-Big Brem",100, 0., 2.7);
00256 TProfile hElectronClassProfile12Higgs("hElectronClassProfile12Higgs","Electron Candidates from Higgs ClassProfile: 2, Endcap-Narrow",100, 0., 2.7);
00257 TProfile hElectronClassProfile13Higgs("hElectronClassProfile13Higgs","Electron Candidates from Higgs ClassProfile: 3, Endcap-Showering",100, 0., 2.7);
00258
00259 TH1F hElectronEtaHiggs("hElectronEtaHiggs","Electron Candidates from Higgs Eta",100, 0., 2.7);
00260 TH1F hElectronEtaClass0Higgs("hElectronEtaClass0Higgs","Electron Candidates from Higgs Eta Class: 0, Barrel-Golden",100, 0., 2.7);
00261 TH1F hElectronEtaClass1Higgs("hElectronEtaClass1Higgs","Electron Candidates from Higgs Eta Class: 1, Barrel-Big Brem",100, 0., 2.7);
00262 TH1F hElectronEtaClass2Higgs("hElectronEtaClass2Higgs","Electron Candidates from Higgs Eta Class: 2, Barrel-Narrow",100, 0., 2.7);
00263 TH1F hElectronEtaClass3Higgs("hElectronEtaClass3Higgs","Electron Candidates from Higgs Eta Class: 3, Barrel-Showering",100, 0., 2.7);
00264 TH1F hElectronEtaClass4Higgs("hElectronEtaClass4Higgs","Electron Candidates from Higgs Eta Class: 4, Barrel-Cracks",100, 0., 2.7);
00265 TH1F hElectronEtaClass10Higgs("hElectronEtaClass10Higgs","Electron Candidates from Higgs Eta Class: 0, Endcap-Golden",100, 0., 2.7);
00266 TH1F hElectronEtaClass11Higgs("hElectronEtaClass11Higgs","Electron Candidates from Higgs Eta Class: 1, Endcap-Big Brem",100, 0., 2.7);
00267 TH1F hElectronEtaClass12Higgs("hElectronEtaClass12Higgs","Electron Candidates from Higgs Eta Class: 2, Endcap-Narrow",100, 0., 2.7);
00268 TH1F hElectronEtaClass13Higgs("hElectronEtaClass13Higgs","Electron Candidates from Higgs Eta Class: 3, Endcap-Showering",100, 0., 2.7);
00269
00270 TH1F hElectronMomentumEOverPHiggs("hElectronMomentumEOverPHiggs","Electron Candidates from Higgs - E/p",200,0.,2.);
00271 TH1F hElectronMomentumEnergyHiggs("hElectronMomentumEnergyHiggs","Electron Candidates from Higgs - E/Etrue",200,0.,2.);
00272 TH1F hElectronMomentumTrackHiggs("hElectronMomentumTrackHiggs","Electron Candidates from Higgs - p/Etrue",200,0.,2.);
00273 TH1F hElectronMomentumFinalHiggs("hElectronMomentumFinalHiggs","Electron Candidates from Higgs - Efinal/Etrue",200,0.,2.);
00274
00275 TH2F hElectronMomentumEnergy2Higgs("hElectronMomentumEnergy2Higgs","Electron Candidates from Higgs - E/Etrue vs E/p",200,0.,2.,200,0.,2.);
00276 TH2F hElectronMomentumTrack2Higgs("hElectronMomentumTrack2Higgs","Electron Candidates from Higgs - p/Etrue vs E/p",200,0.,2.,200,0.,2.5);
00277 TH2F hElectronMomentumFinal2Higgs("hElectronMomentumFinal2Higgs","Electron Candidates from Higgs - Efinal/Etrue vs E/p",200,0.,2.,200,0.,2.);
00278
00279 // Initialisation of some variables (mainly counters and temporary)
00280 Int_t nb = 0;
00281 Int_t nevent = 0, ievent = 0;
00282 int ne = 0, nTrigAccL1=0, nTrigAccHLT=0, nTrigAccHLTelec=0,nPresAcc = 0, nIsoAcc = 0, nAmbBef = 0, nAmbAft = 0, n2e2pAcc = 0;
00283 int nCutsIso = 0, nCutsPt = 0, nCutsRealZMass = 0, nCutsVirtualZMass = 0, nCutsHiggsMass = 0;
00284 int mc4e = 0, mc4t = 0, mc2e2t = 0, mc4eCuts = 0, mc4tCuts = 0, mc2e2tCuts = 0;
00285 int nconsidered=0,niso=0;
00286 TObjArray *electrons = 0, *superClusters = 0;
00287 TObjArray *tracks=0;
00288 HepLorentzVector realZ, virtualZ, higgs, realZOrig, virtualZOrig, higgsOrig;
00289 HepLorentzVector tmpElRec, tmpElGen, tmpElGenOrigB, tmpElGenOrigA;
00290 XANAElectronCandidate *tmpElectronCandidate;
00291
00292 // define XANA event
00293 XANAEsdEvent *event = new XANAEsdEvent();
00294 XANAGeneratorEvent *generatorEvent = new XANAGeneratorEvent();
00295 XANAGeantEvent *geantEvent = new XANAGeantEvent();
00296
00297 // version TChain
00298 Short_t numberOfFiles = 0;
00299 TChain *T=new TChain("ESD");
00300 char tmp[512];
00301
00302 // files to load, dependant on the key word given at runtime
00303 // to be completed by other datasets
00304 if (strcmp(argv[1], "hzz150") == 0) {
00305 cout<<"hzz4e 150"<<endl;
00306 numberOfFiles = 18;
00307 for (Int_t m_fi = 0; m_fi < numberOfFiles ; m_fi++) {
00308 sprintf(tmp,"/sps/cms/baffioni/XanadooTrees/hg03_hzz_4e_150_XANA_153_pix/hg03_hzz_4e_150.root.%03d",m_fi);
00309 TFile F(tmp);
00310 if (!(F.GetNkeys())) { continue; }
00311 T->Add(tmp);
00312 }
00313 }
00314 else if (strcmp(argv[1], "hzz120") == 0) {
00315 cout<<"hzz4e 120"<<endl;
00316 numberOfFiles = 20;
00317 for (Int_t m_fi = 0; m_fi < numberOfFiles ; m_fi++) {
00318 sprintf(tmp,"/sps/cms/baffioni/XanadooTrees/hg03_hzz_4e_120_XANA_153_pix/hg03_hzz_4e_120.root.%03d",m_fi);
00319 TFile F(tmp);
00320 if (!(F.GetNkeys())) { continue; }
00321 T->Add(tmp);
00322 }
00323 }
00324 else if(strcmp(argv[1], "tt4e") == 0){
00325 cout<<"ttbar->4e"<<endl;
00326 numberOfFiles = 160;
00327 for (Int_t m_fi = 0; m_fi < numberOfFiles ; m_fi++) {
00328 if (m_fi == 59 || m_fi == 106 ||m_fi == 158 || m_fi == 159 || m_fi == 149 || m_fi == 19|| m_fi == 31) {continue;}
00329 sprintf(tmp,"/sps/cms/baffioni/XanadooTrees/eg03_tt_4e_XANA_153_pix/eg03_tt_4e.root.%03d",m_fi);
00330 TFile F(tmp);
00331 if (!(F.GetNkeys())) {continue; }
00332 T->Add(tmp);
00333 }
00334 }
00335 else{
00336 cout << "no such data yet"<<endl;
00337 return 0;
00338 }
00339
00340 // init settings
00341 initSettings(T, &event, &generatorEvent, &geantEvent);
00342
00343 // end of initialization phase
00344
00345 //************************************************************************************************************
00346
00347 // Start main loop on all events
00348 nevent = (Int_t)T->GetEntries();
00349 cout << "================== Total number of events: " << nevent << " ================" << endl;
00350
00351 // for (Int_t i=0; i<T->GetEntries(); i++) {
00352 for (Int_t i=0; i<100; i++) {
00353 nb += T->GetEntry(i);
00354 ievent++;
00355
00356 if (i%100 == 0) cout << "================ Event: " << i << " ================" << endl;
00357 //if (i%1 == 0) cout << "================ Event: " << i << " ================" << endl;
00358
00359 ne = event->getNumberOfElectronCandidates(); electrons = event->getElectronCandidates();
00360 tracks= event->getTracks();
00361 hElectronNumber.Fill((Axis_t) ne);
00362
00363 /*
00364 // cross check
00365 for (int i3 = 0; i3 < ne; i3++) {
00366 tmpElectronCandidate = (XANAElectronCandidate *) electrons->At(i3);
00367 cout << "Super Cluster Pointer = " << tmpElectronCandidate->getSeedCluster() << endl;
00368 }
00369 */
00370
00371 //**********************************************************************************************************
00372
00373 // FIXME: problem with Super Clusters in endcap
00374 // this part of code should do nothing with data from XANADOO_1_5_3
00375 // to be deleted as soon as we don't use old data anymore
00376
00377 //cout<<"fix superclusters"<<endl;
00378 superClusters = event->getSuperClusters();
00379
00380 for (int i3 = 0; i3 < ne; i3++) {
00381 tmpElectronCandidate = (XANAElectronCandidate *) electrons->At(i3);
00382
00383 // cout << "Super Cluster Pointer = " << tmpElectronCandidate->getSuperCluster() << endl;
00384 for (int isc = 0 ; isc < superClusters->GetLast() + 1; isc ++) {
00385 XANASuperCluster * sc = (XANASuperCluster *) superClusters->At(isc);
00386 if (tmpElectronCandidate->getSuperCluster()==0 &&
00387 tmpElectronCandidate->getSuperClusterPosition() == sc->getPosition()) {
00388 cout<<" !!! entering the super cluster bug correction "<<endl;
00389 cout<<" !!! should not be the case if data >= XANA1_5_3 are used"<<endl;
00390 // cout << " SC: Cand en, SC en, SS pos, brems: " << tmpElectronCandidate->getSuperClusterEnergy() << ", " <<
00391 // sc->getEnergy() << ", " << sc->getPosition().eta() << " "
00392 // << sc->getNumberOfBrems() << endl;
00393 float oldEnergy = tmpElectronCandidate->getSuperClusterEnergy();
00394 tmpElectronCandidate->setSuperClusterEnergy(sc->getEnergy());
00395 //cout << "new supercluster endcap: old energy " << oldEnergy << " new energy " << sc->getEnergy() << endl;
00396 tmpElectronCandidate->setSuperCluster(sc);
00397 //CC also update MomentumAtVertex
00398 Hep3Vector ptrack;
00399 if (tmpElectronCandidate->getElectronTrack()->getAlgoName() == "GSF") {
00400 XANAElectronGSFTrack *gsfetr = (XANAElectronGSFTrack *)(tmpElectronCandidate->getElectronTrack());
00401 ptrack = gsfetr->getGSFMomentumAtVertexMode();
00402 } else {
00403 ptrack = tmpElectronCandidate->getElectronTrack()->getMomentumAtVertex();
00404 }
00405 HepLorentzVector momentum(ptrack.x(), ptrack.y(), ptrack.z(), 0.);
00406 momentum *= sc->getEnergy() / ptrack.mag();
00407 momentum.setE(sc->getEnergy());
00408 //cout << "new electron endcap: old momentum " << tmpElectronCandidate->getMomentumAtVertex() << " new momentum " << momentum << endl;
00409 tmpElectronCandidate->setMomentumAtVertex(momentum);
00410 //CC also update Eoverp
00411 //cout << "new electron endcap: old E/p " << tmpElectronCandidate->getESuperClusterOverP() << " new E/p " << sc->getEnergy() / ptrack.mag() << endl;
00412 tmpElectronCandidate->setMomentumAtVertex(momentum);
00413 tmpElectronCandidate->setESuperClusterOverP(sc->getEnergy()/ptrack.mag());
00414 //CC also update had/E
00415 //cout << "new electron endcap: old h/e " << tmpElectronCandidate->getHadronicOverEm() << " new h/e " << sc->getHadronicOverEm()*oldEnergy/sc->getEnergy() << endl;
00416 tmpElectronCandidate->setHadronicOverEm(sc->getHadronicOverEm()*oldEnergy/sc->getEnergy());
00417 }
00418 }
00419 }
00420
00421
00422 //****************************************************************************************
00423 // part of code to unescale electron candidates
00424 // CC in case we do want to remove all corrections
00425
00426 if (!egammaEscaleCorrection) {
00427 //cout<< "unEscale electron candidates"<<endl;
00428 for (int i3 = 0; i3 < ne; i3++) {
00429 tmpElectronCandidate = (XANAElectronCandidate *) electrons->At(i3);
00430 // from federico
00431 float newEnergy = tmpElectronCandidate->getSuperCluster()->getEnergy()/
00432 tmpElectronCandidate->getSuperCluster()->getEnergyScaleFactor();
00433 HepLorentzVector oldMomentum = tmpElectronCandidate->getMomentumAtVertex();
00434 HepLorentzVector newMomentum(oldMomentum.x()*newEnergy/oldMomentum.t(),
00435 oldMomentum.y()*newEnergy/oldMomentum.t(),
00436 oldMomentum.z()*newEnergy/oldMomentum.t(),
00437 newEnergy);
00438 tmpElectronCandidate->setMomentumAtVertex(newMomentum);
00439 }
00440 }
00441
00442
00443
00444 //**********************************************************************************
00445
00446 // select only 4e decay for signal
00447 // -------------------------------------------------------
00448
00449 //cout<<"XANASelectMc"<<endl;
00450
00451 bool issignal=false;
00452 XANASelectMcHZZ selectorMc(generatorEvent);
00453
00454 //only for signal be careful to add a new signal if there's one
00455 if (strcmp(argv[1], "hzz150") == 0 || strcmp(argv[1], "hzz120") == 0) {
00456 issignal=true;
00457 // Analysing H-->ZZ decay tree
00458 if (selectorMc.isFourElectrons()) mc4e++;
00459 if (selectorMc.isFourTaus()) mc4t++;
00460 if (selectorMc.isTwoElectronsTwoTaus()) mc2e2t++;
00461
00462 if (!selectorMc.isFourElectrons()) {clean(event, generatorEvent, geantEvent);continue;}
00463 }
00464 nconsidered++;
00465
00466 //**********************************************************************************
00467
00469 // Trigger //
00471
00472 // cout<<"XANAnalysisTrigger"<<endl;
00473 XANATriggerInfo triggerInfo = event->getTriggerInfo();
00474 XANAnalysisTrigger trigger(&triggerInfo);
00475
00476 if (trigger.getL1Decision() == true) nTrigAccL1++;
00477 else {clean(event, generatorEvent, geantEvent);continue; }
00478
00479 if (trigger.getHLTDecision() == true) nTrigAccHLT++;
00480 else {clean(event, generatorEvent, geantEvent);continue; }
00481
00482 if (trigger.getHLTElecDecision() == true) nTrigAccHLTelec++;
00483 else {clean(event, generatorEvent, geantEvent);continue; }
00484
00485 //**********************************************************************************
00486
00488 // Preselection //
00490
00491 // 1. Simple Cuts
00492 // --------------
00493
00494 // Ask for 2 electrons and 2 positrons all with pT > 3 GeV with no isolation requirement
00495 //cout<<"XANAEventSelec"<<endl;
00496 XANAEventSelect eventSelect(event);
00497 vector<float> ptCut;
00498 ptCut.push_back(3.);
00499 ptCut.push_back(3.);
00500 if ((eventSelect.NumberOfElCandSelect( 2, -1, 0, ptCut)) &&
00501 (eventSelect.NumberOfElCandSelect( 2, 1, 0, ptCut))) {
00502 nPresAcc++;
00503 } else {clean(event, generatorEvent, geantEvent);continue;}
00504
00505
00506 // 2. Loose isolation
00507 // ------------------
00508
00509 // Calculation isolation on all Electron candidate in the event
00510 //cout<<"XANAEvtIsolation"<<endl;
00511 XANAEvtIsolation isolation(event);
00512 std::vector<XANAIsolation> isolationVector = isolation.getIsolationVector();
00513 std::vector<XANAIsolation>::const_iterator isoIter;
00514
00515 // Fill some histograms
00516 for (isoIter = isolationVector.begin(); isoIter != isolationVector.end(); isoIter++){
00517 hIsoNumberTracks.Fill((Axis_t) isoIter->getNumberTracks());
00518 hIsoPtTracks.Fill((Axis_t) isoIter->getPtTracks());
00519 hIsoNumberECinsideCone.Fill((Axis_t) isoIter->getNumberECinsideCone());
00520 hIsoPtMax.Fill((Axis_t) isoIter->getPtMax());
00521 }
00522
00523
00524 // Very loose isolation (CASE 1): Ask for four isolated electrons with pTsum < 10 GeV inside a cone
00525 if (isolation.getIsolationPtSum(10.,4)) nIsoAcc++;
00526 else {clean(event, generatorEvent, geantEvent);continue;}
00527
00528 // Very loose isolation (CASE 2): Ask for four isolated electrons with pTmax < 5 GeV inside a cone
00529 //if (isolation.getIsolationPtMax(5.,4)) nIsoAcc++;
00530
00531
00532
00533 //**********************************************************************************
00534
00536 // Electron reconstruction //
00538
00539 // Associate reconstructed and generated electrons and make some analysis
00540 //cout<<"XANAElectronRecGen"<<endl;
00541 XANAElectronRecGen anaMcEvent(event,generatorEvent);
00542 map<XANAElectronCandidate *, XANAGeneratorParticle *> electronsRecGen = anaMcEvent.getElectronsRecGen();
00543 map<XANAElectronCandidate *, XANAGeneratorParticle *>::const_iterator iterEleRecGen;
00544
00545 for (iterEleRecGen = electronsRecGen.begin(); iterEleRecGen != electronsRecGen.end(); iterEleRecGen++){
00546 tmpElRec = (iterEleRecGen->first)->getMomentumAtVertex();
00547 tmpElGen = (iterEleRecGen->second)->getMomentum();
00548 hElectronRecGenPt.Fill(tmpElRec.perp()/tmpElGen.perp());
00549 hElectronRecGenDr.Fill(tmpElRec.deltaR(tmpElGen));
00550 }
00551
00552
00553 //**********************************************************************************
00554
00556 // Electron preselection //
00558
00559 //cout<<"XANAElectronSelect"<<endl;
00560 XANAElectronSelect electronSelector(electrons);
00561 electronSelector.ElSelect_dphi(0.1);
00562 TObjArray *selectedElectrons = electronSelector.getSelectedElectrons();
00563
00564
00565
00567 // Ambiguity resolving //
00569
00570 // single resolving
00571 //cout<<"XANAmbiguityResolve"<<endl;
00572 XANAmbiguityResolve ambResolve(selectedElectrons,selectedElectrons->GetLast() + 1);
00573 ambResolve.ResolveByEoverP();
00574
00575 // double resolving
00576
00577 //XANAmbiguityResolve ambResolve1(selectedElectrons,selectedElectrons->GetLast() + 1);
00578 //ambResolve1.ResolveByEoverP();
00579
00580 // XANAmbiguityResolve ambResolve(ambResolve1.getResolvedElcands(),ambResolve1.getResolvedElcands()->GetLast() + 1);
00581 // ambResolve.ResolveByEoverP();
00582
00583 // Prepare for calculation of average number of Electron Candidates before and after ambiguity resolving
00584 nAmbBef+=ne;
00585 nAmbAft+=ambResolve.getResolvedElcands()->GetLast()+1;
00586 // Fill some histograms
00587 hElectronNumberAfterAmb.Fill((Axis_t) ambResolve.getResolvedElcands()->GetLast() + 1);
00588 hElectronNumberBeforeAfter.Fill((Axis_t) ne - ambResolve.getResolvedElcands()->GetLast() - 1);
00589
00590
00591 //**********************************************************************************
00592
00594 // Electron classification //
00596
00597 //cout<<" XANAElectronClassification"<<endl;
00598 XANAElectronClassification electronClasses(ambResolve.getResolvedElcands());
00599 map<XANAElectronCandidate *, int> elClassMap = electronClasses.getElectronClass();
00600 map<XANAElectronCandidate *, int>::const_iterator elClassMapIter;
00601 float uncorrEnergy = 0;
00602
00603 // some checks
00604 // cout << "------- Electron classification ------- " << ambResolve.getResolvedElcands()->GetLast() + 1 << endl;
00605 for (elClassMapIter = elClassMap.begin(); elClassMapIter != elClassMap.end(); elClassMapIter++) {
00606 // cout << "El. candidate pointer and class : " << elClassMapIter->first << " " << elClassMapIter->second << endl;
00607 hElectronClasses.Fill((Axis_t) elClassMapIter->second);
00608 tmpElectronCandidate = elClassMapIter->first;
00609 float myEta = tmpElectronCandidate->getSuperClusterPosition().eta();
00610 uncorrEnergy = tmpElectronCandidate->getSuperClusterEnergy();
00611 if (elClassMapIter->second >= 0 && elClassMapIter->second < 5 && tmpElectronCandidate->getSuperCluster())
00612 uncorrEnergy = uncorrEnergy/tmpElectronCandidate->getSuperCluster()->getEnergyScaleFactor();
00613 tmpElGen = ((electronsRecGen.find(tmpElectronCandidate))->second)->getMomentum();
00614 Axis_t ratio = uncorrEnergy/tmpElGen.e();
00615 hElectronEta.Fill(myEta);
00616 if (elClassMapIter->second == 0) {
00617 hElectronClass0.Fill(ratio); if (ratio > 0.8 && ratio < 2.) hElectronClassProfile0.Fill(myEta, ratio);
00618 hElectronEtaClass0.Fill(myEta);
00619 }
00620 else if (elClassMapIter->second == 1) {
00621 hElectronClass1.Fill(ratio); if (ratio > 0.8 && ratio < 2.) hElectronClassProfile1.Fill(myEta, ratio);
00622 hElectronEtaClass1.Fill(myEta);
00623 }
00624 else if (elClassMapIter->second == 2) {
00625 hElectronClass2.Fill(ratio); if (ratio > 0.8 && ratio < 2.) hElectronClassProfile2.Fill(myEta, ratio);
00626 hElectronEtaClass2.Fill(myEta);
00627 }
00628 else if (elClassMapIter->second == 3) {
00629 hElectronClass3.Fill(ratio); if (ratio > 0.8 && ratio < 2.) hElectronClassProfile3.Fill(myEta, ratio);
00630 hElectronEtaClass3.Fill(myEta);
00631 }
00632 else if (elClassMapIter->second == 4) {
00633 hElectronClass4.Fill(ratio); if (ratio > 0.8 && ratio < 2.) hElectronClassProfile4.Fill(myEta, ratio);
00634 hElectronEtaClass4.Fill(myEta);
00635 }
00636 else if (elClassMapIter->second == 10) {
00637 hElectronClass10.Fill(ratio); if (ratio > 0.8 && ratio < 2.) hElectronClassProfile10.Fill(myEta, ratio);
00638 hElectronEtaClass10.Fill(myEta);
00639 }
00640 else if (elClassMapIter->second == 11) {
00641 hElectronClass11.Fill(ratio); if (ratio > 0.8 && ratio < 2.) hElectronClassProfile11.Fill(myEta, ratio);
00642 hElectronEtaClass11.Fill(myEta);
00643 }
00644 else if (elClassMapIter->second == 12) {
00645 hElectronClass12.Fill(ratio); if (ratio > 0.8 && ratio < 2.) hElectronClassProfile12.Fill(myEta, ratio);
00646 hElectronEtaClass12.Fill(myEta);
00647 }
00648 else if (elClassMapIter->second == 13) {
00649 hElectronClass13.Fill(ratio); if (ratio > 0.8 && ratio < 2.)hElectronClassProfile13.Fill(myEta, ratio);
00650 hElectronEtaClass13.Fill(myEta);
00651 }
00652 }
00653
00654
00655 //**********************************************************************************
00656
00658 // Electron corrections //
00660
00661 float corrEnergy;
00662
00663 if (federicoEscaleCorrection) {
00664
00665 XANAElectronCorrection electronCorrection(elClassMap);
00666 corrEnergy = 0.;
00667
00668 // cross checks
00669
00670 for (elClassMapIter = elClassMap.begin(); elClassMapIter != elClassMap.end(); elClassMapIter++) {
00671 tmpElectronCandidate = elClassMapIter->first;
00672 float myEta = tmpElectronCandidate->getSuperClusterPosition().eta();
00673 tmpElGen = ((electronsRecGen.find(tmpElectronCandidate))->second)->getMomentum();
00674 corrEnergy = tmpElectronCandidate->getSuperClusterEnergy();
00675 Axis_t ratio = corrEnergy/tmpElGen.e();
00676 if (elClassMapIter->second == 0) { hElectronClassCorr0.Fill(ratio); if (ratio > 0.8 && ratio < 2.) hElectronClassProfileCorr0.Fill(myEta, ratio); }
00677 else if (elClassMapIter->second == 1) { hElectronClassCorr1.Fill(ratio); if (ratio > 0.8 && ratio < 2.) hElectronClassProfileCorr1.Fill(myEta, ratio); }
00678 else if (elClassMapIter->second == 2) { hElectronClassCorr2.Fill(ratio); if (ratio > 0.8 && ratio < 2.) hElectronClassProfileCorr2.Fill(myEta, ratio); }
00679 else if (elClassMapIter->second == 3) { hElectronClassCorr3.Fill(ratio); if (ratio > 0.8 && ratio < 2.) hElectronClassProfileCorr3.Fill(myEta, ratio); }
00680 else if (elClassMapIter->second == 4) { hElectronClassCorr4.Fill(ratio); if (ratio > 0.8 && ratio < 2.) hElectronClassProfileCorr4.Fill(myEta, ratio); }
00681 else if (elClassMapIter->second == 10) { hElectronClassCorr10.Fill(ratio); if (ratio > 0.8 && ratio < 2.) hElectronClassProfileCorr10.Fill(myEta, ratio); }
00682 else if (elClassMapIter->second == 11) { hElectronClassCorr11.Fill(ratio); if (ratio > 0.8 && ratio < 2.) hElectronClassProfileCorr11.Fill(myEta, ratio); }
00683 else if (elClassMapIter->second == 12) { hElectronClassCorr12.Fill(ratio); if (ratio > 0.8 && ratio < 2.) hElectronClassProfileCorr12.Fill(myEta, ratio); }
00684 else if (elClassMapIter->second == 13) { hElectronClassCorr13.Fill(ratio); if (ratio > 0.8 && ratio < 2.) hElectronClassProfileCorr13.Fill(myEta, ratio); }
00685 }
00686
00687 }
00688
00689 //**********************************************************************************
00690
00692 // Electron momentum calculation //
00694 if (momentumCorrection) {
00695
00696 XANAElectronMomentum electronMomentum(elClassMap);
00697
00698 // cross checks
00699 for (elClassMapIter = elClassMap.begin(); elClassMapIter != elClassMap.end(); elClassMapIter++) {
00700 tmpElectronCandidate = elClassMapIter->first;
00701 XANAElectronGSFTrack *gsfTrack = (XANAElectronGSFTrack *) tmpElectronCandidate->getElectronTrack();
00702 HepVector3D tmpTrackRec = gsfTrack->getGSFMomentumAtVertexMode();
00703 tmpElGen = ((electronsRecGen.find(tmpElectronCandidate))->second)->getMomentum();
00704 Axis_t eOverP = tmpElectronCandidate->getSuperClusterEnergy()/tmpTrackRec.mag();
00705 Axis_t eOverETrue = tmpElectronCandidate->getSuperClusterEnergy()/tmpElGen.t();
00706 Axis_t pOverETrue = tmpTrackRec.mag()/tmpElGen.t();
00707 HepLorentzVector final = tmpElectronCandidate->getMomentumAtVertex();
00708 Axis_t finalOverETrue = final.t()/tmpElGen.t();
00709
00710 hElectronMomentumEOverP.Fill(eOverP);
00711 hElectronMomentumEnergy.Fill(eOverETrue); hElectronMomentumEnergy2.Fill(eOverP, eOverETrue);
00712 hElectronMomentumTrack.Fill(pOverETrue); hElectronMomentumTrack2.Fill(eOverP, pOverETrue);
00713 hElectronMomentumFinal.Fill(finalOverETrue); hElectronMomentumFinal2.Fill(eOverP, finalOverETrue);
00714 }
00715
00716 }
00717
00718 //********************************************************************
00719
00721 // Higgs Analysis //
00723
00724 // Isolation
00725 // select events with >=2e+ and >=2e- isolated electrons with Pt>5GeV
00726
00727
00728 XANAEventSelect eventSelect2(ambResolve.getResolvedElcands(),tracks);
00729 vector<float> ptCut2;
00730 ptCut2.push_back(5.);
00731 ptCut2.push_back(5.);
00732
00733 if ((eventSelect2.NumberOfElCandSelect( 2, -1, 1, ptCut2)) &&
00734 (eventSelect2.NumberOfElCandSelect( 2, 1, 1, ptCut2))) {
00735 niso++;
00736 } else {clean(event, generatorEvent, geantEvent);continue;}
00737 TObjArray * isoel=eventSelect2.getGoodElCands(0,1,5);
00738
00739
00740 // Reconstruct Z, Z*, Higgs (using unambigous isolated electrons only)
00741 // ----------------------------------------------------------
00742 // XANASelectHZZ selector(electrons,ne);
00743 XANASelectHZZ selector(isoel,isoel->GetLast()+1);
00744
00745 // internal consisteny, asks for 2e- and 2e+
00746 if (selector.isValid()) n2e2pAcc++;
00747 else {clean(event, generatorEvent, geantEvent);continue; }
00748
00749 // Fill temporary array for easier navigation
00750 XANAElectronCandidate *electronsHiggs[4] = {selector.getElectronFromRealZ(),
00751 selector.getPositronFromRealZ(),
00752 selector.getElectronFromVirtualZ(),
00753 selector.getPositronFromVirtualZ()};
00754
00755
00756
00757 realZ = selector.getRealZ(); hRealZMass.Fill(realZ.m());
00758 virtualZ = selector.getVirtualZ(); hVirtualZMass.Fill(virtualZ.m());
00759 higgs = selector.getHiggs(); hHiggsMass.Fill(higgs.m());
00760
00761
00762
00763 // Fill temporary array for easier navigation
00764 XANAGeneratorParticle *electronsHiggsOrig[4];
00765 XANAGeneratorParticle *electronsHiggsOrigAftPh[4];
00766 if(issignal){
00767 // original electrons before PHOTOS
00768
00769 electronsHiggsOrig [0]= selectorMc.getElectronFromRealZ();
00770 electronsHiggsOrig [1]= selectorMc.getPositronFromRealZ();
00771 electronsHiggsOrig [2]= selectorMc.getElectronFromVirtualZ();
00772 electronsHiggsOrig [3]= selectorMc.getPositronFromVirtualZ() ;
00773
00774 // original electrons after PHOTOS
00775 electronsHiggsOrigAftPh[0] = selectorMc.getElectronFromRealZAftPh();
00776 electronsHiggsOrigAftPh[1] = selectorMc.getPositronFromRealZAftPh();
00777 electronsHiggsOrigAftPh[2] = selectorMc.getElectronFromVirtualZAftPh();
00778 electronsHiggsOrigAftPh[3] = selectorMc.getPositronFromVirtualZAftPh() ;
00779
00780 // Some analysis
00781 realZOrig = selectorMc.getRealZ()->getMomentum(); hRealZMassRecoOrig.Fill(realZ.m()/realZOrig.m());
00782 virtualZOrig = selectorMc.getVirtualZ()->getMomentum(); hVirtualZMassRecoOrig.Fill(virtualZ.m()/virtualZOrig.m());
00783 higgsOrig = selectorMc.getHiggs()->getMomentum(); hHiggsMassRecoOrig.Fill(higgs.m()/higgsOrig.m());
00784 }
00785
00786
00787 // Electrons from Higgs - Analysis
00788 // Reconstructed versus generated (for Higgs electrons) and reconstructed versus original electrons
00789 for (int ieh = 0; ieh < 4; ieh++){
00790 hElectronPt.Fill(electronsHiggs[ieh]->getPt());
00791 tmpElRec = electronsHiggs[ieh]->getMomentumAtVertex();
00792 tmpElGen = ((electronsRecGen.find(electronsHiggs[ieh]))->second)->getMomentum();
00793 hElectronRecGenPtHiggs.Fill(tmpElRec.perp()/tmpElGen.perp());
00794 hElectronRecGenDrHiggs.Fill(tmpElRec.deltaR(tmpElGen));
00795
00796 if (issignal && selectorMc.isFourElectrons()) {
00797 tmpElGenOrigB = electronsHiggsOrig[ieh]->getMomentum();
00798 hElectronRecGenPtHiggsOrig.Fill(tmpElRec.perp()/tmpElGenOrigB.perp());
00799 hElectronRecGenDrHiggsOrig.Fill(tmpElRec.deltaR(tmpElGenOrigB));
00800 tmpElGenOrigA = electronsHiggsOrigAftPh[ieh]->getMomentum();
00801
00802 hElectronRecGenPtHiggsOrigAftPh.Fill(tmpElRec.perp()/tmpElGenOrigA.perp());
00803 hElectronRecGenDrHiggsOrigAftPh.Fill(tmpElRec.deltaR(tmpElGenOrigA));
00804 hElectronPhotosAftBefPt.Fill(tmpElGenOrigA.perp()/tmpElGenOrigB.perp());
00805 hElectronPhotosAftBefE.Fill(tmpElGenOrigA.t()/tmpElGenOrigB.t());
00806 // if (tmpElGenOrigA.perp()/tmpElGenOrigB.perp()> 1)
00807 // cout << tmpElGenOrigA.perp()/tmpElGenOrigB.perp() << " | " << tmpElGenOrigA.t()/tmpElGenOrigB.t() << endl;
00808 }
00809
00810
00811 // find to which classes electrons belongs
00812 int elHiggsClass = (elClassMap.find(electronsHiggs[ieh]))->second;
00813 hElectronClassesHiggs.Fill((Axis_t) elHiggsClass);
00814
00815 // use the same analysis as above
00816 tmpElectronCandidate = electronsHiggs[ieh];
00817 XANAElectronGSFTrack *gsfTrack = (XANAElectronGSFTrack *) tmpElectronCandidate->getElectronTrack();
00818 HepVector3D tmpTrackRec = gsfTrack->getGSFMomentumAtVertexMode();
00819 Axis_t eOverP = tmpElectronCandidate->getSuperClusterEnergy()/tmpTrackRec.mag();
00820 Axis_t eOverETrue = tmpElectronCandidate->getSuperClusterEnergy()/tmpElGen.t();
00821 Axis_t pOverETrue = tmpTrackRec.mag()/tmpElGen.t();
00822 HepLorentzVector final = tmpElectronCandidate->getMomentumAtVertex();
00823 //std::cout << "tmpElGen " << tmpElGen << std::endl;
00824 //std::cout << "tmpElRec " << tmpElRec << std::endl;
00825 //std::cout << "final " << final << std::endl;
00826 Axis_t finalOverETrue = final.t()/tmpElGen.t();
00827
00828 hElectronMomentumEOverPHiggs.Fill(eOverP);
00829 hElectronMomentumEnergyHiggs.Fill(eOverETrue); hElectronMomentumEnergy2Higgs.Fill(eOverP, eOverETrue);
00830 hElectronMomentumTrackHiggs.Fill(pOverETrue); hElectronMomentumTrack2Higgs.Fill(eOverP, pOverETrue);
00831 hElectronMomentumFinalHiggs.Fill(finalOverETrue); hElectronMomentumFinal2Higgs.Fill(eOverP, finalOverETrue);
00832
00833 // Classes analysis
00834 int electronClass = (elClassMap.find(electronsHiggs[ieh]))->second;
00835 float myEta = electronsHiggs[ieh]->getSuperClusterPosition().eta();
00836 tmpElGen = ((electronsRecGen.find(electronsHiggs[ieh]))->second)->getMomentum();
00837 corrEnergy = electronsHiggs[ieh]->getSuperClusterEnergy();
00838 Axis_t ratio = corrEnergy/tmpElGen.e();
00839 hElectronEtaHiggs.Fill(myEta);
00840 if (electronClass == 0) {
00841 hElectronClass0Higgs.Fill(ratio); if (ratio > 0.8 && ratio < 2.) hElectronClassProfile0Higgs.Fill(myEta, ratio);
00842 hElectronEtaClass0Higgs.Fill(myEta);
00843 }
00844 else if (electronClass == 1) {
00845 hElectronClass1Higgs.Fill(ratio); if (ratio > 0.8 && ratio < 2.) hElectronClassProfile1Higgs.Fill(myEta, ratio);
00846 hElectronEtaClass1Higgs.Fill(myEta);
00847 }
00848 else if (electronClass == 2) {
00849 hElectronClass2Higgs.Fill(ratio); if (ratio > 0.8 && ratio < 2.) hElectronClassProfile2Higgs.Fill(myEta, ratio);
00850 hElectronEtaClass2Higgs.Fill(myEta);
00851 }
00852 else if (electronClass == 3) {
00853 hElectronClass3Higgs.Fill(ratio); if (ratio > 0.8 && ratio < 2.) hElectronClassProfile3Higgs.Fill(myEta, ratio);
00854 hElectronEtaClass3Higgs.Fill(myEta);
00855 }
00856 else if (electronClass == 4) {
00857 hElectronClass4Higgs.Fill(ratio); if (ratio > 0.8 && ratio < 2.) hElectronClassProfile4Higgs.Fill(myEta, ratio);
00858 hElectronEtaClass4Higgs.Fill(myEta);
00859 }
00860 else if (electronClass == 10) {
00861 hElectronClass10Higgs.Fill(ratio); if (ratio > 0.8 && ratio < 2.) hElectronClassProfile10Higgs.Fill(myEta, ratio);
00862 hElectronEtaClass10Higgs.Fill(myEta);
00863 }
00864 else if (electronClass == 11) {
00865 hElectronClass11Higgs.Fill(ratio); if (ratio > 0.8 && ratio < 2.) hElectronClassProfile11Higgs.Fill(myEta, ratio);
00866 hElectronEtaClass11Higgs.Fill(myEta);
00867 }
00868 else if (electronClass == 12) {
00869 hElectronClass12Higgs.Fill(ratio); if (ratio > 0.8 && ratio < 2.) hElectronClassProfile12Higgs.Fill(myEta, ratio);
00870 hElectronEtaClass12Higgs.Fill(myEta);
00871 }
00872 else if (electronClass == 13) {
00873 hElectronClass13Higgs.Fill(ratio); if (ratio > 0.8 && ratio < 2.) hElectronClassProfile13Higgs.Fill(myEta, ratio);
00874 hElectronEtaClass13Higgs.Fill(myEta);
00875 }
00876 }
00877
00878 //***********************************************************************
00879
00880
00881 //this part is commented because the isolation is now done before for the moment
00882
00883 // // Calculate isolation in the region around four selected electrons
00884 // // -----------------------------------------------------------------
00885 // map<XANAElectronCandidate *,bool> mapSelectedElectrons;
00886 // for (int ieh = 0; ieh < 4; ieh++) mapSelectedElectrons[electronsHiggs[ieh]] = true;
00887 // XANAEvtRegionIsolation regionIsolation(event,mapSelectedElectrons);
00888
00889 // // Fill some histograms
00890 // hIsoRegionNumberTracks.Fill((Axis_t) regionIsolation.getNumberTracks());
00891 // hIsoRegionPtTracks.Fill((Axis_t) regionIsolation.getPtTracks());
00892 // hIsoRegionPtMax.Fill((Axis_t) regionIsolation.getPtMax());
00893
00894
00895
00896 //**********************************************************************
00897
00898 // Cuts
00899 // ----
00900
00901 // 1. Isolation cut
00902 // ................
00903 // already done
00904 // if (regionIsolation.getPtTracks() < 20.) nCutsIso++; // will be replaced by isIsolated() method
00905 // else {clean(event, generatorEvent, geantEvent);continue;}
00906
00907
00908
00909 // 2. pT electrons cut
00910 // ...................
00911 if (selector.isPtCut(5.,10.,15.,20.)) nCutsPt++;
00912 else {clean(event, generatorEvent, geantEvent);continue;}
00913
00914 // 3. Real Z mass cut
00915 // ..................
00916 if (selector.isRealZMassCut(XANASelectHZZ::ZMASSPDG-15.,XANASelectHZZ::ZMASSPDG+15.)) nCutsRealZMass++;
00917 //if (selector.isRealZMassCut(40.,120.)) nCutsRealZMass++;
00918 else {clean(event, generatorEvent, geantEvent);continue; }
00919
00920 // 4. Virtual Z mass cut
00921 // .....................
00922 if (selector.isVirtualZMassCut(10.,60.)) nCutsVirtualZMass++;
00923 //if (selector.isVirtualZMassCut(0.,100.)) nCutsVirtualZMass++;
00924 else {clean(event, generatorEvent, geantEvent);continue; }
00925
00926 // 5. Higgs mass cut
00927 // .....................
00928 if (selector.isHiggsMassCut(140.,160.)) nCutsHiggsMass++;
00929 //if (selector.isHiggsMassCut(80.,180.)) nCutsHiggsMass++;
00930 else {clean(event, generatorEvent, geantEvent);continue; }
00931
00932
00933 // some analysis after all cuts
00934
00935 hRealZMassCuts.Fill(realZ.m());
00936 hVirtualZMassCuts.Fill(virtualZ.m());
00937 hHiggsMassCuts.Fill(higgs.m());
00938
00939 for (int ieh = 0; ieh < 4; ieh++){
00940 hElectronPtCuts.Fill(electronsHiggs[ieh]->getPt());
00941 tmpElRec = electronsHiggs[ieh]->getMomentumAtVertex();
00942 tmpElGen = ((electronsRecGen.find(electronsHiggs[ieh]))->second)->getMomentum();
00943 hElectronRecGenPtHiggsCuts.Fill(tmpElRec.perp()/tmpElGen.perp());
00944 hElectronRecGenDrHiggsCuts.Fill(tmpElRec.deltaR(tmpElGen));
00945 }
00946
00947 // Surviving probability for different classes
00948 if (issignal&&selectorMc.isFourElectrons()) mc4eCuts++;
00949 if (issignal&&selectorMc.isFourTaus()) mc4tCuts++;
00950 if (issignal&&selectorMc.isTwoElectronsTwoTaus()) mc2e2tCuts++;
00951
00952
00953 event->clear();
00954 generatorEvent->clear();
00955 geantEvent->clear();
00956
00957 }
00958
00959 // cout << "XANASelectMcHZZ: " << " 4E = " << 100.*mc4e/nevent << "% | 4T = " <<
00960 // 100.*mc4t/nevent << "% | 2E2T = " << 100.*mc2e2t/nevent << "% | OTHER(2E2M?) = " <<
00961 // 100.*(nevent-mc4e-mc4t-mc2e2t)/nevent << "%" << endl << endl;
00962 cout << " -------------------------------------------------------------------------------" << endl;
00963 cout << " ACCEPTANCES " << endl;
00964 cout << " TOTAL NUMBER OF EVENTS = " << nevent << endl;
00965 cout << " TOTAL NUMBER OF EVENTS read = " << ievent << endl;
00966 cout << " NUMBER of 4e events (signal) = " <<nconsidered<<endl;
00967 cout << " -------------------------------------------------------------------------------" << endl;
00968 cout << " #EVENTS RELATIVE ACC. (%) APSOLUTE ACC. (%) " << endl;
00969 cout << " -------------------------------------------------------------------------------" << endl;
00970 cout << " 1. TRIGGER (L1) = \t" << nTrigAccL1 << "\t\t" << 100.*nTrigAccL1/nconsidered
00971 << "\t\t\t" << 100.*nTrigAccL1/nconsidered << endl;
00972 cout << " 1. TRIGGER (HLT) = \t" << nTrigAccHLT << "\t\t" << 100.*nTrigAccHLT/nTrigAccL1
00973 << "\t\t\t" << 100.*nTrigAccHLT/nconsidered << endl;
00974 cout << " 1. TRIGGER (HLT elec) = \t" << nTrigAccHLTelec << "\t\t" << 100.*nTrigAccHLTelec/nTrigAccHLT
00975 << "\t\t\t" << 100.*nTrigAccHLTelec/nconsidered << endl;
00976 cout << " 2. SIMPLE CUTS = \t" << nPresAcc << "\t\t" << 100.*nPresAcc/nTrigAccHLTelec
00977 << "\t\t\t" << 100.*nPresAcc/nconsidered << endl;
00978 cout << " 3. LOOSE ISOLATION = \t" << nIsoAcc << "\t\t" << 100.*nIsoAcc/nPresAcc
00979 << "\t\t\t" << 100.*nIsoAcc/nconsidered << endl;
00980 cout << " 3. TIGH ISOLATION = \t" << niso << "\t\t" << 100.*niso/nIsoAcc
00981 << "\t\t\t" << 100.*niso/nconsidered << endl;
00982
00983
00984 cout << " 4. HZZ RECO = \t" << n2e2pAcc << "\t\t" << 100.*n2e2pAcc/niso
00985 << "\t\t\t" << 100.*n2e2pAcc/nconsidered << endl;
00986 // cout << " 5. ISOLATION = \t" << nCutsIso << "\t\t" << 100.*nCutsIso/n2e2pAcc
00987 // << "\t\t\t" << 100.*nCutsIso/nconsidered << endl;
00988 cout << " 6. ELECTRON PT = \t" << nCutsPt << "\t\t" << 100.*nCutsPt/n2e2pAcc
00989 << "\t\t\t" << 100.*nCutsPt/nconsidered << endl;
00990 cout << " 7. REAL Z MASS = \t" << nCutsRealZMass << "\t\t" << 100.*nCutsRealZMass/nCutsPt
00991 << "\t\t\t" << 100.*nCutsRealZMass/nconsidered << endl;
00992 cout << " 8. VIRTUAL Z MASS = \t" << nCutsVirtualZMass << "\t\t" << 100.*nCutsVirtualZMass/nCutsRealZMass
00993 << "\t\t\t" << 100.*nCutsVirtualZMass/nconsidered << endl;
00994 cout << " 9. HIGGS MASS = \t" << nCutsHiggsMass << "\t\t" << 100.*nCutsHiggsMass/nCutsVirtualZMass
00995 << "\t\t\t" << 100.*nCutsHiggsMass/nconsidered << endl;
00996 cout << " -------------------------------------------------------------------------------" << endl;
00997 cout << " TOTAL ACCEPTANCE = " << 100.*nCutsHiggsMass/nconsidered << "%" << endl;
00998 cout << " -------------------------------------------------------------------------------" << endl;
00999
01000 cout << nevent << " events and " << nb << " bytes read " << endl;
01001 cout << endl;
01002
01003 // Exporting histograms to file
01004
01005 f.Write();
01006 f.Close();
01007
01008 event->Delete();
01009 generatorEvent->Delete();
01010 geantEvent->Delete();
01011 stopWatch.Stop();
01012 cout << "Time elapsed - CPU: " << stopWatch.CpuTime() <<
01013 " real: " << stopWatch.RealTime() << endl;
01014
01015 }
|
1.3.5