G4SynchrotronRadiationInMat Class Reference

#include <G4SynchrotronRadiationInMat.hh>

Inheritance diagram for G4SynchrotronRadiationInMat:

G4VDiscreteProcess G4VProcess

Public Member Functions

 G4SynchrotronRadiationInMat (const G4String &processName="SynchrotronRadiation", G4ProcessType type=fElectromagnetic)
virtual ~G4SynchrotronRadiationInMat ()
G4double GetMeanFreePath (const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)
G4VParticleChangePostStepDoIt (const G4Track &track, const G4Step &Step)
G4double GetPhotonEnergy (const G4Track &trackData, const G4Step &stepData)
G4double GetRandomEnergySR (G4double, G4double)
G4double GetProbSpectrumSRforInt (G4double)
G4double GetIntProbSR (G4double)
G4double GetProbSpectrumSRforEnergy (G4double)
G4double GetEnergyProbSR (G4double)
G4double GetIntegrandForAngleK (G4double)
G4double GetAngleK (G4double)
G4double GetAngleNumberAtGammaKsi (G4double)
G4bool IsApplicable (const G4ParticleDefinition &)
void SetRootNumber (G4int rn)
void SetVerboseLevel (G4int v)
void SetKsi (G4double ksi)
void SetEta (G4double eta)
void SetPsiGamma (G4double psg)
void SetOrderAngleK (G4double ord)

Static Public Member Functions

static G4double GetLambdaConst ()
static G4double GetEnergyConst ()

Detailed Description

Definition at line 68 of file G4SynchrotronRadiationInMat.hh.


Constructor & Destructor Documentation

G4SynchrotronRadiationInMat::G4SynchrotronRadiationInMat ( const G4String processName = "SynchrotronRadiation",
G4ProcessType  type = fElectromagnetic 
)

Definition at line 123 of file G4SynchrotronRadiationInMat.cc.

References fSynchrotronRadiation, G4TransportationManager::GetPropagatorInField(), G4TransportationManager::GetTransportationManager(), and G4VProcess::SetProcessSubType().

00124                      :G4VDiscreteProcess (processName, type),
00125   LowestKineticEnergy (10.*keV),
00126   HighestKineticEnergy (100.*TeV),
00127   TotBin(200),
00128   theGamma (G4Gamma::Gamma() ),
00129   theElectron ( G4Electron::Electron() ),
00130   thePositron ( G4Positron::Positron() ), 
00131   GammaCutInKineticEnergy(0),
00132   ElectronCutInKineticEnergy(0),
00133   PositronCutInKineticEnergy(0),
00134   ParticleCutInKineticEnergy(0),
00135   fAlpha(0.0), fRootNumber(80),
00136   fVerboseLevel( verboseLevel )
00137 {
00138   G4TransportationManager* transportMgr = G4TransportationManager::GetTransportationManager();
00139 
00140   fFieldPropagator = transportMgr->GetPropagatorInField();
00141   SetProcessSubType(fSynchrotronRadiation);
00142   CutInRange = GammaCutInKineticEnergyNow = ElectronCutInKineticEnergyNow = 
00143     PositronCutInKineticEnergyNow =  ParticleCutInKineticEnergyNow = fKsi = 
00144     fPsiGamma = fEta = fOrderAngleK = 0.0;
00145 }

G4SynchrotronRadiationInMat::~G4SynchrotronRadiationInMat (  )  [virtual]

Definition at line 152 of file G4SynchrotronRadiationInMat.cc.

00153 {}


Member Function Documentation

G4double G4SynchrotronRadiationInMat::GetAngleK ( G4double   ) 

Definition at line 636 of file G4SynchrotronRadiationInMat.cc.

References GetIntegrandForAngleK(), G4Integrator< T, F >::Laguerre(), and CLHEP::detail::n.

Referenced by GetAngleNumberAtGammaKsi().

00637 {
00638   fEta = eta; // should be > 0. !
00639   G4int n;
00640   G4double result, a;
00641 
00642   a = fAlpha;        // always = 0.
00643   n = fRootNumber;   // around default = 80
00644 
00645   G4Integrator<G4SynchrotronRadiationInMat, G4double(G4SynchrotronRadiationInMat::*)(G4double)> integral;
00646 
00647   result = integral.Laguerre(this, 
00648            &G4SynchrotronRadiationInMat::GetIntegrandForAngleK, a, n);
00649 
00650   return result;
00651 }

