| Geant4 Cross Reference |
1 // 1 //
2 // ******************************************* 2 // ********************************************************************
3 // * License and Disclaimer 3 // * License and Disclaimer *
4 // * 4 // * *
5 // * The Geant4 software is copyright of th 5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided 6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License 7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/ 8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. 9 // * include a list of copyright holders. *
10 // * 10 // * *
11 // * Neither the authors of this software syst 11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing fin 12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warran 13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assum 14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file 15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitatio 16 // * for the full disclaimer and the limitation of liability. *
17 // * 17 // * *
18 // * This code implementation is the result 18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboratio 19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distri 20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you ag 21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publicati 22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Sof 23 // * acceptance of all terms of the Geant4 Software license. *
24 // ******************************************* 24 // ********************************************************************
25 // 25 //
>> 26 // $Id: G4LossTableBuilder.cc,v 1.22 2007/01/16 14:30:59 vnivanch Exp $
>> 27 // GEANT4 tag $Name: geant4-08-02-patch-01 $
>> 28 //
26 // ------------------------------------------- 29 // -------------------------------------------------------------------
27 // 30 //
28 // GEANT4 Class file 31 // GEANT4 Class file
29 // 32 //
30 // 33 //
31 // File name: G4LossTableBuilder 34 // File name: G4LossTableBuilder
32 // 35 //
33 // Author: Vladimir Ivanchenko 36 // Author: Vladimir Ivanchenko
34 // 37 //
35 // Creation date: 03.01.2002 38 // Creation date: 03.01.2002
36 // 39 //
37 // Modifications: 40 // Modifications:
38 // 41 //
39 // 23-01-03 V.Ivanchenko Cut per region 42 // 23-01-03 V.Ivanchenko Cut per region
40 // 21-07-04 V.Ivanchenko Fix problem of range 43 // 21-07-04 V.Ivanchenko Fix problem of range for dedx=0
41 // 08-11-04 Migration to new interface of Stor 44 // 08-11-04 Migration to new interface of Store/Retrieve tables (V.Ivanchenko)
42 // 07-12-04 Fix of BuildDEDX table (V.Ivanchen 45 // 07-12-04 Fix of BuildDEDX table (V.Ivanchenko)
43 // 27-03-06 Add bool options isIonisation (V.I 46 // 27-03-06 Add bool options isIonisation (V.Ivanchenko)
44 // 16-01-07 Fill new (not old) DEDX table (V.I 47 // 16-01-07 Fill new (not old) DEDX table (V.Ivanchenko)
45 // 12-02-07 Use G4LPhysicsFreeVector for the i <<
46 // 24-06-09 Removed hidden bin in G4PhysicsVec <<
47 // 48 //
48 // Class Description: 49 // Class Description:
49 // 50 //
50 // ------------------------------------------- 51 // -------------------------------------------------------------------
51 // 52 //
52 //....oooOO0OOooo........oooOO0OOooo........oo 53 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
53 //....oooOO0OOooo........oooOO0OOooo........oo 54 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
54 55
55 #include "G4LossTableBuilder.hh" 56 #include "G4LossTableBuilder.hh"
56 #include "G4SystemOfUnits.hh" <<
57 #include "G4PhysicsTable.hh" 57 #include "G4PhysicsTable.hh"
58 #include "G4PhysicsLogVector.hh" 58 #include "G4PhysicsLogVector.hh"
59 #include "G4PhysicsTableHelper.hh" 59 #include "G4PhysicsTableHelper.hh"
60 #include "G4PhysicsFreeVector.hh" <<
61 #include "G4ProductionCutsTable.hh" <<
62 #include "G4MaterialCutsCouple.hh" <<
63 #include "G4Material.hh" <<
64 #include "G4VEmModel.hh" <<
65 #include "G4ParticleDefinition.hh" <<
66 #include "G4LossTableManager.hh" <<
67 #include "G4EmParameters.hh" <<
68 <<
69 G4bool G4LossTableBuilder::baseMatFlag = false <<
70 std::vector<G4double>* G4LossTableBuilder::the <<
71 std::vector<G4int>* G4LossTableBuilder::the <<
72 std::vector<G4bool>* G4LossTableBuilder::the <<
73 60
74 //....oooOO0OOooo........oooOO0OOooo........oo 61 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
75 62
76 G4LossTableBuilder::G4LossTableBuilder(G4bool << 63 G4LossTableBuilder::G4LossTableBuilder()
77 : isInitializer(master) << 64 {}
78 { <<
79 theParameters = G4EmParameters::Instance(); <<
80 if (nullptr == theFlag) { <<
81 theDensityFactor = new std::vector<G4doubl <<
82 theDensityIdx = new std::vector<G4int>; <<
83 theFlag = new std::vector<G4bool>; <<
84 } <<
85 } <<
86 65
87 //....oooOO0OOooo........oooOO0OOooo........oo 66 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
88 67
89 G4LossTableBuilder::~G4LossTableBuilder() 68 G4LossTableBuilder::~G4LossTableBuilder()
90 { << 69 {}
91 if (isInitializer) { <<
92 delete theDensityFactor; <<
93 delete theDensityIdx; <<
94 delete theFlag; <<
95 theDensityFactor = nullptr; <<
96 theDensityIdx = nullptr; <<
97 theFlag = nullptr; <<
98 } <<
99 } <<
100 <<
101 //....oooOO0OOooo........oooOO0OOooo........oo <<
102 <<
103 const std::vector<G4int>* G4LossTableBuilder:: <<
104 { <<
105 return theDensityIdx; <<
106 } <<
107 <<
108 //....oooOO0OOooo........oooOO0OOooo........oo <<
109 <<
110 const std::vector<G4double>* G4LossTableBuilde <<
111 { <<
112 return theDensityFactor; <<
113 } <<
114 70
115 //....oooOO0OOooo........oooOO0OOooo........oo 71 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
116 72
117 G4bool G4LossTableBuilder::GetFlag(std::size_t << 73 void G4LossTableBuilder::BuildDEDXTable(G4PhysicsTable* dedxTable,
>> 74 const std::vector<G4PhysicsTable*>& list)
118 { 75 {
119 return (idx < theFlag->size()) ? (*theFlag)[ << 76 size_t n_processes = list.size();
120 } << 77 if(1 >= n_processes) return;
121 78
122 //....oooOO0OOooo........oooOO0OOooo........oo << 79 size_t n_vectors = (list[0])->length();
>> 80 if(0 >= n_vectors) return;
123 81
124 G4bool G4LossTableBuilder::GetBaseMaterialFlag << 82 G4bool b;
125 { <<
126 return baseMatFlag; <<
127 } <<
128 83
129 //....oooOO0OOooo........oooOO0OOooo........oo << 84 G4PhysicsVector* pv = (*(list[0]))[0];
>> 85 size_t nbins = pv->GetVectorLength();
>> 86 G4double elow = pv->GetLowEdgeEnergy(0);
>> 87 G4double ehigh = pv->GetLowEdgeEnergy(nbins);
>> 88 for (size_t i=0; i<n_vectors; i++) {
130 89
131 void << 90 pv = new G4PhysicsLogVector(elow, ehigh, nbins);
132 G4LossTableBuilder::BuildDEDXTable(G4PhysicsTa << 91 for (size_t j=0; j<nbins; j++) {
133 const std:: <<
134 { <<
135 InitialiseBaseMaterials(dedxTable); <<
136 std::size_t n_processes = list.size(); <<
137 if(1 >= n_processes) { return; } <<
138 <<
139 std::size_t nCouples = dedxTable->size(); <<
140 //G4cout << "Nproc= " << n_processes << " nC <<
141 // << dedxTable->size() << G4endl; <<
142 if(0 >= nCouples) { return; } <<
143 <<
144 for (std::size_t i=0; i<nCouples; ++i) { <<
145 auto pv0 = static_cast<G4PhysicsLogVector* <<
146 //if (0 == i) G4cout << i << ". pv0=" << p <<
147 if(pv0 == nullptr) { continue; } <<
148 std::size_t npoints = pv0->GetVectorLength <<
149 auto pv = new G4PhysicsLogVector(*pv0); <<
150 for (std::size_t j=0; j<npoints; ++j) { <<
151 G4double dedx = 0.0; 92 G4double dedx = 0.0;
152 for (std::size_t k=0; k<n_processes; ++k << 93 G4double energy = pv->GetLowEdgeEnergy(j);
153 const G4PhysicsVector* pv1 = (*(list[k]))[i] << 94
154 //if (0 == i) G4cout << " " << k << ". p << 95 for (size_t k=0; k<n_processes; k++) {
155 dedx += (*pv1)[j]; << 96 dedx += ((*(list[k]))[i])->GetValue(energy, b);
156 } 97 }
157 pv->PutValue(j, dedx); 98 pv->PutValue(j, dedx);
>> 99 G4PhysicsTableHelper::SetPhysicsVector(dedxTable, i, pv);
158 } 100 }
159 if(splineFlag) { pv->FillSecondDerivatives <<
160 G4PhysicsTableHelper::SetPhysicsVector(ded <<
161 } 101 }
162 //G4cout << "### G4LossTableBuilder::BuildDE <<
163 //G4cout << *dedxTable << G4endl; <<
164 } 102 }
165 103
166 //....oooOO0OOooo........oooOO0OOooo........oo 104 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
167 105
168 void G4LossTableBuilder::BuildRangeTable(const 106 void G4LossTableBuilder::BuildRangeTable(const G4PhysicsTable* dedxTable,
169 G4Phy << 107 G4PhysicsTable* rangeTable,
>> 108 G4bool isIonisation)
170 // Build range table from the energy loss tabl 109 // Build range table from the energy loss table
171 { 110 {
172 //G4cout << "### G4LossTableBuilder::BuildRa << 111 size_t n_vectors = dedxTable->length();
173 //G4cout << *const_cast<G4PhysicsTable*>(ded << 112 if(!n_vectors) return;
174 const std::size_t nCouples = dedxTable->size <<
175 if(0 >= nCouples) { return; } <<
176 <<
177 const std::size_t n = 100; <<
178 const G4double del = 1.0/(G4double)n; <<
179 <<
180 for (std::size_t i=0; i<nCouples; ++i) { <<
181 auto pv = static_cast<G4PhysicsLogVector*> <<
182 if((pv == nullptr) || (isBaseMatActive && <<
183 std::size_t npoints = pv->GetVectorLength( <<
184 std::size_t bin0 = 0; <<
185 G4double elow = pv->Energy(0); <<
186 G4double ehigh = pv->Energy(npoints-1); <<
187 G4double dedx1 = (*pv)[0]; <<
188 <<
189 // protection for specific cases dedx=0 <<
190 if(dedx1 == 0.0) { <<
191 for (std::size_t k=1; k<npoints; ++k) { <<
192 ++bin0; <<
193 elow = pv->Energy(k); <<
194 dedx1 = (*pv)[k]; <<
195 if(dedx1 > 0.0) { break; } <<
196 } <<
197 npoints -= bin0; <<
198 } <<
199 113
200 // initialisation of a new vector << 114 G4bool b;
201 if(npoints < 3) { npoints = 3; } << 115 size_t n = 100;
202 << 116 G4double del = 1.0/(G4double)n;
203 delete (*rangeTable)[i]; << 117
204 G4PhysicsLogVector* v; << 118 for (size_t i=0; i<n_vectors; i++) {
205 if(0 == bin0) { v = new G4PhysicsLogVector << 119
206 else { v = new G4PhysicsLogVector(elow, eh << 120 if (rangeTable->GetFlag(i) || !isIonisation) {
207 << 121 G4PhysicsVector* pv = (*dedxTable)[i];
208 // assumed dedx proportional to beta << 122 size_t nbins = pv->GetVectorLength();
209 G4double energy1 = v->Energy(0); << 123 size_t bin0 = 0;
210 G4double range = 2.*energy1/dedx1; << 124 G4double elow = pv->GetLowEdgeEnergy(0);
211 /* << 125 G4double ehigh = pv->GetLowEdgeEnergy(nbins);
212 G4cout << "New Range vector Npoints=" << v << 126 G4double dedx1 = pv->GetValue(elow, b);
213 << " coupleIdx=" << i << " spline=" << v- << 127
214 << " Elow=" << v->GetMinEnergy() <<" Ehig << 128 if(dedx1 == 0.0) {
215 << " DEDX(Elow)=" << dedx1 << " R(Elow)=" << 129 for (size_t k=1; k<nbins; k++) {
216 */ << 130 bin0++;
217 v->PutValue(0,range); << 131 elow = pv->GetLowEdgeEnergy(k);
218 << 132 dedx1 = pv->GetValue(elow, b);
219 for (std::size_t j=1; j<npoints; ++j) { << 133 if(dedx1 > 0.0) break;
220 << 134 }
221 G4double energy2 = v->Energy(j); << 135 nbins -= bin0;
222 G4double de = (energy2 - energy1) * <<
223 G4double energy = energy2 + de*0.5; <<
224 G4double sum = 0.0; <<
225 std::size_t idx = j - 1; <<
226 for (std::size_t k=0; k<n; ++k) { <<
227 energy -= de; <<
228 dedx1 = pv->Value(energy, idx); <<
229 if(dedx1 > 0.0) { sum += de/dedx1; } <<
230 } 136 }
231 range += sum; <<
232 /* <<
233 if(energy < 10.) <<
234 G4cout << "j= " << j << " e1= " << energy1 < <<
235 << " n= " << n << " range=" << range< <<
236 */ <<
237 v->PutValue(j,range); <<
238 energy1 = energy2; <<
239 } <<
240 if(splineFlag) { v->FillSecondDerivatives( <<
241 G4PhysicsTableHelper::SetPhysicsVector(ran <<
242 } <<
243 //G4cout << "### Range table" << G4endl; <<
244 //G4cout << *rangeTable << G4endl; <<
245 } <<
246 137
247 //....oooOO0OOooo........oooOO0OOooo........oo << 138 G4PhysicsLogVector* v = new G4PhysicsLogVector(elow, ehigh, nbins);
248 139
249 void << 140 G4double range = 2.*elow/dedx1;
250 G4LossTableBuilder::BuildInverseRangeTable(con << 141 //G4double range = elow/dedx1;
251 G4P << 142 v->PutValue(0,range);
252 // Build inverse range table from the energy l << 143 G4double energy1 = elow;
253 { <<
254 std::size_t nCouples = rangeTable->size(); <<
255 if(0 >= nCouples) { return; } <<
256 <<
257 for (std::size_t i=0; i<nCouples; ++i) { <<
258 G4PhysicsVector* pv = (*rangeTable)[i]; <<
259 if((pv == nullptr) || (isBaseMatActive && <<
260 std::size_t npoints = pv->GetVectorLength( <<
261 <<
262 delete (*invRangeTable)[i]; <<
263 auto v = new G4PhysicsFreeVector(npoints, <<
264 <<
265 for (std::size_t j=0; j<npoints; ++j) { <<
266 G4double e = pv->Energy(j); <<
267 G4double r = (*pv)[j]; <<
268 v->PutValues(j,r,e); <<
269 } <<
270 if (splineFlag) { v->FillSecondDerivatives <<
271 v->EnableLogBinSearch(theParameters->Numbe <<
272 144
273 G4PhysicsTableHelper::SetPhysicsVector(inv << 145 for (size_t j=1; j<nbins; j++) {
274 } <<
275 //G4cout << "### Inverse range table" << G4e <<
276 //G4cout << *invRangeTable << G4endl; <<
277 } <<
278 146
279 //....oooOO0OOooo........oooOO0OOooo........oo << 147 G4double energy2 = pv->GetLowEdgeEnergy(j+bin0);
>> 148 G4double dedx2 = pv->GetValue(energy2, b);
>> 149 G4double de = (energy2 - energy1) * del;
>> 150 G4double energy = energy1 - de*0.5;
280 151
281 void G4LossTableBuilder::InitialiseBaseMateria << 152 G4bool yes = true;
282 { << 153 if(dedx1 < DBL_MIN || dedx2 < DBL_MIN) yes = false;
283 if(!isInitializer) { return; } <<
284 const G4ProductionCutsTable* theCoupleTable= <<
285 G4ProductionCutsTable::GetProductionCutsTa <<
286 std::size_t nCouples = theCoupleTable->GetTa <<
287 std::size_t nFlags = theFlag->size(); <<
288 /* <<
289 G4cout << "### InitialiseBaseMaterials: nCou <<
290 << " nFlags=" << nFlags << " isInit:" << is <<
291 << " baseMat:" << baseMatFlag << G4endl; <<
292 */ <<
293 // define base material flag <<
294 if(isBaseMatActive && !baseMatFlag) { <<
295 for(G4int i=0; i<(G4int)nCouples; ++i) { <<
296 if(nullptr != theCoupleTable->GetMateria <<
297 baseMatFlag = true; <<
298 isInitialized = false; <<
299 break; <<
300 } <<
301 } <<
302 } <<
303 154
304 if(nFlags != nCouples) { isInitialized = fal << 155 G4double fac, f;
305 if(isInitialized) { return; } <<
306 156
307 // reserve memory << 157 if(yes) fac = std::log(dedx2/dedx1)/std::log(energy2/energy1);
308 theFlag->resize(nCouples, true); << 158 else fac = (dedx2 - dedx1)/(energy2 - energy1);
309 theDensityFactor->resize(nCouples,1.0); <<
310 theDensityIdx->resize(nCouples, 0); <<
311 <<
312 // define default flag and index of used mat <<
313 for (G4int i=0; i<(G4int)nCouples; ++i) { <<
314 (*theFlag)[i] = (nullptr == table) ? true <<
315 (*theDensityIdx)[i] = i; <<
316 } <<
317 isInitialized = true; <<
318 if (!baseMatFlag) { return; } <<
319 159
320 // use base materials << 160 for (size_t k=0; k<n; k++) {
321 for (G4int i=0; i<(G4int)nCouples; ++i) { << 161 energy += de;
322 // base material is needed only for a coup << 162 if(yes) f = dedx1*std::exp(fac*std::log(energy/energy1));
323 // initialised and for which tables will b << 163 else f = dedx1 + fac*(energy - energy1);
324 auto couple = theCoupleTable->GetMaterialC << 164 if(f > DBL_MIN) range += de/f;
325 auto pcuts = couple->GetProductionCuts(); << 165 }
326 auto mat = couple->GetMaterial(); << 166 // G4cout << "Range i= " <<i << " j= " << j << G4endl;
327 auto bmat = mat->GetBaseMaterial(); << 167 v->PutValue(j,range);
328 << 168 energy1 = energy2;
329 // base material exists - find it and chec << 169 dedx1 = dedx2;
330 if(nullptr != bmat) { <<
331 for(G4int j=0; j<(G4int)nCouples; ++j) { <<
332 if(j == i) { continue; } <<
333 auto bcouple = theCoupleTable->GetMaterialCu <<
334 <<
335 if(bcouple->GetMaterial() == bmat && <<
336 bcouple->GetProductionCuts() == pcuts) { <<
337 <<
338 // based couple exist in the same region <<
339 (*theDensityFactor)[i] = mat->GetDensity() <<
340 (*theDensityIdx)[i] = j; <<
341 (*theFlag)[i] = false; <<
342 <<
343 // ensure that there will no double initia <<
344 (*theDensityFactor)[j] = 1.0; <<
345 (*theDensityIdx)[j] = j; <<
346 (*theFlag)[j] = true; <<
347 break; <<
348 } <<
349 } 170 }
>> 171 G4PhysicsTableHelper::SetPhysicsVector(rangeTable, i, v);
350 } 172 }
351 } 173 }
352 } 174 }
353 175
354 //....oooOO0OOooo........oooOO0OOooo........oo 176 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
355 177
356 G4PhysicsTable* << 178 void G4LossTableBuilder::BuildInverseRangeTable(const G4PhysicsTable* rangeTable,
357 G4LossTableBuilder::BuildTableForModel(G4Physi << 179 G4PhysicsTable* invRangeTable,
358 G4VEmMo << 180 G4bool isIonisation)
359 const G << 181 // Build inverse range table from the energy loss table
360 G4doubl <<
361 G4bool <<
362 { 182 {
363 // check input << 183 size_t n_vectors = rangeTable->length();
364 G4PhysicsTable* table = G4PhysicsTableHelper << 184 if(!n_vectors) return;
365 if (nullptr == table) { return table; } << 185 G4bool b;
366 if (aTable != nullptr && aTable != table) { << 186
367 aTable->clearAndDestroy(); << 187 for (size_t i=0; i<n_vectors; i++) {
368 delete aTable; << 188
369 } << 189 if (invRangeTable->GetFlag(i) || !isIonisation) {
370 << 190 G4PhysicsVector* pv = (*rangeTable)[i];
371 InitialiseBaseMaterials(table); << 191 size_t nbins = pv->GetVectorLength();
372 G4int nbins = theParameters->NumberOfBinsPer << 192 G4double elow = pv->GetLowEdgeEnergy(0);
>> 193 G4double ehigh = pv->GetLowEdgeEnergy(nbins-1);
>> 194 G4double rlow = pv->GetValue(elow, b);
>> 195 G4double rhigh = pv->GetValue(ehigh, b);
>> 196
>> 197 rhigh *= std::exp(std::log(rhigh/rlow)/((G4double)(nbins-1)));
>> 198
>> 199 G4PhysicsLogVector* v = new G4PhysicsLogVector(rlow, rhigh, nbins);
>> 200
>> 201 v->PutValue(0,elow);
>> 202 G4double energy1 = elow;
>> 203 G4double range1 = rlow;
>> 204 G4double energy2 = elow;
>> 205 G4double range2 = rlow;
>> 206 size_t ilow = 0;
>> 207 size_t ihigh;
>> 208
>> 209 for (size_t j=1; j<nbins; j++) {
>> 210
>> 211 G4double range = v->GetLowEdgeEnergy(j);
>> 212
>> 213 for (ihigh=ilow+1; ihigh<nbins; ihigh++) {
>> 214 energy2 = pv->GetLowEdgeEnergy(ihigh);
>> 215 range2 = pv->GetValue(energy2, b);
>> 216 if(range2 >= range || ihigh == nbins-1) {
>> 217 ilow = ihigh - 1;
>> 218 energy1 = pv->GetLowEdgeEnergy(ilow);
>> 219 range1 = pv->GetValue(energy1, b);
>> 220 break;
>> 221 }
>> 222 }
373 223
374 // Access to materials << 224 G4double e = std::log(energy1) +
375 const G4ProductionCutsTable* theCoupleTable= << 225 std::log(energy2/energy1)*std::log(range/range1)/std::log(range2/range1);
376 G4ProductionCutsTable::GetProductionCu <<
377 std::size_t numOfCouples = theCoupleTable->G <<
378 /* <<
379 G4cout << " G4LossTableBuilder::BuildTable <<
380 << " isMaster=" << isInitializer << " model <<
381 << " " << part->GetParticleName() << G4end <<
382 */ <<
383 G4PhysicsLogVector* aVector = nullptr; <<
384 <<
385 for(G4int i=0; i<(G4int)numOfCouples; ++i) { <<
386 //G4cout << i << ". " << (*theFlag)[i] << <<
387 if (table->GetFlag(i)) { <<
388 <<
389 // create physics vector and fill it <<
390 auto couple = theCoupleTable->GetMateria <<
391 delete (*table)[i]; <<
392 <<
393 // if start from zero then change the sc <<
394 const G4Material* mat = couple->GetMater <<
395 <<
396 G4double tmin = std::max(emin, model->Mi <<
397 if(0.0 >= tmin) { tmin = CLHEP::eV; } <<
398 G4int n = nbins; <<
399 <<
400 if(tmin >= emax) { <<
401 aVector = nullptr; <<
402 } else { <<
403 n *= G4lrint(std::log10(emax/tmin)); <<
404 n = std::max(n, 3); <<
405 aVector = new G4PhysicsLogVector(tmin, <<
406 } <<
407 226
408 if(nullptr != aVector) { << 227 v->PutValue(j,std::exp(e));
409 //G4cout << part->GetParticleName() << <<
410 // << " emin= " << tmin << " emax=" << e <<
411 for(G4int j=0; j<=n; ++j) { <<
412 G4double e = aVector->Energy(j); <<
413 G4double y = model->Value(couple, part, e) <<
414 //G4cout << " " << j << ") E=" << e < <<
415 aVector->PutValue(j, y); <<
416 } <<
417 if(spline) { aVector->FillSecondDeriva <<
418 } 228 }
419 G4PhysicsTableHelper::SetPhysicsVector(t << 229 G4PhysicsTableHelper::SetPhysicsVector(invRangeTable, i, v);
420 } 230 }
421 } 231 }
422 /* <<
423 G4cout << "G4LossTableBuilder::BuildTableFor <<
424 << part->GetParticleName() << " and " <<
425 << " " << table << G4endl; <<
426 */ <<
427 //G4cout << *table << G4endl; <<
428 return table; <<
429 } 232 }
430 233
431 //....oooOO0OOooo........oooOO0OOooo........oo 234 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
>> 235
432 236