MA5SandBox: cms_exo_16_022.cpp

File cms_exo_16_022.cpp, 10.1 KB (added by Benjamin Fuks, 7 years ago)
Line 
1#include "SampleAnalyzer/User/Analyzer/cms_exo_16_022.h"
2using namespace MA5;
3using namespace std;
4
5double calcSumPt(const RecLeptonFormat* mylepton, double coneSize)
6{
7 double sumPt_ = 0;
8 for(unsigned int c=0; c<mylepton->isolCones().size(); c++)
9 {
10 if(!(fabs(mylepton->isolCones()[c].deltaR() - coneSize)<0.0001)) continue;
11 sumPt_ = mylepton->isolCones()[c].sumPT();
12 }
13 return sumPt_;
14}
15
16
17
18// -----------------------------------------------------------------------------
19// Initialize
20// function called one time at the beginning of the analysis
21// cms_exo_16_022.cpp version 3.
22// -----------------------------------------------------------------------------
23bool cms_exo_16_022::Initialize(const MA5::Configuration& cfg, const std::map<std::string,std::string>& parameters)
24{
25 INFO << " <><><><><><><><><><><><><><><><><><><><><><><><><><><>" << endmsg;
26 INFO << " <> Analysis: CMS-EXO-16-022, <>" << endmsg;
27 INFO << " <> (Long lived particles) <>" << endmsg;
28 INFO << " <> Recasted by: J.Chang, D.Kang, P.Wu and S.Yang <>" << endmsg;
29 INFO << " <> Contact: lovejesus99wwjd@gmail.com <>" << endmsg;
30 INFO << " <> dayou17@gmail.com <>" << endmsg;
31 INFO << " <> peiwen.wu123@gmail.com <>" << endmsg;
32 INFO << " <> slowmoyang@gmail.com <>" << endmsg;
33 INFO << " <> Based on MadAnalysis 5 v1.6.25 <>" << endmsg;
34 INFO << " <> For more information, see <>" << endmsg;
35 INFO << " <> http://madanalysis.irmp.ucl.ac.be/wiki/PhysicsAnalysisDatabase" << endmsg;
36 INFO << " <><><><><><><><><><><><><><><><><><><><><><><><>" << endmsg;
37
38
39 // Declaration of the signal regions
40 Manager()->AddRegionSelection("SR1");
41 Manager()->AddRegionSelection("SR2");
42 Manager()->AddRegionSelection("SR3");
43
44 // Declaration of the preselection cuts
45 Manager()->AddCut("one e, one m");
46 Manager()->AddCut("e v0, vz");
47 Manager()->AddCut("m v0, vz");
48
49 Manager()->AddCut("opposite charge");
50 Manager()->AddCut("dR_emu>0.5");
51
52 string SR1Cut[] = {
53 "SR1"};
54 Manager()->AddCut("200<|d0|<500",SR1Cut);
55
56 string SR2Cut[] = {
57 "SR2"};
58 Manager()->AddCut("500<|d0|<1000",SR2Cut);
59
60 string SR3Cut[] = {
61 "SR3"};
62 Manager()->AddCut("1000<|d0|<100000",SR3Cut);
63
64 // Histogram declaration
65 Manager()->AddHisto("electron_d0", 40,0,1000);
66 Manager()->AddHisto("muon_d0", 40,0,1000);
67 Manager()->AddHisto("dR(e,mu)", 50,0,5);
68 Manager()->AddHisto("eta(e)", 50,-5,5);
69 Manager()->AddHisto("eta(mu)", 50,-5,5);
70 Manager()->AddHisto("PT(e)", 50,0,1000);
71 Manager()->AddHisto("PT(mu)", 50,0,1000);
72
73
74 // Histogram declaration
75 Manager()->AddHisto("electron_d0_c", 40,0,1000);
76 Manager()->AddHisto("muon_d0_c", 40,0,1000);
77 Manager()->AddHisto("dR(e,mu)_c", 50,0,5);
78 Manager()->AddHisto("eta(e)_c", 50,-5,5);
79 Manager()->AddHisto("eta(mu)_c", 50,-5,5);
80 Manager()->AddHisto("PT(e)_c", 50,0,1000);
81 Manager()->AddHisto("PT(mu)_c", 50,0,1000);
82
83 return true;
84}
85
86// -----------------------------------------------------------------------------
87// Finalize
88// function called one time at the end of the analysis
89// -----------------------------------------------------------------------------
90void cms_exo_16_022::Finalize(const SampleFormat& summary, const std::vector<SampleFormat>& files)
91{
92 cout << "BEGIN Finalization" << endl;
93 cout << "END Finalization" << endl;
94}
95
96// -----------------------------------------------------------------------------
97// Execute
98// function called each time one event is read
99// -----------------------------------------------------------------------------
100bool cms_exo_16_022::Execute(SampleFormat& sample, const EventFormat& event)
101{
102
103 double myWeight;
104 if(Configuration().IsNoEventWeight()) myWeight=1.;
105 else if(event.mc()->weight()!=0.) myWeight=event.mc()->weight();
106 else
107 {
108 WARNING << "Found one event with a zero weight. Skipping..." << endmsg;
109 return false;
110 }
111 Manager()->InitializeForNewEvent(myWeight);
112
113 if (event.rec()!=0)
114 {
115 std::vector<const RecLeptonFormat*> SignalElectrons,SignalElectronsSR1, SignalElectronsSR2, SignalElectronsSR3;
116 std::vector<const RecLeptonFormat*> SignalMuons, SignalMuonsSR1, SignalMuonsSR2, SignalMuonsSR3;
117
118
119
120 // ------- Fill the electron container ------------------------------------------------
121
122
123 for (MAuint32 i=0;i<event.rec()->electrons().size();i++)
124 {
125 const RecLeptonFormat *myElec = &(event.rec()->electrons()[i]);
126 double eta = fabs(myElec->eta());
127 double pt = myElec->pt();
128 double sumpt = calcSumPt(myElec, 0.3);
129
130 // ------- isolation condition ------- //
131 if(eta > 1.57 && eta < 2.4){if(!(sumpt/pt < 0.065)) continue;}
132 if(eta < 1.44) {if(!(sumpt/pt < 0.035)) continue; }
133
134 // ------- preselection ------- //
135 if( (eta > 1.44) and (eta<1.57) ) continue;
136 if(eta > 2.4) continue;
137 if(pt < 42) continue;
138
139 SignalElectrons.push_back(myElec);
140 }
141
142
143 // --------------------------------------------------------------------------------
144
145
146
147
148 // ------- Fill the muon container ------------------------------------------------
149
150 for (MAuint32 i=0;i<event.rec()->muons().size();i++)
151 {
152
153 const RecLeptonFormat *myMuon = &(event.rec()->muons()[i]);
154 double eta = fabs(myMuon->eta());
155 double pt = myMuon->pt();
156 double sumpt = calcSumPt(myMuon, 0.4);
157
158 // ------- isolation ------- //
159
160 if(!(sumpt/pt < 0.15)) continue;
161
162 // ------- preselection ------- //
163
164 if(eta > 2.4) continue;
165 if(pt < 40) continue;
166
167 SignalMuons.push_back(myMuon);
168 }
169
170
171 // --------------------------------------------------------------------------------
172
173 int ne = SignalElectrons.size(); // ---- # of e pass preselection
174 int nm = SignalMuons.size(); // ---- # of mu pass preselection
175 int ce = 0; // ---- charge of e
176 int cm = 0; // ---- charge of mu
177 double d0e = 0; // ---- d0 of e
178 double d0m = 0; // ---- d0 of mu
179 double pte = 0; // ---- pt of e
180 double ptm = 0; // ---- pt of mu
181 double etae = 0; // ---- eta of e
182 double etam = 0; // ---- eta of mu
183 double dr_em=0; // ---- DeltaR (e,mu)
184
185 double exd = 0;
186 double eyd = 0;
187 double ezd = 0;
188 double ev0 = 0;
189 double mxd = 0;
190 double myd = 0;
191 double mzd = 0;
192 double mv0 = 0;
193
194 // ------- only 1 e and 1 muon is allowed, no for interatin is needed ---
195
196 if(!Manager()->ApplyCut((ne == 1 && nm == 1),"one e, one m"))
197 return true;
198
199 int i=0;
200 int j=0;
201 d0e = fabs(SignalElectrons[j]->d0());
202 pte = SignalElectrons[j]->pt();
203 etae = SignalElectrons[j]->eta();
204 ce = SignalElectrons[j]->charge();
205
206 d0m = fabs(SignalMuons[i]->d0());
207 ptm = SignalMuons[i]->pt();
208 etam = SignalMuons[i]->eta();
209 cm = SignalMuons[i]->charge();
210
211 dr_em = SignalMuons[i]->dr(SignalElectrons[j]);
212
213 exd = fabs(SignalElectrons[j]->closestPoint().X());
214 eyd = fabs(SignalElectrons[j]->closestPoint().Y());
215 ezd = fabs(SignalElectrons[j]->closestPoint().Z());
216 ev0 = sqrt( exd * exd + eyd * eyd);
217 mxd = fabs(SignalMuons[i]->closestPoint().X());
218 myd = fabs(SignalMuons[i]->closestPoint().Y());
219 mzd = fabs(SignalMuons[i]->closestPoint().Z());
220 mv0 = sqrt( mxd * mxd + myd * myd);
221
222 // ------- Fill Histograms before dr, charge, and number of e, mu selection. ------
223
224 Manager()->FillHisto("electron_d0", d0e);
225 Manager()->FillHisto("muon_d0", d0m);
226 Manager()->FillHisto("dR(e,mu)", dr_em);
227 Manager()->FillHisto("eta(e)", etae);
228 Manager()->FillHisto("eta(mu)", etam);
229 Manager()->FillHisto("PT(e)", pte);
230 Manager()->FillHisto("PT(mu)", ptm);
231
232 // ------- Apply Cuts -------------------------------------------------------------
233
234 if(!Manager()->ApplyCut((ezd < 300 && ev0 < 40),"e v0, vz"))
235 return true;
236 if(!Manager()->ApplyCut((mzd < 300 && mv0 < 40),"m v0, vz"))
237 return true;
238
239 if(!Manager()->ApplyCut((ce!=cm),"opposite charge"))
240 return true;
241
242 if(!Manager()->ApplyCut((dr_em > 0.5),"dR_emu>0.5"))
243 return true;
244 Manager()->FillHisto("dR(e,mu)", dr_em);
245
246 if(!Manager()->ApplyCut((d0e>0.2 && d0m > 0.2) && (d0e < 100. && d0m < 100.) ,"200<|d0|<500"))
247 return true;
248
249
250 if(!Manager()->ApplyCut((d0e>0.5 && d0m > 0.5) && (d0e < 100. && d0m < 100.) ,"500<|d0|<1000"))
251 return true;
252
253
254 if(!Manager()->ApplyCut((d0e>1. && d0m > 1.) && (d0e < 100. && d0m < 100.) ,"1000<|d0|<100000"))
255 return true;
256
257
258
259 // --------------------------------------------------------------------------------
260
261
262 // ------- Fill Histograms before dr, charge, and number of e, mu selection. ------
263
264 Manager()->FillHisto("electron_d0_c", d0e);
265 Manager()->FillHisto("muon_d0_c", d0m);
266 Manager()->FillHisto("dR(e,mu)_c", dr_em);
267 Manager()->FillHisto("eta(e)_c", etae);
268 Manager()->FillHisto("eta(mu)_c", etam);
269 Manager()->FillHisto("PT(e)_c", pte);
270 Manager()->FillHisto("PT(mu)_c", ptm);
271
272
273 }
274
275 return true;
276}
277