G4double G4SynchrotronRadiationInMat::GetAngleNumberAtGammaKsi ( G4double   ) 

Definition at line 657 of file G4SynchrotronRadiationInMat.cc.

References GetAngleK().

00658 {
00659   G4double result, funK, funK2, gpsi2 = gpsi*gpsi;
00660 
00661   fPsiGamma    = gpsi;
00662   fEta         = 0.5*fKsi*(1 + gpsi2)*std::sqrt(1 + gpsi2);
00663  
00664   fOrderAngleK = 1./3.;
00665   funK         = GetAngleK(fEta); 
00666   funK2        = funK*funK;
00667 
00668   result       = gpsi2*funK2/(1 + gpsi2);
00669 
00670   fOrderAngleK = 2./3.;
00671   funK         = GetAngleK(fEta); 
00672   funK2        = funK*funK;
00673 
00674   result      += funK2;
00675   result      *= (1 + gpsi2)*fKsi;
00676      
00677   return result;
00678 }

static G4double G4SynchrotronRadiationInMat::GetEnergyConst (  )  [inline, static]

Definition at line 110 of file G4SynchrotronRadiationInMat.hh.

00110 { return fEnergyConst; };

G4double G4SynchrotronRadiationInMat::GetEnergyProbSR ( G4double   ) 

Definition at line 599 of file G4SynchrotronRadiationInMat.cc.

References GetProbSpectrumSRforEnergy(), G4Integrator< T, F >::Laguerre(), CLHEP::detail::n, and G4INCL::Math::pi.

00600 {
00601   if (ksi <= 0.) return 1.0;
00602   fKsi = ksi; // should be > 0. !
00603   G4int n;
00604   G4double result, a;
00605 
00606   a = fAlpha;        // always = 0.
00607   n = fRootNumber;   // around default = 80
00608 
00609   G4Integrator<G4SynchrotronRadiationInMat, G4double(G4SynchrotronRadiationInMat::*)(G4double)> integral;
00610 
00611   result = integral.Laguerre(this, 
00612            &G4SynchrotronRadiationInMat::GetProbSpectrumSRforEnergy, a, n);
00613 
00614   result *= 9.*std::sqrt(3.)*ksi/8./pi;
00615 
00616   return result;
00617 }

G4double G4SynchrotronRadiationInMat::GetIntegrandForAngleK ( G4double   ) 

Definition at line 623 of file G4SynchrotronRadiationInMat.cc.

Referenced by GetAngleK().

00624 {
00625   G4double result, hypCos=std::cosh(t);
00626   
00627   result  = std::cosh(fOrderAngleK*t)*std::exp(t - fEta*hypCos); // fEta > 0. !
00628   result /= hypCos;
00629   return result;
00630 }

G4double G4SynchrotronRadiationInMat::GetIntProbSR ( G4double   ) 

Definition at line 560 of file G4SynchrotronRadiationInMat.cc.

References GetProbSpectrumSRforInt(), G4Integrator< T, F >::Laguerre(), CLHEP::detail::n, and G4INCL::Math::pi.

00561 {
00562   if (ksi <= 0.) return 1.0;
00563   fKsi = ksi; // should be > 0. !
00564   G4int n;
00565   G4double result, a;
00566 
00567   a = fAlpha;        // always = 0.
00568   n = fRootNumber;   // around default = 80
00569 
00570   G4Integrator<G4SynchrotronRadiationInMat, G4double(G4SynchrotronRadiationInMat::*)(G4double)> integral;
00571 
00572   result = integral.Laguerre(this, 
00573            &G4SynchrotronRadiationInMat::GetProbSpectrumSRforInt, a, n);
00574 
00575   result *= 3./5./pi;
00576 
00577   return result;
00578 }

static G4double G4SynchrotronRadiationInMat::GetLambdaConst (  )  [inline, static]

