Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4DNAMolecularDissociation Class Reference

#include <G4DNAMolecularDissociation.hh>

Inheritance diagram for G4DNAMolecularDissociation:
Collaboration diagram for G4DNAMolecularDissociation:

Public Member Functions

 G4DNAMolecularDissociation (const G4String &processName="DNAMolecularDecay", G4ProcessType type=fDecay)
 
virtual ~G4DNAMolecularDissociation ()
 
virtual G4bool IsApplicable (const G4ParticleDefinition &)
 
virtual G4double PostStepGetPhysicalInteractionLength (const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)
 
G4double AtRestGetPhysicalInteractionLength (const G4Track &track, G4ForceCondition *condition)
 
G4VParticleChangeAtRestDoIt (const G4Track &track, const G4Step &step)
 
G4VParticleChangePostStepDoIt (const G4Track &track, const G4Step &step)
 
void SetVerbose (G4int)
 
void SetDecayDisplacer (const G4ParticleDefinition *, G4VMolecularDecayDisplacer *)
 
G4VMolecularDecayDisplacerGetDecayDisplacer (const G4ParticleDefinition *)
 
void SetDisplacer (const G4ParticleDefinition *, G4VMolecularDecayDisplacer *)
 
G4VMolecularDecayDisplacerGetDisplacer (const G4ParticleDefinition *)
 
- Public Member Functions inherited from G4VITRestDiscreteProcess
 G4VITRestDiscreteProcess (const G4String &, G4ProcessType aType=fNotDefined)
 
 G4VITRestDiscreteProcess (const G4VITRestDiscreteProcess &)
 
virtual ~G4VITRestDiscreteProcess ()
 
virtual G4double AlongStepGetPhysicalInteractionLength (const G4Track &, G4double, G4double, G4double &, G4GPILSelection *)
 
virtual G4VParticleChangeAlongStepDoIt (const G4Track &, const G4Step &)
 
- Public Member Functions inherited from G4VITProcess
 G4VITProcess (const G4String &name, G4ProcessType type=fNotDefined)
 
virtual ~G4VITProcess ()
 
 G4VITProcess (const G4VITProcess &other)
 
G4VITProcessoperator= (const G4VITProcess &other)
 
G4int operator== (const G4VITProcess &right) const
 
G4int operator!= (const G4VITProcess &right) const
 
size_t GetProcessID () const
 
G4shared_ptr< G4ProcessState_LockGetProcessState ()
 
void SetProcessState (G4shared_ptr< G4ProcessState_Lock > aProcInfo)
 
void ResetProcessState ()
 
virtual void StartTracking (G4Track *)
 
virtual void BuildPhysicsTable (const G4ParticleDefinition &)
 
G4double GetInteractionTimeLeft ()
 
virtual void ResetNumberOfInteractionLengthLeft ()
 
G4bool ProposesTimeStep () const
 
- Public Member Functions inherited from G4VProcess
 G4VProcess (const G4String &aName="NoName", G4ProcessType aType=fNotDefined)
 
 G4VProcess (const G4VProcess &right)
 
virtual ~G4VProcess ()
 
G4int operator== (const G4VProcess &right) const
 
G4int operator!= (const G4VProcess &right) const
 
G4double GetCurrentInteractionLength () const
 
void SetPILfactor (G4double value)
 
G4double GetPILfactor () const
 
G4double AlongStepGPIL (const G4Track &track, G4double previousStepSize, G4double currentMinimumStep, G4double &proposedSafety, G4GPILSelection *selection)
 
G4double AtRestGPIL (const G4Track &track, G4ForceCondition *condition)
 
G4double PostStepGPIL (const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)
 
virtual void PreparePhysicsTable (const G4ParticleDefinition &)
 
virtual G4bool StorePhysicsTable (const G4ParticleDefinition *, const G4String &, G4bool)
 
virtual G4bool RetrievePhysicsTable (const G4ParticleDefinition *, const G4String &, G4bool)
 
const G4StringGetPhysicsTableFileName (const G4ParticleDefinition *, const G4String &directory, const G4String &tableName, G4bool ascii=false)
 
const G4StringGetProcessName () const
 
G4ProcessType GetProcessType () const
 
