Skip to content

Commit 23aa7a2

Browse files
author
Markus
committed
first commit
0 parents  commit 23aa7a2

3 files changed

Lines changed: 820 additions & 0 deletions

File tree

‎SUSYLooperHists.C‎

Lines changed: 229 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,229 @@
1+
#define SUSYLooperHists_cxx
2+
#include "SUSYLooperHists.h"
3+
#include <TH2.h>
4+
#include <TStyle.h>
5+
#include <TCanvas.h>
6+
#include <vector>
7+
8+
#include <TF1.h>
9+
#include <TMath.h>
10+
#include <HTransformToCS.c>
11+
#include <TMatrixF.h>
12+
#include <TVectorF.h>
13+
14+
TFile* SUSYLooperHists::Loop()
15+
{
16+
// In a ROOT session, you can do:
17+
// Root > .L SUSYLooperHists.C
18+
// Root > SUSYLooperHists t
19+
// Root > t.GetEntry(12); // Fill t data members with entry number 12
20+
// Root > t.Show(); // Show values of entry 12
21+
// Root > t.Show(16); // Read and show values of entry 16
22+
// Root > t.Loop(); // Loop on all entries
23+
//
24+
25+
// This is the loop skeleton where:
26+
// jentry is the global entry number in the chain
27+
// ientry is the entry number in the current Tree
28+
// Note that the argument to GetEntry must be:
29+
// jentry for TChain::GetEntry
30+
// ientry for TTree::GetEntry and TBranch::GetEntry
31+
//
32+
// To read only selected branches, Insert statements like:
33+
// METHOD1:
34+
// fChain->SetBranchStatus("*",0); // disable all branches
35+
// fChain->SetBranchStatus("branchname",1); // activate branchname
36+
// METHOD2: replace line
37+
// fChain->GetEntry(jentry); //read all branches
38+
//by b_branchname->GetEntry(ientry); //read only this branch
39+
40+
41+
if (fChain == 0) return NULL;
42+
Long64_t nentries = fChain->GetEntriesFast();
43+
HTransformToCS mtCS;
44+
mtCS.SetCmsEnergy(8000);
45+
46+
47+
TFile* myFile = new TFile(outputFileName,"recreate");
48+
TH1D* LepTwoZmass= new TH1D("LepTwoZmass",";approximate Z(#tau#tau) mass [GeV];",100,0,2000);
49+
TH2D* scan = new TH2D("scan","scan",8,112.5,312.5,40,112.5,312.5);
50+
TH2D* scanA = new TH2D("scanA","scan",8,112.5,312.5,40,112.5,312.5);
51+
TH2D* scanB = new TH2D("scanB","scan",8,112.5,312.5,40,112.5,312.5);
52+
TH2D* scanC = new TH2D("scanC","scan",8,112.5,312.5,40,112.5,312.5);
53+
// baseline plots +++++++++++++++++++++
54+
TH1D* CutFlowAntonis = new TH1D("CutFlowAntonis",";P(l);",20,-0.5,19.5);
55+
CutFlowAntonis->Sumw2();
56+
57+
Long64_t nbytes = 0, nb = 0;
58+
for (Long64_t jentry=0; jentry<nentries;jentry++) {
59+
60+
Long64_t ientry = LoadTree(jentry);
61+
if (ientry < 0) break;
62+
nb = fChain->GetEntry(jentry); nbytes += nb;
63+
64+
if(jentry%10000==0) cout << "Event: "<<jentry <<" weight: " << weight <<endl;
65+
scan->Fill(GenSusyMStop,GenSusyMNeutralino);
66+
// select lepton and put it into TLorentzvector
67+
68+
TLorentzVector lep;
69+
TLorentzVector met;
70+
vector <TLorentzVector> lepVec;
71+
for(int i =0;i<nLepGood;i++)
72+
{
73+
TLorentzVector aLep;
74+
aLep.SetPtEtaPhiM(LepGood_pt[i],LepGood_eta[i],LepGood_phi[i],0.1);
75+
lepVec.push_back(aLep);
76+
}
77+
if(nLepGood>0) lep.SetPtEtaPhiM(LepGood_pt[0],LepGood_eta[0],LepGood_phi[0],0.1);
78+
else lep.SetPtEtaPhiM(0,0,0,0.1);
79+
met.SetPtEtaPhiM(met_pt,0.,met_phi,0);
80+
float MT = sqrt(lep.Pt()*2*met.Pt()*(1-cos(lep.DeltaPhi(met))));
81+
TLorentzVector jetSum;
82+
vector <TLorentzVector> myVecJets;
83+
jetSum.SetPtEtaPhiM(0,0,0,0);
84+
unsigned int nb40=0;
85+
unsigned int nb30L=0;
86+
87+
float HT30 = 0;
88+
for (int i=0;i<nJet;i++)
89+
{
90+
TLorentzVector myjet;
91+
//cout <<" pt: " << Jet_pt[i] << " "<< Jet_eta[i]<< " "<< Jet_phi[i]<< " QG "<< Jet_quarkGluonID[i] <<endl;
92+
myjet.SetPtEtaPhiM(Jet_pt[i],Jet_eta[i],Jet_phi[i],Jet_mass[i]);
93+
jetSum=jetSum+myjet;
94+
myVecJets.push_back(myjet);
95+
if(Jet_btagCSV[i]>0.679) {
96+
97+
nb40++;
98+
}
99+
100+
if(Jet_btagCSV[i]>0.244) {
101+
102+
nb30L++;
103+
}
104+
105+
if(fabs(Jet_eta[i])<4.5&&Jet_pt[i]>30) HT30=HT30+Jet_pt[i];
106+
107+
// QG->Fill(Jet_quarkGluonID[i]);
108+
}
109+
110+
// preselection syncronizd to Antonis ++++++++++++++++++++
111+
if( nLepGood == 0 ) continue;
112+
CutFlowAntonis->Fill(0.);
113+
if(met.Pt()<200) continue;
114+
if(nLepGood<2) continue;
115+
if(nJet==0) continue;
116+
if(nTauGood!=0) continue;
117+
CutFlowAntonis->Fill(1.);
118+
if(fabs(LepGood_pdgId[0]* LepGood_pdgId[1])<122) continue;
119+
CutFlowAntonis->Fill(2.);
120+
if(lep.Pt()<5||fabs(lep.Eta())>1.5||lep.Pt()>60) continue;
121+
CutFlowAntonis->Fill(3.);
122+
if(fabs(LepGood_pdgId[0])==11&&lep.Pt()<7) continue;
123+
CutFlowAntonis->Fill(4.);
124+
// check that a cut was ahead in nLepton !!!!!!!!
125+
TLorentzVector secondLep;
126+
secondLep.SetPtEtaPhiM(LepGood_pt[1], LepGood_eta[1], LepGood_phi[1],0.1);
127+
if(secondLep.Pt()<3||fabs(secondLep.Eta())>1.5||secondLep.Pt()>60) continue;
128+
CutFlowAntonis->Fill(5.);
129+
// lepton cuts
130+
if(LepGood_relIso[0]*lep.Pt()>10.) continue;
131+
if(LepGood_relIso[1]*secondLep.Pt()>5.) continue;
132+
if(LepGood_relIso[1]>.5) continue;
133+
CutFlowAntonis->Fill(6.);
134+
// JET CUTS
135+
if(nJet>0){
136+
if (Jet_pt[0] < 150) continue;
137+
if(fabs(Jet_eta[0]) > 2.4) continue;
138+
}
139+
if(nJet>2){
140+
if (Jet_pt[2] > 60) continue;
141+
}
142+
CutFlowAntonis->Fill(7.);
143+
if(met.Pt()/HT30<2./3.) continue;
144+
if(nb40!=0) continue;
145+
CutFlowAntonis->Fill(9.);
146+
// preselection syncronizd to Antonis ++++++++++++++++++++
147+
148+
if(MT>60&&MT<100) continue;
149+
150+
151+
float pairmass = DiTau_InvMass(met,lep,secondLep,0);
152+
if(lep.Pt()>25||secondLep.Pt()>15) continue;
153+
if (pairmass>2000) pairmass=1999.;
154+
if (pairmass<0) pairmass=0.;
155+
156+
LepTwoZmass->Fill(pairmass,weight);
157+
158+
159+
if((LepGood_pdgId[0]* LepGood_pdgId[1]<-121))
160+
{
161+
162+
if(pairmass>160||pairmass<20){
163+
CutFlowAntonis->Fill(10., weight*puWeight);
164+
scanB->Fill(GenSusyMStop,GenSusyMNeutralino);
165+
166+
if(nb30L==0)
167+
{
168+
scanC->Fill(GenSusyMStop,GenSusyMNeutralino);
169+
}
170+
for (int i=0;i<nJet;i++)
171+
{
172+
if(Jet_btagCSV[i]>0.244 ) {
173+
// LooseBPt->Fill(Jet_pt[i],weight );
174+
}
175+
// QG->Fill(Jet_quarkGluonID[i]);
176+
}
177+
178+
}
179+
180+
181+
else // same sign leptons
182+
{
183+
184+
}
185+
}
186+
187+
}
188+
189+
190+
scanA->Divide(scan);
191+
scanB->Divide(scan);
192+
scanC->Divide(scan);
193+
194+
CutFlowAntonis->DrawCopy();
195+
196+
myFile->Write();
197+
return (myFile);
198+
// if (Cut(ientry) < 0) continue;
199+
}
200+
201+
double SUSYLooperHists::DiTau_InvMass( TLorentzVector Met, TLorentzVector L1, TLorentzVector L2, float Al ){
202+
TLorentzVector T1,T2;
203+
double DiTauMass;
204+
TMatrixF A(2,2);
205+
TVectorF C(2),X(2);
206+
A(0,0)=L1.Px();
207+
A(0,1)=L2.Px();
208+
A(1,0)=L1.Py();
209+
A(1,1)=L2.Py();
210+
A=A.Invert();
211+
C(0)=(Met+L1+L2).Px();
212+
C(1)=(Met+L1+L2).Py();
213+
X=A*C;// double X0i=X(0), X1i=X(1);
214+
//---------------[ MET ReAlignement subsection ]------------------------------
215+
if(X(0)<0||X(1)<0){
216+
if ( fabs(L1.DeltaPhi(Met))>Al && fabs(L2.DeltaPhi(Met))>Al ) {/*DO NOTHING just normaly a non-Z event!*/}
217+
else if( fabs(L1.DeltaPhi(Met))<Al && fabs(L2.DeltaPhi(Met))>Al ) Met.SetPtEtaPhiM(Met.Pt(),0,L1.Phi(),0);
218+
else if( fabs(L2.DeltaPhi(Met))<Al && fabs(L1.DeltaPhi(Met))>Al ) Met.SetPtEtaPhiM(Met.Pt(),0,L2.Phi(),0);
219+
else if( fabs(L1.DeltaPhi(Met))<Al && fabs(L2.DeltaPhi(Met))<Al && fabs(L1.DeltaPhi(Met))<fabs(L2.DeltaPhi(Met)) ) Met.SetPtEtaPhiM(Met.Pt(),0,L1.Phi(),0);
220+
else if( fabs(L1.DeltaPhi(Met))<Al && fabs(L2.DeltaPhi(Met))<Al && fabs(L1.DeltaPhi(Met))>fabs(L2.DeltaPhi(Met)) ) Met.SetPtEtaPhiM(Met.Pt(),0,L2.Phi(),0);
221+
}//---------------------------------------------------------------------------
222+
C(0)=(Met+L1+L2).Px(); C(1)=(Met+L1+L2).Py(); X=A*C;
223+
T1.SetPxPyPzE( L1.Px()*X(0), L1.Py()*X(0), L1.Pz()*X(0), sqrt( 3.1571 +L1.P()*L1.P()*X(0)*X(0) ) );
224+
T2.SetPxPyPzE( L2.Px()*X(1), L2.Py()*X(1), L2.Pz()*X(1), sqrt( 3.1571 +L2.P()*L2.P()*X(1)*X(1) ) );
225+
if( X(0)>0 && X(1)>0 ) DiTauMass=(T1+T2).M(); else DiTauMass=-(T1+T2).M(); return DiTauMass;
226+
//if((X(0)!=X0i||X(1)!=X1i))std::cout<<X(0)<<" "<<X(1)<<" <--"<<X0i<<" "<<X1i<<" RMETal.phi="<<(T1+T2-L1-L2).Phi()<<" RMETal.eta"<<(T1+T2-L1-L2).Eta()<<" MZ="<<DiTauMass<<endl;
227+
}
228+
229+

0 commit comments

Comments
 (0)