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