G4ParticleChangeForLoss Class Reference

#include <G4ParticleChangeForLoss.hh>

Inheritance diagram for G4ParticleChangeForLoss:

G4VParticleChange

Public Member Functions

 G4ParticleChangeForLoss ()
virtual ~G4ParticleChangeForLoss ()
G4StepUpdateStepForAlongStep (G4Step *Step)
G4StepUpdateStepForPostStep (G4Step *Step)
void InitializeForAlongStep (const G4Track &)
void InitializeForPostStep (const G4Track &)
void AddSecondary (G4DynamicParticle *aParticle)
G4double GetProposedCharge () const
void SetProposedCharge (G4double theCharge)
G4double GetCharge () const
void ProposeCharge (G4double finalCharge)
G4double GetProposedKineticEnergy () const
void SetProposedKineticEnergy (G4double proposedKinEnergy)
const G4ThreeVectorGetProposedMomentumDirection () const
void SetProposedMomentumDirection (const G4ThreeVector &dir)
const G4ThreeVectorGetMomentumDirection () const
void ProposeMomentumDirection (G4double Px, G4double Py, G4double Pz)
void ProposeMomentumDirection (const G4ThreeVector &Pfinal)
const G4ThreeVectorGetProposedPolarization () const
void ProposePolarization (const G4ThreeVector &dir)
void ProposePolarization (G4double Px, G4double Py, G4double Pz)
const G4TrackGetCurrentTrack () const
void SetLowEnergyLimit (G4double elimit)
virtual void DumpInfo () const
virtual G4bool CheckIt (const G4Track &)

Protected Member Functions

 G4ParticleChangeForLoss (const G4ParticleChangeForLoss &right)
G4ParticleChangeForLossoperator= (const G4ParticleChangeForLoss &right)

Detailed Description

Definition at line 60 of file G4ParticleChangeForLoss.hh.


Constructor & Destructor Documentation

G4ParticleChangeForLoss::G4ParticleChangeForLoss (  ) 

Definition at line 52 of file G4ParticleChangeForLoss.cc.

References G4VParticleChange::debugFlag, G4cout, G4endl, NormalCondition, G4VParticleChange::theSteppingControlFlag, and G4VParticleChange::verboseLevel.

00053   : G4VParticleChange(), currentTrack(0), proposedKinEnergy(0.),
00054     lowEnergyLimit(1.0*eV), currentCharge(0.)
00055 {
00056   theSteppingControlFlag = NormalCondition;
00057   debugFlag = false;
00058 #ifdef G4VERBOSE
00059   if (verboseLevel>2) {
00060     G4cout << "G4ParticleChangeForLoss::G4ParticleChangeForLoss() " << G4endl;
00061   }
00062 #endif
00063 }

G4ParticleChangeForLoss::~G4ParticleChangeForLoss (  )  [virtual]

Definition at line 65 of file G4ParticleChangeForLoss.cc.

References G4cout, G4endl, and G4VParticleChange::verboseLevel.

00066 {
00067 #ifdef G4VERBOSE
00068   if (verboseLevel>2) {
00069     G4cout << "G4ParticleChangeForLoss::~G4ParticleChangeForLoss() " << G4endl;
00070   }
00071 #endif
00072 }

G4ParticleChangeForLoss::G4ParticleChangeForLoss ( const G4ParticleChangeForLoss right  )  [protected]

Definition at line 75 of file G4ParticleChangeForLoss.cc.

References currentCharge, currentTrack, G4cout, G4endl, lowEnergyLimit, proposedKinEnergy, proposedMomentumDirection, and G4VParticleChange::verboseLevel.

00076   : G4VParticleChange(right)
00077 {
00078   if (verboseLevel>1) {
00079     G4cout << "G4ParticleChangeForLoss::  copy constructor is called " << G4endl;
00080   }
00081   currentTrack = right.currentTrack;
00082   proposedKinEnergy = right.proposedKinEnergy;
00083   lowEnergyLimit = right.lowEnergyLimit;
00084   currentCharge = right.currentCharge;
00085   proposedMomentumDirection = right.proposedMomentumDirection;
00086 }


Member Function Documentation

void G4ParticleChangeForLoss::AddSecondary ( G4DynamicParticle aParticle  )  [inline]

