CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
AnalysisErsatz.cc
Go to the documentation of this file.
1 /*#include "ElectroWeakAnalysis/ZEE/interface/AnalysisErsatz.h"
2 AnalysisErsatz::AnalysisErsatz(const edm::ParameterSet& ps)
3 {
4  MCTruthCollection_ = ps.getParameter<edm::InputTag>("MCTruthCollection");
5  ElectronCollection_ = ps.getParameter<edm::InputTag>("ElectronCollection");
6  GenMEtCollection_ = ps.getParameter<edm::InputTag>("GenMEtCollection");
7  //T1MEtCollection_ = ps.getParameter<edm::InputTag>("T1MEtCollection");
8  PfMEtCollection_ = ps.getParameter<edm::InputTag>("PfMEtCollection");
9  TcMEtCollection_ = ps.getParameter<edm::InputTag>("TcMEtCollection");
10  CaloMEtCollection_ = ps.getParameter<edm::InputTag>("CaloMEtCollection");
11  TriggerEvent_ = ps.getParameter<edm::InputTag>("TriggerEvent");
12  TriggerPath_ = ps.getParameter<edm::InputTag>("TriggerPath");
13  TriggerResults_ = ps.getParameter<edm::InputTag>("TriggerResults");
14  TriggerName_ = ps.getParameter<std::string>("TriggerName");
15  ErsatzEvent_ = ps.getParameter<bool>("ErsatzEvent");
16  double CElecPtMin = ps.getParameter<double>("CElecPtMin");
17  double CEB_siEiE = ps.getParameter<double>("CEB_sigmaIEtaIEta");
18  double CEB_dPhiIn = ps.getParameter<double>("CEB_deltaPhiIn");
19  double CEB_dEtaIn = ps.getParameter<double>("CEB_deltaEtaIn");
20  double CEB_EcalIso = ps.getParameter<double>("CEB_EcalIso");
21  double CEB_HcalIso = ps.getParameter<double>("CEB_HcalIso");
22  double CEB_TrckIso = ps.getParameter<double>("CEB_TrckIso");
23  double CEE_siEiE = ps.getParameter<double>("CEE_sigmaIEtaIEta");
24  double CEE_dPhiIn = ps.getParameter<double>("CEE_deltaPhiIn");
25  double CEE_dEtaIn = ps.getParameter<double>("CEE_deltaEtaIn");
26  double CEE_EcalIso = ps.getParameter<double>("CEE_EcalIso");
27  double CEE_HcalIso = ps.getParameter<double>("CEE_HcalIso");
28  double CEE_TrckIso = ps.getParameter<double>("CEE_TrckIso");
29  mW_ = ps.getParameter<double>("mW");
30  mZ_ = ps.getParameter<double>("mZ");
31 
32  CutVector_.resize(13);
33  CutVector_[EtCut_] = CElecPtMin;
34  CutVector_[EB_sIhIh_] = CEB_siEiE;
35  CutVector_[EB_dPhiIn_] = CEB_dPhiIn;
36  CutVector_[EB_dEtaIn_] = CEB_dEtaIn;
37  CutVector_[EB_TrckIso_] = CEB_TrckIso;
38  CutVector_[EB_EcalIso_] = CEB_EcalIso;
39  CutVector_[EB_HcalIso_] = CEB_HcalIso;
40  CutVector_[EE_sIhIh_] = CEE_siEiE;
41  CutVector_[EE_dPhiIn_] = CEE_dPhiIn;
42  CutVector_[EE_dEtaIn_] = CEE_dEtaIn;
43  CutVector_[EE_TrckIso_] = CEE_TrckIso;
44  CutVector_[EE_EcalIso_] = CEE_EcalIso;
45  CutVector_[EE_HcalIso_] = CEE_HcalIso;
46 
47 }
48 
49 
50 AnalysisErsatz::~AnalysisErsatz()
51 {
52 }
53 
54 
55 // ------------ method called once each job just before starting event loop ------------
56 void AnalysisErsatz::beginJob()
57 {
58  edm::Service<TFileService> fs;
59  t_ = fs->make<TTree>("AnalysisData", "Analysis of Ersatz MEt Properties");
60 
61  t_->Branch("Boson_pt", &Boson_pt_,"Boson_pt/D");
62  t_->Branch("Boson_y", &Boson_y_,"Boson_y/D");
63  t_->Branch("Boson_phi", &Boson_phi_,"Boson_phi/D");
64  t_->Branch("Boson_m", &Boson_m_,"Boson_m/D");
65  t_->Branch("McElec1_pt", McElec1_pt_,"McElec1_pt[4]/D");
66  t_->Branch("McElec1_eta", McElec1_eta_,"McElec1_eta[4]/D");
67  t_->Branch("McElec3_pt", McElec3_pt_,"McElec3_pt[4]/D");
68  t_->Branch("McElec3_eta", McElec3_eta_,"McElec3_eta[4]/D");
69  t_->Branch("RndmInt", &RndmInt_, "RndmInt/I");
70  t_->Branch("RndmTrig", &RndmTrig_, "RndmTrig/I");
71  t_->Branch("RndmReco", &RndmReco_, "RndmReco/I");
72  t_->Branch("OthrTrig", &OthrTrig_, "OthrTrig/I");
73  t_->Branch("OthrReco", &OthrReco_, "OthrReco/I");
74  t_->Branch("RndmMcElec_pt", &RndmMcElec_pt_,"RndmMcElec_pt/D");
75  t_->Branch("RndmMcElec_eta", &RndmMcElec_eta_,"RndmMcElec_eta/D");
76  t_->Branch("RndmMcElec_phi", &RndmMcElec_phi_,"RndmMcElec_phi/D");
77  t_->Branch("RndmMcElec_Rescaled_pt", &RndmMcElec_Rescaled_pt_,"RndmMcElec_Rescaled_pt/D");
78  t_->Branch("RndmMcElec_Rescaled_eta", &RndmMcElec_Rescaled_eta_,"RndmMcElec_Rescaled_eta/D");
79  t_->Branch("RndmMcElec_Rescaled_phi", &RndmMcElec_Rescaled_phi_,"RndmMcElec_Rescaled_phi/D");
80  t_->Branch("RndmMcElecTRIG_pt", &RndmMcElecTRIG_pt_,"RndmMcElecTRIG_pt/D");
81  t_->Branch("RndmMcElecTRIG_eta", &RndmMcElecTRIG_eta_,"RndmMcElecTRIG_eta/D");
82  t_->Branch("RndmMcElecRECO_pt", &RndmMcElecRECO_pt_,"RndmMcElecRECO_pt/D");
83  t_->Branch("RndmMcElecRECO_eta", &RndmMcElecRECO_eta_,"RndmMcElecRECO_eta/D");
84  t_->Branch("OthrMcElec_pt", &OthrMcElec_pt_,"OthrMcElec_pt/D");
85  t_->Branch("OthrMcElec_eta", &OthrMcElec_eta_,"OthrMcElec_eta/D");
86  t_->Branch("OthrMcElec_phi", &OthrMcElec_phi_,"OthrMcElec_phi/D");
87  t_->Branch("OthrMcElec_Rescaled_pt", &OthrMcElec_Rescaled_pt_,"OthrMcElec_Rescaled_pt/D");
88  t_->Branch("OthrMcElec_Rescaled_eta", &OthrMcElec_Rescaled_eta_,"OthrMcElec_Rescaled_eta/D");
89  t_->Branch("OthrMcElec_Rescaled_phi", &OthrMcElec_Rescaled_phi_,"OthrMcElec_Rescaled_phi/D");
90  t_->Branch("OthrMcElecTRIG_pt", &OthrMcElecTRIG_pt_,"OthrMcElecTRIG_pt/D");
91  t_->Branch("OthrMcElecTRIG_eta", &OthrMcElecTRIG_eta_,"OthrMcElecTRIG_eta/D");
92  t_->Branch("OthrMcElecRECO_pt", &OthrMcElecRECO_pt_,"OthrMcElecRECO_pt/D");
93  t_->Branch("OthrMcElecRECO_eta", &OthrMcElecRECO_eta_,"OthrMcElecRECO_eta/D");
94  t_->Branch("McNu_pt", &McNu_pt_,"McNu_pt/D");
95  t_->Branch("McNu_eta", &McNu_eta_,"McNu_eta/D");
96  t_->Branch("McNu_phi", &McNu_phi_,"McNu_phi/D");
97  t_->Branch("McNu_vx", &McNu_vx_,"McNu_vx/D");
98  t_->Branch("McNu_vy", &McNu_vy_,"McNu_vy/D");
99  t_->Branch("McNu_vz", &McNu_vz_,"McNu_vz/D");
100  t_->Branch("McLeptons_dPhi", &McLeptons_dPhi_,"McLeptons_dPhi/D");
101  t_->Branch("McLeptons_dEta", &McLeptons_dEta_,"McLeptons_dEta/D");
102  t_->Branch("McLeptons_dR", &McLeptons_dR_,"McLeptons_dR/D");
103  t_->Branch("nSelElecs", &nSelElecs_,"nSelElecs/I");
104  t_->Branch("elec_q", elec_q_,"elec_q[4]/D");
105  t_->Branch("elec_pt", elec_pt_,"elec_pt[4]/D");
106  t_->Branch("elec_eta", elec_eta_,"elec_eta[4]/D");
107  t_->Branch("elec_phi", elec_phi_,"elec_phi[4]/D");
108  t_->Branch("elec_pt25", &elec_pt25_,"elec_pt25/D");
109  t_->Branch("elec_eta25", &elec_eta25_,"elec_eta25/D");
110  t_->Branch("elec_phi25", &elec_phi25_,"elec_phi25/D");
111  t_->Branch("elec_sIhIh", elec_sIhIh_, "elec_sIhIh[4]/D");
112  t_->Branch("elec_dPhiIn", elec_dPhiIn_, "elec_dPhiIn[4]/D");
113  t_->Branch("elec_dEtaIn", elec_dEtaIn_, "elec_dEtaIn[4]/D");
114  t_->Branch("elec_trckIso", elec_trckIso_,"elec_trckIso[4]/D");
115  t_->Branch("elec_ecalIso", elec_ecalIso_,"elec_ecalIso[4]/D");
116  t_->Branch("elec_hcalIso", elec_hcalIso_,"elec_hcalIso[4]/D");
117  t_->Branch("elec_e2x5Max", elec_e2x5Max_,"elec_e2x5Max[4]/D");
118  t_->Branch("elec_e1x5Max", elec_e1x5Max_,"elec_e1x5Max[4]/D");
119  t_->Branch("elec_e5x5", elec_e5x5_,"elec_e5x5[4]/D");
120  t_->Branch("elec_hoe", elec_hoe_,"elec_hoe[4]/D");
121  t_->Branch("elec_eop", elec_eop_,"elec_eop[4]/D");
122  t_->Branch("elec_pin", elec_pin_,"elec_pin[4]/D");
123  t_->Branch("elec_pout", elec_pout_,"elec_pout[4]/D");
124 
125  t_->Branch("Selected_nuPt", Selected_nuPt_, "Selected_nuPt[4]/D");
126  t_->Branch("Selected_nuEta", Selected_nuEta_, "Selected_nuEta[4]/D");
127  t_->Branch("Selected_nuPhi", Selected_nuPhi_, "Selected_nuPhi[4]/D");
128  t_->Branch("caloMEt", &caloMEt_,"caloMEt/D");
129  t_->Branch("t1MEt", &t1MEt_,"t1MEt/D");
130  t_->Branch("t1MEtPhi", &t1MEtPhi_,"t1MEtPhi/D");
131  t_->Branch("t1SumEt", &t1SumEt_,"t1SumEt/D");
132  t_->Branch("pfMEt", &pfMEt_,"pfMEt/D");
133  t_->Branch("pfMEtPhi", &pfMEtPhi_,"pfMEtPhi/D");
134  t_->Branch("pfSumEt", &pfSumEt_,"pfSumEt/D");
135  t_->Branch("tcMEt", &tcMEt_,"tcMEt/D");
136  t_->Branch("tcMEtPhi", &tcMEtPhi_,"tcMEtPhi/D");
137  t_->Branch("tcSumEt", &tcSumEt_,"tcSumEt/D");
138  t_->Branch("caloSumEt", &caloSumEt_,"caloSumEt/D");
139  t_->Branch("caloMEt25", &caloMEt25_,"caloMEt25/D");
140  t_->Branch("caloMEt30", &caloMEt30_,"caloMEt30/D");
141  t_->Branch("caloMEtECAL25", &caloMEtECAL25_,"caloMEtECAL25/D");
142  t_->Branch("caloMEtECAL30", &caloMEtECAL30_,"caloMEtECAL30/D");
143  t_->Branch("caloMEtPhi", &caloMEtPhi_,"caloMEtPhi/D");
144  t_->Branch("caloMEtPhi25", &caloMEtPhi25_,"caloMEtPhi25/D");
145  t_->Branch("caloMEtPhi30", &caloMEtPhi30_,"caloMEtPhi30/D");
146  t_->Branch("caloMEtPhiECAL25", &caloMEtPhiECAL25_,"caloMEtPhiECAL25/D");
147  t_->Branch("caloMEtPhiECAL30", &caloMEtPhiECAL30_,"caloMEtPhiECAL30/D");
148  t_->Branch("caloMt", caloMt_,"caloMt[4]/D");
149  t_->Branch("caloMt25", caloMt25_,"caloMt25[4]/D");
150  t_->Branch("caloMt30", caloMt30_,"caloMt30[4]/D");
151  t_->Branch("caloUESumEt", &caloUESumEt_, "caloUESumEt_/D");
152  t_->Branch("nHltObj", &nHltObj_, "nHltObj/I");
153  t_->Branch("HltObj_pt_", HltObj_pt_,"HltObj_pt_[4]/D");
154  t_->Branch("HltObj_eta_", HltObj_eta_,"HltObj_eta_[4]/D");
155  t_->Branch("genMEt", &genMEt_,"genMEt/D");
156  t_->Branch("genUESumEt", &genUESumEt_, "genUESumEt_/D");
157 
158 }
159 // ------------ method called to for each event ------------
160 void AnalysisErsatz::analyze(const edm::Event& evt, const edm::EventSetup& es)
161 {
162  caloMEt_ = -99.; caloSumEt_ = -99.; caloUESumEt_ = -99.;
163  caloMEt25_ = -99.; caloMEt30_ = -99.;
164  caloMEtECAL25_ = -99.; caloMEtECAL30_ =- 99.;
165  caloMEtPhi_ = -99.; caloMEtPhi25_ = -99.; caloMEtPhi30_ = -99.;
166  caloMEtPhiECAL25_ = -99.; caloMEtPhiECAL30_ =- 99.;
167  genMEt_ = -99.; genUESumEt_ = -99.; genMEt25_ = -99.;
168  t1MEt_ = -99.; t1MEtPhi_ = -99.; t1SumEt_ = -99.;
169  pfMEt_ = -99.; pfMEtPhi_ = -99.; pfSumEt_ = -99.;
170  tcMEt_ = -99.; tcMEtPhi_ = -99.; tcSumEt_ = -99.;
171  nHltObj_ = -99; nSelElecs_ = -99;
172  Boson_pt_ = -99.; Boson_y_ = -99.; Boson_m_ = -99.; Boson_mt_ = -99.; Boson_phi_ = -99.;
173  McNu_pt_ = -99.; McNu_eta_ = -99.; McNu_phi_ = -99.;
174  McNu_vx_= -99.; McNu_vy_= -99.; McNu_vz_ = -99.;
175  McLeptons_dPhi_ = -99.; McLeptons_dEta_ = -99.; McLeptons_dR_ = -99.;
176  RndmMcElec_pt_ = -99.; RndmMcElec_eta_ = -99.; RndmMcElec_phi_ = -99.;
177  RndmMcElec_Rescaled_pt_ = -99.; RndmMcElec_Rescaled_eta_ = -99.; RndmMcElec_Rescaled_phi_ = -99.;
178  RndmMcElecTRIG_pt_ = -99.; RndmMcElecTRIG_eta_ = -99.; RndmMcElecRECO_pt_ = -99.; RndmMcElecRECO_eta_ = -99.;
179  OthrMcElec_pt_ = -99.; OthrMcElec_eta_ = -99.; OthrMcElec_phi_ = -99.;
180  OthrMcElec_Rescaled_pt_ = -99.; OthrMcElec_Rescaled_eta_ = -99.; OthrMcElec_Rescaled_phi_ = -99.;
181  OthrMcElecTRIG_pt_ = -99.; OthrMcElecTRIG_eta_ = -99.; OthrMcElecRECO_pt_ = -99.; OthrMcElecRECO_eta_ = -99.;
182  RndmTrig_ = 0; RndmReco_ = 0; OthrTrig_ = 0; OthrReco_ = 0;
183  elec_pt25_ = -99.; elec_eta25_= -99.; elec_eta25_= -99.;
184  for(int i = 0; i < nEntries_arr_; ++i)
185  {
186  McElec1_pt_[i] = -99.; McElec1_eta_[i] = -99.;
187  McElec3_pt_[i] = -99.; McElec3_eta_[i] = -99.;
188  elec_q_[i] = -99.; elec_pt_[i] = -99.; elec_eta_[i]= -99.; elec_phi_[i]= -99.;
189  elec_trckIso_[i] = -99.; elec_ecalIso_[i] = -99.; elec_hcalIso_[i] = -99.;
190  elec_sIhIh_[i] = -99.; elec_dPhiIn_[i] = -99.; elec_dEtaIn_[i] = -99.;
191  elec_e5x5_[i] = -99.; elec_e2x5Max_[i] = -99.; elec_e1x5Max_[i] = -99.;
192  elec_hoe_[i] = -99.; elec_eop_[i] = -99.; elec_pin_[i] = -99.; elec_pout_[i] = -99.;
193  Selected_nuPt_[i] = -99.; Selected_nuEta_[i] = -99.; Selected_nuPhi_[i] = -99.;
194  caloMt_[i] = -99.; caloMt25_[i] = -99.; caloMt30_[i] = -99.;
195  genMt_[i] = -99.;
196  }
197  edm::Handle<reco::GenParticleCollection> pGenPart;
198  try
199  {
200  evt.getByLabel(MCTruthCollection_, pGenPart);
201  }catch(cms::Exception& ex)
202  {
203  edm::LogError("")<<"Error! Can't get collection with label "<< MCTruthCollection_.label();
204  }
205  edm::Handle<reco::GsfElectronCollection> pElectrons;
206  try
207  {
208  evt.getByLabel(ElectronCollection_, pElectrons);
209  }catch(cms::Exception &ex)
210  {
211  edm::LogError("analyze") <<"Can't get collection with label "<< ElectronCollection_.label();
212  }
213  edm::Handle<reco::CaloMETCollection> pCaloMEt;
214  try
215  {
216  evt.getByLabel(CaloMEtCollection_, pCaloMEt);
217  }catch(cms::Exception &ex)
218  {
219  edm::LogError("analyze")<<"Can't get collection with label "<< CaloMEtCollection_.label();
220  }
221  //edm::Handle<reco::METCollection> pT1MEt;
222  //try
223  //{
224  // evt.getByLabel(T1MEtCollection_, pT1MEt);
225  //}catch(cms::Exception &ex)
226  //{
227  // edm::LogError("analyze")<<"Can't get collection with label "<< T1MEtCollection_.label();
228  //}
229  edm::Handle<reco::PFMETCollection> pPfMEt;
230  try
231  {
232  evt.getByLabel(PfMEtCollection_, pPfMEt);
233  }catch(cms::Exception &ex)
234  {
235  edm::LogError("analyze")<<"Can't get collection with label "<< PfMEtCollection_.label();
236  }
237  edm::Handle<reco::METCollection> pTcMEt;
238  try
239  {
240  evt.getByLabel(TcMEtCollection_, pTcMEt);
241  }catch(cms::Exception &ex)
242  {
243  edm::LogError("analyze")<<"Can't get collection with label "<< TcMEtCollection_.label();
244  }
245  edm::Handle<reco::GenMETCollection> pGenMEt;
246  try
247  {
248  evt.getByLabel(GenMEtCollection_, pGenMEt);
249  }catch(cms::Exception &ex)
250  {
251  edm::LogError("analyze")<<"Can't get collection with label "<< GenMEtCollection_.label();
252  }
253  edm::Handle<edm::TriggerResults> pTriggerResults;
254  try
255  {
256  evt.getByLabel(TriggerResults_, pTriggerResults);
257  }catch(cms::Exception &ex)
258  {
259  edm::LogError("analyze")<<"Cant get collection with label "<< TriggerResults_.label();
260  }
261  const edm::TriggerResults* HltRes = pTriggerResults.product();
262  edm::Handle<trigger::TriggerEvent> pHLT;
263  if(HltRes->accept(34) !=0)
264  {
265  try
266  {
267  evt.getByLabel(TriggerEvent_, pHLT);
268  }catch(cms::Exception &ex)
269  {
270  edm::LogError("analyze")<<"Can't get collection with label "<< TriggerEvent_.label();
271  }
272  }
273  edm::LogDebug_("analyse","", 143)<<"Have obtained collections."<<std::endl;
274  int m = 0; int n = 0; int k = 0;
275  bool BosonAnalysed = false; bool nuInEta25 = false; bool nuInEta30 = false;
276  double elecEt = 0;
277  std::vector<reco::GenParticleCollection::const_iterator> Leptons;
278  const reco::GenParticleCollection *McCand = pGenPart.product();
279  math::XYZTLorentzVector RndmMcElec, OthrMcElec;
280  edm::Service<edm::RandomNumberGenerator> rng;
281  math::XYZTLorentzVector RndmMcElec_alt, OthrMcElec_alt, Zboson;
282  CLHEP::RandFlat flatDistribution(rng->getEngine(), 0, 2);
283  double RandomNum = flatDistribution.fire();
284  RndmInt_ = int(RandomNum);
285 // std::cout<<"Random Number = "<< RandomNum <<"\t int = "<< int(RandomNum)<<std::endl;
286  for(reco::GenParticleCollection::const_iterator McP = McCand->begin(); McP != McCand->end(); ++McP)
287  {
288  const reco::Candidate* mum = McP->mother();
289  if(McP->pdgId()==11)
290  {
291  edm::LogDebug_("", "", 216)<<"Found electron, mother = "<< mum->pdgId() <<"\t status = "<< McP->status()
292  <<"\tpt = "<< McP->pt() <<"\teta = "<< McP->eta();
293  }
294  if(McP->pdgId()==-11)
295  {
296  edm::LogDebug_("", "", 221)<<"Found positron, mother = "<< mum->pdgId() <<"\t status = "<< McP->status()
297  <<"\tpt = "<< McP->pt() <<"\teta = "<< McP->eta();
298  }
299  if(abs(McP->pdgId())==12)
300  {
301  edm::LogDebug_("", "", 216)<<"Found neutrino, mother = "<< mum->pdgId() <<"\t status = "<< McP->status()
302  <<"\tpt = "<< McP->pt() <<"\teta = "<< McP->eta();
303  }
304  if(abs(McP->pdgId())==11 && (abs(mum->pdgId()) == 24|| abs(mum->pdgId()) == 23))
305  {
306 // Leptons.push_back(McP);
307  if(!BosonAnalysed)
308  {
309  Boson_pt_ = mum->pt();
310  Boson_y_ = mum->y();
311  Boson_phi_ = mum->phi();
312  Boson_m_ = mum->mass();
313  Boson_mt_ = mum->mt();
314  if(abs(mum->pdgId() == 23)) Zboson = mum->p4();
315  BosonAnalysed = true;
316  }
317  McElec3_pt_[k] = McP->pt();
318  McElec3_eta_[k] = McP->eta();
319  if(k == int(RandomNum))
320  {
321  RndmMcElec_alt = McP->p4();
322  }else OthrMcElec_alt = McP->p4();
323  ++k;
324  }
325  if(abs(McP->pdgId())==12 && abs(mum->pdgId()) == 24)
326  {
327  Leptons.push_back(McP);
328  edm::LogDebug_("","",328)<<"Pushed neutrino back into Leptons. Leptons.size() = "<< Leptons.size();
329  McNu_pt_ = McP->pt();
330  McNu_eta_ = McP->eta();
331  edm::LogDebug_("","",332)<<"ECAL eta = "<< McNu_ECALeta_;
332  McNu_phi_ = McP->phi();
333  McNu_vx_ = McP->vx();
334  McNu_vy_ = McP->vy();
335  McNu_vz_ = McP->vz();
336  if(fabs(McNu_eta_) < 2.5) nuInEta25 = true;
337  if(fabs(McNu_eta_) < 3.0) nuInEta30 = true;
338  ++n;
339  }
340  if(abs(McP->pdgId())==11 && McP->status() == 1 && (abs(mum->pdgId()) == 11))
341  {
342  Leptons.push_back(McP);
343  edm::LogDebug_("","",344)<<"Pushed electron back into Leptons. Leptons.size() = "<< Leptons.size();
344  McElec1_pt_[m] = McP->pt();
345  McElec1_eta_[m] = McP->eta();
346  if(m == int(RandomNum))
347  {
348  RndmMcElec = McP->p4();
349  RndmMcElec_pt_ = McElec1_pt_[m];
350  RndmMcElec_eta_ = McElec1_eta_[m];
351  RndmMcElec_phi_ = McP->phi();
352  }
353  else{
354  OthrMcElec = McP->p4();
355  OthrMcElec_pt_ = McElec1_pt_[m];
356  OthrMcElec_eta_ = McElec1_eta_[m];
357  OthrMcElec_phi_ = McP->phi();
358  }
359 
360  elecEt += McP->pt();
361  ++m;
362  }
363  }
364  edm::LogDebug_("", "", 362)<<"Size of Leptons = "<< Leptons.size();
365  McLeptons_dPhi_ = reco::deltaPhi(Leptons[0]->phi(), Leptons[1]->phi());
366  McLeptons_dEta_ = Leptons[0]->eta() - Leptons[1]->eta();
367  McLeptons_dR_ = reco::deltaR(*Leptons[0], *Leptons[1]);
368  edm::LogDebug_("","",369)<<"McLeptons_dR_ = "<< McLeptons_dR_;
369 
370  math::XYZTLorentzVector Wboson;
371  if(McNu_pt_ < -98.)
372  {
373  edm::LogDebug_("", "", 303)<<"Z pt = "<< Zboson.Pt() << "Z boson mass = " << Zboson.M();
374  edm::LogDebug_("","", 307)<<"Z boson in lab frame = ("<<Zboson.Px()<<", "<<Zboson.Py()<<", "
375  <<Zboson.Pz()<<", "<<Zboson.E()<<")";
376 
377  Wboson.SetCoordinates(Zboson.Px(), Zboson.Py(), Zboson.Pz(), sqrt(Zboson.P2()+(mW_*mW_*Zboson.M2())/(mZ_*mZ_)));
378  edm::LogDebug_("","", 307)<<"W boson in lab frame = ("<<Wboson.Px()<<", "<<Wboson.Py()<<", "
379  <<Wboson.Pz()<<", "<<Wboson.E()<<")";
380  ROOT::Math::Boost CoMBoost(Zboson.BoostToCM());
381  edm::LogDebug_("","", 307)<<"RndmElec in lab frame = ("<<RndmMcElec_alt.Px()<<", "<<RndmMcElec_alt.Py()<<", "
382  <<RndmMcElec_alt.Pz()<<", "<<RndmMcElec_alt.E()<<")";
383  edm::LogDebug_("","", 307)<<"OthrElec in lab frame = ("<<OthrMcElec_alt.Px()<<", "<<OthrMcElec_alt.Py()<<", "
384  <<OthrMcElec_alt.Pz()<<", "<<OthrMcElec_alt.E()<<")";
385  RndmMcElec_alt = CoMBoost(RndmMcElec_alt);
386  OthrMcElec_alt = CoMBoost(OthrMcElec_alt);
387  edm::LogDebug_("","", 307)<<"RndmElec in Z rest frame = ("<<RndmMcElec_alt.Px()<<", "<<RndmMcElec_alt.Py()<<", "
388  <<RndmMcElec_alt.Pz()<<", "<<RndmMcElec_alt.E()<<")";
389  edm::LogDebug_("","", 307)<<"OthrElec in Z rest frame = ("<<OthrMcElec_alt.Px()<<", "<<OthrMcElec_alt.Py()<<", "
390  <<OthrMcElec_alt.Pz()<<", "<<OthrMcElec_alt.E()<<")";
391  RndmMcElec_alt *= mW_/mZ_;
392  OthrMcElec_alt *= mW_/mZ_;
393 
394  double E_W = Wboson.E();
395  ROOT::Math::Boost BackToLab(Wboson.Px()/E_W, Wboson.Py()/E_W, Wboson.Pz()/E_W);
396 
397  RndmMcElec_alt = BackToLab(RndmMcElec_alt);
398  RndmMcElec_Rescaled_pt_ = RndmMcElec_alt.Pt();
399  RndmMcElec_Rescaled_eta_ = RndmMcElec_alt.Eta();
400  RndmMcElec_Rescaled_phi_ = RndmMcElec_alt.Phi();
401 
402  OthrMcElec_alt = BackToLab(OthrMcElec_alt);
403  OthrMcElec_Rescaled_pt_ = OthrMcElec_alt.Pt();
404  OthrMcElec_Rescaled_eta_ = OthrMcElec_alt.Eta();
405  OthrMcElec_Rescaled_phi_ = OthrMcElec_alt.Phi();
406 
407  math::XYZTLorentzVector sum = OthrMcElec_alt+RndmMcElec_alt;
408  edm::LogDebug_("","", 307)<<"RndmElec back in lab frame = ("<<RndmMcElec_alt.Px()<<", "<<RndmMcElec_alt.Py()<<", "
409  <<RndmMcElec_alt.Pz()<<", "<<RndmMcElec_alt.E()<<")";
410  edm::LogDebug_("","", 307)<<"OthrElec back in lab frame = ("<<OthrMcElec_alt.Px()<<", "<<OthrMcElec_alt.Py()<<", "
411  <<OthrMcElec_alt.Pz()<<", "<<OthrMcElec_alt.E()<<")";
412  edm::LogDebug_("","", 307)<<"OthrElec +RndmElec = ("<<sum.Px()<<", "<<sum.Py()<<", "
413  <<sum.Pz()<<", "<<sum.E()<<")";
414  }else{
415  edm::LogDebug_("","", 416)<<"McNu_pt_ = "<<McNu_pt_;
416  RndmMcElec_Rescaled_pt_ = RndmMcElec_pt_;
417  edm::LogDebug_("","",416)<<" RndmMcElec_Rescaled_pt_ = "<< RndmMcElec_Rescaled_pt_;
418  OthrMcElec_Rescaled_pt_ = OthrMcElec_pt_;
419  }
420  //TriggerNames_.init(*HltRes);
421 
422  edm::LogDebug_("","", 420)<<"HltRes->accept() = "<< HltRes->accept(34);
423  if(HltRes->accept(34) ==0) edm::LogError("")<<"Event did not pass HLT path 34, assumed to be "<< TriggerName_ <<"!";
424  unsigned int fId = 999;
425  if(HltRes->accept(34) !=0)
426  {
427  fId = pHLT->filterIndex(TriggerPath_); // something wrong with this step
428  edm::LogDebug_("","",426)<<"fId = pHLT->filterIndex("<< TriggerPath_<<") = "<< fId;
429  const trigger::Keys& ring = pHLT->filterKeys(fId);
430  const trigger::TriggerObjectCollection& HltObjColl = pHLT->getObjects();
431  nHltObj_ = ring.size();
432  for(int k = 0; k < nHltObj_; ++k)
433  {
434  const trigger::TriggerObject& HltObj = HltObjColl[ring[k]];
435  if(reco::deltaR(RndmMcElec, HltObj) < 0.1) RndmTrig_ = 1;
436  if(reco::deltaR(OthrMcElec, HltObj) < 0.1) OthrTrig_ = 1;
437  if(k < 4)
438  {
439  HltObj_pt_[k] = HltObj.pt();
440  HltObj_eta_[k] = HltObj.eta();
441  }
442  }
443  }
444  if(RndmTrig_ != 0)
445  {
446  RndmMcElecTRIG_pt_ = RndmMcElec_pt_;
447  RndmMcElecTRIG_eta_ = RndmMcElec_eta_;
448  }
449  if(OthrTrig_ != 0)
450  {
451  OthrMcElecTRIG_pt_ = OthrMcElec_pt_;
452  OthrMcElecTRIG_eta_ = OthrMcElec_eta_;
453  }
454  const reco::GenMETCollection* genMEtCollection = pGenMEt.product();
455  reco::GenMETCollection::const_iterator genMEt = genMEtCollection->begin();
456  genMEt_ = genMEt->pt();
457  genUESumEt_ = genMEt->sumEt() - elecEt;
458 // std::cout<<"genMEt->sumEt() - elecEt = "<< genMEt->sumEt()<<" - "<< elecEt <<" = "<< genUESumEt_;
459 
460  const reco::CaloMETCollection* caloMEtCollection = pCaloMEt.product();
461  reco::CaloMETCollection::const_iterator met = caloMEtCollection->begin();
462  caloMEt_ = met->pt();
463  edm::LogDebug_("","",462)<<"caloMEt_ = "<< caloMEt_;
464  caloMEtPhi_ = met->phi();
465  caloSumEt_ = met->sumEt();
466 
467  //const reco::METCollection* t1MEtCollection = pT1MEt.product();
468  //reco::METCollection::const_iterator t1met = t1MEtCollection->begin();
469  //t1MEt_ = t1met->pt();
470  //edm::LogDebug_("","",462)<<"t1MEt_ = "<< t1MEt_;
471  //t1MEtPhi_ = t1met->phi();
472  //t1SumEt_ = t1met->sumEt();
473 
474  const reco::PFMETCollection* pfMEtCollection = pPfMEt.product();
475  reco::PFMETCollection::const_iterator pfmet = pfMEtCollection->begin();
476  pfMEt_ = pfmet->pt();
477  edm::LogDebug_("","",462)<<"pfMEt_ = "<< pfMEt_;
478  pfMEtPhi_ = pfmet->phi();
479  pfSumEt_ = pfmet->sumEt();
480 
481  const reco::METCollection* tcMEtCollection = pTcMEt.product();
482  reco::METCollection::const_iterator tcmet = tcMEtCollection->begin();
483  tcMEt_ = tcmet->pt();
484  edm::LogDebug_("","",462)<<"tcMEt_ = "<< tcMEt_;
485  tcMEtPhi_ = tcmet->phi();
486  tcSumEt_ = tcmet->sumEt();
487 
488  if(fabs(McNu_ECALeta_) < 2.5)
489  {
490  caloMEtECAL25_ = met->pt();
491  caloMEtPhiECAL25_ = met->phi();
492  }
493  if(fabs(McNu_ECALeta_) < 3.0)
494  {
495  caloMEtECAL30_ = met->pt();
496  caloMEtPhiECAL30_ = met->phi();
497  }
498  if(nuInEta25)
499  {
500  genMEt25_ = genMEt->pt();
501  caloMEt25_ = met->pt();
502  caloMEtPhi25_ = met->phi();
503  }
504  if(nuInEta30){
505  caloMEt30_ = met->pt();
506  caloMEtPhi30_ = met->phi();
507  }
508 
509  std::vector<reco::GsfElectronRef> UniqueElectrons = uniqueElectronFinder(pElectrons);
510  caloUESumEt_ = met->sumEt();
511  for(std::vector<reco::GsfElectronRef>::const_iterator Relec = UniqueElectrons.begin(); Relec != UniqueElectrons.end(); ++Relec)
512  {
513  reco::GsfElectronRef elec = *Relec;
514  math::XYZTLorentzVector sc = PhysicsVector(met->vertex(), *(elec->superCluster()));
515  if(reco::deltaR(RndmMcElec, *elec) < 0.1) RndmReco_ = 1;
516  if(reco::deltaR(OthrMcElec, *elec) < 0.1) OthrReco_ = 1;
517  caloUESumEt_ -= sc.Pt();
518  }
519  if(RndmReco_ != 0)
520  {
521  RndmMcElecRECO_pt_ = RndmMcElec_pt_;
522  RndmMcElecRECO_eta_ = RndmMcElec_eta_;
523  }
524  if(OthrReco_ != 0)
525  {
526  OthrMcElecRECO_pt_ = OthrMcElec_pt_;
527  OthrMcElecRECO_eta_ = OthrMcElec_eta_;
528  }
529  edm::LogDebug_("analyse","", 230)<<"Analysed UE information"<<std::endl;
530  if(HltRes->accept(34) != 0)
531  {
532  std::vector<reco::GsfElectronRef> SelectedElectrons = electronSelector(UniqueElectrons, pHLT, fId, CutVector_);
533  nSelElecs_ = SelectedElectrons.size();
534  m = 0;
535  for(std::vector<reco::GsfElectronRef>::const_iterator Relec = SelectedElectrons.begin(); Relec != SelectedElectrons.end(); ++Relec)
536  {
537  reco::GsfElectronRef elec = *Relec;
538  if(elec->pt() > CutVector_[0])
539  {
540  elec_q_[m] = elec->charge();
541  elec_pt_[m] = elec->pt();
542  elec_eta_[m] = elec->eta();
543  elec_phi_[m] = elec->phi();
544  elec_trckIso_[m] = elec->isolationVariables03().tkSumPt;
545  elec_ecalIso_[m] = elec->isolationVariables04().ecalRecHitSumEt;
546  elec_hcalIso_[m] = elec->isolationVariables04().hcalDepth1TowerSumEt
547  + elec->isolationVariables04().hcalDepth2TowerSumEt;
548  elec_sIhIh_[m] = elec->scSigmaIEtaIEta();
549  elec_dPhiIn_[m] = elec->deltaPhiSuperClusterTrackAtVtx();
550  elec_dEtaIn_[m] = elec->deltaEtaSuperClusterTrackAtVtx();
551  elec_e5x5_[m] = elec->scE5x5();
552  elec_e2x5Max_[m] = elec->scE2x5Max();
553  elec_e2x5Max_[m] = elec->scE1x5();
554  elec_hoe_[m] = elec->hadronicOverEm();
555  elec_eop_[m] = elec->eSuperClusterOverP();
556  elec_pin_[m] = elec->trackMomentumAtVtx().R();
557  elec_pout_[m] = elec->trackMomentumOut().R();
558 
559  Selected_nuPt_[m] = McNu_pt_;
560  Selected_nuEta_[m] = McNu_eta_;
561  Selected_nuPhi_[m] = McNu_phi_;
562  caloMt_[m] = sqrt(2.*elec->pt()*met->pt()*(1-cos(reco::deltaPhi(elec->phi(), met->phi()))));
563  if(nuInEta25)
564  {
565  caloMt25_[m] = caloMt_[m];
566  elec_pt25_ = elec->pt();
567  elec_eta25_ = elec->eta();
568  elec_phi25_ = elec->phi();
569  }
570  if(nuInEta30) caloMt30_[m] = caloMt_[m];
571  ++m;
572  }
573  }
574  }
575 
576  edm::LogDebug_("analyse","", 248)<<"Analysed final selection information"<<std::endl;
577  t_->Fill();
578 }
579 
580 
581 
582 // ------------ method called once each job just after ending the event loop ------------
583 void AnalysisErsatz::endJob() {
584 }
585 
586 //define this as a plug-in
587 DEFINE_FWK_MODULE(AnalysisErsatz);
588 */