BOSS 7.0.6
BESIII Offline Software System
Loading...
Searching...
No Matches
LambdaCReconstruction.cxx
Go to the documentation of this file.
1//
2// All D+ decay modes Reconstruction
3//
4#include <fstream>
5#include "GaudiKernel/MsgStream.h"
6#include "GaudiKernel/AlgFactory.h"
7#include "GaudiKernel/ISvcLocator.h"
8#include "GaudiKernel/SmartDataPtr.h"
9#include "GaudiKernel/IDataProviderSvc.h"
10
12#include "EventModel/Event.h"
14
21
26#include "BesDChain/CDPi0List.h"
27#include "BesDChain/CDEtaList.h"
28#include "BesDChain/CDKsList.h"
31
32#include "McTruth/McParticle.h"
34#include "VertexFit/VertexFit.h"
37
53
54#include "DTagAlg/utility.h"
56
57
58
59using namespace Event;
60
61//*******************************************************************************************
62LambdaCReconstruction::LambdaCReconstruction(const std::string& name, ISvcLocator* pSvcLocator):
63 Algorithm(name, pSvcLocator) {
64 //Declare the properties
65 declareProperty( "debug", m_debug = false );
66 declareProperty( "ReadBeamEFromDB", m_ReadBeamEFromDB = false );
67 declareProperty( "UseCalibBeamE", m_usecalibBeamE = false );
68 declareProperty( "UseVertexfit", m_usevertexfit = false );
69 declareProperty( "BeamE", m_beamE = 2.3);
70 declareProperty( "LcList", m_decaylist = "test.txt" );
71 }
72
73//******************************************************************************************
75 MsgStream log(msgSvc(), name());
76 log << MSG::INFO << "in initialize()" <<endreq;
77
78 m_irun=-100;
79 m_beta.setX(0.011);
80 m_beta.setY(0);
81 m_beta.setZ(0);
82 chanlist=getlist(m_decaylist);
83
84 return StatusCode::SUCCESS;
85}
86
87//********************************************************************************************
89 MsgStream log(msgSvc(), name());
90 log << MSG::INFO << "in finalize()" << endreq;
91
92 chanlist.clear();
93
94 return StatusCode::SUCCESS;
95}
96
97StatusCode LambdaCReconstruction::registerEvtRecDTagCol(
98 EvtRecDTagCol* aNewEvtRecDTagCol, MsgStream& log) {
99 StatusCode sc = eventSvc()->registerObject("/Event/EvtRec/EvtRecDTagCol",
100 aNewEvtRecDTagCol);
101 if (sc != StatusCode::SUCCESS) {
102 log << MSG::FATAL << "Could not register EvtRecDTagCol in TDS!" << endreq;
103 }
104 return sc;
105}
106
107//*********************************************************************************************
109 MsgStream log(msgSvc(), name());
110 log << MSG::INFO << "in execute()" << endreq;
111
112 StatusCode sc;
113
114 //////////////////
115 // Read REC data
116 /////////////////
117 SmartDataPtr<Event::EventHeader> eventHeader(eventSvc(),"/Event/EventHeader");
118 int event= eventHeader->eventNumber();
119 // if ( m_debug || ( (event & 0x3FF) == 0 ) )
120 //std::cout << "event: " << event << std::endl;
121
122 SmartDataPtr<EvtRecEvent> recEvent(eventSvc(), EventModel::EvtRec::EvtRecEvent);
123 SmartDataPtr<EvtRecTrackCol> recTrackCol(eventSvc(), EventModel::EvtRec::EvtRecTrackCol);
124 log << MSG::DEBUG << "run and event = " << eventHeader->runNumber()
125 << " " << eventHeader->eventNumber() << endreq;
126 log << MSG::DEBUG <<"ncharg, nneu, tottks = "
127 << recEvent->totalCharged() << " , "
128 << recEvent->totalNeutral() << " , "
129 << recEvent->totalTracks() <<endreq;
130
131 EvtRecTrackIterator charged_begin = recTrackCol->begin();
132 EvtRecTrackIterator charged_end = charged_begin + recEvent->totalCharged();
133
134 EvtRecTrackIterator neutral_begin = recTrackCol->begin()+recEvent->totalCharged();
135 EvtRecTrackIterator neutral_end = recTrackCol->begin()+recEvent->totalTracks();
136
137
138 SmartDataPtr<EvtRecPi0Col> recPi0Col(eventSvc(), "/Event/EvtRec/EvtRecPi0Col");
139 if ( ! recPi0Col ) {
140 log << MSG::FATAL << "Could not find EvtRecPi0Col" << endreq;
141 return StatusCode::FAILURE;
142 }
143
144
145 SmartDataPtr<EvtRecEtaToGGCol> recEtaToGGCol(eventSvc(), "/Event/EvtRec/EvtRecEtaToGGCol");
146 if ( ! recEtaToGGCol ) {
147 log << MSG::FATAL << "Could not find EvtRecEtaToGGCol" << endreq;
148 return StatusCode::FAILURE;
149 }
150
151
152 SmartDataPtr<EvtRecVeeVertexCol> recVeeVertexCol(eventSvc(), "/Event/EvtRec/EvtRecVeeVertexCol");
153 if ( ! recVeeVertexCol ) {
154 log << MSG::FATAL << "Could not find EvtRecVeeVertexCol" << endreq;
155 return StatusCode::FAILURE;
156 }
157
158
159 SmartDataPtr<EvtRecDTagCol> recDTagCol(eventSvc(), EventModel::EvtRec::EvtRecDTagCol);
160 if (!recDTagCol) {
161 log << MSG::FATAL << "EvtRecDTagCol is not registered yet" << endreq;
162 return StatusCode::FAILURE;
163 }
164
165
166
167 //get primary vertex from db
168 Hep3Vector xorigin(0,0,0);
169
170 IVertexDbSvc* vtxsvc;
171 Gaudi::svcLocator()->service("VertexDbSvc", vtxsvc);
172 if (vtxsvc->isVertexValid()) {
173
174 //vertex[0] = vx; vertex[1]= vy; vertex[2] = vz;
175 double* vertex = vtxsvc->PrimaryVertex();
176 xorigin.setX(vertex[0]);
177 xorigin.setY(vertex[1]);
178 xorigin.setZ(vertex[2]);
179 }
180
181 utility util;
182
183 /////////////////////////////
184 //reconstruct particle lists
185 /////////////////////////////
189
190 CDChargedPionList pionList(charged_begin, charged_end, pionSelector);
191 CDChargedKaonList kaonList(charged_begin, charged_end, kaonSelector);
192 CDProtonList protonList(charged_begin, charged_end, protonSelector);
193
194 CDPhotonList photonList(neutral_begin, neutral_end, photonSelector);
195
196 CDKsList ksList(ksSelector);
197 dc_fill(ksList, recVeeVertexCol->begin(), recVeeVertexCol->end());
198
199 // do a secondary vertex fit and cut on the results
200 map<EvtRecVeeVertex*, vector< double > > fitinfo_Ks;
201 for( CDKsList::iterator ksit = ksList.particle_begin(); ksit != ksList.particle_end(); ++ksit ){
202 EvtRecVeeVertex* ks = const_cast<EvtRecVeeVertex*>( (*ksit).particle().navKshort() );
203 fitinfo_Ks[ks] = util.SecondaryVFit(ks, vtxsvc);
204 }
205
206 //CDLambdaList lambdaList(lambdaSelector);
207 //dc_fill(lambdaList, recVeeVertexCol->begin(), recVeeVertexCol->end());
208 CDLambdaList lambdaList( recVeeVertexCol->begin(), recVeeVertexCol->end(), lambdaSelector);
209 map<EvtRecVeeVertex*, vector< double > > fitinfo_Lambda;
210 for( CDLambdaList::iterator lambdait = lambdaList.particle_begin(); lambdait != lambdaList.particle_end(); ++lambdait ){
211 EvtRecVeeVertex* lambda = const_cast<EvtRecVeeVertex*>( (*lambdait).particle().navLambda() );
212 fitinfo_Lambda[lambda] = util.SecondaryVFit_Lambda(lambda, vtxsvc);
213 }
214
215 CDPi0List pi0List(pi0Selector);
216 dc_fill(pi0List, recPi0Col->begin(), recPi0Col->end());
217
218 CDEtaList etaList(etatoGGSelector);
219 dc_fill(etaList, recEtaToGGCol->begin(), recEtaToGGCol->end());
220
221 //pion/kaon list with PID
225 CDChargedPionList pionList_tight(charged_begin, charged_end, pionSelector);
226 CDChargedKaonList kaonList_tight(charged_begin, charged_end, kaonSelector);
227 CDProtonList protonList_tight(charged_begin, charged_end, protonSelector);
228 int run = eventHeader->runNumber();
229 m_ievt = eventHeader->eventNumber();
230 m_nChrg = recEvent->totalCharged();
231 m_nNeu = recEvent->totalNeutral();
232 m_nPion = pionList.size();
233 m_nKaon = kaonList.size();
234 m_nProton = protonList.size();
235 m_nPi0 = pi0List.size();
236 m_nKs = ksList.size();
237 m_nLambda = lambdaList.size();
238
239
240 ///////////////////////
241 // get beam energy and beta
242 ///////////////////////
243
244
245 if(m_ReadBeamEFromDB && m_irun!=run){
246 m_irun=run;
247 if(m_usecalibBeamE)
248 m_readDb.setcalib(true);
249 m_beamE=m_readDb.getbeamE(m_irun,m_beamE);
250 if(run>0)
251 m_beta=m_readDb.getbeta();
252 // cout<<"use beam E from data base:"<<m_beamE<<endl;
253 }
254 double ebeam=m_beamE;
255
256
257 //////////////////////////////
258 //reconstruct decay lists
259 /////////////////////////////
260
261
262 for(int list=0;list<chanlist.size();list++){
263
264 string channel=chanlist[list];
265 vector<int> numchan;
267 lambdaCSelector.setbeta(m_beta);
268 CDDecayList decaylist(lambdaCSelector);
269
270 //K+/-: 1, Pi+/-:2, Pi0:3, Eta: 4, Ks:5, epTopipieta:6, epTorhogam:7, proton:8, Lambda:9, Sigma0:10, SigmaP:11, Omega:12
271 //the fist element of the vector stands for decay mode,
272 //the rest will be particles, and size of the vector minus 1 will be number of daughers.
273
274 if(channel=="LambdacPtoKsP") {
275 numchan.push_back( EvtRecDTag::kLambdacPtoKsP);
276 numchan.push_back(8);
277 numchan.push_back(5);
278 decaylist=protonList.plus()* ksList;
279 }
280 else if(channel=="LambdacPtoKPiP") {
281 numchan.push_back( EvtRecDTag::kLambdacPtoKPiP );
282 numchan.push_back(8);
283 numchan.push_back(1);
284 numchan.push_back(2);
285 decaylist=protonList.plus()* kaonList.minus()* pionList.plus();
286 }
287 else if(channel=="LambdacPtoKsPi0P") {
288 numchan.push_back( EvtRecDTag::kLambdacPtoKsPi0P);
289 numchan.push_back(8);
290 numchan.push_back(5);
291 numchan.push_back(3);
292 decaylist=protonList.plus()* ksList* pi0List;
293 }
294 else if(channel=="LambdacPtoKsPiPiP") {
295 numchan.push_back( EvtRecDTag::kLambdacPtoKsPiPiP);
296 numchan.push_back(8);
297 numchan.push_back(5);
298 numchan.push_back(2);
299 numchan.push_back(2);
300 decaylist=protonList.plus()* ksList* pionList.plus()* pionList.minus();
301 }
302 else if(channel=="LambdacPtoKPiPi0P") {
303 numchan.push_back( EvtRecDTag::kLambdacPtoKPiPi0P);
304 numchan.push_back(8);
305 numchan.push_back(1);
306 numchan.push_back(2);
307 numchan.push_back(3);
308 decaylist=protonList.plus()* kaonList.minus()* pionList.plus()* pi0List;
309 }
310 else if(channel=="LambdacPtoPiPiP") {
311 numchan.push_back( EvtRecDTag::kLambdacPtoPiPiP );
312 numchan.push_back(8);
313 numchan.push_back(2);
314 numchan.push_back(2);
315 decaylist=protonList.plus()* pionList.plus()* pionList.minus();
316 }
317 else if(channel=="LambdacPtoLambdaPi") {
318 numchan.push_back( EvtRecDTag::kLambdacPtoLambdaPi );
319 numchan.push_back(9);
320 numchan.push_back(2);
321 decaylist=lambdaList()* pionList.plus();
322 }
323 else if(channel=="LambdacPtoLambdaPiPi0") {
324 numchan.push_back( EvtRecDTag::kLambdacPtoLambdaPiPi0 );
325 numchan.push_back(9);
326 numchan.push_back(2);
327 numchan.push_back(3);
328 decaylist=lambdaList()* pionList.plus()* pi0List;
329 }
330 else if(channel=="LambdacPtoLambdaPiEta") {
331 numchan.push_back( EvtRecDTag::kLambdacPtoLambdaPiEta );
332 numchan.push_back(9);
333 numchan.push_back(2);
334 numchan.push_back(4);
335 decaylist=lambdaList()* pionList.plus()* etaList;;
336 }
337 else if(channel=="LambdacPtoLambdaPiPiPi") {
338 numchan.push_back( EvtRecDTag::kLambdacPtoLambdaPiPiPi );
339 numchan.push_back(9);
340 numchan.push_back(2);
341 numchan.push_back(2);
342 numchan.push_back(2);
343 decaylist=lambdaList()* pionList.plus()* pionList.minus() *pionList.plus();
344 }
345 else if(channel=="LambdacPtoLambdaPiOmega") {
346 numchan.push_back( EvtRecDTag::kLambdacPtoLambdaPiOmega );
347 numchan.push_back(9);
348 numchan.push_back(2);
349 numchan.push_back(12);
351 omgList=pionList.plus()* pionList.minus() *pi0List;
352 decaylist=lambdaList()* pionList.plus()* omgList;
353 }
354 else if(channel=="LambdacPtoPiSIGMA0LambdaGam") {
355 numchan.push_back( EvtRecDTag::kLambdacPtoPiSIGMA0LambdaGam );
356 numchan.push_back(2);
357 numchan.push_back(10);
359 sgmList=lambdaList()* photonList;
360 decaylist=pionList.plus()* sgmList;
361 }
362 else if(channel=="LambdacPtoPiPi0SIGMA0LambdaGam") {
364 numchan.push_back(2);
365 numchan.push_back(3);
366 numchan.push_back(10);
368 sgmList=lambdaList()* photonList;
369 decaylist=pionList.plus()* pi0List* sgmList;
370 }
371 else if(channel=="LambdacPtoPi0SIGMAPi0P") {
372 numchan.push_back( EvtRecDTag::kLambdacPtoPi0SIGMAPi0P );
373 numchan.push_back(3);
374 numchan.push_back(11);
376 sgmList=pi0List* protonList.plus();
377 decaylist=pi0List* sgmList;
378 }
379 else if(channel=="LambdacPtoPiPiSIGMAPi0P") {
380 numchan.push_back( EvtRecDTag::kLambdacPtoPiPiSIGMAPi0P );
381 numchan.push_back(2);
382 numchan.push_back(2);
383 numchan.push_back(11);
385 sgmList=pi0List* protonList.plus();
386 decaylist=pionList.plus()* pionList.minus()* sgmList;
387 }
388 else if(channel=="LambdacPtoOmegaSIGMAPi0P") {
389 numchan.push_back( EvtRecDTag::kLambdacPtoOmegaSIGMAPi0P );
390 numchan.push_back(12);
391 numchan.push_back(11);
393 omgList=pionList.plus()* pionList.minus()* pi0List;
395 sgmList=pi0List* protonList.plus();
396 decaylist=omgList* sgmList;
397 }
398 else if(channel=="LambdacPtoPiPiPi0SIGMAPi0P") {
399 numchan.push_back( EvtRecDTag::kLambdacPtoPiPiPi0SIGMAPi0P );
400 numchan.push_back(2);
401 numchan.push_back(2);
402 numchan.push_back(3);
403 numchan.push_back(11);
405 sgmList=pi0List* protonList.plus();
406 decaylist=pionList.plus()* pionList.minus()*pi0List* sgmList;
407 }
408
409
410 CDDecayList::iterator Lc_begin =decaylist.particle_begin();
411 CDDecayList::iterator Lc_end =decaylist.particle_end();
412
413 for ( CDDecayList::iterator it = Lc_begin; it != Lc_end; it++ ) {
414
415 EvtRecDTag* recDTag = new EvtRecDTag;
416 recDTag->setdecayMode( (EvtRecDTag::DecayMode)numchan[0] );
417
418 vector<int> trackid, showerid;
419 vector<int> kaonid, pionid, protonid;
420 int numofchildren=numchan.size()-1;
421
422 for(int i=0; i< numofchildren;i++){
423
424 const CDCandidate& daughter=(*it).particle().child(i);
425
426 if(numchan[i+1]==1){
427 const EvtRecTrack* track=daughter.track();
428 trackid.push_back(track->trackId());
429 kaonid.push_back(track->trackId());
430 }
431 else if(numchan[i+1]==2){
432 const EvtRecTrack* track=daughter.track();
433 trackid.push_back(track->trackId());
434 pionid.push_back(track->trackId());
435 }
436 else if(numchan[i+1]==8){
437 const EvtRecTrack* track=daughter.track();
438 trackid.push_back(track->trackId());
439 protonid.push_back(track->trackId());
440 }
441 else if ( numchan[i+1]==3){
442 const EvtRecTrack* hiEnGamma=daughter.navPi0()->hiEnGamma();
443 const EvtRecTrack* loEnGamma=daughter.navPi0()->loEnGamma();
444 showerid.push_back(hiEnGamma->trackId());
445 showerid.push_back(loEnGamma->trackId());
446 }
447 else if ( numchan[i+1]==4){
448 const EvtRecTrack* hiEnGamma=daughter.navEta()->hiEnGamma();
449 const EvtRecTrack* loEnGamma=daughter.navEta()->loEnGamma();
450 showerid.push_back(hiEnGamma->trackId());
451 showerid.push_back(loEnGamma->trackId());
452 }
453 else if ( numchan[i+1]==5){
454 EvtRecVeeVertex* aKsCand = const_cast<EvtRecVeeVertex*>( daughter.navKshort() );
455 recDTag->addToFitInfo(aKsCand->mass(),fitinfo_Ks[aKsCand][0],fitinfo_Ks[aKsCand][1],fitinfo_Ks[aKsCand][2]);
456 EvtRecTrack* pion1Trk = aKsCand->daughter(0);
457 EvtRecTrack* pion2Trk = aKsCand->daughter(1);
458 trackid.push_back(pion1Trk->trackId());
459 trackid.push_back(pion2Trk->trackId());
460 }
461 else if ( numchan[i+1]==9){
462 EvtRecVeeVertex* aLambdaCand = const_cast<EvtRecVeeVertex*>( daughter.navLambda() );
463 recDTag->addToFitInfo(aLambdaCand->mass(),fitinfo_Lambda[aLambdaCand][0],fitinfo_Lambda[aLambdaCand][1],fitinfo_Lambda[aLambdaCand][2]);
464
465 int index[2] = {0, 1};
466 if (aLambdaCand->vertexId() < 0){
467 index[0] = 1;
468 index[1] = 0;
469 }
470
471 EvtRecTrack* protonTrk = aLambdaCand->daughter(index[0]);
472 EvtRecTrack* pionTrk = aLambdaCand->daughter(index[1]);
473 trackid.push_back(protonTrk->trackId());
474 trackid.push_back(pionTrk->trackId());
475 }
476 else if (numchan[i+1]==6){
477 const CDCandidate& apion = daughter.decay().child(0);
478 const CDCandidate& spion = daughter.decay().child(1);
479 const CDCandidate& eta = daughter.decay().child(2);
480 const EvtRecTrack* apiontrk = apion.track();
481 const EvtRecTrack* spiontrk = spion.track();
482 const EvtRecTrack* hiEnGamma=eta.navEta()->hiEnGamma();
483 const EvtRecTrack* loEnGamma=eta.navEta()->loEnGamma();
484
485 trackid.push_back(apiontrk->trackId());
486 trackid.push_back(spiontrk->trackId());
487 showerid.push_back(hiEnGamma->trackId());
488 showerid.push_back(loEnGamma->trackId());
489
490 }
491 else if (numchan[i+1]==7){
492 const CDCandidate& rho = daughter.decay().child(0);
493 const CDCandidate& apion = rho.decay().child(0);
494 const CDCandidate& spion = rho.decay().child(1);
495 const CDCandidate& gamma = daughter.decay().child(1);
496
497 const EvtRecTrack* apiontrk = apion.track();
498 const EvtRecTrack* spiontrk = spion.track();
499 const EvtRecTrack* gammatrk = gamma.photon();
500
501
502 trackid.push_back(apiontrk->trackId());
503 trackid.push_back(spiontrk->trackId());
504 showerid.push_back(gammatrk->trackId());
505
506 }
507 else if (numchan[i+1]==10){
508 const CDCandidate& gamma = daughter.decay().child(1);
509 const EvtRecTrack* gammatrk = gamma.photon();
510 showerid.push_back(gammatrk->trackId());
511
512 EvtRecVeeVertex* aLambdaCand = const_cast<EvtRecVeeVertex*>( daughter.decay().child(0).navLambda() );
513 recDTag->addToFitInfo(aLambdaCand->mass(),fitinfo_Lambda[aLambdaCand][0],fitinfo_Lambda[aLambdaCand][1],fitinfo_Lambda[aLambdaCand][2]);
514
515 int index[2] = {0, 1};
516 if (aLambdaCand->vertexId() < 0){
517 index[0] = 1;
518 index[1] = 0;
519 }
520
521 EvtRecTrack* protonTrk = aLambdaCand->daughter(index[0]);
522 EvtRecTrack* pionTrk = aLambdaCand->daughter(index[1]);
523 trackid.push_back(protonTrk->trackId());
524 trackid.push_back(pionTrk->trackId());
525
526 }
527 else if (numchan[i+1]==11){
528 const CDCandidate& pi0 = daughter.decay().child(0);
529 const CDCandidate& agam = pi0.decay().child(0);
530 const CDCandidate& sgam = pi0.decay().child(1);
531 const CDCandidate& P = daughter.decay().child(1);
532
533 const EvtRecTrack* agamtrk = agam.photon();
534 const EvtRecTrack* sgamtrk = sgam.photon();
535 const EvtRecTrack* Ptrk = P.track();
536
537
538 showerid.push_back(agamtrk->trackId());
539 showerid.push_back(sgamtrk->trackId());
540 trackid.push_back(Ptrk->trackId());
541 }
542 else if (numchan[i+1]==12){
543 const CDCandidate& apion = daughter.decay().child(0);
544 const CDCandidate& spion = daughter.decay().child(1);
545 const CDCandidate& pi0 = daughter.decay().child(2);
546 const EvtRecTrack* apiontrk = apion.track();
547 const EvtRecTrack* spiontrk = spion.track();
548 const EvtRecTrack* hiEnGamma=pi0.navPi0()->hiEnGamma();
549 const EvtRecTrack* loEnGamma=pi0.navPi0()->loEnGamma();
550
551 trackid.push_back(apiontrk->trackId());
552 trackid.push_back(spiontrk->trackId());
553 showerid.push_back(hiEnGamma->trackId());
554 showerid.push_back(loEnGamma->trackId());
555
556 }
557
558 }//end of filling track and shower ids
559
560
561 saveLcInfo(it, ebeam, numofchildren, recDTag);
562 savetrack(trackid,showerid,charged_begin,charged_end,neutral_begin,neutral_end,recDTag);
563 pidtag(kaonid,pionid,protonid, kaonList_tight, pionList_tight, protonList_tight, recDTag);
564
565 if(m_usevertexfit){
566 HepLorentzVector p4change_vfit=util.vfit(channel, kaonid, pionid, protonid, xorigin, charged_begin);
567 recDTag->setp4(recDTag->p4()+p4change_vfit);
568 }
569
570
571 trackid.clear();
572 showerid.clear();
573 kaonid.clear();
574 pionid.clear();
575 protonid.clear();
576
577 //write dtag object out
578 recDTagCol->push_back(recDTag);
579
580 }//end of decaylist iterator
581
582 numchan.clear();
583
584 }//end of reconstrucing all D+ decay lists
585
586 return StatusCode::SUCCESS;
587}
588
589
590void LambdaCReconstruction::saveLcInfo(CDDecayList::iterator it, double ebeam, int numofchildren, EvtRecDTag* recDTag){
591
592 double mass = (*it).particle().mass();
593 int charge= (*it).particle().charge();
594 HepLorentzVector p4((*it).particle().momentum(), (*it).particle().energy());
595 recDTag->setp4(p4);
596
597 p4.boost(-m_beta);
598 double mbc2_CMS = ebeam*ebeam - p4.v().mag2();
599 double mbc_CMS = mbc2_CMS > 0 ? sqrt( mbc2_CMS ) : -10;
600 double deltaE_CMS = p4.t() - ebeam;
601
602 if((*it).particle().userTag()==1)
603 recDTag->settype( EvtRecDTag::Tight );
604 else
605 recDTag->settype( EvtRecDTag::Loose );
606 recDTag->setcharge(charge);
607 recDTag->setcharm(charge);
608 recDTag->setnumOfChildren(numofchildren);
609 recDTag->setmass(mass);
610 recDTag->setmBC(mbc_CMS);
611 recDTag->setbeamE(ebeam);
612 recDTag->setdeltaE(deltaE_CMS);
613
614}
615
616void LambdaCReconstruction::savetrack(vector<int> trackid, vector<int> showerid, EvtRecTrackIterator charged_begin, EvtRecTrackIterator charged_end,
617 EvtRecTrackIterator neutral_begin, EvtRecTrackIterator neutral_end,EvtRecDTag* recDTag){
618
619 vector<EvtRecTrackIterator> trktemp;
620 vector<EvtRecTrackIterator> shrtemp;
621
622 //fill tracks
623 for(EvtRecTrackIterator trk=charged_begin; trk<charged_end;trk++){
624
625 bool isothertrack=true;
626 for(int i=0; i<trackid.size(); i++){
627 if( (*trk)->trackId()==trackid[i] ){
628 trktemp.push_back(trk);
629 isothertrack=false;
630 break;
631 }
632 }
633 if(isothertrack)
634 recDTag->addOtherTrack(*trk);
635 }
636 for(int i=0; i<trackid.size();i++){
637 for(int j=0; j<trktemp.size(); j++){
638 EvtRecTrackIterator trk=trktemp[j];
639 if( (*trk)->trackId()==trackid[i])
640 recDTag->addTrack(*trktemp[j]);
641 }
642 }
643
644
645 //fill showers
646 for(EvtRecTrackIterator shr=neutral_begin; shr<neutral_end;shr++){
647 bool isothershower=true;
648 for(int i=0; i<showerid.size(); i++){
649 if( (*shr)->trackId()==showerid[i] ){
650 shrtemp.push_back(shr);
651 isothershower=false;
652 break;
653 }
654 }
655 if(isothershower)
656 recDTag->addOtherShower(*shr);
657 }
658
659 for(int i=0; i<showerid.size();i++){
660 for(int j=0; j<shrtemp.size(); j++){
661 EvtRecTrackIterator shr=shrtemp[j];
662 if( (*shr)->trackId()==showerid[i])
663 recDTag->addShower(*shrtemp[j]);
664 }
665 }
666
667
668}
669
670void LambdaCReconstruction::pidtag(vector<int> kaonid, vector<int> pionid, vector<int> protonid, CDChargedKaonList& kaonList, CDChargedPionList& pionList, CDProtonList& protonList,EvtRecDTag* recDTag){
671
672 bool iskaon=false,ispion=false,isproton=false;
673
674
675 // save track ids which passed pion/kaon cuts
676
677 for (CDChargedKaonList::iterator kit = kaonList.particle_begin(); kit != kaonList.particle_end(); kit++) {
678 recDTag->addKaonId( (*kit).particle().track() );
679 }
680
681 for (CDChargedPionList::iterator pit = pionList.particle_begin(); pit != pionList.particle_end(); pit++) {
682 recDTag->addPionId( (*pit).particle().track() );
683 }
684 for (CDProtonList::iterator pt = protonList.particle_begin(); pt != protonList.particle_end(); pt++) {
685 recDTag->addProtonId( (*pt).particle().track() );
686 }
687
688
689 /*
690 for(int i=0; i<kaonid.size(); i++){
691 bool ithkaon=false;
692 for (CDChargedKaonList::iterator kit = kaonList.particle_begin(); kit != kaonList.particle_end(); kit++) {
693 if((*kit).particle().track()->trackId()==kaonid[i]){
694 ithkaon=true;
695 break;
696 }
697 }
698 if(!ithkaon) break;
699 if(i==kaonid.size()-1)
700 iskaon=true;
701 }
702
703 for(int i=0; i<pionid.size(); i++){
704 bool ithpion=false;
705 for (CDChargedPionList::iterator pit = pionList.particle_begin(); pit != pionList.particle_end(); pit++) {
706 if((*pit).particle().track()->trackId()==pionid[i]){
707 ithpion=true;
708 break;
709 }
710 }
711 if(!ithpion) break;
712 if(i==pionid.size()-1)
713 ispion=true;
714 }
715
716
717 if( iskaon && ispion)
718 recDTag->settype( EvtRecDTag::Tight );
719 else if( (kaonid.size()==0 && ispion) || (pionid.size()==0 && iskaon))
720 recDTag->settype( EvtRecDTag::Tight );
721 else if( kaonid.size()==0 && pionid.size()==0 )
722 recDTag->settype( EvtRecDTag::Tight );
723 */
724
725}
726
727
728
729vector<string> LambdaCReconstruction::getlist(string& filename ){
730
731 string channel;
732 vector<string> temp;
733
734 ifstream inFile;
735
736 inFile.open(filename.c_str());
737 if (!inFile) {
738 cout << "Unable to open decay list file";
739 exit(1); // terminate with error
740 }
741
742 while (inFile >> channel) {
743 temp.push_back(channel);
744 }
745
746 inFile.close();
747
748 return temp;
749
750}
double P(RecMdcKalTrack *trk)
double mass
void dc_fill(DCFillableChargedList< Charged > &aFillableList, WitnessIterator first, WitnessIterator last)
double ebeam
ObjectVector< EvtRecDTag > EvtRecDTagCol
Definition: EvtRecDTag.h:286
EvtRecTrackCol::iterator EvtRecTrackIterator
Definition: EvtRecTrack.h:111
LambdaCSelector lambdaCSelector
LocalChargedSigmaSelector chargedsigmaSelector
LocalEtatoGGSelector etatoGGSelector
LocalKaonSelector kaonSelector
LocalKsSelector ksSelector
LocalLambdaSelector lambdaSelector
LocalPhotonSelector photonSelector
LocalPi0Selector pi0Selector
LocalPionSelector pionSelector
LocalProtonSelector protonSelector
LocalSigma0Selector sigma0Selector
LocalomegatoPiPiPi0Selector omegatoPiPiPi0Selector
IMessageSvc * msgSvc()
virtual const DecayEvidence & decay() const
virtual const EvtRecTrack * photon() const
virtual const EvtRecVeeVertex * navLambda() const
virtual const EvtRecTrack * track() const
virtual const EvtRecVeeVertex * navKshort() const
virtual const EvtRecPi0 * navPi0() const
virtual const EvtRecEtaToGG * navEta() const
const CDCandidate & child(unsigned int aPosition) const
Definition: CDDecay.cxx:247
@ kLambdacPtoKPiP
Definition: EvtRecDTag.h:161
@ kLambdacPtoLambdaPiOmega
Definition: EvtRecDTag.h:171
@ kLambdacPtoKsP
Definition: EvtRecDTag.h:160
@ kLambdacPtoPiSIGMA0LambdaGam
Definition: EvtRecDTag.h:173
@ kLambdacPtoKPiPi0P
Definition: EvtRecDTag.h:164
@ kLambdacPtoOmegaSIGMAPi0P
Definition: EvtRecDTag.h:177
@ kLambdacPtoKsPi0P
Definition: EvtRecDTag.h:162
@ kLambdacPtoLambdaPi
Definition: EvtRecDTag.h:167
@ kLambdacPtoPi0SIGMAPi0P
Definition: EvtRecDTag.h:175
@ kLambdacPtoLambdaPiPi0
Definition: EvtRecDTag.h:168
@ kLambdacPtoPiPiSIGMAPi0P
Definition: EvtRecDTag.h:176
@ kLambdacPtoLambdaPiEta
Definition: EvtRecDTag.h:169
@ kLambdacPtoKsPiPiP
Definition: EvtRecDTag.h:163
@ kLambdacPtoPiPiP
Definition: EvtRecDTag.h:165
@ kLambdacPtoPiPiPi0SIGMAPi0P
Definition: EvtRecDTag.h:178
@ kLambdacPtoPiPi0SIGMA0LambdaGam
Definition: EvtRecDTag.h:174
@ kLambdacPtoLambdaPiPiPi
Definition: EvtRecDTag.h:170
void addOtherTrack(const SmartRef< EvtRecTrack > track)
Definition: EvtRecDTag.h:232
void settype(SelType type)
Definition: EvtRecDTag.h:211
void setdecayMode(DecayMode decayMode)
Definition: EvtRecDTag.h:210
HepLorentzVector p4() const
Definition: EvtRecDTag.h:194
void setp4(HepLorentzVector p4)
Definition: EvtRecDTag.h:219
void setmass(double mass)
Definition: EvtRecDTag.h:213
void addToFitInfo(double ksmass, double chi2, double length, double error)
Definition: EvtRecDTag.h:221
void addOtherShower(const SmartRef< EvtRecTrack > shower)
Definition: EvtRecDTag.h:234
void addKaonId(const SmartRef< EvtRecTrack > kaonId)
Definition: EvtRecDTag.h:238
void setdeltaE(double deltaE)
Definition: EvtRecDTag.h:215
void setcharm(int charm)
Definition: EvtRecDTag.h:217
void setmBC(double mBC)
Definition: EvtRecDTag.h:214
void setcharge(int charge)
Definition: EvtRecDTag.h:216
void addPionId(const SmartRef< EvtRecTrack > pionId)
Definition: EvtRecDTag.h:236
void setbeamE(double beamE)
Definition: EvtRecDTag.h:212
void addProtonId(const SmartRef< EvtRecTrack > protonId)
Definition: EvtRecDTag.h:240
void addShower(const SmartRef< EvtRecTrack > shower)
Definition: EvtRecDTag.h:230
void setnumOfChildren(int numOfChildren)
Definition: EvtRecDTag.h:218
void addTrack(const SmartRef< EvtRecTrack > track)
Definition: EvtRecDTag.h:228
const EvtRecTrack * hiEnGamma() const
Definition: EvtRecEtaToGG.h:30
const EvtRecTrack * loEnGamma() const
Definition: EvtRecEtaToGG.h:31
const EvtRecTrack * loEnGamma() const
Definition: EvtRecPi0.h:31
const EvtRecTrack * hiEnGamma() const
Definition: EvtRecPi0.h:30
int trackId() const
Definition: EvtRecTrack.h:32
SmartRef< EvtRecTrack > & daughter(int i)
int vertexId() const
double mass() const
virtual bool isVertexValid()=0
virtual double * PrimaryVertex()=0
void savetrack(vector< int >, vector< int >, EvtRecTrackIterator, EvtRecTrackIterator, EvtRecTrackIterator, EvtRecTrackIterator, EvtRecDTag *)
void pidtag(vector< int >, vector< int >, vector< int >, CDChargedKaonList &, CDChargedPionList &, CDProtonList &, EvtRecDTag *)
void saveLcInfo(CDDecayList::iterator, double, int, EvtRecDTag *)
vector< string > getlist(string &filename)
LambdaCReconstruction(const std::string &name, ISvcLocator *pSvcLocator)
void setebeam(double ebeam)
void setbeta(Hep3Vector beta)
void setpidtype(int type)
void setpidtype(int type)
void setpidtype(int type)
virtual int size() const
ChosenChargeList< Charged, CandidateClass > & plus() const
ChosenChargeList< Charged, CandidateClass > & minus() const
virtual iterator particle_begin()
Definition: DecayList.cc:169
virtual iterator particle_end()
Definition: DecayList.cc:176
vector< double > SecondaryVFit_Lambda(EvtRecVeeVertex *lambda, IVertexDbSvc *vtxsvc)
Definition: utility.cxx:295
HepLorentzVector vfit(string channel, vector< int > kaonid, vector< int > pionid, HepPoint3D vx, EvtRecTrackIterator charged_begin)
Definition: utility.cxx:59
vector< double > SecondaryVFit(EvtRecVeeVertex *ks, IVertexDbSvc *vtxsvc)
Definition: utility.cxx:195
_EXTERN_ std::string EvtRecEvent
Definition: EventModel.h:116
_EXTERN_ std::string EvtRecDTagCol
Definition: EventModel.h:122
_EXTERN_ std::string EvtRecTrackCol
Definition: EventModel.h:117
Definition: Event.h:21
float charge