Definition at line 260 of file G4ParticleChangeForLoss.hh.

References G4VParticleChange::AddSecondary(), G4Track::GetGlobalTime(), G4Track::GetPosition(), G4Track::GetTouchableHandle(), and G4Track::SetTouchableHandle().

00261 {
00262   //  create track
00263   G4Track* aTrack = new G4Track(aParticle, currentTrack->GetGlobalTime(),
00264                                            currentTrack->GetPosition());
00265 
00266   //   Touchable handle is copied to keep the pointer
00267   aTrack->SetTouchableHandle(currentTrack->GetTouchableHandle());
00268 
00269   //  add a secondary
00270   G4VParticleChange::AddSecondary(aTrack);
00271 }

G4bool G4ParticleChangeForLoss::CheckIt ( const G4Track  )  [virtual]

Reimplemented from G4VParticleChange.

Definition at line 160 of file G4ParticleChangeForLoss.cc.

References G4VParticleChange::accuracyForException, G4VParticleChange::accuracyForWarning, G4VParticleChange::CheckIt(), DumpInfo(), EventMustBeAborted, G4cout, G4endl, G4Exception(), G4Track::GetDefinition(), G4Track::GetKineticEnergy(), G4ParticleDefinition::GetParticleName(), and G4Track::GetPosition().

00161 {
00162   G4bool    itsOK = true;
00163   G4bool    exitWithError = false;
00164 
00165   G4double  accuracy;
00166 
00167   // Energy should not be lager than initial value
00168   accuracy = ( proposedKinEnergy - aTrack.GetKineticEnergy())/MeV;
00169   if (accuracy > accuracyForWarning) {
00170     itsOK = false;
00171     exitWithError = (accuracy > accuracyForException);
00172 #ifdef G4VERBOSE
00173     G4cout << "G4ParticleChangeForLoss::CheckIt: ";
00174     G4cout << "KinEnergy become larger than the initial value!" 
00175            << "  Difference:  " << accuracy  << "[MeV] " <<G4endl;
00176     G4cout << aTrack.GetDefinition()->GetParticleName()
00177            << " E=" << aTrack.GetKineticEnergy()/MeV
00178            << " pos=" << aTrack.GetPosition().x()/m
00179            << ", " << aTrack.GetPosition().y()/m
00180            << ", " << aTrack.GetPosition().z()/m
00181            <<G4endl;
00182 #endif
00183   }
00184 
00185   // dump out information of this particle change
00186 #ifdef G4VERBOSE
00187   if (!itsOK) DumpInfo();
00188 #endif
00189 
00190   // Exit with error
00191   if (exitWithError) {
00192     G4Exception("G4ParticleChangeForLoss::CheckIt",
00193                 "TRACK004", EventMustBeAborted,
00194                 "energy was  illegal");
00195   }
00196 
00197   //correction
00198   if (!itsOK) {
00199     proposedKinEnergy = aTrack.GetKineticEnergy();
00200   }
00201 
00202   itsOK = (itsOK) && G4VParticleChange::CheckIt(aTrack);
00203   return itsOK;
00204 }

void G4ParticleChangeForLoss::DumpInfo (  )  const [virtual]

Reimplemented from G4VParticleChange.

Definition at line 136 of file G4ParticleChangeForLoss.cc.

References G4VParticleChange::DumpInfo(), G4cout, and G4endl.

Referenced by CheckIt().

00137 {
00138 // use base-class DumpInfo
00139   G4VParticleChange::DumpInfo();
00140 
00141   G4int oldprc = G4cout.precision(3);
00142   G4cout << "        Charge (eplus)   : "
00143        << std::setw(20) << currentCharge/eplus
00144        << G4endl;
00145   G4cout << "        Kinetic Energy (MeV): "
00146        << std::setw(20) << proposedKinEnergy/MeV
00147        << G4endl;
00148   G4cout << "        Momentum Direct - x : "
00149        << std::setw(20) << proposedMomentumDirection.x()
00150        << G4endl;
00151   G4cout << "        Momentum Direct - y : "
00152        << std::setw(20) << proposedMomentumDirection.y()
00153        << G4endl;
00154   G4cout << "        Momentum Direct - z : "
00155        << std::setw(20) << proposedMomentumDirection.z()
00156        << G4endl;
00157   G4cout.precision(oldprc);
00158 }

