BOSS 6.6.4.p01
BESIII Offline Software System
Loading...
Searching...
No Matches
BesMcTruthWriter Class Reference

#include <BesMcTruthWriter.hh>

Public Member Functions

 BesMcTruthWriter ()
 
 ~BesMcTruthWriter ()
 
void SaveMcTruth ()
 
void SaveMcParticle ()
 
void AddMother (Event::McParticle *, Event::McParticleCol *)
 
void SaveDecayMode ()
 
void SaveMdcTruth ()
 
void SaveTofTruth ()
 
void SaveEmcTruth ()
 
void SaveMucTruth ()
 

Detailed Description

Definition at line 23 of file BesMcTruthWriter.hh.

Constructor & Destructor Documentation

◆ BesMcTruthWriter()

BesMcTruthWriter::BesMcTruthWriter ( )

Definition at line 52 of file BesMcTruthWriter.cc.

53{
54 m_DigiMan = G4DigiManager::GetDMpointer();
55 //mdcGeoPointer=BesMdcGeoParameter::GetGeo();
56}

◆ ~BesMcTruthWriter()

BesMcTruthWriter::~BesMcTruthWriter ( )

Definition at line 58 of file BesMcTruthWriter.cc.

59{
60}

Member Function Documentation

◆ AddMother()

void BesMcTruthWriter::AddMother ( Event::McParticle currentParticle,
Event::McParticleCol particleCol 
)

Definition at line 291 of file BesMcTruthWriter.cc.

292{
293 if (currentParticle->leafParticle()) {
294 return;
295 }
296 Event::McParticleCol::iterator iter;
297 bool found = false;
298 for ( iter = particleCol->begin(); iter != particleCol->end(); iter++) {
299 if (currentParticle->vertexIndex1() == (*iter)->vertexIndex0()) {
300 found = true;
301 (*iter)->setMother(currentParticle);
302 currentParticle->addDaughter(*iter);
303 AddMother((*iter), particleCol);
304 }
305 }
306 if (!found) std::cout << "AddMother: inconsistency was found!" << std::endl;
307
308}
EvtStreamInputIterator< typename Generator::result_type > iter(Generator gen, int N=0)
void AddMother(Event::McParticle *, Event::McParticleCol *)
int vertexIndex1() const
Definition: McParticle.h:125
bool leafParticle() const
Retrieve whether this is a leaf particle.
Definition: McParticle.cxx:19
void addDaughter(const SmartRef< McParticle > d)
add a daugther particle to this particle
Definition: McParticle.h:147

Referenced by AddMother(), and SaveMcParticle().

◆ SaveDecayMode()

void BesMcTruthWriter::SaveDecayMode ( )

Definition at line 310 of file BesMcTruthWriter.cc.

311{
312 SmartDataPtr<DecayMode> decayMode(m_evtSvc,"/Event/MC/DecayMode");
313 if(!decayMode)
314 {
315 //G4cout<<"Could not retrieve DecayMode"<<G4endl;
316 DecayMode* aDecayMode = new DecayMode;
317 //register Decay Mode to TDS
318 int dm[10]={0,0,0,0,0,0,0,0,0,0};
319 aDecayMode->putData(dm,10);
320 StatusCode scDM = m_evtSvc->registerObject("/Event/MC/DecayMode", aDecayMode);
321 if(scDM!=StatusCode::SUCCESS)
322 G4cout<< "Could not register DecayMode" <<G4endl;
323 }
324}
void putData(int *data, unsigned int size)

Referenced by SaveMcTruth().

◆ SaveEmcTruth()

void BesMcTruthWriter::SaveEmcTruth ( )

Definition at line 497 of file BesMcTruthWriter.cc.

