Geant4  10.00.p02
XrayFluoDataSet.cc
Go to the documentation of this file.
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 //
27 // $Id: XrayFluoDAtaSet.cc
28 // GEANT4 tag $Name: xray_fluo-V03-02-00
29 //
30 // Author: Elena Guardincerri (Elena.Guardincerri@ge.infn.it)
31 //
32 // History:
33 // -----------
34 // 28 Nov 2001 Elena Guardincerri Created
35 //
36 // -------------------------------------------------------------------
37 
38 #include "XrayFluoDataSet.hh"
39 #include <fstream>
40 #include "G4VDataSetAlgorithm.hh"
41 
43  G4DataVector* points,
44  G4DataVector* values,
45  const G4VDataSetAlgorithm* interpolation,
46  G4double unitE, G4double unitData)
47  :energies(points), data(values), algorithm(interpolation)
48 {
49  numberOfBins = energies->size();
50  unit1 = unitE;
51  unit2 = unitData;
52 
53  G4cout << "XrayFluo FluoDataSet created" << G4endl;
54 }
55 
57  const G4String& dataFile,
58  const G4VDataSetAlgorithm* interpolation,
59  G4double unitE, G4double unitData)
60  : algorithm(interpolation)
61 {
62  energies = new G4DataVector;
63  data = new G4DataVector;
64  unit1 = unitE;
65  unit2 = unitData;
66  LoadData(dataFile);
67  numberOfBins = energies->size();
68 
69  G4cout << "XrayFluo FluoDataSet created" << G4endl;
70 
71 }
72 
73 
74 // Destructor
75 
77 {
78  delete energies;
79  delete data;
80  G4cout << "XrayFluo FluoDataSet deleted" << G4endl;
81 }
82 
83 
85 {
86  G4double value;
87  G4double e0 = (*energies)[0];
88  // Protections
89  size_t bin = FindBinLocation(e);
90  if (bin == numberOfBins)
91  {
92 
93  value = (*data)[bin];
94  }
95  else if (e <= e0)
96  {
97 
98  value = (*data)[0];
99  }
100  else
101  {
102  value = algorithm->Calculate(e,bin,*energies,*data);
103  }
104 
105  return value;
106 }
107 
109 {
110  // Protection against call outside allowed range
111  G4double e0 = (*energies)[0];
112  if (energy < e0)
113  {
114 
115  energy = e0;
116  }
117 
118  size_t lowerBound = 0;
119  size_t upperBound = numberOfBins - 1;
120 
121  // Binary search
122  while (lowerBound <= upperBound)
123  {
124  size_t midBin = (lowerBound + upperBound)/2;
125  if ( energy < (*energies)[midBin] ) upperBound = midBin-1;
126  else lowerBound = midBin+1;
127  }
128 
129  return upperBound;
130 }
131 
132 
134 {
135  // Build the complete string identifying the file with the data set
136 
137 
138  G4String dirFile = "";
139 
140  char* path;
141 
142 #ifndef XRAYDATA
143 #define XRAYDATA PWD
144 #endif
145 
146  path = getenv("XRAYDATA");
147 
148  G4cout << path << G4endl;
149  G4cout << fileName << G4endl;
150 
151 
152  G4String pathString(path);
153  dirFile = pathString + "/" + fileName + ".dat";
154 
155  std::ifstream file(dirFile);
156  std::filebuf* lsdp = file.rdbuf();
157 
158  if (! (lsdp->is_open()) )
159  {
161  execp << "XrayFluoDataSet - data file: " + dirFile + " not found"<<G4endl;
162  G4Exception("XrayFluoDataSet::LoadData()","example-xray_fluorescence01",
163  FatalException, execp);
164  }
165  G4double a = 0;
166  G4int k = 1;
167 
168  do
169  {
170  file >> a;
171  G4int nColumns = 2;
172  // The file is organized into two columns:
173  // 1st column is the energy
174  // 2nd column is the corresponding value
175  // The file terminates with the pattern: -1 -1
176  // -2 -2
177  if (a == -1 || a == -2)
178  {
179 
180  }
181  else
182  {
183  if (k%nColumns != 0)
184  {
185  G4double e = a * unit1;
186  energies->push_back(e);
187 
188  k++;
189 
190  }
191  else if (k%nColumns == 0)
192  {
193  G4double value = a * unit2;
194  data->push_back(value);
195 
196  k = 1;
197  }
198  }
199 
200  } while (a != -2); // end of file
201 
202  file.close();
203  return true;
204 
205 }
207 {
208  size_t size = numberOfBins;
209  for (size_t i=0; i<size; i++)
210  {
211  G4double e = (*energies)[i] / unit1;
212  G4double sigma = (*data)[i] / unit2 ;
213  G4cout << "Point: "
214  << e
215  << " - Data value : "
216  << sigma
217  << G4endl;
218  }
219 }
220 
221 
222 
XrayFluoDataSet(G4int Z, G4DataVector *points, G4DataVector *values, const G4VDataSetAlgorithm *interpolation, G4double unitE=CLHEP::MeV, G4double unitData=CLHEP::barn)
std::ostringstream G4ExceptionDescription
Definition: globals.hh:76
G4double a
Definition: TRTMaterials.hh:39
int G4int
Definition: G4Types.hh:78
virtual G4double Calculate(G4double point, G4int bin, const G4DataVector &energies, const G4DataVector &data) const =0
void PrintData() const
G4GLOB_DLL std::ostream G4cout
bool G4bool
Definition: G4Types.hh:79
G4bool LoadData(const G4String &dataFile)
G4int FindBinLocation(G4double energy) const
const G4VDataSetAlgorithm * algorithm
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
G4double energy(const ThreeVector &p, const G4double m)
G4double FindValue(G4double e, G4int) const
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76
G4DataVector * data
G4DataVector * energies