G4PolarizedPEEffectModel Class Reference

#include <G4PolarizedPEEffectModel.hh>

Inheritance diagram for G4PolarizedPEEffectModel:

G4PEEffectFluoModel G4VEmModel

Public Member Functions

 G4PolarizedPEEffectModel (const G4ParticleDefinition *p=0, const G4String &nam="Polarized-PhotoElectric")
void Initialise (const G4ParticleDefinition *pd, const G4DataVector &dv)
virtual ~G4PolarizedPEEffectModel ()
virtual void SampleSecondaries (std::vector< G4DynamicParticle * > *, const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double tmin, G4double maxEnergy)

Detailed Description

Definition at line 64 of file G4PolarizedPEEffectModel.hh.


Constructor & Destructor Documentation

G4PolarizedPEEffectModel::G4PolarizedPEEffectModel ( const G4ParticleDefinition p = 0,
const G4String nam = "Polarized-PhotoElectric" 
)

Definition at line 68 of file G4PolarizedPEEffectModel.cc.

00070   : G4PEEffectFluoModel(nam),crossSectionCalculator(0),verboseLevel(0)
00071 {
00072 }

G4PolarizedPEEffectModel::~G4PolarizedPEEffectModel (  )  [virtual]

Definition at line 76 of file G4PolarizedPEEffectModel.cc.

00077 {
00078   if (crossSectionCalculator) delete crossSectionCalculator;
00079 }


Member Function Documentation

void G4PolarizedPEEffectModel::Initialise ( const G4ParticleDefinition pd,
const G4DataVector dv 
) [virtual]

Reimplemented from G4PEEffectFluoModel.

Definition at line 83 of file G4PolarizedPEEffectModel.cc.

References G4PEEffectFluoModel::Initialise().

00085 {
00086   G4PEEffectFluoModel::Initialise(pd,dv);
00087   if (!crossSectionCalculator)
00088     crossSectionCalculator = new G4PolarizedPEEffectCrossSection();
00089 }

void G4PolarizedPEEffectModel::SampleSecondaries ( std::vector< G4DynamicParticle * > *  ,
const G4MaterialCutsCouple ,
const G4DynamicParticle ,
G4double  tmin,
G4double  maxEnergy 
) [virtual]

Reimplemented from G4PEEffectFluoModel.

Definition at line 93 of file G4PolarizedPEEffectModel.cc.

References G4cout, G4VEmModel::GetCurrentElement(), G4PolarizationHelper::GetFrame(), G4DynamicParticle::GetKineticEnergy(), G4DynamicParticle::GetMomentumDirection(), G4PolarizedPEEffectCrossSection::GetPol2(), G4DynamicParticle::GetPolarization(), G4PolarizationHelper::GetRandomFrame(), G4PolarizedPEEffectCrossSection::Initialize(), G4StokesVector::InvRotateAz(), G4StokesVector::p1(), G4StokesVector::p2(), G4StokesVector::p3(), G4StokesVector::RotateAz(), G4PEEffectFluoModel::SampleSecondaries(), G4VPolarizedCrossSection::SetMaterial(), G4StokesVector::SetPhoton(), and G4StokesVector::ZERO.

00098 {
00099   //  std::vector<G4DynamicParticle*>* vdp = 
00100   G4PEEffectFluoModel::SampleSecondaries(vdp,couple, dp, tmin, maxEnergy);
00101 
00102   if(vdp && vdp->size()>0) {
00103     G4double gamEnergy0 = dp->GetKineticEnergy();
00104     G4double lepEnergy1 = (*vdp)[0]->GetKineticEnergy();
00105     G4double sintheta   = dp->GetMomentumDirection().cross((*vdp)[0]->GetMomentumDirection()).mag();
00106     if (sintheta>1.) sintheta=1.;
00107 
00108     G4StokesVector beamPol = dp->GetPolarization();
00109     beamPol.SetPhoton();
00110     //    G4cout<<" beamPol "<<beamPol<<G4endl;
00111 
00112     // determine interaction plane
00113     G4ThreeVector  nInteractionFrame = 
00114       G4PolarizationHelper::GetFrame(dp->GetMomentumDirection(), 
00115                                      (*vdp)[0]->GetMomentumDirection());
00116     //    G4cout<<" nInteractionFrame = "<<nInteractionFrame<<G4endl;
00117     if (dp->GetMomentumDirection().cross((*vdp)[0]->GetMomentumDirection()).mag()<1.e-10) {
00118       //      G4cout<<" nInteractionFrame not well defined "<<G4endl;
00119       //      G4cout<<" choosing random interaction frame "<<G4endl;
00120 
00121       nInteractionFrame = G4PolarizationHelper::GetRandomFrame(dp->GetMomentumDirection());
00122       //      G4cout<<"new nInteractionFrame = "<<nInteractionFrame<<G4endl;
00123     }
00124     /*
00125     else {
00126       G4ThreeVector mom1=dp->GetMomentumDirection();
00127       G4ThreeVector mom2=(*vdp)[0]->GetMomentumDirection();
00128       G4cout<<" mom1 = "<<mom1<<G4endl;
00129       G4cout<<" mom2 = "<<mom2<<G4endl;
00130       G4ThreeVector x=mom1.cross(mom2);
00131       G4cout<<" mom1 x mom2 = "<<x<<" "<<x.mag()<<G4endl;
00132       G4cout<<"  norm       = "<<(1./x.mag()*x)<<" "<<G4endl;
00133     }
00134     */
00135 
00136     // transform polarization into interaction frame
00137     beamPol.InvRotateAz(nInteractionFrame,dp->GetMomentumDirection());
00138 
00139     // calulcate polarization transfer
00140     crossSectionCalculator->SetMaterial(GetCurrentElement()->GetN(), // number of nucleons
00141                                         GetCurrentElement()->GetZ(), 
00142                                         GetCurrentElement()->GetfCoulomb());
00143     crossSectionCalculator->Initialize(gamEnergy0, lepEnergy1, sintheta,
00144                                        beamPol, G4StokesVector::ZERO);
00145 
00146     // deterimine final state polarization
00147     G4StokesVector lep1Pol = crossSectionCalculator->GetPol2();
00148     //    G4cout<<" lepPol "<<lep1Pol<<G4endl;
00149     lep1Pol.RotateAz(nInteractionFrame,(*vdp)[0]->GetMomentumDirection());
00150     (*vdp)[0]->SetPolarization(lep1Pol.p1(),
00151                                lep1Pol.p2(),
00152                                lep1Pol.p3());
00153 
00154     //    G4cout<<" lepPol "<<lep1Pol<<G4endl;
00155     size_t num = vdp->size();
00156     if (num!=1) G4cout<<" WARNING "<<num<<" secondaries in polarized photo electric effect not supported!\n"; 
00157   }
00158 }


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