TopJetCombination/src/TtHadSimpleBestJetComb.cc

0001 #include "TopQuarkAnalysis/TopJetCombination/interface/TtHadSimpleBestJetComb.h"
0002 //
0003 // adapted Id: TtSemiSimpleBestJetComb.cc,v 1.2 2007/06/09 01:17:40 lowette Exp
0004 // for fully hadronic channel
0005
0006 /**
0007   \class    TtHadSimpleBestJetComb
0008   \brief    Based on the TtSemiSimpleBestJetComb.by: Jan Heyninck
0009    version: TtSemiSimpleBestJetComb.cc,v 1.2 2007/06/09 01:17:40 lowette Exp
0010
0011    Simple method to get the correct jet combination in hadronic ttbar events.
0012    Returns the solution with the highest probChi^2 value, starting from TtHadEvtSolutions.
0013    If more than one solution has a low enough Chi^2, the chosen best solution is that which
0014    either have both b-jets closer to a particular W-jet or the one in which the two angles added
0015    together. WARNING WARNING WARNING, last option for selection best solution by taking the sum
0016    of angles needs to be checked/approved by experts!!!!
0017 */
0018
0019 TtHadSimpleBestJetComb::TtHadSimpleBestJetComb() {}
0020
0021 TtHadSimpleBestJetComb::~TtHadSimpleBestJetComb() {}
0022
0023 int TtHadSimpleBestJetComb::operator()(std::vector<TtHadEvtSolution>& sols) {
0024   // search the highest probChi^2 value in the among the different jet combination solutions
0025   double maxProbChi2 = 0.;
0026   for (unsigned int s = 0; s < sols.size(); s++) {
0027     maxProbChi2 = std::max(maxProbChi2, sols[s].getProbChi2());
0028   }
0029
0030   //search indices of original solutions with highest probChi2 value and select those solutionsthat are close to the highest Chi^2
0031   std::vector<unsigned int> indices;
0032   indices.clear();
0033   for (unsigned int s = 0; s < sols.size(); s++) {
0034     if (fabs(sols[s].getProbChi2() - maxProbChi2) < 0.0001)
0035       indices.push_back(s);
0036   }
0037
0038   // TtHadSolutionMaker takes light jet combinations into account, but not b-jet ambiguity...
0039   int bestSol = -999;
0040   double prev_W1b = 999.;
0041   double prev_W2b = 999.;
0042   if (maxProbChi2 > 0.) {
0043     if (indices.size() == 1)
0044       bestSol = indices[0];
0045     if (indices.size() > 1) {  //for more than one solution...
0046       for (unsigned int i = 0; i != indices.size(); i++) {
0047         double DPhi_W1b0 = fabs(sols[indices[i]].getFitHadW_plus().phi() - sols[indices[i]].getFitHadb().phi());
0048         double DPhi_W1b1 = fabs(sols[indices[i]].getFitHadW_plus().phi() - sols[indices[i]].getFitHadbbar().phi());
0049         double DPhi_W2b0 = fabs(sols[indices[i]].getFitHadW_minus().phi() - sols[indices[i]].getFitHadb().phi());
0050         double DPhi_W2b1 = fabs(sols[indices[i]].getFitHadW_minus().phi() - sols[indices[i]].getFitHadbbar().phi());
0051
0052         if (DPhi_W1b0 > 3.1415)
0053           DPhi_W1b0 = 2. * 3.1415 - DPhi_W1b0;
0054         if (DPhi_W1b1 > 3.1415)
0055           DPhi_W1b1 = 2. * 3.1415 - DPhi_W1b1;
0056         if (DPhi_W2b0 > 3.1415)
0057           DPhi_W2b0 = 2. * 3.1415 - DPhi_W2b0;
0058         if (DPhi_W2b1 > 3.1415)
0059           DPhi_W2b1 = 2. * 3.1415 - DPhi_W2b1;
0060         // Select as best solution the one which either has both b-jets closer to a particular W-jet
0061         // or the one in which the two angles added together are lower than the other.....FIXME!!!
0062         // W1b0 and W1b1 is a pair, W2b0 and W2b1
0063         if (DPhi_W1b0 < DPhi_W2b0 && DPhi_W1b1 < DPhi_W2b1) {
0064           if (DPhi_W1b0 < prev_W1b && DPhi_W1b1 < prev_W2b) {
0065             bestSol = indices[i];
0066           }
0067         }
0068         if (DPhi_W1b0 > DPhi_W2b0 && DPhi_W1b1 > DPhi_W2b1) {
0069           if (DPhi_W2b0 < prev_W1b && DPhi_W2b1 < prev_W2b) {
0070             bestSol = indices[i];
0071           }
0072         }
0073         if ((DPhi_W1b0 < DPhi_W2b0 && DPhi_W1b1 > DPhi_W2b1) || (DPhi_W1b0 > DPhi_W2b0 && DPhi_W1b1 < DPhi_W2b1)) {
0074           if ((DPhi_W1b0 + DPhi_W1b1) < (DPhi_W2b0 + DPhi_W2b1)) {
0075             if (DPhi_W1b0 < prev_W1b && DPhi_W1b1 < prev_W2b) {
0076               bestSol = indices[i];
0077             }
0078           } else {
0079             if (DPhi_W2b0 < prev_W1b && DPhi_W2b1 < prev_W2b) {
0080               bestSol = indices[i];
0081             }
0082           }
0083         }
0084       }
0085     }
0086   }
0087   return bestSol;
0088 }