Geant4.10
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Public Member Functions
G4PolarizedAnnihilationModel Class Reference

#include <G4PolarizedAnnihilationModel.hh>

Inheritance diagram for G4PolarizedAnnihilationModel:
G4eeToTwoGammaModel G4VEmModel

Public Member Functions

 G4PolarizedAnnihilationModel (const G4ParticleDefinition *p=0, const G4String &nam="Polarized-Annihilation")
 
virtual ~G4PolarizedAnnihilationModel ()
 
virtual void Initialise (const G4ParticleDefinition *, const G4DataVector &)
 
virtual G4double ComputeCrossSectionPerElectron (const G4ParticleDefinition *, G4double kinEnergy, G4double cut, G4double emax)
 
void ComputeAsymmetriesPerElectron (G4double gammaEnergy, G4double &valueX, G4double &valueA, G4double &valueT)
 
virtual void SampleSecondaries (std::vector< G4DynamicParticle * > *, const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double tmin, G4double maxEnergy)
 
void SetTargetPolarization (const G4ThreeVector &pTarget)
 
void SetBeamPolarization (const G4ThreeVector &pBeam)
 
const G4ThreeVectorGetTargetPolarization () const
 
const G4ThreeVectorGetBeamPolarization () const
 
const G4ThreeVectorGetFinalGamma1Polarization () const
 
const G4ThreeVectorGetFinalGamma2Polarization () const
 
- Public Member Functions inherited from G4eeToTwoGammaModel
 G4eeToTwoGammaModel (const G4ParticleDefinition *p=0, const G4String &nam="eplus2gg")
 
virtual ~G4eeToTwoGammaModel ()
 
virtual G4double ComputeCrossSectionPerAtom (const G4ParticleDefinition *, G4double kinEnergy, G4double Z, G4double A=0., G4double cutEnergy=0., G4double maxEnergy=DBL_MAX)
 
virtual G4double CrossSectionPerVolume (const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy, G4double maxEnergy)
 
- Public Member Functions inherited from G4VEmModel
 G4VEmModel (const G4String &nam)
 
virtual ~G4VEmModel ()
 
virtual void InitialiseLocal (const G4ParticleDefinition *, G4VEmModel *masterModel)
 
virtual void InitialiseForMaterial (const G4ParticleDefinition *, const G4Material *)
 
virtual void InitialiseForElement (const G4ParticleDefinition *, G4int Z)
 
virtual G4double ComputeDEDXPerVolume (const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=DBL_MAX)
 
virtual G4double ChargeSquareRatio (const G4Track &)
 
