Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/processes/electromagnetic/utils/src/G4LossTableBuilder.cc

Version: [ ReleaseNotes ] [ 1.0 ] [ 1.1 ] [ 2.0 ] [ 3.0 ] [ 3.1 ] [ 3.2 ] [ 4.0 ] [ 4.0.p1 ] [ 4.0.p2 ] [ 4.1 ] [ 4.1.p1 ] [ 5.0 ] [ 5.0.p1 ] [ 5.1 ] [ 5.1.p1 ] [ 5.2 ] [ 5.2.p1 ] [ 5.2.p2 ] [ 6.0 ] [ 6.0.p1 ] [ 6.1 ] [ 6.2 ] [ 6.2.p1 ] [ 6.2.p2 ] [ 7.0 ] [ 7.0.p1 ] [ 7.1 ] [ 7.1.p1 ] [ 8.0 ] [ 8.0.p1 ] [ 8.1 ] [ 8.1.p1 ] [ 8.1.p2 ] [ 8.2 ] [ 8.2.p1 ] [ 8.3 ] [ 8.3.p1 ] [ 8.3.p2 ] [ 9.0 ] [ 9.0.p1 ] [ 9.0.p2 ] [ 9.1 ] [ 9.1.p1 ] [ 9.1.p2 ] [ 9.1.p3 ] [ 9.2 ] [ 9.2.p1 ] [ 9.2.p2 ] [ 9.2.p3 ] [ 9.2.p4 ] [ 9.3 ] [ 9.3.p1 ] [ 9.3.p2 ] [ 9.4 ] [ 9.4.p1 ] [ 9.4.p2 ] [ 9.4.p3 ] [ 9.4.p4 ] [ 9.5 ] [ 9.5.p1 ] [ 9.5.p2 ] [ 9.6 ] [ 9.6.p1 ] [ 9.6.p2 ] [ 9.6.p3 ] [ 9.6.p4 ] [ 10.0 ] [ 10.0.p1 ] [ 10.0.p2 ] [ 10.0.p3 ] [ 10.0.p4 ] [ 10.1 ] [ 10.1.p1 ] [ 10.1.p2 ] [ 10.1.p3 ] [ 10.2 ] [ 10.2.p1 ] [ 10.2.p2 ] [ 10.2.p3 ] [ 10.3 ] [ 10.3.p1 ] [ 10.3.p2 ] [ 10.3.p3 ] [ 10.4 ] [ 10.4.p1 ] [ 10.4.p2 ] [ 10.4.p3 ] [ 10.5 ] [ 10.5.p1 ] [ 10.6 ] [ 10.6.p1 ] [ 10.6.p2 ] [ 10.6.p3 ] [ 10.7 ] [ 10.7.p1 ] [ 10.7.p2 ] [ 10.7.p3 ] [ 10.7.p4 ] [ 11.0 ] [ 11.0.p1 ] [ 11.0.p2 ] [ 11.0.p3, ] [ 11.0.p4 ] [ 11.1 ] [ 11.1.1 ] [ 11.1.2 ] [ 11.1.3 ] [ 11.2 ] [ 11.2.1 ] [ 11.2.2 ] [ 11.3.0 ]

Diff markup