Definition at line 109 of file G4SynchrotronRadiationInMat.hh.

00109 { return fLambdaConst; };

G4double G4SynchrotronRadiationInMat::GetMeanFreePath ( const G4Track track,
G4double  previousStepSize,
G4ForceCondition condition 
) [virtual]

Implements G4VDiscreteProcess.

Definition at line 174 of file G4SynchrotronRadiationInMat.cc.

References DBL_MAX, G4PropagatorInField::FindAndSetFieldManager(), G4cout, G4endl, G4DynamicParticle::GetDefinition(), G4FieldManager::GetDetectorField(), G4Track::GetDynamicParticle(), G4Field::GetFieldValue(), G4Track::GetGlobalTime(), G4DynamicParticle::GetKineticEnergy(), G4DynamicParticle::GetMass(), G4DynamicParticle::GetMomentumDirection(), G4ParticleDefinition::GetPDGCharge(), G4Track::GetPosition(), G4DynamicParticle::GetTotalEnergy(), G4DynamicParticle::GetTotalMomentum(), G4Track::GetVolume(), and NotForced.

00177 {
00178   // gives the MeanFreePath in GEANT4 internal units
00179   G4double MeanFreePath;
00180 
00181   const G4DynamicParticle* aDynamicParticle = trackData.GetDynamicParticle();
00182   // G4Material* aMaterial = trackData.GetMaterial();
00183 
00184   //G4bool isOutRange ;
00185  
00186   *condition = NotForced ;
00187 
00188   G4double gamma = aDynamicParticle->GetTotalEnergy()/
00189                    aDynamicParticle->GetMass();
00190 
00191   G4double particleCharge = aDynamicParticle->GetDefinition()->GetPDGCharge();
00192 
00193   G4double KineticEnergy = aDynamicParticle->GetKineticEnergy();
00194 
00195   if ( KineticEnergy <  LowestKineticEnergy || gamma < 1.0e3 )  MeanFreePath = DBL_MAX; 
00196   else
00197   {
00198 
00199     G4ThreeVector  FieldValue;
00200     const G4Field*   pField = 0;
00201 
00202     G4FieldManager* fieldMgr=0;
00203     G4bool          fieldExertsForce = false;
00204 
00205     if( (particleCharge != 0.0) )
00206     {
00207       fieldMgr = fFieldPropagator->FindAndSetFieldManager( trackData.GetVolume() );
00208 
00209       if ( fieldMgr != 0 ) 
00210       {
00211         // If the field manager has no field, there is no field !
00212 
00213         fieldExertsForce = ( fieldMgr->GetDetectorField() != 0 );
00214       }
00215     }
00216     if ( fieldExertsForce )
00217     {  
00218       pField = fieldMgr->GetDetectorField() ;
00219       G4ThreeVector  globPosition = trackData.GetPosition();
00220 
00221       G4double  globPosVec[4], FieldValueVec[6];
00222 
00223       globPosVec[0] = globPosition.x();
00224       globPosVec[1] = globPosition.y();
00225       globPosVec[2] = globPosition.z();
00226       globPosVec[3] = trackData.GetGlobalTime();
00227 
00228       pField->GetFieldValue( globPosVec, FieldValueVec );
00229 
00230       FieldValue = G4ThreeVector( FieldValueVec[0], 
00231                                    FieldValueVec[1], 
00232                                    FieldValueVec[2]   );
00233 
00234        
00235 
00236       G4ThreeVector unitMomentum = aDynamicParticle->GetMomentumDirection(); 
00237       G4ThreeVector unitMcrossB  = FieldValue.cross(unitMomentum) ;
00238       G4double perpB             = unitMcrossB.mag() ;
00239       G4double beta              = aDynamicParticle->GetTotalMomentum()/
00240                                    (aDynamicParticle->GetTotalEnergy()    );
00241 
00242       if( perpB > 0.0 ) MeanFreePath = fLambdaConst*beta/perpB;
00243       else              MeanFreePath = DBL_MAX;       
00244     }
00245     else  MeanFreePath = DBL_MAX;    
00246   }
00247   if(fVerboseLevel > 0)
00248   {
00249     G4cout<<"G4SynchrotronRadiationInMat::MeanFreePath = "<<MeanFreePath/m<<" m"<<G4endl;
00250   }
00251   return MeanFreePath; 
00252 } 

