BOSS 8.0.0
BESIII Offline Software System
Loading...
Searching...
No Matches
EventPreSelectAlg.cxx
Go to the documentation of this file.
1#include <CLHEP/Vector/LorentzVector.h>
2#include <Gaudi/PluginServiceV2.h>
3#include <GaudiKernel/NTuple.h>
4#include <GaudiKernel/StatusCode.h>
5#include <string>
6#include <vector>
7
8#include "EmcRecEventModel/RecEmcShower.h"
9#include "EventModel/EventHeader.h"
10#include "EvtRecEvent/EvtRecEvent.h"
11#include "EvtRecEvent/EvtRecTrack.h"
12#include "MdcRawEvent/MdcDigi.h"
13#include "MdcRecEvent/RecMdcTrack.h"
14
15#include "EventPreSelectAlg.h"
16
17typedef std::vector<int> Vint;
18typedef std::vector<HepLorentzVector> Vp4;
19
20const double mpi = 0.13957;
21
22DECLARE_COMPONENT( EventPreSelectAlg )
23
24std::array<int, 43> EventPreSelectAlg::idmax = {
25 40, 44, 48, 56, 64, 72, 80, 80, 76, 76, 88, 88, 100, 100, 112,
26 112, 128, 128, 140, 140, 160, 160, 160, 160, 176, 176, 176, 176, 208, 208,
27 208, 208, 240, 240, 240, 240, 256, 256, 256, 256, 288, 288, 288 };
28
29StatusCode EventPreSelectAlg::initialize() {
30 info() << "In initialize()" << endmsg;
31
32 m_bhabhaEmcECut = 0.7 * m_ecm;
33 m_bhabhaMaxECut = 0.3 * m_ecm;
34 m_bhabhaSecECut = 0.1 * m_ecm;
35 m_dimuEHighCut = 0.27 * m_ecm;
36 m_dimuELowCut = 0.027 * m_ecm;
37 m_diphotonEmcECut = 0.7 * m_ecm;
38 m_diphotonSecECut = 0.3 * m_ecm;
39 m_hadronChaECut = 0.3 * m_ecm;
40 m_hadronNeuECut = 0.3 * m_ecm;
41
42 if ( m_ecm == 2.4 )
43 {
44 m_bhabhaEmcECut = 0.7 * m_ecm;
45 m_bhabhaMaxECut = 0.15 * m_ecm;
46 m_bhabhaSecECut = 0.05 * m_ecm;
47 m_diphotonEmcECut = 0.7 * m_ecm;
48 m_diphotonSecECut = 0.2 * m_ecm;
49 }
50
51 info() << "Creating tools..." << endmsg;
52
53 // retrieve IEvtWriterTool
54 std::vector<std::string> tools_to_retrieve;
55 if ( m_selectBhabha )
56 {
57 tools_to_retrieve.push_back( "SelectBarrelBhabha" );
58 tools_to_retrieve.push_back( "SelectEndcapBhabha" );
59 }
60
61 if ( m_selectDimu )
62 {
63 tools_to_retrieve.push_back( "SelectBarrelDimu" );
64 tools_to_retrieve.push_back( "SelectEndcapDimu" );
65 }
66
67 if ( m_selectDiphoton )
68 {
69 tools_to_retrieve.push_back( "SelectBarrelDiphoton" );
70 tools_to_retrieve.push_back( "SelectEndcapDiphoton" );
71 }
72
73 if ( m_selectHadron ) { tools_to_retrieve.push_back( "SelectHadron" ); }
74
75 for ( const auto& tool_name : tools_to_retrieve )
76 {
77 IEventWriterTool* tool = nullptr;
78 StatusCode sc = toolSvc()->retrieveTool( "EventWriterTool", tool_name, tool, this );
79 if ( sc.isFailure() )
80 {
81 error() << "Error retrieving tool " << tool_name << endmsg;
82 return sc;
83 }
84 else info() << "Successfully retrieved tool " << tool_name << endmsg;
85
86 m_preSelTools[tool_name] = tool;
87 }
88
89 // retrieve IDimuPreSelectTool
90 if ( m_selectDimu )
91 {
92 StatusCode sc = toolSvc()->retrieveTool( "DimuPreSelectTool", m_dimuTool, this );
93 if ( sc.isFailure() )
94 {
95 error() << "Error retrieving tool DimuPreSelectTool" << endmsg;
96 return sc;
97 }
98 else info() << "Successfully retrieved tool DimuPreSelectTool" << endmsg;
99 }
100
101 if ( m_output )
102 {
103 StatusCode status;
104 NTuplePtr nt0( ntupleSvc(), "FILE1/hadron" );
105 if ( nt0 ) m_tuple0 = nt0;
106 else
107 {
108 m_tuple0 = ntupleSvc()->book( "FILE1/hadron", CLID_ColumnWiseTuple, "N-Tuple example" );
109 if ( m_tuple0 )
110 {
111 status = m_tuple0->addItem( "esum", m_esum );
112 status = m_tuple0->addItem( "eemc", m_eemc );
113 status = m_tuple0->addItem( "etot", m_etot );
114 status = m_tuple0->addItem( "nGood", m_nGood );
115 status = m_tuple0->addItem( "nCharge", m_nCharge );
116 status = m_tuple0->addItem( "nGam", m_nGam );
117 status = m_tuple0->addItem( "ptot", m_ptot );
118 status = m_tuple0->addItem( "pp", m_pp );
119 status = m_tuple0->addItem( "pm", m_pm );
120 status = m_tuple0->addItem( "run", m_runnb );
121 status = m_tuple0->addItem( "event", m_evtnb );
122 status = m_tuple0->addItem( "maxE", m_maxE );
123 status = m_tuple0->addItem( "secE", m_secE );
124 status = m_tuple0->addItem( "dthe", m_dThe );
125 status = m_tuple0->addItem( "dphi", m_dPhi );
126 status = m_tuple0->addItem( "mdcHit1", m_mdcHit1 );
127 status = m_tuple0->addItem( "mdcHit2", m_mdcHit2 );
128 }
129 else
130 {
131 error() << " Cannot book N-tuple:" << long( m_tuple0 ) << endmsg;
132 return StatusCode::FAILURE;
133 }
134 }
135
136 NTuplePtr nt1( ntupleSvc(), "FILE1/vxyz" );
137 if ( nt1 ) m_tuple1 = nt1;
138 else
139 {
140 m_tuple1 = ntupleSvc()->book( "FILE1/vxyz", CLID_ColumnWiseTuple, "ks N-Tuple example" );
141 if ( m_tuple1 )
142 {
143 status = m_tuple1->addItem( "vx0", m_vx0 );
144 status = m_tuple1->addItem( "vy0", m_vy0 );
145 status = m_tuple1->addItem( "vz0", m_vz0 );
146 status = m_tuple1->addItem( "vr0", m_vr0 );
147 status = m_tuple1->addItem( "theta0", m_theta0 );
148 status = m_tuple1->addItem( "p0", m_p0 );
149 status = m_tuple1->addItem( "pt0", m_pt0 );
150 }
151 else
152 {
153 error() << " Cannot book N-tuple:" << long( m_tuple1 ) << endmsg;
154 return StatusCode::FAILURE;
155 }
156 }
157
158 NTuplePtr nt2( ntupleSvc(), "FILE1/photon" );
159 if ( nt2 ) m_tuple2 = nt2;
160 else
161 {
162 m_tuple2 =
163 ntupleSvc()->book( "FILE1/photon", CLID_ColumnWiseTuple, "ks N-Tuple example" );
164 if ( m_tuple2 )
165 {
166 status = m_tuple2->addItem( "dthe", m_dthe );
167 status = m_tuple2->addItem( "dphi", m_dphi );
168 status = m_tuple2->addItem( "dang", m_dang );
169 status = m_tuple2->addItem( "eraw", m_eraw );
170 }
171 else
172 {
173 error() << " Cannot book N-tuple:" << long( m_tuple2 ) << endmsg;
174 return StatusCode::FAILURE;
175 }
176 }
177
178 if ( m_selectDimu )
179 {
180 NTuplePtr nt3( ntupleSvc(), "FILE1/dimu" );
181 if ( nt3 ) m_tuple3 = nt3;
182 else
183 {
184 m_tuple3 =
185 ntupleSvc()->book( "FILE1/dimu", CLID_ColumnWiseTuple, "ks N-Tuple example" );
186 if ( m_tuple3 ) { m_dimuTool->BookNtuple( m_tuple3 ); }
187 else
188 {
189 error() << " Cannot book N-tuple:" << long( m_tuple3 ) << endmsg;
190 return StatusCode::FAILURE;
191 }
192 }
193 }
194 }
195
196 m_events = 0;
197 m_barrelBhabhaNumber = 0;
198 m_endcapBhabhaNumber = 0;
199 m_barrelDimuNumber = 0;
200 m_endcapDimuNumber = 0;
201 m_hadronNumber = 0;
202 m_barrelDiphotonNumber = 0;
203 m_endcapDiphotonNumber = 0;
204
205 info() << "successfully return from initialize()" << endmsg;
206 return StatusCode::SUCCESS;
207}
208
209StatusCode EventPreSelectAlg::execute() {
210 info() << "In execute" << endmsg;
211
212 bool isBarrelBhabha = false;
213 bool isEndcapBhabha = false;
214 bool isBarrelDimu = false;
215 bool isEndcapDimu = false;
216 bool isHadron = false;
217 bool isBarrelDiphoton = false;
218 bool isEndcapDiphoton = false;
219
220 SmartDataPtr<Event::EventHeader> eventHeader( eventSvc(), "/Event/EventHeader" );
221 if ( !eventHeader )
222 {
223 error() << "Could not retrieve event header" << endmsg;
224 return StatusCode::FAILURE;
225 }
226
227 int run = eventHeader->runNumber();
228 int event = eventHeader->eventNumber();
229 m_events++;
230
231 SmartDataPtr<EvtRecEvent> evtRecEvent( eventSvc(), EventModel::EvtRec::EvtRecEvent );
232 if ( !evtRecEvent )
233 {
234 error() << "Could not retrieve EvtRecEvent" << endmsg;
235 return StatusCode::FAILURE;
236 }
237
238 debug() << "ncharg, nneu, tottks = " << evtRecEvent->totalCharged() << " , "
239 << evtRecEvent->totalNeutral() << " , " << evtRecEvent->totalTracks() << endmsg;
240
241 SmartDataPtr<EvtRecTrackCol> evtRecTrkCol( eventSvc(), EventModel::EvtRec::EvtRecTrackCol );
242 if ( !evtRecTrkCol )
243 {
244 error() << "Could not retrieve EvtRecTrackCol" << endmsg;
245 return StatusCode::FAILURE;
246 }
247
248 if ( evtRecEvent->totalTracks() != evtRecTrkCol->size() ) return StatusCode::SUCCESS;
249
250 // -------- Good Charged Track Selection
251 std::vector<EvtRecTrack*> goodTracks;
252
253 int nCharge = 0;
254 for ( int i = 0; i < evtRecEvent->totalCharged(); i++ )
255 {
256 EvtRecTrack* recTrk = ( *evtRecTrkCol )[i];
257 if ( !recTrk->isMdcTrackValid() ) continue;
258
259 RecMdcTrack* mdcTrk = recTrk->mdcTrack();
260 double vx0 = mdcTrk->x();
261 double vy0 = mdcTrk->y();
262 double vz0 = mdcTrk->z();
263 double vr0 = mdcTrk->r();
264 double theta0 = mdcTrk->theta();
265 double p0 = mdcTrk->p();
266 double pt0 = mdcTrk->pxy();
267
268 if ( m_output )
269 {
270 m_vx0 = vx0;
271 m_vy0 = vy0;
272 m_vz0 = vz0;
273 m_vr0 = vr0;
274 m_theta0 = theta0;
275 m_p0 = p0;
276 m_pt0 = pt0;
277 m_tuple1->write().ignore();
278 }
279
280 if ( fabs( vz0 ) >= m_vz0cut ) continue;
281 if ( vr0 >= m_vr0cut ) continue;
282 if ( pt0 >= m_pt0HighCut ) continue;
283 if ( pt0 <= m_pt0LowCut ) continue;
284
285 goodTracks.push_back( recTrk );
286 nCharge += mdcTrk->charge();
287 }
288 int nGood = goodTracks.size();
289
290 // -------- Good Photon Selection
291 std::vector<EvtRecTrack*> goodGammas;
292 for ( int i = evtRecEvent->totalCharged(); i < evtRecEvent->totalTracks(); i++ )
293 {
294 EvtRecTrack* recTrk = ( *evtRecTrkCol )[i];
295 if ( !recTrk->isEmcShowerValid() ) continue;
296
297 RecEmcShower* emcTrk = recTrk->emcShower();
298 double eraw = emcTrk->energy();
299 if ( m_output )
300 {
301 m_eraw = eraw;
302 m_tuple2->write().ignore();
303 }
304
305 if ( eraw < m_energyThreshold ) continue;
306 goodGammas.push_back( recTrk );
307 }
308 int nGam = goodGammas.size();
309
310 // -------- Assign 4-momentum to each charged track
311 Vp4 ppip{ 0 }, ppim{ 0 };
312
313 double echarge = 0.; // total energy of charged track
314 double ptot = 0.; // total momentum of charged track
315 double etot = 0.; // total energy in MDC and EMC
316 double eemc = 0.; // total energy in EMC
317 double pp = 0.;
318 double pm = 0.;
319
320 for ( const auto& goodTrk : goodTracks )
321 {
322 if ( goodTrk->isMdcTrackValid() )
323 {
324 // RecMdcKalTrack* mdcKalTrk = (*itTrk)->mdcKalTrack();//After ParticleID,
325 // use RecMdcKalTrack substitute RecMdcTrack
326 // RecMdcKalTrack::setPidType(RecMdcKalTrack::pion);//PID can set to
327 // electron, muon, pion, kaon and proton;The default setting is pion If
328 // use this algorithm in reconstruction job, mdcKalTrk->p()=0!!! I don't
329 // know why.
330
331 RecMdcTrack* mdcTrk = goodTrk->mdcTrack();
332
333 ptot += mdcTrk->p();
334
335 HepLorentzVector ptrk;
336 ptrk.setPx( mdcTrk->px() );
337 ptrk.setPy( mdcTrk->py() );
338 ptrk.setPz( mdcTrk->pz() );
339 double p3 = ptrk.mag();
340 ptrk.setE( sqrt( p3 * p3 + mpi * mpi ) );
341 ptrk = ptrk.boost( -0.011, 0, 0 ); // boost to cms
342
343 echarge += ptrk.e();
344 etot += ptrk.e();
345
346 if ( mdcTrk->charge() > 0 )
347 {
348 ppip.push_back( ptrk );
349 pp = ptrk.rho();
350 }
351 else
352 {
353 ppim.push_back( ptrk );
354 pm = ptrk.rho();
355 }
356 }
357
358 if ( goodTrk->isEmcShowerValid() )
359 {
360 RecEmcShower* emcTrk = goodTrk->emcShower();
361 double eraw = emcTrk->energy();
362 double phi = emcTrk->phi();
363 double the = emcTrk->theta();
364
365 HepLorentzVector ptrk;
366 ptrk.setPx( eraw * sin( the ) * cos( phi ) );
367 ptrk.setPy( eraw * sin( the ) * sin( phi ) );
368 ptrk.setPz( eraw * cos( the ) );
369 ptrk.setE( eraw );
370 ptrk = ptrk.boost( -0.011, 0, 0 ); // boost to cms
371
372 eemc += ptrk.e();
373 etot += ptrk.e();
374 }
375 }
376
377 // -------- Assign 4-momentum to each photon
378 Vp4 pGam{ 0 };
379
380 double eneu = 0; // total energy of neutral track
381 for ( const auto& goodGamTrk : goodGammas )
382 {
383 RecEmcShower* emcTrk = goodGamTrk->emcShower();
384 double eraw = emcTrk->energy();
385 double phi = emcTrk->phi();
386 double the = emcTrk->theta();
387 HepLorentzVector ptrk;
388 ptrk.setPx( eraw * sin( the ) * cos( phi ) );
389 ptrk.setPy( eraw * sin( the ) * sin( phi ) );
390 ptrk.setPz( eraw * cos( the ) );
391 ptrk.setE( eraw );
392 ptrk = ptrk.boost( -0.011, 0, 0 ); // boost to cms
393 pGam.push_back( ptrk );
394 eneu += ptrk.e();
395 etot += ptrk.e();
396 eemc += ptrk.e();
397 }
398
399 double esum = echarge + eneu;
400
401 // -------- Use EMC shower information only
402 double maxE = 0.;
403 double secE = 0.;
404 double maxThe = 999.;
405 double maxPhi = 999.;
406 double secThe = 999.;
407 double secPhi = 999.;
408 int npart = 999.;
409 double dphi = 999.;
410 double dthe = 999.;
411 int mdcHit1 = 0.;
412 int mdcHit2 = 0.;
413
414 SmartDataPtr<RecEmcShowerCol> aShowerCol( eventSvc(), "/Event/Recon/RecEmcShowerCol" );
415 if ( !aShowerCol )
416 {
417 warning() << "Could not find RecEmcShowerCol" << endmsg;
418 return StatusCode::SUCCESS;
419 }
420
421 int ishower = 0;
422 RecEmcShowerCol::iterator iShowerCol;
423 for ( iShowerCol = aShowerCol->begin(); iShowerCol != aShowerCol->end(); iShowerCol++ )
424 {
425
426 if ( ishower == 0 )
427 {
428 maxE = ( *iShowerCol )->e5x5();
429 maxThe = ( *iShowerCol )->theta();
430 maxPhi = ( *iShowerCol )->phi();
431 npart = ( *iShowerCol )->module();
432 }
433 else if ( ishower == 1 )
434 {
435 secE = ( *iShowerCol )->e5x5();
436 secThe = ( *iShowerCol )->theta();
437 secPhi = ( *iShowerCol )->phi();
438 }
439 else if ( ishower == 2 ) { break; }
440
441 ishower++;
442 }
443
444 if ( aShowerCol->size() >= 2 )
445 {
446 dphi = ( fabs( maxPhi - secPhi ) - CLHEP::pi ) * 180. / CLHEP::pi;
447 dthe = ( fabs( maxThe + secThe ) - CLHEP::pi ) * 180. / CLHEP::pi;
448
449 double phi1 = maxPhi < 0 ? maxPhi + CLHEP::twopi : maxPhi;
450 double phi2 = secPhi < 0 ? secPhi + CLHEP::twopi : secPhi;
451
452 // Define sector (phi11,phi12) and (phi21,phi22)
453 double phi11 = min( phi1, phi2 );
454 double phi22 = max( phi1, phi2 );
455 double phi12 = ( phi11 + phi22 - CLHEP::pi ) * 0.5;
456 double phi21 = ( phi11 + phi22 + CLHEP::pi ) * 0.5;
457 if ( phi12 < 0. ) phi12 += CLHEP::twopi;
458 if ( phi21 > CLHEP::twopi ) phi21 -= CLHEP::twopi;
459
460 SmartDataPtr<MdcDigiCol> mdcDigiCol( evtSvc(), "/Event/Digi/MdcDigiCol" );
461 if ( !mdcDigiCol )
462 {
463 error() << "Could not find MdcDigiCol!" << endmsg;
464 return StatusCode::FAILURE;
465 }
466 int hitnum = mdcDigiCol->size();
467 for ( int i = 0; i < hitnum; i++ )
468 {
469 MdcDigi* digi = dynamic_cast<MdcDigi*>( mdcDigiCol->containedObject( i ) );
470 double time = digi->getTimeChannel();
471 double charge = digi->getChargeChannel();
472
473 if ( time == 0x7FFFFFFF || charge == 0x7FFFFFFF ) continue;
474
475 Identifier id( digi->identify() );
476 unsigned int iphi = MdcID::wire( id );
477 unsigned int ilayer = MdcID::layer( id );
478
479 if ( ilayer >= 43 ) error() << "MDC(" << ilayer << "," << iphi << ")" << endmsg;
480
481 double phi = CLHEP::twopi * iphi / idmax[ilayer];
482
483 if ( WhetherSector( phi, phi11, phi12 ) ) mdcHit1++;
484 else if ( WhetherSector( phi, phi21, phi22 ) ) mdcHit2++;
485 }
486 }
487
488 // If it's bhabha, return.
489 // if(nGood==2 && pp+pm>3.4) return StatusCode::SUCCESS;
490
491 // -------- Select each event
492 // Select bhabha
493 if ( eemc > m_bhabhaEmcECut && maxE > m_bhabhaMaxECut && secE > m_bhabhaSecECut &&
494 abs( dthe ) < m_bhabhaDTheCut &&
495 ( ( dphi > m_bhabhaDPhiCut1 && dphi < m_bhabhaDPhiCut2 ) ||
496 ( dphi > m_bhabhaDPhiCut3 && dphi < m_bhabhaDPhiCut4 ) ) )
497 {
498 if ( npart == 1 && mdcHit1 > m_bhabhaMdcHitCutB && mdcHit2 > m_bhabhaMdcHitCutB )
499 {
500 isBarrelBhabha = true;
501 m_barrelBhabhaNumber++;
502 }
503 else if ( ( npart == 0 || npart == 2 ) && mdcHit1 > m_bhabhaMdcHitCutE &&
504 mdcHit2 > m_bhabhaMdcHitCutE )
505 {
506 isEndcapBhabha = true;
507 m_endcapBhabhaNumber++;
508 }
509 }
510
511 // Select dimu
512 if ( m_selectDimu )
513 {
514 if ( m_dimuTool->IsDimu() == 1 )
515 {
516 isBarrelDimu = true;
517 m_barrelDimuNumber++;
518 }
519 else if ( m_dimuTool->IsDimu() == 2 )
520 {
521 isEndcapDimu = true;
522 m_endcapDimuNumber++;
523 }
524 }
525
526 // Select hadron
527 if ( ( nGood >= 1 && esum > m_hadronChaECut ) || ( nGood == 0 && esum > m_hadronNeuECut ) )
528 {
529 isHadron = true;
530 m_hadronNumber++;
531 }
532
533 // Select diphoton
534 if ( nGood == 0 && eemc > m_diphotonEmcECut && secE > m_diphotonSecECut &&
535 fabs( dthe ) < m_diphotonDTheCut && dphi > m_diphotonDPhiCut1 &&
536 dphi < m_diphotonDPhiCut2 )
537 {
538 if ( npart == 1 )
539 {
540 isBarrelDiphoton = true;
541 m_barrelDiphotonNumber++;
542 }
543 else if ( npart == 0 || npart == 2 )
544 {
545 isEndcapDiphoton = true;
546 m_endcapDiphotonNumber++;
547 }
548 }
549
550 // -------- Write to root file
551 if ( m_selectBhabha && isBarrelBhabha )
552 m_preSelTools["SelectBarrelBhabha"]->write().ignore();
553 if ( m_selectBhabha && isEndcapBhabha )
554 m_preSelTools["SelectEndcapBhabha"]->write().ignore();
555
556 if ( m_selectDimu && isBarrelDimu ) m_preSelTools["SelectBarrelDimu"]->write().ignore();
557 if ( m_selectDimu && isEndcapDimu ) m_preSelTools["SelectEndcapDimu"]->write().ignore();
558
559 if ( m_selectHadron && isHadron ) m_preSelTools["SelectHadron"]->write().ignore();
560
561 if ( m_selectDiphoton && isBarrelDiphoton )
562 m_preSelTools["SelectBarrelDiphoton"]->write().ignore();
563 if ( m_selectDiphoton && isEndcapDiphoton )
564 m_preSelTools["SelectEndcapDiphoton"]->write().ignore();
565
566 if ( m_output )
567 {
568 m_runnb = run;
569 m_evtnb = event;
570 m_esum = esum;
571 m_eemc = eemc;
572 m_etot = etot;
573 m_nCharge = nCharge;
574 m_nGood = nGood;
575 m_nGam = nGam;
576 m_ptot = ptot;
577 m_pp = pp;
578 m_pm = pm;
579 m_maxE = maxE;
580 m_secE = secE;
581 m_dThe = dthe;
582 m_dPhi = dphi;
583 m_mdcHit1 = mdcHit1;
584 m_mdcHit2 = mdcHit2;
585 m_tuple0->write().ignore();
586 }
587
588 return StatusCode::SUCCESS;
589}
590
591StatusCode EventPreSelectAlg::finalize() {
592 info() << "In finalize()" << endmsg;
593
594 if ( m_selectDimu && m_output ) { m_dimuTool->Print(); }
595
596 always() << "Total events: " << m_events << endmsg;
597
598 if ( m_selectBhabha )
599 {
600 always() << "Selected barrel bhabha: " << m_barrelBhabhaNumber << endmsg;
601 always() << "Selected endcap bhabha: " << m_endcapBhabhaNumber << endmsg;
602 }
603
604 if ( m_selectDimu )
605 {
606 always() << "Selected barrel dimu: " << m_barrelDimuNumber << endmsg;
607 always() << "Selected endcap dimu: " << m_endcapDimuNumber << endmsg;
608 }
609
610 if ( m_selectHadron ) { always() << "Selected hadron: " << m_hadronNumber << endmsg; }
611
612 if ( m_selectDiphoton )
613 {
614 always() << "Selected barrel diphoton: " << m_barrelDiphotonNumber << endmsg;
615 always() << "Selected endcap diphoton: " << m_endcapDiphotonNumber << endmsg;
616 }
617
618 return StatusCode::SUCCESS;
619}
620
621bool EventPreSelectAlg::WhetherSector( double ph, double ph1, double ph2 ) {
622 double phi1 = min( ph1, ph2 );
623 double phi2 = max( ph1, ph2 );
624 double delta = 0.0610865; // 3.5*3.1415926/180.
625 if ( ( phi2 - phi1 ) < CLHEP::pi )
626 {
627 phi1 -= delta;
628 phi2 += delta;
629 if ( phi1 < 0. ) phi1 += CLHEP::twopi;
630 if ( phi2 > CLHEP::twopi ) phi2 -= CLHEP::twopi;
631 double tmp1 = min( phi1, phi2 );
632 double tmp2 = max( phi1, phi2 );
633 phi1 = tmp1;
634 phi2 = tmp2;
635 }
636 else
637 {
638 phi1 += delta;
639 phi2 -= delta;
640 }
641 if ( ( phi2 - phi1 ) < CLHEP::pi )
642 {
643 if ( ph <= phi2 && ph >= phi1 ) return true;
644 else return false;
645 }
646 else
647 {
648 if ( ph >= phi2 || ph <= phi1 ) return true;
649 else return false;
650 }
651}
DECLARE_COMPONENT
DECLARE_COMPONENT(BesBdkRc)
phi2
Double_t phi2
Definition DataBase/tau_mode.c:8
etot
Double_t etot
Definition DataBase/tau_mode.c:8
phi1
Double_t phi1
Definition DataBase/tau_mode.c:8
time
Double_t time
Definition DataBase/tau_mode.c:7
min
#define min(a, b)
Definition Eepipi/src/ee2eepp/basesv5.1/f2c.h:157
max
#define max(a, b)
Definition Eepipi/src/ee2eepp/basesv5.1/f2c.h:158
mpi
double mpi
Definition EvtPhokharaDef.hh:17
Vp4
std::vector< HepLorentzVector > Vp4
Definition Gam4pikp.cxx:38
Vint
std::vector< int > Vint
Definition Gam4pikp.cxx:37
sin
double sin(const BesAngle a)
Definition InstallArea/x86_64-el9-gcc13-dbg/include/MdcGeom/BesAngle.h:185
cos
double cos(const BesAngle a)
Definition InstallArea/x86_64-el9-gcc13-dbg/include/MdcGeom/BesAngle.h:187
ntupleSvc
INTupleSvc * ntupleSvc()
Definition ServiceAccessor.h:51
EventPreSelectAlg::initialize
StatusCode initialize() override
Definition EventPreSelectAlg.cxx:29
EventPreSelectAlg::finalize
StatusCode finalize() override
Definition EventPreSelectAlg.cxx:591
EventPreSelectAlg::execute
StatusCode execute() override
Definition EventPreSelectAlg.cxx:209
MdcID::layer
static int layer(const Identifier &id)
Values of different levels (failure returns 0).
Definition MdcID.cxx:47
MdcID::wire
static int wire(const Identifier &id)
Definition MdcID.cxx:52
RawData::identify
virtual Identifier identify() const
Definition RawData.cxx:15
RawData::getChargeChannel
unsigned int getChargeChannel() const
Definition RawData.cxx:35
RawData::getTimeChannel
unsigned int getTimeChannel() const
Definition RawData.cxx:32