add: user action, edep statistic;

This commit is contained in:
liuyihui 2022-05-19 11:47:35 +08:00
parent 3183f8e6d2
commit a6dcf70dbf
43 changed files with 675 additions and 25 deletions

View File

@ -20,7 +20,6 @@ add_executable(DESCSS main.cpp ${sources} ${headers})
target_link_libraries(DESCSS ${Geant4_LIBRARIES})
set(DESCSS_SCRIPTS
electron.mac
vis.mac
assets/model.txt
)

View File

@ -120,7 +120,18 @@
"vmanip": "cpp",
"element": "cpp",
"tree": "cpp",
"vertices": "cpp"
"vertices": "cpp",
"c1d": "cpp",
"dps": "cpp",
"rotf": "cpp",
"vec2": "cpp",
"vec3": "cpp",
"vec4": "cpp",
"bsf": "cpp",
"senum": "cpp",
"value": "cpp",
"ntuple": "cpp",
"free_seg": "cpp"
}
}
}

View File

@ -23,6 +23,7 @@ public:
void SetParticleType(G4String s) { particleType = s; };
void ConstructSDandField();
G4String GetParticleType() const { return particleType; };
std::vector<G4String> GetSensitiveList() const { return sensitiveList; };
private:
void ConstructSectionSphere(G4LogicalVolume* fMotherLogical, G4double zBias);
@ -35,7 +36,8 @@ private:
private:
G4String sex = "Male";
G4String model = "MIRD";
G4String particleType = "TE";
G4String particleType = "GCR_H";
std::vector<G4String> sensitiveList;
std::map<std::string, G4bool> sensitivities;
G4HumanPhantomMaterial* material;
G4HumanPhantomMessenger* messenger;

22
include/EventAction.h Normal file
View File

@ -0,0 +1,22 @@
#ifndef DESCSS_EventAction_h
#define DESCSS_EventAction_h
#include "G4UserEventAction.hh"
#include "globals.hh"
class RunAction;
class EventAction : public G4UserEventAction {
public:
EventAction(RunAction* runAction);
~EventAction() override;
void BeginOfEventAction(const G4Event*) override;
void EndOfEventAction(const G4Event*) override;
private:
RunAction* fRunAction = nullptr;
std::map<G4String, G4double> fEdep;
};
#endif

24
include/RunAction.h Normal file
View File

@ -0,0 +1,24 @@
#ifndef DESCSS_RunAction_h
#define DESCSS_RunAction_h
#include "G4Accumulable.hh"
#include "G4UserRunAction.hh"
#include "globals.hh"
class G4Run;
class RunAction : public G4UserRunAction {
public:
RunAction();
~RunAction() override;
void BeginOfRunAction(const G4Run*) override;
void EndOfRunAction(const G4Run*) override;
void AddEdep(std::map<G4String, G4double> edep);
private:
std::vector<G4String> sensitiveList;
};
#endif

15
macro/GCR_Al.mac Normal file
View File

@ -0,0 +1,15 @@
# 多线程设置
/run/numberOfThreads 12
# verbose
/control/verbose 1
/run/verbose 0
/event/verbose 0
/tracking/verbose 0
/vis/verbose errors
/type/set GCR_Al
# 初始化
/run/initialize
/run/beamOn 100000

15
macro/GCR_Ar.mac Normal file
View File

@ -0,0 +1,15 @@
# 多线程设置
/run/numberOfThreads 12
# verbose
/control/verbose 1
/run/verbose 0
/event/verbose 0
/tracking/verbose 0
/vis/verbose errors
/type/set GCR_Ar
# 初始化
/run/initialize
/run/beamOn 100000

15
macro/GCR_B.mac Normal file
View File

@ -0,0 +1,15 @@
# 多线程设置
/run/numberOfThreads 12
# verbose
/control/verbose 1
/run/verbose 0
/event/verbose 0
/tracking/verbose 0
/vis/verbose errors
/type/set GCR_B
# 初始化
/run/initialize
/run/beamOn 100000

15
macro/GCR_Be.mac Normal file
View File

