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

#include <G4VPhononProcess.hh>

Inheritance diagram for G4VPhononProcess:
Collaboration diagram for G4VPhononProcess:

Public Member Functions

 G4VPhononProcess (const G4String &processName)
 
virtual ~G4VPhononProcess ()
 
virtual G4bool IsApplicable (const G4ParticleDefinition &aPD)
 
virtual void StartTracking (G4Track *track)
 
virtual void EndTracking ()
 
- Public Member Functions inherited from G4VDiscreteProcess
 G4VDiscreteProcess (const G4String &, G4ProcessType aType=fNotDefined)
 
 G4VDiscreteProcess (G4VDiscreteProcess &)
 
virtual ~G4VDiscreteProcess ()
 
virtual G4double PostStepGetPhysicalInteractionLength (const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)
 
virtual G4VParticleChangePostStepDoIt (const G4Track &, const G4Step &)
 
virtual G4double AlongStepGetPhysicalInteractionLength (const G4Track &, G4double, G4double, G4double &, G4GPILSelection *)
 
virtual G4double AtRestGetPhysicalInteractionLength (const G4Track &, G4ForceCondition *)
 
virtual G4VParticleChangeAtRestDoIt (const G4Track &, const G4Step &)
 
virtual G4VParticleChangeAlongStepDoIt (const G4Track &, const G4Step &)
 
- 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 BuildPhysicsTable (const G4ParticleDefinition &)
 
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 SetProcessManager (const G4ProcessManager *)
 
virtual const G4ProcessManagerGetProcessManager ()
 
virtual void ResetNumberOfInteractionLengthLeft ()
 
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 G4int GetPolarization (const G4Track &track) const
 
virtual G4int GetPolarization (const G4Track *track) const
 
virtual G4int ChoosePolarization (G4double Ldos, G4double STdos, G4double FTdos) const
 
virtual G4TrackCreateSecondary (G4int polarization, const G4ThreeVector &K, G4double energy) const
 
- Protected Member Functions inherited from G4VDiscreteProcess
virtual G4double GetMeanFreePath (const G4Track &aTrack, G4double previousStepSize, G4ForceCondition *condition)=0
 
- Protected Member Functions inherited from G4VProcess
void SubtractNumberOfInteractionLengthLeft (G4double previousStepSize)
 
void ClearNumberOfInteractionLengthLeft ()
 

Protected Attributes

G4PhononTrackMaptrackKmap
 
const G4LatticePhysicaltheLattice
 
- 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
 

Additional Inherited Members

- Static Public Member Functions inherited from G4VProcess
static const G4StringGetProcessTypeName (G4ProcessType)
 

Detailed Description

Definition at line 42 of file G4VPhononProcess.hh.

Constructor & Destructor Documentation

G4VPhononProcess::G4VPhononProcess ( const G4String processName)

Definition at line 54 of file G4VPhononProcess.cc.

55  : G4VDiscreteProcess(processName, fPhonon),
57  currentTrack(0) {}
static G4PhononTrackMap * GetInstance()
const G4LatticePhysical * theLattice
G4PhononTrackMap * trackKmap
G4VPhononProcess::~G4VPhononProcess ( )
virtual

Definition at line 59 of file G4VPhononProcess.cc.

59 {;}

Member Function Documentation

G4int G4VPhononProcess::ChoosePolarization ( G4double  Ldos,
G4double  STdos,
G4double  FTdos 
) const
protectedvirtual

Definition at line 104 of file G4VPhononProcess.cc.

105  {
106  G4double norm = Ldos + STdos + FTdos;
107  G4double cProbST = STdos/norm;
108  G4double cProbFT = FTdos/norm + cProbST;
109 
110  // NOTE: Order of selection done to match previous random sequences
111  G4double modeMixer = G4UniformRand();
112  if (modeMixer<cProbST) return G4PhononPolarization::TransSlow;
113  if (modeMixer<cProbFT) return G4PhononPolarization::TransFast;
115 }
#define G4UniformRand()
Definition: Randomize.hh:97
double G4double
Definition: G4Types.hh:76

Here is the caller graph for this function:

G4Track * G4VPhononProcess::CreateSecondary ( G4int  polarization,
const G4ThreeVector K,
G4double  energy 
) const
protectedvirtual

Definition at line 120 of file G4VPhononProcess.cc.