void SetProcessType (G4ProcessType)
 
G4int GetProcessSubType () const
 
void SetProcessSubType (G4int)
 
virtual void EndTracking ()
 
virtual void SetProcessManager (const G4ProcessManager *)
 
virtual const G4ProcessManagerGetProcessManager ()
 
G4double GetNumberOfInteractionLengthLeft () const
 
G4double GetTotalNumberOfInteractionLengthTraversed () const
 
G4bool isAtRestDoItIsEnabled () const
 
G4bool isAlongStepDoItIsEnabled () const
 
G4bool isPostStepDoItIsEnabled () const
 
virtual void DumpInfo () const
 
void SetVerboseLevel (G4int value)
 
G4int GetVerboseLevel () const
 
virtual void SetMasterProcess (G4VProcess *masterP)
 
const G4VProcessGetMasterProcess () const
 
virtual void BuildWorkerPhysicsTable (const G4ParticleDefinition &part)
 
virtual void PrepareWorkerPhysicsTable (const G4ParticleDefinition &)
 

Protected Member Functions

virtual G4VParticleChangeDecayIt (const G4Track &, const G4Step &)
 
virtual G4double GetMeanLifeTime (const G4Track &, G4ForceCondition *)
 
virtual G4double GetMeanFreePath (const G4Track &, G4double, G4ForceCondition *)
 
- Protected Member Functions inherited from G4VITRestDiscreteProcess
 G4VITRestDiscreteProcess ()
 
G4VITRestDiscreteProcessoperator= (const G4VITRestDiscreteProcess &right)
 
- Protected Member Functions inherited from G4VITProcess
void RetrieveProcessInfo ()
 
void CreateInfo ()
 
template<typename T >
T * GetState ()
 
virtual void SubtractNumberOfInteractionLengthLeft (G4double previousStepSize)
 
virtual void ClearInteractionTimeLeft ()
 
virtual void ClearNumberOfInteractionLengthLeft ()
 
void SetInstantiateProcessState (G4bool flag)
 
G4bool InstantiateProcessState ()
 
- Protected Member Functions inherited from G4VProcess
void SubtractNumberOfInteractionLengthLeft (G4double previousStepSize)
 
void ClearNumberOfInteractionLengthLeft ()
 

Additional Inherited Members

- Static Public Member Functions inherited from G4VITProcess
static const size_t & GetMaxProcessIndex ()
 
- Static Public Member Functions inherited from G4VProcess
static const G4StringGetProcessTypeName (G4ProcessType)
 
- Protected Attributes inherited from G4VITProcess
G4shared_ptr< G4ProcessStatefpState
 
G4bool fProposesTimeStep
 
- Protected Attributes inherited from G4VProcess
const G4ProcessManageraProcessManager
 
G4VParticleChangepParticleChange
 
G4ParticleChange aParticleChange
 
G4double theNumberOfInteractionLengthLeft
 
G4double currentInteractionLength
 
G4double theInitialNumberOfInteractionLength
 
G4String theProcessName
 
G4String thePhysicsTableFileName
 
G4ProcessType theProcessType
 
G4int theProcessSubType
 
G4double thePILfactor
 
G4bool enableAtRestDoIt
 
G4bool enableAlongStepDoIt
 
G4bool enablePostStepDoIt
 
G4int verboseLevel
 

Detailed Description

G4DNAMolecularDissociation should be called only for molecules. It will dissociate the molecules using the decay associated to this molecule and if a displacement scheme has been registered, it will place the products to the expected position.

Definition at line 66 of file G4DNAMolecularDissociation.hh.

Constructor & Destructor Documentation

G4DNAMolecularDissociation::G4DNAMolecularDissociation ( const G4String processName = "DNAMolecularDecay",
G4ProcessType  type = fDecay 
)

Definition at line 52 of file G4DNAMolecularDissociation.cc.