G4double G4ParticleChangeForLoss::GetCharge (  )  const [inline]

Definition at line 160 of file G4ParticleChangeForLoss.hh.

00161 {
00162   return currentCharge;
00163 }

const G4Track * G4ParticleChangeForLoss::GetCurrentTrack (  )  const [inline]

Definition at line 207 of file G4ParticleChangeForLoss.hh.

Referenced by G4PolarizedMollerBhabhaModel::SampleSecondaries().

00208 {
00209   return currentTrack;
00210 }

const G4ThreeVector & G4ParticleChangeForLoss::GetMomentumDirection (  )  const [inline]

Definition at line 182 of file G4ParticleChangeForLoss.hh.

00183 {
00184   return proposedMomentumDirection;
00185 }

G4double G4ParticleChangeForLoss::GetProposedCharge (  )  const [inline]

Definition at line 155 of file G4ParticleChangeForLoss.hh.

00156 {
00157   return currentCharge;
00158 }

G4double G4ParticleChangeForLoss::GetProposedKineticEnergy (  )  const [inline]

Definition at line 145 of file G4ParticleChangeForLoss.hh.

Referenced by G4EmBiasingManager::ApplySecondaryBiasing(), and G4VEnergyLossProcess::PostStepDoIt().

00146 {
00147   return proposedKinEnergy;
00148 }

const G4ThreeVector & G4ParticleChangeForLoss::GetProposedMomentumDirection (  )  const [inline]

Definition at line 176 of file G4ParticleChangeForLoss.hh.

Referenced by G4EmBiasingManager::ApplySecondaryBiasing(), and G4ePolarizedBremsstrahlungModel::SampleSecondaries().

00177 {
00178   return proposedMomentumDirection;
00179 }

const G4ThreeVector & G4ParticleChangeForLoss::GetProposedPolarization (  )  const [inline]

Definition at line 213 of file G4ParticleChangeForLoss.hh.

00214 {
00215   return proposedPolarization;
00216 }

void G4ParticleChangeForLoss::InitializeForAlongStep ( const G4Track  )  [inline]

Definition at line 232 of file G4ParticleChangeForLoss.hh.

References G4DynamicParticle::GetCharge(), G4Track::GetDynamicParticle(), G4Track::GetKineticEnergy(), G4Track::GetTrackStatus(), G4Track::GetWeight(), G4VParticleChange::InitializeSecondaries(), G4VParticleChange::isParentWeightProposed, G4VParticleChange::theLocalEnergyDeposit, G4VParticleChange::theNonIonizingEnergyDeposit, G4VParticleChange::theParentWeight, and G4VParticleChange::theStatusChange.

Referenced by G4VEnergyLossProcess::AlongStepDoIt(), and G4NuclearStopping::AlongStepDoIt().

00233 {
00234   theStatusChange = track.GetTrackStatus();
00235   theLocalEnergyDeposit = 0.0;
00236   theNonIonizingEnergyDeposit = 0.0;
00237   InitializeSecondaries(track);
00238   theParentWeight = track.GetWeight();
00239   isParentWeightProposed = false;
00240   proposedKinEnergy = track.GetKineticEnergy();
00241   currentCharge = track.GetDynamicParticle()->GetCharge();
00242 }

void G4ParticleChangeForLoss::InitializeForPostStep ( const G4Track  )  [inline]

Definition at line 244 of file G4ParticleChangeForLoss.hh.

References G4DynamicParticle::GetCharge(), G4Track::GetDynamicParticle(), G4Track::GetKineticEnergy(), G4Track::GetMomentumDirection(), G4Track::GetPolarization(), G4Track::GetTrackStatus(), G4Track::GetWeight(), G4VParticleChange::InitializeSecondaries(), G4VParticleChange::isParentWeightProposed, G4VParticleChange::theLocalEnergyDeposit, G4VParticleChange::theNonIonizingEnergyDeposit, G4VParticleChange::theParentWeight, and G4VParticleChange::theStatusChange.

Referenced by G4VEnergyLossProcess::PostStepDoIt().

