Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/processes/electromagnetic/adjoint/src/G4AdjointPhotoElectricModel.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/adjoint/src/G4AdjointPhotoElectricModel.cc (Version 11.3.0) and /processes/electromagnetic/adjoint/src/G4AdjointPhotoElectricModel.cc (Version 9.2.p1)


  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                                                << 
 27 #include "G4AdjointPhotoElectricModel.hh"          26 #include "G4AdjointPhotoElectricModel.hh"
 28                                                << 
 29 #include "G4AdjointCSManager.hh"                   27 #include "G4AdjointCSManager.hh"
                                                   >>  28 
                                                   >>  29 
                                                   >>  30 #include "G4Integrator.hh"
                                                   >>  31 #include "G4TrackStatus.hh"
                                                   >>  32 #include "G4ParticleChange.hh"
 30 #include "G4AdjointElectron.hh"                    33 #include "G4AdjointElectron.hh"
                                                   >>  34 #include  "G4Gamma.hh"
 31 #include "G4AdjointGamma.hh"                       35 #include "G4AdjointGamma.hh"
 32 #include "G4Gamma.hh"                          << 
 33 #include "G4ParticleChange.hh"                 << 
 34 #include "G4PEEffectFluoModel.hh"              << 
 35 #include "G4PhysicalConstants.hh"              << 
 36 #include "G4TrackStatus.hh"                    << 
 37                                                    36 
 38 //////////////////////////////////////////////     37 ////////////////////////////////////////////////////////////////////////////////
 39 G4AdjointPhotoElectricModel::G4AdjointPhotoEle <<  38 //
 40   : G4VEmAdjointModel("AdjointPEEffect")       <<  39 G4AdjointPhotoElectricModel::G4AdjointPhotoElectricModel():
 41                                                <<  40  G4VEmAdjointModel("AdjointPEEffect")
 42 {                                              <<  41 
 43   SetUseMatrix(false);                         <<  42 { SetUseMatrix(false);
 44   SetApplyCutInRange(false);                   <<  43   current_eEnergy =0.;
 45                                                <<  44   totAdjointCS=0.;
 46   fAdjEquivDirectPrimPart   = G4AdjointGamma:: << 
 47   fAdjEquivDirectSecondPart = G4AdjointElectro << 
 48   fDirectPrimaryPart        = G4Gamma::Gamma() << 
 49   fSecondPartSameType       = false;           << 
 50   fDirectModel              = new G4PEEffectFl << 
 51 }                                                  45 }
 52                                                << 
 53 ////////////////////////////////////////////// << 
 54 G4AdjointPhotoElectricModel::~G4AdjointPhotoEl << 
 55                                                << 
 56 //////////////////////////////////////////////     46 ////////////////////////////////////////////////////////////////////////////////
 57 void G4AdjointPhotoElectricModel::SampleSecond <<  47 //
 58   const G4Track& aTrack, G4bool isScatProjToPr <<  48 G4AdjointPhotoElectricModel::~G4AdjointPhotoElectricModel()
 59   G4ParticleChange* fParticleChange)           <<  49 {;}
 60 {                                              <<  50 
 61   if(isScatProjToProj)                         <<  51 ////////////////////////////////////////////////////////////////////////////////
 62     return;                                    <<  52 //
 63                                                <<  53 void G4AdjointPhotoElectricModel::SampleSecondaries(const G4Track& aTrack,
 64   // Compute the fTotAdjointCS vectors if not  <<  54                                 G4bool IsScatProjToProjCase,
 65   // couple and electron energy                <<  55         G4ParticleChange* fParticleChange)
 66   const G4DynamicParticle* aDynPart = aTrack.G <<  56 { if (IsScatProjToProjCase) return ;
 67   G4double electronEnergy           = aDynPart <<  57 
 68   G4ThreeVector electronDirection   = aDynPart <<  58   //Compute the totAdjointCS vectors if not already done for the current couple and electron energy
 69   fPreStepAdjointCS =                          <<  59   const G4MaterialCutsCouple* aCouple = aTrack.GetMaterialCutsCouple();
 70     fTotAdjointCS;  // The last computed CS wa <<  60   const G4DynamicParticle* aDynPart =  aTrack.GetDynamicParticle() ;
 71   AdjointCrossSection(aTrack.GetMaterialCutsCo <<  61   G4double electronEnergy = aDynPart->GetKineticEnergy();
 72                       isScatProjToProj);       <<  62   G4ThreeVector electronDirection= aDynPart->GetMomentumDirection() ;
 73   fPostStepAdjointCS = fTotAdjointCS;          <<  63   totAdjointCS = AdjointCrossSection(aCouple, electronEnergy,IsScatProjToProjCase);
 74                                                <<  64         
 75   // Sample element                            <<  65  
 76   const G4ElementVector* theElementVector =    <<  66   //Sample gamma energy
 77     fCurrentMaterial->GetElementVector();      <<  67   //-------------
 78   size_t nelm      = fCurrentMaterial->GetNumb <<  68     /////////////////////////////////////////////////////////////////////////////////
 79   G4double rand_CS = G4UniformRand() * fXsec[n <<  69 //      Module:   G4ContinuousGainOfEnergy.hh
 80   for(fIndexElement = 0; fIndexElement < nelm  <<  70 //  Author:         L. Desorgher
 81   {                                            <<  71 //  Date:   1 September 2007
 82     if(rand_CS < fXsec[fIndexElement])         <<  72 //  Organisation:   SpaceIT GmbH
 83       break;                                   <<  73 //  Customer:       ESA/ESTEC
 84   }                                            <<  74 /////////////////////////////////////////////////////////////////////////////////
 85                                                <<  75 //
 86   // Sample shell and binding energy           <<  76 // CHANGE HISTORY
 87   G4int nShells = (*theElementVector)[fIndexEl <<  77 // --------------
 88   rand_CS       = fShellProb[fIndexElement][nS <<  78 //      ChangeHistory: 
 89   G4int i;                                     <<  79 //    1 September 2007 creation by L. Desorgher     
 90   for(i = 0; i < nShells - 1; ++i)             <<  80 //
 91   {                                            <<  81 //-------------------------------------------------------------
 92     if(rand_CS < fShellProb[fIndexElement][i]) <<  82 //  Documentation:
 93       break;                                   <<  83 //    Modell for the adjoint compton scattering
 94   }                                            <<  84 //
 95   G4double gammaEnergy =                       <<  85 
 96     electronEnergy + (*theElementVector)[fInde <<  86   //Sample element
 97                                                <<  87   //-------------
 98   // Sample cos theta                          <<  88    const G4ElementVector* theElementVector = currentMaterial->GetElementVector();
 99   // Copy of the G4PEEfectFluoModel cos theta  <<  89    const G4double* theAtomNumDensityVector = currentMaterial->GetVecNbOfAtomsPerVolume();
100   // ElecCosThetaDistribution. This method can <<  90    size_t nelm =  currentMaterial->GetNumberOfElements();
101   // G4PEEffectFluoModel because it is a frien <<  91    G4double rand_CS= totAdjointCS*G4UniformRand();
102   G4double cos_theta = 1.;                     <<  92    for (index_element=0; index_element<nelm-1; index_element++){
103   G4double gamma     = 1. + electronEnergy / e <<  93   if (rand_CS<xsec[index_element]) break;
104   if(gamma <= 5.)                              <<  94    }
105   {                                            <<  95   
106     G4double beta = std::sqrt(gamma * gamma -  <<  96    //Sample shell and binding energy
107     G4double b    = 0.5 * gamma * (gamma - 1.) <<  97    //-------------
108                                                <<  98    rand_CS= totAdjointCS*G4UniformRand()/theAtomNumDensityVector[index_element];
109     G4double rndm, term, greject, grejsup;     <<  99    G4int nShells = (*theElementVector)[index_element]->GetNbOfAtomicShells();
110     if(gamma < 2.)                             << 100    G4int i  = 0;  
111       grejsup = gamma * gamma * (1. + b - beta << 101    for (i=0; i<nShells-1; i++){
112     else                                       << 102   if (rand_CS<shell_prob[index_element][i]) break;
113       grejsup = gamma * gamma * (1. + b + beta << 103    }
114                                                << 104    G4double gammaEnergy= electronEnergy+(*theElementVector)[index_element]->GetAtomicShell(i);
115     do                                         << 105   
116     {                                          << 106   //Sample cos theta
117       rndm      = 1. - 2. * G4UniformRand();   << 107   //Copy of the G4PEEffectModel cos theta sampling method ElecCosThetaDistribution.   
118       cos_theta = (rndm + beta) / (rndm * beta << 108   //This method cannot be used directly from G4PEEffectModel because it is a friend method. I should ask Vladimir to change that  
119       term      = 1. - beta * cos_theta;       << 109   //------------------------------------------------------------------------------------------------  
120       greject = (1. - cos_theta * cos_theta) * << 110   //G4double cos_theta = theDirectPEEffectModel->ElecCosThetaDistribution(electronEnergy);
121       // Loop checking, 07-Aug-2015, Vladimir  << 111   
122     } while(greject < G4UniformRand() * grejsu << 112    G4double  cos_theta = 1.;
                                                   >> 113    G4double gamma   = 1. + electronEnergy/electron_mass_c2;
                                                   >> 114    if (gamma <= 5.) {
                                                   >> 115     G4double beta  = std::sqrt(gamma*gamma-1.)/gamma;
                                                   >> 116   G4double b     = 0.5*gamma*(gamma-1.)*(gamma-2);
                                                   >> 117 
                                                   >> 118     G4double rndm,term,greject,grejsup;
                                                   >> 119     if (gamma < 2.) grejsup = gamma*gamma*(1.+b-beta*b);
                                                   >> 120     else            grejsup = gamma*gamma*(1.+b+beta*b);
                                                   >> 121 
                                                   >> 122     do {  rndm = 1.-2*G4UniformRand();
                                                   >> 123           cos_theta = (rndm+beta)/(rndm*beta+1.);
                                                   >> 124           term = 1.-beta*cos_theta;
                                                   >> 125           greject = (1.-cos_theta*cos_theta)*(1.+b*term)/(term*term);
                                                   >> 126     } while(greject < G4UniformRand()*grejsup);
123   }                                               127   }
124                                                << 128   
125   // direction of the adjoint gamma electron      129   // direction of the adjoint gamma electron
126   G4double sin_theta = std::sqrt(1. - cos_thet << 130   //---------------------------------------
127   G4double phi       = twopi * G4UniformRand() << 131   
128   G4double dirx      = sin_theta * std::cos(ph << 132      
129   G4double diry      = sin_theta * std::sin(ph << 133   G4double sin_theta = std::sqrt(1.-cos_theta*cos_theta);
130   G4double dirz      = cos_theta;              << 134   G4double Phi     = twopi * G4UniformRand();
131   G4ThreeVector adjoint_gammaDirection(dirx, d << 135   G4double dirx = sin_theta*std::cos(Phi),diry = sin_theta*std::sin(Phi),dirz = cos_theta;
                                                   >> 136   G4ThreeVector adjoint_gammaDirection(dirx,diry,dirz);
132   adjoint_gammaDirection.rotateUz(electronDire    137   adjoint_gammaDirection.rotateUz(electronDirection);
133                                                << 138   
134   // Weight correction                         << 139   
135   CorrectPostStepWeight(fParticleChange, aTrac << 140   
136                         gammaEnergy, isScatPro << 141   //Weight correction
137                                                << 142  //-----------------------             
138   // Create secondary and modify fParticleChan << 143   CorrectPostStepWeight(fParticleChange, aTrack.GetWeight(), electronEnergy,gammaEnergy); 
139   G4DynamicParticle* anAdjointGamma = new G4Dy << 144  
140     G4AdjointGamma::AdjointGamma(), adjoint_ga << 145   
141                                                << 146   
                                                   >> 147   //Create secondary and modify fParticleChange 
                                                   >> 148   //--------------------------------------------
                                                   >> 149   G4DynamicParticle* anAdjointGamma = new G4DynamicParticle (
                                                   >> 150                        G4AdjointGamma::AdjointGamma(),adjoint_gammaDirection, gammaEnergy);
142   fParticleChange->ProposeTrackStatus(fStopAnd    151   fParticleChange->ProposeTrackStatus(fStopAndKill);
143   fParticleChange->AddSecondary(anAdjointGamma    152   fParticleChange->AddSecondary(anAdjointGamma);
144 }                                              << 153      
                                                   >> 154       
145                                                   155 
146 ////////////////////////////////////////////// << 156   
147 void G4AdjointPhotoElectricModel::CorrectPostS << 
148   G4ParticleChange* fParticleChange, G4double  << 
149   G4double adjointPrimKinEnergy, G4double proj << 
150 {                                              << 
151   G4double new_weight = old_weight;            << 
152                                                << 
153   G4double w_corr =                            << 
154     G4AdjointCSManager::GetAdjointCSManager()- << 
155     fFactorCSBiasing;                          << 
156   w_corr *= fPostStepAdjointCS / fPreStepAdjoi << 
157                                                << 
158   new_weight *= w_corr * projectileKinEnergy / << 
159   fParticleChange->SetParentWeightByProcess(fa << 
160   fParticleChange->SetSecondaryWeightByProcess << 
161   fParticleChange->ProposeParentWeight(new_wei << 
162 }                                              << 
163                                                   157 
164 ////////////////////////////////////////////// << 158 } 
165 G4double G4AdjointPhotoElectricModel::AdjointC << 
166   const G4MaterialCutsCouple* aCouple, G4doubl << 
167   G4bool isScatProjToProj)                     << 
168 {                                              << 
169   if(isScatProjToProj)                         << 
170     return 0.;                                 << 
171                                                << 
172   G4double totBiasedAdjointCS = 0.;            << 
173   if(aCouple != fCurrentCouple || fCurrenteEne << 
174   {                                            << 
175     fTotAdjointCS = 0.;                        << 
176     DefineCurrentMaterialAndElectronEnergy(aCo << 
177     const G4ElementVector* theElementVector =  << 
178       fCurrentMaterial->GetElementVector();    << 
179     const G4double* theAtomNumDensityVector =  << 
180       fCurrentMaterial->GetVecNbOfAtomsPerVolu << 
181     size_t nelm = fCurrentMaterial->GetNumberO << 
182     for(fIndexElement = 0; fIndexElement < nel << 
183     {                                          << 
184       fTotAdjointCS += AdjointCrossSectionPerA << 
185                          (*theElementVector)[f << 
186                        theAtomNumDensityVector << 
187       fXsec[fIndexElement] = fTotAdjointCS;    << 
188     }                                          << 
189                                                   159 
190     totBiasedAdjointCS = std::min(fTotAdjointC << 160 ////////////////////////////////////////////////////////////////////////////////
191     fFactorCSBiasing   = totBiasedAdjointCS /  << 161 //      
192   }                                            << 162 
193   return totBiasedAdjointCS;                   << 163 G4double G4AdjointPhotoElectricModel::AdjointCrossSection(const G4MaterialCutsCouple* aCouple,
                                                   >> 164         G4double electronEnergy,
                                                   >> 165         G4bool IsScatProjToProjCase)
                                                   >> 166 { if (IsScatProjToProjCase) return 0.;
                                                   >> 167   if (aCouple !=currentCouple || current_eEnergy !=electronEnergy) {
                                                   >> 168     totAdjointCS = 0.;
                                                   >> 169   DefineCurrentMaterialAndElectronEnergy(aCouple, electronEnergy);
                                                   >> 170     const G4ElementVector* theElementVector = currentMaterial->GetElementVector();
                                                   >> 171     const G4double* theAtomNumDensityVector = currentMaterial->GetVecNbOfAtomsPerVolume();
                                                   >> 172     size_t nelm =  currentMaterial->GetNumberOfElements();
                                                   >> 173     for (index_element=0;index_element<nelm;index_element++){
                                                   >> 174     
                                                   >> 175     totAdjointCS +=AdjointCrossSectionPerAtom((*theElementVector)[index_element],electronEnergy)*theAtomNumDensityVector[index_element];
                                                   >> 176     xsec[index_element] = totAdjointCS;
                                                   >> 177   } 
                                                   >> 178   } 
                                                   >> 179   return totAdjointCS;
                                                   >> 180   
                                                   >> 181   
194 }                                                 182 }
195                                                << 
196 //////////////////////////////////////////////    183 ////////////////////////////////////////////////////////////////////////////////
197 G4double G4AdjointPhotoElectricModel::AdjointC << 184 //      
198   const G4Element* anElement, G4double electro << 185 
199 {                                              << 186 G4double G4AdjointPhotoElectricModel::AdjointCrossSectionPerAtom(const G4Element* anElement,G4double electronEnergy)
200   G4int nShells        = anElement->GetNbOfAto << 187 { 
201   G4double Z           = anElement->GetZ();    << 188   G4int nShells = anElement->GetNbOfAtomicShells();
202   G4double gammaEnergy = electronEnergy + anEl << 189   G4double Z= anElement->GetZ();
203   G4double CS          = fDirectModel->Compute << 190   G4double N= anElement->GetN();
204     G4Gamma::Gamma(), gammaEnergy, Z, 0., 0.,  << 191   G4int i  = 0;  
205   G4double adjointCS = 0.;                     << 192   G4double B0=anElement->GetAtomicShell(0);
206   if(CS > 0.)                                  << 193   G4double gammaEnergy = electronEnergy+B0;
207     adjointCS += CS / gammaEnergy;             << 194   G4double adjointCS = theDirectPEEffectModel->ComputeCrossSectionPerAtom(G4Gamma::Gamma(),gammaEnergy,Z,N,0.,0.)*electronEnergy/gammaEnergy; 
208   fShellProb[fIndexElement][0] = adjointCS;    << 195   shell_prob[index_element][0] = adjointCS;                                          
209   for(G4int i = 1; i < nShells; ++i)           << 196   for (i=1;i<nShells;i++){
210   {                                            << 197     //G4cout<<i<<std::endl;
211     G4double Bi1 = anElement->GetAtomicShell(i << 198     G4double Bi_= anElement->GetAtomicShell(i-1);
212     G4double Bi  = anElement->GetAtomicShell(i << 199   G4double Bi = anElement->GetAtomicShell(i);
213     if(electronEnergy < Bi1 - Bi)              << 200   //G4cout<<Bi_<<'\t'<<Bi<<std::endl;
214     {                                          << 201   if (electronEnergy <Bi_-Bi) {
215       gammaEnergy = electronEnergy + Bi;       << 202     gammaEnergy = electronEnergy+Bi;
216       CS          = fDirectModel->ComputeCross << 203     adjointCS +=theDirectPEEffectModel->ComputeCrossSectionPerAtom(G4Gamma::Gamma(),gammaEnergy,anElement->GetZ(),N,0.,0.)*electronEnergy/gammaEnergy;
217                                                << 204   }
218       if(CS > 0.)                              << 205   shell_prob[index_element][i] = adjointCS; 
219         adjointCS += CS / gammaEnergy;         << 206   
220     }                                          << 
221     fShellProb[fIndexElement][i] = adjointCS;  << 
222   }                                               207   }
223   adjointCS *= electronEnergy;                 << 208   
224   return adjointCS;                               209   return adjointCS;
225 }                                              << 210   
226                                                << 211 }       
227 //////////////////////////////////////////////    212 ////////////////////////////////////////////////////////////////////////////////
228 void G4AdjointPhotoElectricModel::DefineCurren << 213 //      
229   const G4MaterialCutsCouple* couple, G4double << 214 
230 {                                              << 215 void G4AdjointPhotoElectricModel::DefineCurrentMaterialAndElectronEnergy(const G4MaterialCutsCouple* couple, G4double anEnergy)
231   fCurrentCouple   = const_cast<G4MaterialCuts << 216 { currentCouple   = const_cast<G4MaterialCutsCouple*> (couple);
232   fCurrentMaterial = const_cast<G4Material*>(c << 217   currentMaterial = const_cast<G4Material*> (couple->GetMaterial());
233   fCurrenteEnergy = anEnergy;                  << 218   currentCoupleIndex = couple->GetIndex();
234   fDirectModel->SetCurrentCouple(couple);      << 219   currentMaterialIndex = currentMaterial->GetIndex();
                                                   >> 220   current_eEnergy = anEnergy; 
235 }                                                 221 }
236                                                   222