53  :
54  G4VITRestDiscreteProcess(processName, type)
55 {
56  // set Process Sub Type
57  SetProcessSubType(59); // DNA sub-type
58  enableAlongStepDoIt = false;
59  enablePostStepDoIt = true;
60  enableAtRestDoIt = true;
61 
62  fVerbose = 0;
63 
64 #ifdef G4VERBOSE
65  if(verboseLevel > 1)
66  {
67  G4cout << "G4MolecularDissociationProcess constructor " << " Name:"
68  << processName << G4endl;
69  }
70 #endif
71 
73 
74  fDecayAtFixedTime = true;
75  fProposesTimeStep = true;
76 }
G4int verboseLevel
Definition: G4VProcess.hh:368
G4bool enableAtRestDoIt
Definition: G4VProcess.hh:350
G4GLOB_DLL std::ostream G4cout
void SetProcessSubType(G4int)
Definition: G4VProcess.hh:432
G4bool enablePostStepDoIt
Definition: G4VProcess.hh:352
G4VParticleChange * pParticleChange
Definition: G4VProcess.hh:283
G4ParticleChange aParticleChange
Definition: G4VProcess.hh:289
G4bool fProposesTimeStep
#define G4endl
Definition: G4ios.hh:61
G4bool enableAlongStepDoIt
Definition: G4VProcess.hh:351

Here is the call graph for this function:

G4DNAMolecularDissociation::~G4DNAMolecularDissociation ( )
virtual

Definition at line 80 of file G4DNAMolecularDissociation.cc.

81 {
82  DisplacementMap::iterator it = fDisplacementMap.begin();
83 
84  for(; it != fDisplacementMap.end(); it++)
85  {
86  if(it->second)
87  {
88  delete it->second;
89  it->second = 0;
90  }
91  }
92  fDisplacementMap.clear();
93 }

Member Function Documentation

G4VParticleChange * G4DNAMolecularDissociation::AtRestDoIt ( const G4Track track,
const G4Step step 
)
inlinevirtual

Reimplemented from G4VITRestDiscreteProcess.

Definition at line 155 of file G4DNAMolecularDissociation.hh.

159 {
162  return DecayIt(track, step);
163 }
virtual void ClearNumberOfInteractionLengthLeft()
virtual void ClearInteractionTimeLeft()
virtual G4VParticleChange * DecayIt(const G4Track &, const G4Step &)

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4DNAMolecularDissociation::AtRestGetPhysicalInteractionLength ( const G4Track track,
G4ForceCondition condition 
)
inlinevirtual

Reimplemented from G4VITRestDiscreteProcess.

Definition at line 142 of file G4DNAMolecularDissociation.hh.

145 {
146  if(fDecayAtFixedTime)
147  {
148  return GetMeanLifeTime(track, condition);
149  }
150 
152 // return G4VITRestProcess::AtRestGetPhysicalInteractionLength(track, condition);
153 }
G4double condition(const G4ErrorSymMatrix &m)
virtual G4double GetMeanLifeTime(const G4Track &, G4ForceCondition *)
virtual G4double AtRestGetPhysicalInteractionLength(const G4Track &, G4ForceCondition *)

Here is the call graph for this function:

G4VParticleChange * G4DNAMolecularDissociation::DecayIt ( const G4Track track,
const G4Step  
)
protectedvirtual

Definition at line 141 of file G4DNAMolecularDissociation.cc.

