CGEM BOSS 6.6.5.f
BESIII Offline Software System
Loading...
Searching...
No Matches
BesCgemSD.cc
Go to the documentation of this file.
1//-------------------------------------------------------------------------------------//
2// BOOST --- BESIII Object_Oriented Simulation Tool //
3//-------------------------------------------------------------------------------------//
4/*
5 * =====================================================================================
6 *
7 * Filename: BesCgemSD.cc
8 * Description:
9 * Conventions:
10 * gv_ : global variable
11 * lv_ : local variable used in function
12 * lvd_ : local variable double
13 * lv3_ : local variable ThreeVector
14 * m_* : normal member of class
15 * sm_* : static member of class
16 * m_M_* : class data member, Material of each layer
17 * m_N_* : class data member, Number of layers (CgemLayer,GemFoil)
18 * m_R_* : class data member, Radius of each (material) layer, and so on
19 * m_L_* : class data member, Length of each layer or hole pitch
20 * m_T_* : class data member, Thickness of each (material) layer
21 * m_W* : calss data member, Width of strip
22 * m_A_* : class data member, Angle of anode VStrip
23 * Version: CgemSim-01-00-00
24 * Created: 01/02/2014 08:59:21 AM
25 * Revision: CgemSim-00-00-01(in CMT version CgemBoss-0.0.1 written by xiuql)
26 * Compiler: gcc
27 * Author: [email protected]
28 * Organization: DG1,EPC,IHEP
29 * History:
30 * <Num> <Author> <Time> <Version> <remark>
31 * 0 juxd 20140102 01-00-00 created,
32 *
33 * =====================================================================================
34 */
35
36//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
37/* Header file: BOSS */
38#include "BesCgemSD.hh"
39#include "ReadBoostRoot.hh"
40
41//#include "GaudiKernel/SmartIF.h"
42//#include "GaudiKernel/Bootstrap.h"
43//#include "GaudiKernel/IProperty.h"
44//#include "GaudiKernel/ServiceLocatorHelper.h"
45//#include "GaudiKernel/ISvcLocator.h"
46
47/* Header file: Geant4 */
48#include "G4RunManager.hh"
49#include "G4SDManager.hh"
50#include "G4HCofThisEvent.hh"
51#include "G4Step.hh"
52#include "G4ThreeVector.hh"
53#include "G4UnitsTable.hh"
54#include "G4ios.hh"
55#include "Randomize.hh"
56#include "globals.hh"
57
58/* Header file: ROOT */
59
60/* Header file: C++ */
61#include <iostream>
62#include <iomanip>
63#include <fstream>
64#include <vector>
65#include <string>
66#include <cfloat>
67
68// Gaudi
69#include "GaudiKernel/IMessageSvc.h"
70#include "GaudiKernel/MsgStream.h"
71#include "GaudiKernel/ISvcLocator.h"
72#include "GaudiKernel/Bootstrap.h"
73#include "GaudiKernel/SmartIF.h"
74#include "GaudiKernel/Property.h"
75#include "GaudiKernel/IJobOptionsSvc.h"
76
77using namespace std;
78
79//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
80 BesCgemSD::BesCgemSD(G4String name)
81:BesSensitiveDetector(name)
82{
83 m_F_printTruth = 0;
84 m_F_outputTruth = 0;
85 myCgemTruth=NULL;
86
87 /* Define hit collection name handled by SD */
88 collectionName.insert("BesCgemHitsCollection");
89 collectionName.insert("BesCgemTruthCollection");
90
91 // get info from JobOptionsSvc.
92 ISvcLocator* iface = Gaudi::svcLocator();
93 SmartIF<IJobOptionsSvc> jos(iface->service("JobOptionsSvc"));
94 const std::vector<const Property*>* ppp = jos->getProperties("BesCgemSD");
95
96 string outputFile;
97 //outputFile += "/dat/outputMcTruth.txt";
98
99 if (ppp) {
100 //std::cout << "BesCgemSD size: " << ppp->size() << std::endl;
101 for (std::vector<const Property*>::const_iterator it = ppp->begin();
102 it != ppp->end(); ++it) {
103 if ( (*it)->name() == "outputFile") {
104 outputFile = (*it)->toString();
105 }
106 }
107 }
108 IMessageSvc* msgSvc;
109 Gaudi::svcLocator() -> service("MessageSvc", msgSvc);
110 MsgStream log(msgSvc, "BesCgemSD::BesCgemSD():BesSensitiveDetector()");
111
112 /* Save McTruth information to compare with Root RawData Information */
113 if (m_F_outputTruth == 1)
114 {
115 log<< MSG::INFO << "CgemSim::BesCgemSD::Initialize(), McTruth Hit will be output to a txt file!"
116 << endreq
117 << "outputFile Path: " << outputFile << endreq;
118
119 m_outputTruth.open(outputFile.c_str(), ios::out);
120
121 if (!m_outputTruth.good())
122 {
123 log<< MSG::ERROR << "CgemSim::BesCgemSD::Initialize(), Fail to open McTruth output file: "
124 << outputFile << endreq;
125 return ;
126 }
127 } /* End of 'if (m_outputTruth == 1)' */
128}
129
130//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
131BesCgemSD::~BesCgemSD()
132{
133 m_outputTruth.close();
134}
135
136//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
137/* Invoked at the beginning of each event */
138void BesCgemSD::Initialize(G4HCofThisEvent* HCE)
139{
140 /* Instantiate hit collection */
141 m_collection_hit = new BesCgemHitsCollection(SensitiveDetectorName,collectionName[0]);
142
143 //Add by sunxh (lvi_ID_hit)
144 lvi_ID_hit = 0;
145 //end add
146 //IMessageSvc* msgSvc;
147 //Gaudi::svcLocator() -> service("MessageSvc", msgSvc);
148 //MsgStream log(msgSvc, "BesCgemSD::Initialize()");
149 //log<< MSG::INFO << "|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||" << endreq;
150 //log<< MSG::INFO << "CgemSim::BesCgemSD::Initialize(), Begin !!!" << endreq;
151}
152
153//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
154/* MC Truth */
155void BesCgemSD::BeginOfTruthEvent(const G4Event* f_event)
156{
157 m_collection_truth = new BesCgemHitsCollection(SensitiveDetectorName,collectionName[1]);
158 lvi_ID_hit=0;
159
160}
161
162//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
163/* MC Truth */
164void BesCgemSD::EndOfTruthEvent(const G4Event* f_event)
165{
166 static G4int HLID=-1;
167 if (HLID < 0)
168 {
169 HLID = G4SDManager::GetSDMpointer()->GetCollectionID(collectionName[1]);
170 }
171
172 G4HCofThisEvent *HCE = f_event->GetHCofThisEvent();
173
174 IMessageSvc* msgSvc;
175 Gaudi::svcLocator() -> service("MessageSvc", msgSvc);
176 MsgStream log(msgSvc, "BesCgemSD::EndOfTruthEvent()");
177 //log<< MSG::INFO << "CgemSim::BesCgemSD::EndOfTruthEvent(), Entries of truth collection is "
178 // << m_collection_truth->entries() << endreq;
179
180 HCE->AddHitsCollection(HLID, m_collection_truth);
181
182 /* Print and output McTruth information */
183 if (m_F_printTruth==1 or m_F_outputTruth==1)
184 {
185 /* Print McTruth information */
186 if (m_F_printTruth == 1)
187 {
188 log<< MSG::INFO << "-------------------------------------------------------------------------------"
189 << endreq;
190 log<< MSG::INFO << "CgemSim::BesCgemSD::EndOfTruthEvent(), McTruth information: " << endreq;
191 log<< MSG::INFO << left << "Track " << left << "Layer " << left << "PDG "
192 << left << setw(12) << "R_Pre "
193 << left << setw(12) << "R_Post "
194 << left << setw(40) << "pos_Pre "
195 << left << setw(40) << "pos_Post "
196 << left << setw(12) << "P_Pre "
197 << left << setw(12) << "P_Post "
198 << left << setw(16) << "StepL "
199 << left << setw(16) <<"Time "
200 << left << "Edep "
201 << endreq;
202 }
203
204
205 for (G4int ii=0; ii<m_collection_truth->entries(); ii++)
206 {
207 BesCgemHit *lv_truth = (*m_collection_truth)[ii];
208
209 /* Print McTruth information */
210 if (m_F_printTruth == 1)
211 {
212 G4cout << left << setw( 6) << lv_truth-> GetTrackID ()
213 << left << setw( 6) << lv_truth-> GetLayerID ()
214 << left << setw( 6) << lv_truth-> GetPDGCode ()
215 << left << setw(12) << lv_truth-> GetPositionOfPrePoint ().rho()
216 << left << setw(12) << lv_truth-> GetPositionOfPostPoint ().rho()
217 << left << lv_truth-> GetPositionOfPrePoint()
218 << left << lv_truth-> GetPositionOfPostPoint()
219 << left << setw(12) << lv_truth-> GetMomentumOfPrePoint().mag()
220 << left << setw(12) << lv_truth-> GetMomentumOfPostPoint().mag()
221 << left << setw(16) << lv_truth-> GetStepLength ()
222 << left << setw(16) << lv_truth-> GetGlobalTime ()
223 << left << setw(16) << lv_truth-> GetTotalEnergyDeposit ()
224 << G4endl;
225 }
226
227 /* Save McTruth information to compare with Root RawData Information */
228 if (m_F_outputTruth == 1)
229 {
230 m_outputTruth << left << setw( 5) << lv_truth->GetTrackID()
231 << left << setw( 3) << lv_truth->GetLayerID()
232 << left << setw(12) << lv_truth->GetTotalEnergyDeposit()
233 << left << setw(12) << lv_truth->GetPositionOfPrePoint().x()
234 << left << setw(12) << lv_truth->GetPositionOfPrePoint().y()
235 << left << setw(12) << lv_truth->GetPositionOfPrePoint().z()
236 << left << setw(12) << lv_truth->GetPositionOfPostPoint().x()
237 << left << setw(12) << lv_truth->GetPositionOfPostPoint().y()
238 << left << setw(12) << lv_truth->GetPositionOfPostPoint().z()
239 << left << setw(12) << lv_truth->GetMomentumOfPrePoint().x()
240 << left << setw(12) << lv_truth->GetMomentumOfPrePoint().y()
241 << left << setw(12) << lv_truth->GetMomentumOfPrePoint().z()
242 << left << setw(12) << lv_truth->GetMomentumOfPostPoint().x()
243 << left << setw(12) << lv_truth->GetMomentumOfPostPoint().y()
244 << left << setw(12) << lv_truth->GetMomentumOfPostPoint().z()
245 << G4endl;
246 }
247 } /* End of 'for (G4int ii=0; ii<m_collection_truth->entries(); ii++)' */
248
249 //m_outputTruth.close();
250
251 } /* End of 'if (m_F_printTruth==1 or m_F_outputTruth==1)' */
252}
253//************************************************************************
254//************************************************************************
255//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
256/* Invoked at the beginning of every step in the SD volume */
257G4bool BesCgemSD::ProcessHits(G4Step* f_step, G4TouchableHistory* f_touchable)
258{
259 bool printFlag=false;
260
261 IMessageSvc* msgSvc;
262 Gaudi::svcLocator() -> service("MessageSvc", msgSvc);
263 MsgStream log(msgSvc, "BesCgemSD::ProcessHits()");
264
265 /* Event -> Track -> Step */
266 G4Track *lv_track = f_step->GetTrack();
267
268 /* Check tracks */
269 /* Skip neutral tracks */
270 G4double charge = lv_track->GetDefinition()->GetPDGCharge();
271 if (charge == 0.)
272 {
273 if(printFlag)
274 log<< MSG::WARNING << "CgemSim::BesCgemSD::ProcessHits(), Track is neutral, skip the track" << endreq;
275 return false;
276 }
277 G4int lvi_pdg_code = lv_track->GetDefinition()->GetPDGEncoding();
278
279 /* Global Time (time since the current event began ) */
280 G4double lvd_global_time = lv_track->GetGlobalTime();
281 if (isnan(lvd_global_time))
282 {
283 if(printFlag)
284 log<< MSG::WARNING << "CgemSim::BesCgemSD::ProcessHits(), Gloabal time is nan" << endreq;
285 return false;
286 }
287 /* MDC T window is 2 microsecond */
288 if (lvd_global_time > 2000)
289 {
290 if(printFlag)
291 log<< MSG::WARNING << "CgemSim::BesCgemSD::ProcessHits(), Time window is >2000" << endreq;
292 return false;
293 }
294
295 /* Skip no energy deposit tracks */
296 /* if (lvd_E_deposit == 0.)
297 {
298 log<< MSG::WARNING << "CgemSim::BesCgemSD::ProcessHits(), Energy deposition is 0,lose track="
299 << lvi_ID_track << endreq;
300 return false;
301 }*/
302
303 /* Skip no length tracks */
304 G4double lvd_L_step = f_step->GetStepLength();
305 if (lvd_L_step == 0.)
306 {
307 if(printFlag)
308 log<< MSG::WARNING << "CgemSim::BesCgemSD::ProcessHits(), Step length is 0, skip track" << endreq;
309 return false;
310 }
311
312 /* A G4Step object consists of two points */
313 G4StepPoint *lv_step_pre = f_step->GetPreStepPoint();
314 G4StepPoint *lv_step_post = f_step->GetPostStepPoint();
315 G4TrackStatus currentTrackStatus = lv_track->GetTrackStatus();
316
317 /* Skip delta electron, added by sunxh */
318 G4int lvi_RdtElectron = 0;
319 G4String process_name="Generator";
320 if (NULL != lv_track->GetCreatorProcess()){
321 process_name = lv_track->GetCreatorProcess()->GetProcessName();
322 if (process_name == "hIoni" || process_name == "ionIoni" || process_name == "conv" || process_name == "compt" || process_name == "muIoni" || process_name == "eIoni" || process_name == "phot" || process_name == "muMinusCaptureAtRest" ) lvi_RdtElectron = 1;
323 //if(lvi_RdtElectron != 1&&lv_track->GetDefinition()->GetPDGEncoding()!=22){
324 // G4cout << "process_name = " << process_name << G4endl;
325 // G4cout << "PDG = " << lv_track->GetDefinition()->GetPDGEncoding()<<"\tTrack ID = "<<lv_track->GetTrackID() <<"\tParent ID = "<<lv_track->GetParentID() << G4endl;
326 //}
327 //cout<<"particle pdg="<<lvi_pdg_code<<" with CreatorProcess "<<process_name<<endl;
328 }
329 if(lvi_RdtElectron) return false;
330
331 //else {
332 // G4cout << "process_name = null" << G4endl;
333 // G4cout << "PDG = " << lv_track->GetDefinition()->GetPDGEncoding()<<"\tTrack ID = "<<lv_track->GetTrackID() <<"\tParent ID = "<<lv_track->GetParentID() << G4endl;
334 //}
335 // G4cout << "PDG = " << lv_track->GetDefinition()->GetPDGEncoding() << G4endl;
336 // else G4cout << "No process_name PDG = " << lv_track->GetDefinition()->GetPDGEncoding() << G4endl;
337 // if (abs(lv_track->GetDefinition()->GetPDGEncoding()) == 11 ) lvi_RdtElectron = 1;
338 // else{
339 // if (abs(lv_track->GetDefinition()->GetPDGEncoding()) == 11) {
340 // G4cout << "No process_name PDG = " << lv_track->GetDefinition()->GetPDGEncoding() << G4endl;
341 // }
342 // }
343 //End Add
344
345 /* Local variables used to get Hit information */
346 /* Get track ID */
347 G4int lvi_ID_track = 0;
348 G4int lvi_ID_truth = 0;
349 G4int lvi_ID_G4Track = 0;
350
351 /* In the following we use 'SetTrackID(lvi_ID_track)' */
352 //lvi_ID_track = lv_track->GetTrackID();
353 //GetCurrentTrackIndex(lvi_ID_truth, lvi_ID_G4Track);
354
355 //Add by sunxh
356 lvi_ID_truth = lv_track->GetTrackID();
357 GetCurrentTrackIndex(lvi_ID_track, lvi_ID_G4Track);
358 // G4cout<<"PDG = "<<lv_track->GetDefinition()->GetPDGEncoding()<<"\tprocess_name = "<<process_name<<"\tlvi_ID_truth = " << lvi_ID_truth<<"\tlvi_ID_track = "<<lvi_ID_track<<"\tlvi_ID_G4Track = "<<lvi_ID_G4Track<<G4endl;
359 //end
360
361 G4int lvi_ID_layer = 0;
362 //G4int lvi_pdg_code = lv_track->GetDefinition()->GetPDGEncoding();
363 G4int lvi_ID_parent = lv_track->GetParentID();
364 G4double lvd_E_deposit = f_step->GetTotalEnergyDeposit();
365 /* Positions in the global coordinate system */
366 G4ThreeVector lv3_XYZ_pre = lv_step_pre ->GetPosition();
367 G4ThreeVector lv3_XYZ_post = lv_step_post->GetPosition();
368
369 G4ThreeVector lv3_P_pre = lv_step_pre ->GetMomentum();
370 G4ThreeVector lv3_P_post = lv_step_post->GetMomentum();
371
372 if(printFlag) {
373 log<< MSG::INFO << "..............................................................................." << endreq;
374 log<< MSG::INFO << "CgemSim::BesCgemSD::ProcessHits(), Process a hit at track="
375 << lvi_ID_track << endreq;
376 }
377
378
379
380 /* Check and Save Hits and MCTruth Hits */
381 /* Get pre step point status and its volume */
382 /* Volume depth: 0-Gap_D, 1-CgemLayer, 2-Cgem, 3-BES */
383 G4StepStatus lv_status_pre = lv_step_pre ->GetStepStatus();
384 G4StepStatus lv_status_post = lv_step_post->GetStepStatus();
385
386 /* Retrieve a touchable by creating a handle for it */
387 const G4VTouchable *lv_TOUCHABLE = lv_step_pre ->GetTouchable();
388 /* Get the current volume(where the step has just gone through) and its logical name */
389 G4String lvs_volume_name = lv_TOUCHABLE ->GetVolume(0)->GetLogicalVolume()->GetName();
390
391 /* Save Hits in SD Geometry only! */
392 //if ((lv_status_pre == fGeomBoundary) && (lvs_volume_name == "Gap_D_Logic"))
393
394
395 //if (lvs_volume_name == "Gap_D_logic")
396 if (lvs_volume_name == "Gap_D_logic0" or lvs_volume_name == "Gap_D_logic1" or lvs_volume_name == "Gap_D_logic2")
397 {
398
399 lvi_ID_layer = lv_TOUCHABLE->GetVolume(1)->GetCopyNo();
400
401 if(printFlag) {
402
403 log<< MSG::INFO << "CgemSim::BesCgemSD::ProcessHits(), Track enters the SD volume: "
404 << lvs_volume_name << endreq;
405
406 /* CgemLayer copy number */
407
408 // KeyType lv_key(lvi_ID_track, lvi_ID_layer);
409
410 /* Create Hits */
411 log<< MSG::INFO << "CgemSim::BesCgemSD::ProcessHits(), Create a hit at layer="
412 << lvi_ID_layer << " track=" << lvi_ID_track << endreq;
413 }
414
415 BesCgemHit *l_new_hit = new BesCgemHit();
416
417 //Add by sunxh
418 l_new_hit->SetRdtEl ( lvi_RdtElectron );
419 l_new_hit->SetHitID ( lvi_ID_hit );
420 //end
421 l_new_hit->SetTrackID ( lvi_ID_track );
422 l_new_hit->SetLayerID ( lvi_ID_layer );
423 l_new_hit->SetPDGCode ( lvi_pdg_code );
424 l_new_hit->SetParentID ( lvi_ID_parent );
425 l_new_hit->SetGlobalTime ( lvd_global_time );
426 l_new_hit->SetTotalEnergyDeposit ( lvd_E_deposit );
427 l_new_hit->SetStepLength ( lvd_L_step );
428 l_new_hit->SetPositionOfPrePoint ( lv3_XYZ_pre );
429 l_new_hit->SetPositionOfPostPoint ( lv3_XYZ_post );
430 l_new_hit->SetMomentumOfPrePoint ( lv3_P_pre );
431 l_new_hit->SetMomentumOfPostPoint ( lv3_P_post );
432 l_new_hit->SetCreatorProcess ( process_name );
433
434 m_collection_hit->insert(l_new_hit);
435
436 //Add by sunxh
437 hit_ID_vector.push_back(lvi_ID_hit);
438 lvi_ID_hit++;
439 //end
440
441
442
443 /* Create Truth Hit */
444 /* MCTruth Hit only contains the point input and output the volume (not step by step) */
445 /* Meaning of 'MCTruth hit' parameters is not the same as 'hit' */
446 /* All steps in SD volume will be merged to 1 step */
447 /* Check over that the particle has just entered in the current volume */
448 /* It is at the first step in the volume; the preStepPoint is at the boundary */
449
450 KeyType lv_key(lvi_ID_track, lvi_ID_layer);
451
452 if (m_collection_truth)
453 {
454 int newVer=1;
455
456 if(newVer==0) {
457 m_map_truth_it = m_map_truth.find(lv_key);
458 if (m_map_truth_it == m_map_truth.end())
459 {
460 log<< MSG::INFO << "CgemSim::BesCgemSD::ProcessHits(), Create a Truth at layer="
461 << lvi_ID_layer << " track=" << lvi_ID_track << endreq;
462
463 BesCgemHit *l_new_truth = new BesCgemHit();
464
465 l_new_truth->SetTrackID ( lvi_ID_track ); /* SAME */
466 l_new_truth->SetLayerID ( lvi_ID_layer ); /* SAME */
467 l_new_truth->SetPDGCode ( lvi_pdg_code ); /* SAME */
468 l_new_truth->SetParentID ( lvi_ID_parent ); /* SAME */
469 l_new_truth->SetGlobalTime ( lvd_global_time ); /* DIFF */
470 l_new_truth->SetTotalEnergyDeposit ( lvd_E_deposit ); /* SAME */
471 l_new_truth->SetStepLength ( lvd_L_step ); /* SAME */
472 l_new_truth->SetPositionOfPrePoint ( lv3_XYZ_pre ); /* SAME */
473 l_new_truth->SetPositionOfPostPoint ( lv3_XYZ_post ); /* SAME */
474 l_new_truth->SetMomentumOfPrePoint ( lv3_P_pre ); /* SAME */
475 l_new_truth->SetMomentumOfPostPoint ( lv3_P_post ); /* SAME */
476
477 //Add by sunxh
478 int tmp[2000];
479 for(int ii =0; ii < hit_ID_vector.size();ii++){
480 tmp[ii] = hit_ID_vector[ii];
481 }
482 l_new_truth->AddIdentifier(tmp,hit_ID_vector.size());
483 hit_ID_vector.clear();
484 //end
485
486
487 m_collection_truth->insert(l_new_truth);
488
489 m_map_truth[lv_key] = m_collection_truth->entries() - 1;
490
491 if (lv_step_pre->GetStepStatus() == fGeomBoundary)
492 {
493 log<< MSG::INFO << "CgemSim::BesCgemSD::ProcessHits(), Truth starts from boundry!!!"
494 << endreq;
495 }
496 else
497 {
498 log<< MSG::INFO << "CgemSim::BesCgemSD::ProcessHits(), Truth doesn't start from boundry!!!"
499 << endreq;
500 }
501 }//if (m_map_truth_it == m_map_truth.end())
502 else
503 {
504 if(printFlag)
505 log<< MSG::INFO << "CgemSim::BesCgemSD::ProcessHits(), Add a step to Truth at layer="
506 << lvi_ID_layer << " track=" << lvi_ID_track << endreq;
507
508 G4int lvi_ID_truth = (*m_map_truth_it).second;
509 G4int lvi_ID_track_pre = (*m_collection_truth)[lvi_ID_truth]->GetTrackID();
510 G4double lvd_E_pre = (*m_collection_truth)[lvi_ID_truth]->GetTotalEnergyDeposit();
511 G4double lvd_L_step_pre = (*m_collection_truth)[lvi_ID_truth]->GetStepLength();
512
513 (*m_collection_truth)[lvi_ID_truth]->SetTotalEnergyDeposit(lvd_E_deposit + lvd_E_pre);
514 (*m_collection_truth)[lvi_ID_truth]->SetStepLength(lvd_L_step + lvd_L_step_pre);
515 (*m_collection_truth)[lvi_ID_truth]->SetPositionOfPostPoint(lv3_XYZ_post);
516 (*m_collection_truth)[lvi_ID_truth]->SetMomentumOfPostPoint(lv3_P_post);
517 }
518
519 /* Check over that the particle is leaving the current volume */
520 /* It is at the last step in the volume; the postStepPoint is at the boundary */
521 if(printFlag)
522 if (lv_step_post->GetStepStatus() == fGeomBoundary)
523 {
524 log<< MSG::INFO << "CgemSim::BesCgemSD::ProcessHits(), Finish a Truth at layer="
525 << lvi_ID_layer << " track=" << lvi_ID_track << endreq;
526 }
527 }// old vjikersion
528 else
529 {
530 if(myCgemTruth==NULL)
531 {
532 myCgemTruth= new BesCgemHit();
533 if(printFlag)
534 cout<<"creat a new myCgemTruth, prePos is Boundary: "<<(lv_status_pre==fGeomBoundary)<<" trkID="<<lvi_ID_track<<", LayerID="<<lvi_ID_layer<<endl;
535 // myCgemTruth->SetTrackID ( lvi_ID_truth); /* SAME */
536 myCgemTruth->SetTrackID ( lvi_ID_track); /* SAME */
537 myCgemTruth->SetLayerID ( lvi_ID_layer ); /* SAME */
538 myCgemTruth->SetPDGCode ( lvi_pdg_code ); /* SAME */
539 myCgemTruth->SetParentID ( lvi_ID_parent ); /* SAME */
540 myCgemTruth->SetGlobalTime ( lvd_global_time ); /* DIFF */
541 myCgemTruth->SetTotalEnergyDeposit ( lvd_E_deposit ); /* SAME */
542 myCgemTruth->SetStepLength ( lvd_L_step ); /* SAME */
543 myCgemTruth->SetPositionOfPrePoint ( lv3_XYZ_pre ); /* SAME */
544 myCgemTruth->SetPositionOfPostPoint ( lv3_XYZ_post ); /* SAME */
545 myCgemTruth->SetMomentumOfPrePoint ( lv3_P_pre ); /* SAME */
546 myCgemTruth->SetMomentumOfPostPoint ( lv3_P_post ); /* SAME */
547 myCgemTruth->SetCreatorProcess ( process_name );
548 double dphi=lv3_XYZ_pre.phi()-lv3_XYZ_post.phi();
549 double dX=dphi*87.5;
550 //if(lvi_ID_layer==0) cout<<"dX="<<dX<<" , nXStrips="<<dX/0.65<<endl;
551 }
552 else
553 {
554 G4double lvd_E_pre = myCgemTruth->GetTotalEnergyDeposit();
555 G4double lvd_L_step_pre = myCgemTruth->GetStepLength();
556 G4ThreeVector pos_pre_last = myCgemTruth->GetPositionOfPrePoint();
557 G4ThreeVector pos_post_last = myCgemTruth->GetPositionOfPostPoint();
558 //cout<<"add CGEM truth, last hit: "<<pos_pre_last<<" to "<<pos_post_last<<endl;
559 //cout<<" new hit: "<<lv3_XYZ_pre <<" to "<<lv3_XYZ_post<<endl;
560 myCgemTruth->SetTotalEnergyDeposit(lvd_E_deposit + lvd_E_pre);
561 myCgemTruth->SetStepLength(lvd_L_step + lvd_L_step_pre);
562 myCgemTruth->SetPositionOfPostPoint(lv3_XYZ_post);
563 myCgemTruth->SetMomentumOfPostPoint(lv3_P_post);
564 if(printFlag)
565 cout<<"update myCgemTruth: prePos is Boundary: "<<(lv_status_pre==fGeomBoundary)<<", postPos is Boundary:"<<(lv_status_post==fGeomBoundary)<<", trkID="<<lvi_ID_track<<", LayerID="<<lvi_ID_layer<<endl;
566 }
567 if( myCgemTruth && (currentTrackStatus==fStopAndKill || lv_status_post==fGeomBoundary) )
568 {
569 //Add by sunxh
570 int tmp[2000];
571 for(int ii =0; ii < hit_ID_vector.size();ii++){
572 tmp[ii] = hit_ID_vector[ii];
573 }
574 myCgemTruth->AddIdentifier(tmp,hit_ID_vector.size());
575 //TArrayI hit_tmp=myCgemTruth->GetIdentifier();
576 //end
577
578 m_collection_truth->insert(myCgemTruth);
579 myCgemTruth=NULL;
580
581 //Add by sunxh
582 hit_ID_vector.clear();
583 //end
584
585 if(printFlag)
586 cout<<"save a myCgemTruth: prePos is Boundary: "<<(lv_status_pre==fGeomBoundary)<<", postPos is Boundary:"<<(lv_status_post==fGeomBoundary)<<", trkID="<<lvi_ID_track<<", LayerID="<<lvi_ID_layer<<endl;
587 }
588 }//if(newVer==0)
589 } /* End of 'if (m_collection_truth)' */
590 } /* End of 'if (lvs_volume_name == "Gap_D_Logic")' */
591
592 return true;
593}
594
595//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
596/* Invoked at the end of processing an event (even the event is aborted) */
597void BesCgemSD::EndOfEvent(G4HCofThisEvent* HCE)
598{
599 IMessageSvc* msgSvc;
600 Gaudi::svcLocator() -> service("MessageSvc", msgSvc);
601 MsgStream log(msgSvc, "BesCgemSD::EndOfEvent()");
602
603 bool printFlag=false;
604 m_map_truth.clear();
605
606 /* Get the unique ID for this collection */
607 /* GetCollectionID() is a heavy operation, shouldn't be called for every event. */
608 /* It's available after SD is constructed and registed to G4SDManager. */
609 /* It can't be invoked in the constructor. */
610 static G4int HCID = -1;
611 if (HCID < 0)
612 {
613 HCID = G4SDManager::GetSDMpointer()->GetCollectionID(collectionName[0]);
614 }
615
616 /* Attach hit collection to G4HCofThisEvent object */
617 if(printFlag)
618 {
619 log<< MSG::INFO << "-------------------------------------------------------------------------------" << endreq;
620 log<< MSG::INFO << "CgemSim::BesCgemSD::EndOfEvent(), Entries of hits collection is "
621 << m_collection_hit->entries() << endreq;
622 }
623
624 HCE->AddHitsCollection(HCID, m_collection_hit);
625
626 if (verboseLevel > 0)
627 {
628 m_collection_hit->PrintAllHits();
629 }
630}
631
632//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
BesCgemHitsCollection
G4THitsCollection< BesCgemHit > BesCgemHitsCollection
Definition: InstallArea/include/CgemSim/CgemSim/BesCgemHit.hh:149
BesCgemHit::SetTotalEnergyDeposit
void SetTotalEnergyDeposit(G4double f_E_deposit)
Definition: InstallArea/include/CgemSim/CgemSim/BesCgemHit.hh:81
BesCgemHit::SetGlobalTime
void SetGlobalTime(G4double f_global_time)
Definition: InstallArea/include/CgemSim/CgemSim/BesCgemHit.hh:80
BesCgemHit::AddIdentifier
void AddIdentifier(G4int f_ID_Identifier[2000], G4int N_dim)
Definition: InstallArea/include/CgemSim/CgemSim/BesCgemHit.hh:92
BesCgemHit::SetMomentumOfPostPoint
void SetMomentumOfPostPoint(G4ThreeVector f_P_post)
Definition: InstallArea/include/CgemSim/CgemSim/BesCgemHit.hh:86
BesCgemHit::SetPositionOfPrePoint
void SetPositionOfPrePoint(G4ThreeVector f_XYZ_pre)
Definition: InstallArea/include/CgemSim/CgemSim/BesCgemHit.hh:83
BesCgemHit::SetMomentumOfPrePoint
void SetMomentumOfPrePoint(G4ThreeVector f_P_pre)
Definition: InstallArea/include/CgemSim/CgemSim/BesCgemHit.hh:85
BesCgemHit::SetPositionOfPostPoint
void SetPositionOfPostPoint(G4ThreeVector f_XYZ_post)
Definition: InstallArea/include/CgemSim/CgemSim/BesCgemHit.hh:84
BesCgemSD::EndOfTruthEvent
void EndOfTruthEvent(const G4Event *)
Definition: BesCgemSD.cc:164
BesCgemSD::EndOfEvent
void EndOfEvent(G4HCofThisEvent *)
Definition: BesCgemSD.cc:597
BesCgemSD::ProcessHits
G4bool ProcessHits(G4Step *, G4TouchableHistory *)
Definition: BesCgemSD.cc:257
BesCgemSD::BeginOfTruthEvent
void BeginOfTruthEvent(const G4Event *)
Definition: BesCgemSD.cc:155
BesCgemSD::BesCgemSD
BesCgemSD(G4String)
Definition: BesCgemSD.cc:80
BesCgemSD::Initialize
void Initialize(G4HCofThisEvent *)
Definition: BesCgemSD.cc:138
BesSensitiveDetector::GetCurrentTrackIndex
void GetCurrentTrackIndex(G4int &trackIndex, G4int &g4TrackId) const
Definition: BesSensitiveDetector.cc:30