G4V3DNucleus.hh

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 //
00027 // $Id$
00028 //
00029 // 20110805  M. Kelsey -- Change nucleon vector to use objects, not pointers
00030 
00031 #ifndef G4V3DNucleus_h
00032 #define G4V3DNucleus_h 1
00033 
00034 class G4Nucleon;
00035 class G4VNuclearDensity;
00036 #include "G4DynamicParticle.hh"
00037 #include "Randomize.hh"
00038 #include <utility>  
00039 #include <vector>
00040 
00041 class G4V3DNucleus 
00042 {
00043 
00044   public:
00045       G4V3DNucleus();
00046       virtual ~G4V3DNucleus();
00047 
00048   private:
00049       G4V3DNucleus(const G4V3DNucleus &right);
00050       const G4V3DNucleus & operator=(const G4V3DNucleus &right);
00051       int operator==(const G4V3DNucleus &right) const;
00052       int operator!=(const G4V3DNucleus &right) const;
00053 
00054   public:
00055       virtual void Init(G4int theA, G4int theZ) = 0;
00056       virtual G4bool StartLoop() = 0;
00057       virtual G4Nucleon * GetNextNucleon() = 0;
00058       virtual const std::vector<G4Nucleon> & GetNucleons() = 0;
00059       virtual G4int GetMassNumber() = 0;
00060       virtual G4double GetMass() = 0;
00061       virtual G4int GetCharge() = 0;
00062       virtual G4double GetNuclearRadius() = 0;
00063       virtual G4double GetNuclearRadius(const G4double maxRelativeDensity) = 0;
00064       virtual G4double GetOuterRadius() = 0;
00065       virtual G4double CoulombBarrier() = 0;
00066       virtual void DoLorentzBoost(const G4LorentzVector & theBoost) = 0;
00067       virtual void DoLorentzBoost(const G4ThreeVector & theBeta) = 0;
00068       virtual void DoLorentzContraction(const G4LorentzVector & theBoost) = 0;
00069       virtual void DoLorentzContraction(const G4ThreeVector & theBeta) = 0;
00070       virtual void DoTranslation(const G4ThreeVector & theShift) = 0;
00071       virtual const G4VNuclearDensity * GetNuclearDensity() const = 0;
00072       virtual void SortNucleonsIncZ() = 0;
00073       virtual void SortNucleonsDecZ() = 0;
00074 
00075   public:
00076       std::pair<G4double, G4double> ChooseImpactXandY(G4double maxImpact);
00077       std::pair<G4double, G4double> RefetchImpactXandY(){return theImpactParameter;}
00078 
00079   private:
00080   
00081     std::pair<G4double, G4double> theImpactParameter;
00082 
00083 };
00084 
00085 inline
00086 std::pair<G4double, G4double> G4V3DNucleus::
00087 ChooseImpactXandY(G4double maxImpact)
00088 {
00089   G4double x,y;
00090   do
00091   {
00092     x = 2*G4UniformRand() - 1;
00093     y = 2*G4UniformRand() - 1;
00094   }
00095   while(x*x + y*y > 1);
00096   G4double impactX = x*(maxImpact); 
00097   G4double impactY = y*(maxImpact);
00098   theImpactParameter.first = impactX;
00099   theImpactParameter.second = impactY;
00100   return theImpactParameter;
00101 }
00102 
00103 
00104 #endif
00105 
00106 

Generated on Mon May 27 17:50:09 2013 for Geant4 by  doxygen 1.4.7