G4double G4SynchrotronRadiationInMat::GetPhotonEnergy ( const G4Track trackData,
const G4Step stepData 
)

Definition at line 421 of file G4SynchrotronRadiationInMat.cc.

References G4PropagatorInField::FindAndSetFieldManager(), G4UniformRand, G4DynamicParticle::GetDefinition(), G4FieldManager::GetDetectorField(), G4Track::GetDynamicParticle(), G4Field::GetFieldValue(), G4DynamicParticle::GetMass(), G4DynamicParticle::GetMomentumDirection(), G4ParticleDefinition::GetPDGCharge(), G4Track::GetPosition(), G4DynamicParticle::GetTotalEnergy(), and G4Track::GetVolume().

00424 {
00425   G4int i ;
00426   G4double energyOfSR = -1.0 ;
00427   //G4Material* aMaterial=trackData.GetMaterial() ;
00428 
00429   const G4DynamicParticle* aDynamicParticle=trackData.GetDynamicParticle();
00430 
00431   G4double gamma = aDynamicParticle->GetTotalEnergy()/
00432                    (aDynamicParticle->GetMass()              ) ;
00433 
00434   G4double particleCharge = aDynamicParticle->GetDefinition()->GetPDGCharge();
00435 
00436   G4ThreeVector  FieldValue;
00437   const G4Field*   pField = 0 ;
00438 
00439   G4FieldManager* fieldMgr=0;
00440   G4bool          fieldExertsForce = false;
00441 
00442   if( (particleCharge != 0.0) )
00443   {
00444     fieldMgr = fFieldPropagator->FindAndSetFieldManager( trackData.GetVolume() );
00445     if ( fieldMgr != 0 ) 
00446     {
00447       // If the field manager has no field, there is no field !
00448 
00449       fieldExertsForce = ( fieldMgr->GetDetectorField() != 0 );
00450     }
00451   }
00452   if ( fieldExertsForce )
00453   {  
00454     pField = fieldMgr->GetDetectorField();
00455     G4ThreeVector  globPosition = trackData.GetPosition();
00456     G4double  globPosVec[3], FieldValueVec[3];
00457 
00458     globPosVec[0] = globPosition.x();
00459     globPosVec[1] = globPosition.y();
00460     globPosVec[2] = globPosition.z();
00461 
00462     pField->GetFieldValue( globPosVec, FieldValueVec );
00463     FieldValue = G4ThreeVector( FieldValueVec[0], 
00464                                    FieldValueVec[1], 
00465                                    FieldValueVec[2]   );
00466        
00467     G4ThreeVector unitMomentum = aDynamicParticle->GetMomentumDirection(); 
00468     G4ThreeVector unitMcrossB  = FieldValue.cross(unitMomentum) ;
00469     G4double perpB             = unitMcrossB.mag();
00470     if( perpB > 0.0 )
00471     {
00472       // M-C of synchrotron photon energy
00473 
00474       G4double random = G4UniformRand() ;
00475       for(i=0;i<200;i++)
00476       {
00477         if(random >= fIntegralProbabilityOfSR[i]) break ;
00478       }
00479       energyOfSR = 0.0001*i*i*fEnergyConst*gamma*gamma*perpB ;
00480 
00481       // check against insufficient energy
00482 
00483       if(energyOfSR <= 0.0)
00484       {
00485         return -1.0 ;
00486       }
00487       //G4double kineticEnergy = aDynamicParticle->GetKineticEnergy();
00488       //G4ParticleMomentum 
00489       //particleDirection = aDynamicParticle->GetMomentumDirection();
00490 
00491       // Gamma production cut in this material
00492       //G4double 
00493       //gammaEnergyCut = (G4Gamma::GetCutsInEnergy())[aMaterial->GetIndex()];
00494 
00495       // SR photon has energy more than the current material cut
00496       // M-C of its direction
00497       
00498       //G4double Teta = G4UniformRand()/gamma ;    // Very roughly
00499 
00500       //G4double Phi  = twopi * G4UniformRand() ;
00501     }       
00502     else
00503     {
00504       return -1.0 ;
00505     }
00506   }     
00507   return energyOfSR ;
00508 }