00245 {
00246   theStatusChange = track.GetTrackStatus();
00247   theLocalEnergyDeposit = 0.0;
00248   theNonIonizingEnergyDeposit = 0.0;
00249   InitializeSecondaries(track);
00250   theParentWeight = track.GetWeight();
00251   isParentWeightProposed = false;
00252   proposedKinEnergy = track.GetKineticEnergy();
00253   currentCharge = track.GetDynamicParticle()->GetCharge();
00254   proposedMomentumDirection = track.GetMomentumDirection();
00255   proposedPolarization = track.GetPolarization();
00256   currentTrack = &track;
00257 }

G4ParticleChangeForLoss & G4ParticleChangeForLoss::operator= ( const G4ParticleChangeForLoss right  )  [protected]

Definition at line 89 of file G4ParticleChangeForLoss.cc.

References currentCharge, currentTrack, G4VParticleChange::fSetSecondaryWeightByProcess, G4cout, G4endl, G4VParticleChange::isParentWeightProposed, proposedKinEnergy, proposedMomentumDirection, G4FastVector< Type, N >::SetElement(), G4VParticleChange::theListOfSecondaries, G4VParticleChange::theLocalEnergyDeposit, G4VParticleChange::theNumberOfSecondaries, G4VParticleChange::theParentWeight, G4VParticleChange::theStatusChange, G4VParticleChange::theSteppingControlFlag, and G4VParticleChange::verboseLevel.

00091 {
00092 #ifdef G4VERBOSE
00093   if (verboseLevel>1) {
00094     G4cout << "G4ParticleChangeForLoss:: assignment operator is called " << G4endl;
00095   }
00096 #endif
00097 
00098   if (this != &right) {
00099     if (theNumberOfSecondaries>0) {
00100 #ifdef G4VERBOSE
00101       if (verboseLevel>0) {
00102         G4cout << "G4ParticleChangeForLoss: assignment operator Warning  ";
00103         G4cout << "theListOfSecondaries is not empty ";
00104       }
00105 #endif
00106        for (G4int index= 0; index<theNumberOfSecondaries; index++){
00107          if ( (*theListOfSecondaries)[index] ) delete (*theListOfSecondaries)[index] ;
00108        }
00109     }
00110     delete theListOfSecondaries; 
00111     theListOfSecondaries =  new G4TrackFastVector();
00112     theNumberOfSecondaries = right.theNumberOfSecondaries;
00113     for (G4int index = 0; index<theNumberOfSecondaries; index++){
00114       G4Track* newTrack =  new G4Track(*((*right.theListOfSecondaries)[index] ));
00115       theListOfSecondaries->SetElement(index, newTrack);                            }
00116 
00117     theStatusChange = right.theStatusChange;
00118     theLocalEnergyDeposit = right.theLocalEnergyDeposit;
00119     theSteppingControlFlag = right.theSteppingControlFlag;
00120     theParentWeight = right.theParentWeight;
00121     isParentWeightProposed = right.isParentWeightProposed;
00122     fSetSecondaryWeightByProcess = right.fSetSecondaryWeightByProcess;
00123 
00124     currentTrack = right.currentTrack;
00125     proposedKinEnergy = right.proposedKinEnergy;
00126     currentCharge = right.currentCharge;
00127     proposedMomentumDirection = right.proposedMomentumDirection;
00128   }
00129   return *this;
00130 }

void G4ParticleChangeForLoss::ProposeCharge ( G4double  finalCharge  )  [inline]

Definition at line 170 of file G4ParticleChangeForLoss.hh.

00171 {
00172   currentCharge = theCharge;
00173 }

void G4ParticleChangeForLoss::ProposeMomentumDirection ( const G4ThreeVector Pfinal  )  [inline]

Definition at line 188 of file G4ParticleChangeForLoss.hh.

00189 {
00190   proposedMomentumDirection = dir;
00191 }

void G4ParticleChangeForLoss::ProposeMomentumDirection ( G4double  Px,
G4double  Py,
G4double  Pz 
) [inline]

Definition at line 200 of file G4ParticleChangeForLoss.hh.

Referenced by G4EmBiasingManager::ApplySecondaryBiasing(), G4PenelopeIonisationModel::SampleSecondaries(), G4PenelopeBremsstrahlungModel::SampleSecondaries(), G4LivermoreIonisationModel::SampleSecondaries(), and G4LivermoreBremsstrahlungModel::SampleSecondaries().