498{
499 //Emc McHit collection defined in BOSS
500 Event::EmcMcHitCol* aEmcMcHitCol = new Event::EmcMcHitCol;
501
502 G4int fullMc = 1;
503 if(fullMc==1) { //full Mc Truth
504 G4int HCID = -1;
505 HCID = m_DigiMan->GetHitsCollectionID("BesEmcTruthHitsList");
506 if(HCID>0)
507 {
509 HC = (BesEmcTruthHitsCollection*) (m_DigiMan->GetHitsCollection(HCID));
510 G4int n_hit = HC->entries();
511 if(n_hit>0)
512 {
513 //arrange hits in hits collection in order of trackIndex
514 BesEmcTruthHit* hit;
515 vector<BesEmcTruthHit*>* vecHC = HC->GetVector();
516 for(int i=0;i<n_hit-1;i++) {
517 for(int j=i+1;j<n_hit;j++) {
518 if((*vecHC)[i]->GetTrackIndex()>(*vecHC)[j]->GetTrackIndex()) {
519 hit = (*vecHC)[i];
520 (*vecHC)[i] = (*vecHC)[j];
521 (*vecHC)[j] = hit;
522 }
523 }
524 }
525
526 for(G4int i=0;i<n_hit;i++)
527 {
528 BesEmcTruthHit* hit;
529 hit = (*HC)[i];
530
531 std::map<Identifier,double> hitMap;
532 std::map<Identifier,double>::const_iterator iHitMap;
533 hitMap.clear();
534 for(iHitMap=hit->Begin();iHitMap!=hit->End();iHitMap++) {
535 hitMap[iHitMap->first]=iHitMap->second;
536 }
537
538 Event::EmcMcHit* emcMcHit = new Event::EmcMcHit(hit->GetIdentify(), hit->GetTrackIndex(),
539 hit->GetPosition().x(), hit->GetPosition().y(), hit->GetPosition().z(),
540 hit->GetMomentum().x(), hit->GetMomentum().y(), hit->GetMomentum().z(),
541 hit->GetEDep() );
542
543 emcMcHit->setHitEmc(hit->GetHitEmc());
544 emcMcHit->setPDGCode(hit->GetPDGCode());
545 emcMcHit->setPDGCharge(hit->GetPDGCharge());
546 emcMcHit->setTime(hit->GetTime());
547 emcMcHit->setHitMap(hitMap);
548
549 aEmcMcHitCol->push_back(emcMcHit);
550 }
551 }
552 }
553 } else { //simple Mc Truth
554 G4int HCID = -1;
555 HCID = m_DigiMan->GetHitsCollectionID("BesEmcHitsList");
556 if(HCID>0)
557 {
558 BesEmcHitsCollection* HC = 0;
559 HC = (BesEmcHitsCollection*) (m_DigiMan->GetHitsCollection(HCID));
560 G4int n_hit = HC->entries();
561 if(n_hit>0)
562 {
563 //arrange hits in hits collection in order of trackIndex
564 BesEmcHit* hit;
565 vector<BesEmcHit*>* vecHC = HC->GetVector();
566 for(int i=0;i<n_hit-1;i++)
567 for(int j=i+1;j<n_hit;j++)
568 if((*vecHC)[i]->GetTrackIndex()>(*vecHC)[j]->GetTrackIndex())
569 {
570 hit = (*vecHC)[i];
571 (*vecHC)[i] = (*vecHC)[j];
572 (*vecHC)[j] = hit;
573 }
574
575 //Emc McHit collection defined in BOSS
576 Event::EmcMcHitCol* aEmcMcHitCol = new Event::EmcMcHitCol;
577
578 for(G4int i=0;i<n_hit;i++)
579 {
580 hit = (*HC)[i];
582
583 std::map<Identifier,double> hitMap;
584 Event::EmcMcHit* emcMcHit = new Event::EmcMcHit(ident, hit->GetTrackIndex(),
585 hit->GetPosCrystal().x(), hit->GetPosCrystal().y(), hit->GetPosCrystal().z(),
586 hit->GetMomentum().x(), hit->GetMomentum().y(), hit->GetMomentum().z(),
587 hit->GetEdepCrystal() );
588
589 aEmcMcHitCol->push_back(emcMcHit);
590 }
591 }
592 }
593 }
594
595 //register EMC McTruth collection to TDS
596 StatusCode scEmc = m_evtSvc->registerObject("/Event/MC/EmcMcHitCol", aEmcMcHitCol);
597 if(scEmc!=StatusCode::SUCCESS)
598 G4cout<< "Could not register EMC McTruth collection" <<G4endl;
599
600 //retrieve EMC McTruth from TDS
601 /*SmartDataPtr<Event::EmcMcHitCol> aMcHitCol(m_evtSvc,"/Event/MC/EmcMcHitCol");
602 if(!aMcHitCol)
603 G4cout<<"Could not retrieve EMC McTruth collection"<<G4endl;
604
605 Event::EmcMcHitCol::iterator iMcHitCol;
606 for(iMcHitCol=aMcHitCol->begin(); iMcHitCol!=aMcHitCol->end(); iMcHitCol++)
607 {
608 const Identifier ident = (*iMcHitCol)->identify();
609 G4cout<<(*iMcHitCol)->getTrackIndex();
610 G4cout<<" "<<EmcID::barrel_ec(ident);
611 G4cout<<" "<<EmcID::theta_module(ident);
612 G4cout<<" "<<EmcID::phi_module(ident);
613 G4cout<<" "<<(*iMcHitCol)->getPositionX();
614 G4cout<<" "<<(*iMcHitCol)->getPositionY();
615 G4cout<<" "<<(*iMcHitCol)->getPositionZ();
616 G4cout<<" "<<(*iMcHitCol)->getPx();
617 G4cout<<" "<<(*iMcHitCol)->getPy();
618 G4cout<<" "<<(*iMcHitCol)->getPz();
619 G4cout<<" "<<(*iMcHitCol)->getDepositEnergy();
620 G4cout<<G4endl;
621 }
622 G4cout<<"end of retrieve EMC McTruth"<<G4endl;
623 */
624}
G4THitsCollection< BesEmcTruthHit > BesEmcTruthHitsCollection
Definition: BesEmcHit.hh:181
G4THitsCollection< BesEmcHit > BesEmcHitsCollection
Definition: BesEmcHit.hh:83
G4int GetNumPhiCrystal()
Definition: BesEmcHit.hh:63
G4ThreeVector GetPosCrystal()
Definition: BesEmcHit.hh:59
G4double GetEdepCrystal()
Definition: BesEmcHit.hh:56
G4int GetTrackIndex()
Definition: BesEmcHit.hh:64
G4ThreeVector GetMomentum()
Definition: BesEmcHit.hh:66
G4int GetNumThetaCrystal()
Definition: BesEmcHit.hh:62
G4int GetPartId()
Definition: BesEmcHit.hh:61
G4int GetTrackIndex() const
Definition: BesEmcHit.hh:136
G4int GetHitEmc() const
Definition: BesEmcHit.hh:138
G4ThreeVector GetPosition() const
Definition: BesEmcHit.hh:145
G4double GetTime() const
Definition: BesEmcHit.hh:143
std::map< Identifier, G4double >::const_iterator End() const
Definition: BesEmcHit.cc:182
G4double GetPDGCharge() const
Definition: BesEmcHit.hh:140
G4double GetEDep() const
Definition: BesEmcHit.hh:142
std::map< Identifier, G4double >::const_iterator Begin() const
Definition: BesEmcHit.cc:177
G4int GetPDGCode() const
Definition: BesEmcHit.hh:139
G4ThreeVector GetMomentum() const
Definition: BesEmcHit.hh:144
Identifier GetIdentify() const
Definition: BesEmcHit.hh:135
static Identifier crystal_id(const unsigned int barrel_ec, const unsigned int theta_module, const unsigned int phi_module)
For a single crystal.
Definition: EmcID.cxx:71
void setPDGCode(int code)
Definition: EmcMcHit.h:75
void setHitEmc(int is)
Definition: EmcMcHit.h:74
void setPDGCharge(double charge)
Definition: EmcMcHit.h:76
void setTime(double time)
Definition: EmcMcHit.h:77
void setHitMap(std::map< Identifier, double > &hitMap)
Definition: EmcMcHit.h:94
ObjectVector< EmcMcHit > EmcMcHitCol
Definition: EmcMcHit.h:132

