G4VPartonStringModel.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 //
00030 //      GEANT 4 class implementation file
00031 //
00032 //      ---------------- G4VPartonStringModel ----------------
00033 //             by Gunter Folger, May 1998.
00034 //      abstract class for all Parton String Models
00035 // ------------------------------------------------------------
00036 // debug switch
00037 //#define debug_PartonStringModel
00038 
00039 
00040 #include "G4VPartonStringModel.hh"
00041 #include "G4ios.hh"
00042 #include "G4ShortLivedConstructor.hh"
00043 
00044 #include "G4ParticleTable.hh"
00045 #include "G4IonTable.hh"
00046 
00047 
00048 G4VPartonStringModel::G4VPartonStringModel(const G4String& modelName)
00049     : G4VHighEnergyGenerator(modelName),
00050       stringFragmentationModel(0),
00051       theThis(0)
00052 {
00053 //  Make shure Shotrylived partyicles are constructed.
00054         G4ShortLivedConstructor ShortLived;
00055         ShortLived.ConstructParticle();
00056 }
00057 
00058 G4VPartonStringModel::~G4VPartonStringModel()
00059 {
00060 }
00061 
00062 G4KineticTrackVector * G4VPartonStringModel::Scatter(const G4Nucleus &theNucleus, 
00063                                                 const G4DynamicParticle &aPrimary)
00064 {  
00065   G4ExcitedStringVector * strings = NULL;
00066 
00067   G4DynamicParticle thePrimary=aPrimary;
00068   
00069   G4LorentzRotation toZ;
00070   G4LorentzVector Ptmp=thePrimary.Get4Momentum();
00071   toZ.rotateZ(-1*Ptmp.phi());
00072   toZ.rotateY(-1*Ptmp.theta());
00073   thePrimary.Set4Momentum(toZ*Ptmp);
00074   G4LorentzRotation toLab(toZ.inverse());
00075 
00076   G4int attempts = 0, maxAttempts=20;
00077   while ( strings  == NULL )
00078   {
00079         if (attempts++ > maxAttempts ) 
00080         {
00081                 throw G4HadronicException(__FILE__, __LINE__, "G4VPartonStringModel::Scatter(): fails to generate strings");
00082         }
00083 
00084         theThis->Init(theNucleus,thePrimary);
00085 
00086         strings = GetStrings();
00087   }
00088   
00089   G4double stringEnergy(0);
00090   G4LorentzVector SumStringMom(0.,0.,0.,0.);
00091 
00092   for ( unsigned int astring=0; astring < strings->size(); astring++)
00093   {
00094 //    rotate string to lab frame, models have it aligned to z
00095     stringEnergy += (*strings)[astring]->GetLeftParton()->Get4Momentum().t();
00096     stringEnergy += (*strings)[astring]->GetRightParton()->Get4Momentum().t();
00097     (*strings)[astring]->LorentzRotate(toLab);
00098     SumStringMom+=(*strings)[astring]->Get4Momentum();
00099   }
00100 
00101   G4double InvMass=SumStringMom.mag();   
00102 
00103 //#define debug_PartonStringModel
00104 #ifdef debug_PartonStringModel
00105 
00106      G4V3DNucleus * fancynucleus=theThis->GetWoundedNucleus();
00107   
00108        // loop over wounded nucleus
00109      G4int hits(0), charged_hits(0);
00110      G4ThreeVector hitNucleonMomentum(0.,0.,0.);
00111      G4Nucleon * theCurrentNucleon = fancynucleus->StartLoop() ? fancynucleus->GetNextNucleon() : NULL;
00112      while( theCurrentNucleon )
00113      {
00114        if(theCurrentNucleon->AreYouHit()) 
00115        {
00116          hits++;
00117          hitNucleonMomentum += theCurrentNucleon->Get4Momentum().vect();
00118          if ( theCurrentNucleon->GetDefinition() == G4Proton::Proton() )  ++charged_hits;
00119        }
00120        theCurrentNucleon = fancynucleus->GetNextNucleon();
00121      }
00122      
00123      G4int initialZ=fancynucleus->GetCharge();
00124      G4int initialA=fancynucleus->GetMassNumber();
00125      G4double initial_mass=G4ParticleTable::GetParticleTable()->GetIonTable()->GetIonMass(initialZ,initialA);
00126      G4double final_mass = G4ParticleTable::GetParticleTable()->GetIonTable()->GetIonMass(initialZ-charged_hits, initialA-hits);
00127      G4cout << "G4VPSM: strE, hit nucleons, Primary, SumStringE, nucleus intial, nucleus final, excitation estimate "
00128             << stringEnergy << " "    << hits << ", " << Ptmp.e() << ", "<< SumStringMom.e() << ", "
00129                 << initial_mass<< ", " << final_mass<< ", "  << Ptmp.e() + initial_mass - final_mass - stringEnergy  << G4endl;
00130      G4cout << "momentum balance " <<  thePrimary.GetMomentum() + hitNucleonMomentum - SumStringMom.vect()<< G4endl;
00131 #endif
00132 
00133 //  Fragment strings
00134 
00135   G4KineticTrackVector * theResult = 0;
00136   G4double SumMass(0.); 
00137   attempts = 0; 
00138   maxAttempts=100;
00139   do 
00140   {         
00141    attempts++;   
00142    if(theResult != 0)
00143    {
00144     std::for_each(theResult->begin(), theResult->end(), DeleteKineticTrack());
00145     delete theResult;
00146    }
00147    theResult = stringFragmentationModel->FragmentStrings(strings);
00148    if(attempts > maxAttempts ) break;
00149 
00150     //G4cout<<"G4endl<<"G4VPartonStringModel:: Final Result, Size "<<theResult->size()<<G4endl;
00151 
00152    SumMass=0.;
00153     //G4LorentzVector SumP(0.,0.,0.,0.);
00154    for ( unsigned int i=0; i < theResult->size(); i++)
00155    {
00156     SumMass+=(*theResult)[i]->GetDefinition()->GetPDGMass();
00157      //SumP+=(*theResult)[i]->Get4Momentum();
00158      //G4cout<<i<<" : "<<(*theResult)[i]->GetDefinition()->GetParticleName();
00159      //G4cout<<"p= " << (*theResult)[i]->Get4Momentum()<<" m= "<<(*theResult)[i]->Get4Momentum().mag()<<G4endl;
00160    }
00161 
00162      //G4cout<<"SumP "<<SumP<<G4endl;
00163   } while(SumMass > InvMass);
00164 
00165   std::for_each(strings->begin(), strings->end(), DeleteString() );
00166   delete strings;
00167 
00168   return theResult;
00169 }
00170 
00171 void G4VPartonStringModel::ModelDescription(std::ostream& outFile) const
00172 {
00173         outFile << GetModelName() << " has no description yet.\n";
00174 }

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