G4PolarizedBremsstrahlungCrossSection.cc

Go to the documentation of this file.
00001 //
00002 // ********************************************************************
00003 // * License and Disclaimer                                           *
00004 // *                                                                  *
00005 // * The  Geant4 software  is  copyright of the Copyright Holders  of *
00006 // * the Geant4 Collaboration.  It is provided  under  the terms  and *
00007 // * conditions of the Geant4 Software License,  included in the file *
00008 // * LICENSE and available at  http://cern.ch/geant4/license .  These *
00009 // * include a list of copyright holders.                             *
00010 // *                                                                  *
00011 // * Neither the authors of this software system, nor their employing *
00012 // * institutes,nor the agencies providing financial support for this *
00013 // * work  make  any representation or  warranty, express or implied, *
00014 // * regarding  this  software system or assume any liability for its *
00015 // * use.  Please see the license in the file  LICENSE  and URL above *
00016 // * for the full disclaimer and the limitation of liability.         *
00017 // *                                                                  *
00018 // * This  code  implementation is the result of  the  scientific and *
00019 // * technical work of the GEANT4 collaboration.                      *
00020 // * By using,  copying,  modifying or  distributing the software (or *
00021 // * any work based  on the software)  you  agree  to acknowledge its *
00022 // * use  in  resulting  scientific  publications,  and indicate your *
00023 // * acceptance of all terms of the Geant4 Software license.          *
00024 // ********************************************************************
00025 //
00026 // $Id: G4PolarizedBremsstrahlungCrossSection.cc 69847 2013-05-16 09:36:18Z gcosmo $
00027 //
00028 // -------------------------------------------------------------------
00029 //
00030 // GEANT4 Class file
00031 //
00032 //
00033 // File name:     G4PolarizedBremsstrahlungCrossSection
00034 //
00035 // Author:        Andreas Schaelicke on the base of Karim Laihems code
00036 //
00037 // Creation date: 16.08.2006
00038 //
00039 
00040 #include "G4PolarizedBremsstrahlungCrossSection.hh"
00041 #include "G4PhysicalConstants.hh"
00042 
00043 G4bool G4PolarizedBremsstrahlungCrossSection::scrnInitialized=false;
00044 G4double G4PolarizedBremsstrahlungCrossSection::SCRN [3][20];  
00045 // screening function lookup table;
00046 
00047 
00048 void G4PolarizedBremsstrahlungCrossSection::InitializeMe()
00049 {
00050   if (!scrnInitialized) {
00051     SCRN [1][1]=  0.5   ; SCRN [2][1] = 0.0145;
00052     SCRN [1][2]=  1.0   ; SCRN [2][2] = 0.0490;
00053     SCRN [1][3]=  2.0   ; SCRN [2][3] = 0.1400;
00054     SCRN [1][4]=  4.0   ; SCRN [2][4] = 0.3312;
00055     SCRN [1][5]=  8.0   ; SCRN [2][5] = 0.6758;
00056     SCRN [1][6]=  15.0  ; SCRN [2][6] = 1.126;
00057     SCRN [1][7]=  20.0  ; SCRN [2][7] = 1.367;
00058     SCRN [1][8]=  25.0  ; SCRN [2][8] = 1.564;
00059     SCRN [1][9]=  30.0  ; SCRN [2][9] = 1.731;
00060     SCRN [1][10]= 35.0  ; SCRN [2][10]= 1.875;
00061     SCRN [1][11]= 40.0  ; SCRN [2][11]= 2.001;
00062     SCRN [1][12]= 45.0  ; SCRN [2][12]= 2.114;
00063     SCRN [1][13]= 50.0  ; SCRN [2][13]= 2.216;
00064     SCRN [1][14]= 60.0  ; SCRN [2][14]= 2.393;
00065     SCRN [1][15]= 70.0  ; SCRN [2][15]= 2.545;
00066     SCRN [1][16]= 80.0  ; SCRN [2][16]= 2.676;
00067     SCRN [1][17]= 90.0  ; SCRN [2][17]= 2.793;
00068     SCRN [1][18]= 100.0 ; SCRN [2][18]= 2.897;
00069     SCRN [1][19]= 120.0 ; SCRN [2][19]= 3.078;   
00070 
00071     scrnInitialized=true;
00072   }
00073 }
00074 
00075 G4PolarizedBremsstrahlungCrossSection::G4PolarizedBremsstrahlungCrossSection()
00076 {
00077   InitializeMe();
00078 }
00079 
00080 
00081 void G4PolarizedBremsstrahlungCrossSection::Initialize(
00082      G4double aLept0E, G4double aGammaE, G4double sintheta,
00083      const G4StokesVector & beamPol,
00084      const G4StokesVector & /*p1*/,
00085      G4int /*flag*/)
00086 {
00087 //   G4cout<<"G4PolarizedBremsstrahlungCrossSection::Initialize \n"
00088 //      <<"lepE = "<<aLept0E
00089 //      <<"gamE = "<<aGammaE
00090 //      <<"sint = "<<sintheta<<"\n"
00091 //      <<"beamPol="<<beamPol<<"\n";
00092 
00093   G4double aLept1E = aLept0E - aGammaE;
00094 
00095   G4double Stokes_S1  = beamPol.x()   ; 
00096   G4double Stokes_S2  = beamPol.y()   ; 
00097   G4double Stokes_S3  = beamPol.z()   ; 
00098   // **************************************************************************  
00099     
00100   G4double m0_c2  = electron_mass_c2; 
00101   G4double Lept0E = aLept0E/m0_c2+1.,   Lept0E2 = Lept0E * Lept0E ;
00102   G4double GammaE = aGammaE/m0_c2,      GammaE2 = GammaE * GammaE ;
00103   G4double Lept1E = aLept1E/m0_c2+1.,   Lept1E2 = Lept1E * Lept1E ;
00104     
00105 
00106   //  const G4Element* theSelectedElement = theModel->SelectedAtom();
00107 
00108   // *******  Gamma Transvers Momentum
00109     
00110   G4double TMom = std::sqrt(Lept0E2 -1.)* sintheta;
00111   G4double u    = TMom       , u2 =u * u ;
00112   G4double Xsi  = 1./(1.+u2)                      , Xsi2 = Xsi * Xsi  ; 
00113 
00114   //  G4double theZ  = theSelectedElement->GetZ();
00115     
00116   //  G4double fCoul = theSelectedElement->GetfCoulomb();
00117   G4double delta = 12. * std::pow(theZ, 1./3.) * 
00118     Lept0E * Lept1E * Xsi / (121. * GammaE); 
00119   G4double GG=0.;
00120 
00121   if(delta < 0.5) {
00122     GG = std::log(2.* Lept0E * Lept1E / GammaE) - 2. - fCoul; 
00123   }
00124   else if ( delta < 120) {
00125     for (G4int j=2; j<=19; j++)  {
00126       if(SCRN[1][j] >= delta)    {
00127         GG =std::log(2 * Lept0E * Lept1E / GammaE) - 2 - fCoul
00128           -(SCRN[2][j-1]+(delta-SCRN[1][j-1])*(SCRN[2][j]-SCRN[2][j-1])
00129             /(SCRN[1][j]-SCRN[1][j-1]));
00130         break;
00131       }
00132     }
00133   }
00134   else  {
00135     G4double alpha_sc  = (111 * std::pow(theZ, -1./3.)) / Xsi;
00136     GG = std::log(alpha_sc)- 2 - fCoul;
00137   }
00138 
00139   if(GG<-1) GG=-1;     // *KL* do we need this ?!
00140 
00141   G4double I_Lept   = (Lept0E2 + Lept1E2) * (3.+2.*GG) - 2 * Lept0E * Lept1E * (1. + 4. * u2 * Xsi2 * GG);
00142   G4double F_Lept   = Lept1E * 4. * GammaE *  u * Xsi * (1. - 2 * Xsi) * GG / I_Lept;
00143   G4double E_Lept   = Lept0E * 4. * GammaE *  u * Xsi * (2. * Xsi - 1.) * GG / I_Lept; 
00144   G4double M_Lept   = 4. * Lept0E * Lept1E * (1. + GG - 2. * Xsi2 * u2 * GG) / I_Lept ;
00145   G4double P_Lept   = GammaE2 * (1. + 8. * GG * (Xsi - 0.5)*(Xsi - 0.5)) / I_Lept ;
00146     
00147   G4double Stokes_SS1 = M_Lept * Stokes_S1 + E_Lept * Stokes_S3;
00148   G4double Stokes_SS2 = M_Lept * Stokes_S2 ;
00149   G4double Stokes_SS3 = (M_Lept + P_Lept) * Stokes_S3 + F_Lept * Stokes_S1; 
00150 
00151   theFinalLeptonPolarization.setX(Stokes_SS1);
00152   theFinalLeptonPolarization.setY(Stokes_SS2);
00153   theFinalLeptonPolarization.setZ(Stokes_SS3);
00154 
00155   if(theFinalLeptonPolarization.mag2()>1) { 
00156     G4cout<<" WARNING in pol-brem theFinalLeptonPolarization \n";
00157     G4cout
00158       <<"\t"<<theFinalLeptonPolarization
00159       <<"\t GG\t"<<GG
00160       <<"\t delta\t"<<delta
00161       <<G4endl;
00162     theFinalLeptonPolarization.setX(0);
00163     theFinalLeptonPolarization.setY(0);
00164     theFinalLeptonPolarization.setZ(Stokes_SS3);
00165     if(Stokes_SS3>1) theFinalLeptonPolarization.setZ(1);
00166   }
00167 
00168 
00169   G4double I_Gamma   = (Lept0E2 + Lept1E2)*(3+2*GG) - 2 * Lept0E * Lept1E * (1 + 4 * u2 * Xsi2 * GG);
00170   G4double D_Gamma   = 8 * Lept0E * Lept1E * u2 * Xsi2 * GG / I_Gamma;
00171   G4double L_Gamma   = GammaE * ((Lept0E + Lept1E) * (3 + 2 * GG) 
00172                                  - 2 * Lept1E * (1 + 4 * u2 * Xsi2 * GG))/I_Gamma;   
00173   G4double T_Gamma   = 4 * GammaE * Lept1E * Xsi * u * (2 * Xsi - 1) * GG / I_Gamma ;
00174     
00175   G4double Stokes_P1 = D_Gamma ;
00176   G4double Stokes_P2 = 0 ;
00177   G4double Stokes_P3 = (Stokes_S3*L_Gamma + Stokes_S1*T_Gamma) ;
00178 
00179   theFinalGammaPolarization.SetPhoton();
00180 
00181   theFinalGammaPolarization.setX(Stokes_P1);
00182   theFinalGammaPolarization.setY(Stokes_P2);
00183   theFinalGammaPolarization.setZ(Stokes_P3);
00184     
00185   if(theFinalGammaPolarization.mag2()>1) {
00186     G4cout<<" WARNING in pol-brem theFinalGammaPolarization \n";
00187     G4cout
00188       <<"\t"<<theFinalGammaPolarization
00189       <<"\t GG\t"<<GG
00190       <<"\t delta\t"<<delta
00191       <<G4endl;
00192   }
00193 }
00194 
00195 G4double G4PolarizedBremsstrahlungCrossSection::XSection(const G4StokesVector & /*pol2*/,
00196                                                          const G4StokesVector & /*pol3*/)
00197 {
00198   G4cout<<"ERROR dummy routine G4PolarizedBremsstrahlungCrossSection::XSection called \n";
00199   return 0;
00200 }
00201 
00202   // return expected mean polarisation
00203 G4StokesVector G4PolarizedBremsstrahlungCrossSection::GetPol2()
00204 {
00205   // electron/positron
00206   return theFinalLeptonPolarization;
00207 }
00208 G4StokesVector G4PolarizedBremsstrahlungCrossSection::GetPol3()
00209 {
00210   // photon
00211   return theFinalGammaPolarization;;
00212 }
00213 
00214 

Generated on Mon May 27 17:49:22 2013 for Geant4 by  doxygen 1.4.7