Referenced by SaveMcTruth().

◆ SaveMcParticle()

void BesMcTruthWriter::SaveMcParticle ( )

Definition at line 87 of file BesMcTruthWriter.cc.

88{
90 vector<BesTruthTrack*>* trackList = sensitiveManager->GetTrackList();
91 G4int nTrack = trackList->size();
92 BesTruthTrack* track;
93
94 vector<BesTruthVertex*>* vertexList = sensitiveManager->GetVertexList();
95 G4int nVertex = vertexList->size();
96 BesTruthVertex* vertex;
97
98 //arrange TruthTrack in trackList in order of trackIndex
99 for(int i=0;i<nTrack-1;i++)
100 for(int j=i+1;j<nTrack;j++)
101 if((*trackList)[i]->GetIndex()>(*trackList)[j]->GetIndex())
102 {
103 track=(*trackList)[i];
104 (*trackList)[i]=(*trackList)[j];
105 (*trackList)[j]=track;
106 }
107
109
110 //loop over tracks
111 for(int i=0;i<nTrack;i++)
112 {
113 track = (*trackList) [i];
114
115 // find out the start point
116 bool isPrimary = false;
117 bool startPointFound = false;
118 BesTruthVertex* startPoint;
119
120 for (int i=0;i<nVertex;i++)
121 {
122 vertex = (*vertexList) [i];
123 if ( track->GetVertex()->GetIndex() == vertex->GetIndex() )
124 {
125 if (!vertex->GetParentTrack()) isPrimary = true;
126 startPointFound = true;
127 startPoint = vertex;
128 break;
129 }
130 }
131
132 if (!startPointFound)
133 std::cout << "error in finding the start point of a track" <<std::endl;
134
135
136 bool endPointFound = false;
137 // find out the end point
138 for (int i=0;i<nVertex;i++)
139 {
140 vertex = (*vertexList) [i];
141 if( track->GetTerminalVertex() )
142 {
143 if (track->GetTerminalVertex()->GetIndex() == vertex->GetIndex() )
144 {
145 //create the mc particle
146 Event::McParticle* mcParticle = new Event::McParticle;
147
148 // initialization
149 HepLorentzVector initialMomentum(track->GetP4().x()/1000., track->GetP4().y()/1000., track->GetP4().z()/1000., track->GetP4().t()/1000.);
150
151 HepLorentzVector initialPosition(startPoint->GetPosition().x()/10., startPoint->GetPosition().y()/10., startPoint->GetPosition().z()/10.,startPoint->GetTime());
152
153 mcParticle->initialize(track->GetPDGCode(), 0, initialMomentum, initialPosition);
154
155 // Set track index
156 mcParticle->setTrackIndex( track->GetIndex() );
157
158 // Adding status flag
159 if (isPrimary ) mcParticle->addStatusFlag(Event::McParticle::PRIMARY);
160
161 if (track->GetDaughterIndexes().size()==0)
163
164 //std::cout<<"pdg:"<<track->GetPDGCode()<<" "
165 // <<"source:"<<track->GetSource()<<" "
166 // <<std::endl;
167 if (track->GetSource()=="FromGenerator")
169 else if(track->GetSource()=="FromG4")
171 else
173
174 //std::cout<<"status:"<<mcParticle->statusFlags()<<" "
175 // <<"FGN:"<<mcParticle->decayFromGenerator()<<" "
176 // <<"FG4:"<<mcParticle->decayInFlight()<<std::endl;
177 // Adding vertex index
178 mcParticle->setVertexIndex0(startPoint->GetIndex());
179 mcParticle->setVertexIndex1(vertex->GetIndex());
180
181 // Set the final position
182 HepLorentzVector finalPosition(vertex->GetPosition().x()/10., vertex->GetPosition().y()/10., vertex->GetPosition().z()/10., vertex->GetTime());
183 mcParticle->finalize(finalPosition);
184
185 // push back
186 particleCol->push_back(mcParticle);
187
188 endPointFound = true;
189 }
190 }
191 }
192
193 if (!endPointFound)
194 {
195 //std::cout << "warning in finding the end point of a track" <<std::endl;
196 if (track->GetDaughterIndexes().size()==0)
197 {
198 // create the mc particle
199 Event::McParticle* mcParticle = new Event::McParticle;
200 // initialization
201 HepLorentzVector initialMomentum( track->GetP4().x()/1000., track->GetP4().y()/1000., track->GetP4().z()/1000., track->GetP4().t()/1000.);
202 HepLorentzVector initialPosition(startPoint->GetPosition().x()/10., startPoint->GetPosition().y()/10., startPoint->GetPosition().z()/10., startPoint->GetTime());
203
204 mcParticle->initialize(track->GetPDGCode(), 0, initialMomentum, initialPosition);
205
206 // Set track index
207 mcParticle->setTrackIndex( track->GetIndex() );
208
209 // Adding status flag
211 if (isPrimary ) mcParticle->addStatusFlag(Event::McParticle::PRIMARY);
212
213 // Adding vertex index
214 mcParticle->setVertexIndex0( startPoint->GetIndex() );
215 mcParticle->setVertexIndex1( -99 );
216
217 //std::cout<<"pdg:"<<track->GetPDGCode()<<" "
218 // <<"source:"<<track->GetSource()<<" "
219 // <<std::endl;
220 if (track->GetSource()=="FromGenerator")
222 else if(track->GetSource()=="FromG4")
224 else
226
227 //std::cout<<"status:"<<mcParticle->statusFlags()<<" "
228 // <<"FGN:"<<mcParticle->decayFromGenerator()<<" "
229 // <<"FG4:"<<mcParticle->decayInFlight()<<std::endl;
230
231 // push back
232 particleCol->push_back(mcParticle);
233 continue;
234 }
235 }
236 }
237
238 // Get primary McParticles
239 SmartRefVector<Event::McParticle> primaryParticleCol;
240 Event::McParticleCol::iterator iter_particle;
241 for ( iter_particle = particleCol->begin();
242 iter_particle != particleCol->end(); iter_particle++) {
243
244 if ((*iter_particle)->primaryParticle()) {
245 Event::McParticle* mcPart = (Event::McParticle*)(*iter_particle);
246 primaryParticleCol.push_back(mcPart);
247 }
248 }
249
250 if (primaryParticleCol.empty())
251 std::cout << "no primary particle found!" << std::endl;
252
253 // Add mother particle recursively
254 SmartRefVector<Event::McParticle>::iterator iter_primary;
255 for ( iter_primary = primaryParticleCol.begin();
256 iter_primary != primaryParticleCol.end(); iter_primary++) {
257 if ( !((*iter_primary)->leafParticle()) )
258 AddMother((*iter_primary), particleCol);
259 }
260
261 //register McParticle collection to TDS
262 StatusCode sc = m_evtSvc->registerObject("/Event/MC/McParticleCol", particleCol);
263 if(sc!=StatusCode::SUCCESS)
264 G4cout<< "Could not register McParticle collection" <<G4endl;
265
266 //retrive McParticle from TDS
267 /*SmartDataPtr<Event::McParticleCol> mcParticleCol(m_evtSvc,"/Event/MC/McParticleCol");
268 if(!mcParticleCol)
269 G4cout<<"Could not retrieve McParticelCol"<<G4endl;
270 Event::McParticleCol::iterator iter_mc = mcParticleCol->begin();
271 for (;iter_mc != mcParticleCol->end(); iter_mc++)
272 {
273 G4cout<<(*iter_mc)->getTrackIndex();
274 G4cout<<" "<<(*iter_mc)->particleProperty();
275 G4cout<<" "<<(*iter_mc)->getVertexIndex0();
276 G4cout<<" "<<(*iter_mc)->getVertexIndex1();
277 G4cout<<" "<<(*iter_mc)->initialFourMomentum().x();
278 G4cout<<" "<<(*iter_mc)->initialFourMomentum().y();
279 G4cout<<" "<<(*iter_mc)->initialFourMomentum().z();
280 G4cout<<" "<<(*iter_mc)->initialFourMomentum().t();
281 G4cout<<" "<<(*iter_mc)->initialPosition().x();
282 G4cout<<" "<<(*iter_mc)->initialPosition().y();
283 G4cout<<" "<<(*iter_mc)->initialPosition().z();
284 G4cout<<G4endl;
285 }
286 G4cout<<"end of retrieve McParticle from TDS"<<G4endl;
287 */
288
289}
std::vector< BesTruthVertex * > * GetVertexList()
std::vector< BesTruthTrack * > * GetTrackList()
static BesSensitiveManager * GetSensitiveManager()
G4int GetPDGCode() const
HepLorentzVector GetP4() const
BesTruthVertex * GetTerminalVertex() const
BesTruthVertex * GetVertex() const
G4String GetSource()
vector< int > GetDaughterIndexes() const
G4int GetIndex() const
BesTruthTrack * GetParentTrack() const
G4double GetTime() const
G4ThreeVector GetPosition() const
G4int GetIndex() const
void setVertexIndex0(int index0)
methods for setting and getting vertex indexes
Definition: McParticle.h:119
@ PRIMARY
Decayed in flight.
Definition: McParticle.h:37
@ ERROR
this particle is a leaf in the particle tree
Definition: McParticle.h:39
@ DECAYFLT
Decayed by generator.
Definition: McParticle.h:36
@ LEAF
primary particle
Definition: McParticle.h:38
void initialize(StdHepId id, unsigned int statusBits, const HepLorentzVector &initialMomentum, const HepLorentzVector &initialPosition, const std::string process="")
Set the initial attributes of the McParticle.
Definition: McParticle.cxx:50
void addStatusFlag(unsigned int flag)
add a new flag to the status flags
Definition: McParticle.h:103
void finalize(const HepLorentzVector &finalPosition)
Set the final attributes of the McParticle.
Definition: McParticle.cxx:85
void setVertexIndex1(int index1)
Definition: McParticle.h:123
void setTrackIndex(int trackIndex)
Definition: McParticle.h:128
ObjectList< McParticle > McParticleCol
Definition: McParticle.h:205

