ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
WLSPrimaryGeneratorAction.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file WLSPrimaryGeneratorAction.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 //
26 //
29 //
30 //
31 #include "G4ios.hh"
32 #include "G4Event.hh"
33 
35 
36 #include "G4Material.hh"
38 
39 #include "G4ParticleTable.hh"
40 #include "G4ParticleDefinition.hh"
41 
42 #include "G4PhysicsTable.hh"
43 
44 #include "Randomize.hh"
45 
47 
50 
51 #include "G4SystemOfUnits.hh"
52 
53 #include "G4AutoLock.hh"
54 
55 namespace {
56  G4Mutex gen_mutex = G4MUTEX_INITIALIZER;
57 }
58 
59 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
60 
62 
65 {
66  fDetector = dc;
67  fIntegralTable = NULL;
68 
70 
72 
73  // G4String particleName;
74  // G4ParticleTable* particleTable = G4ParticleTable::GetParticleTable();
75 
76  fTimeConstant = 0.;
77 
78 // fParticleGun->SetParticleDefinition(particleTable->
79 // FindParticle(particleName="opticalphoton"));
80 }
81 
82 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
83 
85 {
86  delete fParticleGun;
87  delete fGunMessenger;
88  if (fIntegralTable) {
90  delete fIntegralTable;
91  }
92 }
93 
94 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
95 
97 {
99 }
100 
101 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
102 
104 {
105  if (fIntegralTable) return;
106 
107  const G4MaterialTable* theMaterialTable = G4Material::GetMaterialTable();
108 
109  G4int numOfMaterials = G4Material::GetNumberOfMaterials();
110 
111  if(!fIntegralTable)fIntegralTable = new G4PhysicsTable(numOfMaterials);
112 
113  for (G4int i=0 ; i < numOfMaterials; i++) {
114 
115  G4PhysicsOrderedFreeVector* aPhysicsOrderedFreeVector =
117 
118  G4Material* aMaterial = (*theMaterialTable)[i];
119 
120  G4MaterialPropertiesTable* aMaterialPropertiesTable =
121  aMaterial->GetMaterialPropertiesTable();
122 
123  if (aMaterialPropertiesTable) {
124  G4MaterialPropertyVector* theWLSVector =
125  aMaterialPropertiesTable->GetProperty("WLSCOMPONENT");
126 
127  if (theWLSVector) {
128  G4double currentIN = (*theWLSVector)[0];
129  if (currentIN >= 0.0) {
130  G4double currentPM = theWLSVector->Energy(0);
131  G4double currentCII = 0.0;
132  aPhysicsOrderedFreeVector->
133  InsertValues(currentPM , currentCII);
134  G4double prevPM = currentPM;
135  G4double prevCII = currentCII;
136  G4double prevIN = currentIN;
137 
138  for (size_t j = 1;
139  j < theWLSVector->GetVectorLength();
140  j++)
141  {
142  currentPM = theWLSVector->Energy(j);
143  currentIN = (*theWLSVector)[j];
144  currentCII = 0.5 * (prevIN + currentIN);
145  currentCII = prevCII + (currentPM - prevPM) * currentCII;
146  aPhysicsOrderedFreeVector->
147  InsertValues(currentPM, currentCII);
148  prevPM = currentPM;
149  prevCII = currentCII;
150  prevIN = currentIN;
151  }
152  }
153  }
154  }
155  fIntegralTable->insertAt(i,aPhysicsOrderedFreeVector);
156  }
157 }
158 
159 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
160 
162 {
163  if (!fFirst) {
164  fFirst = true;
166  }
167 
168 #ifdef use_sampledEnergy
169  const G4MaterialTable* theMaterialTable = G4Material::GetMaterialTable();
170 
171  G4double sampledEnergy = 3*eV;
172 
173  for (size_t j=0 ; j<theMaterialTable->size() ; j++) {
174  G4Material* fMaterial = (*theMaterialTable)[j];
175  if (fMaterial->GetName() == "PMMA" ) {
176  G4MaterialPropertiesTable* aMaterialPropertiesTable =
177  fMaterial->GetMaterialPropertiesTable();
178  const G4MaterialPropertyVector* WLSIntensity =
179  aMaterialPropertiesTable->GetProperty("WLSCOMPONENT");
180 
181  if (WLSIntensity) {
182  G4int MaterialIndex = fMaterial->GetIndex();
183  G4PhysicsOrderedFreeVector* WLSIntegral =
184  (G4PhysicsOrderedFreeVector*)((*fIntegralTable)(MaterialIndex));
185 
186  G4double CIImax = WLSIntegral->GetMaxValue();
187  G4double CIIvalue = G4UniformRand()*CIImax;
188 
189  sampledEnergy = WLSIntegral->GetEnergy(CIIvalue);
190  }
191  }
192  }
193 
194  //fParticleGun->SetParticleEnergy(sampledEnergy);
195 #endif
196 
197  //The code behing this line is not thread safe because polarization
198  //and time are randomly selected and GPS properties are global
199  G4AutoLock l(&gen_mutex);
200  if(fParticleGun->GetParticleDefinition()->GetParticleName()=="opticalphoton"){
203  }
204 
206 }
207 
208 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
209 
211 {
212  G4double angle = G4UniformRand() * 360.0*deg;
213  SetOptPhotonPolar(angle);
214 }
215 
216 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
217 
219 {
220  if (fParticleGun->GetParticleDefinition()->GetParticleName()!="opticalphoton")
221  {
222  G4cout << "-> warning from WLSPrimaryGeneratorAction::SetOptPhotonPolar()"
223  << ": the ParticleGun is not an opticalphoton" << G4endl;
224  return;
225  }
226 
227  G4ThreeVector normal (1., 0., 0.);
229  G4ThreeVector product = normal.cross(kphoton);
230  G4double modul2 = product*product;
231 
232  G4ThreeVector e_perpend (0., 0., 1.);
233  if (modul2 > 0.) e_perpend = (1./std::sqrt(modul2))*product;
234  G4ThreeVector e_paralle = e_perpend.cross(kphoton);
235 
236  G4ThreeVector polar = std::cos(angle)*e_paralle + std::sin(angle)*e_perpend;
238 
239 }
240 
241 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
242 
244 {
245  G4double time = -std::log(G4UniformRand())*fTimeConstant;
247 }