G4FTFParticipants.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 //
00027 // $Id$
00028 // GEANT4 tag $Name:  $
00029 //
00030 // ------------------------------------------------------------
00031 //      GEANT 4 class implementation file
00032 //
00033 //      ---------------- G4FTFParticipants----------------
00034 //             by Gunter Folger, June 1998.
00035 //       class finding colliding particles in FTFPartonStringModel
00036 //  Changed in a part by V. Uzhinsky in oder to put in correcpondence
00037 //        with original FRITIOF mode. November - December 2006.
00038 //  Ajusted for (anti) nucleus - nucleus interactions by V. Uzhinsky.
00039 //                    (February 2011)
00040 // ------------------------------------------------------------
00041 
00042 #include <utility>
00043 #include <vector>
00044 #include <algorithm>
00045 
00046 #include "G4FTFParticipants.hh"
00047 #include "G4ios.hh"
00048 #include "Randomize.hh"
00049 #include "G4SystemOfUnits.hh"
00050 #include "G4FTFParameters.hh"                            // Uzhi 29.03.08
00051 #include "G4DiffractiveSplitableHadron.hh"
00052 #include "G4VSplitableHadron.hh"
00053 
00054 G4FTFParticipants::G4FTFParticipants() :
00055   theProjectileNucleus(0),
00056   currentInteraction(-1)
00057 {
00058 }
00059 
00060 G4FTFParticipants::G4FTFParticipants(const G4FTFParticipants &): G4VParticipants()
00061   , theProjectileNucleus(0), currentInteraction(-1)   //A.R. 14-Aug-2012 Coverity fix.
00062 {
00063         G4Exception("G4FTFParticipants::G4FTFParticipants()","HAD_FTF_001",
00064                 FatalException," Must not use copy ctor()");
00065 }
00066 
00067 
00068 G4FTFParticipants::~G4FTFParticipants()
00069 {
00070         if ( theProjectileNucleus != NULL ) delete theProjectileNucleus;
00071 }
00072 
00073 //-------------------------------------------------------------------------
00074 
00075 void G4FTFParticipants::SetProjectileNucleus(G4V3DNucleus * aNucleus)
00076 {
00077   if (theProjectileNucleus) delete theProjectileNucleus;
00078   
00079   theProjectileNucleus = aNucleus;
00080 }
00081 
00082 G4V3DNucleus * G4FTFParticipants::GetProjectileNucleus()
00083 {
00084   return theProjectileNucleus;
00085 }
00086 
00087 void G4FTFParticipants::InitProjectileNucleus(G4int theA, G4int theZ)
00088 {
00089         if ( theProjectileNucleus == NULL ) theProjectileNucleus = new G4Fancy3DNucleus();
00090         theProjectileNucleus->Init(theA, theZ);
00091         theProjectileNucleus->SortNucleonsDecZ();
00092 }
00093 //-------------------------------------------------------------------------
00094 void G4FTFParticipants::GetList(const G4ReactionProduct  &thePrimary,
00095                                       G4FTFParameters    *theParameters) 
00096 { 
00097 //G4cout<<"Participants::GetList"<<G4endl;
00098 //G4cout<<"thePrimary "<<thePrimary.GetMomentum()<<G4endl;
00099     StartLoop();  // reset Loop over Interactions
00100 
00101     for(unsigned int i=0; i<theInteractions.size(); i++) delete theInteractions[i];
00102     theInteractions.clear();
00103 
00104     G4double deltaxy=2 * fermi;                       // Extra nuclear radius
00105 //G4cout<<"theProjectileNucleus "<<theProjectileNucleus<<G4endl;
00106     if(theProjectileNucleus == 0)
00107     { // Hadron-nucleus or anti-baryon-nucleus interactions
00108 //G4cout<<"Hadron-nucleus or anti-baryon-nucleus interactions"<<G4endl;
00109 
00110      G4double impactX(0.), impactY(0.);
00111 
00112      G4VSplitableHadron * primarySplitable=new G4DiffractiveSplitableHadron(thePrimary);
00113 //G4cout<<"Prim in Part "<<primarySplitable->Get4Momentum()<<G4endl;
00114      G4double xyradius;                          
00115      xyradius =theNucleus->GetOuterRadius() + deltaxy; // Impact parameter sampling
00116                                                     
00117 //    G4bool nucleusNeedsShift = true;                // Uzhi 20 July 2009
00118     
00119      do 
00120      {
00121          std::pair<G4double, G4double> theImpactParameter;
00122          theImpactParameter = theNucleus->ChooseImpactXandY(xyradius);
00123          impactX = theImpactParameter.first; 
00124          impactY = theImpactParameter.second;
00125 
00126          G4ThreeVector thePosition(impactX, impactY, -DBL_MAX);     
00127          primarySplitable->SetPosition(thePosition);                
00128 
00129          theNucleus->StartLoop();
00130          G4Nucleon * nucleon;
00131 
00132 G4int TrN(0);
00133          while ( (nucleon=theNucleus->GetNextNucleon()) ) 
00134          {
00135            G4double impact2= sqr(impactX - nucleon->GetPosition().x()) +
00136                              sqr(impactY - nucleon->GetPosition().y());
00137 
00138            if ( theParameters->GetProbabilityOfInteraction(impact2/fermi/fermi) 
00139                 > G4UniformRand() )
00140            {
00141                 primarySplitable->SetStatus(1);        // It takes part in the interaction
00142 
00143                 G4VSplitableHadron * targetSplitable=0;
00144                 if ( ! nucleon->AreYouHit() )
00145                 {
00146                     targetSplitable= new G4DiffractiveSplitableHadron(*nucleon);
00147                     nucleon->Hit(targetSplitable);
00148                     nucleon->SetBindingEnergy(3.*nucleon->GetBindingEnergy()); 
00149 //G4cout<<" Part nucl "<<TrN<<" "<<nucleon->Get4Momentum()<<G4endl;
00150 //G4cout<<" Part nucl "<<G4endl;
00151                     targetSplitable->SetStatus(1);     // It takes part in the interaction
00152                 }
00153                 G4InteractionContent * aInteraction = 
00154                                        new G4InteractionContent(primarySplitable);
00155                 aInteraction->SetTarget(targetSplitable);
00156                 aInteraction->SetTargetNucleon(nucleon);     // Uzhi 16.07.09
00157                 aInteraction->SetStatus(1);                  // Uzhi Feb26
00158                 theInteractions.push_back(aInteraction);
00159            }
00160 TrN++;
00161          } 
00162      } while ( theInteractions.size() == 0 );
00163 
00164      //if ( theInteractions.size() == 0 ) delete primarySplitable; //A.R. 14-Aug-2012 Coverity fix
00165 
00166 //      G4cout << "Number of Hit nucleons " << theInteractions.size()
00167 //              << "\t" << impactX/fermi << "\t"<<impactY/fermi
00168 //              << "\t" << std::sqrt(sqr(impactX)+sqr(impactY))/fermi <<G4endl;
00169 
00170      return;
00171     }       // end of if(theProjectileNucleus == 0)
00172 
00173 //-------------------------------------------------------------------
00174 //                Projectile and target are nuclei
00175 //-------------------------------------------------------------------
00176 //VU    G4VSplitableHadron * primarySplitable=new G4DiffractiveSplitableHadron(thePrimary);
00177 //G4cout<<"Prim in Part "<<primarySplitable->Get4Momentum()<<G4endl;
00178 //G4cout<<"Projectile and target are nuclei"<<G4endl;
00179 //G4cout<<thePrimary.GetMomentum()<<G4endl;
00180 //G4cout<<"Part Pr Tr "<<theProjectileNucleus<<" "<<theNucleus<<G4endl;
00181 
00182 
00183     G4double xyradius;                          
00184     xyradius =theProjectileNucleus->GetOuterRadius() +  // Impact parameter sampling
00185                         theNucleus->GetOuterRadius() + deltaxy;
00186 
00187     G4double impactX(0.), impactY(0.);
00188 
00189     do
00190     {
00191 //G4cout<<"New interaction list"<<G4endl;
00192          std::pair<G4double, G4double> theImpactParameter;
00193          theImpactParameter = theNucleus->ChooseImpactXandY(xyradius);
00194          impactX = theImpactParameter.first; 
00195          impactY = theImpactParameter.second;
00196 //G4cout<<"B "<<std::sqrt(sqr(impactX)+sqr(impactY))/fermi<<G4endl;
00197 
00198          G4ThreeVector thePosition(impactX, impactY, -DBL_MAX);     
00199 //VU         primarySplitable->SetPosition(thePosition);                
00200 
00201          theProjectileNucleus->StartLoop();
00202          G4Nucleon * ProjectileNucleon;
00203 G4int PrNuclN(0);
00204 
00205          while ( (ProjectileNucleon=theProjectileNucleus->GetNextNucleon()) ) 
00206          {
00207            G4VSplitableHadron * ProjectileSplitable=0;
00208 //G4cout<<G4endl<<"Prj N mom "<<ProjectileNucleon->Get4Momentum()<<"-------------"<<G4endl;
00209            theNucleus->StartLoop();
00210            G4Nucleon * TargetNucleon;
00211 
00212 G4int TrNuclN(0);
00213            while ( (TargetNucleon=theNucleus->GetNextNucleon()) )
00214            {
00215 //G4cout<<"Trg N mom "<<TargetNucleon->Get4Momentum()<<G4endl;
00216             G4double impact2=
00217             sqr(impactX+ProjectileNucleon->GetPosition().x()-TargetNucleon->GetPosition().x())+
00218             sqr(impactY+ProjectileNucleon->GetPosition().y()-TargetNucleon->GetPosition().y());
00219 
00220             G4VSplitableHadron * TargetSplitable=0;
00221 
00222             if ( theParameters->GetProbabilityOfInteraction(impact2/fermi/fermi) 
00223                  > G4UniformRand() )
00224             { // An Interaction has happend!
00225 //G4cout<<"An Interaction has happend"<<G4endl;
00226 //G4cout<<"PrN TrN "<<PrNuclN<<" "<<TrNuclN<<" "<<ProjectileNucleon->GetPosition().z()/fermi<<" "<<TargetNucleon->GetPosition().z()/fermi<<" "<<ProjectileNucleon->GetPosition().z()/fermi + TargetNucleon->GetPosition().z()/fermi <<G4endl;
00227 
00228              if ( ! ProjectileNucleon->AreYouHit() )
00229              { // Projectile nucleon was not involved until now.
00230               ProjectileSplitable= new G4DiffractiveSplitableHadron(*ProjectileNucleon);
00231               ProjectileNucleon->Hit(ProjectileSplitable);
00232               ProjectileNucleon->SetBindingEnergy(3.*ProjectileNucleon->GetBindingEnergy());
00233               ProjectileSplitable->SetStatus(1);     // It takes part in the interaction
00234              }
00235              else
00236              {  // Projectile nucleon was involved before.
00237               ProjectileSplitable=ProjectileNucleon->GetSplitableHadron();
00238              } // End of if ( ! Projectileucleon->AreYouHit() )
00239 
00240              if ( ! TargetNucleon->AreYouHit() )
00241              {  // Target nucleon was not involved until now
00242               TargetSplitable= new G4DiffractiveSplitableHadron(*TargetNucleon);
00243               TargetNucleon->Hit(TargetSplitable);
00244               TargetNucleon->SetBindingEnergy(3.*ProjectileNucleon->GetBindingEnergy());
00245               TargetSplitable->SetStatus(1);     // It takes part in the interaction
00246              }
00247              else
00248              {  // Target nucleon was involved before.
00249               TargetSplitable=TargetNucleon->GetSplitableHadron();
00250              } // End of if ( ! TargetNeucleon->AreYouHit() )
00251 
00252              G4InteractionContent * anInteraction = 
00253                                    new G4InteractionContent(ProjectileSplitable);
00254              anInteraction->SetTarget(TargetSplitable);
00255              anInteraction->SetTargetNucleon(TargetNucleon);
00256              anInteraction->SetStatus(1);                      // Uzhi Feb26
00257 //             anInteraction->SetInteractionTime(ProjectileNucleon->GetPosition().z()+
00258 //                                                   TargetNucleon->GetPosition().z());
00259 //G4cout<<"Z's pr tr "<<ProjectileNucleon->GetPosition().z()/fermi<<" "<<TargetNucleon->GetPosition().z()/fermi<<" "<<ProjectileNucleon->GetPosition().z()/fermi + TargetNucleon->GetPosition().z()/fermi <<G4endl;
00260              theInteractions.push_back(anInteraction);
00261 //G4cout<<"Ppr tr "<<ProjectileSplitable<<" "<<TargetSplitable<<G4endl;
00262             } // End of An Interaction has happend!
00263 TrNuclN++;
00264            } // End of while ( (TargetNucleon=theNucleus->GetNextNucleon()) )
00265 PrNuclN++;
00266          } // End of   while ( (ProjectileNucleon=theProjectileNucleus->GetNextNucleon()) )
00267     }  while ( theInteractions.size() == 0 );  // end of while ( theInteractions.size() == 0 )
00268 
00269 //std::sort(theInteractions.begin(),theInteractions.end()); // ????
00270 
00271 //      G4cout << "Number of primary collisions " << theInteractions.size() 
00272 //              << "\t" << impactX/fermi << "\t"<<impactY/fermi
00273 //              << "\t" << std::sqrt(sqr(impactX)+sqr(impactY))/fermi <<G4endl;
00274 //G4int Uzhi; G4cin >> Uzhi;
00275     return;
00276 }
00277 //--------------------------------------------------------------
00278 
00279 // Implementation (private) methods

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