Referenced by SaveMcTruth().

◆ SaveMcTruth()

void BesMcTruthWriter::SaveMcTruth ( )

Definition at line 62 of file BesMcTruthWriter.cc.

63{
64 //interface to event data service
65 ISvcLocator* svcLocator = Gaudi::svcLocator();
66 StatusCode sc=svcLocator->service("EventDataSvc", m_evtSvc);
67 if (sc.isFailure())
68 G4cout<<"Could not accesss EventDataSvc!"<<G4endl;
69
70 DataObject *aMcEvent;
71 m_evtSvc->findObject("/Event/MC",aMcEvent);
72 if(aMcEvent==NULL) {
73 McEvent* aMcEvent = new McEvent;
74 sc = m_evtSvc->registerObject("/Event/MC",aMcEvent);
75 if(sc!=StatusCode::SUCCESS)
76 G4cout<< "Could not register McEvent" <<G4endl;
77 }
78
85}
Definition: McEvent.h:8

Referenced by BesTDSWriter::SaveAll().

◆ SaveMdcTruth()

void BesMcTruthWriter::SaveMdcTruth ( )

Definition at line 326 of file BesMcTruthWriter.cc.

327{
328 //Mdc McHit collection defined in BOSS
329 Event::MdcMcHitCol* aMdcMcHitCol = new Event::MdcMcHitCol;
330
331 G4int HCID = -1;
332 HCID = m_DigiMan->GetHitsCollectionID("BesMdcTruthCollection");
333 if(HCID>0)
334 {
335 BesMdcHitsCollection* HC = 0;
336 HC = (BesMdcHitsCollection*) (m_DigiMan->GetHitsCollection(HCID));
337 G4int n_hit = HC->entries();
338 if(n_hit>0)
339 {
340 //arrange hits in hits collection in order of trackIndex
341 BesMdcHit* hit;
342 vector<BesMdcHit*>* vecHC = HC->GetVector();
343 for(int i=0;i<n_hit-1;i++)
344 for(int j=i+1;j<n_hit;j++)
345 if((*vecHC)[i]->GetTrackID()>(*vecHC)[j]->GetTrackID())
346 {
347 hit = (*vecHC)[i];
348 (*vecHC)[i] = (*vecHC)[j];
349 (*vecHC)[j] = hit;
350 }
351
352 //push back MdcMcHits to MdcMcHitCol in TDS
353 for(G4int i=0;i<n_hit;i++)
354 {
355 hit = (*HC)[i];
356 const Identifier ident = MdcID::wire_id ( hit->GetLayerNo(), hit->GetCellNo() );
357 Event::MdcMcHit* mdcMcHit = new Event::MdcMcHit(ident, hit->GetTrackID(),
358 hit->GetPos().x(), hit->GetPos().y(), hit->GetPos().z(),
359 hit->GetDriftD(), hit->GetEdep(),hit->GetPosFlag() );
360 aMdcMcHitCol->push_back(mdcMcHit);
361 }
362
363 }
364 }
365
366 //register MDC McTruth collection to TDS
367 StatusCode scMdc = m_evtSvc->registerObject("/Event/MC/MdcMcHitCol", aMdcMcHitCol);
368 if(scMdc!=StatusCode::SUCCESS)
369 G4cout<< "Could not register MDC McTruth collection" <<G4endl;
370
371 //retrieve MDC McTruth from TDS
372 /*SmartDataPtr<Event::MdcMcHitCol> aMcHitCol(m_evtSvc,"/Event/MC/MdcMcHitCol");
373 if(!aMcHitCol)
374 G4cout<<"Could not retrieve MDC McTruth collection"<<G4endl;
375
376 Event::MdcMcHitCol::iterator iMcHitCol;
377 for(iMcHitCol=aMcHitCol->begin(); iMcHitCol!=aMcHitCol->end(); iMcHitCol++)
378 {
379 const Identifier ident = (*iMcHitCol)->identify();
380 G4cout<<(*iMcHitCol)->getTrackIndex();
381 G4cout<<" "<<MdcID::layer(ident);
382 G4cout<<" "<<MdcID::wire(ident);
383 G4cout<<" "<<(*iMcHitCol)->getDepositEnergy();
384 G4cout<<" "<<(*iMcHitCol)->getDriftDistance();
385 G4cout<<" "<<(*iMcHitCol)->getPositionX();
386 G4cout<<" "<<(*iMcHitCol)->getPositionY();
387 G4cout<<" "<<(*iMcHitCol)->getPositionZ();
388 G4cout<<G4endl;
389 }
390 G4cout<<"end of retrieve MDC McTruth collection"<<G4endl;
391 */
392}
G4THitsCollection< BesMdcHit > BesMdcHitsCollection
Definition: BesMdcHit.hh:78
G4double GetEdep()
Definition: BesMdcHit.hh:53
G4int GetTrackID()
Definition: BesMdcHit.hh:50
G4int GetPosFlag()
Definition: BesMdcHit.hh:60
G4int GetLayerNo()
Definition: BesMdcHit.hh:51
G4double GetDriftD()
Definition: BesMdcHit.hh:55
G4ThreeVector GetPos()
Definition: BesMdcHit.hh:54
G4int GetCellNo()
Definition: BesMdcHit.hh:52
static Identifier wire_id(int wireType, int layer, int wire)
For a single wire.
Definition: MdcID.cxx:77
ObjectVector< MdcMcHit > MdcMcHitCol
Definition: MdcMcHit.h:88