G4double G4SynchrotronRadiationInMat::GetProbSpectrumSRforEnergy ( G4double   ) 

Definition at line 584 of file G4SynchrotronRadiationInMat.cc.

Referenced by GetEnergyProbSR().

00585 {
00586   G4double result, hypCos=std::cosh(t);
00587   
00588   result = std::cosh(5.*t/3.)*std::exp(t - fKsi*hypCos); // fKsi > 0. !
00589   result /= hypCos;
00590   return result;
00591 }

G4double G4SynchrotronRadiationInMat::GetProbSpectrumSRforInt ( G4double   ) 

Definition at line 544 of file G4SynchrotronRadiationInMat.cc.

Referenced by GetIntProbSR().

00545 {
00546   G4double result, hypCos2, hypCos=std::cosh(t);
00547 
00548   hypCos2 = hypCos*hypCos;  
00549   result = std::cosh(5.*t/3.)*std::exp(t-fKsi*hypCos); // fKsi > 0. !
00550   result /= hypCos2;
00551   return result;
00552 }

G4double G4SynchrotronRadiationInMat::GetRandomEnergySR ( G4double  ,
G4double   
)

Definition at line 514 of file G4SynchrotronRadiationInMat.cc.

References G4UniformRand, and position.

Referenced by PostStepDoIt().

00515 {
00516   G4int i, iMax;
00517   G4double energySR, random, position;
00518 
00519   iMax = 200;
00520   random = G4UniformRand();
00521 
00522   for( i = 0; i < iMax; i++ )
00523   {
00524     if( random >= fIntegralProbabilityOfSR[i] ) break;
00525   }
00526   if(i <= 0 )        position = G4UniformRand(); // 0.
00527   else if( i>= iMax) position = G4double(iMax);
00528   else position = i + G4UniformRand(); // -1
00529   // 
00530   // it was in initial implementation:
00531   //       energyOfSR = 0.0001*i*i*fEnergyConst*gamma*gamma*perpB ;
00532 
00533   energySR = 0.0001*position*position*fEnergyConst*gamma*gamma*perpB;
00534 
00535   if( energySR  < 0. ) energySR = 0.;
00536 
00537   return energySR;
00538 }

G4bool G4SynchrotronRadiationInMat::IsApplicable ( const G4ParticleDefinition  )  [virtual]

Reimplemented from G4VProcess.

Definition at line 157 of file G4SynchrotronRadiationInMat.cc.

00158 {
00159 
00160   return ( ( &particle == (const G4ParticleDefinition *)theElectron ) ||
00161            ( &particle == (const G4ParticleDefinition *)thePositron )    );
00162 
00163 }

G4VParticleChange * G4SynchrotronRadiationInMat::PostStepDoIt ( const G4Track track,
const G4Step Step 
) [virtual]

Reimplemented from G4VDiscreteProcess.

Definition at line 259 of file G4SynchrotronRadiationInMat.cc.

References G4ParticleChange::AddSecondary(), G4VProcess::aParticleChange, fcos(), G4PropagatorInField::FindAndSetFieldManager(), fStopAndKill, fStopButAlive, G4cout, G4endl, G4UniformRand, G4Gamma::Gamma(), G4DynamicParticle::GetDefinition(), G4FieldManager::GetDetectorField(), G4Track::GetDynamicParticle(), G4Field::GetFieldValue(), G4Track::GetGlobalTime(), G4DynamicParticle::GetKineticEnergy(), G4DynamicParticle::GetMass(), G4DynamicParticle::GetMomentumDirection(), G4ParticleDefinition::GetPDGCharge(), G4Track::GetPosition(), GetRandomEnergySR(), G4DynamicParticle::GetTotalEnergy(), G4Track::GetVolume(), G4ParticleChange::Initialize(), G4VDiscreteProcess::PostStepDoIt(), G4ParticleChange::ProposeEnergy(), G4VParticleChange::ProposeLocalEnergyDeposit(), G4ParticleChange::ProposeMomentumDirection(), G4VParticleChange::ProposeTrackStatus(), G4VParticleChange::SetNumberOfSecondaries(), and G4DynamicParticle::SetPolarization().

