/ElectroWeakAnalysis/ZEE/src/ErsatzMEt.cc
C++ | 948 lines | 777 code | 68 blank | 103 comment | 60 complexity | 2313a34fb4c823d43c26c4f86f2bc9f0 MD5 | raw file
1#include "ElectroWeakAnalysis/ZEE/interface/ErsatzMEt.h" 2#include "FWCore/Common/interface/TriggerNames.h" 3 4ErsatzMEt::ErsatzMEt(const edm::ParameterSet& ps) 5{ 6 MCTruthCollection_ = ps.getParameter<edm::InputTag>("MCTruthCollection"); 7 ElectronCollection_ = ps.getParameter<edm::InputTag>("ElectronCollection"); 8 HybridScCollection_ = ps.getParameter<edm::InputTag>("HybridScCollection"); 9 M5x5ScCollection_ = ps.getParameter<edm::InputTag>("M5x5ScCollection"); 10 GenMEtCollection_ = ps.getParameter<edm::InputTag>("GenMEtCollection"); 11 CaloMEtCollection_ = ps.getParameter<edm::InputTag>("CaloMEtCollection"); 12 //T1MEtCollection_ = ps.getParameter<edm::InputTag>("T1MEtCollection"); 13 PfMEtCollection_ = ps.getParameter<edm::InputTag>("PfMEtCollection"); 14 TcMEtCollection_ = ps.getParameter<edm::InputTag>("TcMEtCollection"); 15 TriggerEvent_ = ps.getParameter<edm::InputTag>("TriggerEvent"); 16 TriggerPath_ = ps.getParameter<edm::InputTag>("TriggerPath"); 17 TriggerResults_ = ps.getParameter<edm::InputTag>("TriggerResults"); 18 TriggerName_ = ps.getParameter<std::string>("TriggerName"); 19 HLTPathCheck_ = ps.getParameter<bool>("HLTPathCheck"); 20 21 Zevent_ = ps.getParameter<bool>("Zevent"); 22 mW_ = ps.getParameter<double>("mW"); 23 mZ_ = ps.getParameter<double>("mZ"); 24 mTPmin_ = ps.getParameter<double>("mTPmin"); 25 mTPmax_ = ps.getParameter<double>("mTPmax"); 26 BarrelEtaMax_ = ps.getParameter<double>("BarrelEtaMax"); 27 EndCapEtaMin_ = ps.getParameter<double>("EndCapEtaMin"); 28 EndCapEtaMax_ = ps.getParameter<double>("EndCapEtaMax"); 29 30 hyb_fCorrPSet_ = ps.getParameter<edm::ParameterSet>("hyb_fCorrPSet"); 31 m5x5_fCorrPSet_ = ps.getParameter<edm::ParameterSet>("m5x5_fCorrPSet"); 32 33 double CElecPtMin = ps.getParameter<double>("CElecPtMin"); 34 double CEB_siEiE = ps.getParameter<double>("CEB_sigmaIEtaIEta"); 35 double CEB_dPhiIn = ps.getParameter<double>("CEB_deltaPhiIn"); 36 double CEB_dEtaIn = ps.getParameter<double>("CEB_deltaEtaIn"); 37 double CEB_EcalIso = ps.getParameter<double>("CEB_EcalIso"); 38 double CEB_HcalIso = ps.getParameter<double>("CEB_HcalIso"); 39 double CEB_TrckIso = ps.getParameter<double>("CEB_TrckIso"); 40 double CEE_siEiE = ps.getParameter<double>("CEE_sigmaIEtaIEta"); 41 double CEE_dPhiIn = ps.getParameter<double>("CEE_deltaPhiIn"); 42 double CEE_dEtaIn = ps.getParameter<double>("CEE_deltaEtaIn"); 43 double CEE_EcalIso = ps.getParameter<double>("CEE_EcalIso"); 44 double CEE_HcalIso = ps.getParameter<double>("CEE_HcalIso"); 45 double CEE_TrckIso = ps.getParameter<double>("CEE_TrckIso"); 46 47 CutVector_.resize(13); 48 CutVector_[EtCut_] = CElecPtMin; 49 CutVector_[EB_sIhIh_] = CEB_siEiE; 50 CutVector_[EB_dPhiIn_] = CEB_dPhiIn; 51 CutVector_[EB_dEtaIn_] = CEB_dEtaIn; 52 CutVector_[EB_TrckIso_] = CEB_TrckIso; 53 CutVector_[EB_EcalIso_] = CEB_EcalIso; 54 CutVector_[EB_HcalIso_] = CEB_HcalIso; 55 CutVector_[EE_sIhIh_] = CEE_siEiE; 56 CutVector_[EE_dPhiIn_] = CEE_dPhiIn; 57 CutVector_[EE_dEtaIn_] = CEE_dEtaIn; 58 CutVector_[EE_TrckIso_] = CEE_TrckIso; 59 CutVector_[EE_EcalIso_] = CEE_EcalIso; 60 CutVector_[EE_HcalIso_] = CEE_HcalIso; 61 62 for(std::vector<double>::const_iterator it = CutVector_.begin(); it != CutVector_.end(); ++it) 63 { 64 edm::LogDebug_("","",101)<<"CutVector_ = "<< *it; 65 } 66} 67 68ErsatzMEt::~ErsatzMEt() 69{ 70} 71 72// ------------ method called once each job just before starting event loop ------------ 73void ErsatzMEt::beginJob() 74{ 75 edm::Service<TFileService> fs; 76 77 t_ = fs->make<TTree>("ErsatzMEt", "Data on ErsatzMEt"); 78 79 edm::LogDebug_("","", 75)<<"Creating Ersatz MEt branches."; 80 81 t_->Branch("nTags", &nTags_, "nTags/I"); 82 t_->Branch("nProbes", &nProbes_, "nProbes/I"); 83 84 t_->Branch("ErsatzV1CaloMEt", ErsatzV1CaloMEt_, "ErsatzV1CaloMEt[4]/D"); 85 t_->Branch("ErsatzV1CaloMt", ErsatzV1CaloMt_, "ErsatzV1CaloMt[4]/D"); 86 t_->Branch("ErsatzV1CaloMEtPhi", ErsatzV1CaloMEtPhi_, "ErsatzV1CaloMEtPhi[4]/D"); 87 t_->Branch("ErsatzV2CaloMEt", ErsatzV2CaloMEt_, "ErsatzV2CaloMEt[4]/D"); 88 t_->Branch("ErsatzV2CaloMEtPhi", ErsatzV2CaloMEtPhi_, "ErsatzV2CaloMEtPhi[4]/D"); 89 t_->Branch("ErsatzV2CaloMt", ErsatzV2CaloMt_, "ErsatzV2CaloMt[4]/D"); 90 t_->Branch("ErsatzV3CaloMEt", ErsatzV3CaloMEt_, "ErsatzV3CaloMEt[4]/D"); 91 t_->Branch("ErsatzV3CaloMEtPhi", ErsatzV3CaloMEtPhi_, "ErsatzV3CaloMEtPhi[4]/D"); 92 t_->Branch("ErsatzV3CaloMt", ErsatzV3CaloMt_, "ErsatzV3CaloMt[4]/D"); 93 t_->Branch("ErsatzV4CaloMEt", ErsatzV4CaloMEt_, "ErsatzV4CaloMEt[4]/D"); 94 t_->Branch("ErsatzV4CaloMEtPhi", ErsatzV4CaloMEtPhi_, "ErsatzV4CaloMEtPhi[4]/D"); 95 t_->Branch("ErsatzV4CaloMt", ErsatzV4CaloMt_, "ErsatzV4CaloMt[4]/D"); 96 t_->Branch("ErsatzV1T1MEt", ErsatzV1T1MEt_, "ErsatzV1T1MEt[4]/D"); 97 t_->Branch("ErsatzV1T1Mt", ErsatzV1T1Mt_, "ErsatzV1T1Mt[4]/D"); 98 t_->Branch("ErsatzV1T1MEtPhi", ErsatzV1T1MEtPhi_, "ErsatzV1T1MEtPhi[4]/D"); 99 t_->Branch("ErsatzV1PfMEt", ErsatzV1PfMEt_, "ErsatzV1PfMEt[4]/D"); 100 t_->Branch("ErsatzV1PfMt", ErsatzV1PfMt_, "ErsatzV1PfMt[4]/D"); 101 t_->Branch("ErsatzV1PfMEtPhi", ErsatzV1TcMEtPhi_, "ErsatzV1PfMEtPhi[4]/D"); 102 t_->Branch("ErsatzV1TcMEt", ErsatzV1TcMEt_, "ErsatzV1TcMEt[4]/D"); 103 t_->Branch("ErsatzV1TcMt", ErsatzV1TcMt_, "ErsatzV1TcMt[4]/D"); 104 t_->Branch("ErsatzV1TcMEtPhi", ErsatzV1TcMEtPhi_, "ErsatzV1TcMEtPhi[4]/D"); 105 106 t_->Branch("CaloMEt", &CaloMEt_, "CaloMEt/D"); 107 t_->Branch("CaloMEtphi", &CaloMEtphi_, "CaloMEtphi/D"); 108 t_->Branch("T1MEt", &T1MEt_, "T1MEt/D"); 109 t_->Branch("T1MEtphi", &T1MEtphi_, "T1MEtphi/D"); 110 t_->Branch("PfMEt", &PfMEt_, "PfMEt/D"); 111 t_->Branch("PfMEtphi", &PfMEtphi_, "PfMEtphi/D"); 112 t_->Branch("TcMEt", &TcMEt_, "TcMEt/D"); 113 t_->Branch("TcMEtphi", &TcMEtphi_, "TcMEtphi/D"); 114 115 edm::LogDebug_("","", 91)<<"Creating electron branches."; 116 t_->Branch("tag_q", tag_q_,"tag_q[4]/I"); 117 t_->Branch("tag_pt", tag_pt_,"tag_pt[4]/D"); 118 t_->Branch("tag_eta", tag_eta_,"tag_eta[4]/D"); 119 t_->Branch("tag_phi", tag_phi_,"tag_phi[4]/D"); 120 t_->Branch("tag_sIhIh", tag_sIhIh_, "tag_sIhIh[4]/D"); 121 t_->Branch("tag_dPhiIn", tag_dPhiIn_, "tag_dPhiIn[4]/D"); 122 t_->Branch("tag_dEtaIn", tag_dEtaIn_, "tag_dEtaIn[4]/D"); 123 t_->Branch("tag_trckIso", tag_trckIso_,"tag_trckIso[4]/D"); 124 t_->Branch("tag_ecalIso", tag_ecalIso_,"tag_ecalIso[4]/D"); 125 t_->Branch("tag_hcalIso", tag_hcalIso_,"tag_hcalIso[4]/D"); 126 t_->Branch("tag_e2x5Max", tag_e2x5Max_,"tag_e2x5Max[4]/D"); 127 t_->Branch("tag_e1x5Max", tag_e1x5Max_,"tag_e1x5Max[4]/D"); 128 t_->Branch("tag_e5x5", tag_e5x5_,"tag_e5x5[4]/D"); 129 t_->Branch("tag_hoe", tag_hoe_,"tag_hoe[4]/D"); 130 t_->Branch("tag_eop", tag_eop_,"tag_eop[4]/D"); 131 t_->Branch("tag_pin", tag_pin_,"tag_pin[4]/D"); 132 t_->Branch("tag_pout", tag_pout_,"tag_pout[4]/D"); 133 t_->Branch("tag_rescPt", tag_rescPt_, "tag_rescPt[4]/D"); 134 t_->Branch("tag_rescEta", tag_rescEta_, "tag_rescEta[4]/D"); 135 t_->Branch("tag_rescPhi", tag_rescPhi_, "tag_rescPhi[4]/D"); 136 137 edm::LogDebug_("","", 103)<<"Creating ersatz neutrino branches."; 138 t_->Branch("probe_q", probe_q_,"probe_q[4]/I"); 139 t_->Branch("probe_pt", probe_pt_,"probe_pt[4]/D"); 140 t_->Branch("probe_eta", probe_eta_,"probe_eta[4]/D"); 141 t_->Branch("probe_phi", probe_phi_,"probe_phi[4]/D"); 142 t_->Branch("probe_sIhIh", probe_sIhIh_, "probe_sIhIh[4]/D"); 143 t_->Branch("probe_dPhiIn", probe_dPhiIn_, "probe_dPhiIn[4]/D"); 144 t_->Branch("probe_dEtaIn", probe_dEtaIn_, "probe_dEtaIn[4]/D"); 145 t_->Branch("probe_trckIso", probe_trckIso_,"probe_trckIso[4]/D"); 146 t_->Branch("probe_ecalIso", probe_ecalIso_,"probe_ecalIso[4]/D"); 147 t_->Branch("probe_hcalIso", probe_hcalIso_,"probe_hcalIso[4]/D"); 148 t_->Branch("probe_e2x5Max", probe_e2x5Max_,"probe_e2x5Max[4]/D"); 149 t_->Branch("probe_e1x5Max", probe_e1x5Max_,"probe_e1x5Max[4]/D"); 150 t_->Branch("probe_e5x5", probe_e5x5_,"probe_e5x5[4]/D"); 151 t_->Branch("probe_hoe", probe_hoe_,"probe_hoe[4]/D"); 152 t_->Branch("probe_eop", probe_eop_,"probe_eop[4]/D"); 153 t_->Branch("probe_pin", probe_pin_,"probe_pin[4]/D"); 154 t_->Branch("probe_pout", probe_pout_,"probe_pout[4]/D"); 155 t_->Branch("probe_rescPt", probe_rescPt_, "probe_rescPt[4]/D"); 156 t_->Branch("probe_rescEta", probe_rescEta_, "probe_rescEta[4]/D"); 157 t_->Branch("probe_rescPhi", probe_rescPhi_, "probe_rescPhi[4]/D"); 158 159 t_->Branch("Z_m", Z_m_, "Z_m[4]/D"); 160 t_->Branch("Z_pt", Z_pt_, "Z_pt[4]/D"); 161 t_->Branch("Z_eta", Z_eta_, "Z_eta[4]/D"); 162 t_->Branch("Z_y", Z_y_, "Z_y[4]/D"); 163 t_->Branch("Z_phi", Z_phi_, "Z_phi[4]/D"); 164 t_->Branch("Z_rescM", Z_rescM_, "Z_rescM[4]/D"); 165 t_->Branch("Z_rescPt", Z_rescPt_, "Z_rescPt[4]/D"); 166 t_->Branch("Z_rescEta", Z_rescEta_, "Z_rescEta[4]/D"); 167 t_->Branch("Z_rescY", Z_rescY_, "Z_rescY[4]/D"); 168 t_->Branch("Z_rescPhi", Z_rescPhi_, "Z_rescPhi[4]/D"); 169 t_->Branch("Z_probe_dPhi",Z_probe_dPhi_,"Z_probe_dPhi[4]/D"); 170 171 t_->Branch("probe_sc_pt",probe_sc_pt_, "probe_sc_pt[4]/D"); 172 t_->Branch("probe_sc_eta",probe_sc_eta_, "probe_sc_eta[4]/D"); 173 t_->Branch("probe_sc_phi", probe_sc_phi_, "probe_sc_phi[4]/D"); 174 t_->Branch("probe_sc_E",probe_sc_E_, "probe_sc_E[4]/D"); 175 t_->Branch("probe_sc_rawE",probe_sc_rawE_, "probe_sc_rawE[4]/D"); 176 t_->Branch("probe_sc_nClus", probe_sc_nClus_, "probe_sc_nClus[4]/D"); 177 t_->Branch("probe_scV2_E",probe_scV2_E_, "probe_scV2_E[4]/D"); 178 t_->Branch("probe_scV3_E",probe_scV3_E_, "probe_scV3_E[4]/D"); 179 t_->Branch("probe_scV4_E",probe_scV4_E_, "probe_scV4_E[4]/D"); 180 t_->Branch("probe_d_MCE_SCE", probe_d_MCE_SCE_, "probe_d_MCE_SCE[4]/D"); 181 182 t_->Branch("ErsatzV1_Mesc", ErsatzV1_Mesc_, "ErsatzV1_Mesc[4]/D"); 183 t_->Branch("ErsatzV1_rescMesc", ErsatzV1_rescMesc_, "ErsatzV1_rescMesc[4]/D"); 184 185 t_->Branch("McElec_nFinal", &McElec_nFinal_, "McElec_nFinal/I"); 186 187 if(Zevent_){ 188 t_->Branch("McZ_m", &McZ_m_, "McZ_m/D"); 189 t_->Branch("McZ_rescM", &McZ_rescM_, "McZ_rescM/D"); 190 t_->Branch("McZ_Pt", &McZ_pt_, "McZ_Pt/D"); 191 t_->Branch("McZ_y", &McZ_y_, "McZ_y/D"); 192 t_->Branch("McZ_Eta", &McZ_eta_, "McZ_Eta/D"); 193 t_->Branch("McZ_Phi", &McZ_phi_, "McZ_Phi/D"); 194 t_->Branch("McZ_rescPt", &McZ_rescPt_, "McZ_Pt/D"); 195 t_->Branch("McZ_rescY", &McZ_rescY_, "McZ_rescY/D"); 196 t_->Branch("McZ_rescEta", &McZ_rescEta_, "McZ_Eta/D"); 197 t_->Branch("McZ_rescPhi", &McZ_rescPhi_, "McZ_Phi/D"); 198 t_->Branch("McElec_nZmum", &McElec_nZmum_, "McElec_nZmum/I"); 199 t_->Branch("McElec_eta", &McElec_eta_, "McElec_eta[4]/D"); 200 t_->Branch("McElec_pt", &McElec_pt_, "McElec_pt[4]/D"); 201 t_->Branch("McElec_phi", &McElec_phi_, "McElec_phi[4]/D"); 202 t_->Branch("McElec_rescEta", &McElec_rescEta_, "McElec_rescEta[4]/D"); 203 t_->Branch("McElec_rescPhi", &McElec_rescPhi_, "McElec_rescPhi[4]/D"); 204 t_->Branch("McElec_rescPt", &McElec_rescPt_, "McElec_rescPt[4]/D"); 205 t_->Branch("McProbe_eta", &McProbe_eta_, "McProbe_eta[4]/D"); 206 t_->Branch("McProbe_pt", &McProbe_pt_, "McProbe_pt[4]/D"); 207 t_->Branch("McProbe_phi", &McProbe_phi_, "McProbe_phi[4]/D"); 208 t_->Branch("McProbe_rescEta", &McProbe_rescEta_, "McProbe_rescEta[4]/D"); 209 t_->Branch("McProbe_rescPt", &McProbe_rescPt_, "McProbe_rescPt[4]/D"); 210 t_->Branch("McProbe_rescPhi", &McProbe_rescPhi_, "McProbe_rescPhi[4]/D"); 211 t_->Branch("McElecProbe_dPhi", &McElecProbe_dPhi_, "McElecProbe_dPhi/D"); 212 t_->Branch("McElecProbe_dEta", &McElecProbe_dEta_, "McElecProbe_dEta/D"); 213 t_->Branch("McElecProbe_dR", &McElecProbe_dR_, "McElecProbe_dR/D"); 214 } 215 216} 217 218void ErsatzMEt::analyze(const edm::Event& evt, const edm::EventSetup& es) 219{ 220 221 es.get<CaloGeometryRecord>().get(geoHandle_); 222 es.get<CaloTopologyRecord>().get(pTopology_); 223 224 edm::LogDebug_("","", 151)<<"Initialising variables."; 225 nTags_ = -99; nProbes_ = -99; 226 CaloMEt_ = -99.; CaloMEtphi_ = -99.; 227 T1MEt_ = -99.; T1MEtphi_ = -99.; 228 PfMEt_ = -99.; PfMEtphi_ = -99.; 229 TcMEt_ = -99.; TcMEtphi_ = -99.; 230 if(Zevent_) 231 { 232 McZ_m_ = -99.; McZ_pt_ = -99.; McZ_y_ = -99.; McZ_eta_ = -99.; McZ_phi_ = -99.; 233 McZ_rescM_ = -99.; McZ_rescPt_ = -99.; McZ_rescY_ = -99.; McZ_rescEta_ = -99.; McZ_rescPhi_ = -99.; 234 McElec_nZmum_ = -99; McElec_nFinal_ = -99; 235 } 236 237 for(int i = 0; i < nEntries_arr_; ++i) 238 { 239 tag_q_[i] = -99; 240 tag_pt_[i] = -99.; tag_eta_[i] = -99.; tag_phi_[i] = -99.; 241 tag_rescPt_[i] = -99.; tag_rescEta_[i] = -99.; tag_rescPhi_[i] = -99.; 242 tag_trckIso_[i] = -99.; tag_ecalIso_[i] = -99.; tag_hcalIso_[i] = -99.; 243 tag_sIhIh_[i] = -99.; tag_dPhiIn_[i] = -99.; tag_dEtaIn_[i] = -99.; 244 tag_e5x5_[i] = -99.; tag_e2x5Max_[i] = -99.; tag_e1x5Max_[i] = -99.; 245 tag_hoe_[i] = -99.; tag_eop_[i] = -99.; tag_pin_[i] = -99.; tag_pout_[i] = -99.; 246 247 probe_q_[i] = -99; 248 probe_pt_[i] = -99.; probe_eta_[i] = -99.; probe_phi_[i] = -99.; 249 probe_rescPt_[i] = -99.; probe_rescEta_[i] = -99.; probe_rescPhi_[i] = -99.; 250 probe_trckIso_[i] = -99.; probe_ecalIso_[i] = -99.; probe_hcalIso_[i] = -99.; 251 probe_sIhIh_[i] = -99.; probe_dPhiIn_[i] = -99.; probe_dEtaIn_[i] = -99.; 252 probe_e5x5_[i] = -99.; probe_e2x5Max_[i] = -99.; probe_e1x5Max_[i] = -99.; 253 probe_hoe_[i] = -99.; probe_eop_[i] = -99.; probe_pin_[i] = -99.; probe_pout_[i] = -99.; 254 255 Z_pt_[i] = -99.; Z_y_[i] = -99.; Z_eta_[i] = -99.; Z_phi_[i] = -99.; Z_m_[i] = -99.; 256 Z_rescPt_[i] = -99.; Z_rescY_[i] = -99.; Z_rescEta_[i] = -99.; Z_rescPhi_[i] = -99.; Z_rescM_[i] = -99.; Z_probe_dPhi_[i] = -99.; 257 258 ErsatzV1_Mesc_[i] = -99.; ErsatzV1_rescMesc_[i] = -99.; 259 ErsatzV2_Mesc_[i] = -99.; ErsatzV2_rescMesc_[i] = -99.; 260 ErsatzV3_Mesc_[i] = -99.; ErsatzV3_rescMesc_[i] = -99.; 261 ErsatzV4_Mesc_[i] = -99.; ErsatzV4_rescMesc_[i] = -99.; 262 ErsatzV1CaloMEt_[i] = -99.; ErsatzV1CaloMt_[i] = -99.; ErsatzV1CaloMEtPhi_[i] = -99.; 263 ErsatzV2CaloMEt_[i] = -99.; ErsatzV2CaloMt_[i] = -99.; ErsatzV2CaloMEtPhi_[i] = -99.; 264 ErsatzV3CaloMEt_[i] = -99.; ErsatzV3CaloMt_[i] = -99.; ErsatzV3CaloMEtPhi_[i] = -99.; 265 ErsatzV4CaloMEt_[i] = -99.; ErsatzV4CaloMt_[i] = -99.; ErsatzV4CaloMEtPhi_[i] = -99.; 266 ErsatzV1T1MEt_[i] = -99.; ErsatzV1T1Mt_[i] = -99.; ErsatzV1T1MEtPhi_[i] = -99.; 267 ErsatzV1PfMEt_[i] = -99.; ErsatzV1PfMt_[i] = -99.; ErsatzV1PfMEtPhi_[i] = -99.; 268 ErsatzV1TcMEt_[i] = -99.; ErsatzV1TcMt_[i] = -99.; ErsatzV1TcMEtPhi_[i] = -99.; 269 270 probe_sc_pt_[i] = -99.; probe_sc_eta_[i] = -99.; probe_sc_phi_[i] = -99.; 271 probe_sc_E_[i] = -99.; probe_sc_rawE_[i] = -99.; 272 probe_scV2_E_[i] = -99.; 273 probe_scV3_E_[i] = -99.; 274 probe_scV4_E_[i] = -99.; 275 276 if(Zevent_) 277 { 278 McElec_pt_[i] = -99.; McElec_eta_[i] = -99.; McElec_phi_[i] = -99.; 279 McElec_rescPt_[i] = -99.; McElec_rescEta_[i] = -99.; McElec_rescPhi_[i] = -99.; 280 McProbe_pt_[i] = -99.; McProbe_eta_[i] = -99.; McProbe_phi_[i] = -99.; 281 McProbe_rescPt_[i] = -99.; McProbe_rescEta_[i] = -99.; McProbe_rescPhi_[i] = -99.; 282 McElecProbe_dPhi_[i] = -99.; McElecProbe_dEta_[i] = -99.; McElecProbe_dR_[i] = -99.; 283 } 284 285 edm::LogDebug_("","",180)<<"Initialisation of array index "<< i <<" completed."; 286 } 287 //Get Collections 288 edm::Handle<reco::GenParticleCollection> pGenPart; 289 try 290 { 291 evt.getByLabel(MCTruthCollection_, pGenPart); 292 }catch(cms::Exception &ex) 293 { 294 edm::LogError("analyze") <<"Can't get collection with label "<< MCTruthCollection_.label(); 295 } 296 edm::Handle<reco::GsfElectronCollection> pElectrons; 297 try 298 { 299 evt.getByLabel(ElectronCollection_, pElectrons); 300 }catch(cms::Exception &ex) 301 { 302 edm::LogError("analyze") <<"Can't get collection with label "<< ElectronCollection_.label(); 303 } 304 edm::Handle<reco::SuperClusterCollection> pHybrid; 305 try 306 { 307 evt.getByLabel(HybridScCollection_, pHybrid); 308 }catch(cms::Exception &ex) 309 { 310 edm::LogError("analyze") <<"Can't get collection with label "<< HybridScCollection_.label(); 311 } 312 edm::Handle<reco::SuperClusterCollection> pM5x5; 313 try 314 { 315 evt.getByLabel(M5x5ScCollection_, pM5x5); 316 }catch(cms::Exception &ex) 317 { 318 edm::LogError("analyze") <<"Can't get collection with label "<< M5x5ScCollection_.label(); 319 } 320 edm::Handle<reco::CaloMETCollection> pCaloMEt; 321 try 322 { 323 evt.getByLabel(CaloMEtCollection_, pCaloMEt); 324 }catch(cms::Exception &ex) 325 { 326 edm::LogError("analyze")<<"Can't get collection with label "<< CaloMEtCollection_.label(); 327 } 328 edm::Handle<reco::METCollection> pT1MEt; 329// try 330// { 331// evt.getByLabel(T1MEtCollection_, pT1MEt); 332// }catch(cms::Exception &ex) 333// { 334// edm::LogError("analyze")<<"Can't get collection with label "<< T1MEtCollection_.label(); 335// } 336 edm::Handle<reco::PFMETCollection> pPfMEt; 337 try 338 { 339 evt.getByLabel(PfMEtCollection_, pPfMEt); 340 }catch(cms::Exception &ex) 341 { 342 edm::LogError("analyze")<<"Can't get collection with label "<< PfMEtCollection_.label(); 343 } 344 edm::Handle<reco::METCollection> pTcMEt; 345 try 346 { 347 evt.getByLabel(TcMEtCollection_, pTcMEt); 348 }catch(cms::Exception &ex) 349 { 350 edm::LogError("analyze")<<"Can't get collection with label "<< TcMEtCollection_.label(); 351 } 352 edm::Handle<reco::GenMETCollection> pGenMEt; 353 try 354 { 355 evt.getByLabel(GenMEtCollection_, pGenMEt); 356 }catch(cms::Exception &ex) 357 { 358 edm::LogError("analyze")<<"Can't get collection with label "<< GenMEtCollection_.label(); 359 } 360 edm::Handle<edm::TriggerResults> pTriggerResults; 361 try 362 { 363 evt.getByLabel(TriggerResults_, pTriggerResults); 364 }catch(cms::Exception &ex) 365 { 366 edm::LogError("analyze")<<"Cant get collection with label "<< TriggerResults_.label(); 367 } 368 edm::Handle<trigger::TriggerEvent> pHLT; 369 try 370 { 371 evt.getByLabel(TriggerEvent_, pHLT); 372 }catch(cms::Exception &ex) 373 { 374 edm::LogError("analyze")<<"Can't get collection with label "<< TriggerEvent_.label(); 375 } 376 377 std::vector<math::XYZTLorentzVector>McElecs,McElecsFinalState; 378 std::vector<math::XYZTLorentzVector> McElecsResc; 379 if(Zevent_) 380 { 381 edm::LogDebug_("","",289)<<"Analysing MC properties."; 382 const reco::GenParticleCollection *McCand = pGenPart.product(); 383 math::XYZTLorentzVector Zboson, RescZboson, McElec1, McElec2; 384 for(reco::GenParticleCollection::const_iterator McP = McCand->begin(); McP != McCand->end(); ++McP) 385 { 386 const reco::Candidate* mum = McP->mother(); 387 if(abs(McP->pdgId())==11 && abs(mum->pdgId()) == 23) 388 { 389 McElecs.push_back(McP->p4()); 390 if(abs(mum->pdgId() == 23)) Zboson = mum->p4(); 391 392 std::cout <<"Found electron, ID = "<< McP->pdgId() <<"\t status = "<< McP->status()<<std::endl; 393 if(McP->status() != 1) 394 { 395 const reco::Candidate* McPD = McP->daughter(0); 396 McPD = McPD->mother(); 397 while(McPD->status() != 1) 398 { 399 int n = McPD->numberOfDaughters(); 400 std::cout<< McPD->pdgId() << " : status = "<<McPD->status() 401 <<"\tNumber of Daughters = "<< n <<std::endl; 402 for(int j = 0; j < n; ++ j) 403 { 404 const reco::Candidate *d = McPD->daughter( j ); 405 std::cout <<"Daughter "<< j <<"\t id = "<< d->pdgId() << std::endl; 406 if(abs(d->pdgId()) == 11) 407 { 408 McPD = d; 409 break; 410 } 411 } 412 } 413 std::cout<< McPD->pdgId() << " : status = "<<McPD->status()<<"\tAdding to vector!"<<std::endl; 414 McElecsFinalState.push_back(McPD->p4()); 415 }else McElecsFinalState.push_back(McP->p4()); 416 } 417 } 418 McZ_m_ = Zboson.M(); McZ_pt_ = Zboson.Pt(); McZ_phi_ = Zboson.Phi(); McZ_eta_ = Zboson.Eta(); McZ_y_ = Zboson.Rapidity(); 419 McElec_nZmum_ =McElecs.size(); 420 McElec_nFinal_ =McElecsFinalState.size(); 421 edm::LogDebug_("","",309)<<"MC electrons with Z mother = "<< McElec_nZmum_ 422 <<"\tFinal state MC electrons = "<< McElec_nFinal_; 423 424 McElecsResc.resize(2); 425// RescZboson.SetCoordinates(Zboson.Px(), Zboson.Py(), Zboson.Pz(), sqrt(Zboson.P2()+(mW_*mW_*Zboson.M2())/(mZ_*mZ_))); 426 RescZboson.SetCoordinates(Zboson.Px()*mW_/mZ_, Zboson.Py()*mW_/mZ_, Zboson.Pz()*mW_/mZ_, Zboson.E()*mW_/mZ_); 427 McZ_rescM_ = RescZboson.M(); McZ_rescPt_ = RescZboson.Pt(); McZ_rescEta_ = RescZboson.Eta(); McZ_rescPhi_ = RescZboson.Phi(); 428 McZ_rescY_ = RescZboson.Rapidity(); 429 ROOT::Math::Boost CoMBoost(Zboson.BoostToCM()); 430 431 math::XYZTLorentzVector RescMcElec0 = CoMBoost(McElecsFinalState[0]); 432 math::XYZTLorentzVector RescMcElec1 = CoMBoost(McElecsFinalState[1]); 433 RescMcElec0 *= mW_/mZ_; 434 RescMcElec1 *= mW_/mZ_; 435 436 double E_W = RescZboson.E(); 437 ROOT::Math::Boost BackToLab(RescZboson.Px()/E_W, RescZboson.Py()/E_W, RescZboson.Pz()/E_W); 438 439 RescMcElec0 = BackToLab(RescMcElec0); 440// RndmMcElec_Rescaled_pt_ = RescMcElec0.Pt(); 441// RndmMcElec_Rescaled_eta_ = RescMcElec0.Eta(); 442// RndmMcElec_Rescaled_phi_ = RescMcElec0.Phi(); 443 444 RescMcElec1 = BackToLab(RescMcElec1); 445// OthrMcElec_Rescaled_pt_ = RescMcElec1.Pt(); 446// OthrMcElec_Rescaled_eta_ = RescMcElec1.Eta(); 447// OthrMcElec_Rescaled_phi_ = RescMcElec1.Phi(); 448 McElecsResc[0] = RescMcElec0; 449 McElecsResc[1] = RescMcElec1; 450 math::XYZTLorentzVector sum = RescMcElec1+RescMcElec0; 451 edm::LogDebug_("","", 307)<<"McElecsResc[0] + McElecsResc[1] = ("<<sum.Px()<<", "<<sum.Py()<<", " 452 <<sum.Pz()<<", "<<sum.E()<<")"; 453 } 454 455 const edm::TriggerResults* HltRes = pTriggerResults.product(); 456 const edm::TriggerNames & triggerNames = evt.triggerNames(*HltRes); 457 if(HLTPathCheck_) 458 { 459 for(unsigned int itrig = 0; itrig < HltRes->size(); ++itrig) 460 { 461 std::string nom = triggerNames.triggerName(itrig); 462 edm::LogInfo("")<< itrig <<" : Name = "<< nom <<"\t Accepted = "<< HltRes->accept(itrig); 463 } 464 } 465 if(HltRes->accept(34) ==0) edm::LogError("")<<"Event did not pass "<< triggerNames.triggerName(34)<<"!"; 466 if(HltRes->accept(34) !=0) 467 { 468 std::vector<reco::GsfElectronRef> UniqueElectrons; 469 UniqueElectrons = uniqueElectronFinder(pElectrons); 470 edm::LogDebug_("","ErsatzMEt",192)<<"Unique electron size = "<<UniqueElectrons.size(); 471 std::vector<reco::GsfElectronRef> SelectedElectrons; 472 const unsigned int fId = pHLT->filterIndex(TriggerPath_); 473 std::cout << "Filter Id = " << fId << std::endl; 474 SelectedElectrons = electronSelector(UniqueElectrons, pHLT, fId, CutVector_); 475 nTags_ = SelectedElectrons.size(); 476 edm::LogDebug_("","ErsatzMEt",197)<<"Selected electron size = "<<nTags_; 477 478 iComb_ = 0; 479 if(Zevent_) 480 { 481 //Match MC electrons to the selected electrons and store some of their properties in the tree. 482 //The properties of the other MC electron (i.e. that not selected) are also stored. 483 for(std::vector<reco::GsfElectronRef>::const_iterator elec = SelectedElectrons.begin(); 484 elec != SelectedElectrons.end(); ++elec) 485 { 486 for(int m = 0; m < 2; ++m) 487 { 488 double dRLimit = 99.; 489 double dR = reco::deltaR(McElecs[m], *(*elec)); 490 if(dR < dRLimit) 491 { 492 dRLimit = dR; 493 McElec_pt_[iComb_] = McElecs[m].pt(); 494 McElec_eta_[iComb_] = McElecs[m].eta(); 495 McElec_rescPt_[iComb_] = McElecsResc[m].pt(); 496 McElec_rescEta_[iComb_] = McElecsResc[m].eta(); 497 } 498 } 499 } 500 } 501 502 std::map<reco::GsfElectronRef, reco::GsfElectronRef> TagProbePairs; 503 TagProbePairs = probeFinder(SelectedElectrons, pElectrons); 504 nProbes_ = TagProbePairs.size(); 505 edm::LogDebug_("", "ErsatzMEt", 209)<<"Number of tag-probe pairs = "<< TagProbePairs.size(); 506 507 if(!TagProbePairs.empty()) 508 { 509 const reco::CaloMETCollection* caloMEtCollection = pCaloMEt.product(); 510 const reco::MET calomet = *(caloMEtCollection->begin()); 511 CaloMEt_ = calomet.pt(); 512 CaloMEtphi_ = calomet.phi(); 513 514 //const reco::METCollection* t1MEtCollection = pT1MEt.product(); 515 //const reco::MET t1met = *(t1MEtCollection->begin()); 516 //T1MEt_ = t1met.pt(); 517 //T1MEtphi_ = t1met.phi(); 518 519 const reco::PFMETCollection* pfMEtCollection = pPfMEt.product(); 520 const reco::MET pfmet = *(pfMEtCollection->begin()); 521 PfMEt_ = pfmet.pt(); 522 PfMEtphi_ = pfmet.phi(); 523 524 const reco::METCollection* tcMEtCollection = pTcMEt.product(); 525 const reco::MET tcmet = *(tcMEtCollection->begin()); 526 TcMEt_ = tcmet.pt(); 527 TcMEtphi_ = tcmet.phi(); 528 529 reco::MET ersatzMEt; 530 531 for(std::map<reco::GsfElectronRef, reco::GsfElectronRef>::const_iterator it = TagProbePairs.begin(); 532 it != TagProbePairs.end(); ++it) 533 { 534 edm::LogDebug_("","DelendumLoop", 293)<<"iComb_ = "<< iComb_; 535 tag_q_[iComb_] = it->first->charge(); 536 edm::LogDebug_("","",360)<<"tag charge = "<< tag_q_[iComb_]; 537 tag_pt_[iComb_] = it->first->pt(); 538 tag_eta_[iComb_] = it->first->eta(); 539 tag_phi_[iComb_] = it->first->phi(); 540 edm::LogDebug_("","ErsatzMEt", 364)<<"tag pt = "<< tag_pt_[iComb_] 541 <<"\teta = "<< tag_eta_[iComb_]<<"\tphi = "<< tag_phi_[iComb_]; 542 tag_trckIso_[iComb_] = it->first->isolationVariables03().tkSumPt; 543 tag_ecalIso_[iComb_] = it->first->isolationVariables04().ecalRecHitSumEt; 544 tag_hcalIso_[iComb_] = it->first->isolationVariables04().hcalDepth1TowerSumEt 545 + it->first->isolationVariables04().hcalDepth2TowerSumEt; 546 edm::LogDebug_("","ErsatzMEt", 370)<<"tag trackiso = "<< tag_trckIso_[iComb_] 547 <<"\tecaliso = "<< tag_ecalIso_[iComb_]<<"\thcaliso = "<< tag_hcalIso_[iComb_]; 548 tag_sIhIh_[iComb_] = it->first->scSigmaIEtaIEta(); 549 tag_dPhiIn_[iComb_] = it->first->deltaPhiSuperClusterTrackAtVtx(); 550 tag_dEtaIn_[iComb_] = it->first->deltaEtaSuperClusterTrackAtVtx(); 551 edm::LogDebug_("","ErsatzMEt", 245)<<"tag sIhIh = "<< tag_sIhIh_[iComb_] 552 <<"\tdPhiIn = "<< tag_dPhiIn_[iComb_]<<"\tdEtaIn = "<< tag_dEtaIn_[iComb_]; 553 tag_e5x5_[iComb_] = it->first->scE5x5(); 554 tag_e2x5Max_[iComb_] = it->first->scE2x5Max(); 555 tag_e2x5Max_[iComb_] = it->first->scE1x5(); 556 edm::LogDebug_("","ErsatzMEt", 245)<<"tag e5x5 = "<< tag_e5x5_[iComb_] 557 <<"\te2x5Max = "<< tag_e2x5Max_[iComb_]<<"\te1x5Max = "<< tag_e1x5Max_[iComb_]; 558 tag_hoe_[iComb_] = it->first->hadronicOverEm(); 559 tag_eop_[iComb_] = it->first->eSuperClusterOverP(); 560 tag_pin_[iComb_] = it->first->trackMomentumAtVtx().R(); 561 tag_pout_[iComb_] = it->first->trackMomentumOut().R(); 562 edm::LogDebug_("","ErsatzMEt", 245)<<"tag hoe = "<<tag_hoe_[iComb_]<<"\tpoe = "<<tag_eop_[iComb_] 563 <<"\tpin = "<< tag_pin_[iComb_]<<"\tpout = "<< tag_pout_[iComb_]; 564 probe_q_[iComb_] = it->first->charge(); 565 edm::LogDebug_("","",360)<<"probe charge = "<< probe_q_[iComb_]; 566 probe_pt_[iComb_] = it->second->pt(); 567 probe_eta_[iComb_] = it->second->eta(); 568 probe_phi_[iComb_] = it->second->phi(); 569 edm::LogDebug_("","ErsatzMEt", 245)<<"probe pt = "<< probe_pt_[iComb_] 570 <<"\teta = "<< probe_eta_[iComb_]<<"\tphi = "<< probe_phi_[iComb_]; 571 probe_trckIso_[iComb_] = it->second->isolationVariables03().tkSumPt; 572 probe_ecalIso_[iComb_] = it->second->isolationVariables04().ecalRecHitSumEt; 573 probe_hcalIso_[iComb_] = it->second->isolationVariables04().hcalDepth1TowerSumEt 574 + it->second->isolationVariables04().hcalDepth2TowerSumEt; 575 edm::LogDebug_("","ErsatzMEt", 245)<<"probe trackiso = "<< probe_trckIso_[iComb_] 576 <<"\tecaliso = "<< probe_ecalIso_[iComb_]<<"\thcaliso = "<< probe_phi_[iComb_]; 577 probe_sIhIh_[iComb_] = it->second->scSigmaIEtaIEta(); 578 probe_dPhiIn_[iComb_] = it->second->deltaPhiSuperClusterTrackAtVtx(); 579 probe_dEtaIn_[iComb_] = it->second->deltaEtaSuperClusterTrackAtVtx(); 580 edm::LogDebug_("","ErsatzMEt", 245)<<"probe sIhIh = "<< probe_sIhIh_[iComb_] 581 <<"\tdPhiIn = "<< probe_dPhiIn_[iComb_]<<"\tdEtaIn = "<< probe_dEtaIn_[iComb_]; 582 probe_e5x5_[iComb_] = it->second->scE5x5(); 583 probe_e2x5Max_[iComb_] = it->second->scE2x5Max(); 584 probe_e2x5Max_[iComb_] = it->second->scE1x5(); 585 edm::LogDebug_("","ErsatzMEt", 245)<<"probe e5x5 = "<< probe_e5x5_[iComb_] 586 <<"\te2x5Max = "<< probe_e2x5Max_[iComb_]<<"\te1x5Max = "<< probe_e1x5Max_[iComb_]; 587 probe_hoe_[iComb_] = it->second->hadronicOverEm(); 588 probe_eop_[iComb_] = it->second->eSuperClusterOverP(); 589 probe_pin_[iComb_] = it->second->trackMomentumAtVtx().R(); 590 probe_pout_[iComb_] = it->second->trackMomentumOut().R(); 591 edm::LogDebug_("","ErsatzMEt", 245)<<"probe hoe = "<<probe_hoe_[iComb_]<<"\tpoe = "<<probe_eop_[iComb_] 592 <<"\tpin = "<< probe_pin_[iComb_]<<"\tpout = "<< probe_pout_[iComb_]; 593 594 double dRLimit = 0.2; 595 for(unsigned int mcEId = 0; mcEId < McElecs.size(); ++mcEId) 596 { 597// double dR = reco::deltaR((*(*mcEl)), probeVec); 598 double dR = reco::deltaR(McElecs[mcEId], it->second->p4()); 599 if(dR < dRLimit) 600 { 601 dRLimit = dR; 602 McProbe_pt_[iComb_] = McElecs[mcEId].pt(); 603 McProbe_eta_[iComb_] = McElecs[mcEId].eta(); 604 McProbe_phi_[iComb_] = McElecs[mcEId].phi(); 605 McProbe_rescPt_[iComb_] = McElecsResc[mcEId].pt(); 606 McProbe_rescEta_[iComb_] = McElecsResc[mcEId].eta(); 607 McProbe_rescPhi_[iComb_] = McElecsResc[mcEId].phi(); 608 probe_d_MCE_SCE_[iComb_] = McElecs[mcEId].energy() - it->second->superCluster()->rawEnergy(); 609 McElecProbe_dPhi_[iComb_] = reco::deltaPhi(McElecs[mcEId].phi(), McElecs[(mcEId+1)%2].phi()); 610 McElecProbe_dEta_[iComb_] = fabs(McElecs[mcEId].eta() - McElecs[(mcEId+1)%2].eta()); 611 McElecProbe_dR_[iComb_] = reco::deltaR(McElecs[mcEId], McElecs[(mcEId+1)%2]); 612 } 613 } 614 615 // Uncorrected supercluster V1 616 reco::SuperCluster scV1 = *(it->second->superCluster()); 617 math::XYZTLorentzVector probe_scV1_detVec = DetectorVector(scV1); 618 probe_sc_pt_[iComb_] = probe_scV1_detVec.pt(); 619 probe_sc_eta_[iComb_] = scV1.eta(); 620 probe_sc_phi_[iComb_] = scV1.phi(); 621 probe_sc_nClus_[iComb_] = scV1.clustersSize(); 622 probe_sc_E_[iComb_] = scV1.energy(); 623 probe_sc_rawE_[iComb_] = scV1.rawEnergy(); 624 625 ersatzMEt = ersatzFabrik(it->first, scV1, calomet, 1); 626 ErsatzV1CaloMEt_[iComb_] = ersatzMEt.pt(); 627 ErsatzV1CaloMEtPhi_[iComb_] = ersatzMEt.phi(); 628 //ersatzMEt = ersatzFabrik(it->first, it->second, t1met); 629 //ErsatzV1T1MEt_[iComb_] = ersatzMEt.pt(); 630 //ErsatzV1T1MEtPhi_[iComb_] = ersatzMEt.phi(); 631 ersatzMEt = ersatzFabrik(it->first, it->second, pfmet); 632 ErsatzV1PfMEt_[iComb_] = ersatzMEt.pt(); 633 ErsatzV1PfMEtPhi_[iComb_] = ersatzMEt.phi(); 634 ersatzMEt = ersatzFabrik(it->first, it->second, tcmet); 635 ErsatzV1TcMEt_[iComb_] = ersatzMEt.pt(); 636 ErsatzV1TcMEtPhi_[iComb_] = ersatzMEt.phi(); 637 638 // fEta corrected supercluster V2 639 reco::SuperCluster scV2; 640 if(fabs(probe_sc_eta_[iComb_]) < 1.479) 641 { 642 scV2 = fEtaScCorr(scV1); 643 }else{ 644 scV2 = scV1; 645 } 646 probe_scV2_E_[iComb_] = scV2.energy(); 647 ersatzMEt = ersatzFabrik(it->first, scV2, calomet, 2); 648 ErsatzV2CaloMEt_[iComb_] = ersatzMEt.pt(); 649 ErsatzV2CaloMEtPhi_[iComb_] = ersatzMEt.phi(); 650 651 // fBrem corrected supercluster V3 652 reco::SuperCluster scV3; 653 if(fabs(probe_sc_eta_[iComb_]) < 1.479) 654 { 655 scV3 = fBremScCorr(scV1, hyb_fCorrPSet_); 656 }else{ 657 scV3 = fBremScCorr(scV1, m5x5_fCorrPSet_); 658 } 659 probe_scV3_E_[iComb_] = scV3.energy(); 660 ersatzMEt = ersatzFabrik(it->first, scV3, calomet, 3); 661 ErsatzV3CaloMEt_[iComb_] = ersatzMEt.pt(); 662 ErsatzV3CaloMEtPhi_[iComb_] = ersatzMEt.phi(); 663 664 // Fully corrected supercluster V4 665 reco::SuperCluster scV4; 666 if(fabs(probe_sc_eta_[iComb_]) < 1.479) 667 { 668 scV4 = fBremScCorr(scV1, hyb_fCorrPSet_); 669 }else{ 670 scV4 = fBremScCorr(scV1, m5x5_fCorrPSet_); 671 } 672 probe_scV4_E_[iComb_] = scV4.energy(); 673 ersatzMEt = ersatzFabrik(it->first, scV4, calomet, 4); 674 ErsatzV4CaloMEt_[iComb_] = ersatzMEt.pt(); 675 ErsatzV4CaloMEtPhi_[iComb_] = ersatzMEt.phi(); 676 677 ++iComb_; 678 } 679 t_->Fill(); 680 } 681 } 682} 683 684std::map<reco::GsfElectronRef, reco::GsfElectronRef> ErsatzMEt::probeFinder(const std::vector<reco::GsfElectronRef>& tags, 685 const edm::Handle<reco::GsfElectronCollection> pElectrons) 686{ 687 const reco::GsfElectronCollection *probeCands = pElectrons.product(); 688 std::map<reco::GsfElectronRef, reco::GsfElectronRef> TagProbes; 689 for(std::vector<reco::GsfElectronRef>::const_iterator tagelec = tags.begin(); tagelec != tags.end(); ++tagelec) 690 { 691 reco::GsfElectronRef tag = *tagelec; 692 std::pair<reco::GsfElectronRef, reco::GsfElectronRef> TagProbePair; 693 int nProbesPerTag = 0; 694 int index = 0; 695 for(reco::GsfElectronCollection::const_iterator probeelec = probeCands->begin(); probeelec != probeCands->end(); ++probeelec) 696 { 697 reco::GsfElectronRef probe(pElectrons, index); 698 double probeScEta = probe->superCluster()->eta(); 699 if(probe->superCluster() != tag->superCluster() && fabs(probeScEta) < 2.5) 700 { 701 if(fabs(probeScEta) < 1.4442 || fabs(probeScEta) > 1.560) 702 { 703 double invmass = ROOT::Math::VectorUtil::InvariantMass(tag->p4(), probe->p4()); 704 if(mTPmin_ <= invmass && invmass <= mTPmax_) 705 { 706 TagProbePair = std::make_pair(tag, probe); 707 ++nProbesPerTag; 708 } 709 } 710 } 711 ++index; 712 } 713 //nGsfElectrons_ = index; 714 if(nProbesPerTag == 1) TagProbes.insert(TagProbePair); 715 } 716 return TagProbes; 717} 718 719reco::MET ErsatzMEt::ersatzFabrik(const reco::GsfElectronRef& elec, 720 const reco::SuperCluster& sc, 721 const reco::MET& met, 722 const int corr) 723{ 724 const math::XYZPoint ZVertex(elec->TrackPositionAtVtx().X(), elec->TrackPositionAtVtx().Y(),elec->TrackPositionAtVtx().Z()); 725 726 math::XYZTLorentzVector nu, boost_nu, ele, boost_ele; 727 reco::SuperCluster elecSc = *(elec->superCluster()); 728 nu = PhysicsVectorRaw(met.vertex(), sc); 729 boost_nu = PhysicsVectorRaw(ZVertex, sc); 730 ele = PhysicsVectorRaw(met.vertex(), elecSc); 731 boost_ele = ele; 732 733 //Should use reco vertex for best Z->ee measurement. 734 edm::LogDebug_("ersatzFabrikV1", "", 569)<<"elec = ("<< elec->p4().Px() << ", "<< elec->p4().Py()<< ", "<< elec->p4().Pz() << ", "<< elec->p4().E()<<")"; 735 math::XYZTLorentzVector Zboson = boost_nu + elec->p4(); 736 edm::LogDebug_("ersatzFabrikV1", "", 569)<<"Z pt = "<< Zboson.Pt() << "Z boson mass = " << Zboson.M(); 737 edm::LogDebug_("ersatzFabrikV1","", 570)<<"Z boson in lab frame = ("<<Zboson.Px()<<", "<<Zboson.Py()<<", " 738 <<Zboson.Pz()<<", "<<Zboson.E()<<")"; 739 math::XYZTLorentzVector RescZboson(Zboson.Px(), Zboson.Py(), Zboson.Pz(), sqrt(Zboson.P2()+(mW_*mW_*Zboson.M2())/(mZ_*mZ_))); 740 edm::LogDebug_("ersatzFabrikV1","", 573)<<"W boson in lab frame = ("<<RescZboson.Px()<<", "<<RescZboson.Py()<<", " 741 <<RescZboson.Pz()<<", "<<RescZboson.E()<<")"; 742 ROOT::Math::Boost BoostToZRestFrame(Zboson.BoostToCM()); 743 edm::LogDebug_("ersatzFabrikV1","", 576)<<"Electron in lab frame = ("<< boost_ele.Px()<<", "<< boost_ele.Py()<<", " 744 << boost_ele.Pz()<<", "<< boost_ele.E()<<")"; 745 edm::LogDebug_("ersatzFabrikV1","", 578)<<"Ersatz Neutrino in lab frame = ("<< boost_nu.Px()<<", "<< boost_nu.Py()<<", " 746 << boost_nu.Pz()<<", "<< boost_nu.E()<<")"; 747 boost_ele = BoostToZRestFrame(boost_ele); 748 boost_nu = BoostToZRestFrame(boost_nu); 749 edm::LogDebug_("ersatzFabrikV1","", 582)<<"Electron in Z rest frame = ("<<boost_ele.Px()<<", "<<boost_ele.Py()<<", " 750 <<boost_ele.Pz()<<", "<<boost_ele.E()<<")"; 751 edm::LogDebug_("ersatzFabrikV1","", 584)<<"Ersatz Neutrino in Z rest frame = ("<<boost_nu.Px()<<", "<<boost_nu.Py()<<", " 752 <<boost_nu.Pz()<<", "<<boost_nu.E()<<")"; 753 boost_ele *= mW_/mZ_; 754 boost_nu *= mW_/mZ_; 755 756 double E_W = RescZboson.E(); 757 ROOT::Math::Boost BackToLab(RescZboson.Px()/E_W, RescZboson.Py()/E_W, RescZboson.Pz()/E_W); 758 math::XYZTLorentzVector metVec(-99999., -99999., -99., -99999.); 759 boost_ele = BackToLab(boost_ele); 760 761 boost_nu = BackToLab(boost_nu); 762 math::XYZTLorentzVector sum = boost_nu+boost_ele; 763 edm::LogDebug_("ersatzFabrikV1","", 597)<<"Electron back in lab frame = ("<<boost_ele.Px()<<", "<<boost_ele.Py()<<", " 764 <<boost_ele.Pz()<<", "<<boost_ele.E()<<")"; 765 edm::LogDebug_("ersatzFabrikV1","", 599)<<"Ersatz Neutrino back in lab frame = ("<<boost_nu.Px()<<", "<<boost_nu.Py()<<", " 766 <<boost_nu.Pz()<<", "<<boost_nu.E()<<")"; 767 edm::LogDebug_("ersatzFabrikV1","", 601)<<"boost_ele + boost_nu = ("<<sum.Px()<<", "<<sum.Py()<<", " 768 <<sum.Pz()<<", "<<sum.E()<<")"; 769 770 nu.SetXYZT(nu.X(), nu.Y(), 0., nu.T()); 771 ele.SetXYZT(ele.X(), ele.Y(), 0., ele.T()); 772 boost_ele.SetXYZT(boost_ele.X(), boost_ele.Y(), 0., boost_ele.T()); 773 metVec = met.p4() + nu + ele - boost_ele; 774 775 reco::MET ersatzMEt(metVec, met.vertex()); 776 if (corr == 1) 777 { 778 //Z_caloV1_m_[iComb_] = Zboson.M(); 779 //Z_caloV1_pt_[iComb_] = Zboson.Pt(); 780 //Z_caloV1_y_[iComb_] = Zboson.Y(); 781 //Z_caloV1_eta_[iComb_] = Zboson.Eta(); 782 //Z_caloV1_phi_[iComb_] = Zboson.Phi(); 783 //Z_caloV1_rescM_[iComb_] = RescZboson.M(); 784 //Z_caloV1_rescPt_[iComb_] = RescZboson.Pt(); 785 //Z_caloV1_rescY_[iComb_] = RescZboson.Y(); 786 //Z_caloV1_rescEta_[iComb_] = RescZboson.Eta(); 787 //Z_caloV1_rescPhi_[iComb_] = RescZboson.Phi(); 788 //Z_caloV1_probe_dPhi_[iComb_] = reco::deltaPhi(Zboson.Phi(), elec->phi()); 789 //tag_caloV1_rescPt_[iComb_] = boost_ele.Pt(); 790 //tag_caloV1_rescEta_[iComb_] = boost_ele.Eta(); 791 //tag_caloV1_rescPhi_[iComb_] = boost_ele.Phi(); 792 //probe_caloV1_rescPt_[iComb_] = boost_nu.Pt(); 793 //probe_caloV1_rescEta_[iComb_] = boost_nu.Eta(); 794 //probe_caloV1_rescPhi_[iComb_] = boost_nu.Phi(); 795 ErsatzV1_Mesc_[iComb_] = ROOT::Math::VectorUtil::InvariantMass(elec->p4(), boost_nu); 796 ErsatzV1_rescMesc_[iComb_] = ROOT::Math::VectorUtil::InvariantMass(ele, nu); 797 ErsatzV1CaloMt_[iComb_] = sqrt(2.*boost_ele.Pt()*ersatzMEt.pt()* 798 (1-cos(reco::deltaPhi(boost_ele.Phi(), ersatzMEt.phi())))); 799 } 800 if (corr == 2) 801 { 802 //Z_caloV2_m_[iComb_] = Zboson.M(); 803 //Z_caloV2_pt_[iComb_] = Zboson.Pt(); 804 //Z_caloV2_y_[iComb_] = Zboson.Y(); 805 //Z_caloV2_eta_[iComb_] = Zboson.Eta(); 806 //Z_caloV2_phi_[iComb_] = Zboson.Phi(); 807 //Z_caloV2_rescM_[iComb_] = RescZboson.M(); 808 //Z_caloV2_rescPt_[iComb_] = RescZboson.Pt(); 809 //Z_caloV2_rescY_[iComb_] = RescZboson.Y(); 810 //Z_caloV2_rescEta_[iComb_] = RescZboson.Eta(); 811 //Z_caloV2_rescPhi_[iComb_] = RescZboson.Phi(); 812 //Z_caloV2_probe_dPhi_[iComb_] = reco::deltaPhi(Zboson.Phi(), boost_elec->phi()); 813 //tag_caloV2_rescPt_[iComb_] = boost_ele.Pt(); 814 //tag_caloV2_rescEta_[iComb_] = boost_ele.Eta(); 815 //tag_caloV2_rescPhi_[iComb_] = boost_ele.Phi(); 816 //probe_caloV2_rescPt_[iComb_] = boost_nu.Pt(); 817 //probe_caloV2_rescEta_[iComb_] = boost_nu.Eta(); 818 //probe_caloV2_rescPhi_[iComb_] = boost_nu.Phi(); 819 ErsatzV2_Mesc_[iComb_] = ROOT::Math::VectorUtil::InvariantMass(elec->p4(), boost_nu); 820 ErsatzV2_rescMesc_[iComb_] = ROOT::Math::VectorUtil::InvariantMass(ele, nu); 821 ErsatzV2CaloMt_[iComb_] = sqrt(2.*boost_ele.Pt()*ersatzMEt.pt()* 822 (1-cos(reco::deltaPhi(boost_ele.Phi(), ersatzMEt.phi())))); 823 } 824 if (corr == 3) 825 { 826 //Z_caloV3_m_[iComb_] = Zboson.M(); 827 //Z_caloV3_pt_[iComb_] = Zboson.Pt(); 828 //Z_caloV3_y_[iComb_] = Zboson.Y(); 829 //Z_caloV3_eta_[iComb_] = Zboson.Eta(); 830 //Z_caloV3_phi_[iComb_] = Zboson.Phi(); 831 //Z_caloV3_rescM_[iComb_] = RescZboson.M(); 832 //Z_caloV3_rescPt_[iComb_] = RescZboson.Pt(); 833 //Z_caloV3_rescY_[iComb_] = RescZboson.Y(); 834 //Z_caloV3_rescEta_[iComb_] = RescZboson.Eta(); 835 //Z_caloV3_rescPhi_[iComb_] = RescZboson.Phi(); 836 //Z_caloV3_probe_dPhi_[iComb_] = reco::deltaPhi(Zboson.Phi(), boost_elec->phi()); 837 //tag_caloV3_rescPt_[iComb_] = boost_ele.Pt(); 838 //tag_caloV3_rescEta_[iComb_] = boost_ele.Eta(); 839 //tag_caloV3_rescPhi_[iComb_] = boost_ele.Phi(); 840 //probe_caloV3_rescPt_[iComb_] = boost_nu.Pt(); 841 //probe_caloV3_rescEta_[iComb_] = boost_nu.Eta(); 842 //probe_caloV3_rescPhi_[iComb_] = boost_nu.Phi(); 843 ErsatzV3_Mesc_[iComb_] = ROOT::Math::VectorUtil::InvariantMass(elec->p4(), boost_nu); 844 ErsatzV3_rescMesc_[iComb_] = ROOT::Math::VectorUtil::InvariantMass(ele, nu); 845 ErsatzV3CaloMt_[iComb_] = sqrt(2.*boost_ele.Pt()*ersatzMEt.pt()* 846 (1-cos(reco::deltaPhi(boost_ele.Phi(), ersatzMEt.phi())))); 847 } 848 if (corr == 4) 849 { 850 //Z_caloV4_m_[iComb_] = Zboson.M(); 851 //Z_caloV4_pt_[iComb_] = Zboson.Pt(); 852 //Z_caloV4_y_[iComb_] = Zboson.Y(); 853 //Z_caloV4_eta_[iComb_] = Zboson.Eta(); 854 //Z_caloV4_phi_[iComb_] = Zboson.Phi(); 855 //Z_caloV4_rescM_[iComb_] = RescZboson.M(); 856 //Z_caloV4_rescPt_[iComb_] = RescZboson.Pt(); 857 //Z_caloV4_rescY_[iComb_] = RescZboson.Y(); 858 //Z_caloV4_rescEta_[iComb_] = RescZboson.Eta(); 859 //Z_caloV4_rescPhi_[iComb_] = RescZboson.Phi(); 860 //Z_caloV4_probe_dPhi_[iComb_] = reco::deltaPhi(Zboson.Phi(), boost_elec->phi()); 861 //tag_caloV4_rescPt_[iComb_] = boost_ele.Pt(); 862 //tag_caloV4_rescEta_[iComb_] = boost_ele.Eta(); 863 //tag_caloV4_rescPhi_[iComb_] = boost_ele.Phi(); 864 //probe_caloV4_rescPt_[iComb_] = boost_nu.Pt(); 865 //probe_caloV4_rescEta_[iComb_] = boost_nu.Eta(); 866 //probe_caloV4_rescPhi_[iComb_] = boost_nu.Phi(); 867 ErsatzV4_Mesc_[iComb_] = ROOT::Math::VectorUtil::InvariantMass(elec->p4(), boost_nu); 868 ErsatzV4_rescMesc_[iComb_] = ROOT::Math::VectorUtil::InvariantMass(ele, nu); 869 ErsatzV4CaloMt_[iComb_] = sqrt(2.*boost_ele.Pt()*ersatzMEt.pt()* 870 (1-cos(reco::deltaPhi(boost_ele.Phi(), ersatzMEt.phi())))); 871 } 872 return ersatzMEt; 873} 874 875reco::MET ErsatzMEt::ersatzFabrik(const reco::GsfElectronRef& tag, 876 const reco::GsfElectronRef& probe, 877 const reco::MET& met) 878{ 879 math::XYZTLorentzVector elec, nu, boost_elec, boost_nu; 880 boost_elec = tag->p4(); 881 edm::LogDebug_("ersatzFabrikV1", "", 858)<<"boost_elec = ("<< boost_elec.Px() << ", "<< boost_elec.Py()<< ", "<< boost_elec.Pz() << ", "<< boost_elec.E()<<")"; 882 boost_nu = probe->p4(); 883 edm::LogDebug_("ersatzFabrikV1", "", 860)<<"boost_nu = ("<< boost_nu.Px() << ", "<< boost_nu.Py()<< ", "<< boost_nu.Pz() << ", "<< boost_nu.E()<<")"; 884 math::XYZTLorentzVector Zboson = boost_elec + boost_nu; 885 edm::LogDebug_("ersatzFabrikV1", "", 862)<<"Zboson = ("<< Zboson.Px() << ", "<< Zboson.Py()<< ", "<< Zboson.Pz() << ", "<< Zboson.E()<<")"; 886 math::XYZTLorentzVector RescZboson(Zboson.Px(), Zboson.Py(), Zboson.Pz(), sqrt(Zboson.P2()+(mW_*mW_*Zboson.M2())/(mZ_*mZ_))); 887 edm::LogDebug_("ersatzFabrikV1", "", 864)<<"RescZboson = ("<< RescZboson.Px() << ", "<< RescZboson.Py()<< ", "<< RescZboson.Pz() << ", "<< RescZboson.E()<<")"; 888 ROOT::Math::Boost BoostToZRestFrame(Zboson.BoostToCM()); 889 elec = BoostToZRestFrame(boost_elec); 890 edm::LogDebug_("ersatzFabrikV1", "", 867)<<"boost_elec (in Z rest frame) = ("<< elec.Px() << ", "<< elec.Py()<< ", "<< elec.Pz() << ", "<< elec.E()<<")"; 891 nu = BoostToZRestFrame(boost_nu); 892 edm::LogDebug_("ersatzFabrikV1", "", 869)<<"boost_nu (in Z rest frame) = ("<< nu.Px() << ", "<< nu.Py()<< ", "<< nu.Pz() << ", "<< nu.E()<<")"; 893 elec *= mW_/mZ_; 894 edm::LogDebug_("ersatzFabrikV1", "", 871)<<"elec (in Z rest frame) = ("<< elec.Px() << ", "<< elec.Py()<< ", "<< elec.Pz() << ", "<< elec.E()<<")"; 895 nu *= mW_/mZ_; 896 edm::LogDebug_("ersatzFabrikV1", "", 873)<<"nu (in Z rest frame) = ("<< nu.Px() << ", "<< nu.Py()<< ", "<< nu.Pz() << ", "<< nu.E()<<")"; 897 ROOT::Math::Boost BoostBackToLab(RescZboson.Px()/RescZboson.E(), RescZboson.Py()/RescZboson.E(), RescZboson.Pz()/RescZboson.E()); 898 math::XYZTLorentzVector metVec(-99999., -99999., -99., -99999.); 899 elec = BoostBackToLab(elec); 900 edm::LogDebug_("ersatzFabrikV1", "", 877)<<"elec = ("<< elec.Px() << ", "<< elec.Py()<< ", "<< elec.Pz() << ", "<< elec.E()<<")"; 901 nu = BoostBackToLab(nu); 902 edm::LogDebug_("ersatzFabrikV1", "", 879)<<"nu = ("<< nu.Px() << ", "<< nu.Py()<< ", "<< nu.Pz() << ", "<< nu.E()<<")"; 903 Z_m_[iComb_] = Zboson.M(); 904 Z_pt_[iComb_] = Zboson.Pt(); 905 Z_y_[iComb_] = Zboson.Y(); 906 Z_eta_[iComb_] = Zboson.Eta(); 907 Z_phi_[iComb_] = Zboson.Phi(); 908 Z_rescM_[iComb_] = RescZboson.M(); 909 Z_rescPt_[iComb_] = RescZboson.Pt(); 910 Z_rescY_[iComb_] = RescZboson.Y(); 911 Z_rescEta_[iComb_] = RescZboson.Eta(); 912 Z_rescPhi_[iComb_] = RescZboson.Phi(); 913 Z_probe_dPhi_[iComb_] = reco::deltaPhi(Zboson.Phi(), boost_elec.phi()); 914 tag_rescPt_[iComb_] = elec.Pt(); 915 tag_rescEta_[iComb_] = elec.Eta(); 916 tag_rescPhi_[iComb_] = elec.Phi(); 917 probe_rescPt_[iComb_] = nu.Pt(); 918 probe_rescEta_[iComb_] = nu.Eta(); 919 probe_rescPhi_[iComb_] = nu.Phi(); 920 elec.SetXYZT(elec.X(), elec.Y(), 0., elec.T()); 921 nu.SetXYZT(nu.X(), nu.Y(), 0., nu.T()); 922 boost_elec.SetXYZT(boost_elec.X(), boost_elec.Y(), 0., boost_elec.T()); 923 metVec = met.p4() + nu + elec - boost_elec; 924 reco::MET ersatzMEt(metVec, met.vertex()); 925 return ersatzMEt; 926} 927 928bool ErsatzMEt::isInBarrel(double eta) 929{ 930 return (fabs(eta) < BarrelEtaMax_); 931} 932 933bool ErsatzMEt::isInEndCap(double eta) 934{ 935 return (fabs(eta) < EndCapEtaMax_ && fabs(eta) > EndCapEtaMin_); 936} 937 938bool ErsatzMEt::isInFiducial(double eta) 939{ 940 return isInBarrel(eta) || isInEndCap(eta); 941} 942 943// ------------ method called once each job just after ending the event loop ------------ 944void ErsatzMEt::endJob() { 945} 946//define this as a plug-in 947DEFINE_FWK_MODULE(ErsatzMEt); 948