00201 {
00202   proposedMomentumDirection.setX(Px);
00203   proposedMomentumDirection.setY(Py);
00204   proposedMomentumDirection.setZ(Pz);
00205 }

void G4ParticleChangeForLoss::ProposePolarization ( G4double  Px,
G4double  Py,
G4double  Pz 
) [inline]

Definition at line 225 of file G4ParticleChangeForLoss.hh.

00226 {
00227   proposedPolarization.setX(Px);
00228   proposedPolarization.setY(Py);
00229   proposedPolarization.setZ(Pz);
00230 }

void G4ParticleChangeForLoss::ProposePolarization ( const G4ThreeVector dir  )  [inline]

Definition at line 219 of file G4ParticleChangeForLoss.hh.

Referenced by G4PolarizedMollerBhabhaModel::SampleSecondaries(), and G4ePolarizedBremsstrahlungModel::SampleSecondaries().

00220 {
00221   proposedPolarization = dir;
00222 }

void G4ParticleChangeForLoss::SetLowEnergyLimit ( G4double  elimit  )  [inline]

Definition at line 273 of file G4ParticleChangeForLoss.hh.

Referenced by G4VEnergyLossProcess::BuildPhysicsTable().

00274 {
00275   lowEnergyLimit = elimit;
00276 }

void G4ParticleChangeForLoss::SetProposedCharge ( G4double  theCharge  )  [inline]

Definition at line 165 of file G4ParticleChangeForLoss.hh.

Referenced by G4VEnergyLossProcess::AlongStepDoIt().

00166 {
00167   currentCharge = theCharge;
00168 }

void G4ParticleChangeForLoss::SetProposedKineticEnergy ( G4double  proposedKinEnergy  )  [inline]

Definition at line 150 of file G4ParticleChangeForLoss.hh.

Referenced by G4VEnergyLossProcess::AlongStepDoIt(), G4NuclearStopping::AlongStepDoIt(), G4EmBiasingManager::ApplySecondaryBiasing(), G4SeltzerBergerModel::SampleSecondaries(), G4PolarizedMollerBhabhaModel::SampleSecondaries(), G4PenelopeIonisationModel::SampleSecondaries(), G4PenelopeBremsstrahlungModel::SampleSecondaries(), G4PAIPhotonModel::SampleSecondaries(), G4PAIModel::SampleSecondaries(), G4MuPairProductionModel::SampleSecondaries(), G4MuBremsstrahlungModel::SampleSecondaries(), G4MollerBhabhaModel::SampleSecondaries(), G4LivermoreIonisationModel::SampleSecondaries(), G4LivermoreBremsstrahlungModel::SampleSecondaries(), G4ICRU73QOModel::SampleSecondaries(), G4eBremsstrahlungRelModel::SampleSecondaries(), G4eBremsstrahlungModel::SampleSecondaries(), G4eBremParametrizedModel::SampleSecondaries(), and G4BraggIonModel::SampleSecondaries().

00151 {
00152   proposedKinEnergy = energy;
00153 }

void G4ParticleChangeForLoss::SetProposedMomentumDirection ( const G4ThreeVector dir  )  [inline]

Definition at line 194 of file G4ParticleChangeForLoss.hh.

Referenced by G4SeltzerBergerModel::SampleSecondaries(), G4PolarizedMollerBhabhaModel::SampleSecondaries(), G4PAIPhotonModel::SampleSecondaries(), G4PAIModel::SampleSecondaries(), G4MuPairProductionModel::SampleSecondaries(), G4MuBremsstrahlungModel::SampleSecondaries(), G4MollerBhabhaModel::SampleSecondaries(), G4ICRU73QOModel::SampleSecondaries(), G4eBremsstrahlungRelModel::SampleSecondaries(), G4eBremsstrahlungModel::SampleSecondaries(), G4eBremParametrizedModel::SampleSecondaries(), and G4BraggIonModel::SampleSecondaries().

00195 {
00196   proposedMomentumDirection = dir;
00197 }

G4Step * G4ParticleChangeForLoss::UpdateStepForAlongStep ( G4Step Step  )  [virtual]

Reimplemented from G4VParticleChange.

Definition at line 210 of file G4ParticleChangeForLoss.cc.

