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 #define INCLXX_IN_GEANT4_MODE 1
00034
00035 #include "globals.hh"
00036
00044 #ifndef G4INCLPROJECTILEREMNANT_HH_
00045 #define G4INCLPROJECTILEREMNANT_HH_
00046
00047 #include "G4INCLCluster.hh"
00048 #include "G4INCLRandom.hh"
00049 #include <vector>
00050 #include <map>
00051 #include <numeric>
00052 #include <functional>
00053
00054 namespace G4INCL {
00055
00057 G4int shuffleComponentsHelper(G4int range);
00058
00059 class ProjectileRemnant : public Cluster {
00060
00061 typedef std::vector<G4double> EnergyLevels;
00062 typedef std::map<long, G4double> EnergyLevelMap;
00063
00064 public:
00065 ProjectileRemnant(ParticleSpecies const species, const G4double kineticEnergy)
00066 : Cluster(species.theZ, species.theA) {
00067
00068
00069 setTableMass();
00070
00071
00072 const G4double projectileMass = getMass();
00073 const G4double energy = kineticEnergy + projectileMass;
00074 const G4double momentumZ = std::sqrt(energy*energy - projectileMass*projectileMass);
00075
00076
00077 initializeParticles();
00078 internalBoostToCM();
00079 putParticlesOffShell();
00080
00081
00082
00083 storeEnergyLevels();
00084
00085
00086 const ThreeVector aBoostVector = ThreeVector(0.0, 0.0, momentumZ / energy);
00087 boost(-aBoostVector);
00088
00089
00090 freezeInternalMotion();
00091
00092
00093 makeProjectileSpectator();
00094 }
00095
00096 ~ProjectileRemnant() {
00097 deleteStoredComponents();
00098 clearEnergyLevels();
00099 }
00100
00102 void reset();
00103
00109 void removeParticle(Particle * const p, const G4double theProjectileCorrection);
00110
00119 ParticleList addDynamicalSpectators(ParticleList pL);
00120
00130 ParticleList addMostDynamicalSpectators(ParticleList pL);
00131
00133 void deleteStoredComponents() {
00134 for(std::map<long,Particle*>::const_iterator p=storedComponents.begin(); p!=storedComponents.end(); ++p)
00135 delete p->second;
00136 clearStoredComponents();
00137 }
00138
00140 void clearStoredComponents() {
00141 storedComponents.clear();
00142 }
00143
00145 void clearEnergyLevels() {
00146 theInitialEnergyLevels.clear();
00147 theGroundStateEnergies.clear();
00148 }
00149
00159 G4double computeExcitationEnergy(const long exceptID) const;
00160
00161 EnergyLevels getPresentEnergyLevels(const long exceptID) const {
00162 EnergyLevels theEnergyLevels;
00163 for(ParticleIter p=particles.begin(); p!=particles.end(); ++p) {
00164 if((*p)->getID()!=exceptID) {
00165 EnergyLevelMap::const_iterator i = theInitialEnergyLevels.find((*p)->getID());
00166
00167 theEnergyLevels.push_back(i->second);
00168 }
00169 }
00170
00171 return theEnergyLevels;
00172 }
00173
00175 void storeComponents() {
00176 for(ParticleIter p=particles.begin(); p!=particles.end(); ++p) {
00177
00178 storedComponents[(*p)->getID()]=new Particle(**p);
00179 }
00180 }
00181
00183 G4int getNumberStoredComponents() const {
00184 return storedComponents.size();
00185 }
00186
00188 void storeEnergyLevels() {
00189 EnergyLevels energies;
00190
00191 for(ParticleIter p=particles.begin(); p!=particles.end(); ++p) {
00192 const G4double theCMEnergy = (*p)->getEnergy();
00193
00194 theInitialEnergyLevels[(*p)->getID()] = theCMEnergy;
00195 energies.push_back(theCMEnergy);
00196 }
00197
00198 std::sort(energies.begin(), energies.end());
00199
00200 theGroundStateEnergies.resize(energies.size());
00201
00202
00203 std::partial_sum(energies.begin(), energies.end(), theGroundStateEnergies.begin());
00204 }
00205
00206 private:
00207
00209 ParticleList shuffleStoredComponents() {
00210 ParticleList pL = getStoredComponents();
00211 std::vector<Particle *> theVector(pL.begin(),pL.end());
00212 std::random_shuffle(theVector.begin(), theVector.end(), shuffleComponentsHelper);
00213 return ParticleList(theVector.begin(),theVector.end());
00214 }
00215
00216 ParticleList getStoredComponents() const {
00217 ParticleList pL;
00218 for(std::map<long,Particle*>::const_iterator p=storedComponents.begin(); p!=storedComponents.end(); ++p)
00219 pL.push_back(p->second);
00220 return pL;
00221 }
00222
00224 ThreeVector const &getStoredMomentum(Particle const * const p) const {
00225 std::map<long,Particle*>::const_iterator i = storedComponents.find(p->getID());
00226 if(i==storedComponents.end()) {
00227 ERROR("Couldn't find particle " << p->getID() << " in the list of projectile components" << std::endl);
00228 return p->getMomentum();
00229 } else {
00230 return i->second->getMomentum();
00231 }
00232 }
00233
00240 G4bool addDynamicalSpectator(Particle * const p);
00241
00243
00244
00245
00246
00247
00248
00249
00250
00251
00252
00254 std::map<long, Particle*> storedComponents;
00255
00257 EnergyLevelMap theInitialEnergyLevels;
00258
00260 EnergyLevels theGroundStateEnergies;
00261
00262 };
00263 }
00264
00265 #endif // G4INCLPROJECTILEREMNANT_HH_
00266