@ -0,0 +1,15 @@
# 多线程设置
/run/numberOfThreads 12
# verbose
/control/verbose 1
/run/verbose 0
/event/verbose 0
/tracking/verbose 0
/vis/verbose errors
/type/set GCR_Be
# 初始化
/run/initialize
/run/beamOn 100000

15
macro/GCR_C.mac Normal file
View File

@ -0,0 +1,15 @@
# 多线程设置
/run/numberOfThreads 12
# verbose
/control/verbose 1
/run/verbose 0
/event/verbose 0
/tracking/verbose 0
/vis/verbose errors
/type/set GCR_C
# 初始化
/run/initialize
/run/beamOn 100000

15
macro/GCR_Ca.mac Normal file
View File

@ -0,0 +1,15 @@
# 多线程设置
/run/numberOfThreads 12
# verbose
/control/verbose 1
/run/verbose 0
/event/verbose 0
/tracking/verbose 0
/vis/verbose errors
/type/set GCR_Ca
# 初始化
/run/initialize
/run/beamOn 100000

15
macro/GCR_Cl.mac Normal file
View File

@ -0,0 +1,15 @@
# 多线程设置
/run/numberOfThreads 12
# verbose
/control/verbose 1
/run/verbose 0
/event/verbose 0
/tracking/verbose 0
/vis/verbose errors
/type/set GCR_Cl
# 初始化
/run/initialize
/run/beamOn 100000

15
macro/GCR_Co.mac Normal file
View File

@ -0,0 +1,15 @@
# 多线程设置
/run/numberOfThreads 12
# verbose
/control/verbose 1
/run/verbose 0
/event/verbose 0
/tracking/verbose 0
/vis/verbose errors
/type/set GCR_Co
# 初始化
/run/initialize
/run/beamOn 100000

15
macro/GCR_Cr.mac Normal file
View File

@ -0,0 +1,15 @@
# 多线程设置
/run/numberOfThreads 12
# verbose
/control/verbose 1
/run/verbose 0
/event/verbose 0
/tracking/verbose 0
/vis/verbose errors
/type/set GCR_Cr
# 初始化
/run/initialize
/run/beamOn 100000

15
macro/GCR_F.mac Normal file
View File

@ -0,0 +1,15 @@
# 多线程设置
/run/numberOfThreads 12
# verbose
/control/verbose 1
/run/verbose 0
/event/verbose 0
/tracking/verbose 0
/vis/verbose errors
/type/set GCR_F
# 初始化
/run/initialize
/run/beamOn 100000

15
macro/GCR_Fe.mac Normal file
View File

@ -0,0 +1,15 @@
# 多线程设置
/run/numberOfThreads 12
# verbose
/control/verbose 1
/run/verbose 0
/event/verbose 0
/tracking/verbose 0
/vis/verbose errors
/type/set GCR_Fe
# 初始化
/run/initialize
/run/beamOn 100000

15
macro/GCR_H.mac Normal file
View File

@ -0,0 +1,15 @@
# 多线程设置
/run/numberOfThreads 12
# verbose
/control/verbose 1
/run/verbose 0
/event/verbose 0
/tracking/verbose 0
/vis/verbose errors
/type/set GCR_H
# 初始化
/run/initialize
/run/beamOn 100000

15
macro/GCR_He.mac Normal file
View File

@ -0,0 +1,15 @@
# 多线程设置
/run/numberOfThreads 12
# verbose
/control/verbose 1
/run/verbose 0
/event/verbose 0
/tracking/verbose 0
/vis/verbose errors
/type/set GCR_He
# 初始化
/run/initialize
/run/beamOn 100000

15
macro/GCR_K.mac Normal file
View File

@ -0,0 +1,15 @@
# 多线程设置
/run/numberOfThreads 12
# verbose
/control/verbose 1
/run/verbose 0
/event/verbose 0
/tracking/verbose 0
/vis/verbose errors
/type/set GCR_K
# 初始化
/run/initialize
/run/beamOn 100000

