CMS 3D CMS Logo
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
List of all members | Classes | Public Member Functions | Private Attributes
TtHadLRSignalSelObservables Class Reference

#include <TtHadLRSignalSelObservables.h>

Classes

struct  CompareBdisc
 
struct  CompareDouble
 
struct  CompareET
 

Public Member Functions

void operator() (TtHadEvtSolution &)
 
 TtHadLRSignalSelObservables ()
 
 ~TtHadLRSignalSelObservables ()
 

Private Attributes

CompareBdisc BdiscComparator
 
CompareDouble dComparator
 
CompareET EtComparator
 
std::vector< std::pair
< unsigned int, double > > 		evtselectVarVal
 

Detailed Description

Definition at line 26 of file TtHadLRSignalSelObservables.h.

Constructor & Destructor Documentation

TtHadLRSignalSelObservables::TtHadLRSignalSelObservables ( )

Definition at line 4 of file TtHadLRSignalSelObservables.cc.

5 {
6 }
TtHadLRSignalSelObservables::~TtHadLRSignalSelObservables ( )

Definition at line 8 of file TtHadLRSignalSelObservables.cc.

9 {
10 }

Member Function Documentation

void TtHadLRSignalSelObservables::operator() ( TtHadEvtSolution &  TS)

Definition at line 12 of file TtHadLRSignalSelObservables.cc.

References BdiscComparator, funct::C, funct::cos(), dComparator, relval_parameters_module::energy, EtComparator, evtselectVarVal, TtHadEvtSolution::getHadb(), TtHadEvtSolution::getHadbbar(), TtHadEvtSolution::getHadj(), TtHadEvtSolution::getHadk(), TtHadEvtSolution::getHadp(), TtHadEvtSolution::getHadq(), class-composition::H, i, edm::isNotFinite(), cmsBatch::log, WZElectronSkims53X_cff::MT, N, PI, funct::pow(), EnergyCorrector::pt, TtHadEvtSolution::setLRSignalEvtObservables(), funct::sin(), python.multivaluedict::sort(), and mathSSE::sqrt().