00262 {
00263   aParticleChange.Initialize(trackData);
00264 
00265   const G4DynamicParticle* aDynamicParticle=trackData.GetDynamicParticle();
00266 
00267   G4double gamma = aDynamicParticle->GetTotalEnergy()/
00268                    (aDynamicParticle->GetMass()              );
00269 
00270   if(gamma <= 1.0e3 ) 
00271   {
00272     return G4VDiscreteProcess::PostStepDoIt(trackData,stepData);
00273   }
00274   G4double particleCharge = aDynamicParticle->GetDefinition()->GetPDGCharge();
00275 
00276   G4ThreeVector  FieldValue;
00277   const G4Field*   pField = 0 ;
00278 
00279   G4FieldManager* fieldMgr=0;
00280   G4bool          fieldExertsForce = false;
00281 
00282   if( (particleCharge != 0.0) )
00283   {
00284     fieldMgr = fFieldPropagator->FindAndSetFieldManager( trackData.GetVolume() );
00285     if ( fieldMgr != 0 ) 
00286     {
00287       // If the field manager has no field, there is no field !
00288 
00289       fieldExertsForce = ( fieldMgr->GetDetectorField() != 0 );
00290     }
00291   }
00292   if ( fieldExertsForce )
00293   {
00294     pField = fieldMgr->GetDetectorField() ;
00295     G4ThreeVector  globPosition = trackData.GetPosition() ;
00296     G4double  globPosVec[4], FieldValueVec[6] ;
00297     globPosVec[0] = globPosition.x() ;
00298     globPosVec[1] = globPosition.y() ;
00299     globPosVec[2] = globPosition.z() ;
00300     globPosVec[3] = trackData.GetGlobalTime();
00301 
00302     pField->GetFieldValue( globPosVec, FieldValueVec ) ;
00303     FieldValue = G4ThreeVector( FieldValueVec[0], 
00304                                    FieldValueVec[1], 
00305                                    FieldValueVec[2]   );
00306        
00307     G4ThreeVector unitMomentum = aDynamicParticle->GetMomentumDirection(); 
00308     G4ThreeVector unitMcrossB = FieldValue.cross(unitMomentum);
00309     G4double perpB = unitMcrossB.mag() ;
00310     if(perpB > 0.0)
00311     {
00312       // M-C of synchrotron photon energy
00313 
00314       G4double energyOfSR = GetRandomEnergySR(gamma,perpB);
00315 
00316       if(fVerboseLevel > 0)
00317       {
00318         G4cout<<"SR photon energy = "<<energyOfSR/keV<<" keV"<<G4endl;
00319       }
00320       // check against insufficient energy
00321 
00322       if( energyOfSR <= 0.0 )
00323       {
00324         return G4VDiscreteProcess::PostStepDoIt(trackData,stepData);
00325       }
00326       G4double kineticEnergy = aDynamicParticle->GetKineticEnergy();
00327       G4ParticleMomentum 
00328       particleDirection = aDynamicParticle->GetMomentumDirection();
00329 
00330       // M-C of its direction, simplified dipole busted approach
00331       
00332       // G4double Teta = G4UniformRand()/gamma ;    // Very roughly
00333 
00334       G4double cosTheta, sinTheta, fcos, beta;
00335 
00336   do
00337   { 
00338     cosTheta = 1. - 2.*G4UniformRand();
00339     fcos     = (1 + cosTheta*cosTheta)*0.5;
00340   }
00341   while( fcos < G4UniformRand() );
00342 
00343   beta = std::sqrt(1. - 1./(gamma*gamma));
00344 
00345   cosTheta = (cosTheta + beta)/(1. + beta*cosTheta);
00346 
00347   if( cosTheta >  1. ) cosTheta =  1.;
00348   if( cosTheta < -1. ) cosTheta = -1.;
00349 
00350   sinTheta = std::sqrt(1. - cosTheta*cosTheta );
00351 
00352       G4double Phi  = twopi * G4UniformRand() ;
00353 
00354       G4double dirx = sinTheta*std::cos(Phi) , 
00355                diry = sinTheta*std::sin(Phi) , 
00356                dirz = cosTheta;
00357 
00358       G4ThreeVector gammaDirection ( dirx, diry, dirz);
00359       gammaDirection.rotateUz(particleDirection);   
00360  
00361       // polarization of new gamma
00362 
00363       // G4double sx =  std::cos(Teta)*std::cos(Phi);
00364       // G4double sy =  std::cos(Teta)*std::sin(Phi);
00365       // G4double sz = -std::sin(Teta);
00366 
00367       G4ThreeVector gammaPolarization = FieldValue.cross(gammaDirection);
00368       gammaPolarization = gammaPolarization.unit();
00369 
00370       // (sx, sy, sz);
00371       // gammaPolarization.rotateUz(particleDirection);
00372 
00373       // create G4DynamicParticle object for the SR photon
00374  
00375       G4DynamicParticle* aGamma= new G4DynamicParticle ( G4Gamma::Gamma(),
00376                                                          gammaDirection, 
00377                                                          energyOfSR      );
00378       aGamma->SetPolarization( gammaPolarization.x(),
00379                                gammaPolarization.y(),
00380                                gammaPolarization.z() );
00381 
00382 
00383       aParticleChange.SetNumberOfSecondaries(1);
00384       aParticleChange.AddSecondary(aGamma); 
00385 
00386       // Update the incident particle 
00387    
00388       G4double newKinEnergy = kineticEnergy - energyOfSR ;
00389       
00390       if (newKinEnergy > 0.)
00391       {
00392         aParticleChange.ProposeMomentumDirection( particleDirection );
00393         aParticleChange.ProposeEnergy( newKinEnergy );
00394         aParticleChange.ProposeLocalEnergyDeposit (0.); 
00395       } 
00396       else
00397       { 
00398         aParticleChange.ProposeEnergy( 0. );
00399         aParticleChange.ProposeLocalEnergyDeposit (0.);
00400         G4double charge = aDynamicParticle->GetDefinition()->GetPDGCharge();   
00401         if (charge<0.)
00402         {
00403            aParticleChange.ProposeTrackStatus(fStopAndKill) ;
00404         }
00405         else      
00406         {
00407           aParticleChange.ProposeTrackStatus(fStopButAlive) ;
00408         }
00409       } 
00410     }       
00411     else
00412     {
00413       return G4VDiscreteProcess::PostStepDoIt(trackData,stepData);
00414     }
00415   }     
00416   return G4VDiscreteProcess::PostStepDoIt(trackData,stepData);
00417 }

