G4ePolarizedBremsstrahlungModel.cc

Go to the documentation of this file.
00001 //
00002 // ********************************************************************
00003 // * License and Disclaimer                                           *
00004 // *                                                                  *
00005 // * The  Geant4 software  is  copyright of the Copyright Holders  of *
00006 // * the Geant4 Collaboration.  It is provided  under  the terms  and *
00007 // * conditions of the Geant4 Software License,  included in the file *
00008 // * LICENSE and available at  http://cern.ch/geant4/license .  These *
00009 // * include a list of copyright holders.                             *
00010 // *                                                                  *
00011 // * Neither the authors of this software system, nor their employing *
00012 // * institutes,nor the agencies providing financial support for this *
00013 // * work  make  any representation or  warranty, express or implied, *
00014 // * regarding  this  software system or assume any liability for its *
00015 // * use.  Please see the license in the file  LICENSE  and URL above *
00016 // * for the full disclaimer and the limitation of liability.         *
00017 // *                                                                  *
00018 // * This  code  implementation is the result of  the  scientific and *
00019 // * technical work of the GEANT4 collaboration.                      *
00020 // * By using,  copying,  modifying or  distributing the software (or *
00021 // * any work based  on the software)  you  agree  to acknowledge its *
00022 // * use  in  resulting  scientific  publications,  and indicate your *
00023 // * acceptance of all terms of the Geant4 Software license.          *
00024 // ********************************************************************
00025 //
00026 // $Id: G4ePolarizedBremsstrahlungModel.cc 69847 2013-05-16 09:36:18Z gcosmo $
00027 //
00028 // -------------------------------------------------------------------
00029 //
00030 // GEANT4 Class file
00031 //
00032 //
00033 // File name:     G4ePolarizedBremsstrahlungModel
00034 //
00035 // Author:        Karim Laihem
00036 //
00037 // Creation date: 12.03.2005
00038 //
00039 // Modifications:
00040 //    19-08-06 addapted to accomodate geant481 structure
00041 //
00042 //
00043 // Class Description:
00044 //
00045 //
00046 // -------------------------------------------------------------------
00047 //
00048 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
00049 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
00050 
00051 #include "G4ePolarizedBremsstrahlungModel.hh"
00052 #include "G4PolarizedBremsstrahlungCrossSection.hh"
00053 #include "G4ParticleChangeForLoss.hh"
00054 #include "G4PolarizationHelper.hh"
00055 
00056 G4ePolarizedBremsstrahlungModel::G4ePolarizedBremsstrahlungModel(const G4ParticleDefinition* p,
00057                                                const G4String& nam)
00058   : G4eBremsstrahlungModel(p,nam),
00059     crossSectionCalculator(0)
00060 {
00061 }
00062 
00063 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
00064 
00065 G4ePolarizedBremsstrahlungModel::~G4ePolarizedBremsstrahlungModel()
00066 {
00067   if (crossSectionCalculator) delete crossSectionCalculator;
00068 }
00069   
00070 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
00071 
00072 void G4ePolarizedBremsstrahlungModel::Initialise(const G4ParticleDefinition* p, 
00073                                                          const G4DataVector& d)
00074 {
00075   G4eBremsstrahlungModel::Initialise(p,d);
00076   if (!crossSectionCalculator)
00077     crossSectionCalculator = new G4PolarizedBremsstrahlungCrossSection();
00078 }
00079 
00080 
00081 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
00082 
00083 
00084 void G4ePolarizedBremsstrahlungModel::SampleSecondaries(std::vector<G4DynamicParticle*>* vdp,
00085                                                         const G4MaterialCutsCouple* couple,
00086                                                         const G4DynamicParticle* dp,
00087                                                         G4double tmin,
00088                                                         G4double maxEnergy)
00089 {
00090   G4eBremsstrahlungModel::SampleSecondaries(vdp,couple,dp,tmin,maxEnergy);
00091   G4int num = vdp->size();
00092 
00093   if(num > 0) {
00094     G4double lepEnergy0 = dp->GetKineticEnergy();
00095     G4double gamEnergy1 = (*vdp)[0]->GetKineticEnergy();
00096     G4double sintheta   = dp->GetMomentumDirection().cross((*vdp)[0]->GetMomentumDirection()).mag();
00097     if (sintheta>1.) sintheta=1.;
00098 
00099 
00100     G4StokesVector beamPol = dp->GetPolarization();
00101 
00102     // determine interaction plane
00103     G4ThreeVector  nInteractionFrame = 
00104       G4PolarizationHelper::GetFrame(dp->GetMomentumDirection(), 
00105                  fParticleChange->GetProposedMomentumDirection());
00106 
00107     // transform polarization into interaction frame
00108      beamPol.InvRotateAz(nInteractionFrame,dp->GetMomentumDirection());
00109 
00110     // calulcate polarization transfer
00111     crossSectionCalculator->SetMaterial(GetCurrentElement()->GetN(),  // number of nucleons
00112                                         GetCurrentElement()->GetZ(),
00113                                         GetCurrentElement()->GetfCoulomb());
00114     crossSectionCalculator->Initialize(lepEnergy0, gamEnergy1, sintheta,
00115                                        beamPol, G4StokesVector::ZERO);
00116 
00117     // deterimine final state polarization
00118     G4StokesVector newBeamPol = crossSectionCalculator->GetPol2();
00119     newBeamPol.RotateAz(nInteractionFrame,
00120         fParticleChange->GetProposedMomentumDirection());
00121     fParticleChange->ProposePolarization(newBeamPol);
00122 
00123     if (num!=1) G4cout<<" WARNING "<<num<<" secondaries in polarized bremsstrahlung not supported!\n"; 
00124     for (G4int i=0; i<num; i++) {
00125       G4StokesVector photonPol = crossSectionCalculator->GetPol3();
00126       photonPol.SetPhoton();
00127       photonPol.RotateAz(nInteractionFrame,(*vdp)[i]->GetMomentumDirection());
00128       (*vdp)[i]->SetPolarization(photonPol.p1(),
00129                                  photonPol.p2(),
00130                                  photonPol.p3());
00131     }
00132   }
00133   return;
00134 }
00135 // The emitted gamma energy is sampled using a parametrized formula 

Generated on Mon May 27 17:48:10 2013 for Geant4 by  doxygen 1.4.7