Referenced by SaveMcTruth().

◆ SaveMucTruth()

void BesMcTruthWriter::SaveMucTruth ( )

Definition at line 626 of file BesMcTruthWriter.cc.

627{
628 //Muc McHit collection defined in BOSS
629 Event::MucMcHitCol* aMucMcHitCol = new Event::MucMcHitCol;
630
631 G4int HCID = -1;
632 HCID = m_DigiMan->GetHitsCollectionID("BesMucHitsList");
633 if(HCID>0)
634 {
635 BesMucHitsCollection* HC = 0;
636 HC = (BesMucHitsCollection*) (m_DigiMan->GetHitsCollection(HCID));
637 G4int n_hit = HC->entries();
638 if(n_hit>0)
639 {
640 //arrange hits in hits collection in order of trackIndex
641 BesMucHit* hit;
642 vector<BesMucHit*>* vecHC = HC->GetVector();
643 for(int i=0;i<n_hit-1;i++)
644 for(int j=i+1;j<n_hit;j++)
645 if((*vecHC)[i]->GetTrackIndex()>(*vecHC)[j]->GetTrackIndex())
646 {
647 hit = (*vecHC)[i];
648 (*vecHC)[i] = (*vecHC)[j];
649 (*vecHC)[j] = hit;
650 }
651
652 for(G4int i=0;i<n_hit;i++)
653 {
654 hit = (*HC)[i];
655 Identifier ident = MucID::channel_id ( hit->GetPart(), hit->GetSeg(), hit->GetGap(),
656 hit->GetStrip() );
657 Event::MucMcHit* mucMcHit = new Event::MucMcHit(ident, hit->GetTrackIndex(),
658 hit->GetPos().x(), hit->GetPos().y(), hit->GetPos().z(), hit->GetMomentum().x(),
659 hit->GetMomentum().y(), hit->GetMomentum().z() );
660 aMucMcHitCol->push_back(mucMcHit);
661 }
662
663 }
664 }
665
666 //register MUC McTruth collection to TDS
667 StatusCode scMuc = m_evtSvc->registerObject("/Event/MC/MucMcHitCol", aMucMcHitCol);
668 if(scMuc!=StatusCode::SUCCESS)
669 G4cout<< "Could not register MUC McTruth collection" <<G4endl;
670
671 //retrieve MUC McTruth from TDS
672 /*SmartDataPtr<Event::MucMcHitCol> aMcHitCol(m_evtSvc,"/Event/MC/MucMcHitCol");
673 if(!aMcHitCol)
674 G4cout<<"Could not retrieve MUC McTruth collection"<<G4endl;
675
676 Event::MucMcHitCol::iterator iMcHitCol;
677 for(iMcHitCol=aMcHitCol->begin(); iMcHitCol!=aMcHitCol->end(); iMcHitCol++)
678 {
679 const Identifier ident = (*iMcHitCol)->identify();
680 G4cout<<(*iMcHitCol)->getTrackIndex();
681 G4cout<<" "<<MucID::part(ident);
682 G4cout<<" "<<MucID::seg(ident);
683 G4cout<<" "<<MucID::gap(ident);
684 G4cout<<" "<<MucID::strip(ident);
685 G4cout<<" "<<(*iMcHitCol)->getPositionX();
686 G4cout<<" "<<(*iMcHitCol)->getPositionY();
687 G4cout<<" "<<(*iMcHitCol)->getPositionZ();
688 G4cout<<" "<<(*iMcHitCol)->getPx();
689 G4cout<<" "<<(*iMcHitCol)->getPy();
690 G4cout<<" "<<(*iMcHitCol)->getPz();
691 G4cout<<G4endl;
692 }
693 G4cout<<"end of retrieve MUC McTruth"<<G4endl;
694 */
695}
G4THitsCollection< BesMucHit > BesMucHitsCollection
Definition: BesMucHit.hh:104
G4int GetTrackIndex()
Definition: BesMucHit.hh:63
G4ThreeVector GetPos()
Definition: BesMucHit.hh:68
G4int GetSeg()
Definition: BesMucHit.hh:75
G4int GetStrip()
Definition: BesMucHit.hh:77
G4int GetPart()
Definition: BesMucHit.hh:74
G4ThreeVector GetMomentum()
Definition: BesMucHit.hh:71
G4int GetGap()
Definition: BesMucHit.hh:76
static Identifier channel_id(int barrel_ec, int segment, int layer, int channel)
For a single crystal.
Definition: MucID.cxx:135
ObjectVector< MucMcHit > MucMcHitCol
Definition: MucMcHit.h:89

