00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032
00033
00034
00035
00036
00037
00038
00039
00040
00041
00042
00043
00044
00045
00046
00047
00048
00049
00050
00051
00052 #ifndef G4PolarizedComptonModel_h
00053 #define G4PolarizedComptonModel_h 1
00054
00055 #include "G4KleinNishinaCompton.hh"
00056 #include "G4StokesVector.hh"
00057
00058 class G4ParticleChangeForGamma;
00059 class G4PolarizedComptonCrossSection;
00060
00061 class G4PolarizedComptonModel : public G4KleinNishinaCompton
00062 {
00063
00064 public:
00065
00066 G4PolarizedComptonModel(const G4ParticleDefinition* p = 0,
00067 const G4String& nam = "Polarized-Compton");
00068
00069 virtual G4double ComputeCrossSectionPerAtom(
00070 const G4ParticleDefinition*,
00071 G4double kinEnergy,
00072 G4double Z,
00073 G4double A,
00074 G4double cut,
00075 G4double emax);
00076 G4double ComputeAsymmetryPerAtom(G4double gammaEnergy, G4double Z);
00077 virtual ~G4PolarizedComptonModel();
00078
00079 virtual void SampleSecondaries(std::vector<G4DynamicParticle*>*,
00080 const G4MaterialCutsCouple*,
00081 const G4DynamicParticle*,
00082 G4double tmin,
00083 G4double maxEnergy);
00084
00085
00086 inline void SetTargetPolarization(const G4ThreeVector & pTarget);
00087 inline void SetBeamPolarization(const G4ThreeVector & pBeam);
00088 inline const G4ThreeVector & GetTargetPolarization() const;
00089 inline const G4ThreeVector & GetBeamPolarization() const;
00090 inline const G4ThreeVector & GetFinalGammaPolarization() const;
00091 inline const G4ThreeVector & GetFinalElectronPolarization() const;
00092 private:
00093
00094
00095 G4PolarizedComptonModel & operator=(const G4PolarizedComptonModel &right);
00096 G4PolarizedComptonModel(const G4PolarizedComptonModel&);
00097
00098 G4PolarizedComptonCrossSection * crossSectionCalculator;
00099
00100 G4StokesVector theBeamPolarization;
00101 G4StokesVector theTargetPolarization;
00102
00103 G4StokesVector finalGammaPolarization;
00104 G4StokesVector finalElectronPolarization;
00105
00106 G4int verboseLevel;
00107 };
00108
00109
00110
00111 inline void G4PolarizedComptonModel::SetTargetPolarization(const G4ThreeVector & pTarget)
00112 {
00113 theTargetPolarization = pTarget;
00114 }
00115 inline void G4PolarizedComptonModel::SetBeamPolarization(const G4ThreeVector & pBeam)
00116 {
00117 theBeamPolarization = pBeam;
00118 }
00119 inline const G4ThreeVector & G4PolarizedComptonModel::GetTargetPolarization() const
00120 {
00121 return theTargetPolarization;
00122 }
00123 inline const G4ThreeVector & G4PolarizedComptonModel::GetBeamPolarization() const
00124 {
00125 return theBeamPolarization;
00126 }
00127 inline const G4ThreeVector & G4PolarizedComptonModel::GetFinalGammaPolarization() const
00128 {
00129 return finalGammaPolarization;
00130 }
00131 inline const G4ThreeVector & G4PolarizedComptonModel::GetFinalElectronPolarization() const
00132 {
00133 return finalElectronPolarization;
00134 }
00135
00136 #endif