diff --git a/.gitlab-ci.yml b/.gitlab-ci.yml new file mode 100644 index 0000000000000000000000000000000000000000..925e5145d8e2950cc1707408dd08ed9258b7f7e4 --- /dev/null +++ b/.gitlab-ci.yml @@ -0,0 +1,32 @@ +############################################################################### +# (c) Copyright 2000-2018 CERN for the benefit of the LHCb Collaboration # +# # +# This software is distributed under the terms of the GNU General Public # +# Licence version 3 (GPL Version 3), copied verbatim in the file "COPYING". # +# # +# In applying this licence, CERN does not waive the privileges and immunities # +# granted to it by virtue of its status as an Intergovernmental Organization # +# or submit itself to any jurisdiction. # +############################################################################### +variables: + TARGET_BRANCH: master +check-copyright: + image: gitlab-registry.cern.ch/ci-tools/ci-worker:cc7 + script: + - curl -o lb-check-copyright "https://gitlab.cern.ch/lhcb-core/LbDevTools/raw/master/LbDevTools/SourceTools.py?inline=false" + - python lb-check-copyright origin/${TARGET_BRANCH} +check-formatting: + image: gitlab-registry.cern.ch/lhcb-docker/style-checker + script: + - if [ ! -e .clang-format ] ; then + - curl -o .clang-format "https://gitlab.cern.ch/lhcb-core/LbDevTools/raw/master/LbDevTools/data/default.clang-format?inline=false" + - echo '.clang-format' >> .gitignore + - git add .gitignore + - fi + - curl -o lb-format "https://gitlab.cern.ch/lhcb-core/LbDevTools/raw/master/LbDevTools/SourceTools.py?inline=false" + - python lb-format --format-patch apply-formatting.patch origin/${TARGET_BRANCH} + artifacts: + paths: + - apply-formatting.patch + when: on_failure + expire_in: 1 week diff --git a/Sim/GiGaMTExamples/FastSimulationExample/CMakeLists.txt b/Sim/GiGaMTExamples/FastSimulationExample/CMakeLists.txt new file mode 100644 index 0000000000000000000000000000000000000000..ca3d1e9b5dc5f53a249bdbe31befe3ee97faaaf1 --- /dev/null +++ b/Sim/GiGaMTExamples/FastSimulationExample/CMakeLists.txt @@ -0,0 +1,46 @@ +############################################################################### +# (c) Copyright 2021 CERN for the benefit of the LHCb Collaboration # +# # +# This software is distributed under the terms of the GNU General Public # +# Licence version 3 (GPL Version 3), copied verbatim in the file "COPYING". # +# # +# In applying this licence, CERN does not waive the privileges and immunities # +# granted to it by virtue of its status as an Intergovernmental Organization # +# or submit itself to any jurisdiction. # +############################################################################### +################################################################################ +# Package: GiGaMTExamples/FastSimulationExample +# +# Package that contains all necessary examples to test experiment-indpendent +# Geant4 hooks for fast simulation. +# +################################################################################# +gaudi_subdir(FastSimulationExample v1r0) + +gaudi_depends_on_subdirs(GaudiAlg + Sim/GiGaMTFactories + Sim/GiGaMTFastSimulation + Sim/GiGaMTCore) +AddHepMC3() +if(${Geant4_config_version} VERSION_LESS "10.06") + add_definitions(-DG4MULTITHREADED) + add_definitions(-DG4USE_STD11) +endif() + +gaudi_add_module(ImmediateDeposit + src/ImmediateDeposit/*.cpp + INCLUDE_DIRS GiGaMTFactories + LINK_LIBRARIES GaudiAlgLib GiGaMTCoreRunLib GiGaMTFastSimulationLib) + +add_dependencies(ImmediateDeposit HepMC3Ext) + +gaudi_add_module(ShowerDeposit + src/ShowerDeposit/*.cpp + INCLUDE_DIRS GiGaMTFactories + LINK_LIBRARIES GaudiAlgLib GiGaMTCoreRunLib GiGaMTFastSimulationLib) + +add_dependencies(ShowerDeposit HepMC3Ext) + +gaudi_env(SET FASTSIMULATIONEXAMPLETESTS \${FASTSIMULATIONEXAMPLEROOT}/tests) + +gaudi_add_test(QMTest QMTEST) diff --git a/Sim/GiGaMTExamples/FastSimulationExample/src/ImmediateDeposit/ImmediateDepositModel.cpp b/Sim/GiGaMTExamples/FastSimulationExample/src/ImmediateDeposit/ImmediateDepositModel.cpp new file mode 100644 index 0000000000000000000000000000000000000000..aa7e517c8975e4b88653544fae1f05bb61f6e049 --- /dev/null +++ b/Sim/GiGaMTExamples/FastSimulationExample/src/ImmediateDeposit/ImmediateDepositModel.cpp @@ -0,0 +1,73 @@ +/*****************************************************************************\ +* (c) Copyright 2021 CERN for the benefit of the LHCb Collaboration * +* * +* This software is distributed under the terms of the GNU General Public * +* Licence version 3 (GPL Version 3), copied verbatim in the file "COPYING". * +* * +* In applying this licence, CERN does not waive the privileges and immunities * +* granted to it by virtue of its status as an Intergovernmental Organization * +* or submit itself to any jurisdiction. * +\*****************************************************************************/ +#include "ImmediateDepositModel.h" +#include "Geant4/G4TransportationManager.hh" +#include "Geant4/G4VSensitiveDetector.hh" +#include "GiGaMTFastSimulation/GiGaMTFastSimModelFAC.h" +#include "GiGaMTFastSimulation/GiGaMTFastHit.h" + +namespace ImmediateDeposit { + +Model::Model( G4String modelName, G4Region* envelope ) : G4VFastSimulationModel( modelName, envelope ), m_sensResponse(new SensitiveResponse) {} + +G4bool Model::IsApplicable( const G4ParticleDefinition& ) { return true; } + +G4bool Model::ModelTrigger( const G4FastTrack& ) { + return true; +} + +void Model::DoIt( const G4FastTrack& aFastTrack, G4FastStep& aFastStep ) { + // kill the track so that it will no longer be propagated by G4 + aFastStep.KillPrimaryTrack(); + double initialEnergy = aFastTrack.GetPrimaryTrack()->GetKineticEnergy(); + auto position = aFastTrack.GetPrimaryTrack()->GetPosition(); + auto rot = G4RotationMatrix(); + rot.rotateY(position.getTheta()); + + if(!(std::fabs(position.x()) < 1e-3 && std::fabs(position.y()) < 1e-3)) { + rot.rotateZ(position.getPhi()); + } + + if ( m_eDepFrac <= 1.0 ) { + double tmpEnergyDeposit = initialEnergy * m_eDepFrac; + auto tmpPosition = position; + while ( tmpEnergyDeposit > m_eKillThr ) { + m_sensResponse->simulate(GiGaMTFastHit(tmpPosition, tmpEnergyDeposit), aFastTrack); + if (m_eDepFrac == 1.0) break; + initialEnergy -= tmpEnergyDeposit; + tmpEnergyDeposit = initialEnergy * m_eDepFrac; + tmpPosition = tmpPosition + rot * G4ThreeVector(0, 0, m_sensResponse->GetStepLength()); + } + } +} + +class ModelFactory : public GiGaMTFastSimModelFAC<Model> { + + Gaudi::Property<double> m_eDepFrac {this, "EnergyDepositFraction", 1.0}; + Gaudi::Property<double> m_eKillThr {this, "KillEnergyThreshold", 50.0 * CLHEP::MeV}; + Gaudi::Property<double> m_stepLength {this, "FakeStepLength", 0.0 * CLHEP::mm}; + + public: + virtual ~ModelFactory() = default; + using gaussino_base_class = GiGaMTFastSimModelFAC<Model>; + using gaussino_base_class::GiGaMTFastSimModelFAC; + + virtual Model* construct() const override { + auto tmp = gaussino_base_class::construct(); + tmp->SetEDepFrac(m_eDepFrac); + tmp->SetEKillThr(m_eKillThr); + tmp->SetStepLength(m_stepLength); + return tmp; + }; +}; +} + +DECLARE_COMPONENT_WITH_ID( ImmediateDeposit::ModelFactory, "ImmediateDepositModelFactory" ) diff --git a/Sim/GiGaMTExamples/FastSimulationExample/src/ImmediateDeposit/ImmediateDepositModel.h b/Sim/GiGaMTExamples/FastSimulationExample/src/ImmediateDeposit/ImmediateDepositModel.h new file mode 100644 index 0000000000000000000000000000000000000000..d851d588a1734ecf51cb0085924a66ecbaa53530 --- /dev/null +++ b/Sim/GiGaMTExamples/FastSimulationExample/src/ImmediateDeposit/ImmediateDepositModel.h @@ -0,0 +1,53 @@ +/*****************************************************************************\ +* (c) Copyright 2021 CERN for the benefit of the LHCb Collaboration * +* * +* This software is distributed under the terms of the GNU General Public * +* Licence version 3 (GPL Version 3), copied verbatim in the file "COPYING". * +* * +* In applying this licence, CERN does not waive the privileges and immunities * +* granted to it by virtue of its status as an Intergovernmental Organization * +* or submit itself to any jurisdiction. * +\*****************************************************************************/ +#pragma once + +#include "Geant4/G4VFastSimulationModel.hh" +#include "GiGaMTCoreMessage/IGiGaMessage.h" +#include "GiGaMTFastSimulation/GiGaMTFastSimSensitiveResponse.h" + +namespace ImmediateDeposit { + +class SensitiveResponse : public GiGaMTFastSimSensitiveResponse { + + double m_stepLength {0.0 * CLHEP::mm}; + + inline virtual void modifyFakeStep(G4Step* fakeStep) override { + // this is a hack, so that the model will work on detectors + // where it the step must be finite + fakeStep->SetStepLength(m_stepLength); + } + + public: + + inline void SetStepLength (double stepLength) { m_stepLength = stepLength; } + inline double GetStepLength () { return m_stepLength; } +}; + +class Model : public G4VFastSimulationModel, public GiGaMessage { + + std::unique_ptr<SensitiveResponse> m_sensResponse; + double m_eDepFrac {1.0}; + double m_eKillThr {50.0 * CLHEP::MeV}; + +public: + Model(G4String modelName, G4Region* envelope); + ~Model() = default; + + G4bool IsApplicable( const G4ParticleDefinition& aParticle ) override; + G4bool ModelTrigger( const G4FastTrack& aFastTrack ) override; + void DoIt( const G4FastTrack& aFastTrack, G4FastStep& aFastStep ) override; + + inline void SetEDepFrac(double eDepFrac) { m_eDepFrac = eDepFrac; }; + inline void SetEKillThr(double eKillThr) { m_eKillThr = eKillThr; }; + inline void SetStepLength(double stepLength) { m_sensResponse->SetStepLength(stepLength); }; +}; +} diff --git a/Sim/GiGaMTExamples/FastSimulationExample/src/ShowerDeposit/ShowerDepositModel.cpp b/Sim/GiGaMTExamples/FastSimulationExample/src/ShowerDeposit/ShowerDepositModel.cpp new file mode 100644 index 0000000000000000000000000000000000000000..a712bc58cbd4ba5feea862395a64cfc484994447 --- /dev/null +++ b/Sim/GiGaMTExamples/FastSimulationExample/src/ShowerDeposit/ShowerDepositModel.cpp @@ -0,0 +1,148 @@ +/*****************************************************************************\ +* (c) Copyright 2021 CERN for the benefit of the LHCb Collaboration * +* * +* This software is distributed under the terms of the GNU General Public * +* Licence version 3 (GPL Version 3), copied verbatim in the file "COPYING". * +* * +* In applying this licence, CERN does not waive the privileges and immunities * +* granted to it by virtue of its status as an Intergovernmental Organization * +* or submit itself to any jurisdiction. * +\*****************************************************************************/ +#include "ShowerDepositModel.h" +#include "Geant4/G4TransportationManager.hh" +#include "Geant4/G4VSensitiveDetector.hh" +#include "GiGaMTFastSimulation/GiGaMTFastSimModelFAC.h" +#include "GiGaMTFastSimulation/GiGaMTFastHit.h" + +#include "Geant4/G4Material.hh" +#include "Geant4/G4Electron.hh" +#include "Geant4/G4Positron.hh" +#include "Geant4/G4Gamma.hh" +#include "Geant4/Randomize.hh" + +namespace ShowerDeposit { + +Model::Model( G4String modelName, G4Region* envelope ) : G4VFastSimulationModel( modelName, envelope ), m_sensResponse(new SensitiveResponse) {} + +G4bool Model::IsApplicable( const G4ParticleDefinition& aParticleType) { + return &aParticleType == G4Electron::ElectronDefinition() || + &aParticleType == G4Positron::PositronDefinition() || + &aParticleType == G4Gamma::GammaDefinition(); +} + +G4bool Model::ModelTrigger( const G4FastTrack& ) { + return true; +} + +void Model::DoIt( const G4FastTrack& aFastTrack, G4FastStep& aFastStep ) { + // kill the track so that it will no longer be propagated by G4 + aFastStep.KillPrimaryTrack(); + + double energy = aFastTrack.GetPrimaryTrack()->GetKineticEnergy(); + auto particlePosition = aFastTrack.GetPrimaryTrack()->GetPosition(); + // Calculate how to create energy deposits + // Following PDG 33.5 chapter + // material calculation assumes homogeneous detector (true for Par03 example) + + always("**********************"); + always("Entered DoIt method"); + auto material = G4Material::GetMaterial(m_material); + if(!material) { + warning("Material " + m_material + " not found. Using the one that the track is currently in."); + material = aFastTrack.GetPrimaryTrack()->GetMaterial(); + } + double materialX0 = material->GetRadlen(); + double materialZ = material->GetZ(); + // EC estimation follows PDG fit to solids in Fig. 33.14 (rms 2.2%) + double materialEc = 610 * CLHEP::MeV / (materialZ + 1.24); + // RM estimation follows PDG Eq. (33.37) (rms 2.2%) + double materialRM = 21.2052 * CLHEP::MeV * materialX0 / materialEc; + double particleY = energy / materialEc; + // Estimate shower maximum and alpha parameter of Gamma distribution + // that describes the longitudinal profile (PDG Eq. (33.35)) + // unless alpha is specified by UI command + if(m_alpha < 0) { + // from PDG Eq. (33.36) + G4double particleTmax = std::log(particleY); + if(aFastTrack.GetPrimaryTrack()->GetParticleDefinition() == G4Gamma::GammaDefinition()) { + particleTmax += 0.5; + } else { + particleTmax -= 0.5; + } + m_alpha = particleTmax * m_beta + 1; + } + // Unless sigma of Gaussian distribution describing the transverse profile + // is specified by UI command, use value calculated from Moliere Radius + if(m_sigma < 0) { + // 90% of shower is contained within 1 * R_M + // 1.645 * std dev of Gaussian contains 90% + m_sigma = materialRM / 1.645; + } + // Calculate rotation matrix along the particle momentum direction + // It will rotate the shower axes to match the incoming particle direction + auto rotMatrix = G4RotationMatrix(); + double particleTheta = particlePosition.getTheta(); + double particlePhi = particlePosition.getPhi(); + double epsilon = 1e-3; + rotMatrix.rotateY(particleTheta); + // do not use (random) phi if x==y==0 + if(!(std::fabs(particlePosition.x()) < epsilon && std::fabs(particlePosition.y()) < epsilon)) { + rotMatrix.rotateZ(particlePhi); + } + // Create hits + // First use rejecton sampling to sample from Gamma distribution + // then get random numbers from uniform distribution for azimuthal angle, and + // from Gaussian for radius + G4ThreeVector position; + double gammaMax = gamma((m_alpha - 1) / m_beta, m_alpha, m_beta); + int generatedHits = 0; + while(generatedHits < m_hitsNo) { + double random1 = G4UniformRand() * m_longMaxDepth; + double random2 = G4UniformRand() * gammaMax; + if(gamma(random1, m_alpha, m_beta) >= random2) { + // Generate corresponding rho (phi) from Gaussian (flat) distribution + double phiPosition = G4UniformRand() * 2 * CLHEP::pi; + double rhoPosition = G4RandGauss::shoot(0, m_sigma); + position = particlePosition + + rotMatrix * G4ThreeVector(rhoPosition * std::sin(phiPosition), + rhoPosition * std::cos(phiPosition), + random1 * materialX0); + // Create energy deposit in the detector + // This will call appropriate sensitive detector class + always("Generating hit x: " + std::to_string(position.x()) + " y: " + std::to_string(position.y()) + " z: " + std::to_string(position.z()) ); + m_sensResponse->simulate(GiGaMTFastHit(position, energy / m_hitsNo), aFastTrack); + generatedHits++; + } + } +} + +class ModelFactory : public GiGaMTFastSimModelFAC<Model> { + + Gaudi::Property<double> m_stepLength {this, "FakeStepLength", 0.0 * CLHEP::mm}; + Gaudi::Property<double> m_sigma {this, "SigmaMoliereRadius", -1}; + Gaudi::Property<double> m_alpha {this, "AlphaGammaDistribution", -1}; + Gaudi::Property<double> m_beta {this, "BetaGammaDistribution", .5}; + Gaudi::Property<int> m_hitsNo {this, "HitsPerShower", 100}; + Gaudi::Property<int> m_longMaxDepth {this, "MaxShowerDepth", 30}; + Gaudi::Property<std::string> m_material {this, "HomogeneousMaterial", std::string()}; + + public: + virtual ~ModelFactory() = default; + using gaussino_base_class = GiGaMTFastSimModelFAC<Model>; + using gaussino_base_class::GiGaMTFastSimModelFAC; + + virtual Model* construct() const override { + auto tmp = gaussino_base_class::construct(); + tmp->SetStepLength(m_stepLength); + tmp->SetAlpha(m_alpha); + tmp->SetBeta(m_beta); + tmp->SetSigma(m_sigma); + tmp->SetHitsNo(m_hitsNo); + tmp->SetLongMaxDepth(m_longMaxDepth); + tmp->SetMaterial(m_material); + return tmp; + }; +}; +} + +DECLARE_COMPONENT_WITH_ID( ShowerDeposit::ModelFactory, "ShowerDepositModelFactory" ) diff --git a/Sim/GiGaMTExamples/FastSimulationExample/src/ShowerDeposit/ShowerDepositModel.h b/Sim/GiGaMTExamples/FastSimulationExample/src/ShowerDeposit/ShowerDepositModel.h new file mode 100644 index 0000000000000000000000000000000000000000..9b41bc25dbfff75014fd64bedc316de5c1c8b3c3 --- /dev/null +++ b/Sim/GiGaMTExamples/FastSimulationExample/src/ShowerDeposit/ShowerDepositModel.h @@ -0,0 +1,72 @@ +/*****************************************************************************\ +* (c) Copyright 2021 CERN for the benefit of the LHCb Collaboration * +* * +* This software is distributed under the terms of the GNU General Public * +* Licence version 3 (GPL Version 3), copied verbatim in the file "COPYING". * +* * +* In applying this licence, CERN does not waive the privileges and immunities * +* granted to it by virtue of its status as an Intergovernmental Organization * +* or submit itself to any jurisdiction. * +\*****************************************************************************/ +#pragma once + +#include "Geant4/G4VFastSimulationModel.hh" +#include "GiGaMTCoreMessage/IGiGaMessage.h" +#include "GiGaMTFastSimulation/GiGaMTFastSimSensitiveResponse.h" + +namespace ShowerDeposit { + +class SensitiveResponse : public GiGaMTFastSimSensitiveResponse { + + double m_stepLength {0.0 * CLHEP::mm}; + + inline virtual void modifyFakeStep(G4Step* fakeStep) override { + // this is a hack, so that the model will work on detectors + // where it the step must be finite + fakeStep->SetStepLength(m_stepLength); + } + + public: + + inline void SetStepLength (double stepLength) { m_stepLength = stepLength; } + inline double GetStepLength () { return m_stepLength; } +}; + +class Model : public G4VFastSimulationModel, public GiGaMessage { + + std::unique_ptr<SensitiveResponse> m_sensResponse; + double m_sigma {-1}; + double m_alpha {-1}; + double m_beta {.5}; + int m_hitsNo {100}; + int m_longMaxDepth {30}; + std::string m_material {std::string()}; + + // Gamma distribution + inline double gamma(double x, double alpha, double beta) { + return (std::pow(beta, alpha) / std::tgamma(alpha) * std::pow(x, alpha - 1) * std::exp(-beta * x)); + } + /// Gaussian distribution + inline double gaussian(double x, double sigma = 1, double x0 = 0) { + double tmp = (x - x0) / sigma; + return (1.0 / (std::sqrt(2 * CLHEP::pi) * sigma)) * std::exp(-tmp * tmp / 2); + } + + +public: + Model(G4String modelName, G4Region* envelope); + ~Model() = default; + + G4bool IsApplicable( const G4ParticleDefinition& aParticle ) override; + G4bool ModelTrigger( const G4FastTrack& aFastTrack ) override; + void DoIt( const G4FastTrack& aFastTrack, G4FastStep& aFastStep ) override; + + inline void SetSigma(double eSigma) { m_sigma = eSigma; }; + inline void SetAlpha(double eAlpha) { m_alpha = eAlpha; }; + inline void SetBeta(double eBeta) { m_beta = eBeta; }; + inline void SetHitsNo(int eHitsNo) { m_hitsNo = eHitsNo; }; + inline void SetLongMaxDepth(double eLongMaxDepth) { m_longMaxDepth = eLongMaxDepth; }; + inline void SetMaterial(std::string eMaterial) { m_material = eMaterial; }; + inline void SetStepLength(double stepLength) { m_sensResponse->SetStepLength(stepLength); }; +}; +} diff --git a/Sim/GiGaMTExamples/FastSimulationExample/tests/options/setup_gamma_gun.py b/Sim/GiGaMTExamples/FastSimulationExample/tests/options/setup_gamma_gun.py new file mode 100644 index 0000000000000000000000000000000000000000..63a3c2e07ec5baeee7e90dd661b8f4cfb280658e --- /dev/null +++ b/Sim/GiGaMTExamples/FastSimulationExample/tests/options/setup_gamma_gun.py @@ -0,0 +1,31 @@ +############################################################################### +# (c) Copyright 2021 CERN for the benefit of the LHCb Collaboration # +# # +# This software is distributed under the terms of the GNU General Public # +# Licence version 3 (GPL Version 3), copied verbatim in the file "COPYING". # +# # +# In applying this licence, CERN does not waive the privileges and immunities # +# granted to it by virtue of its status as an Intergovernmental Organization # +# or submit itself to any jurisdiction. # +############################################################################### +from Gaussino.Generation import GenPhase +GenPhase().ParticleGun = True +GenPhase().ParticleGunUseDefault = False + +from Configurables import ParticleGun +pgun = ParticleGun("ParticleGun") + +from Configurables import FixedMomentum +pgun.ParticleGunTool = "FixedMomentum" +pgun.addTool(FixedMomentum, name="FixedMomentum") +from GaudiKernel.SystemOfUnits import GeV +pgun.FixedMomentum.px = 0. * GeV +pgun.FixedMomentum.py = 0. * GeV +pgun.FixedMomentum.pz = 1. * GeV +pgun.FixedMomentum.PdgCodes = [22] + +from Configurables import FlatNParticles +pgun.NumberOfParticlesTool = "FlatNParticles" +pgun.addTool(FlatNParticles, name="FlatNParticles") +pgun.FlatNParticles.MinNParticles = 1 +pgun.FlatNParticles.MaxNParticles = 1 diff --git a/Sim/GiGaMTExamples/FastSimulationExample/tests/options/setup_immediate_deposit_model.py b/Sim/GiGaMTExamples/FastSimulationExample/tests/options/setup_immediate_deposit_model.py new file mode 100644 index 0000000000000000000000000000000000000000..eb2a60a35a50391edbdd81c1c0b740ee62788f99 --- /dev/null +++ b/Sim/GiGaMTExamples/FastSimulationExample/tests/options/setup_immediate_deposit_model.py @@ -0,0 +1,51 @@ +############################################################################### +# (c) Copyright 2021 CERN for the benefit of the LHCb Collaboration # +# # +# This software is distributed under the terms of the GNU General Public # +# Licence version 3 (GPL Version 3), copied verbatim in the file "COPYING". # +# # +# In applying this licence, CERN does not waive the privileges and immunities # +# granted to it by virtue of its status as an Intergovernmental Organization # +# or submit itself to any jurisdiction. # +############################################################################### +from Configurables import GiGaMT +from Configurables import GiGaMTModularPhysListFAC +giga = GiGaMT() +gmpl = giga.addTool(GiGaMTModularPhysListFAC("ModularPL"), name="ModularPL") +giga.PhysicsListFactory = "GiGaMTModularPhysListFAC/ModularPL" +gmpl = giga.ModularPL +gmpl.PhysicsConstructors = [] + +from Configurables import GiGaMT_G4EmStandardPhysics +gmpl.addTool(GiGaMT_G4EmStandardPhysics(), name='EmPhysics') +gmpl.PhysicsConstructors.append('GiGaMT_G4EmStandardPhysics/EmPhysics') + +from Configurables import GiGaMTFastSimPhysFAC +gmpl.addTool(GiGaMTFastSimPhysFAC(), name='FastSimPhys') +gmpl.PhysicsConstructors.append('GiGaMTFastSimPhysFAC/FastSimPhys') + +from Configurables import GiGaMTDetectorConstructionFAC +giga.DetectorConstruction = "GiGaMTDetectorConstructionFAC" +dettool = giga.addTool(GiGaMTDetectorConstructionFAC, + "GiGaMTDetectorConstructionFAC") +dettool.OutputLevel = -10 + +from Configurables import SimpleCaloGeo +SimpleCaloGeo().OutputLevel = 10 +dettool.GiGaMTGeoSvc = "SimpleCaloGeo" + +from Configurables import ApplicationMgr +ApplicationMgr().ExtSvc += [SimpleCaloGeo()] + +from Configurables import FastG4RegionFAC +fast_region = FastG4RegionFAC() +fast_region.DetName="Absorber" +dettool.addTool(fast_region, name='FastRegionAbsorber') +dettool.FastRegionFactories.append("FastG4RegionFAC/FastRegionAbsorber") + +from Configurables import ImmediateDepositModelFactory +fast_model = ImmediateDepositModelFactory() +fast_model.ModelName="ImmediateDepositModel" +fast_model.RegionName="AbsorberSDetFastRegion" +dettool.addTool(fast_model, name='FastModelAbsorber') +dettool.FastModelFactories.append('ImmediateDepositModelFactory/FastModelAbsorber') diff --git a/Sim/GiGaMTExamples/FastSimulationExample/tests/qmtest/fast_simulation.qmt b/Sim/GiGaMTExamples/FastSimulationExample/tests/qmtest/fast_simulation.qmt new file mode 100644 index 0000000000000000000000000000000000000000..6f94a0cf16cfca1228ac17b8204c7628d60855eb --- /dev/null +++ b/Sim/GiGaMTExamples/FastSimulationExample/tests/qmtest/fast_simulation.qmt @@ -0,0 +1,40 @@ +<?xml version="1.0" encoding="UTF-8"?><!DOCTYPE extension PUBLIC '-//QM/2.3/Extension//EN' 'http://www.codesourcery.com/qm/dtds/2.3/-//qm/2.3/extension//en.dtd'> +<!-- + (c) Copyright 2021 CERN for the benefit of the LHCb Collaboration + + This software is distributed under the terms of the GNU General Public + Licence version 3 (GPL Version 3), copied verbatim in the file "COPYING". + + In applying this licence, CERN does not waive the privileges and immunities + granted to it by virtue of its status as an Intergovernmental Organization + or submit itself to any jurisdiction. +--> + +<extension class="GaudiTest.GaudiExeTest" kind="test"> +<argument name="program"><text>gaudirun.py</text></argument> +<argument name="args"><set> + <text>$SIMPLECALOGEOROOT/options/setup_simple_calo.py</text> + <text>$FASTSIMULATIONEXAMPLETESTS/options/setup_gamma_gun.py</text> + <text>$FASTSIMULATIONEXAMPLETESTS/options/setup_immediate_deposit_model.py</text> +</set></argument> +<argument name="timeout"><integer>60</integer></argument> +<argument name="use_temp_dir"><enumeral>true</enumeral></argument> +<argument name="exit_code"><integer>0</integer></argument> +<argument name="environment"><set/></argument> +<argument name="unsupported_platforms"><set/></argument> +<argument name="workdir"><text/></argument> +<argument name="stderr"><text/></argument> +<argument name="options"><text></text></argument> +<argument name="validator"><text> +results = ["SimpleCaloHits total energy: 100000", + "SimpleCaloHits energy mean: 1000", + "SimpleCaloHits energy mean error: 0"] +for result in results: + if stdout.find(result) == -1: + causes.append('missing string') + result['GaudiTest.expected_string'] = result.Quote(result) +</text></argument> +<argument name="resources"><set/></argument> +<argument name="stdin"><text/></argument> +</extension> + diff --git a/Sim/GiGaMTExamples/SimpleCaloGeo/CMakeLists.txt b/Sim/GiGaMTExamples/SimpleCaloGeo/CMakeLists.txt new file mode 100644 index 0000000000000000000000000000000000000000..bb4949ff310450723b91cb9fd17249df0a52a623 --- /dev/null +++ b/Sim/GiGaMTExamples/SimpleCaloGeo/CMakeLists.txt @@ -0,0 +1,23 @@ +################################################################################ +# Package: GiGaMTExamples/SimpleCaloGeo +# +# Package that sets a simplified version of electromagnetic calorimeter. +# +################################################################################# +gaudi_subdir(SimpleCaloGeo v1r0) + +gaudi_depends_on_subdirs(GaudiAlg + Sim/GiGaMTFactories + Sim/GiGaMTCore) +AddHepMC3() +if(${Geant4_config_version} VERSION_LESS "10.06") + add_definitions(-DG4MULTITHREADED) + add_definitions(-DG4USE_STD11) +endif() + +gaudi_add_module(SimpleCaloGeo + src/*.cpp + INCLUDE_DIRS GiGaMTFactories + LINK_LIBRARIES GaudiAlgLib GiGaMTCoreRunLib) + +add_dependencies(SimpleCaloGeo HepMC3Ext) diff --git a/Sim/GiGaMTExamples/SimpleCaloGeo/options/setup_simple_calo.py b/Sim/GiGaMTExamples/SimpleCaloGeo/options/setup_simple_calo.py new file mode 100644 index 0000000000000000000000000000000000000000..8e8f5aa85d652d5f7f6a8fa5ab2932cc6e8d3d35 --- /dev/null +++ b/Sim/GiGaMTExamples/SimpleCaloGeo/options/setup_simple_calo.py @@ -0,0 +1,24 @@ +############################################################################### +# (c) Copyright 2021 CERN for the benefit of the LHCb Collaboration # +# # +# This software is distributed under the terms of the GNU General Public # +# Licence version 3 (GPL Version 3), copied verbatim in the file "COPYING". # +# # +# In applying this licence, CERN does not waive the privileges and immunities # +# granted to it by virtue of its status as an Intergovernmental Organization # +# or submit itself to any jurisdiction. # +############################################################################### +from multiprocessing import cpu_count +from Configurables import GiGaMT +giga = GiGaMT() +giga.NumberOfWorkerThreads = cpu_count() + +from Configurables import Gaussino +Gaussino().EvtMax = 100 +Gaussino().EnableHive = True +Gaussino().ThreadPoolSize = 20 +Gaussino().EventSlots = 20 + +from Configurables import SimpleCaloSaveHits +monitool = giga.addTool(SimpleCaloSaveHits()) +giga.MonitorTools = ["SimpleCaloSaveHits"] diff --git a/Sim/GiGaMTExamples/SimpleCaloGeo/src/SimpleCaloGeo.cpp b/Sim/GiGaMTExamples/SimpleCaloGeo/src/SimpleCaloGeo.cpp new file mode 100755 index 0000000000000000000000000000000000000000..e368e105d53faf2d77204df84311514dcd54b6d1 --- /dev/null +++ b/Sim/GiGaMTExamples/SimpleCaloGeo/src/SimpleCaloGeo.cpp @@ -0,0 +1,243 @@ +/*****************************************************************************\ +* (c) Copyright 2021 CERN for the benefit of the LHCb Collaboration * +* * +* This software is distributed under the terms of the GNU General Public * +* Licence version 3 (GPL Version 3), copied verbatim in the file "COPYING". * +* * +* In applying this licence, CERN does not waive the privileges and immunities * +* granted to it by virtue of its status as an Intergovernmental Organization * +* or submit itself to any jurisdiction. * +\*****************************************************************************/ +#include "Geant4/G4AutoDelete.hh" +#include "Geant4/G4Box.hh" +#include "Geant4/G4Colour.hh" +#include "Geant4/G4GlobalMagFieldMessenger.hh" +#include "Geant4/G4LogicalVolume.hh" +#include "Geant4/G4MultiFunctionalDetector.hh" +#include "Geant4/G4NistManager.hh" +#include "Geant4/G4PSEnergyDeposit.hh" +#include "Geant4/G4PSTrackLength.hh" +#include "Geant4/G4PVPlacement.hh" +#include "Geant4/G4PVReplica.hh" +#include "Geant4/G4PhysicalConstants.hh" +#include "Geant4/G4SDChargedFilter.hh" +#include "Geant4/G4SDManager.hh" +#include "Geant4/G4SystemOfUnits.hh" +#include "Geant4/G4VPrimitiveScorer.hh" +#include "Geant4/G4VisAttributes.hh" + +#include "SimpleCaloGeo.h" +#include "SimpleCaloSD.h" + +DECLARE_COMPONENT( SimpleCaloGeo ) + +StatusCode SimpleCaloGeo::queryInterface( const InterfaceID& id, void** ppI ) { + if ( 0 == ppI ) { + return StatusCode::FAILURE; // RETURN !!! + } else if ( IGiGaMTGeoSvc::interfaceID() == id ) { + *ppI = static_cast<IGiGaMTGeoSvc*>( this ); + } else { + return Service::queryInterface( id, ppI ); // RETURN !!! + } + + addRef(); + + return StatusCode::SUCCESS; +} + +void SimpleCaloGeo::defineMaterials() { + // Lead material defined using NIST Manager + auto nistManager = G4NistManager::Instance(); + nistManager->FindOrBuildMaterial( "G4_Pb" ); + + // Liquid argon material + G4double a; // mass of a mole; + G4double z; // z=mean number of protons; + G4double density; + new G4Material( "liquidArgon", z = 18., a = 39.95 * g / mole, density = 1.390 * g / cm3 ); + // The argon by NIST Manager is a gas with a different density + // Vacuum + new G4Material( "Galactic", z = 1., a = 1.01 * g / mole, density = universe_mean_density, kStateGas, 2.73 * kelvin, + 3.e-18 * pascal ); + // Print materials + G4cout << *( G4Material::GetMaterialTable() ) << G4endl; +} + +G4VPhysicalVolume* SimpleCaloGeo::constructWorld() { + // Lead material defined using NIST Manager + auto nistManager = G4NistManager::Instance(); + nistManager->FindOrBuildMaterial( "G4_Pb" ); + + // Liquid argon material + G4double a; // mass of a mole; + G4double z; // z=mean number of protons; + G4double density; + new G4Material( "liquidArgon", z = 18., a = 39.95 * g / mole, density = 1.390 * g / cm3 ); + // The argon by NIST Manager is a gas with a different density + + // Vacuum + new G4Material( "Galactic", z = 1., a = 1.01 * g / mole, density = universe_mean_density, kStateGas, 2.73 * kelvin, + 3.e-18 * pascal ); + + // Print materials + G4cout << *( G4Material::GetMaterialTable() ) << G4endl; + + // Geometry parameters + G4double absoThickness = 10. * mm; + G4double gapThickness = 5. * mm; + G4double calorSizeXY = 10. * cm; + auto layerThickness = absoThickness + gapThickness; + auto calorThickness = m_nofLayers * layerThickness; + auto worldSizeXY = 1.2 * calorSizeXY; + auto worldSizeZ = 1.2 * calorThickness; + + // Get materials + auto defaultMaterial = G4Material::GetMaterial( "Galactic" ); + auto absorberMaterial = G4Material::GetMaterial( "G4_Pb" ); + auto gapMaterial = G4Material::GetMaterial( "liquidArgon" ); + + if ( !defaultMaterial || !absorberMaterial || !gapMaterial ) { + G4ExceptionDescription msg; + msg << "Cannot retrieve materials already defined."; + G4Exception( "B4DetectorConstruction::DefineVolumes()", "MyCode0001", FatalException, msg ); + } + + // + // World + // + auto worldS = new G4Box( "World", // its name + worldSizeXY / 2, worldSizeXY / 2, worldSizeZ / 2 ); // its size + + auto worldLV = new G4LogicalVolume( worldS, // its solid + defaultMaterial, // its material + "World" ); // its name + + auto worldPV = new G4PVPlacement( 0, + G4ThreeVector(), // at (0,0,0) + worldLV, // its logical volume + "World", // its name + 0, // its mother volume + false, // no boolean operation + 0, // copy number + m_checkOverlaps ); // checking overlaps + + // + // Calorimeter + // + auto calorimeterS = new G4Box( "Calorimeter", // its name + calorSizeXY / 2, calorSizeXY / 2, calorThickness / 2 ); // its size + + auto calorLV = new G4LogicalVolume( calorimeterS, // its solid + defaultMaterial, // its material + "Calorimeter" ); // its name + + new G4PVPlacement( 0, // no rotation + G4ThreeVector(), // at (0,0,0) + calorLV, // its logical volume + "Calorimeter", // its name + worldLV, // its mother volume + false, // no boolean operation + 0, // copy number + m_checkOverlaps ); // checking overlaps + + // + // Layer + // + auto layerS = new G4Box( "Layer", // its name + calorSizeXY / 2, calorSizeXY / 2, layerThickness / 2 ); // its size + + auto layerLV = new G4LogicalVolume( layerS, // its solid + defaultMaterial, // its material + "Layer" ); // its name + new G4PVReplica( "Layer", // its name + layerLV, // its logical volume + calorLV, // its mother + kZAxis, // axis of replication + m_nofLayers, // number of replica + layerThickness ); // witdth of replica + + // + // Absorber + // + auto absorberS = new G4Box( "Abso", // its name + calorSizeXY / 2, calorSizeXY / 2, absoThickness / 2 ); // its size + + m_absorberLV = new G4LogicalVolume( absorberS, // its solid + absorberMaterial, // its material + "AbsoLV" ); // its name + + new G4PVPlacement( 0, // no rotation + G4ThreeVector( 0., 0., -gapThickness / 2 ), // its position + m_absorberLV, // its logical volume + "Abso", // its name + layerLV, // its mother volume + false, // no boolean operation + 0, // copy number + m_checkOverlaps ); // checking overlaps + // + // Gap + // + auto gapS = new G4Box( "Gap", // its name + calorSizeXY / 2, calorSizeXY / 2, gapThickness / 2 ); // its size + + m_gapLV = new G4LogicalVolume( gapS, // its solid + gapMaterial, // its material + "GapLV" ); // its name + + new G4PVPlacement( 0, // no rotation + G4ThreeVector( 0., 0., absoThickness / 2 ), // its position + m_gapLV, // its logical volume + "Gap", // its name + layerLV, // its mother volume + false, // no boolean operation + 0, // copy number + m_checkOverlaps ); // checking overlaps + + // + // print parameters + // + G4cout << G4endl << "------------------------------------------------------------" << G4endl + << "---> The calorimeter is " << m_nofLayers << " layers of: [ " << absoThickness / mm << "mm of " + << absorberMaterial->GetName() << " + " << gapThickness / mm << "mm of " << gapMaterial->GetName() << " ] " + << G4endl << "------------------------------------------------------------" << G4endl; + + // + // Visualization attributes + // + worldLV->SetVisAttributes( G4VisAttributes::GetInvisible() ); + auto simpleBoxVisAtt = new G4VisAttributes( G4Colour( 1.0, 1.0, 1.0 ) ); + simpleBoxVisAtt->SetVisibility( true ); + calorLV->SetVisAttributes( simpleBoxVisAtt ); + // + // Always return the physical World + // + return worldPV; +} + +void SimpleCaloGeo::constructSDandField() { + // G4SDManager::GetSDMpointer()->SetVerboseLevel(1); + // + // Sensitive detectors + // + + auto absoSD = m_absorberFAC->construct(); + G4SDManager::GetSDMpointer()->AddNewDetector( absoSD ); + m_absorberLV->SetSensitiveDetector( absoSD ); + + auto gapSD = m_gapFAC->construct(); + G4SDManager::GetSDMpointer()->AddNewDetector( gapSD ); + m_gapLV->SetSensitiveDetector( gapSD ); + + // + // Magnetic field + // + // Create global magnetic field messenger. + // Uniform magnetic field is then created automatically if + // the field value is not zero. + // G4ThreeVector fieldValue; + // fMagFieldMessenger = new G4GlobalMagFieldMessenger(fieldValue); + // fMagFieldMessenger->SetVerboseLevel(1); + + // Register the field messenger for deleting + // G4AutoDelete::Register(fMagFieldMessenger); +} diff --git a/Sim/GiGaMTExamples/SimpleCaloGeo/src/SimpleCaloGeo.h b/Sim/GiGaMTExamples/SimpleCaloGeo/src/SimpleCaloGeo.h new file mode 100755 index 0000000000000000000000000000000000000000..8a0e0874aa9a678ad7834cd94d2b4a6f44fe50c3 --- /dev/null +++ b/Sim/GiGaMTExamples/SimpleCaloGeo/src/SimpleCaloGeo.h @@ -0,0 +1,45 @@ +/*****************************************************************************\ +* (c) Copyright 2021 CERN for the benefit of the LHCb Collaboration * +* * +* This software is distributed under the terms of the GNU General Public * +* Licence version 3 (GPL Version 3), copied verbatim in the file "COPYING". * +* * +* In applying this licence, CERN does not waive the privileges and immunities * +* granted to it by virtue of its status as an Intergovernmental Organization * +* or submit itself to any jurisdiction. * +\*****************************************************************************/ +#pragma once + +#include "GaudiKernel/Kernel.h" +#include "GaudiKernel/Service.h" +#include "GaudiKernel/StatusCode.h" + +#include "GiGaMTGeo/IGiGaMTGeoSvc.h" +#include "SimpleCaloSD.h" + +class G4VPhysicalVolume; +class G4LogicalVolume; + +class SimpleCaloGeo : public Service, virtual public IGiGaMTGeoSvc { + G4LogicalVolume* m_gapLV; + G4LogicalVolume* m_absorberLV; + + Gaudi::Property<bool> m_checkOverlaps{this, "CheckOverlaps", true}; + Gaudi::Property<int> m_nofLayers{this, "LayersNumber", 10}; + + ToolHandle<SimpleCaloSDFAC> m_absorberFAC{"AbsorberSDet", this}; + ToolHandle<SimpleCaloSDFAC> m_gapFAC{"GapSDet", this}; + + void defineMaterials(); + +protected: + using Service::Service; + +public: + StatusCode initialize() override { return Service::initialize(); } + StatusCode finalize() override { return Service::finalize(); } + + virtual G4VPhysicalVolume* constructWorld() override; + virtual void constructSDandField() override; + virtual StatusCode queryInterface( const InterfaceID& iid, void** pI ) override; +}; diff --git a/Sim/GiGaMTExamples/SimpleCaloGeo/src/SimpleCaloHit.cpp b/Sim/GiGaMTExamples/SimpleCaloGeo/src/SimpleCaloHit.cpp new file mode 100644 index 0000000000000000000000000000000000000000..d25466d108c13069224813e6777d4deb33b163e1 --- /dev/null +++ b/Sim/GiGaMTExamples/SimpleCaloGeo/src/SimpleCaloHit.cpp @@ -0,0 +1,15 @@ +/*****************************************************************************\ +* (c) Copyright 2021 CERN for the benefit of the LHCb Collaboration * +* * +* This software is distributed under the terms of the GNU General Public * +* Licence version 3 (GPL Version 3), copied verbatim in the file "COPYING". * +* * +* In applying this licence, CERN does not waive the privileges and immunities * +* granted to it by virtue of its status as an Intergovernmental Organization * +* or submit itself to any jurisdiction. * +\*****************************************************************************/ + +#include "SimpleCaloHit.h" + +G4ThreadLocal G4Allocator<SimpleCaloHit>* SimpleCaloHitAllocator = 0; + diff --git a/Sim/GiGaMTExamples/SimpleCaloGeo/src/SimpleCaloHit.h b/Sim/GiGaMTExamples/SimpleCaloGeo/src/SimpleCaloHit.h new file mode 100644 index 0000000000000000000000000000000000000000..6dff6931e848300f1e8ec18ff2f767b54aa32ba0 --- /dev/null +++ b/Sim/GiGaMTExamples/SimpleCaloGeo/src/SimpleCaloHit.h @@ -0,0 +1,49 @@ +/*****************************************************************************\ +* (c) Copyright 2021 CERN for the benefit of the LHCb Collaboration * +* * +* This software is distributed under the terms of the GNU General Public * +* Licence version 3 (GPL Version 3), copied verbatim in the file "COPYING". * +* * +* In applying this licence, CERN does not waive the privileges and immunities * +* granted to it by virtue of its status as an Intergovernmental Organization * +* or submit itself to any jurisdiction. * +\*****************************************************************************/ +#pragma once + +#include "Geant4/G4Allocator.hh" +#include "Geant4/G4THitsCollection.hh" +#include "Geant4/G4Threading.hh" +#include "Geant4/G4ThreeVector.hh" +#include "Geant4/G4VHit.hh" + + +class SimpleCaloHit : public G4VHit { + + double m_energy; + +public: + SimpleCaloHit() = default; + inline SimpleCaloHit( const SimpleCaloHit& ) : G4VHit(), m_energy( 0. ) {}; + virtual ~SimpleCaloHit() = default; + + inline void* operator new( size_t ); + inline void operator delete( void* ); + + inline void Add( double de ) { m_energy += de; } + + inline double GetEnergy() const { return m_energy; } +}; + +using SimpleCaloHitsCollection = G4THitsCollection<SimpleCaloHit>; + +extern G4ThreadLocal G4Allocator<SimpleCaloHit>* SimpleCaloHitAllocator; + +inline void* SimpleCaloHit::operator new( size_t ) { + if ( !SimpleCaloHitAllocator ) { SimpleCaloHitAllocator = new G4Allocator<SimpleCaloHit>; } + return (void*)SimpleCaloHitAllocator->MallocSingle(); +} + +inline void SimpleCaloHit::operator delete( void* hit ) { + if ( !SimpleCaloHitAllocator ) { SimpleCaloHitAllocator = new G4Allocator<SimpleCaloHit>; } + SimpleCaloHitAllocator->FreeSingle( (SimpleCaloHit*)hit ); +} diff --git a/Sim/GiGaMTExamples/SimpleCaloGeo/src/SimpleCaloSD.cpp b/Sim/GiGaMTExamples/SimpleCaloGeo/src/SimpleCaloSD.cpp new file mode 100644 index 0000000000000000000000000000000000000000..6d6bc90d78c82e8888a351a651e480decb137875 --- /dev/null +++ b/Sim/GiGaMTExamples/SimpleCaloGeo/src/SimpleCaloSD.cpp @@ -0,0 +1,64 @@ +/*****************************************************************************\ +* (c) Copyright 2021 CERN for the benefit of the LHCb Collaboration * +* * +* This software is distributed under the terms of the GNU General Public * +* Licence version 3 (GPL Version 3), copied verbatim in the file "COPYING". * +* * +* In applying this licence, CERN does not waive the privileges and immunities * +* granted to it by virtue of its status as an Intergovernmental Organization * +* or submit itself to any jurisdiction. * +\*****************************************************************************/ +#include "SimpleCaloSD.h" + +#include "Geant4/G4HCofThisEvent.hh" +#include "Geant4/G4SDManager.hh" +#include "Geant4/G4Step.hh" + +SimpleCaloSD::SimpleCaloSD( const G4String& name ) : G4VSensitiveDetector( name ) { + collectionName.insert( SensitiveDetectorName + "HitsCollection" ); +} + +void SimpleCaloSD::Initialize( G4HCofThisEvent* hce ) { + m_hitsCollection = new SimpleCaloHitsCollection( SensitiveDetectorName, collectionName[0] ); + auto hcID = G4SDManager::GetSDMpointer()->GetCollectionID( collectionName[0] ); + hce->AddHitsCollection( hcID, m_hitsCollection ); + for ( G4int i = 0; i < m_cellsNo; i++ ) { m_hitsCollection->insert( new SimpleCaloHit() ); } +} + +G4bool SimpleCaloSD::Hit( G4Step* aStep ) { + G4TouchableHistory* ROhis = 0; + if ( !isActive() ) return false; + if ( filter ) { + if ( !( filter->Accept( aStep ) ) ) return false; + } + if ( ROgeometry ) { + if ( !( ROgeometry->CheckROVolume( aStep, ROhis ) ) ) return false; + } + return ProcessHits( aStep, ROhis ); +} + +bool SimpleCaloSD::ProcessHits( G4Step* step, G4TouchableHistory* ) { + auto edep = step->GetTotalEnergyDeposit(); + if ( edep == 0. ) return false; + auto touchable = ( step->GetPreStepPoint()->GetTouchable() ); + auto layerNumber = touchable->GetReplicaNumber( 1 ); + auto hit = ( *m_hitsCollection )[layerNumber]; + if ( !hit ) { + G4ExceptionDescription msg; + msg << "Cannot access hit " << layerNumber; + G4Exception( "SimpleCaloSD::ProcessHits()", "MyCode0004", FatalException, msg ); + } + + hit->Add( edep ); + return true; +} + +void SimpleCaloSD::EndOfEvent( G4HCofThisEvent* ) { + if ( !m_hitsCollection ) { + warning( "There are no hits in the hits collection." ); + return; + } +} + +DECLARE_COMPONENT_WITH_ID( SimpleCaloSDFAC, "GapSDet" ) +DECLARE_COMPONENT_WITH_ID( SimpleCaloSDFAC, "AbsorberSDet" ) diff --git a/Sim/GiGaMTExamples/SimpleCaloGeo/src/SimpleCaloSD.h b/Sim/GiGaMTExamples/SimpleCaloGeo/src/SimpleCaloSD.h new file mode 100644 index 0000000000000000000000000000000000000000..349ebdbf578dbdc6a0ee58ad6b806dd935d33aaa --- /dev/null +++ b/Sim/GiGaMTExamples/SimpleCaloGeo/src/SimpleCaloSD.h @@ -0,0 +1,58 @@ +/*****************************************************************************\ +* (c) Copyright 2021 CERN for the benefit of the LHCb Collaboration * +* * +* This software is distributed under the terms of the GNU General Public * +* Licence version 3 (GPL Version 3), copied verbatim in the file "COPYING". * +* * +* In applying this licence, CERN does not waive the privileges and immunities * +* granted to it by virtue of its status as an Intergovernmental Organization * +* or submit itself to any jurisdiction. * +\*****************************************************************************/ +#pragma once + +#include "GiGaMTCoreMessage/IGiGaMessage.h" +#include "GiGaMTFactories/GiGaMTG4SensDetFactory.h" + +#include "Geant4/G4VSensitiveDetector.hh" + +#include "SimpleCaloHit.h" + +class G4Step; +class G4HCofThisEvent; + +class SimpleCaloSD : public G4VSensitiveDetector, public virtual GiGaMessage { +public: + SimpleCaloSD( const G4String& name ); + virtual ~SimpleCaloSD() = default; + + // methods from base class + void Initialize( G4HCofThisEvent* hitCollection ) override; + bool ProcessHits( G4Step* step, G4TouchableHistory* history ) override; + void EndOfEvent( G4HCofThisEvent* hitCollection ) override; + G4bool Hit( G4Step* aStep ); + + SimpleCaloHitsCollection* m_hitsCollection = nullptr; + int m_cellsNo; +}; + +class SimpleCaloSDFAC : public GiGaMTG4SensDetFactory<SimpleCaloSD> { + + // TODO: this is the same property as in SimpleCaloGeo + Gaudi::Property<int> m_cellsNo{this, "CellsNo", 10}; + +public: + using base_class = GiGaMTG4SensDetFactory<SimpleCaloSD>; + using base_class::GiGaMTG4SensDetFactory; + + virtual ~SimpleCaloSDFAC() = default; + virtual SimpleCaloSD* construct() const override { + auto tmp = base_class::construct(); + tmp->SetVerboseLevel( MSG::DEBUG ); + tmp->SetMessageInterface( message_interface() ); + tmp->m_cellsNo = m_cellsNo; + // tmp->m_nhits = &m_nhits; + // tmp->m_energy = &m_energy; + return tmp; + } +}; + diff --git a/Sim/GiGaMTExamples/SimpleCaloGeo/src/SimpleCaloSaveHits.cpp b/Sim/GiGaMTExamples/SimpleCaloGeo/src/SimpleCaloSaveHits.cpp new file mode 100644 index 0000000000000000000000000000000000000000..1b8340f389f631309ed8c3da1066e6681aba873c --- /dev/null +++ b/Sim/GiGaMTExamples/SimpleCaloGeo/src/SimpleCaloSaveHits.cpp @@ -0,0 +1,45 @@ +/*****************************************************************************\ +* (c) Copyright 2021 CERN for the benefit of the LHCb Collaboration * +* * +* This software is distributed under the terms of the GNU General Public * +* Licence version 3 (GPL Version 3), copied verbatim in the file "COPYING". * +* * +* In applying this licence, CERN does not waive the privileges and immunities * +* granted to it by virtue of its status as an Intergovernmental Organization * +* or submit itself to any jurisdiction. * +\*****************************************************************************/ +// Immediate +#include "SimpleCaloSaveHits.h" +#include "SimpleCaloHit.h" +// Gaudi +#include "GaudiKernel/ITHistSvc.h" +// Geant4 +#include "Geant4/G4Event.hh" + +DECLARE_COMPONENT( SimpleCaloSaveHits ) + +StatusCode SimpleCaloSaveHits::finalize() { + if ( msgLevel( MSG::INFO ) ) { + info() << "SimpleCaloHits total energy: " << m_energy.sum() << endmsg; + info() << "SimpleCaloHits energy mean: " << m_energy.mean() << endmsg; + info() << "SimpleCaloHits energy mean error: " << m_energy.meanErr() << endmsg; + } + return extends::finalize(); +} + +StatusCode SimpleCaloSaveHits::monitor( const G4Event& aEvent ) { + G4HCofThisEvent* collections = aEvent.GetHCofThisEvent(); + + for ( int iter_coll = 0; iter_coll < collections->GetNumberOfCollections(); iter_coll++ ) { + G4VHitsCollection* collection = collections->GetHC( iter_coll ); + int n_hit = collection->GetSize(); + if ( msgLevel( MSG::DEBUG ) ) { + debug() << "Spotted #" << n_hit << " hits stored in " << collection->GetName() << " collection." << endmsg; + } + for ( int iter_hit = 0; iter_hit < n_hit; iter_hit++ ) { + SimpleCaloHit* hit = dynamic_cast<SimpleCaloHit*>( collection->GetHit( iter_hit ) ); + if ( hit->GetEnergy() > 0. ) m_energy += hit->GetEnergy(); + } + } + return StatusCode::SUCCESS; +} diff --git a/Sim/GiGaMTExamples/SimpleCaloGeo/src/SimpleCaloSaveHits.h b/Sim/GiGaMTExamples/SimpleCaloGeo/src/SimpleCaloSaveHits.h new file mode 100644 index 0000000000000000000000000000000000000000..e73bf8b6e96d0396afa7228c4d17253d5c5191e2 --- /dev/null +++ b/Sim/GiGaMTExamples/SimpleCaloGeo/src/SimpleCaloSaveHits.h @@ -0,0 +1,25 @@ +/*****************************************************************************\ +* (c) Copyright 2021 CERN for the benefit of the LHCb Collaboration * +* * +* This software is distributed under the terms of the GNU General Public * +* Licence version 3 (GPL Version 3), copied verbatim in the file "COPYING". * +* * +* In applying this licence, CERN does not waive the privileges and immunities * +* granted to it by virtue of its status as an Intergovernmental Organization * +* or submit itself to any jurisdiction. * +\*****************************************************************************/ +#pragma once +// Gaudi +#include "GaudiAlg/GaudiTool.h" +#include "SimInterfaces/IG4MonitoringTool.h" + +class SimpleCaloSaveHits : public extends<GaudiTool, IG4MonitoringTool> { + + mutable Gaudi::Accumulators::StatCounter<double> m_energy{this, "energy"}; + +public: + using extends::extends; + + virtual StatusCode finalize() override; + virtual StatusCode monitor( const G4Event& aEvent ) override; +}; diff --git a/Sim/GiGaMTFactories/src/det/GiGaMTDetectorConstructionFAC.cpp b/Sim/GiGaMTFactories/src/det/GiGaMTDetectorConstructionFAC.cpp index 1d4d6cf3d72508aa1c3e00b686eaf12d5acc15b0..add270b42f1e1a665fb2acd20d190ecc42c1ef34 100644 --- a/Sim/GiGaMTFactories/src/det/GiGaMTDetectorConstructionFAC.cpp +++ b/Sim/GiGaMTFactories/src/det/GiGaMTDetectorConstructionFAC.cpp @@ -1,3 +1,13 @@ +/*****************************************************************************\ +* (c) Copyright 2021 CERN for the benefit of the LHCb Collaboration * +* * +* This software is distributed under the terms of the GNU General Public * +* Licence version 3 (GPL Version 3), copied verbatim in the file "COPYING". * +* * +* In applying this licence, CERN does not waive the privileges and immunities * +* granted to it by virtue of its status as an Intergovernmental Organization * +* or submit itself to any jurisdiction. * +\*****************************************************************************/ #include "GiGaMTDetectorConstructionFAC.h" #include "GiGaMTCoreDet/GiGaMTDetectorConstruction.h" #include "GiGaMTGeo/IGiGaMTGeoSvc.h" @@ -10,6 +20,9 @@ StatusCode GiGaMTDetectorConstructionFAC::initialize() { // Retrieve the factory tools here to avoid the retrieval happening in multiple // threads for ( auto& keypairs : m_sens_dets ) { sc &= keypairs.second.retrieve(); } + for ( auto& fac : m_fast_region_factories ) { sc &= fac.retrieve(); } + for ( auto& fac : m_fast_model_factories ) { sc &= fac.retrieve(); } + return sc; } @@ -26,6 +39,8 @@ G4VUserDetectorConstruction* GiGaMTDetectorConstructionFAC::construct() const { debug() << "Calling SD and Field constructor" << endmsg; m_geoSvc->constructSDandField(); DressVolumes(); + BuildFastRegions(); + BuildFastModels(); } ); return detconst; @@ -52,6 +67,21 @@ void GiGaMTDetectorConstructionFAC::DressVolumes() const { } } +void GiGaMTDetectorConstructionFAC::BuildFastRegions() const { + for ( auto& fast_region_factory : m_fast_region_factories ) { + if ( msgLevel( MSG::DEBUG ) ) { debug() << "Running fast region constructor " << fast_region_factory->name() << endmsg; } + fast_region_factory->construct(); + } +} + +void GiGaMTDetectorConstructionFAC::BuildFastModels() const { + for ( auto& fast_model_factory : m_fast_model_factories ) { + if ( msgLevel( MSG::DEBUG ) ) { debug() << "Running fast model constructor " << fast_model_factory->name() << endmsg; } + fast_model_factory->construct(); + } +} + + #include "Geant4/G4GDMLParser.hh" void GiGaMTDetectorConstructionFAC::SaveGDML( G4LogicalVolume* world ) const { diff --git a/Sim/GiGaMTFactories/src/det/GiGaMTDetectorConstructionFAC.h b/Sim/GiGaMTFactories/src/det/GiGaMTDetectorConstructionFAC.h index faab00ad35c50a7b9e544c5d2a48b0bb4771b249..72e2a38f9dba1bf9468788c197bccd0597fd5d20 100644 --- a/Sim/GiGaMTFactories/src/det/GiGaMTDetectorConstructionFAC.h +++ b/Sim/GiGaMTFactories/src/det/GiGaMTDetectorConstructionFAC.h @@ -1,3 +1,13 @@ +/*****************************************************************************\ +* (c) Copyright 2021 CERN for the benefit of the LHCb Collaboration * +* * +* This software is distributed under the terms of the GNU General Public * +* Licence version 3 (GPL Version 3), copied verbatim in the file "COPYING". * +* * +* In applying this licence, CERN does not waive the privileges and immunities * +* granted to it by virtue of its status as an Intergovernmental Organization * +* or submit itself to any jurisdiction. * +\*****************************************************************************/ #include "Geant4/G4VUserDetectorConstruction.hh" #include "GiGaMTFactories/GiGaFactoryBase.h" #include "GiGaMTFactories/GiGaTool.h" @@ -5,6 +15,7 @@ #include "GaudiAlg/FunctionalDetails.h" #include "GaudiAlg/FunctionalUtilities.h" #include "Utils/ToolProperty.h" +#include "Geant4/G4VFastSimulationModel.hh" class IGiGaMTGeoSvc; class IGaussinoTool; @@ -27,7 +38,13 @@ protected: typedef ToolHandle<GiGaFactoryBase<G4VSensitiveDetector>> SensDetFac; typedef std::map<std::string, SensDetFac> SensDetVolumeMap; + using FastModelFactory = GiGaFactoryBase<G4VFastSimulationModel>; + using FastRegionFactory = GiGaFactoryBase<G4Region>; + void DressVolumes() const; + void BuildFastRegions() const; + void BuildFastModels() const; + void SaveGDML(G4LogicalVolume*) const; ServiceHandle<IGiGaMTGeoSvc> m_geoSvc{this, "GiGaMTGeoSvc", "GiGaMTGeo"}; ToolHandleArray<IGaussinoTool> m_afterGeo{this}; @@ -36,8 +53,18 @@ protected: tool_array_setter(m_afterGeo, m_afterGeoNames), Gaudi::Details::Property::ImmediatelyInvokeHandler{true}}; - Gaudi::Property<std::string> m_schema{this, "Schema", "$GDML_base/src/GDMLSchema/gdml.xsd"}; - Gaudi::Property<std::string> m_outfile{this, "Output", ""}; + Gaudi::Property<std::string> m_schema{this, "Schema", "$GDML_base/src/GDMLSchema/gdml.xsd"}; + Gaudi::Property<std::string> m_outfile{this, "Output", ""}; + ToolHandleArray<FastModelFactory> m_fast_model_factories{this}; + ToolHandleArray<FastRegionFactory> m_fast_region_factories{this}; + Gaudi::Property<std::vector<std::string>> m_fast_model_factories_names{this, "FastModelFactories", {}, + tool_array_setter( m_fast_model_factories, m_fast_model_factories_names), + Gaudi::Details::Property::ImmediatelyInvokeHandler{true}}; + + Gaudi::Property<std::vector<std::string>> m_fast_region_factories_names{this, "FastRegionFactories", {}, + tool_array_setter( m_fast_region_factories, m_fast_region_factories_names), + Gaudi::Details::Property::ImmediatelyInvokeHandler{true}}; + private: SensDetVolumeMap m_sens_dets; Gaudi::Property<SensDetNameVolumesMap> m_namemap{this, "SensDetVolumeMap", {},[this]( Gaudi::Details::PropertyBase& ){ diff --git a/Sim/GiGaMTFastSimulation/CMakeLists.txt b/Sim/GiGaMTFastSimulation/CMakeLists.txt new file mode 100644 index 0000000000000000000000000000000000000000..6c3a3a05d5a08943fc661c92de0c1761346bd30b --- /dev/null +++ b/Sim/GiGaMTFastSimulation/CMakeLists.txt @@ -0,0 +1,29 @@ +################################################################################ +# Package: GiGaMTFastSimulation +# +# Package that ... +# +################################################################################# +gaudi_subdir(GiGaMTFastSimulation v1r0) + +gaudi_depends_on_subdirs(GaudiAlg + Sim/GiGaMTFactories + Sim/GiGaMTCore) +AddHepMC3() + +if(${Geant4_config_version} VERSION_LESS "10.06") + add_definitions(-DG4MULTITHREADED) + add_definitions(-DG4USE_STD11) +endif() + +gaudi_add_library(GiGaMTFastSimulationLib + src/lib/*.cpp + PUBLIC_HEADERS GiGaMTFastSimulation + INCLUDE_DIRS GiGaMTFactories + LINK_LIBRARIES GaudiAlgLib GiGaMTCoreRunLib) + +gaudi_add_module(GiGaMTFastSimulation + src/components/*.cpp + INCLUDE_DIRS GiGaMTFactories + LINK_LIBRARIES GaudiAlgLib GiGaMTCoreRunLib) + diff --git a/Sim/GiGaMTFastSimulation/GiGaMTFastSimulation/GiGaMTFastHit.h b/Sim/GiGaMTFastSimulation/GiGaMTFastSimulation/GiGaMTFastHit.h new file mode 100644 index 0000000000000000000000000000000000000000..188a4c98753b25073842ae053a19b988612e967a --- /dev/null +++ b/Sim/GiGaMTFastSimulation/GiGaMTFastSimulation/GiGaMTFastHit.h @@ -0,0 +1,19 @@ +#pragma once + +#include "Geant4/G4ThreeVector.hh" + +class GiGaMTFastHit { + + double m_Energy = 0; + G4ThreeVector m_Position = G4ThreeVector(); + + public: + inline GiGaMTFastHit() : m_Energy(), m_Position(G4ThreeVector()) {}; + inline GiGaMTFastHit(const G4ThreeVector& aPosition, double aEnergy) : m_Energy(aEnergy), m_Position(aPosition) {}; + virtual ~GiGaMTFastHit() = default; + + inline void SetEnergy(const double& aEnergy) { m_Energy = aEnergy; } + inline double GetEnergy() const { return m_Energy; } + inline void SetPosition(const G4ThreeVector& aPosition) { m_Position = aPosition; } + inline G4ThreeVector GetPosition() const { return m_Position; } +}; diff --git a/Sim/GiGaMTFastSimulation/GiGaMTFastSimulation/GiGaMTFastSimModelFAC.h b/Sim/GiGaMTFastSimulation/GiGaMTFastSimulation/GiGaMTFastSimModelFAC.h new file mode 100644 index 0000000000000000000000000000000000000000..a42c727237dc351da53c7696564fb2879c8fbc60 --- /dev/null +++ b/Sim/GiGaMTFastSimulation/GiGaMTFastSimulation/GiGaMTFastSimModelFAC.h @@ -0,0 +1,48 @@ +/*****************************************************************************\ +* (c) Copyright 2021 CERN for the benefit of the LHCb Collaboration * +* * +* This software is distributed under the terms of the GNU General Public * +* Licence version 3 (GPL Version 3), copied verbatim in the file "COPYING". * +* * +* In applying this licence, CERN does not waive the privileges and immunities * +* granted to it by virtue of its status as an Intergovernmental Organization * +* or submit itself to any jurisdiction. * +\*****************************************************************************/ +#pragma once + +#include "Geant4/G4RegionStore.hh" +#include "Geant4/G4VFastSimulationModel.hh" +#include "GiGaMTFactories/GiGaFactoryBase.h" +#include "GiGaMTFactories/GiGaTool.h" + +template <typename FastSimModel> +class GiGaMTFastSimModelFAC : public extends<GiGaTool, GiGaFactoryBase<G4VFastSimulationModel>> { + static_assert( std::is_base_of<G4VFastSimulationModel, FastSimModel>::value ); + + + Gaudi::Property<std::string> m_region{this, "RegionName", std::string()}; + Gaudi::Property<std::string> m_model{this, "ModelName", std::string()}; +public: + using extends::extends; + FastSimModel* construct() const override { + + if (m_model.value().empty()) { + error() << "Fast model name was not provided." << endmsg; + } + + if (m_region.value().empty()) { + error() << "Fast region name was not provided." << endmsg; + } + + if ( msgLevel( MSG::DEBUG ) ) { debug() << "Loading fast simulation region " << m_region.value() << endmsg; } + G4Region* region = G4RegionStore::GetInstance()->GetRegion( m_region.value() ); + + if ( msgLevel( MSG::DEBUG ) ) { debug() << "Loading Fast Sim Model " << m_model.value() << endmsg; } + auto tmp = new FastSimModel{m_model.value(), region}; + tmp->SetMessageInterface( message_interface() ); + + if ( msgLevel( MSG::DEBUG ) ) { debug() << "Loaded Fast Sim Model " << m_model.value() << endmsg; } + return tmp; + } + +}; diff --git a/Sim/GiGaMTFastSimulation/GiGaMTFastSimulation/GiGaMTFastSimSensitiveResponse.h b/Sim/GiGaMTFastSimulation/GiGaMTFastSimulation/GiGaMTFastSimSensitiveResponse.h new file mode 100644 index 0000000000000000000000000000000000000000..f893c9004fa7abe515c05258be11aeca3711621c --- /dev/null +++ b/Sim/GiGaMTFastSimulation/GiGaMTFastSimulation/GiGaMTFastSimSensitiveResponse.h @@ -0,0 +1,29 @@ +#pragma once + +#include "Geant4/G4TouchableHandle.hh" +#include "Geant4/G4Navigator.hh" +#include "Geant4/G4FastTrack.hh" + +#include "GiGaMTCoreMessage/IGiGaMessage.h" +#include "GiGaMTFastSimulation/GiGaMTFastHit.h" + + +class GiGaMTFastSimSensitiveResponse : public GiGaMessage { + + G4Navigator* m_navigator; + G4TouchableHandle m_touchableHandle; + bool m_navigatorOn; + std::string m_worldWithSdName; + +public: + + GiGaMTFastSimSensitiveResponse(); + inline virtual ~GiGaMTFastSimSensitiveResponse() { delete m_navigator; }; + + void simulate(const GiGaMTFastHit& aHit, const G4FastTrack& aTrack); + inline void SetNameOfWorldWithSD(const std::string& aName) { m_worldWithSdName = aName; } + +protected: + + virtual void modifyFakeStep(G4Step* /*fakeStep*/) {}; +}; diff --git a/Sim/GiGaMTFastSimulation/src/components/GiGaMTFastGiGaRegionFAC.cpp b/Sim/GiGaMTFastSimulation/src/components/GiGaMTFastGiGaRegionFAC.cpp new file mode 100644 index 0000000000000000000000000000000000000000..e043cc741f732a99a9468e0b67d8b3f64dfae30f --- /dev/null +++ b/Sim/GiGaMTFastSimulation/src/components/GiGaMTFastGiGaRegionFAC.cpp @@ -0,0 +1,62 @@ +#include "GiGaMTFastGiGaRegionFAC.h" + +G4Region* GiGaMTFastG4RegionFAC::construct() const { + if ( m_region_name.value().empty() && m_det_name.value().empty() ) { + error() << "No name specified for the new fast region." << endmsg; + } + + auto sens_det_name = m_det_name.value() + m_sens_det_name_suffix.value(); + auto region_name = m_region_name.value().empty() ? sens_det_name + "FastRegion" : m_region_name.value() ; + + auto region = G4RegionStore::GetInstance()->GetRegion( region_name ); + if ( region ) { + warning() << "Fast Region '" + region_name + "' already exists " << endmsg; + return region; + } + + auto sdmanager = G4SDManager::GetSDMpointer(); + std::vector<G4LogicalVolume*> lvolumes; + auto volume_store = G4LogicalVolumeStore::GetInstance(); + + if (!volume_store) { + error() << "Cannot access volume store." << endmsg; + } + if ( !m_det_name.value().empty() ) { + auto sens_det = sdmanager->FindSensitiveDetector( sens_det_name ); + if (!sens_det) { + error() << "Cannot access sensitive detector " << sens_det_name << endmsg; + } + for( auto ivol = volume_store->begin(); volume_store->end() != ivol ; ++ivol ) { + auto vol_sens_det = (*ivol)->GetSensitiveDetector(); + if (vol_sens_det && sens_det == vol_sens_det) { + if (msgLevel( MSG::DEBUG )) { + debug() << "Found " << (*ivol)->GetName() << " in " << sens_det->GetName() << endmsg; + } + lvolumes.push_back(*ivol); + } + } + } else if ( !m_volumes.value().empty() ) { + for ( auto& ivolume : m_volumes.value() ) { + lvolumes.push_back(volume_store->GetVolume( ivolume )); + } + } else { + error() << "No G4LogicalVolumes provided nor found for the region " << region_name << endmsg; + } + + region = new G4Region( region_name ); + + for ( auto& volume : lvolumes ) { + if ( !volume ) { + error() << "G4LogicalVolume '" + volume->GetName() + "' is invalid "; + } + if ( volume->GetRegion() ) { + error() << "G4LogicalVolume '" + volume->GetName() + "' already belongs to another region '" + volume->GetRegion()->GetName(); + } + volume->SetRegion( region ); + region->AddRootLogicalVolume( volume ); + } + + return region; +} + +DECLARE_COMPONENT_WITH_ID(GiGaMTFastG4RegionFAC, "FastG4RegionFAC") diff --git a/Sim/GiGaMTFastSimulation/src/components/GiGaMTFastGiGaRegionFAC.h b/Sim/GiGaMTFastSimulation/src/components/GiGaMTFastGiGaRegionFAC.h new file mode 100644 index 0000000000000000000000000000000000000000..494345ef728ba63a609fb6f7c948bc5b4bcff494 --- /dev/null +++ b/Sim/GiGaMTFastSimulation/src/components/GiGaMTFastGiGaRegionFAC.h @@ -0,0 +1,23 @@ +#pragma once + +#include "GiGaMTFactories/GiGaTool.h" +#include "GiGaMTFactories/GiGaFactoryBase.h" + +#include "Geant4/G4Region.hh" +#include "Geant4/G4RegionStore.hh" +#include "Geant4/G4LogicalVolume.hh" +#include "Geant4/G4LogicalVolumeStore.hh" +#include "Geant4/G4SDManager.hh" + +class GiGaMTFastG4RegionFAC : public extends<GiGaTool, GiGaFactoryBase<G4Region>> { + + Gaudi::Property<std::string> m_region_name{this, "RegionName", std::string()}; + Gaudi::Property<std::string> m_det_name{this, "DetName", std::string()}; + Gaudi::Property<std::string> m_sens_det_name_suffix{this, "SensDetNameSuffix", "SDet"}; + Gaudi::Property<std::vector<std::string>> m_volumes{this, "Volumes", {}}; + +public: + using extends::extends; + virtual G4Region* construct() const override; +}; + diff --git a/Sim/GiGaMTFastSimulation/src/components/GiGaMTFastSimPhysFAC.cpp b/Sim/GiGaMTFastSimulation/src/components/GiGaMTFastSimPhysFAC.cpp new file mode 100644 index 0000000000000000000000000000000000000000..c668112f758d23f3b8d6c983e47445e210a5a13c --- /dev/null +++ b/Sim/GiGaMTFastSimulation/src/components/GiGaMTFastSimPhysFAC.cpp @@ -0,0 +1,49 @@ +/*****************************************************************************\ +* (c) Copyright 2021 CERN for the benefit of the LHCb Collaboration * +* * +* This software is distributed under the terms of the GNU General Public * +* Licence version 3 (GPL Version 3), copied verbatim in the file "COPYING". * +* * +* In applying this licence, CERN does not waive the privileges and immunities * +* granted to it by virtue of its status as an Intergovernmental Organization * +* or submit itself to any jurisdiction. * +\*****************************************************************************/ +#include "Geant4/G4FastSimulationManagerProcess.hh" +#include "Geant4/G4ProcessManager.hh" +#include "Geant4/G4VPhysicsConstructor.hh" +#include "GiGaMTCoreMessage/IGiGaMessage.h" +#include "GiGaMTFactories/GiGaMTG4PhysicsConstrFAC.h" + +class GiGaMTFastSimPhys : public G4VPhysicsConstructor, public GiGaMessage { +public: + GiGaMTFastSimPhys() = default; + void ConstructParticle() override{}; + void ConstructProcess() override { + // TODO: give a possibility to choose a set of particles + // TODO: an option to set parallel geometry + G4FastSimulationManagerProcess* fastSimProcess = new G4FastSimulationManagerProcess( "G4FSMP" ); + auto theParticleIterator = GetParticleIterator(); + theParticleIterator->reset(); + // Fast simulation manager process is available for all the particles + while ( ( *theParticleIterator )() ) { + G4ParticleDefinition* particle = theParticleIterator->value(); + G4ProcessManager* process_manager = particle->GetProcessManager(); + process_manager->AddDiscreteProcess( fastSimProcess ); + } + } +}; + +class GiGaMTFastSimPhysFAC : public extends<GiGaMTPhysConstr, GiGaFactoryBase<G4VPhysicsConstructor>> { +public: + using extends::extends; + GiGaMTFastSimPhys* construct() const override { + if ( msgLevel( MSG::DEBUG ) ) { debug() << "Constructing fast simulation physics" << endmsg; } + auto tmp = new GiGaMTFastSimPhys{}; + tmp->SetMessageInterface( message_interface() ); + tmp->SetVerboseLevel( verbosity() ); + tmp->SetPhysicsName( name() ); + if ( msgLevel( MSG::DEBUG ) ) { debug() << "Constructed fast simulation physics" << endmsg; } + return tmp; + } +}; +DECLARE_COMPONENT_WITH_ID( GiGaMTFastSimPhysFAC, "GiGaMTFastSimPhysFAC" ) diff --git a/Sim/GiGaMTFastSimulation/src/lib/GiGaMTFastSimSensitiveResponse.cpp b/Sim/GiGaMTFastSimulation/src/lib/GiGaMTFastSimSensitiveResponse.cpp new file mode 100644 index 0000000000000000000000000000000000000000..cbaa1f470611ecb6f4605b056a83448821b498f3 --- /dev/null +++ b/Sim/GiGaMTFastSimulation/src/lib/GiGaMTFastSimSensitiveResponse.cpp @@ -0,0 +1,57 @@ + +#include "GiGaMTFastSimulation/GiGaMTFastSimSensitiveResponse.h" + +#include "Geant4/G4TransportationManager.hh" +#include "Geant4/G4VSensitiveDetector.hh" +#include "Geant4/G4TouchableHandle.hh" + +GiGaMTFastSimSensitiveResponse::GiGaMTFastSimSensitiveResponse() { + m_touchableHandle = new G4TouchableHistory(); + m_navigator = new G4Navigator(); + m_navigatorOn = false; + m_worldWithSdName = std::string(); +} + +void GiGaMTFastSimSensitiveResponse::simulate(const GiGaMTFastHit& aHit, const G4FastTrack& aTrack) { + if(aHit.GetEnergy() <= 0) return; + + auto position = aTrack.GetPrimaryTrack()->GetPosition(); + if(!m_navigatorOn) { + G4VPhysicalVolume* worldWithSD = nullptr; + if(m_worldWithSdName.empty()) { + worldWithSD = G4TransportationManager::GetTransportationManager()->GetNavigatorForTracking()->GetWorldVolume(); + } + else { + worldWithSD = G4TransportationManager::GetTransportationManager()->GetParallelWorld(m_worldWithSdName); + } + m_navigator->SetWorldVolume(worldWithSD); + m_navigator->LocateGlobalPointAndUpdateTouchable(position, m_touchableHandle(), false); + m_navigatorOn = true; + } + else { + //m_navigator->LocateGlobalPointAndUpdateTouchable(aTrack.GetInverseAffineTransformation()->TransformPoint(aHit.GetPosition()), G4ThreeVector( 0., 0., 0. ) ,m_touchableHandle()); + m_navigator->LocateGlobalPointAndUpdateTouchable(position, m_touchableHandle()); + } + G4VPhysicalVolume* currentVolume = m_touchableHandle()->GetVolume(); + if(currentVolume != 0) { + auto lvol = currentVolume->GetLogicalVolume(); + auto sensitive = lvol->GetSensitiveDetector(); + // TODO: check whether this is sufficient to make sure that the right volume is chosen + // in principle here the problem is when an ancestor volume is linked to the sensitive det + if(auto mvol = currentVolume->GetMotherLogical(); !sensitive && mvol->GetRegion() == lvol->GetRegion()) { + sensitive = mvol->GetSensitiveDetector(); + } + if(sensitive) { // && currentVolume->GetLogicalVolume()->GetFastSimulationManager()) { + auto fakeStep = new G4Step(); + auto fakePreStepPoint = fakeStep->GetPreStepPoint(); + fakePreStepPoint->SetTouchableHandle( m_touchableHandle ); + fakePreStepPoint->SetMaterialCutsCouple( aTrack.GetPrimaryTrack()->GetMaterialCutsCouple() ); + fakeStep->SetTotalEnergyDeposit(aHit.GetEnergy()); + fakeStep->SetTrack(const_cast<G4Track*>(aTrack.GetPrimaryTrack())); + modifyFakeStep(fakeStep); // if any other custom modifications have to be applied + sensitive->Hit(fakeStep); + } else { + error("Could not locate sensitive detector for the logical volume: " + lvol->GetName()); + } + } +}