Skip to content

Commit

Permalink
[pre-commit.ci] Automatic python and c++ formatting
Browse files Browse the repository at this point in the history
  • Loading branch information
@pre-commit-ci
pre-commit-ci[bot] committed May 6, 2024
1 parent fece4e7 commit fddb8d9
Show file tree
Hide file tree
Showing 5 changed files with 137 additions and 127 deletions.
2 changes: 1 addition & 1 deletion core/opengate_core/opengate_lib/GateChemistryActor.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -40,7 +40,7 @@ GateChemistryActor::GateChemistryActor(pybind11::dict &user_info)
_timeStepModelStr = DictGetStr(user_info, "timestep_model");
_endTime = DictGetDouble(user_info, "end_time");
_reactions = getReactionInputs(user_info, "reactions");
_moleculeCounterVerbose = DictGetInt(user_info, "molecule_counter_verbose");
_moleculeCounterVerbose = DictGetInt(user_info, "molecule_counter_verbose");

setTimeBinsCount(DictGetInt(user_info, "time_bins_count"));
}
Expand Down
207 changes: 107 additions & 100 deletions core/opengate_core/opengate_lib/GateChemistryLongTimeActor.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -8,16 +8,16 @@
#include "GateChemistryLongTimeActor.h"

#include <CLHEP/Units/SystemOfUnits.h>
#include <G4DNABoundingBox.hh>
#include <G4DNAChemistryManager.hh>
#include <G4DNAMolecularReactionTable.hh>
#include <G4DNABoundingBox.hh>
#include <G4MoleculeTable.hh>
#include <G4EmDNAChemistry_option3.hh>
#include <G4EmDNAPhysics_option3.hh>
#include <G4EmParameters.hh>
#include <G4EventManager.hh>
#include <G4H2O.hh>
#include <G4MoleculeCounter.hh>
#include <G4MoleculeTable.hh>
#include <G4RunManager.hh>
#include <G4Scheduler.hh>
#include <G4THitsMap.hh>
Expand All @@ -26,85 +26,88 @@
#include <G4UnitsTable.hh>
#include <G4VChemistryWorld.hh>
#include <G4ios.hh>
#include <memory>
#include <initializer_list>
#include <memory>
#include <pybind11/pybind11.h>
#include <pybind11/pytypes.h>

#include "GateHelpersDict.h"
#include <string>
#include <vector>
#include "GateHelpersDict.h"

#include "../g4_bindings/chemistryadaptator.h"
#include "GateVActor.h"

GateChemistryLongTimeActor::ChemistryWorld::ChemistryWorld(GateChemistryLongTimeActor *actor):
_actor(actor)
{
}
GateChemistryLongTimeActor::ChemistryWorld::ChemistryWorld(
GateChemistryLongTimeActor *actor)
: _actor(actor) {}

void GateChemistryLongTimeActor::ChemistryWorld::ConstructChemistryBoundary() {
auto const &size = _actor->_boundarySize;
std::initializer_list<double> boundingBox{
// high, low
size[0], -size[0], // x
size[1], -size[1], // y
size[2], -size[2] // z
};
fpChemistryBoundary = std::make_unique<G4DNABoundingBox>(boundingBox);
auto const &size = _actor->_boundarySize;
std::initializer_list<double> boundingBox{
// high, low
size[0], -size[0], // x
size[1], -size[1], // y
size[2], -size[2] // z
};
fpChemistryBoundary = std::make_unique<G4DNABoundingBox>(boundingBox);
}

