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