106 G4SDManager::GetSDMpointer()->GetHCtable();
117 G4SDManager::GetSDMpointer()->GetCollectionID(
"physicalCellsFullSim");
122 G4SDManager::GetSDMpointer()->GetCollectionID(
"physicalCellsFastSim");
125 auto hitsCollection =
127 auto physicalFullHitsCollection =
130 auto physicalFastHitsCollection =
134 if(hitsCollection ==
nullptr)
136 G4ExceptionDescription msg;
138 G4Exception(
"Par04EventAction::GetHitsCollection()",
"MyCode0001", FatalException, msg);
140 if(physicalFullHitsCollection ==
nullptr)
142 G4ExceptionDescription msg;
144 G4Exception(
"Par04EventAction::GetHitsCollection()",
"MyCode0001", FatalException, msg);
146 if(physicalFastHitsCollection ==
nullptr)
148 G4ExceptionDescription msg;
150 G4Exception(
"Par04EventAction::GetHitsCollection()",
"MyCode0001", FatalException, msg);
153 auto analysisManager = G4AnalysisManager::Instance();
174 G4EventManager::GetEventManager()->GetConstCurrentEvent()->GetPrimaryVertex();
175 auto primaryParticle = primaryVertex->GetPrimary(0);
176 G4double primaryEnergy = primaryParticle->GetTotalEnergy();
179 auto primaryDirection = primaryParticle->GetMomentumDirection();
180 auto primaryEntrance =
181 primaryVertex->GetPosition() - primaryVertex->GetPosition().z() * primaryDirection;
196 G4double totalEnergy = 0;
203 G4int numNonZeroThresholdCells = 0;
204 G4double tDistance = 0., rDistance = 0., phiDistance = 0.;
205 G4double tFirstMoment = 0., tSecondMoment = 0.;
206 G4double rFirstMoment = 0., rSecondMoment = 0.;
207 G4double phiMean = 0.;
208 for(
size_t iHit = 0; iHit < hitsCollection->entries(); iHit++)
210 hit =
static_cast<Par04Hit*
>(hitsCollection->GetHit(iHit));
219 totalEnergy += hitEn;
224 tFirstMoment += hitEn * tDistance;
225 rFirstMoment += hitEn * rDistance;
226 phiMean += hitEn * phiDistance;
227 analysisManager->FillH1(4, tDistance, hitEn);
228 analysisManager->FillH1(5, rDistance, hitEn);
229 analysisManager->FillH1(10, hitType);
232 fCalEdep[numNonZeroThresholdCells] = hitEn;
233 fCalRho[numNonZeroThresholdCells] = hitRho;
234 fCalPhi[numNonZeroThresholdCells] = hitPhi;
235 fCalZ[numNonZeroThresholdCells] = hitZ;
236 numNonZeroThresholdCells++;
237 analysisManager->FillH1(13, std::log10(hitEn));
238 analysisManager->FillH1(15, hitNum);
242 tFirstMoment /= totalEnergy;
243 rFirstMoment /= totalEnergy;
244 phiMean /= totalEnergy;
245 analysisManager->FillH1(0, primaryEnergy / GeV);
246 analysisManager->FillH1(1, totalEnergy / GeV);
247 analysisManager->FillH1(2, totalEnergy / primaryEnergy);
248 analysisManager->FillH1(3,
fTimer.GetRealElapsed());
249 analysisManager->FillH1(6, tFirstMoment);
250 analysisManager->FillH1(7, rFirstMoment);
251 analysisManager->FillH1(12, numNonZeroThresholdCells);
252 analysisManager->FillH1(14, totalNum);
254 fCalEdep.resize(numNonZeroThresholdCells);
255 fCalRho.resize(numNonZeroThresholdCells);
256 fCalPhi.resize(numNonZeroThresholdCells);
257 fCalZ.resize(numNonZeroThresholdCells);
258 analysisManager->FillNtupleDColumn(0, 0, primaryEnergy);
259 analysisManager->FillNtupleDColumn(0, 1,
fTimer.GetRealElapsed());
261 for(
size_t iHit = 0; iHit < hitsCollection->entries(); iHit++)
263 hit =
static_cast<Par04Hit*
>(hitsCollection->GetHit(iHit));
273 tSecondMoment += hitEn * std::pow(tDistance - tFirstMoment, 2);
274 rSecondMoment += hitEn * std::pow(rDistance - rFirstMoment, 2);
275 analysisManager->FillH1(11, phiDistance - phiMean, hitEn);
278 tSecondMoment /= totalEnergy;
279 rSecondMoment /= totalEnergy;
280 analysisManager->FillH1(8, tSecondMoment);
281 analysisManager->FillH1(9, rSecondMoment);
284 G4double totalPhysicalEnergy = 0;
291 numNonZeroThresholdCells = 0;
292 for(
size_t iHit = 0; iHit < physicalFullHitsCollection->entries(); iHit++)
294 hit =
static_cast<Par04Hit*
>(physicalFullHitsCollection->GetHit(iHit));
302 totalPhysicalEnergy += hitEn;
310 numNonZeroThresholdCells++;
311 analysisManager->FillH1(19, std::log10(hitEn));
312 analysisManager->FillH1(21, hitNum);
315 }
for(
size_t iHit = 0; iHit < physicalFastHitsCollection->entries(); iHit++)
317 hit =
static_cast<Par04Hit*
>(physicalFastHitsCollection->GetHit(iHit));
325 totalPhysicalEnergy += hitEn;
333 numNonZeroThresholdCells++;
334 analysisManager->FillH1(19, std::log10(hitEn));
335 analysisManager->FillH1(21, hitNum);
339 analysisManager->FillH1(16, totalPhysicalEnergy / GeV);
340 analysisManager->FillH1(17, totalPhysicalEnergy / primaryEnergy);
341 analysisManager->FillH1(18, numNonZeroThresholdCells);
342 analysisManager->FillH1(20, totalNum);
347 analysisManager->AddNtupleRow(0);
348 analysisManager->AddNtupleRow(1);
349 analysisManager->AddNtupleRow(2);