Skip to content
Snippets Groups Projects
Commit 57252652 authored by John Chapman's avatar John Chapman
Browse files

Patches from Makoto Asai to allow particles with pre-defined decays and zero lifetime.

parent ec226de0
No related merge requests found
......@@ -24,7 +24,7 @@
// ********************************************************************
//
//
// $Id: G4PrimaryTransformer.cc 72252 2013-07-12 09:04:11Z gcosmo $
// $Id: G4PrimaryTransformer.cc 101916 2016-12-08 08:04:03Z gcosmo $
//
#include "G4PrimaryTransformer.hh"
......@@ -42,8 +42,8 @@
G4PrimaryTransformer::G4PrimaryTransformer()
:verboseLevel(0),trackID(0),
unknown(0),unknownParticleDefined(false),
opticalphoton(0),opticalphotonDefined(false),
unknown(nullptr),unknownParticleDefined(false),
opticalphoton(nullptr),opticalphotonDefined(false),
nWarn(0)
{
particleTable = G4ParticleTable::GetParticleTable();
......@@ -72,11 +72,12 @@ G4TrackVector* G4PrimaryTransformer::GimmePrimaries(G4Event* anEvent,G4int track
trackID = trackIDCounter;
//TV.clearAndDestroy();
for( size_t ii=0; ii<TV.size();ii++)
{ delete TV[ii]; }
for(auto tr : TV) delete tr;
TV.clear();
//Loop over vertices
G4PrimaryVertex* nextVertex = anEvent->GetPrimaryVertex();
while(nextVertex)
while(nextVertex) // Loop checking 12.28.2015 M.Asai
{
GenerateTracks(nextVertex);
nextVertex = nextVertex->GetNext();
......@@ -105,7 +106,7 @@ void G4PrimaryTransformer::GenerateTracks(G4PrimaryVertex* primaryVertex)
#endif
G4PrimaryParticle* primaryParticle = primaryVertex->GetPrimary();
while( primaryParticle != 0 )
while( primaryParticle != 0 ) // Loop checking 12.28.2015 M.Asai
{
GenerateSingleTrack( primaryParticle, X0, Y0, Z0, T0, WV );
primaryParticle = primaryParticle->GetNext();
......@@ -128,7 +129,7 @@ void G4PrimaryTransformer::GenerateSingleTrack
}
#endif
G4PrimaryParticle* daughter = primaryParticle->GetDaughter();
while(daughter)
while(daughter) // Loop checking 12.28.2015 M.Asai
{
GenerateSingleTrack(daughter,x0,y0,z0,t0,wv);
daughter = daughter->GetNext();
......@@ -180,7 +181,7 @@ void G4PrimaryTransformer::GenerateSingleTrack
primaryParticle->GetPolY(),
primaryParticle->GetPolZ());
}
if(primaryParticle->GetProperTime()>0.0)
if(primaryParticle->GetProperTime()>=0.0)
{ DP->SetPreAssignedDecayProperTime(primaryParticle->GetProperTime()); }
// Set Mass if it is specified
......@@ -270,7 +271,7 @@ void G4PrimaryTransformer::SetDecayProducts
= new G4DynamicParticle(partDef,daughter->GetMomentum());
DP->SetPrimaryParticle(daughter);
// Decay proper time for daughter
if(daughter->GetProperTime()>0.0)
if(daughter->GetProperTime()>=0.0)
{ DP->SetPreAssignedDecayProperTime(daughter->GetProperTime()); }
// Set Charge is specified
if (daughter->GetCharge()<DBL_MAX) {
......@@ -279,6 +280,7 @@ void G4PrimaryTransformer::SetDecayProducts
G4double pmas = daughter->GetMass();
if(pmas>=0.)
{ DP->SetMass(pmas); }
DP->SetPolarization(daughter->GetPolX(),daughter->GetPolY(),daughter->GetPolZ());
decayProducts->PushProducts(DP);
SetDecayProducts(daughter,DP);
// Check the particle is properly constructed
......
......@@ -41,7 +41,7 @@ G4PrimaryParticle::G4PrimaryParticle()
direction(0.,0.,1.),kinE(0.),
nextParticle(0),daughterParticle(0),trackID(-1),
mass(-1.),charge(0.),polX(0.),polY(0.),polZ(0.),
Weight0(1.0),properTime(0.0),userInfo(0)
Weight0(1.0),properTime(-1.0),userInfo(0)
{;}
G4PrimaryParticle::G4PrimaryParticle(G4int Pcode)
......@@ -49,7 +49,7 @@ G4PrimaryParticle::G4PrimaryParticle(G4int Pcode)
direction(0.,0.,1.),kinE(0.),
nextParticle(0),daughterParticle(0),trackID(-1),
mass(-1.),charge(0.),polX(0.),polY(0.),polZ(0.),
Weight0(1.0),properTime(0.0),userInfo(0)
Weight0(1.0),properTime(-1.0),userInfo(0)
{
G4code = G4ParticleTable::GetParticleTable()->FindParticle(Pcode);
if (G4code !=0) {
......@@ -64,7 +64,7 @@ G4PrimaryParticle::G4PrimaryParticle(G4int Pcode,
direction(0.,0.,1.),kinE(0.),
nextParticle(0),daughterParticle(0),trackID(-1),
mass(-1.),charge(0.),polX(0.),polY(0.),polZ(0.),
Weight0(1.0),properTime(0.0),userInfo(0)
Weight0(1.0),properTime(-1.0),userInfo(0)
{
G4code = G4ParticleTable::GetParticleTable()->FindParticle(Pcode);
if (G4code !=0) {
......@@ -80,7 +80,7 @@ G4PrimaryParticle::G4PrimaryParticle(G4int Pcode,
direction(0.,0.,1.),kinE(0.),
nextParticle(0),daughterParticle(0),trackID(-1),
mass(-1.),charge(0.),polX(0.),polY(0.),polZ(0.),
Weight0(1.0),properTime(0.0),userInfo(0)
Weight0(1.0),properTime(-1.0),userInfo(0)
{
G4code = G4ParticleTable::GetParticleTable()->FindParticle(Pcode);
if (G4code !=0) {
......@@ -95,7 +95,7 @@ G4PrimaryParticle::G4PrimaryParticle(const G4ParticleDefinition* Gcode)
direction(0.,0.,1.),kinE(0.),
nextParticle(0),daughterParticle(0),trackID(-1),
mass(-1.),charge(0.),polX(0.),polY(0.),polZ(0.),
Weight0(1.0),properTime(0.0),userInfo(0)
Weight0(1.0),properTime(-1.0),userInfo(0)
{
if (G4code !=0) {
PDGcode = Gcode->GetPDGEncoding();
......@@ -110,7 +110,7 @@ G4PrimaryParticle::G4PrimaryParticle(const G4ParticleDefinition* Gcode,
direction(0.,0.,1.),kinE(0.),
nextParticle(0),daughterParticle(0),trackID(-1),
mass(-1.),charge(0.),polX(0.),polY(0.),polZ(0.),
Weight0(1.0),properTime(0.0),userInfo(0)
Weight0(1.0),properTime(-1.0),userInfo(0)
{
if (G4code !=0) {
PDGcode = Gcode->GetPDGEncoding();
......@@ -126,7 +126,7 @@ G4PrimaryParticle::G4PrimaryParticle(const G4ParticleDefinition* Gcode,
direction(0.,0.,1.),kinE(0.),
nextParticle(0),daughterParticle(0),trackID(-1),
mass(-1.),charge(0.),polX(0.),polY(0.),polZ(0.),
Weight0(1.0),properTime(0.0),userInfo(0)
Weight0(1.0),properTime(-1.0),userInfo(0)
{
if (G4code !=0) {
PDGcode = Gcode->GetPDGEncoding();
......@@ -141,7 +141,7 @@ G4PrimaryParticle::G4PrimaryParticle(const G4PrimaryParticle& right)
direction(0.,0.,1.),kinE(0.),
nextParticle(0),daughterParticle(0),trackID(-1),
mass(-1.),charge(0.),polX(0.),polY(0.),polZ(0.),
Weight0(1.0),properTime(0.0),userInfo(0)
Weight0(1.0),properTime(-1.0),userInfo(0)
{
*this = right;
}
......@@ -281,7 +281,7 @@ void G4PrimaryParticle::Print() const
<< polZ << " )"
<< G4endl;
G4cout << " Weight : " << Weight0 << G4endl;
if(properTime>0.0) {
if(properTime>=0.0) {
G4cout << " PreAssigned proper decay time : " << properTime/ns << " [ns] " << G4endl;
}
if(userInfo != 0) { userInfo->Print(); }
......
......@@ -24,7 +24,6 @@
// ********************************************************************
//
//
// $Id: G4Decay.cc 79343 2014-02-24 11:44:06Z gcosmo $
//
//
// --------------------------------------------------------------
......@@ -65,7 +64,7 @@ G4Decay::G4Decay(const G4String& processName)
verboseLevel(1),
HighestValue(20.0),
fRemainderLifeTime(-1.0),
pExtDecayer(0)
pExtDecayer(nullptr)
{
// set Process Sub Type
SetProcessSubType(static_cast<int>(DECAY));
......@@ -81,7 +80,7 @@ G4Decay::G4Decay(const G4String& processName)
G4Decay::~G4Decay()
{
if (pExtDecayer) {
if (pExtDecayer != nullptr) {
delete pExtDecayer;
}
}
......@@ -198,17 +197,17 @@ G4VParticleChange* G4Decay::DecayIt(const G4Track& aTrack, const G4Step& )
//check if thePreAssignedDecayProducts exists
const G4DecayProducts* o_products = (aParticle->GetPreAssignedDecayProducts());
G4bool isPreAssigned = (o_products != 0);
G4DecayProducts* products = 0;
G4bool isPreAssigned = (o_products != nullptr);
G4DecayProducts* products = nullptr;
// decay table
G4DecayTable *decaytable = aParticleDef->GetDecayTable();
// check if external decayer exists
G4bool isExtDecayer = (decaytable == 0) && (pExtDecayer !=0);
G4bool isExtDecayer = (decaytable == nullptr) && (pExtDecayer != nullptr);
// Error due to NO Decay Table
if ( (decaytable == 0) && !isExtDecayer &&!isPreAssigned ){
if ( (decaytable == nullptr) && !isExtDecayer && !isPreAssigned ){
if (GetVerboseLevel()>0) {
G4cout << "G4Decay::DoIt : decay table not defined for ";
G4cout << aParticle->GetDefinition()->GetParticleName()<< G4endl;
......@@ -417,7 +416,8 @@ G4double G4Decay::PostStepGetPhysicalInteractionLength(
//pre-assigned Decay time case
// reminder proper time
fRemainderLifeTime = pTime - track.GetProperTime();
if (fRemainderLifeTime <= 0.0) fRemainderLifeTime = DBL_MIN;
//############## if (fRemainderLifeTime <= 0.0) fRemainderLifeTime = DBL_MIN;
if (fRemainderLifeTime <= 0.0) fRemainderLifeTime = 0.0;
G4double rvalue=0.0;
// use pre-assigned Decay time to determine PIL
......
......@@ -146,6 +146,7 @@ G4RunManagerKernel::G4RunManagerKernel()
<< " * G4AtlasRK4 stepper" << G4endl
<< " * Patch to G4Track https://its.cern.ch/jira/browse/ATLASSIM-2409" << G4endl
<< " * Patch to Geant4MakeRules_cxx.cmake https://its.cern.ch/jira/browse/ATLASSIM-3144" << G4endl
<< " * Patch to G4Decay.cc, G4PrimaryTransformer.cc and G4PrimaryParticle.cc https://its.cern.ch/jira/browse/ATLASSIM-3313" << G4endl
<< "*************************************************************" << G4endl
<< G4endl;
}
......@@ -235,6 +236,7 @@ numberOfParallelWorld(0),geometryNeedsToBeClosed(true),
<< " * G4AtlasRK4 stepper" << G4endl
<< " * Patch to G4Track https://its.cern.ch/jira/browse/ATLASSIM-2409" << G4endl
<< " * Patch to Geant4MakeRules_cxx.cmake https://its.cern.ch/jira/browse/ATLASSIM-3144" << G4endl
<< " * Patch to G4Decay.cc, G4PrimaryTransformer.cc and G4PrimaryParticle.cc https://its.cern.ch/jira/browse/ATLASSIM-3313" << G4endl
<< "*************************************************************" << G4endl
<< G4endl;
break;
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment