Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/processes/hadronic/models/parton_string/diffraction/src/G4DiffractiveSplitableHadron.cc

Version: [ ReleaseNotes ] [ 1.0 ] [ 1.1 ] [ 2.0 ] [ 3.0 ] [ 3.1 ] [ 3.2 ] [ 4.0 ] [ 4.0.p1 ] [ 4.0.p2 ] [ 4.1 ] [ 4.1.p1 ] [ 5.0 ] [ 5.0.p1 ] [ 5.1 ] [ 5.1.p1 ] [ 5.2 ] [ 5.2.p1 ] [ 5.2.p2 ] [ 6.0 ] [ 6.0.p1 ] [ 6.1 ] [ 6.2 ] [ 6.2.p1 ] [ 6.2.p2 ] [ 7.0 ] [ 7.0.p1 ] [ 7.1 ] [ 7.1.p1 ] [ 8.0 ] [ 8.0.p1 ] [ 8.1 ] [ 8.1.p1 ] [ 8.1.p2 ] [ 8.2 ] [ 8.2.p1 ] [ 8.3 ] [ 8.3.p1 ] [ 8.3.p2 ] [ 9.0 ] [ 9.0.p1 ] [ 9.0.p2 ] [ 9.1 ] [ 9.1.p1 ] [ 9.1.p2 ] [ 9.1.p3 ] [ 9.2 ] [ 9.2.p1 ] [ 9.2.p2 ] [ 9.2.p3 ] [ 9.2.p4 ] [ 9.3 ] [ 9.3.p1 ] [ 9.3.p2 ] [ 9.4 ] [ 9.4.p1 ] [ 9.4.p2 ] [ 9.4.p3 ] [ 9.4.p4 ] [ 9.5 ] [ 9.5.p1 ] [ 9.5.p2 ] [ 9.6 ] [ 9.6.p1 ] [ 9.6.p2 ] [ 9.6.p3 ] [ 9.6.p4 ] [ 10.0 ] [ 10.0.p1 ] [ 10.0.p2 ] [ 10.0.p3 ] [ 10.0.p4 ] [ 10.1 ] [ 10.1.p1 ] [ 10.1.p2 ] [ 10.1.p3 ] [ 10.2 ] [ 10.2.p1 ] [ 10.2.p2 ] [ 10.2.p3 ] [ 10.3 ] [ 10.3.p1 ] [ 10.3.p2 ] [ 10.3.p3 ] [ 10.4 ] [ 10.4.p1 ] [ 10.4.p2 ] [ 10.4.p3 ] [ 10.5 ] [ 10.5.p1 ] [ 10.6 ] [ 10.6.p1 ] [ 10.6.p2 ] [ 10.6.p3 ] [ 10.7 ] [ 10.7.p1 ] [ 10.7.p2 ] [ 10.7.p3 ] [ 10.7.p4 ] [ 11.0 ] [ 11.0.p1 ] [ 11.0.p2 ] [ 11.0.p3, ] [ 11.0.p4 ] [ 11.1 ] [ 11.1.1 ] [ 11.1.2 ] [ 11.1.3 ] [ 11.2 ] [ 11.2.1 ] [ 11.2.2 ] [ 11.3.0 ]

  1 //
  2 // ********************************************************************
  3 // * License and Disclaimer                                           *
  4 // *                                                                  *
  5 // * The  Geant4 software  is  copyright of the Copyright Holders  of *
  6 // * the Geant4 Collaboration.  It is provided  under  the terms  and *
  7 // * conditions of the Geant4 Software License,  included in the file *
  8 // * LICENSE and available at  http://cern.ch/geant4/license .  These *
  9 // * include a list of copyright holders.                             *
 10 // *                                                                  *
 11 // * Neither the authors of this software system, nor their employing *
 12 // * institutes,nor the agencies providing financial support for this *
 13 // * work  make  any representation or  warranty, express or implied, *
 14 // * regarding  this  software system or assume any liability for its *
 15 // * use.  Please see the license in the file  LICENSE  and URL above *
 16 // * for the full disclaimer and the limitation of liability.         *
 17 // *                                                                  *
 18 // * This  code  implementation is the result of  the  scientific and *
 19 // * technical work of the GEANT4 collaboration.                      *
 20 // * By using,  copying,  modifying or  distributing the software (or *
 21 // * any work based  on the software)  you  agree  to acknowledge its *
 22 // * use  in  resulting  scientific  publications,  and indicate your *
 23 // * acceptance of all terms of the Geant4 Software license.          *
 24 // ********************************************************************
 25 //
 26 //
 27 //
 28 
 29 // ------------------------------------------------------------
 30 //      GEANT 4 class implementation file
 31 //
 32 //      ---------------- G4DiffractiveSplitableHadron----------------
 33 //             by Gunter Folger, August 1998.
 34 //       class splitting an interacting particle. Used by FTF String Model.
 35 // ------------------------------------------------------------
 36 
 37 #include "G4DiffractiveSplitableHadron.hh"
 38 
 39 #include "G4ParticleDefinition.hh"
 40 #include "Randomize.hh"
 41 
 42 
 43 //============================================================================
 44 
 45 G4DiffractiveSplitableHadron::G4DiffractiveSplitableHadron() 
 46 {
 47   PartonIndex = -1;
 48   G4LorentzVector tmp=G4LorentzVector(0.,0.,0.,0.);
 49   Parton[0] = new G4Parton( 1 );
 50   Parton[1] = new G4Parton(-1 );
 51 
 52   Parton[0]->Set4Momentum(tmp); Parton[1]->Set4Momentum(tmp);
 53 }
 54 
 55 
 56 //============================================================================
 57 
 58 G4DiffractiveSplitableHadron::G4DiffractiveSplitableHadron( const G4ReactionProduct& aPrimary ) :
 59   G4VSplitableHadron( aPrimary )
 60 {
 61   PartonIndex = -1;
 62   Parton[0] = nullptr;
 63   Parton[1] = nullptr;
 64 }
 65 
 66 
 67 //============================================================================
 68 
 69 G4DiffractiveSplitableHadron::G4DiffractiveSplitableHadron( const G4Nucleon& aNucleon ) :
 70   G4VSplitableHadron( aNucleon )
 71 {
 72   PartonIndex = -1;
 73   Parton[0] = nullptr;
 74   Parton[1] = nullptr;
 75 }
 76 
 77 
 78 //============================================================================
 79 
 80 G4DiffractiveSplitableHadron::G4DiffractiveSplitableHadron( const G4VKineticNucleon* aNucleon ) :
 81   G4VSplitableHadron( aNucleon )
 82 {
 83   PartonIndex = -1;
 84   Parton[0] = nullptr;
 85   Parton[1] = nullptr;
 86 }
 87 
 88 
 89 //============================================================================
 90 
 91 G4DiffractiveSplitableHadron::~G4DiffractiveSplitableHadron() {}
 92 
 93 
 94 //============================================================================
 95 
 96 void G4DiffractiveSplitableHadron::SplitUp() {
 97 
 98   if ( IsSplit() ) return;
 99   Splitting();
100   // Split once only...
101   if ( Parton[0] != nullptr ) return;
102 
103   // flavours of quark ends
104   G4int PDGcode = GetDefinition()->GetPDGEncoding();
105   G4int stringStart, stringEnd;
106   ChooseStringEnds( PDGcode, &stringStart, &stringEnd );
107 
108   Parton[0] = new G4Parton( stringStart );
109   Parton[1] = new G4Parton( stringEnd );
110 
111   G4LorentzVector tmp=G4LorentzVector(0.,0.,0.,0.);
112   Parton[0]->Set4Momentum(tmp); Parton[1]->Set4Momentum(tmp);
113 
114   /*                                        // Inversion of a string
115   if ( G4UniformRand() < 1.75 ) {  //0.75
116     Parton[0] = new G4Parton( stringStart );
117     Parton[1] = new G4Parton( stringEnd );
118   } else {
119     Parton[0] = new G4Parton( stringEnd );
120     Parton[1] = new G4Parton( stringStart );
121   }
122   */
123 
124   PartonIndex = -1;
125 }
126 
127 
128 //============================================================================
129 
130 G4Parton* G4DiffractiveSplitableHadron::GetNextParton() {
131   ++PartonIndex;
132   if ( PartonIndex > 1  ||  PartonIndex < 0 ) return nullptr;
133   G4int PartonInd( PartonIndex );
134   if ( PartonIndex == 1 ) PartonIndex = -1;
135   return Parton[ PartonInd ];
136 }
137 
138 
139 //============================================================================
140 
141 G4Parton* G4DiffractiveSplitableHadron::GetNextAntiParton() {
142   ++PartonIndex;
143   if ( PartonIndex > 1  ||  PartonIndex < 0 ) return nullptr;
144   G4int PartonInd( PartonIndex );
145   if ( PartonIndex == 1 ) PartonIndex = -1;
146   return Parton[ PartonInd ];
147 }
148 
149 
150 //============================================================================
151 
152 void G4DiffractiveSplitableHadron::SetFirstParton( G4int PDGcode ) {
153   delete Parton[0];
154   Parton[0] = new G4Parton( PDGcode );
155   G4LorentzVector tmp=G4LorentzVector(0.,0.,0.,0.);
156   Parton[0]->Set4Momentum(tmp);
157 }
158 
159 
160 //============================================================================
161 
162 void G4DiffractiveSplitableHadron::SetSecondParton( G4int PDGcode ) {
163   delete Parton[1];
164   Parton[1] = new G4Parton( PDGcode );
165   G4LorentzVector tmp=G4LorentzVector(0.,0.,0.,0.);
166   Parton[1]->Set4Momentum(tmp);
167 }
168 
169 
170 //============================================================================
171 
172 void G4DiffractiveSplitableHadron::ChooseStringEnds( G4int PDGcode, G4int* aEnd,
173                                                      G4int* bEnd ) const {
174   G4int absPDGcode = std::abs( PDGcode );
175 
176   if ( absPDGcode < 1000 ) {  //--------------------  Meson -------------
177     G4int heavy(0), light(0);
178     if (!((absPDGcode == 111)||(absPDGcode == 221)||(absPDGcode == 331)))
179     {                          // Ordinary mesons =======================
180      heavy = absPDGcode/100;
181      light = (absPDGcode % 100)/10;
182      //G4int anti = std::pow( -1 , std::max( heavy, light ) );
183      G4int anti = 1 - 2*( std::max( heavy, light ) % 2 );
184      if (PDGcode < 0 ) anti *= -1;
185      heavy *= anti;
186      light *= -1 * anti;
187     } 
188     else 
189     {                         // Pi0, Eta, Eta' =======================
190      if ( G4UniformRand() < 0.5 ) {heavy = 1; light = -1;}
191      else                         {heavy = 2; light = -2;}
192     }
193     if ( G4UniformRand() < 0.5 ) {
194       *aEnd = heavy;
195       *bEnd = light;
196     } else {
197       *aEnd = light;
198       *bEnd = heavy;
199     }
200   } else {                    //-------------------- Baryon --------------
201     G4int j1000 = PDGcode/1000;
202     G4int j100  = (PDGcode % 1000)/100;
203     G4int j10   = (PDGcode % 100)/10;
204 
205     if ( absPDGcode > 4000 ) {
206       *aEnd = j10;
207       if ( G4UniformRand() > 0.25 ) {
208         *bEnd = Diquark( j1000, j100, 0 );
209       } else {
210         *bEnd = Diquark( j1000, j100, 1 );
211       }
212       return;
213     }
214 
215     G4double SuppresUUDDSS=1.0/2.0;
216     if ((j1000 == j100) && (j1000 == j10)) SuppresUUDDSS=1.; 
217 
218     const G4int maxNumberOfLoops = 1000;
219     G4int loopCounter = 0;
220     do
221     {
222       G4double random = G4UniformRand();
223 
224       if (random < 0.33333)
225       {
226         if (( j100 == j10 ) && ( G4UniformRand() > SuppresUUDDSS )) continue;
227         *aEnd = j1000;
228         if ( j100 == j10 )             {*bEnd = Diquark( j100, j10, 1 );}
229         else
230           if ( G4UniformRand() > 0.25) {*bEnd = Diquark( j100, j10, 0 );}
231           else                        {*bEnd = Diquark( j100, j10, 1 );}
232         break;
233        }
234        else if (random < 0.66667)
235        {
236         if (( j1000 == j10 ) && ( G4UniformRand() > SuppresUUDDSS )) continue;
237         *aEnd = j100;
238         if ( j1000 == j10 )            {*bEnd = Diquark( j1000, j10, 1 );}
239         else
240           if ( G4UniformRand() > 0.25) {*bEnd = Diquark( j1000, j10, 0 );}
241           else                        {*bEnd = Diquark( j1000, j10, 1 );}
242         break;
243        }
244        else
245        {
246         if (( j1000 == j100 ) && ( G4UniformRand() > SuppresUUDDSS )) continue;
247         *aEnd = j10;
248         if ( j1000 == j100 )           {*bEnd = Diquark( j1000, j100, 1 );}
249         else
250           if ( G4UniformRand() > 0.25) {*bEnd = Diquark( j1000, j100, 0 );}
251           else                        {*bEnd = Diquark( j1000, j100, 1 );}
252         break;
253        }
254     } while ( (true) && 
255               ++loopCounter < maxNumberOfLoops );  /* Loop checking, 10.08.2015, A.Ribon */
256     if ( loopCounter >= maxNumberOfLoops ) {
257       *aEnd = j10; *bEnd = Diquark( j1000, j100, 1 );  // Just something acceptable, without any physics consideration.
258     }
259 
260   }
261 }
262 
263 
264 //============================================================================
265 
266 G4int G4DiffractiveSplitableHadron::Diquark( G4int aquark, G4int bquark, G4int Spin) const {
267   G4int diquarkPDG = std::max( std::abs( aquark ), std::abs( bquark ) ) * 1000 +
268                      std::min( std::abs( aquark ), std::abs( bquark ) ) * 100  +
269                      2*Spin + 1;
270   return ( aquark > 0  &&  bquark > 0 ) ? diquarkPDG : -1*diquarkPDG;
271 }
272 
273