References G4Step::AddNonIonizingEnergyDeposit(), G4Step::AddTotalEnergyDeposit(), G4Track::CalculateVelocity(), G4StepPoint::GetKineticEnergy(), G4Step::GetPostStepPoint(), G4Step::GetPreStepPoint(), G4Step::GetTrack(), G4VParticleChange::isParentWeightProposed, G4StepPoint::SetCharge(), G4StepPoint::SetKineticEnergy(), G4Track::SetKineticEnergy(), G4StepPoint::SetVelocity(), G4StepPoint::SetWeight(), G4VParticleChange::theLocalEnergyDeposit, G4VParticleChange::theNonIonizingEnergyDeposit, and G4VParticleChange::theParentWeight.

00211 {
00212   G4StepPoint* pPostStepPoint = pStep->GetPostStepPoint();
00213   G4StepPoint* pPreStepPoint  = pStep->GetPreStepPoint();
00214   G4Track* pTrack = pStep->GetTrack();
00215 
00216   // accumulate change of the kinetic energy
00217   G4double preKinEnergy = pPreStepPoint->GetKineticEnergy();
00218   G4double kinEnergy = pPostStepPoint->GetKineticEnergy() +
00219     (proposedKinEnergy - preKinEnergy);
00220 
00221   // update kinetic energy and charge
00222   if (kinEnergy < lowEnergyLimit) {
00223     theLocalEnergyDeposit += kinEnergy;
00224     kinEnergy = 0.0;
00225     pPostStepPoint->SetVelocity(0.0);
00226   } else {
00227     pPostStepPoint->SetCharge( currentCharge );
00228     // calculate velocity
00229     pTrack->SetKineticEnergy(kinEnergy); 
00230     pPostStepPoint->SetVelocity(pTrack->CalculateVelocity());
00231     pTrack->SetKineticEnergy(preKinEnergy); 
00232   }
00233   pPostStepPoint->SetKineticEnergy( kinEnergy );
00234 
00235   if (isParentWeightProposed ){
00236     pPostStepPoint->SetWeight( theParentWeight );
00237   }
00238 
00239   pStep->AddTotalEnergyDeposit( theLocalEnergyDeposit );
00240   pStep->AddNonIonizingEnergyDeposit( theNonIonizingEnergyDeposit );
00241   return pStep;
00242 }

G4Step * G4ParticleChangeForLoss::UpdateStepForPostStep ( G4Step Step  )  [virtual]

Reimplemented from G4VParticleChange.

Definition at line 244 of file G4ParticleChangeForLoss.cc.

References G4Step::AddNonIonizingEnergyDeposit(), G4Step::AddTotalEnergyDeposit(), G4Track::CalculateVelocity(), G4Step::GetPostStepPoint(), G4Step::GetTrack(), G4VParticleChange::isParentWeightProposed, G4StepPoint::SetCharge(), G4Track::SetKineticEnergy(), G4StepPoint::SetKineticEnergy(), G4StepPoint::SetMomentumDirection(), G4StepPoint::SetPolarization(), G4StepPoint::SetVelocity(), G4StepPoint::SetWeight(), G4VParticleChange::theLocalEnergyDeposit, G4VParticleChange::theNonIonizingEnergyDeposit, and G4VParticleChange::theParentWeight.

00245 {
00246   G4StepPoint* pPostStepPoint = pStep->GetPostStepPoint();
00247   G4Track* pTrack = pStep->GetTrack();
00248 
00249   pPostStepPoint->SetCharge( currentCharge );
00250   pPostStepPoint->SetMomentumDirection( proposedMomentumDirection );
00251   pPostStepPoint->SetKineticEnergy( proposedKinEnergy );
00252   pTrack->SetKineticEnergy( proposedKinEnergy );
00253   if(proposedKinEnergy > 0.0) {
00254     pPostStepPoint->SetVelocity(pTrack->CalculateVelocity());
00255   } else {
00256     pPostStepPoint->SetVelocity(0.0);
00257   }
00258   pPostStepPoint->SetPolarization( proposedPolarization );
00259 
00260   if (isParentWeightProposed ){
00261     pPostStepPoint->SetWeight( theParentWeight );
00262   }
00263 
00264   pStep->AddTotalEnergyDeposit( theLocalEnergyDeposit );
00265   pStep->AddNonIonizingEnergyDeposit( theNonIonizingEnergyDeposit );
00266   return pStep;
00267 }


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