BOSS 8.0.0
BESIII Offline Software System
Loading...
Searching...
No Matches
DimuPreSelectTool.cxx
Go to the documentation of this file.
1#include <GaudiKernel/NTuple.h>
2#include <GaudiKernel/SmartDataPtr.h>
3
4#include "DstEvent/TofHitStatus.h"
5#include "EmcRecEventModel/RecEmcShower.h"
6#include "MdcRecEvent/RecMdcTrack.h"
7#include "MucRawEvent/MucDigi.h"
8#include "MucRecEvent/RecMucTrack.h"
9#include "TofRecEvent/RecTofTrack.h"
10
11#include "DimuPreSelectTool.h"
12
13DECLARE_COMPONENT( DimuPreSelectTool )
14
15void DimuPreSelectTool::BookNtuple( NTuple::Tuple*& tuple ) {
16 m_output = true;
17 m_passtuple = tuple;
18 if ( !m_passtuple ) error() << "Invalid ntuple in DimuPreSelectTool!" << endmsg;
19 else
20 {
21 m_passtuple->addItem( "run", m_run ).ignore();
22 m_passtuple->addItem( "event", m_event ).ignore();
23 m_passtuple->addItem( "part", m_part ).ignore();
24 m_passtuple->addItem( "c1", m_c1 ).ignore();
25 m_passtuple->addItem( "c2", m_c2 ).ignore();
26 m_passtuple->addItem( "r1", m_r1 ).ignore();
27 m_passtuple->addItem( "r2", m_r2 ).ignore();
28 m_passtuple->addItem( "z1", m_z1 ).ignore();
29 m_passtuple->addItem( "z2", m_z2 ).ignore();
30 m_passtuple->addItem( "p1", m_p1 ).ignore();
31 m_passtuple->addItem( "p2", m_p2 ).ignore();
32 m_passtuple->addItem( "t1", m_t1 ).ignore();
33 m_passtuple->addItem( "t2", m_t2 ).ignore();
34 m_passtuple->addItem( "e1", m_e1 ).ignore();
35 m_passtuple->addItem( "e2", m_e2 ).ignore();
36 m_passtuple->addItem( "theta1", m_theta1 ).ignore();
37 m_passtuple->addItem( "theta2", m_theta2 ).ignore();
38 m_passtuple->addItem( "phi1", m_phi1 ).ignore();
39 m_passtuple->addItem( "phi2", m_phi2 ).ignore();
40 m_passtuple->addItem( "mudigi", m_mudigi ).ignore();
41 m_passtuple->addItem( "mutrk", m_mutrk ).ignore();
42
43 m_passtuple->addItem( "zeroC", m_zeroCFlag ).ignore();
44 m_passtuple->addItem( "vtRZ", m_vtRZFlag ).ignore();
45 m_passtuple->addItem( "pLim", m_pLimFlag ).ignore();
46 m_passtuple->addItem( "pSym", m_pSymFlag ).ignore();
47 m_passtuple->addItem( "tLim", m_tLimFlag ).ignore();
48 m_passtuple->addItem( "eLim", m_eLimFlag ).ignore();
49 m_passtuple->addItem( "eBB", m_eBBFlag ).ignore();
50 m_passtuple->addItem( "partSlt", m_partFlag ).ignore();
51 m_passtuple->addItem( "muDigiN", m_mudigiFlag ).ignore();
52 m_passtuple->addItem( "muTrkN", m_mutrkFlag ).ignore();
53
54 m_passtuple->addItem( "mdc", m_mdcFlag ).ignore();
55 m_passtuple->addItem( "tof", m_tofFlag ).ignore();
56 m_passtuple->addItem( "emc", m_emcFlag ).ignore();
57 m_passtuple->addItem( "muc", m_mucFlag ).ignore();
58 }
59}
60
61int DimuPreSelectTool::IsDimu() {
62 m_mdcPass = false;
63 m_tofPass = false;
64 m_emcPass = false;
65 m_mucPass = false;
66
67 /* Select by MDC info */
68 SmartDataPtr<RecMdcTrackCol> mdcTrkCol( evtSvc(), "/Event/Recon/RecMdcTrackCol" );
69 if ( !mdcTrkCol )
70 {
71 error() << "Could not find RecMdcTrackCol!" << endmsg;
72 return 0;
73 }
74 info() << "MDC tracks: " << mdcTrkCol->size() << endmsg;
75
76 if ( mdcTrkCol->size() != 2 ) return 0;
77
78 m_cutpass[0] += 1;
79
80 int c1 = ( *mdcTrkCol )[0]->charge();
81 double r1 = ( *mdcTrkCol )[0]->r();
82 double z1 = ( *mdcTrkCol )[0]->z();
83 double p1 = ( *mdcTrkCol )[0]->p();
84
85 int c2 = ( *mdcTrkCol )[1]->charge();
86 double r2 = ( *mdcTrkCol )[1]->r();
87 double z2 = ( *mdcTrkCol )[1]->z();
88 double p2 = ( *mdcTrkCol )[1]->p();
89
90 bool mdcflag1 = c1 + c2 == 0;
91 bool mdcflag2 = fabs( r1 ) <= m_vr0cut && fabs( z1 ) < m_vz0cut;
92 bool mdcflag3 = fabs( r2 ) <= m_vr0cut && fabs( z2 ) < m_vz0cut;
93 bool mdcflag4 = p1 < m_pcut_up && p2 < m_pcut_up;
94 // bool mdcflag4 = p1<2 && p1>1;
95 // bool mdcflag5 = p2<2 && p2>1;
96 bool mdcflag5 = fabs( p1 - p2 ) / ( p1 + p2 ) < m_psymcut;
97
98 info() << "r1:\t" << r1 << "\tz1:" << z1 << endmsg;
99 info() << "r2:\t" << r2 << "\tz2:" << z2 << endmsg;
100 info() << "p1:\t" << p1 << "\tp2:" << p2 << endmsg;
101
102 double zeroCFlag, vtRZFlag, pLimFlag, pSymFlag;
103 if ( mdcflag1 )
104 {
105 m_cutpass[1] += 1;
106 zeroCFlag = 1;
107 }
108 if ( mdcflag2 && mdcflag3 )
109 {
110 m_cutpass[2] += 1;
111 vtRZFlag = 1;
112 }
113 if ( mdcflag4 )
114 {
115 m_cutpass[3] += 1;
116 pLimFlag = 1;
117 }
118 if ( mdcflag5 )
119 {
120 m_cutpass[4] += 1;
121 pSymFlag = 1;
122 }
123
124 if ( mdcflag1 && mdcflag2 && mdcflag3 && mdcflag4 && mdcflag5 )
125 {
126 m_mdcPass = true;
127 m_subpass[0] += 1;
128 }
129
130 info() << "MDC selection done!" << endmsg;
131
132 /* Select by TOF info */
133 SmartDataPtr<RecTofTrackCol> tofTrackCol( evtSvc(), "/Event/Recon/RecTofTrackCol" );
134 if ( !tofTrackCol )
135 {
136 error() << "Could not find RecTofTrackCol!" << endmsg;
137 return 0;
138 }
139 info() << "TOF tracks: " << tofTrackCol->size() << endmsg;
140
141 double t1 = 0.;
142 double t2 = 1000;
143 if ( tofTrackCol->size() > 7 && tofTrackCol->size() < 21 )
144 {
145 int goodtof = 0;
146 for ( int itof = 0; itof < tofTrackCol->size(); itof++ )
147 {
148 TofHitStatus* status = new TofHitStatus;
149 status->setStatus( ( *tofTrackCol )[itof]->status() );
150
151 if ( !( status->is_cluster() ) )
152 {
153 delete status;
154 continue;
155 }
156 if ( goodtof == 0 ) t1 = ( *tofTrackCol )[itof]->tof();
157 if ( goodtof == 1 ) t2 = ( *tofTrackCol )[itof]->tof();
158
159 goodtof++;
160 delete status;
161 }
162 }
163
164 double tLimFlag;
165 bool tofflag1 = fabs( t1 - t2 ) < m_tcut; // ns
166 info() << "t1:\t" << t1 << "\tt2:" << t2 << "dt:\t" << fabs( t1 - t2 ) << endmsg;
167 if ( tofflag1 )
168 {
169 m_cutpass[5] += 1;
170 tLimFlag = 1;
171 m_tofPass = true;
172 m_subpass[1] += 1;
173 }
174 info() << "TOF selection done!" << endmsg;
175
176 /* Select by EMC info */
177 SmartDataPtr<RecEmcShowerCol> emcShowerCol( evtSvc(), "/Event/Recon/RecEmcShowerCol" );
178 if ( !emcShowerCol )
179 {
180 error() << "Could not find RecEmcShowerCol!" << endmsg;
181 return 0;
182 }
183 info() << "EMC showers:\t" << emcShowerCol->size() << endmsg;
184
185 if ( emcShowerCol->size() < 2 ) return 0;
186
187 double e1, e2, theta1, theta2, phi1, phi2;
188 int part;
189
190 e1 = ( *emcShowerCol )[0]->energy();
191 e2 = ( *emcShowerCol )[1]->energy();
192 theta1 = ( *emcShowerCol )[0]->theta();
193 theta2 = ( *emcShowerCol )[1]->theta();
194 phi1 = ( *emcShowerCol )[0]->phi();
195 phi2 = ( *emcShowerCol )[1]->phi();
196 part = ( *emcShowerCol )[0]->module();
197
198 bool emcFlag1 = e1 < m_ecut_up && e1 > m_ecut_down;
199 bool emcFlag2 = e2 < m_ecut_up && e2 > m_ecut_down;
200 bool emcFlag3 = fabs( theta1 + theta2 - CLHEP::pi ) < m_dthetacut;
201 bool emcFlag4 = fabs( phi1 - phi2 ) - CLHEP::pi < m_dphicut;
202 bool emcFlag5 = !m_partselect || ( m_partselect == 1 && part == 1 ) ||
203 ( m_partselect == 2 && part != 1 );
204
205 info() << "e1:\t" << e1 << "\te2:\t" << e2 << endmsg;
206 info() << "theta1:\t" << theta1 << "\ttheta2:\t" << theta2 << endmsg;
207 info() << "phi1:\t" << phi1 << "\tphi2:\t" << phi2 << endmsg;
208 info() << "part:\t" << part << "\tpartFlag:\t" << emcFlag5 << endmsg;
209
210 double eLimFlag, eBBFlag, partFlag;
211 if ( emcFlag1 && emcFlag2 )
212 {
213 m_cutpass[6] += 1;
214 eLimFlag = 1;
215 }
216 if ( emcFlag3 && emcFlag4 )
217 {
218 m_cutpass[7] += 1;
219 eBBFlag = 1;
220 }
221 if ( emcFlag5 )
222 {
223 m_cutpass[8] += 1;
224 partFlag = 1;
225 }
226 if ( emcFlag1 && emcFlag2 && emcFlag3 && emcFlag4 && emcFlag5 )
227 {
228 m_emcPass = true;
229 m_subpass[2] += 1;
230 }
231 info() << "EMC selection done!" << endmsg;
232
233 /* Select by MUC info */
234 // Select by MUC Info
235 SmartDataPtr<MucDigiCol> mucDigiCol( evtSvc(), "/Event/Digi/MucDigiCol" );
236 if ( !mucDigiCol )
237 {
238 error() << "Could not find MucDigiCol!" << endmsg;
239 return 0;
240 }
241 SmartDataPtr<RecMucTrackCol> mucTrackCol( evtSvc(), "/Event/Recon/RecMucTrackCol" );
242 if ( !mucTrackCol )
243 {
244 error() << "Could not find RecMucTrackCol" << endmsg;
245 return 0;
246 }
247
248 int mudigiNum, mutrkNum;
249 mudigiNum = mutrkNum = 0;
250 mudigiNum = mucDigiCol->size();
251 mutrkNum = mucTrackCol->size();
252
253 bool mucflag1 = mudigiNum >= m_mudigicut;
254 bool mucflag2 = mutrkNum >= m_mutrkcut;
255
256 info() << "MUC digi:\t" << mudigiNum << "\tMUC track:\t" << mutrkNum << endmsg;
257
258 double mudigiFlag, mutrkFlag;
259 if ( mucflag1 )
260 {
261 m_cutpass[9] += 1;
262 mudigiFlag = 1;
263 }
264 if ( mucflag2 )
265 {
266 m_cutpass[10] += 1;
267 mutrkFlag = 1;
268 }
269 if ( mucflag1 && mucflag2 )
270 {
271 m_mucPass = true;
272 m_subpass[3] += 1;
273 }
274 info() << "MUC selection done!" << endmsg;
275
276 if ( m_output )
277 {
278 m_c1 = c1;
279 m_c2 = c2;
280 m_r1 = r1;
281 m_r2 = r2;
282 m_z1 = z1;
283 m_z2 = z2;
284 m_p1 = p1;
285 m_p2 = p2;
286 m_t1 = t1;
287 m_t2 = t2;
288 m_e1 = e1;
289 m_e2 = e2;
290 m_theta1 = theta1;
291 m_theta2 = theta2;
292 m_phi1 = phi1;
293 m_phi2 = phi2;
294 m_part = part;
295 m_mudigi = mudigiNum;
296 m_mutrk = mutrkNum;
297 m_zeroCFlag = zeroCFlag;
298 m_vtRZFlag = vtRZFlag;
299 m_pLimFlag = pLimFlag;
300 m_pSymFlag = pSymFlag;
301 m_tLimFlag = tLimFlag;
302 m_eLimFlag = eLimFlag;
303 m_eBBFlag = eBBFlag;
304 m_partFlag = partFlag;
305 m_mudigiFlag = mudigiFlag;
306 m_mutrkFlag = mutrkFlag;
307
308 m_mdcFlag = m_mdcPass;
309 m_tofFlag = m_tofPass;
310 m_emcFlag = m_emcPass;
311 m_mucFlag = m_mucPass;
312
313 m_passtuple->write().ignore();
314 }
315
316 /* All selections */
317 if ( m_mdcPass && m_tofPass && m_emcPass && m_mucPass )
318 {
319 m_totalpass += 1;
320 if ( part == 1 ) return 1; // barrel
321 else return 2; // endcap
322 // setFilterPassed(true);
323 }
324
325 return 0;
326}
327
328void DimuPreSelectTool::Print() {
329 const string str[3] = { "all", "barrel", "endcap" };
330
331 cout << "pass 0 - 2 MDC tracks : " << m_cutpass[0] << endl;
332 cout << "pass 1 - 0 charge : " << m_cutpass[1] << endl;
333 cout << "pass 2 - |r|<" << m_vr0cut << " && |z|<" << m_vz0cut << " : "
334 << m_cutpass[2] << endl;
335 cout << "pass 3 - p < " << m_pcut_up << " GeV/c : " << m_cutpass[3] << endl;
336 cout << "pass 4 - p_sym < " << m_psymcut << " : " << m_cutpass[4] << endl;
337 cout << "pass 5 - |t1-t1| < " << m_tcut << " ns : " << m_cutpass[5] << endl;
338 cout << "pass 6 - " << m_ecut_down << " < e < " << m_ecut_up << " : "
339 << m_cutpass[6] << endl;
340 cout << "pass 7 - |dth|<" << m_dthetacut << " && |dphi|<" << m_dphicut << ": "
341 << m_cutpass[7] << endl;
342 cout << "pass 8 - " << str[(int)m_partselect] << " part is selected : " << m_cutpass[8]
343 << endl;
344 cout << "pass 9 - mudigi >= " << m_mudigicut << " : " << m_cutpass[9] << endl;
345 cout << "pass 10- mutrk >= " << m_mutrkcut << " : " << m_cutpass[10] << endl;
346
347 cout << "pass MDC : " << m_subpass[0] << endl;
348 cout << "pass TOF : " << m_subpass[1] << endl;
349 cout << "pass EMC : " << m_subpass[2] << endl;
350 cout << "pass MUC : " << m_subpass[3] << endl;
351
352 cout << "Total event: " << m_totevent << endl;
353 cout << "Dimu event: " << m_totalpass << endl;
354}
p2
double p2[4]
Definition BabayagaNLO.cxx:41
p1
double p1[4]
Definition BabayagaNLO.cxx:41
DECLARE_COMPONENT
DECLARE_COMPONENT(BesBdkRc)
phi2
Double_t phi2
Definition DataBase/tau_mode.c:8
phi1
Double_t phi1
Definition DataBase/tau_mode.c:8
e1
Double_t e1
Definition DataBase/tau_mode.c:7
e2
Double_t e2
Definition DataBase/tau_mode.c:7
theta2
@ theta2
Definition InstallArea/x86_64-el9-gcc13-dbg/include/TrkBase/TrkKalDeriv.h:24
theta1
@ theta1
Definition InstallArea/x86_64-el9-gcc13-dbg/include/TrkBase/TrkKalDeriv.h:24
DimuPreSelectTool::BookNtuple
void BookNtuple(NTuple::Tuple *&tuple) override
Definition DimuPreSelectTool.cxx:15
TofHitStatus::setStatus
void setStatus(unsigned int status)
Definition TofHitStatus.cxx:115