virtual G4double GetChargeSquareRatio (const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
 
virtual G4double GetParticleCharge (const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
 
virtual void StartTracking (G4Track *)
 
virtual void CorrectionsAlongStep (const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double &eloss, G4double &niel, G4double length)
 
virtual G4double Value (const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy)
 
virtual G4double MinPrimaryEnergy (const G4Material *, const G4ParticleDefinition *, G4double cut=0.0)
 
virtual G4double MinEnergyCut (const G4ParticleDefinition *, const G4MaterialCutsCouple *)
 
virtual void SetupForMaterial (const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
 
virtual void DefineForRegion (const G4Region *)
 
void InitialiseElementSelectors (const G4ParticleDefinition *, const G4DataVector &)
 
std::vector
< G4EmElementSelector * > * 
GetElementSelectors ()
 
void SetElementSelectors (std::vector< G4EmElementSelector * > *)
 
G4double ComputeDEDX (const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=DBL_MAX)
 
G4double CrossSection (const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
G4double ComputeMeanFreePath (const G4ParticleDefinition *, G4double kineticEnergy, const G4Material *, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
G4double ComputeCrossSectionPerAtom (const G4ParticleDefinition *, const G4Element *, G4double kinEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
G4int SelectIsotopeNumber (const G4Element *)
 
const G4ElementSelectRandomAtom (const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
const G4ElementSelectRandomAtom (const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
G4int SelectRandomAtomNumber (const G4Material *)
 
void SetParticleChange (G4VParticleChange *, G4VEmFluctuationModel *f=0)
 
void SetCrossSectionTable (G4PhysicsTable *, G4bool isLocal)
 
G4ElementDataGetElementData ()
 
G4PhysicsTableGetCrossSectionTable ()
 
G4VEmFluctuationModelGetModelOfFluctuations ()
 
G4VEmAngularDistributionGetAngularDistribution ()
 
void SetAngularDistribution (G4VEmAngularDistribution *)
 
G4double HighEnergyLimit () const
 
G4double LowEnergyLimit () const
 
G4double HighEnergyActivationLimit () const
 
G4double LowEnergyActivationLimit () const
 
G4double PolarAngleLimit () const
 
G4double SecondaryThreshold () const
 
G4bool LPMFlag () const
 
G4bool DeexcitationFlag () const
 
G4bool ForceBuildTableFlag () const
 
G4bool UseAngularGeneratorFlag () const
 
void SetAngularGeneratorFlag (G4bool)
 
void SetHighEnergyLimit (G4double)
 
void SetLowEnergyLimit (G4double)
 
void SetActivationHighEnergyLimit (G4double)
 
void SetActivationLowEnergyLimit (G4double)
 
G4bool IsActive (G4double kinEnergy)
 
void SetPolarAngleLimit (G4double)
 
void SetSecondaryThreshold (G4double)
 
void SetLPMFlag (G4bool val)
 
void SetDeexcitationFlag (G4bool val)
 
void SetForceBuildTable (G4bool val)
 
void SetMasterThread (G4bool val)
 
G4bool IsMaster () const
 
G4double MaxSecondaryKinEnergy (const G4DynamicParticle *dynParticle)
 
const G4StringGetName () const
 
void SetCurrentCouple (const G4MaterialCutsCouple *)
 
const G4ElementGetCurrentElement () const
 

Additional Inherited Members

- Protected Member Functions inherited from G4VEmModel
G4ParticleChangeForLossGetParticleChangeForLoss ()
 
G4ParticleChangeForGammaGetParticleChangeForGamma ()
 
virtual G4double MaxSecondaryEnergy (const G4ParticleDefinition *, G4double kineticEnergy)
 
const G4MaterialCutsCoupleCurrentCouple () const
 
void SetCurrentElement (const G4Element *)
 
- Protected Attributes inherited from G4VEmModel
G4ElementDatafElementData
 
G4VParticleChangepParticleChange
 
G4PhysicsTablexSectionTable
 
const std::vector< G4double > * theDensityFactor
 
const std::vector< G4int > * theDensityIdx
 
size_t idxTable
 

Detailed Description

Definition at line 63 of file G4PolarizedAnnihilationModel.hh.

Constructor & Destructor Documentation

G4PolarizedAnnihilationModel::G4PolarizedAnnihilationModel ( const G4ParticleDefinition p = 0,
const G4String nam = "Polarized-Annihilation" 
)

Definition at line 64 of file G4PolarizedAnnihilationModel.cc.

66  : G4eeToTwoGammaModel(p,nam),crossSectionCalculator(0),verboseLevel(0),gParticleChange(0),
67  gIsInitialised(false)
68 {
69  crossSectionCalculator=new G4PolarizedAnnihilationCrossSection();
70 }
G4eeToTwoGammaModel(const G4ParticleDefinition *p=0, const G4String &nam="eplus2gg")
G4PolarizedAnnihilationModel::~G4PolarizedAnnihilationModel ( )
virtual

Definition at line 72 of file G4PolarizedAnnihilationModel.cc.

73 {
74  if (crossSectionCalculator) delete crossSectionCalculator;
75 }

Member Function Documentation

void G4PolarizedAnnihilationModel::ComputeAsymmetriesPerElectron ( G4double  gammaEnergy,
G4double valueX,
G4double valueA,
G4double valueT 
)

Definition at line 110 of file G4PolarizedAnnihilationModel.cc.

References python.hepunit::electron_mass_c2, G4cout, G4endl, G4InuclParticleNames::gam, G4StokesVector::P1, G4StokesVector::P2, G4StokesVector::P3, G4PolarizedAnnihilationCrossSection::TotalXSection(), and G4StokesVector::ZERO.

Referenced by ComputeCrossSectionPerElectron().

114 {
115  // *** calculate asymmetries
116  G4double gam = 1. + ene/electron_mass_c2;
117  G4double xs0=crossSectionCalculator->TotalXSection(0.,1.,gam,
120  G4double xsA=crossSectionCalculator->TotalXSection(0.,1.,gam,
123  G4double xsT1=crossSectionCalculator->TotalXSection(0.,1.,gam,
126  G4double xsT2=crossSectionCalculator->TotalXSection(0.,1.,gam,
129  G4double xsT=0.5*(xsT1+xsT2);
130 
131  valueX=xs0;
132  valueA=xsA/xs0-1.;
133  valueT=xsT/xs0-1.;
134  // G4cout<<valueX<<"\t"<<valueA<<"\t"<<valueT<<" energy = "<<gam<<G4endl;
135  if ( (valueA < -1) || (1 < valueA)) {
136  G4cout<< " ERROR PolarizedAnnihilationPS::ComputeAsymmetries \n";
137  G4cout<< " something wrong in total cross section calculation (valueA)\n";
138  G4cout<<"*********** LONG "<<valueX<<"\t"<<valueA<<"\t"<<valueT<<" energy = "<<gam<<G4endl;
139  }
140  if ( (valueT < -1) || (1 < valueT)) {
141  G4cout<< " ERROR PolarizedAnnihilationPS::ComputeAsymmetries \n";
142  G4cout<< " something wrong in total cross section calculation (valueT)\n";
143  G4cout<<"****** TRAN "<<valueX<<"\t"<<valueA<<"\t"<<valueT<<" energy = "<<gam<<G4endl;
144  }
145 }
virtual G4double TotalXSection(G4double xmin, G4double xmax, G4double y, const G4StokesVector &pol0, const G4StokesVector &pol1)
static const G4StokesVector P3
G4GLOB_DLL std::ostream G4cout
static const G4StokesVector P2
static const G4StokesVector P1
float electron_mass_c2
Definition: hepunit.py:274
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76
static const G4StokesVector ZERO
G4double G4PolarizedAnnihilationModel::ComputeCrossSectionPerElectron ( const G4ParticleDefinition pd,
G4double  kinEnergy,
G4double  cut,
G4double  emax 
)
virtual

Reimplemented from G4eeToTwoGammaModel.

Definition at line 89 of file G4PolarizedAnnihilationModel.cc.

References ComputeAsymmetriesPerElectron(), G4eeToTwoGammaModel::ComputeCrossSectionPerElectron(), CLHEP::Hep3Vector::x(), CLHEP::Hep3Vector::y(), and CLHEP::Hep3Vector::z().

94 {
96  cut,emax);
97 
98  G4double polzz = theBeamPolarization.z()*theTargetPolarization.z();
99  G4double poltt = theBeamPolarization.x()*theTargetPolarization.x()
100  + theBeamPolarization.y()*theTargetPolarization.y();
101  if (polzz!=0 || poltt!=0) {
102  G4double xval,lasym,tasym;
103  ComputeAsymmetriesPerElectron(kinEnergy,xval,lasym,tasym);
104  xs*=(1.+polzz*lasym+poltt*tasym);
105  }
106 
107  return xs;
108 }
double x() const
double z() const
virtual G4double ComputeCrossSectionPerElectron(const G4ParticleDefinition *, G4double kinEnergy, G4double cutEnergy=0., G4double maxEnergy=DBL_MAX)
void ComputeAsymmetriesPerElectron(G4double gammaEnergy, G4double &valueX, G4double &valueA, G4double &valueT)
double y() const
double G4double
Definition: G4Types.hh:76
const G4ThreeVector & G4PolarizedAnnihilationModel::GetBeamPolarization ( ) const
inline

Definition at line 132 of file G4PolarizedAnnihilationModel.hh.

133 {
134  return theBeamPolarization;
135 }
const G4ThreeVector & G4PolarizedAnnihilationModel::GetFinalGamma1Polarization ( ) const
inline

Definition at line 136 of file G4PolarizedAnnihilationModel.hh.

137 {
138  return finalGamma1Polarization;
139 }
const G4ThreeVector & G4PolarizedAnnihilationModel::GetFinalGamma2Polarization ( ) const
inline

Definition at line 140 of file G4PolarizedAnnihilationModel.hh.

141 {
142  return finalGamma2Polarization;
143 }
const G4ThreeVector & G4PolarizedAnnihilationModel::GetTargetPolarization ( ) const
inline

Definition at line 128 of file G4PolarizedAnnihilationModel.hh.

129 {
130  return theTargetPolarization;
131 }
void G4PolarizedAnnihilationModel::Initialise ( const G4ParticleDefinition ,
const G4DataVector  
)
virtual

Reimplemented from G4eeToTwoGammaModel.

Definition at line 80 of file G4PolarizedAnnihilationModel.cc.

References G4VEmModel::GetParticleChangeForGamma().

82 {
83  // G4eeToTwoGammaModel::Initialise(part,dv);
84  if(gIsInitialised) return;
85  gParticleChange = GetParticleChangeForGamma();
86  gIsInitialised = true;
87 }
G4ParticleChangeForGamma * GetParticleChangeForGamma()
Definition: G4VEmModel.cc:121
void G4PolarizedAnnihilationModel::SampleSecondaries ( std::vector< G4DynamicParticle * > *  fvect,
const G4MaterialCutsCouple ,
const G4DynamicParticle dp,
G4double  tmin,
G4double  maxEnergy 
)
virtual

Reimplemented from G4eeToTwoGammaModel.

Definition at line 148 of file G4PolarizedAnnihilationModel.cc.

References DBL_MIN, G4PolarizedAnnihilationCrossSection::DiceEpsilon(), python.hepunit::electron_mass_c2, fStopAndKill, G4cout, G4endl, G4UniformRand, G4Gamma::Gamma(), G4ParticleChangeForGamma::GetCurrentTrack(), G4PolarizationHelper::GetFrame(), G4PolarizationManager::GetInstance(), G4DynamicParticle::GetKineticEnergy(), G4VPhysicalVolume::GetLogicalVolume(), G4DynamicParticle::GetMomentumDirection(), G4PolarizedAnnihilationCrossSection::GetPol2(), G4PolarizedAnnihilationCrossSection::GetPol3(), G4Track::GetPolarization(), G4PolarizedAnnihilationCrossSection::getVar(), G4Track::GetVolume(), G4PolarizationManager::GetVolumePolarization(), G4PolarizedAnnihilationCrossSection::Initialize(), G4StokesVector::InvRotateAz(), G4PolarizationManager::IsPolarized(), CLHEP::Hep3Vector::mag2(), G4StokesVector::p1(), G4StokesVector::p2(), G4StokesVector::p3(), G4VParticleChange::ProposeTrackStatus(), G4StokesVector::RotateAz(), CLHEP::Hep3Vector::rotateUz(), G4StokesVector::SetPhoton(), G4DynamicParticle::SetPolarization(), G4ParticleChangeForGamma::SetProposedKineticEnergy(), sqr(), and python.hepunit::twopi.

153 {
154 // G4ParticleChangeForGamma* gParticleChange
155 // = dynamic_cast<G4ParticleChangeForGamma*>(pParticleChange);
156  const G4Track * aTrack = gParticleChange->GetCurrentTrack();
157 
158  // kill primary
159  gParticleChange->SetProposedKineticEnergy(0.);
160  gParticleChange->ProposeTrackStatus(fStopAndKill);
161 
162  // V.Ivanchenko add protection against zero kin energy
163  G4double PositKinEnergy = dp->GetKineticEnergy();
164 
165  if(PositKinEnergy < DBL_MIN) {
166 
167  G4double cosTeta = 2.*G4UniformRand()-1.;
168  G4double sinTeta = std::sqrt((1.0 - cosTeta)*(1.0 + cosTeta));
169  G4double phi = twopi * G4UniformRand();
170  G4ThreeVector dir(sinTeta*std::cos(phi), sinTeta*std::sin(phi), cosTeta);
171  fvect->push_back( new G4DynamicParticle(G4Gamma::Gamma(), dir, electron_mass_c2));
172  fvect->push_back( new G4DynamicParticle(G4Gamma::Gamma(),-dir, electron_mass_c2));
173  return;
174  }
175 
176  // *** obtain and save target and beam polarization ***
178 
179  // obtain polarization of the beam
180  theBeamPolarization = aTrack->GetPolarization();
181 
182  // obtain polarization of the media
183  G4VPhysicalVolume* aPVolume = aTrack->GetVolume();
184  G4LogicalVolume* aLVolume = aPVolume->GetLogicalVolume();
185  const G4bool targetIsPolarized = polarizationManager->IsPolarized(aLVolume);
186  theTargetPolarization = polarizationManager->GetVolumePolarization(aLVolume);
187 
188  // transfer target electron polarization in frame of positron
189  if (targetIsPolarized)
190  theTargetPolarization.rotateUz(dp->GetMomentumDirection());
191 
192  G4ParticleMomentum PositDirection = dp->GetMomentumDirection();
193 
194  // polar asymmetry:
195  G4double polarization = theBeamPolarization.p3()*theTargetPolarization.p3();
196 
197  G4double gamam1 = PositKinEnergy/electron_mass_c2;
198  G4double gama = gamam1+1. , gamap1 = gamam1+2.;
199  G4double sqgrate = std::sqrt(gamam1/gamap1)/2. , sqg2m1 = std::sqrt(gamam1*gamap1);
200 
201  // limits of the energy sampling
202  G4double epsilmin = 0.5 - sqgrate , epsilmax = 0.5 + sqgrate;
203  G4double epsilqot = epsilmax/epsilmin;
204 
205  //
206  // sample the energy rate of the created gammas
207  // note: for polarized partices, the actual dicing strategy
208  // will depend on the energy, and the degree of polarization !!
209  //
210  G4double epsil;
211  G4double gmax=1. + std::fabs(polarization); // crude estimate
212 
213  //G4bool check_range=true;
214 
215  crossSectionCalculator->Initialize(epsilmin, gama, 0., theBeamPolarization, theTargetPolarization);
216  if (crossSectionCalculator->DiceEpsilon()<0) {
217  G4cout<<"ERROR in PolarizedAnnihilationPS::PostStepDoIt\n"
218  <<"epsilmin DiceRoutine not appropriate ! "<<crossSectionCalculator->DiceEpsilon()<<G4endl;
219  //check_range=false;
220  }
221 
222  crossSectionCalculator->Initialize(epsilmax, gama, 0., theBeamPolarization, theTargetPolarization);
223  if (crossSectionCalculator->DiceEpsilon()<0) {
224  G4cout<<"ERROR in PolarizedAnnihilationPS::PostStepDoIt\n"
225  <<"epsilmax DiceRoutine not appropriate ! "<<crossSectionCalculator->DiceEpsilon()<<G4endl;
226  //check_range=false;
227  }
228 
229  G4int ncount=0;
230  G4double trejectmax=0.;
231  G4double treject;
232 
233 
234  do {
235  //
236  epsil = epsilmin*std::pow(epsilqot,G4UniformRand());
237 
238  crossSectionCalculator->Initialize(epsil, gama, 0., theBeamPolarization, theTargetPolarization,1);
239 
240  treject = crossSectionCalculator->DiceEpsilon();
241  treject*=epsil;
242 
243  if (treject>gmax || treject<0.)
244  G4cout<<"ERROR in PolarizedAnnihilationPS::PostStepDoIt\n"
245  <<" eps ("<<epsil<<") rejection does not work properly: "<<treject<<G4endl;
246  ++ncount;
247  if (treject>trejectmax) trejectmax=treject;
248  if (ncount>1000) {
249  G4cout<<"WARNING in PolarizedAnnihilationPS::PostStepDoIt\n"
250  <<"eps dicing very inefficient ="<<trejectmax/gmax
251  <<", "<<treject/gmax<<". For secondary energy = "<<epsil<<" "<<ncount<<G4endl;
252  break;
253  }
254 
255  } while( treject < gmax*G4UniformRand() );
256 
257  //
258  // scattered Gamma angles. ( Z - axis along the parent positron)
259  //
260 
261  G4double cost = (epsil*gamap1-1.)/(epsil*sqg2m1);
262  G4double sint = std::sqrt((1.+cost)*(1.-cost));
263  G4double phi = 0.;
264  G4double beamTrans = std::sqrt(sqr(theBeamPolarization.p1()) + sqr(theBeamPolarization.p2()));
265  G4double targetTrans = std::sqrt(sqr(theTargetPolarization.p1()) + sqr(theTargetPolarization.p2()));
266 
267  // G4cout<<"phi dicing START"<<G4endl;
268  do{
269  phi = twopi * G4UniformRand();
270  crossSectionCalculator->Initialize(epsil, gama, 0., theBeamPolarization, theTargetPolarization,2);
271 
272  G4double gdiced =crossSectionCalculator->getVar(0);
273  gdiced += crossSectionCalculator->getVar(3)*theBeamPolarization.p3()*theTargetPolarization.p3();
274  gdiced += 1.*(std::fabs(crossSectionCalculator->getVar(1))
275  + std::fabs(crossSectionCalculator->getVar(2)))*beamTrans*targetTrans;
276  gdiced += 1.*std::fabs(crossSectionCalculator->getVar(4))
277  *(std::fabs(theBeamPolarization.p3())*targetTrans + std::fabs(theTargetPolarization.p3())*beamTrans);
278 
279  G4double gdist = crossSectionCalculator->getVar(0);
280  gdist += crossSectionCalculator->getVar(3)*theBeamPolarization.p3()*theTargetPolarization.p3();
281  gdist += crossSectionCalculator->getVar(1)*(std::cos(phi)*theBeamPolarization.p1()
282  + std::sin(phi)*theBeamPolarization.p2())
283  *(std::cos(phi)*theTargetPolarization.p1()
284  + std::sin(phi)*theTargetPolarization.p2());
285  gdist += crossSectionCalculator->getVar(2)*(std::cos(phi)*theBeamPolarization.p2()
286  - std::sin(phi)*theBeamPolarization.p1())
287  *(std::cos(phi)*theTargetPolarization.p2()
288  - std::sin(phi)*theTargetPolarization.p1());
289  gdist += crossSectionCalculator->getVar(4)
290  *(std::cos(phi)*theBeamPolarization.p3()*theTargetPolarization.p1()
291  + std::cos(phi)*theBeamPolarization.p1()*theTargetPolarization.p3()
292  + std::sin(phi)*theBeamPolarization.p3()*theTargetPolarization.p2()
293  + std::sin(phi)*theBeamPolarization.p2()*theTargetPolarization.p3());
294 
295  treject = gdist/gdiced;
296  //G4cout<<" treject = "<<treject<<" at phi = "<<phi<<G4endl;
297  if (treject>1.+1.e-10 || treject<0){
298  G4cout<<"!!!ERROR in PolarizedAnnihilationPS::PostStepDoIt\n"
299  <<" phi rejection does not work properly: "<<treject<<G4endl;
300  G4cout<<" gdiced = "<<gdiced<<G4endl;
301  G4cout<<" gdist = "<<gdist<<G4endl;
302  G4cout<<" epsil = "<<epsil<<G4endl;
303  }
304 
305  if (treject<1.e-3) {
306  G4cout<<"!!!ERROR in PolarizedAnnihilationPS::PostStepDoIt\n"
307  <<" phi rejection does not work properly: "<<treject<<"\n";
308  G4cout<<" gdiced="<<gdiced<<" gdist="<<gdist<<"\n";
309  G4cout<<" epsil = "<<epsil<<G4endl;
310  }
311 
312  } while( treject < G4UniformRand() );
313  // G4cout<<"phi dicing END"<<G4endl;
314 
315  G4double dirx = sint*std::cos(phi) , diry = sint*std::sin(phi) , dirz = cost;
316 
317  //
318  // kinematic of the created pair
319  //
320  G4double TotalAvailableEnergy = PositKinEnergy + 2*electron_mass_c2;
321  G4double Phot1Energy = epsil*TotalAvailableEnergy;
322  G4double Phot2Energy =(1.-epsil)*TotalAvailableEnergy;
323 
324  // *** prepare calculation of polarization transfer ***
325  G4ThreeVector Phot1Direction (dirx, diry, dirz);
326 
327  // get interaction frame
328  G4ThreeVector nInteractionFrame =
329  G4PolarizationHelper::GetFrame(PositDirection,Phot1Direction);
330 
331  // define proper in-plane and out-of-plane component of initial spins
332  theBeamPolarization.InvRotateAz(nInteractionFrame,PositDirection);
333  theTargetPolarization.InvRotateAz(nInteractionFrame,PositDirection);
334 
335  // calculate spin transfere matrix
336 
337  crossSectionCalculator->Initialize(epsil,gama,phi,theBeamPolarization,theTargetPolarization,2);
338 
339  // **********************************************************************
340 
341  Phot1Direction.rotateUz(PositDirection);
342  // create G4DynamicParticle object for the particle1
344  Phot1Direction, Phot1Energy);
345  finalGamma1Polarization=crossSectionCalculator->GetPol2();
346  G4double n1=finalGamma1Polarization.mag2();
347  if (n1>1) {
348  G4cout<<"ERROR: PolarizedAnnihilation Polarization Vector at epsil = "
349  <<epsil<<" is too large!!! \n"
350  <<"annihi pol1= "<<finalGamma1Polarization<<", ("<<n1<<")\n";
351  finalGamma1Polarization*=1./std::sqrt(n1);
352  }
353 
354  // define polarization of first final state photon
355  finalGamma1Polarization.SetPhoton();
356  finalGamma1Polarization.RotateAz(nInteractionFrame,Phot1Direction);
357  aParticle1->SetPolarization(finalGamma1Polarization.p1(),
358  finalGamma1Polarization.p2(),
359  finalGamma1Polarization.p3());
360 
361  fvect->push_back(aParticle1);
362 
363 
364  // **********************************************************************
365 
366  G4double Eratio= Phot1Energy/Phot2Energy;
367  G4double PositP= std::sqrt(PositKinEnergy*(PositKinEnergy+2.*electron_mass_c2));
368  G4ThreeVector Phot2Direction (-dirx*Eratio, -diry*Eratio,
369  (PositP-dirz*Phot1Energy)/Phot2Energy);
370  Phot2Direction.rotateUz(PositDirection);
371  // create G4DynamicParticle object for the particle2
373  Phot2Direction, Phot2Energy);
374 
375  // define polarization of second final state photon
376  finalGamma2Polarization=crossSectionCalculator->GetPol3();
377  G4double n2=finalGamma2Polarization.mag2();
378  if (n2>1) {
379  G4cout<<"ERROR: PolarizedAnnihilation Polarization Vector at epsil = "<<epsil<<" is too large!!! \n";
380  G4cout<<"annihi pol2= "<<finalGamma2Polarization<<", ("<<n2<<")\n";
381 
382  finalGamma2Polarization*=1./std::sqrt(n2);
383  }
384  finalGamma2Polarization.SetPhoton();
385  finalGamma2Polarization.RotateAz(nInteractionFrame,Phot2Direction);
386  aParticle2->SetPolarization(finalGamma2Polarization.p1(),
387  finalGamma2Polarization.p2(),
388  finalGamma2Polarization.p3());
389 
390  fvect->push_back(aParticle2);
391 }
const G4ThreeVector & GetPolarization() const
G4double GetKineticEnergy() const
G4double p2() const
static G4PolarizationManager * GetInstance()
int G4int
Definition: G4Types.hh:78
G4double p3() const
#define G4UniformRand()
Definition: Randomize.hh:87
G4GLOB_DLL std::ostream G4cout
bool G4bool
Definition: G4Types.hh:79
const G4ThreeVector & GetMomentumDirection() const
Hep3Vector & rotateUz(const Hep3Vector &)
Definition: ThreeVector.cc:72
static G4ThreeVector GetFrame(const G4ThreeVector &, const G4ThreeVector &)
G4double p1() const
float electron_mass_c2
Definition: hepunit.py:274
void SetPolarization(G4double polX, G4double polY, G4double polZ)
static G4Gamma * Gamma()
Definition: G4Gamma.cc:86
virtual void Initialize(G4double eps, G4double gamma, G4double phi, const G4StokesVector &p0, const G4StokesVector &p1, G4int flag=0)
void InvRotateAz(G4ThreeVector nInteractionFrame, G4ThreeVector particleDirection)
G4LogicalVolume * GetLogicalVolume() const
bool IsPolarized(G4LogicalVolume *lVol) const
#define DBL_MIN
Definition: templates.hh:75
double mag2() const
G4VPhysicalVolume * GetVolume() const
const G4Track * GetCurrentTrack() const
void SetProposedKineticEnergy(G4double proposedKinEnergy)
#define G4endl
Definition: G4ios.hh:61
T sqr(const T &x)
Definition: templates.hh:145
double G4double
Definition: G4Types.hh:76
void ProposeTrackStatus(G4TrackStatus status)
const G4ThreeVector & GetVolumePolarization(G4LogicalVolume *lVol) const
void RotateAz(G4ThreeVector nInteractionFrame, G4ThreeVector particleDirection)
void G4PolarizedAnnihilationModel::SetBeamPolarization ( const G4ThreeVector pBeam)
inline

Definition at line 124 of file G4PolarizedAnnihilationModel.hh.

Referenced by G4eplusPolarizedAnnihilation::ComputeAsymmetry().

125 {
126  theBeamPolarization = pBeam;
127 }
void G4PolarizedAnnihilationModel::SetTargetPolarization ( const G4ThreeVector pTarget)
inline

Definition at line 120 of file G4PolarizedAnnihilationModel.hh.

Referenced by G4eplusPolarizedAnnihilation::ComputeAsymmetry().

121 {
122  theTargetPolarization = pTarget;
123 }

The documentation for this class was generated from the following files: