MA5SandBox: cms_sus_16_048.cpp

File cms_sus_16_048.cpp, 12.4 KB (added by Benjamin Fuks, 5 years ago)
Line 
1#include "SampleAnalyzer/User/Analyzer/cms_sus_16_048.h"
2using namespace MA5;
3using namespace std;
4
5// -----------------------------------------------------------------------------
6// Initialize
7// function called one time at the beginning of the analysis
8// -----------------------------------------------------------------------------
9bool cms_sus_16_048::Initialize(const MA5::Configuration& cfg,
10 const std::map<std::string,std::string>& parameters)
11{
12 INFO << " <><><><><><><><><><><><><><><><><><><><><><><><><>" << endmsg;
13 INFO << " <> Analysis: CMS_SUS_16_048, arXiv:1801.01846 <>" << endmsg;
14 INFO << " <> (soft dilepton LHC13, 35.9 fb^-1)<>" << endmsg;
15 INFO << " <> Recasted by: B. Fuks <>" << endmsg;
16 INFO << " <> Contacts: fuks@lpthe.jussieu.fr <>" << endmsg;
17 INFO << " <> Based on MadAnalysis 5 v1.8 <>" << endmsg;
18 INFO << " <><><><><><><><><><><><><><><><><><><><><><><><><>" << endmsg;
19
20 // Signal regions
21 Manager()->AddRegionSelection("Ewkino_lowMET_M_4to9");
22 Manager()->AddRegionSelection("Ewkino_lowMET_M_10to20");
23 Manager()->AddRegionSelection("Ewkino_lowMET_M_20to30");
24 Manager()->AddRegionSelection("Ewkino_lowMET_M_30to50");
25 Manager()->AddRegionSelection("Ewkino_medMET_M_4to9");
26 Manager()->AddRegionSelection("Ewkino_medMET_M_10to20");
27 Manager()->AddRegionSelection("Ewkino_medMET_M_20to30");
28 Manager()->AddRegionSelection("Ewkino_medMET_M_30to50");
29 Manager()->AddRegionSelection("Ewkino_highMET_M_4to9");
30 Manager()->AddRegionSelection("Ewkino_highMET_M_10to20");
31 Manager()->AddRegionSelection("Ewkino_highMET_M_20to30");
32 Manager()->AddRegionSelection("Ewkino_highMET_M_30to50");
33
34 Manager()->AddRegionSelection("stop_lowMET_PT_5to12");
35 Manager()->AddRegionSelection("stop_lowMET_PT_12to20");
36 Manager()->AddRegionSelection("stop_lowMET_PT_20to30");
37 Manager()->AddRegionSelection("stop_medMET_PT_5to12");
38 Manager()->AddRegionSelection("stop_medMET_PT_12to20");
39 Manager()->AddRegionSelection("stop_medMET_PT_20to30");
40 Manager()->AddRegionSelection("stop_highMET_PT_5to12");
41 Manager()->AddRegionSelection("stop_highMET_PT_12to20");
42 Manager()->AddRegionSelection("stop_highMET_PT_20to30");
43
44 // Preselection cuts
45 std::string ewkinos[] = {"Ewkino_lowMET_M_4to9", "Ewkino_lowMET_M_10to20",
46 "Ewkino_lowMET_M_20to30", "Ewkino_lowMET_M_30to50", "Ewkino_medMET_M_4to9",
47 "Ewkino_medMET_M_10to20", "Ewkino_medMET_M_20to30", "Ewkino_medMET_M_30to50",
48 "Ewkino_highMET_M_4to9", "Ewkino_highMET_M_10to20","Ewkino_highMET_M_20to30",
49 "Ewkino_highMET_M_30to50"};
50 std::string non_stpHM[] = {"Ewkino_lowMET_M_4to9", "Ewkino_lowMET_M_10to20",
51 "Ewkino_lowMET_M_20to30", "Ewkino_lowMET_M_30to50", "Ewkino_medMET_M_4to9",
52 "Ewkino_medMET_M_10to20", "Ewkino_medMET_M_20to30", "Ewkino_medMET_M_30to50",
53 "Ewkino_highMET_M_4to9", "Ewkino_highMET_M_10to20","Ewkino_highMET_M_20to30",
54 "Ewkino_highMET_M_30to50", "stop_lowMET_PT_5to12","stop_lowMET_PT_12to20",
55 "stop_lowMET_PT_20to30", "stop_medMET_PT_5to12", "stop_medMET_PT_12to20",
56 "stop_medMET_PT_20to30"};
57 std::string lowMET[] = {"Ewkino_lowMET_M_4to9", "Ewkino_lowMET_M_10to20",
58 "Ewkino_lowMET_M_20to30", "Ewkino_lowMET_M_30to50", "stop_lowMET_PT_5to12",
59 "stop_lowMET_PT_12to20", "stop_lowMET_PT_20to30"};
60 Manager()->AddCut("Lepton_1_init");
61 Manager()->AddCut("Lepton_1");
62 Manager()->AddCut("Lepton_2");
63 Manager()->AddCut("Lepton_2_nonStopHM", non_stpHM);
64 Manager()->AddCut("Dimuon",lowMET);
65 Manager()->AddCut("SameFlavour",ewkinos);
66 Manager()->AddCut("OSdilepton");
67 Manager()->AddCut("pt(ll)");
68 Manager()->AddCut("M(ll)_bounds", ewkinos);
69 Manager()->AddCut("M(ll)_hadronic_veto", ewkinos);
70
71 // MET cuts
72 std::string ewk_medMET[] = {"Ewkino_medMET_M_4to9", "Ewkino_medMET_M_10to20",
73 "Ewkino_medMET_M_20to30", "Ewkino_medMET_M_30to50"};
74 std::string ewk_highMET[] = {"Ewkino_highMET_M_4to9",
75 "Ewkino_highMET_M_10to20","Ewkino_highMET_M_20to30",
76 "Ewkino_highMET_M_30to50"};
77 std::string stp_medMET[] = {"stop_medMET_PT_5to12", "stop_medMET_PT_12to20",
78 "stop_medMET_PT_20to30"};
79 std::string stp_highMET[] = {"stop_highMET_PT_5to12", "stop_highMET_PT_12to20",
80 "stop_highMET_PT_20to30"};
81 Manager()->AddCut("MET");
82 Manager()->AddCut("lowMET",lowMET);
83 Manager()->AddCut("MET_corrected");
84 Manager()->AddCut("lowMET_corrected",lowMET);
85 Manager()->AddCut("ewk_medMET", ewk_medMET);
86 Manager()->AddCut("ewk_highMET",ewk_highMET);
87 Manager()->AddCut("stp_medMET", stp_medMET);
88 Manager()->AddCut("stp_highMET",stp_highMET);
89
90 // Follow up preselection
91 Manager()->AddCut("ISR_jet");
92 Manager()->AddCut("HT");
93 Manager()->AddCut("MET/HT");
94 Manager()->AddCut("b-veto");
95 Manager()->AddCut("Mtautau");
96 Manager()->AddCut("MT_MET",ewkinos);
97
98 // Mll cuts
99 std::string MLL4_9[] = {"Ewkino_lowMET_M_4to9", "Ewkino_medMET_M_4to9",
100 "Ewkino_highMET_M_4to9"};
101 std::string MLL10_20[] = {"Ewkino_lowMET_M_10to20", "Ewkino_medMET_M_10to20",
102 "Ewkino_highMET_M_10to20"};
103 std::string MLL20_30[] = {"Ewkino_lowMET_M_20to30", "Ewkino_medMET_M_20to30",
104 "Ewkino_highMET_M_20to30"};
105 std::string MLL30_50[] = {"Ewkino_lowMET_M_30to50", "Ewkino_medMET_M_30to50",
106 "Ewkino_highMET_M_30to50"};
107 Manager()->AddCut("MLL4_9",MLL4_9);
108 Manager()->AddCut("MLL10_20",MLL10_20);
109 Manager()->AddCut("MLL20_30",MLL20_30);
110 Manager()->AddCut("MLL30_50",MLL30_50);
111
112 // PT(lepton1)
113 std::string PTL5_12[] = {"stop_lowMET_PT_5to12", "stop_medMET_PT_5to12",
114 "stop_highMET_PT_5to12"};
115 std::string PTL12_20[] = {"stop_lowMET_PT_12to20", "stop_medMET_PT_12to20",
116 "stop_highMET_PT_12to20"};
117 std::string PTL20_30[] = {"stop_lowMET_PT_20to30", "stop_medMET_PT_20to30",
118 "stop_highMET_PT_20to30"};
119 Manager()->AddCut("PTL5_12", PTL5_12);
120 Manager()->AddCut("PTL12_20",PTL12_20);
121 Manager()->AddCut("PTL20_30",PTL20_30);
122
123
124 return true;
125}
126
127// -----------------------------------------------------------------------------
128// Finalize
129// function called one time at the end of the analysis
130// -----------------------------------------------------------------------------
131void cms_sus_16_048::Finalize(const SampleFormat& summary,
132 const std::vector<SampleFormat>& files)
133{
134}
135
136// -----------------------------------------------------------------------------
137// Execute
138// function called each time one event is read
139// -----------------------------------------------------------------------------
140bool cms_sus_16_048::Execute(SampleFormat& sample, const EventFormat& event)
141{
142 // Security for empty events
143 if (event.rec()==0) return true;
144
145 // Event weight
146 double myEventWeight;
147 if(Configuration().IsNoEventWeight()) myEventWeight=1.;
148 else if(event.mc()->weight()!=0.) myEventWeight = event.mc()->weight();
149 else return false;
150 Manager()->InitializeForNewEvent(myEventWeight);
151
152 // Leptons
153 std::vector<const RecLeptonFormat*> Leptons;
154 for(unsigned int ii=0; ii<event.rec()->electrons().size(); ii++)
155 {
156 const RecLeptonFormat *myElec = &(event.rec()->electrons()[ii]);
157 if(myElec->pt()<5. || myElec->abseta()>2.5) continue;
158 Leptons.push_back(myElec);
159 }
160 for(unsigned int ii=0; ii<event.rec()->muons().size(); ii++)
161 {
162 const RecLeptonFormat *myMuon = &(event.rec()->muons()[ii]);
163 if(myMuon->pt()<3.5 || myMuon->abseta()>2.4) continue;
164 Leptons.push_back(myMuon);
165 }
166 SORTER->sort(Leptons);
167 unsigned int Nl = Leptons.size();
168 if(!Manager()->ApplyCut((event.rec()->electrons().size()+event.rec()->muons().size())>0, "Lepton_1_init")) return true;
169
170 // Jets
171 std::vector<const RecJetFormat*> SignalJets;
172 for(unsigned int ii=0; ii<event.rec()->jets().size(); ii++)
173 {
174 const RecJetFormat * myJet = &(event.rec()->jets()[ii]);
175 if(myJet->pt()<25. || myJet->abseta()>2.4) continue;
176 SignalJets.push_back(myJet);
177 }
178 SignalJets = PHYSICS->Isol->JetCleaning(SignalJets, Leptons, 0.2);
179 double HT=0.;
180 unsigned int Nb=0, Nj = SignalJets.size();
181 for(unsigned int ii=0; ii<SignalJets.size(); ii++)
182 {
183 HT+=SignalJets[ii]->pt();
184 if(SignalJets[ii]->btag()) Nb++;
185 }
186
187 // Requirement on the first lepton
188 bool firstlpt = (Nl>=1) && (Leptons[0]->pt()>5.) && (Leptons[0]->pt()<30.);
189 if(!Manager()->ApplyCut(firstlpt, "Lepton_1")) return true;
190
191 // Requirement on the second lepton
192 bool seclpt_general = (Nl>=2) && (Leptons[1]->pt()<30.);
193 if(!Manager()->ApplyCut(seclpt_general, "Lepton_2")) return true;
194 bool seclpt_non_stpHM = (Leptons[1]->pt()>5.);
195 if(!Manager()->ApplyCut(seclpt_non_stpHM, "Lepton_2_nonStopHM")) return true;
196
197 // 1 OS dilepton / flavour requirements
198 bool osdilep = (Nl==2) && (Leptons[0]->charge()*Leptons[1]->charge()<0);
199 bool is_SF = (Leptons[0]->isMuon()==Leptons[1]->isMuon());
200 bool is_dimuon = is_SF && Leptons[0]->isMuon();
201 if(!Manager()->ApplyCut(is_dimuon, "Dimuon")) return true;
202 if(!Manager()->ApplyCut(is_SF, "SameFlavour")) return true;
203 if(!Manager()->ApplyCut(osdilep, "OSdilepton")) return true;
204
205 // Missing Transverse Energy
206 MALorentzVector ptmiss = event.rec()->MET().momentum();
207 double MET = ptmiss.Pt();
208 double MET_cor = (ptmiss + Leptons[0]->momentum()+Leptons[1]->momentum()).Pt();
209
210 // PT(ll) and Mll
211 double ptll = (Leptons[0]->momentum()+Leptons[1]->momentum()).Pt();
212 double mll = (Leptons[0]->momentum()+Leptons[1]->momentum()).M();
213 bool mll_lims = is_SF && (mll>4.) && (mll<50.);
214 bool mll_veto = is_SF && ((mll<9.) || (mll>10.5));
215 if(!Manager()->ApplyCut(ptll>3., "pt(ll)")) return true;
216 if(!Manager()->ApplyCut(mll_lims, "M(ll)_bounds")) return true;
217 if(!Manager()->ApplyCut(mll_veto, "M(ll)_hadronic_veto")) return true;
218
219 // MET cut
220 if(!Manager()->ApplyCut(MET>125., "MET")) return true;
221 if(!Manager()->ApplyCut(is_dimuon && MET<200., "lowMET")) return true;
222 // Trigger efficiency
223 if(is_dimuon && MET < 200.)
224 { Manager()->SetCurrentEventWeight(0.65*myEventWeight); }
225 // Corrected MET cuts
226 if(!Manager()->ApplyCut(is_dimuon && MET_cor>125., "MET_corrected")) return true;
227 if(!Manager()->ApplyCut(is_dimuon && MET_cor<200., "lowMET_corrected"))
228 return true;
229 if(!Manager()->ApplyCut(MET>200. && MET<250., "ewk_medMET")) return true;
230 if(!Manager()->ApplyCut(MET>200. && MET<300., "stp_medMET")) return true;
231 if(!Manager()->ApplyCut(MET>250., "ewk_highMET")) return true;
232 if(!Manager()->ApplyCut(MET>300., "stp_highMET")) return true;
233
234 // ISR jet
235 if(!Manager()->ApplyCut(Nj>0, "ISR_jet")) return true;
236 if(!Manager()->ApplyCut(HT>100., "HT")) return true;
237
238 // MET to HT
239 if(!Manager()->ApplyCut(MET/HT<1.4 && MET/HT>0.6, "MET/HT")) return true;
240
241 // b-veto
242 if(!Manager()->ApplyCut(Nb==0, "b-veto")) return true;
243
244 // M(tau tau)
245 double AA = -(ptmiss.Px()*Leptons[1]->py() - ptmiss.Py()*Leptons[1]->px())/
246 (Leptons[0]->py()*Leptons[1]->px() - Leptons[0]->px()*Leptons[1]->py());
247 double BB = (ptmiss.Px()*Leptons[0]->py() - ptmiss.Py()*Leptons[0]->px())/
248 (Leptons[0]->py()*Leptons[1]->px() - Leptons[0]->px()*Leptons[1]->py());
249 double mtautau = ((AA+1.)*Leptons[0]->momentum() + (BB+1.)*Leptons[1]->momentum()).M();
250 if(!Manager()->ApplyCut(mtautau<0. || mtautau>160., "Mtautau")) return true;
251
252 // MT(li, met)
253 double mt1 = sqrt(2.*Leptons[0]->pt()*MET*(1.-cos(Leptons[0]->dphi_0_pi(ptmiss))));
254 double mt2 = sqrt(2.*Leptons[1]->pt()*MET*(1.-cos(Leptons[1]->dphi_0_pi(ptmiss))));
255 if(!Manager()->ApplyCut(std::max(mt1,mt2)<70., "MT_MET")) return true;
256
257 // Mll
258 if(!Manager()->ApplyCut(mll<9., "MLL4_9")) return true;
259 if(!Manager()->ApplyCut(mll>10.5 && mll<20., "MLL10_20")) return true;
260 if(!Manager()->ApplyCut(mll>20. && mll<30., "MLL20_30")) return true;
261 if(!Manager()->ApplyCut(mll>30., "MLL30_50")) return true;
262
263 // pt(l1)
264 double pt1=Leptons[0]->pt();
265 if(!Manager()->ApplyCut(pt1<12., "PTL5_12")) return true;
266 if(!Manager()->ApplyCut(pt1>12. && pt1<20., "PTL12_20")) return true;
267 if(!Manager()->ApplyCut(pt1>20. && pt1<30., "PTL20_30")) return true;
268
269 // exit
270 return true;
271}
272
273// Overlap Removal
274template<typename T1, typename T2> std::vector<const T1*> cms_sus_16_048::Removal(std::vector<const T1*> &v1,
275 std::vector<const T2*> &v2, const double &drmin)
276{
277 // Determining with objects should be removed
278 std::vector<bool> mask(v1.size(),false);
279 for (unsigned int j=0;j<v1.size();j++)
280 for (unsigned int i=0;i<v2.size();i++)
281 if (v2[i]->dr(v1[j]) < drmin)
282 {
283 mask[j]=true;
284 break;
285 }
286
287 // Building the cleaned container
288 std::vector<const T1*> cleaned_v1;
289 for (unsigned int i=0;i<v1.size();i++)
290 if (!mask[i]) cleaned_v1.push_back(v1[i]);
291
292 return cleaned_v1;
293}
294