Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-02-14 14:31:22

0001 //
0002 // Author:  Jan Heyninck
0003 // Created: Tue Apr  3 17:33:23 PDT 2007
0004 //
0005 //
0006 #include "TopQuarkAnalysis/TopJetCombination/interface/TtSemiLRJetCombObservables.h"
0007 #include "DataFormats/PatCandidates/interface/Jet.h"
0008 #include "DataFormats/Math/interface/deltaR.h"
0009 #include "DataFormats/Math/interface/deltaPhi.h"
0010 #include "DataFormats/Common/interface/Handle.h"
0011 #include "AnalysisDataFormats/TopObjects/interface/TtGenEvent.h"
0012 
0013 // constructor with path; default should not be used
0014 TtSemiLRJetCombObservables::TtSemiLRJetCombObservables(edm::ConsumesCollector &&iC,
0015                                                        const edm::EDGetTokenT<std::vector<pat::Jet> > &jetSourceToken)
0016     : jetSourceToken_(jetSourceToken), genEvtToken_(iC.consumes<TtGenEvent>(edm::InputTag("genEvt"))) {}
0017 
0018 // destructor
0019 TtSemiLRJetCombObservables::~TtSemiLRJetCombObservables() {}
0020 
0021 std::vector<TtSemiLRJetCombObservables::IntBoolPair> TtSemiLRJetCombObservables::operator()(TtSemiEvtSolution &solution,
0022                                                                                             const edm::Event &iEvent,
0023                                                                                             bool matchOnly) {
0024   bool debug = false;
0025 
0026   evtselectVarVal.clear();
0027   evtselectVarMatch.clear();
0028 
0029   // Check whether the objects are matched:
0030   //  bool matchHadt = false; //commented out since gcc461 complained that this variable was set but unused
0031   //  bool matchLept = false; //commented out since gcc461 complained that this variable was set but unused
0032   //  bool matchHadW = false; //commented out since gcc461 complained that this variable was set but unused
0033   //  bool matchLepW = false; //commented out since gcc461 complained that this variable was set but unused
0034   bool matchHadb = false;
0035   bool matchLepb = false;
0036   bool matchHadp = false;
0037   bool matchHadq = false;
0038   bool matchHadpq = false;
0039   bool matchHadqp = false;
0040   bool matchLepl = false;
0041   bool matchLepn = false;
0042 
0043   if (debug)
0044     std::cout << "== start matching the objects " << std::endl;
0045 
0046   double drLepl = 0, drLepn = 0, drHadb = 0, drLepb = 0, drHadp = 0, drHadq = 0, drHadpq = 0,
0047          drHadqp = 0;  // drHadt=0, drLept=0, drLepW=0, drHadW=0;
0048 
0049   edm::Handle<TtGenEvent> genEvent;
0050   iEvent.getByToken(genEvtToken_, genEvent);
0051 
0052   if (genEvent.failedToGet()) {
0053     if (debug)
0054       std::cout << "== not found genEvent " << std::endl;
0055   } else if (genEvent.isValid()) {
0056     if (debug)
0057       std::cout << "== found genEvent " << std::endl;
0058 
0059     if (genEvent->isSemiLeptonic() && genEvent->numberOfBQuarks() == 2) {
0060       //if(debug) std::cout << "== genEvent->quarkFromAntiTop() " << genEvent->quarkFromAntiTop()->pt() << std::endl;
0061       if (debug)
0062         std::cout << "== genEvent->isSemiLeptonic() && genEvent->numberOfBQuarks() == 2 " << std::endl;
0063       if (debug)
0064         std::cout << "== solution.getDecay() == " << solution.getDecay() << std::endl;
0065       if (debug)
0066         std::cout << "== solution.getRecLepm().pt() = " << solution.getRecLepm().pt() << std::endl;
0067       //if(debug) if(solution.getGenLepl() == 0) std::cout << "solution.getGenLepl() == NULL" << std::endl;
0068       if (debug)
0069         std::cout << "== *(solution.getGenLept())" << solution.getGenLept()->pt() << std::endl;
0070       if (debug)
0071         std::cout << "== *(solution.getGenLepl())" << solution.getGenLepl()->pt() << std::endl;
0072       // std::cout << "Semilepton:\n";
0073       // Match the lepton by deltaR
0074       if (solution.getDecay() == "muon")
0075         drLepl = DeltaR<pat::Particle, reco::GenParticle>()(solution.getRecLepm(), *(solution.getGenLepl()));
0076       if (debug)
0077         std::cout << "== matching lepton " << std::endl;
0078       if (solution.getDecay() == "electron")
0079         drLepl = DeltaR<pat::Particle, reco::GenParticle>()(solution.getRecLepe(), *(solution.getGenLepl()));
0080       matchLepl = (drLepl < 0.3);
0081 
0082       if (debug)
0083         std::cout << "== lepton is matched " << std::endl;
0084       // Match the neutrino by deltaR
0085       drLepn = DeltaR<pat::Particle, reco::GenParticle>()(solution.getRecLepn(), *(solution.getGenLepn()));
0086       matchLepn = (drLepn < 0.3);
0087 
0088       // Match the hadronic b by deltaR
0089       drHadb = DeltaR<pat::Jet, reco::GenParticle>()(solution.getRecHadb(), *(solution.getGenHadb()));
0090       matchHadb = (drHadb < 0.3);
0091 
0092       // Match the hadronicleptonic b by deltaR
0093       drLepb = DeltaR<pat::Jet, reco::GenParticle>()(solution.getRecLepb(), *(solution.getGenLepb()));
0094       matchLepb = (drLepb < 0.3);
0095 
0096       // Match the hadronic p by deltaR
0097       drHadp = DeltaR<pat::Jet, reco::GenParticle>()(solution.getRecHadp(), *(solution.getGenHadp()));
0098       matchHadp = (drHadp < 0.3);
0099 
0100       // Match the hadronic pq by deltaR
0101       drHadpq = DeltaR<pat::Jet, reco::GenParticle>()(solution.getRecHadp(), *(solution.getGenHadq()));
0102       matchHadpq = (drHadpq < 0.3);
0103 
0104       // Match the hadronic q by deltaR
0105       drHadq = DeltaR<pat::Jet, reco::GenParticle>()(solution.getRecHadq(), *(solution.getGenHadq()));
0106       matchHadq = (drHadq < 0.3);
0107 
0108       // Match the hadronic qp by deltaR
0109       drHadqp = DeltaR<pat::Jet, reco::GenParticle>()(solution.getRecHadq(), *(solution.getGenHadp()));
0110       matchHadqp = (drHadqp < 0.3);
0111 
0112       //commented out since gcc461 complained that these variables were set but unused!
0113       /*
0114       // Match the hadronic W by deltaR
0115       drHadW = DeltaR<reco::Particle, reco::GenParticle>()(solution.getRecHadW(), *(solution.getGenHadW()));
0116       matchHadW = (drHadW < 0.3);
0117 
0118       // Match the leptonic W by deltaR
0119       drLepW = DeltaR<reco::Particle, reco::GenParticle>()(solution.getRecLepW(), *(solution.getGenLepW()));
0120       matchLepW = (drLepW < 0.3);
0121 
0122       // Match the hadronic t by deltaR
0123       drHadt = DeltaR<reco::Particle, reco::GenParticle>()(solution.getRecHadt(), *(solution.getGenHadt()));
0124       matchHadt = (drHadt < 0.3);
0125 
0126       // Match the leptonic t by deltaR
0127       drLept = DeltaR<reco::Particle, reco::GenParticle>()(solution.getRecLept(), *(solution.getGenLept()));
0128       matchLept = (drLept < 0.3);
0129       */
0130     }
0131   } else {
0132     if (debug)
0133       std::cout << "== nothing " << std::endl;
0134   }
0135 
0136   if (debug)
0137     std::cout << "== objects matched" << std::endl;
0138 
0139   edm::Handle<std::vector<pat::Jet> > jets;
0140   iEvent.getByToken(jetSourceToken_, jets);
0141 
0142   if (debug)
0143     std::cout << "== start calculating observables" << std::endl;
0144 
0145   //obs1 : pt(had top)
0146   double AverageTop = ((solution.getHadb().p4() + solution.getHadq().p4() + solution.getHadp().p4()).pt() +
0147                        (solution.getLepb().p4() + solution.getHadq().p4() + solution.getHadp().p4()).pt() +
0148                        (solution.getHadb().p4() + solution.getLepb().p4() + solution.getHadp().p4()).pt() +
0149                        (solution.getHadb().p4() + solution.getHadq().p4() + solution.getLepb().p4()).pt()) /
0150                       4.;
0151   double Obs1 = ((solution.getHadb().p4() + solution.getHadq().p4() + solution.getHadp().p4()).pt()) / AverageTop;
0152   evtselectVarVal.push_back(IntDblPair(1, Obs1));
0153   evtselectVarMatch.push_back(IntBoolPair(1, ((matchHadq && matchHadp) || (matchHadpq && matchHadqp)) && matchHadb));
0154 
0155   if (debug)
0156     std::cout << "== observable 1 " << Obs1 << std::endl;
0157 
0158   //obs2 : (pt_b1 + pt_b2)/(sum jetpt)
0159   double Obs2 =
0160       (solution.getHadb().pt() + solution.getLepb().pt()) / (solution.getHadp().pt() + solution.getHadq().pt());
0161   evtselectVarVal.push_back(IntDblPair(2, Obs2));
0162   evtselectVarMatch.push_back(
0163       IntBoolPair(2, ((matchHadp && matchHadq) || (matchHadpq && matchHadqp)) && matchHadb && matchLepb));
0164 
0165   if (debug)
0166     std::cout << "== observable 2 " << Obs2 << std::endl;
0167 
0168   //obs3: delta R between lep b and lepton
0169   double Obs3 = -999;
0170   if (solution.getDecay() == "muon")
0171     Obs3 = ROOT::Math::VectorUtil::DeltaR(solution.getLepb().p4(), solution.getRecLepm().p4());
0172   if (solution.getDecay() == "electron")
0173     Obs3 = ROOT::Math::VectorUtil::DeltaR(solution.getLepb().p4(), solution.getRecLepe().p4());
0174   evtselectVarVal.push_back(IntDblPair(3, Obs3));
0175   evtselectVarMatch.push_back(IntBoolPair(3, matchLepb && matchLepl));
0176 
0177   if (debug)
0178     std::cout << "== observable 3 " << Obs3 << std::endl;
0179 
0180   //obs4 : del R ( had b, had W)
0181   double Obs4 =
0182       ROOT::Math::VectorUtil::DeltaR(solution.getHadb().p4(), solution.getHadq().p4() + solution.getHadp().p4());
0183   evtselectVarVal.push_back(IntDblPair(4, Obs4));
0184   evtselectVarMatch.push_back(IntBoolPair(4, matchHadb && ((matchHadp && matchHadp) || (matchHadpq && matchHadqp))));
0185 
0186   if (debug)
0187     std::cout << "== observable 4 " << Obs4 << std::endl;
0188 
0189   //obs5 : del R between light quarkssolution.getHadp().phi(
0190   double Obs5 = ROOT::Math::VectorUtil::DeltaR(solution.getHadq().p4(), solution.getHadp().p4());
0191   evtselectVarVal.push_back(IntDblPair(5, Obs5));
0192   evtselectVarMatch.push_back(IntBoolPair(5, (matchHadp && matchHadq) || (matchHadpq && matchHadqp)));
0193 
0194   if (debug)
0195     std::cout << "== observable 5 " << Obs5 << std::endl;
0196 
0197   //obs6 : b-tagging information
0198   double Obs6 = 0;
0199   if (fabs(solution.getHadb().bDiscriminator("trackCountingJetTags") + 10) < 0.0001 ||
0200       fabs(solution.getLepb().bDiscriminator("trackCountingJetTags") + 10) < 0.0001) {
0201     Obs6 = -10.;
0202   } else {
0203     Obs6 = (solution.getHadb().bDiscriminator("trackCountingJetTags") +
0204             solution.getLepb().bDiscriminator("trackCountingJetTags"));
0205   }
0206   evtselectVarVal.push_back(IntDblPair(6, Obs6));
0207   evtselectVarMatch.push_back(IntBoolPair(6, 1));
0208 
0209   if (debug)
0210     std::cout << "== observable 6 " << Obs6 << std::endl;
0211 
0212   //obs7 : chi2 value of kinematical fit with W-mass constraint
0213   double Obs7 = 0;
0214   if (solution.getProbChi2() < 0) {
0215     Obs7 = -0;
0216   } else {
0217     Obs7 = log10(solution.getProbChi2() + .00001);
0218   }
0219   evtselectVarVal.push_back(IntDblPair(7, Obs7));
0220   evtselectVarMatch.push_back(IntBoolPair(7, ((matchHadp && matchHadq) || (matchHadpq && matchHadqp))));
0221 
0222   if (debug)
0223     std::cout << "== observable 7 " << Obs7 << std::endl;
0224 
0225   //obs8(=7+1)
0226   double Obs8 = solution.getCalHadt().p4().pt();
0227   evtselectVarVal.push_back(IntDblPair(8, Obs8));
0228   evtselectVarMatch.push_back(IntBoolPair(8, matchHadb && ((matchHadp && matchHadq) || (matchHadpq && matchHadqp))));
0229 
0230   if (debug)
0231     std::cout << "== observable 8 " << Obs8 << std::endl;
0232 
0233   //obs9
0234   double Obs9 = fabs(solution.getCalHadt().p4().eta());
0235   evtselectVarVal.push_back(IntDblPair(9, Obs9));
0236   evtselectVarMatch.push_back(IntBoolPair(9, matchHadb && ((matchHadp && matchHadq) || (matchHadpq && matchHadqp))));
0237 
0238   if (debug)
0239     std::cout << "== observable 9 " << Obs9 << std::endl;
0240 
0241   //obs10
0242   double Obs10 = solution.getCalHadt().p4().theta();
0243   evtselectVarVal.push_back(IntDblPair(10, Obs10));
0244   evtselectVarMatch.push_back(IntBoolPair(10, matchHadb && ((matchHadp && matchHadq) || (matchHadpq && matchHadqp))));
0245 
0246   if (debug)
0247     std::cout << "== observable 10 " << Obs10 << std::endl;
0248 
0249   //obs11
0250   double Obs11 = solution.getCalHadW().p4().pt();
0251   evtselectVarVal.push_back(IntDblPair(11, Obs11));
0252   evtselectVarMatch.push_back(IntBoolPair(11, (matchHadp && matchHadq) || (matchHadpq && matchHadqp)));
0253 
0254   if (debug)
0255     std::cout << "== observable 11 " << Obs11 << std::endl;
0256 
0257   //obs12
0258   double Obs12 = fabs(solution.getCalHadW().p4().eta());
0259   evtselectVarVal.push_back(IntDblPair(12, Obs12));
0260   evtselectVarMatch.push_back(IntBoolPair(12, (matchHadp && matchHadq) || (matchHadpq && matchHadqp)));
0261 
0262   if (debug)
0263     std::cout << "== observable 12 " << Obs12 << std::endl;
0264 
0265   //obs13
0266   double Obs13 = solution.getCalHadW().p4().theta();
0267   evtselectVarVal.push_back(IntDblPair(13, Obs13));
0268   evtselectVarMatch.push_back(IntBoolPair(13, (matchHadp && matchHadq) || (matchHadpq && matchHadqp)));
0269 
0270   if (debug)
0271     std::cout << "== observable 13 " << Obs13 << std::endl;
0272 
0273   //obs14
0274   double Obs14 = solution.getCalHadb().p4().pt();
0275   evtselectVarVal.push_back(IntDblPair(14, Obs14));
0276   evtselectVarMatch.push_back(IntBoolPair(14, matchHadb));
0277 
0278   if (debug)
0279     std::cout << "== observable 14 " << Obs14 << std::endl;
0280 
0281   //obs15
0282   double Obs15 = fabs(solution.getCalHadb().p4().eta());
0283   evtselectVarVal.push_back(IntDblPair(15, Obs15));
0284   evtselectVarMatch.push_back(IntBoolPair(15, matchHadb));
0285 
0286   if (debug)
0287     std::cout << "== observable 15 " << Obs15 << std::endl;
0288 
0289   //obs16
0290   double Obs16 = solution.getCalHadb().p4().theta();
0291   evtselectVarVal.push_back(IntDblPair(16, Obs16));
0292   evtselectVarMatch.push_back(IntBoolPair(16, matchHadb));
0293 
0294   if (debug)
0295     std::cout << "== observable 16 " << Obs16 << std::endl;
0296 
0297   //obs17
0298   double Obs17 = solution.getCalHadp().p4().pt();
0299   evtselectVarVal.push_back(IntDblPair(17, Obs17));
0300   evtselectVarMatch.push_back(IntBoolPair(17, matchHadp));
0301 
0302   if (debug)
0303     std::cout << "== observable 17 " << Obs17 << std::endl;
0304 
0305   //obs18
0306   double Obs18 = fabs(solution.getCalHadp().p4().eta());
0307   evtselectVarVal.push_back(IntDblPair(18, Obs18));
0308   evtselectVarMatch.push_back(IntBoolPair(18, matchHadp));
0309 
0310   if (debug)
0311     std::cout << "== observable 18 " << Obs18 << std::endl;
0312 
0313   //obs19
0314   double Obs19 = solution.getCalHadp().p4().theta();
0315   evtselectVarVal.push_back(IntDblPair(19, Obs19));
0316   evtselectVarMatch.push_back(IntBoolPair(19, matchHadp));
0317 
0318   if (debug)
0319     std::cout << "== observable 19 " << Obs19 << std::endl;
0320 
0321   //obs20
0322   double Obs20 = solution.getCalHadq().p4().pt();
0323   evtselectVarVal.push_back(IntDblPair(20, Obs20));
0324   evtselectVarMatch.push_back(IntBoolPair(20, matchHadq));
0325 
0326   if (debug)
0327     std::cout << "== observable 20 " << Obs20 << std::endl;
0328 
0329   //obs21
0330   double Obs21 = fabs(solution.getCalHadq().p4().eta());
0331   evtselectVarVal.push_back(IntDblPair(21, Obs21));
0332   evtselectVarMatch.push_back(IntBoolPair(21, matchHadq));
0333 
0334   if (debug)
0335     std::cout << "== observable 21 " << Obs21 << std::endl;
0336 
0337   //obs22
0338   double Obs22 = solution.getCalHadq().p4().theta();
0339   evtselectVarVal.push_back(IntDblPair(22, Obs22));
0340   evtselectVarMatch.push_back(IntBoolPair(22, matchHadq));
0341 
0342   if (debug)
0343     std::cout << "== observable 22 " << Obs22 << std::endl;
0344 
0345   //obs23
0346   double Obs23 = solution.getCalLept().p4().pt();
0347   evtselectVarVal.push_back(IntDblPair(23, Obs23));
0348   evtselectVarMatch.push_back(IntBoolPair(23, matchLepl && matchLepn && matchLepb));
0349 
0350   if (debug)
0351     std::cout << "== observable 23 " << Obs23 << std::endl;
0352 
0353   //obs24
0354   double Obs24 = fabs(solution.getCalLept().p4().eta());
0355   evtselectVarVal.push_back(IntDblPair(24, Obs24));
0356   evtselectVarMatch.push_back(IntBoolPair(24, matchLepl && matchLepn && matchLepb));
0357 
0358   if (debug)
0359     std::cout << "== observable 24 " << Obs24 << std::endl;
0360 
0361   //obs25
0362   double Obs25 = solution.getCalLept().p4().theta();
0363   evtselectVarVal.push_back(IntDblPair(25, Obs25));
0364   evtselectVarMatch.push_back(IntBoolPair(25, matchLepl && matchLepn && matchLepb));
0365 
0366   if (debug)
0367     std::cout << "== observable 25 " << Obs25 << std::endl;
0368 
0369   //obs26
0370   double Obs26 = solution.getRecLepW().p4().pt();
0371   evtselectVarVal.push_back(IntDblPair(26, Obs26));
0372   evtselectVarMatch.push_back(IntBoolPair(26, matchLepl && matchLepn));
0373 
0374   if (debug)
0375     std::cout << "== observable 26 " << Obs26 << std::endl;
0376 
0377   //obs27
0378   double Obs27 = fabs(solution.getRecLepW().p4().eta());
0379   evtselectVarVal.push_back(IntDblPair(27, Obs27));
0380   evtselectVarMatch.push_back(IntBoolPair(27, matchLepl && matchLepn));
0381 
0382   if (debug)
0383     std::cout << "== observable 27 " << Obs27 << std::endl;
0384 
0385   //obs28
0386   double Obs28 = solution.getRecLepW().p4().theta();
0387   evtselectVarVal.push_back(IntDblPair(28, Obs28));
0388   evtselectVarMatch.push_back(IntBoolPair(28, matchLepl && matchLepn));
0389 
0390   if (debug)
0391     std::cout << "== observable 28 " << Obs28 << std::endl;
0392 
0393   //obs29
0394   double Obs29 = solution.getCalLepb().p4().pt();
0395   evtselectVarVal.push_back(IntDblPair(29, Obs29));
0396   evtselectVarMatch.push_back(IntBoolPair(29, matchLepb));
0397 
0398   if (debug)
0399     std::cout << "== observable 29 " << Obs29 << std::endl;
0400 
0401   //obs30
0402   double Obs30 = fabs(solution.getCalLepb().p4().eta());
0403   evtselectVarVal.push_back(IntDblPair(30, Obs30));
0404   evtselectVarMatch.push_back(IntBoolPair(30, matchLepb));
0405 
0406   if (debug)
0407     std::cout << "== observable 30 " << Obs30 << std::endl;
0408 
0409   //obs31
0410   double Obs31 = solution.getCalLepb().p4().theta();
0411   evtselectVarVal.push_back(IntDblPair(31, Obs31));
0412   evtselectVarMatch.push_back(IntBoolPair(31, matchLepb));
0413 
0414   if (debug)
0415     std::cout << "== observable 31 " << Obs31 << std::endl;
0416 
0417   //obs32
0418   double Obs32 = -999.;
0419   if (solution.getDecay() == "muon")
0420     Obs32 = solution.getRecLepm().p4().pt();
0421   if (solution.getDecay() == "electron")
0422     Obs32 = solution.getRecLepe().p4().pt();
0423   evtselectVarVal.push_back(IntDblPair(32, Obs32));
0424   evtselectVarMatch.push_back(IntBoolPair(32, matchLepl));
0425 
0426   if (debug)
0427     std::cout << "== observable 32 " << Obs32 << std::endl;
0428 
0429   //obs33
0430   double Obs33 = -999.;
0431   if (solution.getDecay() == "muon")
0432     Obs33 = fabs(solution.getRecLepm().p4().eta());
0433   if (solution.getDecay() == "electron")
0434     Obs33 = fabs(solution.getRecLepe().p4().eta());
0435   evtselectVarVal.push_back(IntDblPair(33, Obs33));
0436   evtselectVarMatch.push_back(IntBoolPair(33, matchLepl));
0437 
0438   if (debug)
0439     std::cout << "== observable 33 " << Obs33 << std::endl;
0440 
0441   //obs34
0442   double Obs34 = -999.;
0443   if (solution.getDecay() == "muon")
0444     Obs34 = fabs(solution.getRecLepm().p4().theta());
0445   if (solution.getDecay() == "electron")
0446     Obs34 = fabs(solution.getRecLepe().p4().theta());
0447   evtselectVarVal.push_back(IntDblPair(34, Obs34));
0448   evtselectVarMatch.push_back(IntBoolPair(34, matchLepl));
0449 
0450   if (debug)
0451     std::cout << "== observable 34 " << Obs34 << std::endl;
0452 
0453   //obs35
0454   double Obs35 = solution.getFitLepn().p4().pt();
0455   evtselectVarVal.push_back(IntDblPair(35, Obs35));
0456   evtselectVarMatch.push_back(IntBoolPair(35, matchLepn));
0457 
0458   if (debug)
0459     std::cout << "== observable 35 " << Obs35 << std::endl;
0460 
0461   //obs36
0462   double Obs36 = fabs(solution.getFitLepn().p4().eta());
0463   evtselectVarVal.push_back(IntDblPair(36, Obs36));
0464   evtselectVarMatch.push_back(IntBoolPair(36, matchLepn));
0465 
0466   if (debug)
0467     std::cout << "== observable 36 " << Obs36 << std::endl;
0468 
0469   //obs37
0470   double Obs37 = solution.getFitLepn().p4().theta();
0471   evtselectVarVal.push_back(IntDblPair(37, Obs37));
0472   evtselectVarMatch.push_back(IntBoolPair(37, matchLepn));
0473 
0474   if (debug)
0475     std::cout << "== observable 37 " << Obs37 << std::endl;
0476 
0477   // 2 particle kinematics
0478   //obs38
0479   double Obs38 = fabs(solution.getCalHadW().p4().phi() - solution.getRecLepW().p4().phi());
0480   if (Obs38 > 3.1415927)
0481     Obs38 = 2 * 3.1415927 - Obs31;
0482   if (Obs38 < -3.1415927)
0483     Obs38 = -2 * 3.1415927 - Obs31;
0484   evtselectVarVal.push_back(IntDblPair(38, Obs38));
0485   evtselectVarMatch.push_back(
0486       IntBoolPair(38, matchLepl && matchLepn && ((matchHadp && matchHadq) || (matchHadpq && matchHadqp))));
0487 
0488   if (debug)
0489     std::cout << "== observable 38 " << Obs38 << std::endl;
0490 
0491   //obs39
0492   double Obs39 = fabs(solution.getCalHadW().p4().eta() - solution.getRecLepW().p4().eta());
0493   evtselectVarVal.push_back(IntDblPair(39, Obs39));
0494   evtselectVarMatch.push_back(
0495       IntBoolPair(39, matchLepl && matchLepn && ((matchHadp && matchHadq) || (matchHadpq && matchHadqp))));
0496 
0497   if (debug)
0498     std::cout << "== observable 39 " << Obs39 << std::endl;
0499 
0500   //obs40
0501   double Obs40 = fabs(solution.getCalHadW().p4().theta() - solution.getRecLepW().p4().theta());
0502   evtselectVarVal.push_back(IntDblPair(40, Obs40));
0503   evtselectVarMatch.push_back(
0504       IntBoolPair(40, matchLepl && matchLepn && ((matchHadp && matchHadq) || (matchHadpq && matchHadqp))));
0505 
0506   if (debug)
0507     std::cout << "== observable 40 " << Obs40 << std::endl;
0508 
0509   //obs41
0510   double Obs41 = ROOT::Math::VectorUtil::DeltaR(solution.getCalHadW().p4(), solution.getRecLepW().p4());
0511   evtselectVarVal.push_back(IntDblPair(41, Obs41));
0512   evtselectVarMatch.push_back(
0513       IntBoolPair(41, matchLepl && matchLepn && ((matchHadp && matchHadq) || (matchHadpq && matchHadqp))));
0514 
0515   if (debug)
0516     std::cout << "== observable 41 " << Obs41 << std::endl;
0517 
0518   //obs42
0519   double Obs42 = fabs(solution.getCalHadb().p4().phi() - solution.getCalLepb().p4().phi());
0520   if (Obs42 > 3.1415927)
0521     Obs42 = 2 * 3.1415927 - Obs42;
0522   if (Obs42 < -3.1415927)
0523     Obs42 = -2 * 3.1415927 - Obs42;
0524   evtselectVarVal.push_back(IntDblPair(42, Obs42));
0525   evtselectVarMatch.push_back(IntBoolPair(42, matchHadb && matchLepb));
0526 
0527   if (debug)
0528     std::cout << "== observable 42 " << Obs42 << std::endl;
0529 
0530   //obs43
0531   double Obs43 = fabs(solution.getCalHadb().p4().eta() - solution.getCalLepb().p4().eta());
0532   evtselectVarVal.push_back(IntDblPair(43, Obs43));
0533   evtselectVarMatch.push_back(IntBoolPair(43, matchHadb && matchLepb));
0534 
0535   if (debug)
0536     std::cout << "== observable 43 " << Obs43 << std::endl;
0537 
0538   //obs44
0539   double Obs44 = fabs(solution.getCalHadb().p4().theta() - solution.getCalLepb().p4().theta());
0540   evtselectVarVal.push_back(IntDblPair(44, Obs44));
0541   evtselectVarMatch.push_back(IntBoolPair(44, matchHadb && matchLepb));
0542 
0543   if (debug)
0544     std::cout << "== observable 44 " << Obs44 << std::endl;
0545 
0546   //obs45
0547   double Obs45 = ROOT::Math::VectorUtil::DeltaR(solution.getCalHadb().p4(), solution.getCalLepb().p4());
0548   evtselectVarVal.push_back(IntDblPair(45, Obs45));
0549   evtselectVarMatch.push_back(IntBoolPair(45, matchHadb && matchLepb));
0550 
0551   if (debug)
0552     std::cout << "== observable 45 " << Obs45 << std::endl;
0553 
0554   //obs46
0555   double Obs46 = fabs(solution.getCalHadb().p4().phi() - solution.getCalHadW().p4().phi());
0556   if (Obs46 > 3.1415927)
0557     Obs46 = 2 * 3.1415927 - Obs46;
0558   if (Obs46 < -3.1415927)
0559     Obs46 = -2 * 3.1415927 - Obs46;
0560   evtselectVarVal.push_back(IntDblPair(46, Obs46));
0561   evtselectVarMatch.push_back(IntBoolPair(46, matchHadb && ((matchHadq && matchHadp) || (matchHadpq && matchHadqp))));
0562 
0563   if (debug)
0564     std::cout << "== observable 46 " << Obs46 << std::endl;
0565 
0566   //obs47
0567   double Obs47 = fabs(solution.getCalHadb().p4().eta() - solution.getCalHadW().p4().eta());
0568   evtselectVarVal.push_back(IntDblPair(47, Obs47));
0569   evtselectVarMatch.push_back(IntBoolPair(47, matchHadb && ((matchHadp && matchHadq) || (matchHadpq && matchHadqp))));
0570 
0571   if (debug)
0572     std::cout << "== observable 47 " << Obs47 << std::endl;
0573 
0574   //obs48
0575   double Obs48 = fabs(solution.getCalHadb().p4().theta() - solution.getCalHadW().p4().theta());
0576   evtselectVarVal.push_back(IntDblPair(48, Obs48));
0577   evtselectVarMatch.push_back(IntBoolPair(48, matchHadb && ((matchHadp && matchHadq) || (matchHadpq && matchHadqp))));
0578 
0579   if (debug)
0580     std::cout << "== observable 48 " << Obs48 << std::endl;
0581 
0582   //obs49
0583   double Obs49 = ROOT::Math::VectorUtil::DeltaR(solution.getCalHadb().p4(), solution.getCalHadW().p4());
0584   evtselectVarVal.push_back(IntDblPair(49, Obs49));
0585   evtselectVarMatch.push_back(IntBoolPair(49, matchHadb && ((matchHadp && matchHadq) || (matchHadpq && matchHadqp))));
0586 
0587   if (debug)
0588     std::cout << "== observable 49 " << Obs49 << std::endl;
0589 
0590   //obs50
0591   double Obs50 = fabs(solution.getCalLepb().p4().phi() - solution.getRecLepW().p4().phi());
0592   if (Obs50 > 3.1415927)
0593     Obs50 = 2 * 3.1415927 - Obs50;
0594   if (Obs50 < -3.1415927)
0595     Obs50 = -2 * 3.1415927 - Obs50;
0596   evtselectVarVal.push_back(IntDblPair(50, Obs50));
0597   evtselectVarMatch.push_back(IntBoolPair(50, matchLepb && matchLepn && matchLepl));
0598 
0599   if (debug)
0600     std::cout << "== observable 50 " << Obs50 << std::endl;
0601 
0602   //obs51
0603   double Obs51 = fabs(solution.getCalLepb().p4().eta() - solution.getRecLepW().p4().eta());
0604   evtselectVarVal.push_back(IntDblPair(51, Obs51));
0605   evtselectVarMatch.push_back(IntBoolPair(51, matchLepb && matchLepn && matchLepl));
0606 
0607   if (debug)
0608     std::cout << "== observable 51 " << Obs51 << std::endl;
0609 
0610   //obs52
0611   double Obs52 = fabs(solution.getCalLepb().p4().theta() - solution.getRecLepW().p4().theta());
0612   evtselectVarVal.push_back(IntDblPair(52, Obs52));
0613   evtselectVarMatch.push_back(IntBoolPair(52, matchLepb && matchLepn && matchLepl));
0614 
0615   if (debug)
0616     std::cout << "== observable 52 " << Obs52 << std::endl;
0617 
0618   //obs53
0619   double Obs53 = ROOT::Math::VectorUtil::DeltaR(solution.getCalLepb().p4(), solution.getRecLepW().p4());
0620   evtselectVarVal.push_back(IntDblPair(53, Obs53));
0621   evtselectVarMatch.push_back(IntBoolPair(53, matchLepb && matchLepn && matchLepl));
0622 
0623   if (debug)
0624     std::cout << "== observable 53 " << Obs53 << std::endl;
0625 
0626   //obs54
0627   double Obs54 = fabs(solution.getCalHadp().p4().phi() - solution.getCalHadq().p4().phi());
0628   if (Obs54 > 3.1415927)
0629     Obs54 = 2 * 3.1415927 - Obs54;
0630   if (Obs54 < -3.1415927)
0631     Obs54 = -2 * 3.1415927 - Obs54;
0632   evtselectVarVal.push_back(IntDblPair(54, Obs54));
0633   evtselectVarMatch.push_back(IntBoolPair(54, (matchHadp && matchHadq) || (matchHadpq && matchHadqp)));
0634 
0635   if (debug)
0636     std::cout << "== observable 54 " << Obs54 << std::endl;
0637 
0638   //obs55
0639   double Obs55 = fabs(solution.getCalHadp().p4().eta() - solution.getCalHadq().p4().eta());
0640   evtselectVarVal.push_back(IntDblPair(55, Obs55));
0641   evtselectVarMatch.push_back(IntBoolPair(55, (matchHadp && matchHadq) || (matchHadpq && matchHadqp)));
0642 
0643   if (debug)
0644     std::cout << "== observable 55 " << Obs55 << std::endl;
0645 
0646   //obs56
0647   double Obs56 = fabs(solution.getCalHadp().p4().theta() - solution.getCalHadq().p4().theta());
0648   evtselectVarVal.push_back(IntDblPair(56, Obs56));
0649   evtselectVarMatch.push_back(IntBoolPair(56, (matchHadp && matchHadq) || (matchHadpq && matchHadqp)));
0650 
0651   if (debug)
0652     std::cout << "== observable 56 " << Obs56 << std::endl;
0653 
0654   //obs57
0655   double Obs57 = ROOT::Math::VectorUtil::DeltaR(solution.getCalHadp().p4(), solution.getCalHadq().p4());
0656   evtselectVarVal.push_back(IntDblPair(57, Obs57));
0657   evtselectVarMatch.push_back(IntBoolPair(57, (matchHadp && matchHadq) || (matchHadpq && matchHadqp)));
0658 
0659   if (debug)
0660     std::cout << "== observable 57 " << Obs57 << std::endl;
0661 
0662   //obs58
0663   double Obs58 = -999.;
0664   if (solution.getDecay() == "muon")
0665     Obs58 = fabs(solution.getRecLepm().p4().phi() - solution.getRecLepn().p4().phi());
0666   if (solution.getDecay() == "electron")
0667     Obs58 = fabs(solution.getRecLepe().p4().phi() - solution.getRecLepn().p4().phi());
0668   if (Obs58 > 3.1415927)
0669     Obs58 = 2 * 3.1415927 - Obs58;
0670   if (Obs58 < -3.1415927)
0671     Obs58 = -2 * 3.1415927 - Obs58;
0672   evtselectVarVal.push_back(IntDblPair(58, Obs58));
0673   evtselectVarMatch.push_back(IntBoolPair(58, matchLepl && matchLepn));
0674 
0675   if (debug)
0676     std::cout << "== observable 58 " << Obs58 << std::endl;
0677 
0678   //obs59
0679   double Obs59 = -999.;
0680   if (solution.getDecay() == "muon")
0681     Obs59 = fabs(solution.getRecLepm().p4().eta() - solution.getRecLepn().p4().eta());
0682   if (solution.getDecay() == "electron")
0683     Obs59 = fabs(solution.getRecLepe().p4().eta() - solution.getRecLepn().p4().eta());
0684   evtselectVarVal.push_back(IntDblPair(59, Obs59));
0685   evtselectVarMatch.push_back(IntBoolPair(59, matchLepl && matchLepn));
0686 
0687   if (debug)
0688     std::cout << "== observable 59 " << Obs59 << std::endl;
0689 
0690   //obs60
0691   double Obs60 = -999.;
0692   if (solution.getDecay() == "muon")
0693     Obs60 = fabs(solution.getRecLepm().p4().theta() - solution.getRecLepn().p4().theta());
0694   if (solution.getDecay() == "electron")
0695     Obs60 = fabs(solution.getRecLepe().p4().theta() - solution.getRecLepn().p4().theta());
0696   evtselectVarVal.push_back(IntDblPair(60, Obs60));
0697   evtselectVarMatch.push_back(IntBoolPair(60, matchLepl && matchLepn));
0698 
0699   if (debug)
0700     std::cout << "== observable 60 " << Obs60 << std::endl;
0701 
0702   //obs61
0703   double Obs61 = -999.;
0704   if (solution.getDecay() == "muon")
0705     Obs61 = ROOT::Math::VectorUtil::DeltaR(solution.getRecLepm().p4(), solution.getRecLepn().p4());
0706   if (solution.getDecay() == "electron")
0707     Obs61 = ROOT::Math::VectorUtil::DeltaR(solution.getRecLepe().p4(), solution.getRecLepn().p4());
0708   evtselectVarVal.push_back(IntDblPair(61, Obs61));
0709   evtselectVarMatch.push_back(IntBoolPair(61, matchLepl && matchLepn));
0710 
0711   if (debug)
0712     std::cout << "== observable 61 " << Obs61 << std::endl;
0713 
0714   // miscellaneous event variables
0715 
0716   //obs62
0717   double Obs62 =
0718       ((jets->size() > 4 && (*jets)[3].p4().Et() > 0.00001) ? (*jets)[4].p4().Et() / (*jets)[3].p4().Et() : 1.0);
0719   //double Obs62 = 1;
0720   evtselectVarVal.push_back(IntDblPair(62, Obs62));
0721   evtselectVarMatch.push_back(
0722       IntBoolPair(62, ((matchHadp && matchHadq) || (matchHadpq && matchHadqp)) && matchHadb && matchLepb));
0723 
0724   if (debug)
0725     std::cout << "== observable 62 " << Obs62 << std::endl;
0726 
0727   float calJetsSumEt = 0;
0728   for (unsigned int i = 4; i < jets->size(); i++) {
0729     calJetsSumEt += (*jets)[i].p4().Et();
0730   }
0731 
0732   //obs63
0733   double Obs63_den = (jets->size() > 4) ? ((*jets)[0].p4().Et() + (*jets)[1].p4().Et() + (*jets)[2].p4().Et() +
0734                                            (*jets)[3].p4().Et() + (*jets)[4].p4().Et())
0735                                         : 0.0;
0736   double Obs63 = (Obs63_den > 0.00001) ? calJetsSumEt / Obs63_den : 1.0;
0737   //double Obs63 =1;
0738   evtselectVarVal.push_back(IntDblPair(63, Obs63));
0739   evtselectVarMatch.push_back(
0740       IntBoolPair(63, ((matchHadp && matchHadq) || (matchHadpq && matchHadqp)) && matchHadb && matchLepb));
0741 
0742   if (debug)
0743     std::cout << "== observable 63 " << Obs63 << std::endl;
0744 
0745   //obs64
0746   double Obs64 = solution.getProbChi2();
0747   evtselectVarVal.push_back(IntDblPair(64, Obs64));
0748   evtselectVarMatch.push_back(
0749       IntBoolPair(64, ((matchHadp && matchHadq) || (matchHadpq && matchHadqp)) && matchHadb && matchLepb));
0750 
0751   if (debug)
0752     std::cout << "== observable 64 " << Obs64 << std::endl;
0753 
0754   //obs65
0755   double Obs65 = solution.getFitHadt().p4().mass() - solution.getCalHadt().p4().mass();
0756   evtselectVarVal.push_back(IntDblPair(65, Obs65));
0757   evtselectVarMatch.push_back(
0758       IntBoolPair(65, ((matchHadp && matchHadq) || (matchHadpq && matchHadqp)) && matchHadb && matchLepb));
0759 
0760   if (debug)
0761     std::cout << "== observable 65 " << Obs65 << std::endl;
0762 
0763   //obs66
0764   double Obs66 = solution.getFitLept().p4().mass() - solution.getCalLept().p4().mass();
0765   evtselectVarVal.push_back(IntDblPair(66, Obs66));
0766   evtselectVarMatch.push_back(
0767       IntBoolPair(66, ((matchHadp && matchHadq) || (matchHadpq && matchHadqp)) && matchHadb && matchLepb));
0768 
0769   if (debug)
0770     std::cout << "observables calculated" << std::endl;
0771 
0772   if (!matchOnly)
0773     solution.setLRJetCombObservables(evtselectVarVal);
0774   return evtselectVarMatch;
0775 }