Differences between /processes/electromagnetic/utils/src/G4LossTableBuilder.cc (Version 11.3.0) and /processes/electromagnetic/utils/src/G4LossTableBuilder.cc (Version 9.2.p3)


  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.27 2008/07/22 15:55:15 vnivanch Exp $
                                                   >>  27 // GEANT4 tag $Name: geant4-09-02-patch-03 $
                                                   >>  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     48 // 12-02-07 Use G4LPhysicsFreeVector for the inverse range table (V.Ivanchenko)
 46 // 24-06-09 Removed hidden bin in G4PhysicsVec << 
 47 //                                                 49 //
 48 // Class Description:                              50 // Class Description:
 49 //                                                 51 //
 50 // -------------------------------------------     52 // -------------------------------------------------------------------
 51 //                                                 53 //
 52 //....oooOO0OOooo........oooOO0OOooo........oo     54 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
 53 //....oooOO0OOooo........oooOO0OOooo........oo     55 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
 54                                                    56 
 55 #include "G4LossTableBuilder.hh"                   57 #include "G4LossTableBuilder.hh"
 56 #include "G4SystemOfUnits.hh"                  << 
 57 #include "G4PhysicsTable.hh"                       58 #include "G4PhysicsTable.hh"
 58 #include "G4PhysicsLogVector.hh"                   59 #include "G4PhysicsLogVector.hh"
 59 #include "G4PhysicsTableHelper.hh"                 60 #include "G4PhysicsTableHelper.hh"
 60 #include "G4PhysicsFreeVector.hh"              <<  61 #include "G4LPhysicsFreeVector.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                                                    62 
 74 //....oooOO0OOooo........oooOO0OOooo........oo     63 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
 75                                                    64 
 76 G4LossTableBuilder::G4LossTableBuilder(G4bool  <<  65 G4LossTableBuilder::G4LossTableBuilder() 
 77   : isInitializer(master)                      << 
 78 {                                                  66 {
 79   theParameters = G4EmParameters::Instance();  <<  67   splineFlag = true;
 80   if (nullptr == theFlag) {                    << 
 81     theDensityFactor = new std::vector<G4doubl << 
 82     theDensityIdx = new std::vector<G4int>;    << 
 83     theFlag = new std::vector<G4bool>;         << 
 84   }                                            << 
 85 }                                                  68 }
 86                                                    69 
 87 //....oooOO0OOooo........oooOO0OOooo........oo     70 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
 88                                                    71 
 89 G4LossTableBuilder::~G4LossTableBuilder()          72 G4LossTableBuilder::~G4LossTableBuilder() 
 90 {                                              <<  73 {}
 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                                                    74 
115 //....oooOO0OOooo........oooOO0OOooo........oo     75 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
116                                                    76 
117 G4bool G4LossTableBuilder::GetFlag(std::size_t <<  77 void 
                                                   >>  78 G4LossTableBuilder::BuildDEDXTable(G4PhysicsTable* dedxTable,
                                                   >>  79            const std::vector<G4PhysicsTable*>& list)
118 {                                                  80 {
119   return (idx < theFlag->size()) ? (*theFlag)[ <<  81   size_t n_processes = list.size();
120 }                                              <<  82   if(1 >= n_processes) return;
121                                                    83 
122 //....oooOO0OOooo........oooOO0OOooo........oo <<  84   size_t n_vectors = (list[0])->length();
                                                   >>  85   if(0 >= n_vectors) return;
123                                                    86 
124 G4bool G4LossTableBuilder::GetBaseMaterialFlag <<  87   G4bool b;
125 {                                              << 
126   return baseMatFlag;                          << 
127 }                                              << 
128                                                    88 
129 //....oooOO0OOooo........oooOO0OOooo........oo <<  89   G4PhysicsVector* pv = (*(list[0]))[0];
                                                   >>  90   size_t nbins = pv->GetVectorLength();
                                                   >>  91   G4double elow = pv->GetLowEdgeEnergy(0);
                                                   >>  92   G4double ehigh = pv->GetLowEdgeEnergy(nbins);
                                                   >>  93   for (size_t i=0; i<n_vectors; i++) {
130                                                    94 
131 void                                           <<  95     pv = new G4PhysicsLogVector(elow, ehigh, nbins);
132 G4LossTableBuilder::BuildDEDXTable(G4PhysicsTa <<  96     pv->SetSpline(splineFlag);
133                                    const std:: <<  97     for (size_t j=0; j<nbins; j++) {
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;                         98       G4double dedx = 0.0;
152       for (std::size_t k=0; k<n_processes; ++k <<  99       G4double energy = pv->GetLowEdgeEnergy(j);
153   const G4PhysicsVector* pv1 = (*(list[k]))[i] << 100 
154   //if (0 == i) G4cout << "     " << k << ". p << 101       for (size_t k=0; k<n_processes; k++) {
155   dedx += (*pv1)[j];                           << 102   dedx += ((*(list[k]))[i])->GetValue(energy, b);
156       }                                           103       }
157       pv->PutValue(j, dedx);                      104       pv->PutValue(j, dedx);
                                                   >> 105       G4PhysicsTableHelper::SetPhysicsVector(dedxTable, i, pv);
158     }                                             106     }
159     if(splineFlag) { pv->FillSecondDerivatives << 
160     G4PhysicsTableHelper::SetPhysicsVector(ded << 
161   }                                               107   }
162   //G4cout << "### G4LossTableBuilder::BuildDE << 
163   //G4cout << *dedxTable << G4endl;            << 
164 }                                                 108 }
165                                                   109 
166 //....oooOO0OOooo........oooOO0OOooo........oo    110 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
167                                                   111 
168 void G4LossTableBuilder::BuildRangeTable(const    112 void G4LossTableBuilder::BuildRangeTable(const G4PhysicsTable* dedxTable,
169                                          G4Phy << 113            G4PhysicsTable* rangeTable,
                                                   >> 114            G4bool isIonisation)
170 // Build range table from the energy loss tabl    115 // Build range table from the energy loss table
171 {                                                 116 {
172   //G4cout << "### G4LossTableBuilder::BuildRa << 117   size_t n_vectors = dedxTable->length();
173   //G4cout << *const_cast<G4PhysicsTable*>(ded << 118   if(!n_vectors) return;
174   const std::size_t nCouples = dedxTable->size << 119 
175   if(0 >= nCouples) { return; }                << 120   G4bool b;
176                                                << 121   size_t n = 100;
177   const std::size_t n = 100;                   << 122   G4double del = 1.0/(G4double)n;
178   const G4double del = 1.0/(G4double)n;        << 123 
179                                                << 124   for (size_t i=0; i<n_vectors; i++) {
180   for (std::size_t i=0; i<nCouples; ++i) {     << 125 
181     auto pv = static_cast<G4PhysicsLogVector*> << 126     if (rangeTable->GetFlag(i) || !isIonisation) {
182     if((pv == nullptr) || (isBaseMatActive &&  << 127       G4PhysicsVector* pv = (*dedxTable)[i];
183     std::size_t npoints = pv->GetVectorLength( << 128       size_t nbins = pv->GetVectorLength();
184     std::size_t bin0    = 0;                   << 129       size_t bin0  = 0;
185     G4double elow  = pv->Energy(0);            << 130       G4double elow = pv->GetLowEdgeEnergy(0);
186     G4double ehigh = pv->Energy(npoints-1);    << 131       G4double ehigh = pv->GetLowEdgeEnergy(nbins);
187     G4double dedx1 = (*pv)[0];                 << 132       G4double dedx1  = pv->GetValue(elow, b);
188                                                << 133 
189     // protection for specific cases dedx=0    << 134       // protection for specific cases dedx=0
190     if(dedx1 == 0.0) {                         << 135       if(dedx1 == 0.0) {
191       for (std::size_t k=1; k<npoints; ++k) {  << 136         for (size_t k=1; k<nbins; k++) {
192   ++bin0;                                      << 137           bin0++;
193   elow  = pv->Energy(k);                       << 138           elow  = pv->GetLowEdgeEnergy(k);
194   dedx1 = (*pv)[k];                            << 139           dedx1 = pv->GetValue(elow, b);
195   if(dedx1 > 0.0) { break; }                   << 140           if(dedx1 > 0.0) break;
                                                   >> 141         }
                                                   >> 142         nbins -= bin0;
196       }                                           143       }
197       npoints -= bin0;                         << 
198     }                                          << 
199                                                   144 
200     // initialisation of a new vector          << 145       // initialisation of a new vector
201     if(npoints < 3) { npoints = 3; }           << 146       G4PhysicsLogVector* v = new G4PhysicsLogVector(elow, ehigh, nbins);
                                                   >> 147       v->SetSpline(splineFlag);
                                                   >> 148 
                                                   >> 149       // assumed dedx proportional to beta
                                                   >> 150       G4double range  = 2.*elow/dedx1;
                                                   >> 151       v->PutValue(0,range);
                                                   >> 152       G4double energy1 = elow;
                                                   >> 153 
                                                   >> 154       for (size_t j=1; j<nbins; j++) {
                                                   >> 155 
                                                   >> 156         G4double energy2 = pv->GetLowEdgeEnergy(j+bin0);
                                                   >> 157         G4double de      = (energy2 - energy1) * del;
                                                   >> 158         G4double energy  = energy2 + de*0.5;
                                                   >> 159         for (size_t k=0; k<n; k++) {
                                                   >> 160           energy -= de;
                                                   >> 161           dedx1 = pv->GetValue(energy, b);
                                                   >> 162           if(dedx1 > 0.0) range += de/dedx1;
                                                   >> 163   }
202                                                   164 
203     delete (*rangeTable)[i];                   << 165   //  G4cout << "Range i= " <<i << " j= " << j << G4endl;
204     G4PhysicsLogVector* v;                     << 166         v->PutValue(j,range);
205     if(0 == bin0) { v = new G4PhysicsLogVector << 167         energy1 = energy2;
206     else { v = new G4PhysicsLogVector(elow, eh << 
207                                                << 
208     // assumed dedx proportional to beta       << 
209     G4double energy1 = v->Energy(0);           << 
210     G4double range   = 2.*energy1/dedx1;       << 
211     /*                                         << 
212     G4cout << "New Range vector Npoints=" << v << 
213      << " coupleIdx=" << i << " spline=" << v- << 
214      << " Elow=" << v->GetMinEnergy() <<" Ehig << 
215      << " DEDX(Elow)=" << dedx1 << " R(Elow)=" << 
216     */                                         << 
217     v->PutValue(0,range);                      << 
218                                                << 
219     for (std::size_t j=1; j<npoints; ++j) {    << 
220                                                << 
221       G4double energy2 = v->Energy(j);         << 
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       }                                           168       }
231       range += sum;                            << 169       G4PhysicsTableHelper::SetPhysicsVector(rangeTable, i, v);
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     }                                             170     }
240     if(splineFlag) { v->FillSecondDerivatives( << 
241     G4PhysicsTableHelper::SetPhysicsVector(ran << 
242   }                                               171   }
243   //G4cout << "### Range table" << G4endl;     << 
244   //G4cout << *rangeTable << G4endl;           << 
245 }                                                 172 }
246                                                   173 
247 //....oooOO0OOooo........oooOO0OOooo........oo    174 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
248                                                   175 
249 void                                           << 176 void G4LossTableBuilder::BuildInverseRangeTable(const G4PhysicsTable* rangeTable,
250 G4LossTableBuilder::BuildInverseRangeTable(con << 177             G4PhysicsTable* invRangeTable,
251                                            G4P << 178             G4bool isIonisation)
252 // Build inverse range table from the energy l    179 // Build inverse range table from the energy loss table
253 {                                                 180 {
254   std::size_t nCouples = rangeTable->size();   << 181   size_t n_vectors = rangeTable->length();
255   if(0 >= nCouples) { return; }                << 182   if(!n_vectors) return;
256                                                << 183   G4bool b;
257   for (std::size_t i=0; i<nCouples; ++i) {     << 184 
258     G4PhysicsVector* pv = (*rangeTable)[i];    << 185   for (size_t i=0; i<n_vectors; i++) {
259     if((pv == nullptr) || (isBaseMatActive &&  << 186 
260     std::size_t npoints = pv->GetVectorLength( << 187     if (invRangeTable->GetFlag(i) || !isIonisation) {
                                                   >> 188       G4PhysicsVector* pv = (*rangeTable)[i];
                                                   >> 189       size_t nbins   = pv->GetVectorLength();
                                                   >> 190       G4double elow  = pv->GetLowEdgeEnergy(0);
                                                   >> 191       G4double ehigh = pv->GetLowEdgeEnergy(nbins-1);
                                                   >> 192       G4double rlow  = pv->GetValue(elow, b);
                                                   >> 193       G4double rhigh = pv->GetValue(ehigh, b);
261                                                   194       
262     delete (*invRangeTable)[i];                << 195       G4LPhysicsFreeVector* v = new G4LPhysicsFreeVector(nbins,rlow,rhigh);
263     auto v = new G4PhysicsFreeVector(npoints,  << 196       v->SetSpline(splineFlag);
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                                                   197 
273     G4PhysicsTableHelper::SetPhysicsVector(inv << 198       for (size_t j=0; j<nbins; j++) {
274   }                                            << 199   G4double e  = pv->GetLowEdgeEnergy(j);
275   //G4cout << "### Inverse range table" << G4e << 200   G4double r  = pv->GetValue(e, b);
276   //G4cout << *invRangeTable << G4endl;        << 201         v->PutValues(j,r,e);
277 }                                              << 
278                                                << 
279 //....oooOO0OOooo........oooOO0OOooo........oo << 
280                                                << 
281 void G4LossTableBuilder::InitialiseBaseMateria << 
282 {                                              << 
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       }                                           202       }
301     }                                          << 203       v->PutValues(nbins,rhigh+rlow,ehigh);
302   }                                            << 
303                                                << 
304   if(nFlags != nCouples) { isInitialized = fal << 
305   if(isInitialized) { return; }                << 
306                                                << 
307   // reserve memory                            << 
308   theFlag->resize(nCouples, true);             << 
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                                                   204 
320   // use base materials                        << 205       G4PhysicsTableHelper::SetPhysicsVector(invRangeTable, i, v);
321   for (G4int i=0; i<(G4int)nCouples; ++i) {    << 
322     // base material is needed only for a coup << 
323     // initialised and for which tables will b << 
324     auto couple = theCoupleTable->GetMaterialC << 
325     auto pcuts = couple->GetProductionCuts();  << 
326     auto mat  = couple->GetMaterial();         << 
327     auto bmat = mat->GetBaseMaterial();        << 
328                                                << 
329     // base material exists - find it and chec << 
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       }                                        << 
350     }                                             206     }
351   }                                               207   }
352 }                                                 208 }
353                                                   209 
354 //....oooOO0OOooo........oooOO0OOooo........oo    210 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
355                                                   211 
356 G4PhysicsTable*                                << 
357 G4LossTableBuilder::BuildTableForModel(G4Physi << 
358                                        G4VEmMo << 
359                                        const G << 
360                                        G4doubl << 
361                                        G4bool  << 
362 {                                              << 
363   // check input                               << 
364   G4PhysicsTable* table = G4PhysicsTableHelper << 
365   if (nullptr == table) { return table; }      << 
366   if (aTable != nullptr && aTable != table) {  << 
367     aTable->clearAndDestroy();                 << 
368     delete aTable;                             << 
369   }                                            << 
370                                                << 
371   InitialiseBaseMaterials(table);              << 
372   G4int nbins = theParameters->NumberOfBinsPer << 
373                                                << 
374   // Access to materials                       << 
375   const G4ProductionCutsTable* theCoupleTable= << 
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                                                << 
408       if(nullptr != aVector) {                 << 
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       }                                        << 
419       G4PhysicsTableHelper::SetPhysicsVector(t << 
420     }                                          << 
421   }                                            << 
422   /*                                           << 
423   G4cout << "G4LossTableBuilder::BuildTableFor << 
424          << part->GetParticleName() << " and " << 
425          << "  " << table << G4endl;           << 
426   */                                           << 
427   //G4cout << *table << G4endl;                << 
428   return table;                                << 
429 }                                              << 
430                                                << 
431 //....oooOO0OOooo........oooOO0OOooo........oo << 
432                                                   212