122  {
123  if (verboseLevel>1) {
124  G4cout << GetProcessName() << " CreateSecondary pol " << polarization
125  << " K " << waveVec << " E " << energy << G4endl;
126  }
127 
128  G4ThreeVector vgroup = theLattice->MapKtoVDir(polarization, waveVec);
129  if (verboseLevel>1) G4cout << " MapKtoVDir returned " << vgroup << G4endl;
130 
131  vgroup = theLattice->RotateToGlobal(vgroup);
132  if (verboseLevel>1) G4cout << " RotateToGlobal returned " << vgroup << G4endl;
133 
134  if (verboseLevel && std::fabs(vgroup.mag()-1.) > 0.01) {
135  G4cout << "WARNING: " << GetProcessName() << " vgroup not a unit vector: "
136  << vgroup << G4endl;
137  }
138 
139  G4ParticleDefinition* thePhonon = G4PhononPolarization::Get(polarization);
140 
141  // Secondaries are created at the current track coordinates
142  G4Track* sec = new G4Track(new G4DynamicParticle(thePhonon, vgroup, energy),
143  currentTrack->GetGlobalTime(),
144  currentTrack->GetPosition());
145 
146  // Store wavevector in lookup table for future tracking
147  trackKmap->SetK(sec, theLattice->RotateToGlobal(waveVec));
148 
149  if (verboseLevel>1) {
150  G4cout << GetProcessName() << " secondary K rotated to "
151  << trackKmap->GetK(sec) << G4endl;
152  }
153 
154  sec->SetVelocity(theLattice->MapKtoV(polarization, waveVec));
155  sec->UseGivenVelocity(true);
156 
157  return sec;
158 }
const G4LatticePhysical * theLattice
void SetVelocity(G4double val)
G4int verboseLevel
Definition: G4VProcess.hh:368
const G4ThreeVector & GetPosition() const
G4double MapKtoV(G4int, G4ThreeVector) const
G4PhononTrackMap * trackKmap
void SetK(const G4Track *track, const G4ThreeVector &K)
G4GLOB_DLL std::ostream G4cout
const G4ThreeVector & GetK(const G4Track *track) const
G4ThreeVector MapKtoVDir(G4int, G4ThreeVector) const
G4int Get(const G4ParticleDefinition *aPD)
G4double GetGlobalTime() const
G4bool UseGivenVelocity() const
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
G4double energy(const ThreeVector &p, const G4double m)
G4ThreeVector RotateToGlobal(const G4ThreeVector &dir) const
#define G4endl
Definition: G4ios.hh:61
double mag() const

Here is the call graph for this function:

Here is the caller graph for this function:

void G4VPhononProcess::EndTracking ( )
virtual

Reimplemented from G4VProcess.

Definition at line 87 of file G4VPhononProcess.cc.

87  {
88  G4VProcess::EndTracking(); // Apply base class actions
89  trackKmap->RemoveTrack(currentTrack);
90  currentTrack = 0;
91  theLattice = 0;
92 }
const G4LatticePhysical * theLattice
void RemoveTrack(const G4Track *track)
G4PhononTrackMap * trackKmap
virtual void EndTracking()
Definition: G4VProcess.cc:113

Here is the call graph for this function:

G4int G4VPhononProcess::GetPolarization ( const G4Track track) const
protectedvirtual

Definition at line 97 of file G4VPhononProcess.cc.

97  {
99 }
const G4ParticleDefinition * GetParticleDefinition() const
G4int Get(const G4ParticleDefinition *aPD)

Here is the call graph for this function:

Here is the caller graph for this function:

virtual G4int G4VPhononProcess::GetPolarization ( const G4Track track) const
inlineprotectedvirtual

Definition at line 57 of file G4VPhononProcess.hh.

57  {
58  return GetPolarization(*track);
59  }
virtual G4int GetPolarization(const G4Track &track) const

Here is the call graph for this function:

G4bool G4VPhononProcess::IsApplicable ( const G4ParticleDefinition aPD)
virtual

Reimplemented from G4VProcess.

Reimplemented in G4PhononDownconversion.

Definition at line 64 of file G4VPhononProcess.cc.

64  {
65  return (&aPD==G4PhononLong::Definition() ||
68 }
static G4PhononLong * Definition()
Definition: G4PhononLong.cc:40
static G4PhononTransFast * Definition()
static G4PhononTransSlow * Definition()

Here is the call graph for this function:

void G4VPhononProcess::StartTracking ( G4Track track)
virtual

Reimplemented from G4VProcess.

Definition at line 73 of file G4VPhononProcess.cc.

73  {
74  G4VProcess::StartTracking(track); // Apply base class actions
75 
76  // FIXME: THE WAVEVECTOR SHOULD BE COMPUTED BY INVERTING THE K/V MAP
77  if (!trackKmap->Find(track))
78  trackKmap->SetK(track, track->GetMomentumDirection());
79 
80  currentTrack = track; // Save for use by EndTracking
81 
82  // Fetch lattice for current track once, use in subsequent steps
84  theLattice = LM->GetLattice(track->GetVolume());
85 }
const G4LatticePhysical * theLattice
static G4LatticeManager * GetLatticeManager()
G4PhononTrackMap * trackKmap
G4bool Find(const G4Track *track) const
virtual void StartTracking(G4Track *)
Definition: G4VProcess.cc:101
void SetK(const G4Track *track, const G4ThreeVector &K)
const G4ThreeVector & GetMomentumDirection() const
G4VPhysicalVolume * GetVolume() const
G4LatticeLogical * GetLattice(G4Material *) const

Here is the call graph for this function:

Member Data Documentation

const G4LatticePhysical* G4VPhononProcess::theLattice
protected

Definition at line 72 of file G4VPhononProcess.hh.

G4PhononTrackMap* G4VPhononProcess::trackKmap
protected

Definition at line 71 of file G4VPhononProcess.hh.


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