15
macro/GCR_Li.mac Normal file
View File

@ -0,0 +1,15 @@
# 多线程设置
/run/numberOfThreads 12
# verbose
/control/verbose 1
/run/verbose 0
/event/verbose 0
/tracking/verbose 0
/vis/verbose errors
/type/set GCR_Li
# 初始化
/run/initialize
/run/beamOn 100000

15
macro/GCR_Mg.mac Normal file
View File

@ -0,0 +1,15 @@
# 多线程设置
/run/numberOfThreads 12
# verbose
/control/verbose 1
/run/verbose 0
/event/verbose 0
/tracking/verbose 0
/vis/verbose errors
/type/set GCR_Mg
# 初始化
/run/initialize
/run/beamOn 100000

15
macro/GCR_Mn.mac Normal file
View File

@ -0,0 +1,15 @@
# 多线程设置
/run/numberOfThreads 12
# verbose
/control/verbose 1
/run/verbose 0
/event/verbose 0
/tracking/verbose 0
/vis/verbose errors
/type/set GCR_Mn
# 初始化
/run/initialize
/run/beamOn 100000

15
macro/GCR_N.mac Normal file
View File

@ -0,0 +1,15 @@
# 多线程设置
/run/numberOfThreads 12
# verbose
/control/verbose 1
/run/verbose 0
/event/verbose 0
/tracking/verbose 0
/vis/verbose errors
/type/set GCR_N
# 初始化
/run/initialize
/run/beamOn 100000

15
macro/GCR_Na.mac Normal file
View File

@ -0,0 +1,15 @@
# 多线程设置
/run/numberOfThreads 12
# verbose
/control/verbose 1
/run/verbose 0
/event/verbose 0
/tracking/verbose 0
/vis/verbose errors
/type/set GCR_Na
# 初始化
/run/initialize
/run/beamOn 100000

15
macro/GCR_Ne.mac Normal file
View File

@ -0,0 +1,15 @@
# 多线程设置
/run/numberOfThreads 12
# verbose
/control/verbose 1
/run/verbose 0
/event/verbose 0
/tracking/verbose 0
/vis/verbose errors
/type/set GCR_Ne
# 初始化
/run/initialize
/run/beamOn 100000

15
macro/GCR_Ni.mac Normal file
View File

@ -0,0 +1,15 @@
# 多线程设置
/run/numberOfThreads 12
# verbose
/control/verbose 1
/run/verbose 0
/event/verbose 0
/tracking/verbose 0
/vis/verbose errors
/type/set GCR_Ni
# 初始化
/run/initialize
/run/beamOn 100000

15
macro/GCR_O.mac Normal file
View File

@ -0,0 +1,15 @@
# 多线程设置
/run/numberOfThreads 12
# verbose
/control/verbose 1
/run/verbose 0
/event/verbose 0
/tracking/verbose 0
/vis/verbose errors
/type/set GCR_O
# 初始化
/run/initialize
/run/beamOn 100000

15
macro/GCR_P.mac Normal file
View File

@ -0,0 +1,15 @@
# 多线程设置
/run/numberOfThreads 12
# verbose
/control/verbose 1
/run/verbose 0
/event/verbose 0
/tracking/verbose 0
/vis/verbose errors
/type/set GCR_P
# 初始化
/run/initialize
/run/beamOn 100000

15
macro/GCR_S.mac Normal file
View File

@ -0,0 +1,15 @@
# 多线程设置
/run/numberOfThreads 12
# verbose
/control/verbose 1
/run/verbose 0
/event/verbose 0
/tracking/verbose 0
/vis/verbose errors
/type/set GCR_S
# 初始化
/run/initialize
/run/beamOn 100000

15
macro/GCR_Sc.mac Normal file
View File

@ -0,0 +1,15 @@
# 多线程设置
/run/numberOfThreads 12
# verbose
/control/verbose 1
/run/verbose 0
/event/verbose 0
/tracking/verbose 0
/vis/verbose errors
/type/set GCR_Sc
# 初始化
/run/initialize
/run/beamOn 100000