Referenced by SaveMcTruth().

◆ SaveTofTruth()

void BesMcTruthWriter::SaveTofTruth ( )

Definition at line 394 of file BesMcTruthWriter.cc.

395{
396 //Tof McHit collection defined in BOSS
397 Event::TofMcHitCol* aTofMcHitCol = new Event::TofMcHitCol;
398
399 G4int HCID = -1;
400 HCID = m_DigiMan->GetHitsCollectionID("BesTofHitsList");
401 if(HCID>0)
402 {
403 BesTofHitsCollection* HC = 0;
404 HC = (BesTofHitsCollection*) (m_DigiMan->GetHitsCollection(HCID));
405 G4int n_hit = HC->entries();
406 if(n_hit>0)
407 {
408 //arrange hits in hits collection in order of trackIndex
409 BesTofHit* hit;
410 vector<BesTofHit*>* vecHC = HC->GetVector();
411 for(int i=0;i<n_hit-1;i++)
412 for(int j=i+1;j<n_hit;j++)
413 if((*vecHC)[i]->GetTrackIndex()>(*vecHC)[j]->GetTrackIndex())
414 {
415 hit = (*vecHC)[i];
416 (*vecHC)[i] = (*vecHC)[j];
417 (*vecHC)[j] = hit;
418 }
419
420 //push back TofMcHit collection to TDS
421 for(G4int i=0;i<n_hit;i++)
422 {
423 hit = (*HC)[i];
424 unsigned int scinNum, barrel_ec, layer = 0;
425 scinNum = hit->GetScinNb();
426 barrel_ec = hit->GetPartId();
427
428 //std::cout << "BesMcTruthWriter PartID | scinNum " << barrel_ec << " | " << scinNum << std::endl;
429 //std::cout << "BesMcTruthWriter is_barrel(barrel_ec) " << TofID::is_barrel(barrel_ec) << std::endl;
430
431
432 if (TofID::is_barrel(barrel_ec) && scinNum > TofID::getPHI_BARREL_MAX()) {
433 layer = 1;
434 scinNum = scinNum - TofID::getPHI_BARREL_MAX() - 1;
435 }
436
437
438
439 Identifier ident;
440
441 if(barrel_ec==0 || barrel_ec==1 || barrel_ec==2 ) //Old-Tof
442 {
443 ident = TofID::cell_id ( barrel_ec, layer, scinNum, 0);
444 }
445 else //MRPC Part
446 {
447 int help = BesTofDigitizerEcV4::Calculate_Readoutstrip_number_continuum(hit->GetPos().x()*mm,hit->GetPos().y()*mm, barrel_ec, scinNum);
449 //std::cout << "BesMcTruthWriter Readoutstrip " << help << std::endl;
450 //std::cout << "BesMcTruthWriter unique Identifier " << scinNum << std::endl;
451 ident = TofID::cell_id_mrpc_mc ( barrel_ec, scinNum);
452 }
453
454
455
456 Event::TofMcHit* tofMcHit = new Event::TofMcHit( ident, hit->GetTrackIndex(),
457 hit->GetPos().x(), hit->GetPos().y(), hit->GetPos().z(), hit->GetMomentum().x(),
458 hit->GetMomentum().y(), hit->GetMomentum().z(), hit->GetTrackL(), hit->GetTime() );
459
460 aTofMcHitCol->push_back(tofMcHit);
461 }
462 }
463 }
464
465 //register TOF McTruth collection to TDS
466 StatusCode scTof = m_evtSvc->registerObject("/Event/MC/TofMcHitCol", aTofMcHitCol);
467 if(scTof!=StatusCode::SUCCESS)
468 G4cout<< "Could not register TOF McTruth collection" <<G4endl;
469
470 //retrieve TOF McTruth from TDS
471 /*SmartDataPtr<Event::TofMcHitCol> aMcHitCol(m_evtSvc,"/Event/MC/TofMcHitCol");
472 if(!aMcHitCol)
473 G4cout<<"Could not retrieve TOF McTruth collection"<<G4endl;
474
475 Event::TofMcHitCol::iterator iMcHitCol;
476 for(iMcHitCol=aMcHitCol->begin(); iMcHitCol!=aMcHitCol->end(); iMcHitCol++)
477 {
478 const Identifier ident = (*iMcHitCol)->identify();
479 G4cout<<(*iMcHitCol)->getTrackIndex();
480 G4cout<<" "<<TofID::barrel_ec(ident);;
481 G4cout<<" "<<TofID::layer(ident);
482 G4cout<<" "<<TofID::phi_module(ident);
483 G4cout<<" "<<(*iMcHitCol)->getPositionX();
484 G4cout<<" "<<(*iMcHitCol)->getPositionY();
485 G4cout<<" "<<(*iMcHitCol)->getPositionZ();
486 G4cout<<" "<<(*iMcHitCol)->getPx();
487 G4cout<<" "<<(*iMcHitCol)->getPy();
488 G4cout<<" "<<(*iMcHitCol)->getPz();
489 G4cout<<" "<<(*iMcHitCol)->getTrackLength();
490 G4cout<<" "<<(*iMcHitCol)->getFlightTime();
491 G4cout<<G4endl;
492 }
493 G4cout<<"end of retrieve TOF McTruth"<<G4endl;
494 */
495}
G4THitsCollection< BesTofHit > BesTofHitsCollection
Definition: BesTofHit.hh:108
std::string help()
Definition: HelloServer.cpp:35
static G4int Produce_unique_identifier(G4int module_mrpc_f, G4int readoutstripnumber_f)
static G4int Calculate_Readoutstrip_number_continuum(G4double x_mm, G4double y_mm, G4int partId_f, G4int module_mrpc_f)
G4double GetTime()
Definition: BesTofHit.hh:66
G4ThreeVector GetPos()
Definition: BesTofHit.hh:65
G4int GetScinNb()
Definition: BesTofHit.hh:61
G4int GetPartId()
Definition: BesTofHit.hh:60
G4double GetTrackL()
Definition: BesTofHit.hh:64
G4ThreeVector GetMomentum()
Definition: BesTofHit.hh:69
G4int GetTrackIndex()
Definition: BesTofHit.hh:58
static Identifier cell_id(int barrel_ec, int layer, int phi_module, int end)
For a single crystal.
Definition: TofID.cxx:156
static value_type getPHI_BARREL_MAX()
Definition: TofID.cxx:237
static bool is_barrel(const Identifier &id)
Test for barrel.
Definition: TofID.cxx:78
static Identifier cell_id_mrpc_mc(int partID, int scinNum)
Definition: TofID.cxx:189
ObjectVector< TofMcHit > TofMcHitCol
Definition: TofMcHit.h:100

Referenced by SaveMcTruth().


The documentation for this class was generated from the following files: