GammaComptonRange.cpp Source File

00001 #include <math.h>
00002 #include <iostream.h>
00003 #include "ReactorConstants.hh"
00004 
00005 double GammaComptonRange(double egamma){
00006   //
00007   // Compute the mean free path in scint. for Compton scattering.
00008   //
00009   static ReactorConstants k;
00010   static double pi = k.pi;
00011   static double m = k.Melectron;
00012   static double alpha = 1/k.alphainv;
00013   static double Units = k.XcMeVtoCmsqrd;
00014   //
00015   // Overall cross section scale
00016   //
00017   static double Thomson = 8*pi*alpha*alpha/3/m/m*Units;
00018   double w = egamma/m;
00019   //
00020   // Klein-Nishina cross section
00021   //
00022   double XC = 3*Thomson/4*
00023     ((1+w)/w/w/w*(2*w*(1+w)/(1+2*w)-log(1+2*w))
00024      +log(1+2*w)/2/w-(1+3*w)/(1+2*w)/(1+2*w));
00025   //
00026   // Need electron density to get mfp. Get H and C contributions.
00027   //
00028   static double density = k.density;
00029   static double Na = k.Navogadro;
00030   static double A[2] = {k.A0,k.A1};
00031   static double Z[2] = {k.Z0,k.Z1};
00032   static double f[2] = {k.f0,k.f1};
00033   static double nel[2];
00034   static bool first = 1;
00035   for(int i=0;i<2;i++){
00036     nel[i]=Na*density*f[i]*Z[i]/A[i];
00037     first = 0;
00038   }
00039   //
00040   //  Compute mean free path.
00041   //  
00042   return 1/((nel[0]+nel[1])*XC);
00043 }
00044