15
macro/GCR_Si.mac Normal file
View File

@ -0,0 +1,15 @@
# 多线程设置
/run/numberOfThreads 12
# verbose
/control/verbose 1
/run/verbose 0
/event/verbose 0
/tracking/verbose 0
/vis/verbose errors
/type/set GCR_Si
# 初始化
/run/initialize
/run/beamOn 100000

15
macro/GCR_Ti.mac Normal file
View File

@ -0,0 +1,15 @@
# 多线程设置
/run/numberOfThreads 12
# verbose
/control/verbose 1
/run/verbose 0
/event/verbose 0
/tracking/verbose 0
/vis/verbose errors
/type/set GCR_Ti
# 初始化
/run/initialize
/run/beamOn 100000

15
macro/GCR_V.mac Normal file
View File

@ -0,0 +1,15 @@
# 多线程设置
/run/numberOfThreads 12
# verbose
/control/verbose 1
/run/verbose 0
/event/verbose 0
/tracking/verbose 0
/vis/verbose errors
/type/set GCR_V
# 初始化
/run/initialize
/run/beamOn 100000

View File

@ -2,14 +2,14 @@
/run/numberOfThreads 12
# verbose
/control/verbose 2
/control/verbose 1
/run/verbose 0
/event/verbose 0
/tracking/verbose 0
/vis/verbose errors
/type/set TE
# 初始化
/run/initialize
/type/set TE
/run/beamOn 50
/run/beamOn 923242

15
macro/TP.mac Normal file
View File

@ -0,0 +1,15 @@
# 多线程设置
/run/numberOfThreads 12
# verbose
/control/verbose 1
/run/verbose 0
/event/verbose 0
/tracking/verbose 0
/vis/verbose errors
/type/set TP
# 初始化
/run/initialize
/run/beamOn 100000

View File

