ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4VPhononProcess.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file G4VPhononProcess.cc
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
28 //
29 //
30 // 20131111 Add verbosity to report creating secondaries
31 
32 #include "G4VPhononProcess.hh"
33 #include "G4DynamicParticle.hh"
34 #include "G4ExceptionSeverity.hh"
35 #include "G4LatticeManager.hh"
36 #include "G4LatticePhysical.hh"
37 #include "G4ParticleDefinition.hh"
38 #include "G4PhononLong.hh"
39 #include "G4PhononPolarization.hh"
40 #include "G4PhononTrackMap.hh"
41 #include "G4PhononTransFast.hh"
42 #include "G4PhononTransSlow.hh"
43 #include "G4ProcessType.hh"
44 #include "G4ThreeVector.hh"
45 #include "G4Track.hh"
46 
47 namespace {
48  const G4ThreeVector nullVec(0.,0.,0.); // For convenience below
49 }
50 
51 // Constructor and destructor
52 
54  : G4VDiscreteProcess(processName, fPhonon),
55  trackKmap(G4PhononTrackMap::GetInstance()), theLattice(0),
56  currentTrack(0) {}
57 
59 
60 
61 // Only applies to the known phonon polarization states
62 
64  return (&aPD==G4PhononLong::Definition() ||
67 }
68 
69 
70 // Initialize wave vectors for currently active track(s)
71 
73  G4VProcess::StartTracking(track); // Apply base class actions
74 
75  // FIXME: THE WAVEVECTOR SHOULD BE COMPUTED BY INVERTING THE K/V MAP
76  if (!trackKmap->Find(track))
77  trackKmap->SetK(track, track->GetMomentumDirection());
78 
79  currentTrack = track; // Save for use by EndTracking
80 
81  // Fetch lattice for current track once, use in subsequent steps
83  theLattice = LM->GetLattice(track->GetVolume());
84 }
85 
87  G4VProcess::EndTracking(); // Apply base class actions
89  currentTrack = 0;
90  theLattice = 0;
91 }
92 
93 
94 // For convenience, map phonon type to polarization code
95 
98 }
99 
100 
101 // Generate random polarization from density of states
102 
104  G4double FTdos) const {
105  G4double norm = Ldos + STdos + FTdos;
106  G4double cProbST = STdos/norm;
107  G4double cProbFT = FTdos/norm + cProbST;
108 
109  // NOTE: Order of selection done to match previous random sequences
110  G4double modeMixer = G4UniformRand();
111  if (modeMixer<cProbST) return G4PhononPolarization::TransSlow;
112  if (modeMixer<cProbFT) return G4PhononPolarization::TransFast;
114 }
115 
116 
117 // Create new secondary track from phonon configuration
118 
120  const G4ThreeVector& waveVec,
121  G4double energy) const {
122  if (verboseLevel>1) {
123  G4cout << GetProcessName() << " CreateSecondary pol " << polarization
124  << " K " << waveVec << " E " << energy << G4endl;
125  }
126 
127  G4ThreeVector vgroup = theLattice->MapKtoVDir(polarization, waveVec);
128  if (verboseLevel>1) G4cout << " MapKtoVDir returned " << vgroup << G4endl;
129 
130  vgroup = theLattice->RotateToGlobal(vgroup);
131  if (verboseLevel>1) G4cout << " RotateToGlobal returned " << vgroup << G4endl;
132 
133  if (verboseLevel && std::fabs(vgroup.mag()-1.) > 0.01) {
134  G4cout << "WARNING: " << GetProcessName() << " vgroup not a unit vector: "
135  << vgroup << G4endl;
136  }
137 
138  G4ParticleDefinition* thePhonon = G4PhononPolarization::Get(polarization);
139 
140  // Secondaries are created at the current track coordinates
141  G4Track* sec = new G4Track(new G4DynamicParticle(thePhonon, vgroup, energy),
144 
145  // Store wavevector in lookup table for future tracking
146  trackKmap->SetK(sec, theLattice->RotateToGlobal(waveVec));
147 
148  if (verboseLevel>1) {
149  G4cout << GetProcessName() << " secondary K rotated to "
150  << trackKmap->GetK(sec) << G4endl;
151  }
152 
153  sec->SetVelocity(theLattice->MapKtoV(polarization, waveVec));
154  sec->UseGivenVelocity(true);
155 
156  return sec;
157 }