void G4SynchrotronRadiationInMat::SetEta ( G4double  eta  )  [inline]

Definition at line 115 of file G4SynchrotronRadiationInMat.hh.

00115 { fEta = eta; };

void G4SynchrotronRadiationInMat::SetKsi ( G4double  ksi  )  [inline]

Definition at line 114 of file G4SynchrotronRadiationInMat.hh.

00114 { fKsi = ksi; };

void G4SynchrotronRadiationInMat::SetOrderAngleK ( G4double  ord  )  [inline]

Definition at line 117 of file G4SynchrotronRadiationInMat.hh.

00117 { fOrderAngleK = ord; }; // should be 1/3 or 2/3

void G4SynchrotronRadiationInMat::SetPsiGamma ( G4double  psg  )  [inline]

Definition at line 116 of file G4SynchrotronRadiationInMat.hh.

00116 { fPsiGamma = psg; };

void G4SynchrotronRadiationInMat::SetRootNumber ( G4int  rn  )  [inline]

Definition at line 112 of file G4SynchrotronRadiationInMat.hh.

00112 { fRootNumber = rn; };

void G4SynchrotronRadiationInMat::SetVerboseLevel ( G4int  v  )  [inline]

Reimplemented from G4VProcess.

Definition at line 113 of file G4SynchrotronRadiationInMat.hh.

00113 { fVerboseLevel = v; };


The documentation for this class was generated from the following files:
Generated on Mon May 27 17:53:28 2013 for Geant4 by  doxygen 1.4.7