Geant4  9.6.p02
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4HadronBuilder.cc
Go to the documentation of this file.
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 // $Id$
28 //
29 // -----------------------------------------------------------------------------
30 // GEANT 4 class implementation file
31 //
32 // History:
33 // Gunter Folger, August/September 2001
34 // Create class; algorithm previously in G4VLongitudinalStringDecay.
35 // -----------------------------------------------------------------------------
36 
37 #include "G4HadronBuilder.hh"
38 #include "G4SystemOfUnits.hh"
39 #include "Randomize.hh"
40 #include "G4HadronicException.hh"
41 #include "G4ParticleTable.hh"
42 
43 G4HadronBuilder::G4HadronBuilder(G4double mesonMix, G4double barionMix,
44  std::vector<double> scalarMesonMix,
45  std::vector<double> vectorMesonMix)
46 {
47  mesonSpinMix=mesonMix;
48  barionSpinMix=barionMix;
49  scalarMesonMixings=scalarMesonMix;
50  vectorMesonMixings=vectorMesonMix;
51 }
52 
53 G4ParticleDefinition * G4HadronBuilder::Build(G4ParticleDefinition * black, G4ParticleDefinition * white)
54 {
55 
56  if (black->GetParticleSubType()== "di_quark" || white->GetParticleSubType()== "di_quark" ) {
57 
58 // Barion
59  Spin spin = (G4UniformRand() < barionSpinMix) ? SpinHalf : SpinThreeHalf;
60  return Barion(black,white,spin);
61 
62  } else {
63 
64 // Meson
65  Spin spin = (G4UniformRand() < mesonSpinMix) ? SpinZero : SpinOne;
66  return Meson(black,white,spin);
67 
68  }
69 }
70 
71 //-------------------------------------------------------------------------
72 
73 G4ParticleDefinition * G4HadronBuilder::BuildLowSpin(G4ParticleDefinition * black, G4ParticleDefinition * white)
74 {
75  if ( black->GetParticleSubType()== "quark" && white->GetParticleSubType()== "quark" ) {
76  return Meson(black,white, SpinZero);
77  } else {
78 // will return a SpinThreeHalf Barion if all quarks the same
79  return Barion(black,white, SpinHalf);
80  }
81 }
82 
83 //-------------------------------------------------------------------------
84 
85 G4ParticleDefinition * G4HadronBuilder::BuildHighSpin(G4ParticleDefinition * black, G4ParticleDefinition * white)
86 {
87  if ( black->GetParticleSubType()== "quark" && white->GetParticleSubType()== "quark" ) {
88  return Meson(black,white, SpinOne);
89  } else {
90  return Barion(black,white,SpinThreeHalf);
91  }
92 }
93 
94 //-------------------------------------------------------------------------
95 
96 G4ParticleDefinition * G4HadronBuilder::Meson(G4ParticleDefinition * black,
97  G4ParticleDefinition * white, Spin theSpin)
98 {
99 #ifdef G4VERBOSE
100 // Verify Input Charge
101 
102  G4double charge = black->GetPDGCharge()
103  + white->GetPDGCharge();
104  if (std::abs(charge) > 2 || std::abs(3.*charge - 3*G4int(charge*1.001)) > perCent ) // 1.001 to avoid int(.9999) -> 0
105  {
106  G4cerr << " G4HadronBuilder::Build()" << G4endl;
107  G4cerr << " Invalid total charge found for on input: "
108  << charge<< G4endl;
109  G4cerr << " PGDcode input quark1/quark2 : " <<
110  black->GetPDGEncoding() << " / "<<
111  white->GetPDGEncoding() << G4endl;
112  G4cerr << G4endl;
113  }
114 #endif
115 
116  G4int id1= black->GetPDGEncoding();
117  G4int id2= white->GetPDGEncoding();
118 // G4int ifl1= std::max(std::abs(id1), std::abs(id2));
119  if ( std::abs(id1) < std::abs(id2) )
120  {
121  G4int xchg = id1;
122  id1 = id2;
123  id2 = xchg;
124  }
125 
126  if (std::abs(id1) > 3 )
127  throw G4HadronicException(__FILE__, __LINE__, "G4HadronBuilder::Meson : Illegal Quark content as input");
128 
129  G4int PDGEncoding=0;
130 
131  if (id1 + id2 == 0) {
132  G4double rmix = G4UniformRand();
133  G4int imix = 2*std::abs(id1) - 1;
134  if(theSpin == SpinZero) {
135  PDGEncoding = 110*(1 + (G4int)(rmix + scalarMesonMixings[imix - 1])
136  + (G4int)(rmix + scalarMesonMixings[imix])
137  ) + theSpin;
138  } else {
139  PDGEncoding = 110*(1 + (G4int)(rmix + vectorMesonMixings[imix - 1])
140  + (G4int)(rmix + vectorMesonMixings[imix])
141  ) + theSpin;
142  }
143  } else {
144  PDGEncoding = 100 * std::abs(id1) + 10 * std::abs(id2) + theSpin;
145  G4bool IsUp = (std::abs(id1)&1) == 0; // quark 1 up type quark (u or c)
146  G4bool IsAnti = id1 < 0; // quark 1 is antiquark?
147  if( (IsUp && IsAnti ) || (!IsUp && !IsAnti ) )
148  PDGEncoding = - PDGEncoding;
149  }
150 
151 
152  G4ParticleDefinition * MesonDef=
153  G4ParticleTable::GetParticleTable()->FindParticle(PDGEncoding);
154 #ifdef G4VERBOSE
155  if (MesonDef == 0 ) {
156  G4cerr << " G4HadronBuilder - Warning: No particle for PDGcode= "
157  << PDGEncoding << G4endl;
158  } else if ( ( black->GetPDGCharge() + white->GetPDGCharge()
159  - MesonDef->GetPDGCharge() ) > perCent ) {
160  G4cerr << " G4HadronBuilder - Warning: Incorrect Charge : "
161  << " Quark1/2 = "
162  << black->GetParticleName() << " / "
163  << white->GetParticleName()
164  << " resulting Hadron " << MesonDef->GetParticleName()
165  << G4endl;
166  }
167 #endif
168 
169  return MesonDef;
170 }
171 
172 
173 G4ParticleDefinition * G4HadronBuilder::Barion(G4ParticleDefinition * black,
174  G4ParticleDefinition * white,Spin theSpin)
175 {
176 
177 #ifdef G4VERBOSE
178 // Verify Input Charge
179  G4double charge = black->GetPDGCharge()
180  + white->GetPDGCharge();
181  if (std::abs(charge) > 2 || std::abs(3.*charge - 3*G4int(charge*1.001)) > perCent )
182  {
183  G4cerr << " G4HadronBuilder::Build()" << G4endl;
184  G4cerr << " Invalid total charge found for on input: "
185  << charge<< G4endl;
186  G4cerr << " PGDcode input quark1/quark2 : " <<
187  black->GetPDGEncoding() << " / "<<
188  white->GetPDGEncoding() << G4endl;
189  G4cerr << G4endl;
190  }
191 #endif
192  G4int id1= black->GetPDGEncoding();
193  G4int id2= white->GetPDGEncoding();
194  if ( std::abs(id1) < std::abs(id2) )
195  {
196  G4int xchg = id1;
197  id1 = id2;
198  id2 = xchg;
199  }
200 
201  if (std::abs(id1) < 1000 || std::abs(id2) > 3 )
202  throw G4HadronicException(__FILE__, __LINE__, "G4HadronBuilder::Barion: Illegal quark content as input");
203 
204  G4int ifl1= std::abs(id1)/1000;
205  G4int ifl2 = (std::abs(id1) - ifl1 * 1000)/100;
206  G4int diquarkSpin = std::abs(id1)%10;
207  G4int ifl3 = id2;
208  if (id1 < 0)
209  {
210  ifl1 = - ifl1;
211  ifl2 = - ifl2;
212  }
213  //... Construct barion, distinguish Lambda and Sigma barions.
214  G4int kfla = std::abs(ifl1);
215  G4int kflb = std::abs(ifl2);
216  G4int kflc = std::abs(ifl3);
217 
218  G4int kfld = std::max(kfla,kflb);
219  kfld = std::max(kfld,kflc);
220  G4int kflf = std::min(kfla,kflb);
221  kflf = std::min(kflf,kflc);
222 
223  G4int kfle = kfla + kflb + kflc - kfld - kflf;
224 
225  //... barion with content uuu or ddd or sss has always spin = 3/2
226  theSpin = (kfla == kflb && kflb == kflc)? SpinThreeHalf : theSpin;
227 
228  G4int kfll = 0;
229  if(theSpin == SpinHalf && kfld > kfle && kfle > kflf) {
230 // Spin J=1/2 and all three quarks different
231 // Two states exist: (uds -> lambda or sigma0)
232 // - lambda: s(ud)0 s : 3122; ie. reverse the two lighter quarks
233 // - sigma0: s(ud)1 s : 3212
234  if(diquarkSpin == 1 ) {
235  if ( kfla == kfld) { // heaviest quark in diquark
236  kfll = 1;
237  } else {
238  kfll = (G4int)(0.25 + G4UniformRand());
239  }
240  }
241  if(diquarkSpin == 3 && kfla != kfld)
242  kfll = (G4int)(0.75 + G4UniformRand());
243  }
244 
245  G4int PDGEncoding;
246  if (kfll == 1)
247  PDGEncoding = 1000 * kfld + 100 * kflf + 10 * kfle + theSpin;
248  else
249  PDGEncoding = 1000 * kfld + 100 * kfle + 10 * kflf + theSpin;
250 
251  if (id1 < 0)
252  PDGEncoding = -PDGEncoding;
253 
254 
255  G4ParticleDefinition * BarionDef=
256  G4ParticleTable::GetParticleTable()->FindParticle(PDGEncoding);
257 #ifdef G4VERBOSE
258  if (BarionDef == 0 ) {
259  G4cerr << " G4HadronBuilder - Warning: No particle for PDGcode= "
260  << PDGEncoding << G4endl;
261  } else if ( ( black->GetPDGCharge() + white->GetPDGCharge()
262  - BarionDef->GetPDGCharge() ) > perCent ) {
263  G4cerr << " G4HadronBuilder - Warning: Incorrect Charge : "
264  << " DiQuark/Quark = "
265  << black->GetParticleName() << " / "
266  << white->GetParticleName()
267  << " resulting Hadron " << BarionDef->GetParticleName()
268  << G4endl;
269  }
270 #endif
271 
272  return BarionDef;
273 }
Generated on Sat May 25 2013 14:34:05 for Geant4 by   doxygen 1.8.4