void GateChemistryLongTimeActor::ChemistryWorld::ConstructChemistryComponents() {
constexpr auto water = 55.3; // from NIST material database
constexpr auto moleLiter = CLHEP::mole * CLHEP::liter;
void GateChemistryLongTimeActor::ChemistryWorld::
ConstructChemistryComponents() {
constexpr auto water = 55.3; // from NIST material database
constexpr auto moleLiter = CLHEP::mole * CLHEP::liter;

constexpr double pKw = 14; // at 25°C pK of water is 14
auto const & pH = _actor->_pH;
auto const &pH = _actor->_pH;

auto* moleculeTable = G4MoleculeTable::Instance();
auto *moleculeTable = G4MoleculeTable::Instance();

auto* mH2O = moleculeTable->GetConfiguration("H2O");
auto *mH2O = moleculeTable->GetConfiguration("H2O");
fpChemicalComponent[mH2O] = water / moleLiter;

auto* mH3Op = moleculeTable->GetConfiguration("H3Op(B)");
fpChemicalComponent[mH3Op] = std::pow(10, - pH) / moleLiter; // pH = 7
auto *mH3Op = moleculeTable->GetConfiguration("H3Op(B)");
fpChemicalComponent[mH3Op] = std::pow(10, -pH) / moleLiter; // pH = 7

auto* mOHm = moleculeTable->GetConfiguration("OHm(B)");
auto *mOHm = moleculeTable->GetConfiguration("OHm(B)");
fpChemicalComponent[mOHm] = std::pow(10, -(pKw - pH)) / moleLiter; // pH = 7

auto* mO2 = moleculeTable->GetConfiguration("O2");
auto *mO2 = moleculeTable->GetConfiguration("O2");
fpChemicalComponent[mO2] = (0. / 100) * 0.0013 / moleLiter;

// apply scavanger configurations
for(auto const &scavengerConfig : _actor->_scavengerConfigs) {
auto* m = moleculeTable->GetConfiguration(scavengerConfig.species);
auto const &concentration = scavengerConfig.concentration;
auto const &unitWithPrefix = scavengerConfig.unit;
if(not unitWithPrefix.empty()) {
char unit = *unitWithPrefix.rbegin();
if(unit == 'M') {
double factor = 1;
if(unitWithPrefix.size() > 1) {
auto prefix = unitWithPrefix.substr(0, unitWithPrefix.size()-1);
if(prefix == "u") factor = 1e-6;
else if(prefix == "m") factor = 1e-3;
else // TODO error
;
}

fpChemicalComponent[m] = factor * concentration / moleLiter;
} else if(unit == '%') {
constexpr auto o2 = 0.0013;
double concentrationInM =
fpChemicalComponent[m] = (concentration / 100) * o2 / moleLiter;
}
}
}
// apply scavanger configurations
for (auto const &scavengerConfig : _actor->_scavengerConfigs) {
auto *m = moleculeTable->GetConfiguration(scavengerConfig.species);
auto const &concentration = scavengerConfig.concentration;
auto const &unitWithPrefix = scavengerConfig.unit;
if (not unitWithPrefix.empty()) {
char unit = *unitWithPrefix.rbegin();
if (unit == 'M') {
double factor = 1;
if (unitWithPrefix.size() > 1) {
auto prefix = unitWithPrefix.substr(0, unitWithPrefix.size() - 1);
if (prefix == "u")
factor = 1e-6;
else if (prefix == "m")
factor = 1e-3;
else // TODO error
;
}

fpChemicalComponent[m] = factor * concentration / moleLiter;
} else if (unit == '%') {
constexpr auto o2 = 0.0013;
double concentrationInM = fpChemicalComponent[m] =
(concentration / 100) * o2 / moleLiter;
}
}
}
}

/* *** */

GateChemistryLongTimeActor::GateChemistryLongTimeActor(pybind11::dict &user_info)
GateChemistryLongTimeActor::GateChemistryLongTimeActor(
pybind11::dict &user_info)
: GateVActor(user_info, true) {
fActions.insert("NewStage");
fActions.insert("EndOfRunAction");
Expand All @@ -115,28 +118,29 @@ GateChemistryLongTimeActor::GateChemistryLongTimeActor(pybind11::dict &user_info
_timeStepModelStr = DictGetStr(user_info, "timestep_model");
_endTime = DictGetDouble(user_info, "end_time");
_reactions = getReactionInputs(user_info, "reactions");
_moleculeCounterVerbose = DictGetInt(user_info, "molecule_counter_verbose");
_moleculeCounterVerbose = DictGetInt(user_info, "molecule_counter_verbose");

setTimeBinsCount(DictGetInt(user_info, "time_bins_count"));

_pH = DictGetDouble(user_info, "pH");
_scavengerConfigs = getScavengerConfigs(user_info, "scavengers");
_pH = DictGetDouble(user_info, "pH");
_scavengerConfigs = getScavengerConfigs(user_info, "scavengers");

_doseCutOff = DictGetDouble(user_info, "dose_cutoff");
_doseCutOff = DictGetDouble(user_info, "dose_cutoff");

_resetScavengerForEachBeam = DictGetBool(user_info, "reset_scavenger_for_each_beam");
_resetScavengerForEachBeam =
DictGetBool(user_info, "reset_scavenger_for_each_beam");

// TODO remove
_boundarySize = DictGetVecDouble(user_info, "boundary_size");
// TODO remove
_boundarySize = DictGetVecDouble(user_info, "boundary_size");
}

void GateChemistryLongTimeActor::Initialize(G4HCofThisEvent *hce) {
GateVActor::Initialize(hce);

_chemistryWorld = std::make_unique<ChemistryWorld>(this);
_chemistryWorld->ConstructChemistryComponents();
_chemistryWorld = std::make_unique<ChemistryWorld>(this);
_chemistryWorld->ConstructChemistryComponents();

G4Scheduler::Instance()->ResetScavenger(_resetScavengerForEachBeam);
G4Scheduler::Instance()->ResetScavenger(_resetScavengerForEachBeam);

G4MoleculeCounter::Instance()->SetVerbose(_moleculeCounterVerbose);
G4MoleculeCounter::Instance()->Use();
Expand Down Expand Up @@ -249,33 +253,36 @@ void GateChemistryLongTimeActor::SteppingAction(G4Step *step) {
edep *= step->GetPreStepPoint()->GetWeight();
_edepSum += edep;

auto *track = step->GetTrack();

if(track->GetParentID() == 0 && track->GetCurrentStepNumber() == 1) {
constexpr auto density = .001;
auto volume = _chemistryWorld->GetChemistryBoundary()->Volume();
double dose = (_edepSum / CLHEP::eV) / (density * volume * 6.242e+18);
if(dose > _doseCutOff) {
auto *eventManager = G4EventManager::GetEventManager();
auto *primary = eventManager->GetConstCurrentEvent()->GetPrimaryVertex()->GetPrimary();
auto const &name = primary->GetParticleDefinition()->GetParticleName();
auto energy = primary->GetKineticEnergy();
G4cout << "[GateChemistryLongTimeActor] stop beam line '" << name << "' (" << energy << " MeV) at dose: " << dose << " Gy" << G4endl;

track->SetTrackStatus(fStopAndKill);
auto const *secondaries = track->GetStep()->GetSecondaryInCurrentStep();
for (auto const *secondary : *secondaries) {
if (secondary != nullptr) {
// FIXME
// from UHDR example
// must find how to get non-const access to secondaries
auto *secondaryMut = const_cast<G4Track*>(secondary);
secondaryMut->SetTrackStatus(fStopAndKill);
}
}
eventManager->GetStackManager()->ClearUrgentStack();
}
}
auto *track = step->GetTrack();

if (track->GetParentID() == 0 && track->GetCurrentStepNumber() == 1) {
constexpr auto density = .001;
auto volume = _chemistryWorld->GetChemistryBoundary()->Volume();
double dose = (_edepSum / CLHEP::eV) / (density * volume * 6.242e+18);
if (dose > _doseCutOff) {
auto *eventManager = G4EventManager::GetEventManager();
auto *primary = eventManager->GetConstCurrentEvent()
->GetPrimaryVertex()
->GetPrimary();
auto const &name = primary->GetParticleDefinition()->GetParticleName();
auto energy = primary->GetKineticEnergy();
G4cout << "[GateChemistryLongTimeActor] stop beam line '" << name << "' ("
<< energy << " MeV) at dose: " << dose << " Gy" << G4endl;

track->SetTrackStatus(fStopAndKill);
auto const *secondaries = track->GetStep()->GetSecondaryInCurrentStep();
for (auto const *secondary : *secondaries) {
if (secondary != nullptr) {
// FIXME
// from UHDR example
// must find how to get non-const access to secondaries
auto *secondaryMut = const_cast<G4Track *>(secondary);
secondaryMut->SetTrackStatus(fStopAndKill);
}
}
eventManager->GetStackManager()->ClearUrgentStack();
}
}
}

void GateChemistryLongTimeActor::NewStage() {
Expand Down Expand Up @@ -335,24 +342,24 @@ pybind11::dict GateChemistryLongTimeActor::getData() const {

GateChemistryLongTimeActor::ScavengerConfigs
GateChemistryLongTimeActor::getScavengerConfigs(pybind11::dict &user_info,
std::string const &key) {
ScavengerConfigs scavengerConfigs;
std::string const &key) {
ScavengerConfigs scavengerConfigs;

auto configs = DictGetVecList(user_info, key);
for (auto const &config : configs) {
ScavengerConfig scavengerConfig;
auto configs = DictGetVecList(user_info, key);
for (auto const &config : configs) {
ScavengerConfig scavengerConfig;

scavengerConfig.species = config[0].cast<std::string>();
scavengerConfig.concentration = config[1].cast<double>();
scavengerConfig.unit = config[2].cast<std::string>();
}
scavengerConfig.species = config[0].cast<std::string>();
scavengerConfig.concentration = config[1].cast<double>();
scavengerConfig.unit = config[2].cast<std::string>();
}

return scavengerConfigs;
return scavengerConfigs;
}

GateChemistryLongTimeActor::ReactionInputs
GateChemistryLongTimeActor::getReactionInputs(pybind11::dict &user_info,
std::string const &key) {
std::string const &key) {
ReactionInputs reactionInputs;

auto reactions = DictGetVecList(user_info, key);
Expand Down
45 changes: 22 additions & 23 deletions core/opengate_core/opengate_lib/GateChemistryLongTimeActor.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -38,11 +38,11 @@ class GateChemistryLongTimeActor : public GateVActor {
[[nodiscard]] pybind11::dict getData() const;

public:
struct ScavengerConfig {
std::string species;
double concentration;
std::string unit;
};
struct ScavengerConfig {
std::string species;
double concentration;
std::string unit;
};

using ScavengerConfigs = std::vector<ScavengerConfig>;

Expand All @@ -67,21 +67,20 @@ class GateChemistryLongTimeActor : public GateVActor {
using SpeciesMap = std::map<SpeciesPtr, InnerSpeciesMap>;

private:
class ChemistryWorld: public G4VChemistryWorld {
public:
ChemistryWorld(GateChemistryLongTimeActor *actor);

void ConstructChemistryBoundary() override;
void ConstructChemistryComponents() override;
class ChemistryWorld : public G4VChemistryWorld {
public:
ChemistryWorld(GateChemistryLongTimeActor *actor);

private:
GateChemistryLongTimeActor *_actor;
};
void ConstructChemistryBoundary() override;
void ConstructChemistryComponents() override;

private:
GateChemistryLongTimeActor *_actor;
};

protected:
static ScavengerConfigs getScavengerConfigs(pybind11::dict &user_info,
std::string const &key);
static ScavengerConfigs getScavengerConfigs(pybind11::dict &user_info,
std::string const &key);

static ReactionInputs getReactionInputs(pybind11::dict &user_info,
std::string const &key);
Expand All @@ -99,17 +98,17 @@ class GateChemistryLongTimeActor : public GateVActor {
double _endTime;
std::vector<ReactionInput> _reactions;

std::unique_ptr<G4VChemistryWorld> _chemistryWorld;
std::unique_ptr<G4VChemistryWorld> _chemistryWorld;

double _pH = 7;
ScavengerConfigs _scavengerConfigs;
double _pH = 7;
ScavengerConfigs _scavengerConfigs;

double _doseCutOff = 0;
double _doseCutOff = 0;

bool _resetScavengerForEachBeam = false;
bool _resetScavengerForEachBeam = false;

// TODO remove
std::vector<double> _boundarySize;
// TODO remove
std::vector<double> _boundarySize;
};

#endif
4 changes: 2 additions & 2 deletions core/opengate_core/opengate_lib/pyGateChemistryLongTimeActor.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -14,8 +14,8 @@ namespace py = pybind11;

void init_GateChemistryLongTimeActor(py::module &m) {
py::class_<GateChemistryLongTimeActor,
std::unique_ptr<GateChemistryLongTimeActor, py::nodelete>, GateVActor>(
m, "GateChemistryLongTimeActor")
std::unique_ptr<GateChemistryLongTimeActor, py::nodelete>,
GateVActor>(m, "GateChemistryLongTimeActor")
.def(py::init<py::dict &>())
.def("get_times", &GateChemistryLongTimeActor::getTimes)
.def("get_data", &GateChemistryLongTimeActor::getData);
Expand Down
Loading

0 comments on commit fddb8d9

Please sign in to comment.