@ -1,5 +1,7 @@
#include "ActionInitialization.h"
#include "EventAction.h"
#include "PrimaryGeneratorAction.h"
#include "RunAction.h"
#include <time.h>
#include "Randomize.hh"
@ -15,8 +17,16 @@ ActionInitialization::ActionInitialization() {
ActionInitialization::~ActionInitialization() {}
void ActionInitialization::BuildForMaster() const {}
void ActionInitialization::BuildForMaster() const {
RunAction* runAction = new RunAction;
SetUserAction(runAction);
}
void ActionInitialization::Build() const {
SetUserAction(new PrimaryGeneratorAction);
RunAction* runAction = new RunAction;
SetUserAction(runAction);
EventAction* eventAction = new EventAction(runAction);
SetUserAction(eventAction);
}

View File

@ -23,6 +23,49 @@
DetectorConstruction::DetectorConstruction() {
material = new G4HumanPhantomMaterial();
messenger = new G4HumanPhantomMessenger(this);
sensitiveList.emplace_back("logicalHead");
sensitiveList.emplace_back("logicalSkull");
sensitiveList.emplace_back("logicalBrain");
sensitiveList.emplace_back("logicalTrunk");
sensitiveList.emplace_back("logicalLeftLeg");
sensitiveList.emplace_back("logicalRightLeg");
sensitiveList.emplace_back("logicalLeftArmBone");
sensitiveList.emplace_back("logicalRightArmBone");
sensitiveList.emplace_back("logicalLeftLegBone");
sensitiveList.emplace_back("logicalRightLegBone");
sensitiveList.emplace_back("logicalUpperSpine");
sensitiveList.emplace_back("logicalLeftScapula");
sensitiveList.emplace_back("logicalRightScapula");
sensitiveList.emplace_back("logicalLeftAdrenal");
sensitiveList.emplace_back("logicalRightAdrenal");
sensitiveList.emplace_back("logicalThymus");
sensitiveList.emplace_back("logicalLeftClavicle");
sensitiveList.emplace_back("logicalRightClavicle");
sensitiveList.emplace_back("logicalSmallIntestine");
sensitiveList.emplace_back("logicalRibCage");
sensitiveList.emplace_back("logicalMiddleLowerSpine");
sensitiveList.emplace_back("logicalPelvis");
sensitiveList.emplace_back("logicalStomach");
sensitiveList.emplace_back("logicalUpperLargeIntestine");
sensitiveList.emplace_back("logicalLowerLargeIntestine");
sensitiveList.emplace_back("logicalSpleen");
sensitiveList.emplace_back("logicalPancreas");
sensitiveList.emplace_back("logicalLiver");
sensitiveList.emplace_back("logicalLeftKidney");
sensitiveList.emplace_back("logicalRightKidney");
sensitiveList.emplace_back("logicalUrinaryBladder");
sensitiveList.emplace_back("logicalLeftLung");
sensitiveList.emplace_back("logicalRightLung");
sensitiveList.emplace_back("logicalThyroid");
sensitiveList.emplace_back("logicalLeftOvary");
sensitiveList.emplace_back("logicalRightOvary");
sensitiveList.emplace_back("logicalUterus");
sensitiveList.emplace_back("logicalLeftBreast");
sensitiveList.emplace_back("logicalRightBreast");
sensitiveList.emplace_back("logicalMaleGenitalia");
sensitiveList.emplace_back("logicalLeftTeste");
sensitiveList.emplace_back("logicalRightTeste");
}
DetectorConstruction::~DetectorConstruction() {
@ -289,7 +332,7 @@ void DetectorConstruction::ConstructHumanPhantom(G4VPhysicalVolume* fMotherPhysi
builder->BuildRightKidney("green", true, sensitivities["RightKidney"]);
builder->BuildUrinaryBladder("green", true, sensitivities["UrinaryBladder"]);
builder->BuildHeart("red", true, sensitivities["Hearth"]); // to do MIRD
// builder->BuildHeart("red", true, sensitivities["Hearth"]); // to do MIRD
builder->BuildLeftLung("blue", true, sensitivities["LeftLung"]);
builder->BuildRightLung("blue", true, sensitivities["RightLung"]);
builder->BuildThyroid("orange", true, sensitivities["Thyroid"]);
@ -299,13 +342,8 @@ void DetectorConstruction::ConstructHumanPhantom(G4VPhysicalVolume* fMotherPhysi
builder->BuildRightOvary("purple", true, sensitivities["RightOvary"]);
builder->BuildUterus("purple", true, sensitivities["Uterus"]);
if (model == "MIRD") {
builder->BuildLeftBreast("purple", true, sensitivities["LeftBreast"]);
builder->BuildRightBreast("purple", true, sensitivities["RightBreast"]);
} else {
builder->BuildVoxelLeftBreast("purple", false, sensitivities["LeftBreast"]);
builder->BuildVoxelRightBreast("purple", false, sensitivities["RightBreast"]);
}
builder->BuildLeftBreast("purple", true, sensitivities["LeftBreast"]);
builder->BuildRightBreast("purple", true, sensitivities["RightBreast"]);
}
if (sex == "Male") {
@ -349,6 +387,8 @@ G4VPhysicalVolume* DetectorConstruction::Construct() {
"Phantom", logicWorld, false, 0, true);
ConstructHumanPhantom(physicsPhantom);
// for (auto i = sensitivities.begin(); i != sensitivities.end(); ++i) sensitiveList.emplace_back(i->first);
return physicsWorld;
}
@ -399,14 +439,19 @@ void DetectorConstruction::ConstructSDandField() {
SetSensitiveDetector("logicalSmallIntestine", SD);
SetSensitiveDetector("logicalRibCage", SD);
SetSensitiveDetector("logicalMiddleLowerSpine", SD);
SetSensitiveDetector("logicalPelvis", SD);
SetSensitiveDetector("logicalStomach", SD);
SetSensitiveDetector("logicalUpperLargeIntestine", SD);
SetSensitiveDetector("logicalLowerLargeIntestine", SD);
SetSensitiveDetector("logicalSpleen", SD);
SetSensitiveDetector("logicalPancreas", SD);
SetSensitiveDetector("logicalLiver", SD);
SetSensitiveDetector("logicalLeftKidney", SD);
SetSensitiveDetector("logicalRightKidney", SD);
SetSensitiveDetector("logicalUrinaryBladder", SD);
SetSensitiveDetector("logicalLeftLung", SD);
SetSensitiveDetector("logicalRightLung", SD);
SetSensitiveDetector("logicalThyroid", SD);
if (sex == "Female") {
SetSensitiveDetector("logicalLeftOvary", SD);

22
src/EventAction.cpp Normal file
View File

@ -0,0 +1,22 @@
#include "EventAction.h"
#include "G4HumanPhantomHit.h"
#include "RunAction.h"
#include "G4Event.hh"
#include "G4VHitsCollection.hh"
EventAction::EventAction(RunAction* runAction) : fRunAction(runAction) {}
EventAction::~EventAction() {}
void EventAction::BeginOfEventAction(const G4Event*) { fEdep.clear(); }
void EventAction::EndOfEventAction(const G4Event* event) {
G4HumanPhantomHitsCollection* hc = (G4HumanPhantomHitsCollection*)event->GetHCofThisEvent()->GetHC(0);
int numHits = hc->entries();
for (int i = 0; i < numHits; i++) {
G4String name = (*hc)[0]->GetBodyPartID();
fEdep[name] += (*hc)[i]->GetEdep();
}
if (numHits > 0) fRunAction->AddEdep(fEdep);
}

View File

@ -36,9 +36,4 @@ G4bool G4HumanPhantomSD::ProcessHits(G4Step* aStep, G4TouchableHistory*) {
return true;
}
void G4HumanPhantomSD::EndOfEvent(G4HCofThisEvent*) {
// G4int NbHits = collection->entries();
// G4cout << "\n-------->Hits Collection: in this event, there are " << NbHits
// << " hits in the tracker chambers: " << G4endl;
// for (G4int i = 0; i < NbHits; i++) (*collection)[i]->Print();
}
void G4HumanPhantomSD::EndOfEvent(G4HCofThisEvent*) {}

View File

@ -58,7 +58,7 @@ G4VPhysicalVolume* G4MIRDLiver::Construct(const G4String& volumeName, G4VPhysica
G4SubtractionSolid* liver = new G4SubtractionSolid("Liver", firstLiver, subtrLiver, rm_relative,
G4ThreeVector(10.0 * cm, 0.0 * cm, 0.0 * cm));
G4LogicalVolume* logicLiver = new G4LogicalVolume(liver, soft, "LiverVolume", 0, 0, 0);
G4LogicalVolume* logicLiver = new G4LogicalVolume(liver, soft, "logical" + volumeName, 0, 0, 0);
// Define rotation and position here!
G4RotationMatrix* rm = new G4RotationMatrix();

View File

@ -60,7 +60,7 @@ G4VPhysicalVolume* G4MIRDThyroid::Construct(const G4String& volumeName, G4VPhysi
G4SubtractionSolid* thyroid = new G4SubtractionSolid(
"SecondThyroid", SecondThyroid, SubtrThyroid, relative_matrix_2, G4ThreeVector(0.0 * cm, 0.0 * cm, -5.40 * cm));
G4LogicalVolume* logicThyroid = new G4LogicalVolume(thyroid, soft, "ThyroidVolume", 0, 0, 0);
G4LogicalVolume* logicThyroid = new G4LogicalVolume(thyroid, soft, "logical" + volumeName, 0, 0, 0);
G4RotationMatrix* rm = new G4RotationMatrix();
rm->rotateZ(180. * degree);

View File

@ -15,6 +15,7 @@ G4double C1, C2, C3, C4, sum;
PrimaryGeneratorAction::PrimaryGeneratorAction() {
G4ParticleDefinition* ion = DefineParticle();
fParticleGun = new G4ParticleGun();
fParticleGun->SetParticleCharge(this->Z * eplus);
fParticleGun->SetParticleDefinition(ion);
@ -52,7 +53,6 @@ G4ParticleDefinition* PrimaryGeneratorAction::DefineParticle() {
}
G4ParticleDefinition* ion;
fParticleGun = new G4ParticleGun();
if (Z == -1)
ion = G4ParticleTable::GetParticleTable()->FindParticle("e-");
else if (Z == 1)

85
src/RunAction.cpp Normal file
View File

@ -0,0 +1,85 @@
#include "DetectorConstruction.h"
#include "RunAction.h"
#include "G4AccumulableManager.hh"
#include "G4Run.hh"
#include "G4RunManager.hh"
#include "G4SystemOfUnits.hh"
RunAction::RunAction() {
G4RunManager::GetRunManager()->SetPrintProgress(10000);
const DetectorConstruction* detConstruction =
static_cast<const DetectorConstruction*>(G4RunManager::GetRunManager()->GetUserDetectorConstruction());
sensitiveList = detConstruction->GetSensitiveList();
G4AccumulableManager* accumulableManager = G4AccumulableManager::Instance();
for (int i = 0; i < sensitiveList.size(); i++) {
G4String name = sensitiveList[i];
accumulableManager->CreateAccumulable<G4double>(name + "Edep", 0);
accumulableManager->CreateAccumulable<G4double>(name + "Edep2", 0);
}
}
RunAction::~RunAction() {}
void RunAction::BeginOfRunAction(const G4Run*) {
G4RunManager::GetRunManager()->SetRandomNumberStore(false);
G4AccumulableManager* accumulableManager = G4AccumulableManager::Instance();
accumulableManager->Reset();
}
void RunAction::EndOfRunAction(const G4Run* run) {
G4int nofEvents = run->GetNumberOfEvent();
if (nofEvents == 0) return;
G4AccumulableManager* accumulableManager = G4AccumulableManager::Instance();
accumulableManager->Merge();
std::map<G4String, G4double> edep;
std::map<G4String, G4double> edep2;
std::map<G4String, G4double> rms;
for (int i = 0; i < sensitiveList.size(); i++) {
G4String name = sensitiveList[i];
G4Accumulable<G4double>* fEdep = accumulableManager->GetAccumulable<G4double>(name + "Edep");
G4Accumulable<G4double>* fEdep2 = accumulableManager->GetAccumulable<G4double>(name + "Edep2");
edep[name] = fEdep->GetValue();
edep2[name] = fEdep2->GetValue();
rms[name] = edep2[name] - edep[name] * edep[name] / nofEvents;
if (rms[name] > 0.)
rms[name] = std::sqrt(rms[name]);
else
rms[name] = 0.;
}
if (IsMaster()) {
const DetectorConstruction* detConstruction =
static_cast<const DetectorConstruction*>(G4RunManager::GetRunManager()->GetUserDetectorConstruction());
std::ofstream edepFile("result/edep_" + detConstruction->GetParticleType() + ".txt");
for (int i = 0; i < sensitiveList.size(); i++) {
G4String name = sensitiveList[i];
edepFile << "bodyPart: " << std::setw(30) << name << " | eDep: " << std::setw(15)
<< edep[name] / MeV << " MeV | eDep2: " << std::setw(15)
<< edep2[name] / MeV << " MeV | rms: " << std::setw(15)
<< rms[name] / MeV << " MeV" << G4endl;
}
G4cout << G4endl << "--------------------End of Global Run-----------------------" << G4endl;
} else
G4cout << G4endl << "--------------------End of Local Run------------------------" << G4endl;
}
void RunAction::AddEdep(std::map<G4String, G4double> edep) {
G4AccumulableManager* accumulableManager = G4AccumulableManager::Instance();
for (int i = 0; i < sensitiveList.size(); i++) {
G4String name = sensitiveList[i];
G4double value = edep[name] * MeV;
G4Accumulable<G4double>* fEdep = accumulableManager->GetAccumulable<G4double>(name + "Edep");
G4Accumulable<G4double>* fEdep2 = accumulableManager->GetAccumulable<G4double>(name + "Edep2");
*(fEdep) += value;
*(fEdep2) += value * value;
}
}