ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4VDNAModel.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file G4VDNAModel.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 // Authors: S. Meylan and C. Villagrasa (IRSN, France)
27 // This class is used to support PTB models that come from
28 // M. Bug et al, Rad. Phys and Chem. 130, 459-479 (2017)
29 //
30 
31 #include "G4VDNAModel.hh"
32 #include "G4SystemOfUnits.hh"
33 #include "G4ParticleTable.hh"
34 
35 G4VDNAModel::G4VDNAModel(const G4String &nam, const G4String &applyToMaterial)
36  : fStringOfMaterials(applyToMaterial), fName(nam)
37 {
38 
39 }
40 
42 {
43  // Clean fTableData
44  std::map<G4String, std::map<G4String,G4DNACrossSectionDataSet*,std::less<G4String> > >::iterator posOuter;
45  std::map<G4String,G4DNACrossSectionDataSet*,std::less<G4String> >::iterator posInner;
46  // iterate on each material
47  for (posOuter = fTableData.begin(); posOuter != fTableData.end(); ++posOuter)
48  {
49  // iterate on each particle
50  for(posInner = posOuter->second.begin(); posInner != posOuter->second.end(); ++posInner)
51  {
52  G4DNACrossSectionDataSet* table = posInner->second;
53  if(table != 0) delete table;
54  }
55  }
56 }
57 
58 void G4VDNAModel::AddCrossSectionData(G4String materialName, G4String particleName, G4String fileCS, G4String fileDiffCS, G4double scaleFactor)
59 {
60  fModelMaterials.push_back(materialName);
61  fModelParticles.push_back(particleName);
62  fModelCSFiles.push_back(fileCS);
63  fModelDiffCSFiles.push_back(fileDiffCS);
64  fModelScaleFactors.push_back(scaleFactor);
65 }
66 
67 void G4VDNAModel::AddCrossSectionData(G4String materialName, G4String particleName, G4String fileCS, G4double scaleFactor)
68 {
69  fModelMaterials.push_back(materialName);
70  fModelParticles.push_back(particleName);
71  fModelCSFiles.push_back(fileCS);
72  fModelScaleFactors.push_back(scaleFactor);
73 }
74 
76 {
77  G4String fileElectron, fileDiffElectron;
78  G4String materialName, modelParticleName;
79  G4double scaleFactor;
80 
81  // construct applyToMatVect with materials specified by the user
82  std::vector<G4String> applyToMatVect = BuildApplyToMatVect(fStringOfMaterials);
83 
84  // iterate on each material contained into the fStringOfMaterials variable (through applyToMatVect)
85  for(unsigned int i=0;i<applyToMatVect.size();++i)
86  {
87  // We have selected a material coming from applyToMatVect
88  // We try to find if this material correspond to a model registered material
89  // If it is, then isMatFound becomes true
90  G4bool isMatFound = false;
91 
92  // We iterate on each model registered materials to load the CS data
93  // We have to do a for loop because of the "all" option
94  // applyToMatVect[i] == "all" implies applyToMatVect.size()=1 and we want to iterate on all registered materials
95  for(unsigned int j=0;j<fModelMaterials.size();++j)
96  {
97  if(applyToMatVect[i] == fModelMaterials[j] || applyToMatVect[i] == "all")
98  {
99  isMatFound = true;
100  materialName = fModelMaterials[j];
101  modelParticleName = fModelParticles[j];
102  fileElectron = fModelCSFiles[j];
103  if(!fModelDiffCSFiles.empty()) fileDiffElectron = fModelDiffCSFiles[j];
104  scaleFactor = fModelScaleFactors[j];
105 
106  ReadAndSaveCSFile(materialName, modelParticleName, fileElectron, scaleFactor);
107 
108  if(!fModelDiffCSFiles.empty()) ReadDiffCSFile(materialName, modelParticleName, fileDiffElectron, scaleFactor);
109 
110  }
111  }
112 
113  // check if we found a correspondance, if not: fatal error
114  if(!isMatFound)
115  {
116  std::ostringstream oss;
117  oss << applyToMatVect[i] << " material was not found. It means the material specified in the UserPhysicsList is not a model material for ";
118  oss << particleName;
119  G4Exception("G4VDNAModel::LoadCrossSectionData","em0003",
120  FatalException, oss.str().c_str());
121  return;
122  }
123  }
124 }
125 
126 void G4VDNAModel::ReadDiffCSFile(const G4String&, const G4String&, const G4String&, const G4double)
127 {
128  G4String text("ReadDiffCSFile must be implemented in the model class using a differential cross section data file");
129 
130  G4Exception("G4VDNAModel::ReadDiffCSFile","em0003",
131  FatalException, text);
132 }
133 
134 void G4VDNAModel::EnableForMaterialAndParticle(const G4String &materialName, const G4String &particleName)
135 {
136  fTableData[materialName][particleName] = 0;
137 }
138 
140 {
141  // output material vector
142  std::vector<G4String> materialVect;
143 
144  // if we don't find any "/" then it means we only have one "material" (could be the "all" option)
145  if(materials.find("/")==std::string::npos)
146  {
147  // we add the material to the output vector
148  materialVect.push_back(materials);
149  }
150  // if we have several materials listed in the string then we must retrieve them
151  else
152  {
153  G4String materialsNonIdentified = materials;
154 
155  while(materialsNonIdentified.find_first_of("/") != std::string::npos)
156  {
157  // we select the first material and stop at the "/" caracter
158  G4String mat = materialsNonIdentified.substr(0, materialsNonIdentified.find_first_of("/"));
159  materialVect.push_back(mat);
160 
161  // we remove the previous material from the materialsNonIdentified string
162  materialsNonIdentified = materialsNonIdentified.substr(materialsNonIdentified.find_first_of("/")+1,
163  materialsNonIdentified.size()-materialsNonIdentified.find_first_of("/"));
164  }
165 
166  // we don't find "/" anymore, it means we only have one material string left
167  // we get it
168  materialVect.push_back(materialsNonIdentified);
169  }
170 
171  return materialVect;
172 }
173 
174 void G4VDNAModel::ReadAndSaveCSFile(const G4String& materialName,
175  const G4String& particleName,
176  const G4String& file, G4double scaleFactor)
177 {
178  fTableData[materialName][particleName] = new G4DNACrossSectionDataSet(new G4LogLogInterpolation, eV, scaleFactor);
179  fTableData[materialName][particleName]->LoadData(file);
180 }
181 
183 {
184  G4int level = 0;
185 
186  TableMapData* tableData = GetTableData();
187 
188  std::map< G4String,G4DNACrossSectionDataSet*,std::less<G4String> >::iterator pos;
189  pos = (*tableData)[materialName].find(particle);
190 
191  if (pos != (*tableData)[materialName].end())
192  {
193  G4DNACrossSectionDataSet* table = pos->second;
194 
195  if (table != 0)
196  {
197  G4double* valuesBuffer = new G4double[table->NumberOfComponents()];
198  const size_t n(table->NumberOfComponents());
199  size_t i(n);
200  G4double value = 0.;
201 
202  while (i>0)
203  {
204  i--;
205  valuesBuffer[i] = table->GetComponent(i)->FindValue(k);
206  value += valuesBuffer[i];
207  }
208 
209  value *= G4UniformRand();
210 
211  i = n;
212 
213  while (i > 0)
214  {
215  i--;
216 
217  if (valuesBuffer[i] > value)
218  {
219  delete[] valuesBuffer;
220  return i;
221  }
222  value -= valuesBuffer[i];
223  }
224 
225  if (valuesBuffer) delete[] valuesBuffer;
226 
227  }
228  }
229  else
230  {
231  G4Exception("G4VDNAModel::RandomSelectShell","em0002",
232  FatalException,"Model not applicable to particle type.");
233  }
234  return level;
235 }
236 
238 {
239  // Check if the given material is defined in the simulation
240 
241  G4bool exist (false);
242 
243  double matTableSize = G4Material::GetMaterialTable()->size();
244 
245  for(int i=0;i<matTableSize;i++)
246  {
247  if(materialName == G4Material::GetMaterialTable()->at(i)->GetName())
248  {
249  exist = true;
250  return exist;
251  }
252  }
253 
254  return exist;
255 }
256 
258 {
259  // Check if the given material is defined in the current model class
260 
261  if (fTableData.find(materialName) == fTableData.end())
262  {
263  return false;
264  }
265  else
266  {
267  return true;
268  }
269 }
270 
272 {
273  // To check two things:
274  // 1- is the material existing in model ?
275  // 2- if yes, is the particle defined for that material ?
276 
277  if(IsMaterialExistingInModel(materialName))
278  {
279  if (fTableData[materialName].find(particleName) == fTableData[materialName].end())
280  {
281  return false;
282  }
283  else return true;
284  }
285  else return false;
286 }