Changeset 962 for trunk/source/processes/hadronic/models/parton_string/diffraction/src/G4FTFModel.cc
- Timestamp:
- Apr 6, 2009, 12:30:29 PM (15 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/source/processes/hadronic/models/parton_string/diffraction/src/G4FTFModel.cc
r819 r962 25 25 // 26 26 // 27 // $Id: G4FTFModel.cc,v 1. 7 2007/04/24 10:32:59 gunterExp $28 // GEANT4 tag $Name: geant4-09-0 1-patch-02 $27 // $Id: G4FTFModel.cc,v 1.13 2008/12/09 10:40:52 vuzhinsk Exp $ 28 // GEANT4 tag $Name: geant4-09-02-ref-02 $ 29 29 // 30 30 … … 38 38 39 39 #include "G4FTFModel.hh" 40 #include "G4FTFParameters.hh" // Uzhi 29.03.08 40 41 #include "G4FTFParticipants.hh" 41 42 #include "G4InteractionContent.hh" … … 43 44 #include "G4ParticleDefinition.hh" 44 45 #include "G4ios.hh" 46 #include <utility> // Uzhi 29.03.08 45 47 46 48 // Class G4FTFModel 47 49 48 G4FTFModel::G4FTFModel():theExcitation(new G4DiffractiveExcitation()) // Uzhi 50 G4FTFModel::G4FTFModel():theExcitation(new G4DiffractiveExcitation()), 51 theElastic(new G4ElasticHNScattering()) // Uzhi 29.03.08 49 52 { 50 53 G4VPartonStringModel::SetThisPointer(this); 51 } 52 53 G4FTFModel::G4FTFModel(G4double a, G4double b, G4double c):theExcitation(new G4DiffractiveExcitation()) 54 theParameters=0; // Uzhi 9.12.08 55 } 56 57 /* 58 G4FTFModel::G4FTFModel(G4double , G4double , G4double ):theExcitation(new // Uzhi 9.12.08 G4DiffractiveExcitation()) 54 59 { 55 60 G4VPartonStringModel::SetThisPointer(this); … … 62 67 G4VPartonStringModel::SetThisPointer(this); 63 68 } 64 69 */ 65 70 66 71 67 72 G4FTFModel::~G4FTFModel() 68 {} 73 { 74 if( theParameters != 0 ) delete theParameters; // Uzhi 5.12.08 75 // Because FTF model can be called for various particles 76 // theParameters must be erased at the end of each call. 77 // Thus the delete is olso in G4FTFModel::GetStrings() method 78 if( theExcitation != 0 ) delete theExcitation; // Uzhi 5.12.08 79 if( theElastic != 0 ) delete theElastic; // Uzhi 5.12.08 80 } 69 81 70 82 … … 86 98 } 87 99 100 // ------------------------------------------------------------ 88 101 void G4FTFModel::Init(const G4Nucleus & aNucleus, const G4DynamicParticle & aProjectile) 89 102 { 90 theParticipants.Init(aNucleus.GetN(),aNucleus.GetZ()); // Uzhi N-mass number Z-charge91 92 103 theProjectile = aProjectile; 93 } 94 104 //G4cout<<"G4FTFModel::Init "<<aNucleus.GetN()<<" "<<aNucleus.GetZ()<<G4endl; 105 theParticipants.Init(aNucleus.GetN(),aNucleus.GetZ()); 106 // Uzhi N-mass number Z-charge ------------------------- Uzhi 29.03.08 107 108 // --- cms energy 109 110 G4double s = sqr( theProjectile.GetMass() ) + 111 sqr( G4Proton::Proton()->GetPDGMass() ) + 112 2*theProjectile.GetTotalEnergy()*G4Proton::Proton()->GetPDGMass(); 113 /* 114 G4cout << " primary Total E (GeV): " << theProjectile.GetTotalEnergy()/GeV << G4endl; 115 G4cout << " primary Mass (GeV): " << theProjectile.GetMass() /GeV << G4endl; 116 G4cout << "cms std::sqrt(s) (GeV) = " << std::sqrt(s) / GeV << G4endl; 117 */ 118 119 if( theParameters != 0 ) delete theParameters; // Uzhi 9.12.08 120 theParameters = new G4FTFParameters(theProjectile.GetDefinition(), 121 aNucleus.GetN(),aNucleus.GetZ(), 122 s);// ------------------------- Uzhi 19.04.08 123 //theParameters->SetProbabilityOfElasticScatt(0.); // To turn on/off (1/0) elastic scattering 124 } 125 126 // ------------------------------------------------------------ 95 127 G4ExcitedStringVector * G4FTFModel::GetStrings() 96 128 { 97 theParticipants.BuildInteractions(theProjectile); 98 129 //G4cout<<"theParticipants.GetList"<<G4endl; 130 theParticipants.GetList(theProjectile,theParameters); 131 //G4cout<<"ExciteParticipants()"<<G4endl; 99 132 if (! ExciteParticipants()) return NULL;; 100 133 //G4cout<<"theStrings = BuildStrings()"<<G4endl; 101 134 G4ExcitedStringVector * theStrings = BuildStrings(); 102 135 //G4cout<<"Return to theStrings "<<G4endl; 136 if( theParameters != 0 ) // Uzhi 9.12.08 137 { // Uzhi 9.12.08 138 delete theParameters; // Uzhi 9.12.08 139 theParameters=0; // Uzhi 9.12.08 140 } // Uzhi 9.12.08 103 141 return theStrings; 104 142 } 105 143 144 // ------------------------------------------------------------ 106 145 struct DeleteVSplitableHadron { void operator()(G4VSplitableHadron * aH){delete aH;} }; 107 146 108 G4ExcitedStringVector * G4FTFModel::BuildStrings() 109 { 110 111 // Loop over all collisions; find all primaries, and all target ( targets may 112 // be duplicate in the List ( to unique G4VSplitableHadrons) 113 114 G4ExcitedStringVector * strings; 115 strings = new G4ExcitedStringVector(); 116 117 std::vector<G4VSplitableHadron *> primaries; 118 std::vector<G4VSplitableHadron *> targets; 119 120 theParticipants.StartLoop(); // restart a loop 121 while ( theParticipants.Next() ) 122 { 123 const G4InteractionContent & interaction=theParticipants.GetInteraction(); 124 // do not allow for duplicates ... 125 if ( primaries.end() == std::find(primaries.begin(), primaries.end(), interaction.GetProjectile()) ) 126 primaries.push_back(interaction.GetProjectile()); 127 128 if ( targets.end() == std::find(targets.begin(), targets.end(),interaction.GetTarget()) ) 129 targets.push_back(interaction.GetTarget()); 130 } 131 132 133 // G4cout << "BuildStrings prim/targ " << primaries.entries() << " , " << 134 // targets.entries() << G4endl; 135 136 137 unsigned int ahadron; 138 for ( ahadron=0; ahadron < primaries.size() ; ahadron++) 139 { 140 G4bool isProjectile=true; 141 strings->push_back(theExcitation->String(primaries[ahadron], isProjectile)); 142 } 143 for ( ahadron=0; ahadron < targets.size() ; ahadron++) 144 { 145 G4bool isProjectile=false; 146 strings->push_back(theExcitation->String(targets[ahadron], isProjectile)); 147 } 148 149 std::for_each(primaries.begin(), primaries.end(), DeleteVSplitableHadron()); 150 primaries.clear(); 151 std::for_each(targets.begin(), targets.end(), DeleteVSplitableHadron()); 152 targets.clear(); 153 154 return strings; 155 } 156 147 // ------------------------------------------------------------ 157 148 G4bool G4FTFModel::ExciteParticipants() 158 149 { 159 150 /* // Uzhi 29.03.08 For elastic Scatt. 151 G4cout<<" In ExciteParticipants() "<<theParticipants.theInteractions.size()<<G4endl; 152 G4cout<<" test Params Tot "<<theParameters->GetTotalCrossSection()<<G4endl; 153 G4cout<<" test Params Ela "<<theParameters->GetElasticCrossSection()<<G4endl; 154 155 G4int counter=0; 156 */ // Uzhi 29.03.08 157 158 159 160 160 while (theParticipants.Next()) 161 161 { 162 162 const G4InteractionContent & collision=theParticipants.GetInteraction(); 163 164 //G4cout << " soft colls : " << collision.GetNumberOfSoftCollisions() << G4endl; // Uzhi no match 163 /* 164 counter++; 165 G4cout<<" Inter # "<<counter<<G4endl; 166 */ 165 167 G4VSplitableHadron * projectile=collision.GetProjectile(); 166 168 G4VSplitableHadron * target=collision.GetTarget(); 167 169 168 if ( ! theExcitation->ExciteParticipants(projectile, target) ) 170 // // Uzhi 29.03.08 171 G4bool Successfull; 172 if(G4UniformRand()< theParameters->GetProbabilityOfElasticScatt()) 173 { 174 //G4cout<<"Elastic"<<G4endl; 175 Successfull=theElastic->ElasticScattering(projectile, target, theParameters); 176 } 177 else 178 { 179 //G4cout<<"Inelastic"<<G4endl; 180 Successfull=theExcitation->ExciteParticipants(projectile, target, theParameters); 181 } 182 // if(!Successfull) 183 // // Uzhi 29.03.08 184 185 // if ( ! theExcitation->ExciteParticipants(projectile, target) ) 186 if(!Successfull) 169 187 { 170 188 // give up, clean up … … 176 194 const G4InteractionContent & interaction=theParticipants.GetInteraction(); 177 195 // do not allow for duplicates ... 178 if ( primaries.end() == std::find(primaries.begin(), primaries.end(), interaction.GetProjectile()) ) 196 if ( primaries.end() == std::find(primaries.begin(), primaries.end(), 197 interaction.GetProjectile()) ) 179 198 primaries.push_back(interaction.GetProjectile()); 180 199 181 if ( targets.end() == std::find(targets.begin(), targets.end(),interaction.GetTarget()) ) 200 if ( targets.end() == std::find(targets.begin(), targets.end(), 201 interaction.GetTarget()) ) 182 202 targets.push_back(interaction.GetTarget()); 183 203 } … … 188 208 189 209 return false; 190 } 191 210 } // End of the loop Uzhi 211 } 212 return true; 213 } 214 // ------------------------------------------------------------ 215 G4ExcitedStringVector * G4FTFModel::BuildStrings() 216 { 217 // Loop over all collisions; find all primaries, and all target ( targets may 218 // be duplicate in the List ( to unique G4VSplitableHadrons) 219 220 G4ExcitedStringVector * strings; 221 strings = new G4ExcitedStringVector(); 222 223 std::vector<G4VSplitableHadron *> primaries; 224 std::vector<G4VSplitableHadron *> targets; 225 226 theParticipants.StartLoop(); // restart a loop 227 while ( theParticipants.Next() ) 228 { 229 const G4InteractionContent & interaction=theParticipants.GetInteraction(); 230 // do not allow for duplicates ... 231 if ( primaries.end() == std::find(primaries.begin(), primaries.end(), 232 interaction.GetProjectile()) ) 233 primaries.push_back(interaction.GetProjectile()); 234 235 if ( targets.end() == std::find(targets.begin(), targets.end(), 236 interaction.GetTarget()) ) 237 targets.push_back(interaction.GetTarget()); 192 238 } 193 return true; 194 } 239 240 241 // G4cout << "BuildStrings prim/targ " << primaries.size() << " , " << 242 // targets.size() << G4endl; 243 244 unsigned int ahadron; 245 // Only for hA-interactions Uzhi ------------------------------------- 246 for ( ahadron=0; ahadron < primaries.size() ; ahadron++) 247 { 248 //G4ThreeVector aPosition=primaries[ahadron]->GetPosition(); 249 //G4cout<<"Proj Build "<<aPosition<<" "<<primaries[ahadron]->GetTimeOfCreation()<<G4endl; 250 G4bool isProjectile=true; 251 strings->push_back(theExcitation->String(primaries[ahadron], isProjectile)); 252 } 253 254 for ( ahadron=0; ahadron < targets.size() ; ahadron++) 255 { 256 //G4ThreeVector aPosition=targets[ahadron]->GetPosition(); 257 //G4cout<<"Targ Build "<<aPosition<<" "<<targets[ahadron]->GetTimeOfCreation()<<G4endl; 258 G4bool isProjectile=false; 259 strings->push_back(theExcitation->String(targets[ahadron], isProjectile)); 260 } 261 262 std::for_each(primaries.begin(), primaries.end(), DeleteVSplitableHadron()); 263 primaries.clear(); 264 std::for_each(targets.begin(), targets.end(), DeleteVSplitableHadron()); 265 targets.clear(); 266 267 return strings; 268 } 269 // ------------------------------------------------------------
Note: See TracChangeset
for help on using the changeset viewer.