Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/processes/electromagnetic/xrays/src/G4XTRTransparentRegRadModel.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/xrays/src/G4XTRTransparentRegRadModel.cc (Version 11.3.0) and /processes/electromagnetic/xrays/src/G4XTRTransparentRegRadModel.cc (Version 9.4.p2)


  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 //
 26                                                    27 
 27 #include "G4XTRTransparentRegRadModel.hh"      <<  28 #include <complex>
 28                                                    29 
 29 #include "G4PhysicalConstants.hh"              <<  30 #include "G4XTRTransparentRegRadModel.hh"
                                                   >>  31 #include "Randomize.hh"
                                                   >>  32 #include "G4Integrator.hh"
                                                   >>  33 #include "G4Gamma.hh"
 30                                                    34 
 31 //////////////////////////////////////////////     35 ////////////////////////////////////////////////////////////////////////////
                                                   >>  36 //
 32 // Constructor, destructor                         37 // Constructor, destructor
 33 G4XTRTransparentRegRadModel::G4XTRTransparentR <<  38 
 34   G4LogicalVolume* anEnvelope, G4Material* foi <<  39 G4XTRTransparentRegRadModel::G4XTRTransparentRegRadModel(G4LogicalVolume *anEnvelope,
 35   G4double a, G4double b, G4int n, const G4Str <<  40            G4Material* foilMat,G4Material* gasMat, 
 36   : G4VXTRenergyLoss(anEnvelope, foilMat, gasM <<  41                                          G4double a, G4double b, G4int n,
                                                   >>  42                                          const G4String& processName) :
                                                   >>  43   G4VXTRenergyLoss(anEnvelope,foilMat,gasMat,a,b,n,processName)
 37 {                                                  44 {
 38   G4cout << "Regular transparent X-ray TR  rad <<  45   G4cout<<"Regular transparent X-ray TR  radiator EM process is called"<<G4endl;
 39          << G4endl;                            << 
 40                                                    46 
                                                   >>  47   // Build energy and angular integral spectra of X-ray TR photons from
                                                   >>  48   // a radiator
 41   fExitFlux   = true;                              49   fExitFlux   = true;
 42   fAlphaPlate = 10000;                             50   fAlphaPlate = 10000;
 43   fAlphaGas   = 1000;                              51   fAlphaGas   = 1000;
                                                   >>  52 
                                                   >>  53   //  BuildTable();
 44 }                                                  54 }
 45                                                    55 
 46 //////////////////////////////////////////////     56 ///////////////////////////////////////////////////////////////////////////
 47 G4XTRTransparentRegRadModel::~G4XTRTransparent << 
 48                                                    57 
 49 ////////////////////////////////////////////// <<  58 G4XTRTransparentRegRadModel::~G4XTRTransparentRegRadModel()
 50 void G4XTRTransparentRegRadModel::ProcessDescr << 
 51 {                                                  59 {
 52   out << "Process describing radiator of X-ray <<  60   ;
 53 }                                                  61 }
 54                                                    62 
 55 //////////////////////////////////////////////     63 ///////////////////////////////////////////////////////////////////////////
                                                   >>  64 //
                                                   >>  65 //
                                                   >>  66 
 56 G4double G4XTRTransparentRegRadModel::Spectral     67 G4double G4XTRTransparentRegRadModel::SpectralXTRdEdx(G4double energy)
 57 {                                                  68 {
 58   static constexpr G4double cofPHC = 4. * pi * <<  69   G4double result, sum = 0., tmp, cof1, cof2, cofMin, cofPHC,aMa, bMb, sigma;
 59   G4double result, sum = 0., tmp, cof1, cof2,  << 
 60   G4int k, kMax, kMin;                             70   G4int k, kMax, kMin;
 61                                                    71 
 62   aMa = GetPlateLinearPhotoAbs(energy);            72   aMa = GetPlateLinearPhotoAbs(energy);
 63   bMb = GetGasLinearPhotoAbs(energy);              73   bMb = GetGasLinearPhotoAbs(energy);
 64                                                    74 
 65   if(fCompton)                                     75   if(fCompton)
 66   {                                                76   {
 67     aMa += GetPlateCompton(energy);                77     aMa += GetPlateCompton(energy);
 68     bMb += GetGasCompton(energy);                  78     bMb += GetGasCompton(energy);
 69   }                                                79   }
 70   aMa *= fPlateThick;                              80   aMa *= fPlateThick;
 71   bMb *= fGasThick;                                81   bMb *= fGasThick;
 72                                                    82 
 73   sigma = aMa + bMb;                               83   sigma = aMa + bMb;
                                                   >>  84    
                                                   >>  85   cofPHC  = 4*pi*hbarc;
                                                   >>  86   tmp     = (fSigma1 - fSigma2)/cofPHC/energy;  
                                                   >>  87   cof1    = fPlateThick*tmp;
                                                   >>  88   cof2    = fGasThick*tmp;
 74                                                    89 
 75   tmp  = (fSigma1 - fSigma2) / cofPHC / energy <<  90   cofMin  =  energy*(fPlateThick + fGasThick)/fGamma/fGamma;
 76   cof1 = fPlateThick * tmp;                    <<  91   cofMin += (fPlateThick*fSigma1 + fGasThick*fSigma2)/energy;
 77   cof2 = fGasThick * tmp;                      << 
 78                                                << 
 79   cofMin = energy * (fPlateThick + fGasThick)  << 
 80   cofMin += (fPlateThick * fSigma1 + fGasThick << 
 81   cofMin /= cofPHC;                                92   cofMin /= cofPHC;
 82                                                    93 
                                                   >>  94   //  if (fGamma < 1200) kMin = G4int(cofMin);  // 1200 ?
                                                   >>  95   // else               kMin = 1;
                                                   >>  96 
                                                   >>  97 
 83   kMin = G4int(cofMin);                            98   kMin = G4int(cofMin);
 84   if(cofMin > kMin)                            <<  99   if (cofMin > kMin) kMin++;
 85     kMin++;                                    << 
 86                                                   100 
 87   kMax = kMin + 19;                            << 101   // tmp  = (fPlateThick + fGasThick)*energy*fMaxThetaTR;
                                                   >> 102   // tmp /= cofPHC;
                                                   >> 103   // kMax = G4int(tmp);
                                                   >> 104   // if(kMax < 0) kMax = 0;
                                                   >> 105   // kMax += kMin;
                                                   >> 106   
                                                   >> 107 
                                                   >> 108   kMax = kMin + 19; // 5; // 9; //   kMin + G4int(tmp);
                                                   >> 109 
                                                   >> 110   // tmp /= fGamma;
                                                   >> 111   // if( G4int(tmp) < kMin ) kMin = G4int(tmp);
                                                   >> 112   // G4cout<<"kMin = "<<kMin<<";    kMax = "<<kMax<<G4endl;
 88                                                   113 
 89   for(k = kMin; k <= kMax; k++)                << 114   for( k = kMin; k <= kMax; k++ )
 90   {                                               115   {
 91     tmp    = pi * fPlateThick * (k + cof2) / ( << 116     tmp    = pi*fPlateThick*(k + cof2)/(fPlateThick + fGasThick);
 92     result = (k - cof1) * (k - cof1) * (k + co << 117     result = (k - cof1)*(k - cof1)*(k + cof2)*(k + cof2);
 93                                                   118 
 94     if(k == kMin && kMin == G4int(cofMin))     << 119     if( k == kMin && kMin == G4int(cofMin) )
 95     {                                             120     {
 96       sum +=                                   << 121       sum   += 0.5*std::sin(tmp)*std::sin(tmp)*std::abs(k-cofMin)/result;
 97         0.5 * std::sin(tmp) * std::sin(tmp) *  << 
 98     }                                             122     }
 99     else                                          123     else
100     {                                             124     {
101       sum += std::sin(tmp) * std::sin(tmp) * s << 125       sum   += std::sin(tmp)*std::sin(tmp)*std::abs(k-cofMin)/result;
102     }                                             126     }
                                                   >> 127     //  G4cout<<"k = "<<k<<";    sum = "<<sum<<G4endl;    
103   }                                               128   }
104   result = 4. * (cof1 + cof2) * (cof1 + cof2)  << 129   result = 4.*( cof1 + cof2 )*( cof1 + cof2 )*sum/energy;
105   result *= (1. - std::exp(-fPlateNumber * sig << 130   result *= ( 1. - std::exp(-fPlateNumber*sigma) )/( 1. - std::exp(-sigma) );  
106   return result;                                  131   return result;
107 }                                                 132 }
108                                                   133 
                                                   >> 134 
109 //////////////////////////////////////////////    135 ///////////////////////////////////////////////////////////////////////////
                                                   >> 136 //
110 // Approximation for radiator interference fac    137 // Approximation for radiator interference factor for the case of
111 // fully Regular radiator. The plate and gas g << 138 // fully Regular radiator. The plate and gas gap thicknesses are fixed .
112 // The mean values of the plate and gas gap th << 139 // The mean values of the plate and gas gap thicknesses 
113 // are supposed to be about XTR formation zone << 140 // are supposed to be about XTR formation zones but much less than 
114 // mean absorption length of XTR photons in co << 141 // mean absorption length of XTR photons in coresponding material.
115 G4double G4XTRTransparentRegRadModel::GetStack << 142 
116                                                << 143 G4double 
117                                                << 144 G4XTRTransparentRegRadModel::GetStackFactor( G4double energy, 
                                                   >> 145                                          G4double gamma, G4double varAngle )
118 {                                                 146 {
119   G4double aZa   = fPlateThick / GetPlateForma << 147   /*
120   G4double bZb   = fGasThick / GetGasFormation << 148   G4double result, Za, Zb, Ma, Mb, sigma;
121   G4double aMa   = fPlateThick * GetPlateLinea << 149   
122   G4double bMb   = fGasThick * GetGasLinearPho << 150   Za = GetPlateFormationZone(energy,gamma,varAngle);
123   G4double sigma = aMa * fPlateThick + bMb * f << 151   Zb = GetGasFormationZone(energy,gamma,varAngle);
124   G4double Qa    = std::exp(-0.5 * aMa);       << 152   Ma = GetPlateLinearPhotoAbs(energy);
125   G4double Qb    = std::exp(-0.5 * bMb);       << 153   Mb = GetGasLinearPhotoAbs(energy);
126   G4double Q     = Qa * Qb;                    << 154   sigma = Ma*fPlateThick + Mb*fGasThick;
127                                                << 155 
128   G4complex Ha(Qa * std::cos(aZa), -Qa * std:: << 156   G4complex Ca(1.0+0.5*fPlateThick*Ma/fAlphaPlate,fPlateThick/Za/fAlphaPlate); 
129   G4complex Hb(Qb * std::cos(bZb), -Qb * std:: << 157   G4complex Cb(1.0+0.5*fGasThick*Mb/fAlphaGas,fGasThick/Zb/fAlphaGas); 
130   G4complex H  = Ha * Hb;                      << 158 
                                                   >> 159   G4complex Ha = std::pow(Ca,-fAlphaPlate);  
                                                   >> 160   G4complex Hb = std::pow(Cb,-fAlphaGas);
                                                   >> 161   G4complex H  = Ha*Hb;
                                                   >> 162   G4complex F1 =   (1.0 - Ha)*(1.0 - Hb )/(1.0 - H)
                                                   >> 163                  * G4double(fPlateNumber) ;
                                                   >> 164   G4complex F2 =   (1.0-Ha)*(1.0-Ha)*Hb/(1.0-H)/(1.0-H)
                                                   >> 165                  * (1.0 - std::exp(-0.5*fPlateNumber*sigma)) ;
                                                   >> 166   //    *(1.0 - std::pow(H,fPlateNumber)) ;
                                                   >> 167     G4complex R  = (F1 + F2)*OneInterfaceXTRdEdx(energy,gamma,varAngle);
                                                   >> 168   // G4complex R  = F2*OneInterfaceXTRdEdx(energy,gamma,varAngle);
                                                   >> 169   result       = 2.0*std::real(R);  
                                                   >> 170   return      result;
                                                   >> 171   */
                                                   >> 172    // numerically unstable result
                                                   >> 173 
                                                   >> 174   G4double result, Qa, Qb, Q, aZa, bZb, aMa, bMb, D, sigma; 
                                                   >> 175  
                                                   >> 176   aZa   = fPlateThick/GetPlateFormationZone(energy,gamma,varAngle);
                                                   >> 177   bZb   = fGasThick/GetGasFormationZone(energy,gamma,varAngle);
                                                   >> 178   aMa   = fPlateThick*GetPlateLinearPhotoAbs(energy);
                                                   >> 179   bMb   = fGasThick*GetGasLinearPhotoAbs(energy);
                                                   >> 180   sigma = aMa*fPlateThick + bMb*fGasThick;
                                                   >> 181   Qa    = std::exp(-0.5*aMa);
                                                   >> 182   Qb    = std::exp(-0.5*bMb);
                                                   >> 183   Q     = Qa*Qb;
                                                   >> 184 
                                                   >> 185   G4complex Ha( Qa*std::cos(aZa), -Qa*std::sin(aZa)   );  
                                                   >> 186   G4complex Hb( Qb*std::cos(bZb), -Qb*std::sin(bZb)    );
                                                   >> 187   G4complex H  = Ha*Hb;
131   G4complex Hs = conj(H);                         188   G4complex Hs = conj(H);
132   G4double D =                                 << 189   D            = 1.0 /( (1 - Q)*(1 - Q) + 
133     1.0 / ((1. - Q) * (1. - Q) +               << 190                   4*Q*std::sin(0.5*(aZa + bZb))*std::sin(0.5*(aZa + bZb)) );
134            4. * Q * std::sin(0.5 * (aZa + bZb) << 191   G4complex F1 = (1.0 - Ha)*(1.0 - Hb)*(1.0 - Hs)
135   G4complex F1 =                               << 192                  * G4double(fPlateNumber)*D;
136     (1.0 - Ha) * (1.0 - Hb) * (1.0 - Hs) * G4d << 193   G4complex F2 = (1.0 - Ha)*(1.0 - Ha)*Hb*(1.0 - Hs)*(1.0 - Hs)
137   G4complex F2 = (1.0 - Ha) * (1.0 - Ha) * Hb  << 194                    // * (1.0 - std::pow(H,fPlateNumber)) * D*D;
138                  (1.0 - std::exp(-0.5 * fPlate << 195                  * (1.0 - std::exp(-0.5*fPlateNumber*sigma)) * D*D;
139   G4complex R = (F1 + F2) * OneInterfaceXTRdEd << 196   G4complex R  = (F1 + F2)*OneInterfaceXTRdEdx(energy,gamma,varAngle);
140   return 2.0 * std::real(R);                   << 197   result       = 2.0*std::real(R); 
                                                   >> 198   return      result;
                                                   >> 199   
141 }                                                 200 }
                                                   >> 201 
                                                   >> 202 
                                                   >> 203 //
                                                   >> 204 //
                                                   >> 205 ////////////////////////////////////////////////////////////////////////////
                                                   >> 206 
                                                   >> 207 
                                                   >> 208 
                                                   >> 209 
                                                   >> 210 
                                                   >> 211 
                                                   >> 212 
                                                   >> 213 
142                                                   214