13 {
14  evtselectVarVal.clear();
15  std::vector<pat::Jet> topJets;
16  topJets.clear();
17  topJets.push_back(TS.getHadp());
18  topJets.push_back(TS.getHadq());
19  topJets.push_back(TS.getHadj());
20  topJets.push_back(TS.getHadk());
21  topJets.push_back(TS.getHadb());
22  topJets.push_back(TS.getHadbbar());
23 
24  //Assume the highest Et jets are the b-jets, the others the jets from the W - but how to work out which two are from which W? FIXME!!! for now don't sort jets in Et before after this...
25  TLorentzVector *Hadp = new TLorentzVector();
26  TLorentzVector *Hadq = new TLorentzVector();
27  TLorentzVector *Hadj = new TLorentzVector();
28  TLorentzVector *Hadk = new TLorentzVector();
29  TLorentzVector *Hadb = new TLorentzVector();
30  TLorentzVector *Hadbbar = new TLorentzVector();
31  Hadp->SetPxPyPzE(topJets[0].px(),topJets[0].py(),topJets[0].pz(),topJets[3].energy());
32  Hadq->SetPxPyPzE(topJets[1].px(),topJets[1].py(),topJets[1].pz(),topJets[2].energy());
33  Hadj->SetPxPyPzE(topJets[2].px(),topJets[2].py(),topJets[2].pz(),topJets[2].energy());
34  Hadk->SetPxPyPzE(topJets[3].px(),topJets[3].py(),topJets[3].pz(),topJets[3].energy());
35  Hadb->SetPxPyPzE(topJets[4].px(),topJets[4].py(),topJets[4].pz(),topJets[1].energy());
36  Hadbbar->SetPxPyPzE(topJets[4].px(),topJets[4].py(),topJets[4].pz(),topJets[1].energy());
37 
38  //sort the topJets in Et
39  std::sort(topJets.begin(),topJets.end(),EtComparator);
40 
41  //Et-Sum of the lightest jets
42  double EtSum = topJets[5].et()+topJets[5].et();
43  double Obs1 = (EtSum>0 ? EtSum : -1);
44  evtselectVarVal.push_back(std::pair<unsigned int,double>(1,Obs1));
45 
46  //Log of the difference in Bdisc between the 2nd and the 3rd jets (ordered in Bdisc) - does that still
47  //make sense for the fully hadronic channel as well?FIXME!!!
48 
49  //sort the jets in Bdiscriminant
50  std::sort(topJets.begin(),topJets.end(),BdiscComparator);
51 
52  double BGap = topJets[1].bDiscriminator("trackCountingJetTags") - topJets[2].bDiscriminator("trackCountingJetTags");
53  double Obs2 = (BGap>0 ? log(BGap) : -1);
54  evtselectVarVal.push_back(std::pair<unsigned int,double>(2,Obs2));
55 
56  //Circularity of the event
57  double N=0,D=0,C_tmp=0,C=1000;
58  double nx,ny,nz;
59 
60  // C = 2min(E(pt.n)^2/E(pt)^2) = 2*N/D but it is theorically preferable to use C'=PI/2*min(E|pt.n|/E|pt|), sum over all jets+lepton+MET (cf PhysRevD 48 R3953(Nov 1993))
61 
62  for(unsigned int i=0;i<6;i++){
63  D += fabs(topJets[i].pt());
64  }
65 
66  if((D>0)){
67  // Loop over all the unit vectors in the transverse plane in order to find the miminum :
68  for(unsigned int i=0; i<360; i++){
69  nx = cos((2*PI/360)*i);
70  ny = sin((2*PI/360)*i);
71  nz = 0;
72  N=0;
73  for(unsigned int i=0;i<4;i++){
74  N += fabs(topJets[i].px()*nx+topJets[i].py()*ny+topJets[i].pz()*nz);
75  }
76  C_tmp = 2*N/D;
77  if(C_tmp<C) C = C_tmp;
78  }
79  }
80 
81  double Obs3 = ( C!=1000 ? C : -1);
82  evtselectVarVal.push_back(std::pair<unsigned int,double>(3,Obs3));
83 
84 
85  //HT variable (Et-sum of the six jets)
86  double HT=0;
87  for(unsigned int i=0;i<6;i++){
88  HT += topJets[i].et();
89  }
90 
91  double Obs4 = ( HT!=0 ? HT : -1);
92  evtselectVarVal.push_back(std::pair<unsigned int,double>(4,Obs4));
93 
94  //Transverse Mass of the system
95  math::XYZTLorentzVector pjets;
96  // for the six jets
97  for(unsigned int i=0;i<6;i++){
98  pjets += topJets[i].p4();
99  }
100  double MT = sqrt(std::pow(pjets.mass(),2));
101  double Obs5 = ( MT>0 ? MT : -1);
102  evtselectVarVal.push_back(std::pair<unsigned int,double>(5,Obs5));
103 
104  //CosTheta(Hadp,Hadq) and CosTheta(Hadj,Hadq)
105  //sort the lightJets in Et
106  std::sort(topJets.begin(),topJets.end(),EtComparator);
107 
108  double px1 = topJets[2].px(); double px2 = topJets[3].px();
109  double py1 = topJets[2].py(); double py2 = topJets[3].py();
110  double pz1 = topJets[2].pz(); double pz2 = topJets[3].pz();
111  double E1 = topJets[2].energy(); double E2 = topJets[3].energy();
112  double px3 = topJets[4].px(); double px4 = topJets[5].px();
113  double py3 = topJets[4].py(); double py4 = topJets[5].py();
114  double pz3 = topJets[4].pz(); double pz4 = topJets[5].pz();
115  double E3 = topJets[4].energy(); double E4 = topJets[5].energy();
116 
117  TLorentzVector *LightJet1 = new TLorentzVector();
118  TLorentzVector *LightJet2 = new TLorentzVector();
119  TLorentzVector *LightJet3 = new TLorentzVector();
120  TLorentzVector *LightJet4 = new TLorentzVector();
121 
122  LightJet1->SetPxPyPzE(px1,py1,pz1,E1);
123  //LightJet1->Boost(BoostBackToCM);
124  LightJet2->SetPxPyPzE(px2,py2,pz2,E2);
125  //LightJet2->Boost(BoostBackToCM);
126  LightJet3->SetPxPyPzE(px3,py3,pz3,E3);
127  //LightJet1->Boost(BoostBackToCM);
128  LightJet4->SetPxPyPzE(px4,py4,pz4,E4);
129  //LightJet2->Boost(BoostBackToCM);
130 
131  double CosTheta1 = cos(LightJet2->Angle(LightJet1->Vect()));
132  double CosTheta2 = cos(LightJet4->Angle(LightJet3->Vect()));
133 
134  double Obs6 = ( -1<CosTheta1 ? CosTheta1 : -2);
135  evtselectVarVal.push_back(std::pair<unsigned int,double>(6,Obs6));
136  double Obs7 = ( -1<CosTheta2 ? CosTheta2 : -2);
137  evtselectVarVal.push_back(std::pair<unsigned int,double>(7,Obs7));
138 
139  delete LightJet1;
140  delete LightJet2;
141  delete LightJet3;
142  delete LightJet4;
143 
144  // try to find out more powerful observables related to the b-disc
145  //sort the jets in Bdiscriminant
146  std::sort(topJets.begin(),topJets.end(),BdiscComparator);
147 
148  double BjetsBdiscSum = topJets[0].bDiscriminator("trackCountingJetTags") + topJets[1].bDiscriminator("trackCountingJetTags");
149  double Ljets1BdiscSum = topJets[2].bDiscriminator("trackCountingJetTags") + topJets[3].bDiscriminator("trackCountingJetTags");
150  double Ljets2BdiscSum = topJets[4].bDiscriminator("trackCountingJetTags") + topJets[5].bDiscriminator("trackCountingJetTags");
151  //std::cout<<"BjetsBdiscSum = "<<BjetsBdiscSum<<std::endl;
152  //std::cout<<"LjetsBdiscSum = "<<LjetsBdiscSum<<std::endl;
153 
154  double Obs8 = (Ljets1BdiscSum !=0 ? log(BjetsBdiscSum/Ljets1BdiscSum) : -1);
155  evtselectVarVal.push_back(std::pair<unsigned int,double>(8,Obs8));
156  double Obs9 = (Ljets2BdiscSum !=0 ? log(BjetsBdiscSum/Ljets2BdiscSum) : -1);
157  evtselectVarVal.push_back(std::pair<unsigned int,double>(9,Obs9));
158  double Obs10 = (BGap>0 ? log(BjetsBdiscSum*BGap) : -1);
159  evtselectVarVal.push_back(std::pair<unsigned int,double>(10,Obs10));
160 
161 
162  // Et-Ratio between the two highest in Et jets and four highest jets
163  double Obs11 = (HT!=0 ? (HT-EtSum)/HT : -1);
164  evtselectVarVal.push_back(std::pair<unsigned int,double>(11,Obs11));
165 
166 
167  //Sphericity and Aplanarity without boosting back the system to CM frame
168 
169  TMatrixDSym Matrix(3);
170  double PX2 = std::pow(Hadp->Px(),2)+std::pow(Hadq->Px(),2)+std::pow(Hadb->Px(),2)+
171  std::pow(Hadj->Px(),2)+std::pow(Hadk->Px(),2)+std::pow(Hadbbar->Px(),2);
172  double PY2 = std::pow(Hadp->Py(),2)+std::pow(Hadq->Py(),2)+std::pow(Hadb->Py(),2)+
173  std::pow(Hadj->Py(),2)+std::pow(Hadk->Py(),2)+std::pow(Hadbbar->Py(),2);
174  double PZ2 = std::pow(Hadp->Pz(),2)+std::pow(Hadq->Pz(),2)+std::pow(Hadb->Pz(),2)+
175  std::pow(Hadj->Pz(),2)+std::pow(Hadk->Pz(),2)+std::pow(Hadbbar->Pz(),2);
176  double P2 = PX2+PY2+PZ2;
177 
178  double PXY = Hadp->Px()*Hadp->Py()+Hadq->Px()*Hadq->Py()+Hadb->Px()*Hadb->Py()+
179  Hadj->Px()*Hadj->Py()+Hadk->Px()*Hadk->Py()+Hadbbar->Px()*Hadbbar->Py();
180  double PXZ = Hadp->Px()*Hadp->Pz()+Hadq->Px()*Hadq->Pz()+Hadb->Px()*Hadb->Pz()+
181  Hadj->Px()*Hadj->Pz()+Hadk->Px()*Hadk->Pz()+Hadbbar->Px()*Hadbbar->Pz();
182  double PYZ = Hadp->Py()*Hadp->Pz()+Hadq->Py()*Hadq->Pz()+Hadb->Py()*Hadb->Pz()+
183  Hadj->Py()*Hadj->Pz()+Hadk->Py()*Hadk->Pz()+Hadbbar->Py()*Hadbbar->Pz();
184 
185  Matrix(0,0) = PX2/P2; Matrix(0,1) = PXY/P2; Matrix(0,2) = PXZ/P2;
186  Matrix(1,0) = PXY/P2; Matrix(1,1) = PY2/P2; Matrix(1,2) = PYZ/P2;
187  Matrix(2,0) = PXZ/P2; Matrix(2,1) = PYZ/P2; Matrix(2,2) = PZ2/P2;
188 
189  TMatrixDSymEigen pTensor(Matrix);
190 
191  std::vector<double> EigValues;
192  EigValues.clear();
193  for(int i=0;i<3;i++){
194  EigValues.push_back(pTensor.GetEigenValues()[i]);
195  }
196 
197  std::sort(EigValues.begin(),EigValues.end(),dComparator);
198 
199  double Sphericity = 1.5*(EigValues[1]+EigValues[2]);
200  double Aplanarity = 1.5*EigValues[2];
201 
202  double Obs12 = (edm::isNotFinite(Sphericity) ? -1 : Sphericity);
203  evtselectVarVal.push_back(std::pair<unsigned int,double>(12,Obs12));
204 
205  double Obs13 = (edm::isNotFinite(Aplanarity) ? -1 : Aplanarity);
206  evtselectVarVal.push_back(std::pair<unsigned int,double>(13,Obs13));
207 
208 
209  //Sphericity and Aplanarity with boosting back the system to CM frame
210  TLorentzVector *TtbarSystem = new TLorentzVector();
211  TtbarSystem->SetPx(Hadp->Px()+Hadq->Px()+Hadb->Px()+Hadj->Px()+Hadk->Px()+Hadbbar->Px());
212  TtbarSystem->SetPy(Hadp->Py()+Hadq->Py()+Hadb->Py()+Hadj->Py()+Hadk->Py()+Hadbbar->Py());
213  TtbarSystem->SetPz(Hadp->Pz()+Hadq->Pz()+Hadb->Pz()+Hadj->Pz()+Hadk->Pz()+Hadbbar->Pz());
214  TtbarSystem->SetE(Hadp->E()+Hadq->E()+Hadb->E()+Hadj->E()+Hadk->E()+Hadbbar->E());
215 
216  TVector3 BoostBackToCM = -(TtbarSystem->BoostVector());
217  Hadp->Boost(BoostBackToCM);
218  Hadq->Boost(BoostBackToCM);
219  Hadb->Boost(BoostBackToCM);
220  Hadj->Boost(BoostBackToCM);
221  Hadk->Boost(BoostBackToCM);
222  Hadbbar->Boost(BoostBackToCM);
223 
224  double BOOST_PX2 = std::pow(Hadp->Px(),2)+std::pow(Hadq->Px(),2)+std::pow(Hadb->Px(),2)+
225  std::pow(Hadj->Px(),2)+std::pow(Hadk->Px(),2)+std::pow(Hadbbar->Px(),2);
226  double BOOST_PY2 = std::pow(Hadp->Py(),2)+std::pow(Hadq->Py(),2)+std::pow(Hadb->Py(),2)+
227  std::pow(Hadj->Py(),2)+std::pow(Hadk->Py(),2)+std::pow(Hadbbar->Py(),2);
228  double BOOST_PZ2 = std::pow(Hadp->Pz(),2)+std::pow(Hadq->Pz(),2)+std::pow(Hadb->Pz(),2)+
229  std::pow(Hadj->Pz(),2)+std::pow(Hadk->Pz(),2)+std::pow(Hadbbar->Pz(),2);
230 
231  double BOOST_P2 = BOOST_PX2+BOOST_PY2+BOOST_PZ2;
232 
233  double BOOST_PXY = Hadp->Px()*Hadp->Py()+Hadq->Px()*Hadq->Py()+Hadb->Px()*Hadb->Py()+
234  Hadj->Px()*Hadj->Py()+Hadk->Px()*Hadk->Py()+Hadbbar->Px()*Hadbbar->Py();
235  double BOOST_PXZ = Hadp->Px()*Hadp->Pz()+Hadq->Px()*Hadq->Pz()+Hadb->Px()*Hadb->Pz()+
236  Hadj->Px()*Hadj->Pz()+Hadk->Px()*Hadk->Pz()+Hadbbar->Px()*Hadbbar->Pz();
237  double BOOST_PYZ = Hadp->Py()*Hadp->Pz()+Hadq->Py()*Hadq->Pz()+Hadb->Py()*Hadb->Pz()+
238  Hadj->Py()*Hadj->Pz()+Hadk->Py()*Hadk->Pz()+Hadbbar->Py()*Hadbbar->Pz();
239 
240  TMatrixDSym BOOST_Matrix(3);
241 
242  BOOST_Matrix(0,0) = BOOST_PX2/BOOST_P2; BOOST_Matrix(0,1) = BOOST_PXY/BOOST_P2; BOOST_Matrix(0,2) = BOOST_PXZ/BOOST_P2;
243  BOOST_Matrix(1,0) = BOOST_PXY/BOOST_P2; BOOST_Matrix(1,1) = BOOST_PY2/BOOST_P2; BOOST_Matrix(1,2) = BOOST_PYZ/BOOST_P2;
244  BOOST_Matrix(2,0) = BOOST_PXZ/BOOST_P2; BOOST_Matrix(2,1) = BOOST_PYZ/BOOST_P2; BOOST_Matrix(2,2) = BOOST_PZ2/BOOST_P2;
245 
246  TMatrixDSymEigen BOOST_pTensor(BOOST_Matrix);
247 
248  std::vector<double> BOOST_EigValues;
249  BOOST_EigValues.clear();
250  for(int i=0;i<3;i++){
251  BOOST_EigValues.push_back(BOOST_pTensor.GetEigenValues()[i]);
252  }
253 
254  std::sort(BOOST_EigValues.begin(),BOOST_EigValues.end(),dComparator);
255 
256  double BOOST_Sphericity = 1.5*(BOOST_EigValues[1]+BOOST_EigValues[2]);
257  double BOOST_Aplanarity = 1.5*BOOST_EigValues[2];
258 
259  double Obs14 = ( edm::isNotFinite(BOOST_Sphericity) ? -1 : BOOST_Sphericity );
260  evtselectVarVal.push_back(std::pair<unsigned int,double>(14,Obs14));
261 
262  double Obs15 = ( edm::isNotFinite(BOOST_Aplanarity) ? -1 : BOOST_Aplanarity );
263  evtselectVarVal.push_back(std::pair<unsigned int,double>(15,Obs15));
264 
265 
266  // Centrality of the six jets
267  double H=0;
268  for(unsigned int i=0;i<6;i++){
269  H += topJets[i].energy();
270  }
271  double Obs16 = ( H != 0 ? HT/H : -1 );
272  evtselectVarVal.push_back(std::pair<unsigned int,double>(16,Obs16));
273 
274 
275  // distance from the origin in the (BjetsBdiscSum, Ljets1BdiscSum) and (BjetsBdiscSum, Ljets2BdiscSum)
276  double Obs17 = ( BjetsBdiscSum != 0 && Ljets1BdiscSum != 0 ? 0.707*(BjetsBdiscSum-Ljets1BdiscSum) : -1 );
277  evtselectVarVal.push_back(std::pair<unsigned int,double>(17,Obs17));
278  double Obs18 = ( BjetsBdiscSum != 0 && Ljets2BdiscSum != 0 ? 0.707*(BjetsBdiscSum-Ljets2BdiscSum) : -1 );
279  evtselectVarVal.push_back(std::pair<unsigned int,double>(18,Obs18));
280 
281  // Ratio of the Et Sum of the two lightest jets over HT=Sum of the Et of the six highest jets in Et
282  double Obs19 = ( HT != 0 ? EtSum/HT : -1 );
283  evtselectVarVal.push_back(std::pair<unsigned int,double>(19,Obs19));
284 
285 
286  // Put the vector in the TtHadEvtSolution
287  TS.setLRSignalEvtObservables(evtselectVarVal);
288 }
i
int i
Definition: DBlmapReader.cc:9
PI
#define PI
Definition: HcalBeamMonitor.cc:14
funct::sin
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
TtHadEvtSolution::getHadbbar
pat::Jet getHadbbar() const
Definition: TtHadEvtSolution.cc:85
TtHadEvtSolution::setLRSignalEvtObservables
void setLRSignalEvtObservables(const std::vector< std::pair< unsigned int, double > > &varval)
Definition: TtHadEvtSolution.cc:262
TtHadEvtSolution::getHadb
pat::Jet getHadb() const
Definition: TtHadEvtSolution.cc:37
TtHadEvtSolution::getHadq
pat::Jet getHadq() const
Definition: TtHadEvtSolution.cc:69
TtHadEvtSolution::getHadj
pat::Jet getHadj() const
Definition: TtHadEvtSolution.cc:101
hitfit::Matrix
CLHEP::HepMatrix Matrix
Definition: matutil.h:65
WZElectronSkims53X_cff.MT
string MT
Definition: WZElectronSkims53X_cff.py:219
math::XYZTLorentzVector
XYZTLorentzVectorD XYZTLorentzVector
Lorentz vector with cylindrical internal representation using pseudorapidity.
Definition: LorentzVector.h:29
edm::isNotFinite
bool isNotFinite(T x)
Definition: isFinite.h:10
mathSSE::sqrt
T sqrt(T t)
Definition: SSEVec.h:48
funct::cos
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
TtHadEvtSolution::getHadp
pat::Jet getHadp() const
Definition: TtHadEvtSolution.cc:53
TtHadEvtSolution::getHadk
pat::Jet getHadk() const
Definition: TtHadEvtSolution.cc:117
N
#define N
Definition: blowfish.cc:9
funct::D
DecomposeProduct< arg, typename Div::arg > D
Definition: Factorize.h:150
relval_parameters_module.energy
string energy
Definition: relval_parameters_module.py:29
TtHadLRSignalSelObservables::evtselectVarVal
std::vector< std::pair< unsigned int, double > > evtselectVarVal
Definition: TtHadLRSignalSelObservables.h:67
HT
Definition: HT.h:20
funct::pow
Power< A, B >::type pow(const A &a, const B &b)
Definition: Power.h:40
cmsBatch.log
tuple log
Definition: cmsBatch.py:347

Member Data Documentation

CompareBdisc TtHadLRSignalSelObservables::BdiscComparator
private

Definition at line 55 of file TtHadLRSignalSelObservables.h.

Referenced by operator()().

CompareDouble TtHadLRSignalSelObservables::dComparator
private

Definition at line 65 of file TtHadLRSignalSelObservables.h.

Referenced by operator()().

CompareET TtHadLRSignalSelObservables::EtComparator
private

Definition at line 45 of file TtHadLRSignalSelObservables.h.

Referenced by operator()().

std::vector<std::pair<unsigned int,double> > TtHadLRSignalSelObservables::evtselectVarVal
private

Definition at line 67 of file TtHadLRSignalSelObservables.h.

Referenced by operator()().

Generated for CMSSW Reference Manual by doxygen 1.8.5