143 {
144  // DEBUG
145  // G4cout << "Is calling G4MolecularDecayProcess::DecayIt" << G4endl;
146 
148  const G4Molecule * theMotherMolecule = GetMolecule(track);
149  const G4MoleculeDefinition* moleculeDefinition = theMotherMolecule
150  ->GetDefinition();
151 
152  // DEBUG
153  // G4cout <<"Calling G4MolecularDecayProcess::DecayIt"<<G4endl;
154  // G4cout << "The mother molecule state : " << G4endl;
155  // theMotherMolecule -> PrintState();
156 
157  if(moleculeDefinition->GetDecayTable())
158  {
159  const vector<const G4MolecularDissociationChannel*>* DecayVector =
160  (theMotherMolecule->GetDecayChannel());
161 
162  if(DecayVector == 0)
163  {
164  G4ExceptionDescription exceptionDescription;
165  theMotherMolecule->PrintState();
166  exceptionDescription << "No decay channel was found for the molecule : "
167  << theMotherMolecule->GetName() << G4endl;
168  G4Exception("G4DNAMolecularDissociation::DecayIt",
169  "G4DNAMolecularDissociation::NoDecayChannel",
171  exceptionDescription);
172  return &aParticleChange;
173  }
174 
175  G4int DecayVectorSize = DecayVector->size();
176  // DEBUG
177  // G4cout<< "Number of decay channels : " << DecayVectorSize <<G4endl;
178  G4double RdmValue = G4UniformRand();
179 
180  const G4MolecularDissociationChannel* decayChannel(0);
181  G4int i = 0;
182  do
183  {
184  decayChannel = (*DecayVector)[i];
185  if(RdmValue < decayChannel->GetProbability()) break;
186  RdmValue -= decayChannel->GetProbability();
187  i++;
188  }
189  while(i < DecayVectorSize);
190 
191  // DEBUG
192  // G4cout << "Selected Decay channel : "
193  // << decayChannel->GetName() << G4endl;
194 
195  G4double decayEnergy = decayChannel->GetEnergy();
196  G4int nbProducts = decayChannel->GetNbProducts();
197 
198  if(decayEnergy)
199  {
200  // DEBUG
201  // G4cout << "Deposit energy :"
202  // << decayChannel->GetEnergy()/eV << " eV" << G4endl;
203 
204  aParticleChange.ProposeLocalEnergyDeposit(decayChannel->GetEnergy());
205  }
206 
207  if(nbProducts)
208  {
209 
210  // DEBUG
211  // G4cout << "Number of products :" << nbProducts << G4endl;
212 
213  vector<G4ThreeVector> ProductsDisplacement(nbProducts);
214  G4ThreeVector theMotherMoleculeDisplacement;
215 
216  DisplacementMap::iterator it =
217  fDisplacementMap.find(moleculeDefinition);
218 
219  if(it != fDisplacementMap.end())
220  {
221  G4VMolecularDecayDisplacer* displacer = it->second;
222  ProductsDisplacement = displacer->GetProductsDisplacement(decayChannel);
223  theMotherMoleculeDisplacement =
224  displacer->GetMotherMoleculeDisplacement(decayChannel);
225  }
226  else
227  {
228  G4ExceptionDescription errMsg;
229  errMsg << "No G4MolecularDecayProcess::theDecayDisplacementMap["
230  << theMotherMolecule->GetName() + "]";
231  G4Exception("G4MolecularDecayProcess::DecayIt",
232  "DNAMolecularDecay001",
234  errMsg);
235  }
236 
238 
239 // G4cout << " nbProducts = " << nbProducts << G4endl;
240 // theMotherMolecule->PrintState();
241 
242 #ifdef G4VERBOSE
243  if(fVerbose)
244  {
245  G4cout << "Decay Process : " << theMotherMolecule->GetName()
246  << " (trackID :" << track.GetTrackID() << ") "
247  << decayChannel->GetName() << G4endl;
248  }
249 #endif
250 
251  G4ITNavigator* navigator =
254 
255  for(G4int j = 0; j < nbProducts; j++)
256  {
257  G4Molecule* product = new G4Molecule(decayChannel->GetProduct(j));
258 
259  G4ThreeVector displacement = theMotherMoleculeDisplacement
260  + ProductsDisplacement[j];
261  double mag_displacement = displacement.mag();
262  G4ThreeVector displacement_direction =
263  displacement/(mag_displacement+1e-30);
264 
265  double prNewSafety = DBL_MAX;
266 
267  //double step =
268  navigator->CheckNextStep(track.GetPosition(),
269  displacement_direction,
270  mag_displacement,
271  prNewSafety);
272 
273  //if(prNewSafety < mag_displacement || step < mag_displacement)
274  mag_displacement = min(prNewSafety*0.8, mag_displacement);
275 
276  G4ThreeVector product_pos = track.GetPosition()
277  + displacement_direction * mag_displacement;
278 
279  const G4AffineTransform& transform = navigator
280  ->GetGlobalToLocalTransform();
281 
282  G4ThreeVector localPoint =
283  transform.TransformPoint(product_pos); //track.GetPosition());
284 
285  if(track.GetTouchable()->GetSolid()->Inside(localPoint) !=
287  {
288  G4cout << "Mother volume: "
289  << track.GetTouchable()->GetVolume()->GetName()
290  << G4endl;
291  G4Exception("G4DNAMolecularDissociation::DecayIt",
292  "OUTSIDE_OF_MOTHER_VOLUME",
294  "Product has been placed outside of the volume "
295  "containing the mother molecule");
296  }
297 
298  G4Track* secondary =
299  product->BuildTrack(track.GetGlobalTime(),
300  product_pos);
301 
302  secondary->SetTrackStatus(fAlive);
303 #ifdef G4VERBOSE
304  if(fVerbose)
305  {
306  G4cout << "Product : " << product->GetName() << G4endl;
307  }
308 #endif
309  // add the secondary track in the List
310  aParticleChange.G4VParticleChange::AddSecondary(secondary);
311  }
312 #ifdef G4VERBOSE
313  if(fVerbose) G4cout << "-------------" << G4endl;
314 #endif
315  }
316  //DEBUG
317  else if(fVerbose && decayEnergy)
318  {
319  G4cout << "No products for this channel" << G4endl;
320  G4cout<<"-------------"<<G4endl;
321  }
322  /*
323  else if(!decayEnergy && !nbProducts)
324  {
325  G4ExceptionDescription errMsg;
326  errMsg << "There is no products and no energy specified in the molecular "
327  "dissociation channel";
328  G4Exception("G4MolecularDissociationProcess::DecayIt",
329  "DNAMolecularDissociation002",
330  FatalErrorInArgument,
331  errMsg);
332  }
333  */
334  }
335 
337 
338  return &aParticleChange;
339 }
void SetTrackStatus(const G4TrackStatus aTrackStatus)
const std::vector< const G4MolecularDissociationChannel * > * GetDecayChannel() const
Definition: G4Molecule.cc:469
virtual G4VSolid * GetSolid(G4int depth=0) const
Definition: G4VTouchable.cc:51
std::ostringstream G4ExceptionDescription
Definition: globals.hh:76
const G4ThreeVector & GetPosition() const
int G4int
Definition: G4Types.hh:78
void ProposeLocalEnergyDeposit(G4double anEnergyPart)
const G4String & GetName() const
Definition: G4Molecule.cc:356
static G4ITTransportationManager * GetTransportationManager()
#define G4UniformRand()
Definition: Randomize.hh:97
G4GLOB_DLL std::ostream G4cout
const G4String & GetName() const
virtual EInside Inside(const G4ThreeVector &p) const =0
virtual std::vector< G4ThreeVector > GetProductsDisplacement(const G4MolecularDissociationChannel *) const =0
tuple navigator
Definition: write_gdml.py:35
G4int GetTrackID() const
G4double GetGlobalTime() const
G4Molecule * GetMolecule(const G4Track &track)
Definition: G4Molecule.cc:76
const G4MoleculeDefinition * GetDefinition() const
Definition: G4Molecule.cc:546
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
G4Track * BuildTrack(G4double globalTime, const G4ThreeVector &Position)
Definition: G4Molecule.cc:391
virtual void Initialize(const G4Track &)
G4ThreeVector TransformPoint(const G4ThreeVector &vec) const
virtual G4VPhysicalVolume * GetVolume(G4int depth=0) const
Definition: G4VTouchable.cc:44
void SetNumberOfSecondaries(G4int totSecondaries)
const G4VTouchable * GetTouchable() const
G4ITNavigator * GetNavigatorForTracking() const
void PrintState() const
Definition: G4Molecule.cc:384
G4ParticleChange aParticleChange
Definition: G4VProcess.hh:289
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
#define G4endl
Definition: G4ios.hh:61
virtual G4ThreeVector GetMotherMoleculeDisplacement(const G4MolecularDissociationChannel *) const =0
double G4double
Definition: G4Types.hh:76
void ProposeTrackStatus(G4TrackStatus status)
double mag() const
#define DBL_MAX
Definition: templates.hh:83
const G4MolecularDissociationTable * GetDecayTable() const

Here is the call graph for this function:

Here is the caller graph for this function:

G4VMolecularDecayDisplacer * G4DNAMolecularDissociation::GetDecayDisplacer ( const G4ParticleDefinition molDef)

Definition at line 354 of file G4DNAMolecularDissociation.cc.

355 {
356  return fDisplacementMap[molDef];
357 }
G4VMolecularDecayDisplacer * G4DNAMolecularDissociation::GetDisplacer ( const G4ParticleDefinition molDef)

Definition at line 372 of file G4DNAMolecularDissociation.cc.

373 {
374  return fDisplacementMap[molDef];
375 }
virtual G4double G4DNAMolecularDissociation::GetMeanFreePath ( const G4Track ,
G4double  ,
G4ForceCondition  
)
inlineprotectedvirtual

Implements G4VITRestDiscreteProcess.

Definition at line 116 of file G4DNAMolecularDissociation.hh.

119  {
120  return 0;
121  }
G4double G4DNAMolecularDissociation::GetMeanLifeTime ( const G4Track track,
G4ForceCondition  
)
protectedvirtual

Implements G4VITRestDiscreteProcess.

Definition at line 132 of file G4DNAMolecularDissociation.cc.

134 {
135  G4double output = GetMolecule(track)->GetDecayTime() - track.GetProperTime();
136  return (output > 0 ? output : 0);
137 }
G4double GetProperTime() const
G4Molecule * GetMolecule(const G4Track &track)
Definition: G4Molecule.cc:76
double G4double
Definition: G4Types.hh:76
G4double GetDecayTime() const
Definition: G4Molecule.cc:497

Here is the call graph for this function:

Here is the caller graph for this function:

G4bool G4DNAMolecularDissociation::IsApplicable ( const G4ParticleDefinition aParticleType)
virtual

Reimplemented from G4VProcess.

Definition at line 109 of file G4DNAMolecularDissociation.cc.

110 {
111  if(aParticleType.GetParticleType() == "Molecule")
112  {
113 #ifdef G4VERBOSE
114 
115  if(fVerbose > 1)
116  {
117  G4cout << "G4MolecularDissociation::IsApplicable(";
118  G4cout << aParticleType.GetParticleName() << ",";
119  G4cout << aParticleType.GetParticleType() << ")" << G4endl;
120  }
121 #endif
122  return (true);
123  }
124  else
125  {
126  return false;
127  }
128 }
const G4String & GetParticleName() const
G4GLOB_DLL std::ostream G4cout
const G4String & GetParticleType() const
#define G4endl
Definition: G4ios.hh:61

Here is the call graph for this function:

G4VParticleChange* G4DNAMolecularDissociation::PostStepDoIt ( const G4Track track,
const G4Step step 
)
inlinevirtual

Reimplemented from G4VITRestDiscreteProcess.

Definition at line 93 of file G4DNAMolecularDissociation.hh.

97  {
98  return AtRestDoIt(track, step);
99  }
G4VParticleChange * AtRestDoIt(const G4Track &track, const G4Step &step)

Here is the call graph for this function:

G4double G4DNAMolecularDissociation::PostStepGetPhysicalInteractionLength ( const G4Track track,
G4double  previousStepSize,
G4ForceCondition condition 
)
virtual

Reimplemented from G4VITRestDiscreteProcess.

Definition at line 380 of file G4DNAMolecularDissociation.cc.

383 {
384  return 0; //-1*(-log(1-G4UniformRand())*100*1e-15*s);
385 }
void G4DNAMolecularDissociation::SetDecayDisplacer ( const G4ParticleDefinition molDef,
G4VMolecularDecayDisplacer aDisplacer 
)

Definition at line 344 of file G4DNAMolecularDissociation.cc.

346 {
347  fDisplacementMap[molDef] = aDisplacer;
348 }

Here is the caller graph for this function:

void G4DNAMolecularDissociation::SetDisplacer ( const G4ParticleDefinition molDef,
G4VMolecularDecayDisplacer aDisplacer 
)

Definition at line 362 of file G4DNAMolecularDissociation.cc.

364 {
365  fDisplacementMap[molDef] = aDisplacer;
366 }
void G4DNAMolecularDissociation::SetVerbose ( G4int  verbose)
inline

Definition at line 137 of file G4DNAMolecularDissociation.hh.

138 {
139  fVerbose = verbose ;
140 }

The documentation for this class was generated from the following files: