Geant4 10.7.0
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
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//
28// -----------------------------------------------------------------------------
29// GEANT 4 class implementation file
30//
31// History:
32// Gunter Folger, August/September 2001
33// Create class; algorithm previously in G4VLongitudinalStringDecay.
34// -----------------------------------------------------------------------------
35
36#include "G4HadronBuilder.hh"
37#include "G4SystemOfUnits.hh"
38#include "Randomize.hh"
39#include "G4HadronicException.hh"
40#include "G4ParticleTable.hh"
41
42//#define debug_Hbuilder
43//#define debug_heavyHadrons
44
45G4HadronBuilder::G4HadronBuilder(G4double mesonMix, G4double barionMix,
46 std::vector<double> scalarMesonMix,
47 std::vector<double> vectorMesonMix,
48 G4double Eta_cProb, G4double Eta_bProb)
49{
50 mesonSpinMix = mesonMix;
51 barionSpinMix = barionMix;
52 scalarMesonMixings = scalarMesonMix;
53 vectorMesonMixings = vectorMesonMix;
54
55 ProbEta_c = Eta_cProb;
56 ProbEta_b = Eta_bProb;
57}
58
59G4ParticleDefinition * G4HadronBuilder::Build(G4ParticleDefinition * black, G4ParticleDefinition * white)
60{
61 if (black->GetParticleSubType()== "di_quark" || white->GetParticleSubType()== "di_quark" ) {
62 // Barion
63 Spin spin = (G4UniformRand() < barionSpinMix) ? SpinHalf : SpinThreeHalf;
64 return Barion(black,white,spin);
65 } else {
66 // Meson
67 Spin spin = (G4UniformRand() < mesonSpinMix) ? SpinZero : SpinOne;
68 return Meson(black,white,spin);
69 }
70}
71
72//-------------------------------------------------------------------------
73
74G4ParticleDefinition * G4HadronBuilder::BuildLowSpin(G4ParticleDefinition * black, G4ParticleDefinition * white)
75{
76 if ( black->GetParticleSubType()== "quark" && white->GetParticleSubType()== "quark" ) {
77 return Meson(black,white, SpinZero);
78 } else {
79 // will return a SpinThreeHalf Barion if all quarks the same
80 return Barion(black,white, SpinHalf);
81 }
82}
83
84//-------------------------------------------------------------------------
85
86G4ParticleDefinition * G4HadronBuilder::BuildHighSpin(G4ParticleDefinition * black, G4ParticleDefinition * white)
87{
88 if ( black->GetParticleSubType()== "quark" && white->GetParticleSubType()== "quark" ) {
89 return Meson(black,white, SpinOne);
90 } else {
91 return Barion(black,white,SpinThreeHalf);
92 }
93}
94
95//-------------------------------------------------------------------------
96
97G4ParticleDefinition * G4HadronBuilder::Meson(G4ParticleDefinition * black,
98 G4ParticleDefinition * white, Spin theSpin)
99{
100 #ifdef debug_Hbuilder
101 // Verify Input Charge
102 G4double charge = black->GetPDGCharge() + white->GetPDGCharge();
103 if (std::abs(charge) > 2 || std::abs(3.*charge - 3*G4int(charge*1.001)) > perCent ) // 1.001 to avoid int(.9999) -> 0
104 {
105 G4cerr << " G4HadronBuilder::Build()" << G4endl;
106 G4cerr << " Invalid total charge found for on input: "
107 << charge<< G4endl;
108 G4cerr << " PGDcode input quark1/quark2 : " <<
109 black->GetPDGEncoding() << " / "<<
110 white->GetPDGEncoding() << G4endl;
111 G4cerr << G4endl;
112 }
113 #endif
114
115 G4int id1 = black->GetPDGEncoding();
116 G4int id2 = white->GetPDGEncoding();
117
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 G4int abs_id1 = std::abs(id1);
127
128 if ( abs_id1 > 5 )
129 throw G4HadronicException(__FILE__, __LINE__, "G4HadronBuilder::Meson : Illegal Quark content as input");
130
131 G4int PDGEncoding=0;
132
133 if (id1 + id2 == 0) {
134 if ( abs_id1 < 4) { // light quarks: u, d or s
135 G4double rmix = G4UniformRand();
136 G4int imix = 2*std::abs(id1) - 1;
137 if (theSpin == SpinZero) {
138 PDGEncoding = 110*(1 + (G4int)(rmix + scalarMesonMixings[imix - 1])
139 + (G4int)(rmix + scalarMesonMixings[imix])
140 ) + theSpin;
141 } else {
142 PDGEncoding = 110*(1 + (G4int)(rmix + vectorMesonMixings[imix - 1])
143 + (G4int)(rmix + vectorMesonMixings[imix])
144 ) + theSpin;
145 }
146 } else { // for c and b quarks
147
148 PDGEncoding = abs_id1*100 + abs_id1*10;
149
150 if (PDGEncoding == 440) {
151 if ( G4UniformRand() < ProbEta_c ) {
152 PDGEncoding +=1;
153 } else {
154 PDGEncoding +=3;
155 }
156 }
157 if (PDGEncoding == 550) {
158 if ( G4UniformRand() < ProbEta_b ) {
159 PDGEncoding +=1;
160 } else {
161 PDGEncoding +=3;
162 }
163 }
164 }
165 } else {
166 PDGEncoding = 100 * std::abs(id1) + 10 * std::abs(id2) + theSpin;
167 G4bool IsUp = (std::abs(id1)&1) == 0; // quark 1 up type quark (u or c)
168 G4bool IsAnti = id1 < 0; // quark 1 is antiquark?
169 if ( (IsUp && IsAnti ) || (!IsUp && !IsAnti ) ) PDGEncoding = - PDGEncoding;
170 }
171
172 // ---------------------------------------------------------------------
173 // Special treatment for charmed and bottom mesons : in Geant4 there are
174 // no excited charmed or bottom mesons, therefore we need to transform these
175 // into existing charmed and bottom mesons in Geant4. Whenever possible,
176 // we use the corresponding ground state mesons with the same quantum numbers;
177 // else, we prefer to conserve the electric charge rather than other flavor numbers.
178 #ifdef debug_heavyHadrons
179 G4int initialPDGEncoding = PDGEncoding;
180 #endif
181 if ( std::abs( PDGEncoding ) == 10411 ) // D*0(2400)+ -> D+
182 ( PDGEncoding > 0 ? PDGEncoding = 411 : PDGEncoding = -411 );
183 else if ( std::abs( PDGEncoding ) == 10421 ) // D*0(2400)0 -> D0
184 ( PDGEncoding > 0 ? PDGEncoding = 421 : PDGEncoding = -421 );
185 else if ( std::abs( PDGEncoding ) == 413 ) // D*(2010)+ -> D+
186 ( PDGEncoding > 0 ? PDGEncoding = 411 : PDGEncoding = -411 );
187 else if ( std::abs( PDGEncoding ) == 423 ) // D*(2007)0 -> D0
188 ( PDGEncoding > 0 ? PDGEncoding = 421 : PDGEncoding = -421 );
189 else if ( std::abs( PDGEncoding ) == 10413 ) // D1(2420)+ -> D+
190 ( PDGEncoding > 0 ? PDGEncoding = 411 : PDGEncoding = -411 );
191 else if ( std::abs( PDGEncoding ) == 10423 ) // D1(2420)0 -> D0
192 ( PDGEncoding > 0 ? PDGEncoding = 421 : PDGEncoding = -421 );
193 else if ( std::abs( PDGEncoding ) == 20413 ) // D1(H)+ -> D+
194 ( PDGEncoding > 0 ? PDGEncoding = 411 : PDGEncoding = -411 );
195 else if ( std::abs( PDGEncoding ) == 20423 ) // D1(2430)0 -> D0
196 ( PDGEncoding > 0 ? PDGEncoding = 421 : PDGEncoding = -421 );
197 else if ( std::abs( PDGEncoding ) == 415 ) // D2*(2460)+ -> D+
198 ( PDGEncoding > 0 ? PDGEncoding = 411 : PDGEncoding = -411 );
199 else if ( std::abs( PDGEncoding ) == 425 ) // D2*(2460)0 -> D0
200 ( PDGEncoding > 0 ? PDGEncoding = 421 : PDGEncoding = -421 );
201 else if ( std::abs( PDGEncoding ) == 10431 ) // Ds0*(2317)+ -> Ds+
202 ( PDGEncoding > 0 ? PDGEncoding = 431 : PDGEncoding = -431 );
203 else if ( std::abs( PDGEncoding ) == 433 ) // Ds*+ -> Ds+
204 ( PDGEncoding > 0 ? PDGEncoding = 431 : PDGEncoding = -431 );
205 else if ( std::abs( PDGEncoding ) == 10433 ) // Ds1(2536)+ -> Ds+
206 ( PDGEncoding > 0 ? PDGEncoding = 431 : PDGEncoding = -431 );
207 else if ( std::abs( PDGEncoding ) == 20433 ) // Ds1(2460)+ -> Ds+
208 ( PDGEncoding > 0 ? PDGEncoding = 431 : PDGEncoding = -431 );
209 else if ( std::abs( PDGEncoding ) == 435 ) // Ds2*(2573)+ -> Ds+
210 ( PDGEncoding > 0 ? PDGEncoding = 431 : PDGEncoding = -431 );
211 else if ( std::abs( PDGEncoding ) == 10441 ) PDGEncoding = 441; // chi_c0(1P) -> eta_c
212 else if ( std::abs( PDGEncoding ) == 100441 ) PDGEncoding = 441; // eta_c(2S) -> eta_c
213 else if ( std::abs( PDGEncoding ) == 10443 ) PDGEncoding = 443; // h_c(1P) -> J/psi
214 else if ( std::abs( PDGEncoding ) == 20443 ) PDGEncoding = 443; // chi_c1(1P) -> J/psi
215 else if ( std::abs( PDGEncoding ) == 100443 ) PDGEncoding = 443; // psi(2S) -> J/psi
216 else if ( std::abs( PDGEncoding ) == 30443 ) PDGEncoding = 443; // psi(3770) -> J/psi
217 else if ( std::abs( PDGEncoding ) == 9000443 ) PDGEncoding = 443; // psi(4040) -> J/psi
218 else if ( std::abs( PDGEncoding ) == 9010443 ) PDGEncoding = 443; // psi(4160) -> J/psi
219 else if ( std::abs( PDGEncoding ) == 9020443 ) PDGEncoding = 443; // psi(4415) -> J/psi
220 else if ( std::abs( PDGEncoding ) == 445 ) PDGEncoding = 443; // chi_c2(1P) -> J/psi
221 else if ( std::abs( PDGEncoding ) == 100445 ) PDGEncoding = 443; // chi_c2(2P) -> J/psi
222 // Bottom mesons
223 else if ( std::abs( PDGEncoding ) == 10511 ) // B0*0 -> B0
224 ( PDGEncoding > 0 ? PDGEncoding = 511 : PDGEncoding = -511 );
225 else if ( std::abs( PDGEncoding ) == 10521 ) // B0*+ -> B+
226 ( PDGEncoding > 0 ? PDGEncoding = 521 : PDGEncoding = -521 );
227 else if ( std::abs( PDGEncoding ) == 513 ) // B*0 -> B0
228 ( PDGEncoding > 0 ? PDGEncoding = 511 : PDGEncoding = -511 );
229 else if ( std::abs( PDGEncoding ) == 523 ) // B*+ -> B+
230 ( PDGEncoding > 0 ? PDGEncoding = 521 : PDGEncoding = -521 );
231 else if ( std::abs( PDGEncoding ) == 10513 ) // B1(L)0 -> B0
232 ( PDGEncoding > 0 ? PDGEncoding = 511 : PDGEncoding = -511 );
233 else if ( std::abs( PDGEncoding ) == 10523 ) // B1(L)+ -> B+
234 ( PDGEncoding > 0 ? PDGEncoding = 521 : PDGEncoding = -521 );
235 else if ( std::abs( PDGEncoding ) == 20513 ) // B1(H)0 -> B0
236 ( PDGEncoding > 0 ? PDGEncoding = 511 : PDGEncoding = -511 );
237 else if ( std::abs( PDGEncoding ) == 20523 ) // B1(H)+ -> B+
238 ( PDGEncoding > 0 ? PDGEncoding = 521 : PDGEncoding = -521 );
239 else if ( std::abs( PDGEncoding ) == 515 ) // B2*0 -> B0
240 ( PDGEncoding > 0 ? PDGEncoding = 511 : PDGEncoding = -511 );
241 else if ( std::abs( PDGEncoding ) == 525 ) // B2*+ -> B+
242 ( PDGEncoding > 0 ? PDGEncoding = 521 : PDGEncoding = -521 );
243 else if ( std::abs( PDGEncoding ) == 10531 ) // Bs0*0 -> Bs0
244 ( PDGEncoding > 0 ? PDGEncoding = 531 : PDGEncoding = -531 );
245 else if ( std::abs( PDGEncoding ) == 533 ) // Bs*0 -> Bs0
246 ( PDGEncoding > 0 ? PDGEncoding = 531 : PDGEncoding = -531 );
247 else if ( std::abs( PDGEncoding ) == 10533 ) // Bs1(L)0 -> Bs0
248 ( PDGEncoding > 0 ? PDGEncoding = 531 : PDGEncoding = -531 );
249 else if ( std::abs( PDGEncoding ) == 20533 ) // Bs1(H)0 -> Bs0
250 ( PDGEncoding > 0 ? PDGEncoding = 531 : PDGEncoding = -531 );
251 else if ( std::abs( PDGEncoding ) == 535 ) // Bs2*0 -> Bs0
252 ( PDGEncoding > 0 ? PDGEncoding = 531 : PDGEncoding = -531 );
253 else if ( std::abs( PDGEncoding ) == 10541 ) // Bc0*+ -> Bc+
254 ( PDGEncoding > 0 ? PDGEncoding = 541 : PDGEncoding = -541 );
255 else if ( std::abs( PDGEncoding ) == 543 ) // Bc*+ -> Bc+
256 ( PDGEncoding > 0 ? PDGEncoding = 541 : PDGEncoding = -541 );
257 else if ( std::abs( PDGEncoding ) == 10543 ) // Bc1(L)+ -> Bc+
258 ( PDGEncoding > 0 ? PDGEncoding = 541 : PDGEncoding = -541 );
259 else if ( std::abs( PDGEncoding ) == 20543 ) // Bc1(H)+ -> Bc+
260 ( PDGEncoding > 0 ? PDGEncoding = 541 : PDGEncoding = -541 );
261 else if ( std::abs( PDGEncoding ) == 545 ) // Bc2*+ -> Bc+
262 ( PDGEncoding > 0 ? PDGEncoding = 541 : PDGEncoding = -541 );
263 else if ( std::abs( PDGEncoding ) == 551 ) PDGEncoding = 553; // eta_b(1S) -> Upsilon
264 else if ( std::abs( PDGEncoding ) == 10551 ) PDGEncoding = 553; // chi_b0(1P) -> Upsilon
265 else if ( std::abs( PDGEncoding ) == 100551 ) PDGEncoding = 553; // eta_b(2S) -> Upsilon
266 else if ( std::abs( PDGEncoding ) == 110551 ) PDGEncoding = 553; // chi_b0(2P) -> Upsilon
267 else if ( std::abs( PDGEncoding ) == 200551 ) PDGEncoding = 553; // eta_b(3S) -> Upsilon
268 else if ( std::abs( PDGEncoding ) == 210551 ) PDGEncoding = 553; // chi_b0(3P) -> Upsilon
269 else if ( std::abs( PDGEncoding ) == 10553 ) PDGEncoding = 553; // h_b(1P) -> Upsilon
270 else if ( std::abs( PDGEncoding ) == 20553 ) PDGEncoding = 553; // chi_b1(1P) -> Upsilon
271 else if ( std::abs( PDGEncoding ) == 30553 ) PDGEncoding = 553; // Upsilon_1(1D) -> Upsilon
272 else if ( std::abs( PDGEncoding ) == 100553 ) PDGEncoding = 553; // Upsilon(2S) -> Upsilon
273 else if ( std::abs( PDGEncoding ) == 110553 ) PDGEncoding = 553; // h_b(2P) -> Upsilon
274 else if ( std::abs( PDGEncoding ) == 120553 ) PDGEncoding = 553; // chi_b1(2P) -> Upsilon
275 else if ( std::abs( PDGEncoding ) == 130553 ) PDGEncoding = 553; // Upsilon_1(2D) -> Upsilon
276 else if ( std::abs( PDGEncoding ) == 200553 ) PDGEncoding = 553; // Upsilon(3S) -> Upsilon
277 else if ( std::abs( PDGEncoding ) == 210553 ) PDGEncoding = 553; // h_b(3P) -> Upsilon
278 else if ( std::abs( PDGEncoding ) == 220553 ) PDGEncoding = 553; // chi_b1(3P) -> Upsilon
279 else if ( std::abs( PDGEncoding ) == 300553 ) PDGEncoding = 553; // Upsilon(4S) -> Upsilon
280 else if ( std::abs( PDGEncoding ) == 9000553 ) PDGEncoding = 553; // Upsilon(10860) -> Upsilon
281 else if ( std::abs( PDGEncoding ) == 9010553 ) PDGEncoding = 553; // Upsilon(11020) -> Upsilon
282 else if ( std::abs( PDGEncoding ) == 555 ) PDGEncoding = 553; // chi_b2(1P) -> Upsilon
283 else if ( std::abs( PDGEncoding ) == 10555 ) PDGEncoding = 553; // eta_b2(1D) -> Upsilon
284 else if ( std::abs( PDGEncoding ) == 20555 ) PDGEncoding = 553; // Upsilon_2(1D) -> Upsilon
285 else if ( std::abs( PDGEncoding ) == 100555 ) PDGEncoding = 553; // chi_b2(2P) -> Upsilon
286 else if ( std::abs( PDGEncoding ) == 110555 ) PDGEncoding = 553; // eta_b2(2D) -> Upsilon
287 else if ( std::abs( PDGEncoding ) == 120555 ) PDGEncoding = 553; // Upsilon_2(2D) -> Upsilon
288 else if ( std::abs( PDGEncoding ) == 200555 ) PDGEncoding = 553; // chi_b2(3P) -> Upsilon
289 else if ( std::abs( PDGEncoding ) == 557 ) PDGEncoding = 553; // Upsilon_3(1D) -> Upsilon
290 else if ( std::abs( PDGEncoding ) == 100557 ) PDGEncoding = 553; // Upsilon_3(2D) -> Upsilon
291 #ifdef debug_heavyHadrons
292 if ( initialPDGEncoding != PDGEncoding ) {
293 G4cout << "G4HadronBuilder::Meson : forcing (inexisting in G4) heavy meson with pdgCode="
294 << initialPDGEncoding << " into pdgCode=" << PDGEncoding << G4endl;
295 }
296 #endif
297 // ---------------------------------------------------------------------
298
299 G4ParticleDefinition * MesonDef=
300 G4ParticleTable::GetParticleTable()->FindParticle(PDGEncoding);
301
302 #ifdef debug_Hbuilder
303 if (MesonDef == 0 ) {
304 G4cerr << " G4HadronBuilder - Warning: No particle for PDGcode= "
305 << PDGEncoding << G4endl;
306
307 } else if ( ( black->GetPDGCharge() + white->GetPDGCharge()
308 - MesonDef->GetPDGCharge() ) > perCent ) {
309 G4cerr << " G4HadronBuilder - Warning: Incorrect Charge : "
310 << " Quark1/2 = "
311 << black->GetParticleName() << " / "
312 << white->GetParticleName()
313 << " resulting Hadron " << MesonDef->GetParticleName()
314 << G4endl;
315 }
316 #endif
317
318 return MesonDef;
319}
320
321
322G4ParticleDefinition * G4HadronBuilder::Barion(G4ParticleDefinition * black,
323 G4ParticleDefinition * white,Spin theSpin)
324{
325 #ifdef debug_Hbuilder
326 // Verify Input Charge
327 G4double charge = black->GetPDGCharge() + white->GetPDGCharge();
328 if (std::abs(charge) > 2 || std::abs(3.*charge - 3*G4int(charge*1.001)) > perCent )
329 {
330 G4cerr << " G4HadronBuilder::Build()" << G4endl;
331 G4cerr << " Invalid total charge found for on input: "
332 << charge<< G4endl;
333 G4cerr << " PGDcode input quark1/quark2 : " <<
334 black->GetPDGEncoding() << " / "<<
335 white->GetPDGEncoding() << G4endl;
336 G4cerr << G4endl;
337 }
338 #endif
339
340 G4int id1 = black->GetPDGEncoding();
341 G4int id2 = white->GetPDGEncoding();
342
343 if ( std::abs(id1) < std::abs(id2) )
344 {
345 G4int xchg = id1;
346 id1 = id2;
347 id2 = xchg;
348 }
349
350 if (std::abs(id1) < 1000 || std::abs(id2) > 5 )
351 throw G4HadronicException(__FILE__, __LINE__, "G4HadronBuilder::Barion: Illegal quark content as input");
352
353 G4int ifl1= std::abs(id1)/1000;
354 G4int ifl2 = (std::abs(id1) - ifl1 * 1000)/100;
355 G4int diquarkSpin = std::abs(id1)%10;
356 G4int ifl3 = id2;
357 if (id1 < 0)
358 {
359 ifl1 = - ifl1;
360 ifl2 = - ifl2;
361 }
362 //... Construct barion, distinguish Lambda and Sigma barions.
363 G4int kfla = std::abs(ifl1);
364 G4int kflb = std::abs(ifl2);
365 G4int kflc = std::abs(ifl3);
366
367 G4int kfld = std::max(kfla,kflb);
368 kfld = std::max(kfld,kflc);
369 G4int kflf = std::min(kfla,kflb);
370 kflf = std::min(kflf,kflc);
371
372 G4int kfle = kfla + kflb + kflc - kfld - kflf;
373
374 //... barion with content uuu or ddd or sss has always spin = 3/2
375 theSpin = (kfla == kflb && kflb == kflc)? SpinThreeHalf : theSpin;
376
377 G4int kfll = 0;
378 if (kfld < 6) {
379 if (theSpin == SpinHalf && kfld > kfle && kfle > kflf) {
380 // Spin J=1/2 and all three quarks different
381 // Two states exist: (uds -> lambda or sigma0)
382 // - lambda: s(ud)0 s : 3122; ie. reverse the two lighter quarks
383 // - sigma0: s(ud)1 s : 3212
384 if (diquarkSpin == 1 ) {
385 if ( kfla == kfld) { // heaviest quark in diquark
386 kfll = 1;
387 } else {
388 kfll = (G4int)(0.25 + G4UniformRand());
389 }
390 }
391 if (diquarkSpin == 3 && kfla != kfld)
392 kfll = (G4int)(0.75 + G4UniformRand());
393 }
394 }
395
396 G4int PDGEncoding;
397 if (kfll == 1)
398 PDGEncoding = 1000 * kfld + 100 * kflf + 10 * kfle + theSpin;
399 else
400 PDGEncoding = 1000 * kfld + 100 * kfle + 10 * kflf + theSpin;
401
402 if (id1 < 0)
403 PDGEncoding = -PDGEncoding;
404
405 // ---------------------------------------------------------------------
406 // Special treatment for charmed and bottom baryons : in Geant4 there are
407 // neither excited charmed or bottom baryons, nor baryons with two or three
408 // heavy (c, b) constitutent quarks:
409 // Sigma_c* , Xi_c' , Xi_c* , Omega_c* ,
410 // Xi_cc , Xi_cc* , Omega_cc , Omega_cc* , Omega_ccc ;
411 // Sigma_b* , Xi_b' , Xi_b* , Omega_b*,
412 // Xi_bc , Xi_bc' , Xi_bc* , Omega_bc , Omega_bc' , Omega_bc* ,
413 // Omega_bcc , Omega_bcc* , Xi_bb, Xi_bb* , Omega_bb, Omega_bb* ,
414 // Omega_bbc , Omega_bbc* , Omega_bbb
415 // therefore we need to transform these into existing charmed and bottom
416 // baryons in Geant4. Whenever possible, we use the corresponding ground state
417 // baryons with the same quantum numbers; else, we prefer to conserve the
418 // electric charge rather than other flavor numbers.
419 #ifdef debug_heavyHadrons
420 G4int charmViolation = 0, bottomViolation = 0; // Only positive
421 G4int initialPDGEncoding = PDGEncoding;
422 #endif
423 if ( std::abs( PDGEncoding ) == 4224 ) { // Sigma_c*++ -> Sigma_c++
424 ( PDGEncoding > 0 ? PDGEncoding = 4222 : PDGEncoding = -4222 );
425 } else if ( std::abs( PDGEncoding ) == 4214 ) { // Sigma_c*+ -> Sigma_c+
426 ( PDGEncoding > 0 ? PDGEncoding = 4212 : PDGEncoding = -4212 );
427 } else if ( std::abs( PDGEncoding ) == 4114 ) { // Sigma_c*0 -> Sigma_c0
428 ( PDGEncoding > 0 ? PDGEncoding = 4112 : PDGEncoding = -4112 );
429 } else if ( std::abs( PDGEncoding ) == 4322 ) { // Xi_c'+ -> Xi_c+
430 ( PDGEncoding > 0 ? PDGEncoding = 4232 : PDGEncoding = -4232 );
431 } else if ( std::abs( PDGEncoding ) == 4312 ) { // Xi_c'0 -> Xi_c0
432 ( PDGEncoding > 0 ? PDGEncoding = 4132 : PDGEncoding = -4132 );
433 } else if ( std::abs( PDGEncoding ) == 4324 ) { // Xi_c*+ -> Xi_c+
434 ( PDGEncoding > 0 ? PDGEncoding = 4232 : PDGEncoding = -4232 );
435 } else if ( std::abs( PDGEncoding ) == 4314 ) { // Xi_c*0 -> Xi_c0
436 ( PDGEncoding > 0 ? PDGEncoding = 4132 : PDGEncoding = -4132 );
437 } else if ( std::abs( PDGEncoding ) == 4334 ) { // Omega_c*0 -> Omega_c0
438 ( PDGEncoding > 0 ? PDGEncoding = 4332 : PDGEncoding = -4332 );
439 } else if ( std::abs( PDGEncoding ) == 4412 ) { // Xi_cc+ -> Xi_c+
440 ( PDGEncoding > 0 ? PDGEncoding = 4232 : PDGEncoding = -4232 );
441 #ifdef debug_heavyHadrons
442 charmViolation = 1;
443 #endif
444 } else if ( std::abs( PDGEncoding ) == 4422 ) { // Xi_cc++ -> Sigma_c++ (use Sigma to conserve charge)
445 ( PDGEncoding > 0 ? PDGEncoding = 4222 : PDGEncoding = -4222 );
446 #ifdef debug_heavyHadrons
447 charmViolation = 1;
448 #endif
449 } else if ( std::abs( PDGEncoding ) == 4414 ) { // Xi_cc*+ -> Xi_c+
450 ( PDGEncoding > 0 ? PDGEncoding = 4232 : PDGEncoding = -4232 );
451 #ifdef debug_heavyHadrons
452 charmViolation = 1;
453 #endif
454 } else if ( std::abs( PDGEncoding ) == 4424 ) { // Xi_cc*++ -> Sigma_c++ (use Sigma to conserve charge)
455 ( PDGEncoding > 0 ? PDGEncoding = 4222 : PDGEncoding = -4222 );
456 #ifdef debug_heavyHadrons
457 charmViolation = 1;
458 #endif
459 } else if ( std::abs( PDGEncoding ) == 4432 ) { // Omega_cc+ -> Xi_c+ (use Xi to conserve charge)
460 ( PDGEncoding > 0 ? PDGEncoding = 4232 : PDGEncoding = -4232 );
461 #ifdef debug_heavyHadrons
462 charmViolation = 1;
463 #endif
464 } else if ( std::abs( PDGEncoding ) == 4434 ) { // Omega_cc*+ -> Xi_c+ (use Xi to conserve charge)
465 ( PDGEncoding > 0 ? PDGEncoding = 4232 : PDGEncoding = -4232 );
466 #ifdef debug_heavyHadrons
467 charmViolation = 1;
468 #endif
469 } else if ( std::abs( PDGEncoding ) == 4444 ) { // Omega_ccc++ -> Sigma_c++ (use Sigma to conserve charge)
470 ( PDGEncoding > 0 ? PDGEncoding = 4222 : PDGEncoding = -4222 );
471 #ifdef debug_heavyHadrons
472 charmViolation = 2;
473 #endif
474 // Bottom baryons
475 } else if ( std::abs( PDGEncoding ) == 5114 ) { // Sigma_b*- -> Sigma_b-
476 ( PDGEncoding > 0 ? PDGEncoding = 5112 : PDGEncoding = -5112 );
477 } else if ( std::abs( PDGEncoding ) == 5214 ) { // Sigma_b*0 -> Sigma_b0
478 ( PDGEncoding > 0 ? PDGEncoding = 5212 : PDGEncoding = -5212 );
479 } else if ( std::abs( PDGEncoding ) == 5224 ) { // Sigma_b*+ -> Sigma_b+
480 ( PDGEncoding > 0 ? PDGEncoding = 5222 : PDGEncoding = -5222 );
481 } else if ( std::abs( PDGEncoding ) == 5312 ) { // Xi_b'- -> Xi_b-
482 ( PDGEncoding > 0 ? PDGEncoding = 5132 : PDGEncoding = -5132 );
483 } else if ( std::abs( PDGEncoding ) == 5322 ) { // Xi_b'0 -> Xi_b0
484 ( PDGEncoding > 0 ? PDGEncoding = 5232 : PDGEncoding = -5232 );
485 } else if ( std::abs( PDGEncoding ) == 5314 ) { // Xi_b*- -> Xi_b-
486 ( PDGEncoding > 0 ? PDGEncoding = 5132 : PDGEncoding = -5132 );
487 } else if ( std::abs( PDGEncoding ) == 5324 ) { // Xi_b*0 -> Xi_b0
488 ( PDGEncoding > 0 ? PDGEncoding = 5232 : PDGEncoding = -5232 );
489 } else if ( std::abs( PDGEncoding ) == 5334 ) { // Omega_b*- -> Omega_b-
490 ( PDGEncoding > 0 ? PDGEncoding = 5332 : PDGEncoding = -5332 );
491 } else if ( std::abs( PDGEncoding ) == 5142 ) { // Xi_bc0 -> Xi_b0
492 ( PDGEncoding > 0 ? PDGEncoding = 5232 : PDGEncoding = -5232 );
493 #ifdef debug_heavyHadrons
494 charmViolation = 1;
495 #endif
496 } else if ( std::abs( PDGEncoding ) == 5242 ) { // Xi_bc+ -> Sigma_b+ (use Sigma to conserve charge)
497 ( PDGEncoding > 0 ? PDGEncoding = 5222 : PDGEncoding = -5222 );
498 #ifdef debug_heavyHadrons
499 charmViolation = 1;
500 #endif
501 } else if ( std::abs( PDGEncoding ) == 5412 ) { // Xi_bc'0 -> Xi_b0
502 ( PDGEncoding > 0 ? PDGEncoding = 5232 : PDGEncoding = -5232 );
503 #ifdef debug_heavyHadrons
504 charmViolation = 1;
505 #endif
506 } else if ( std::abs( PDGEncoding ) == 5422 ) { // Xi_bc'+ -> Sigma_b+ (use Sigma to conserve charge)
507 ( PDGEncoding > 0 ? PDGEncoding = 5222 : PDGEncoding = -5222 );
508 #ifdef debug_heavyHadrons
509 charmViolation = 1;
510 #endif
511 } else if ( std::abs( PDGEncoding ) == 5414 ) { // Xi_bc*0 -> Xi_b0
512 ( PDGEncoding > 0 ? PDGEncoding = 5232 : PDGEncoding = -5232 );
513 #ifdef debug_heavyHadrons
514 charmViolation = 1;
515 #endif
516 } else if ( std::abs( PDGEncoding ) == 5424 ) { // Xi_bc*+ -> Sigma_b+ (use Sigma to conserve charge)
517 ( PDGEncoding > 0 ? PDGEncoding = 5222 : PDGEncoding = -5222 );
518 #ifdef debug_heavyHadrons
519 charmViolation = 1;
520 #endif
521 } else if ( std::abs( PDGEncoding ) == 5342 ) { // Omega_bc0 -> Xi_b0 (use Xi to conserve charge)
522 ( PDGEncoding > 0 ? PDGEncoding = 5232 : PDGEncoding = -5232 );
523 #ifdef debug_heavyHadrons
524 charmViolation = 1;
525 #endif
526 } else if ( std::abs( PDGEncoding ) == 5432 ) { // Omega_bc'0 -> Xi_b0 (use Xi to conserve charge)
527 ( PDGEncoding > 0 ? PDGEncoding = 5232 : PDGEncoding = -5232 );
528 #ifdef debug_heavyHadrons
529 charmViolation = 1;
530 #endif
531 } else if ( std::abs( PDGEncoding ) == 5434 ) { // Omega_bc*0 -> Xi_b0 (use Xi to conserve charge)
532 ( PDGEncoding > 0 ? PDGEncoding = 5232 : PDGEncoding = -5232 );
533 #ifdef debug_heavyHadrons
534 charmViolation = 1;
535 #endif
536 } else if ( std::abs( PDGEncoding ) == 5442 ) { // Omega_bcc+ -> Sigma_b+ (use Sigma to conserve charge)
537 ( PDGEncoding > 0 ? PDGEncoding = 5222 : PDGEncoding = -5222 );
538 #ifdef debug_heavyHadrons
539 charmViolation = 2;
540 #endif
541 } else if ( std::abs( PDGEncoding ) == 5444 ) { // Omega_bcc*+ -> Sigma_b+ (use Sigma to conserve charge)
542 ( PDGEncoding > 0 ? PDGEncoding = 5222 : PDGEncoding = -5222 );
543 #ifdef debug_heavyHadrons
544 charmViolation = 2;
545 #endif
546 } else if ( std::abs( PDGEncoding ) == 5512 ) { // Xi_bb- -> Xi_b-
547 ( PDGEncoding > 0 ? PDGEncoding = 5132 : PDGEncoding = -5132 );
548 #ifdef debug_heavyHadrons
549 bottomViolation = 1;
550 #endif
551 } else if ( std::abs( PDGEncoding ) == 5522 ) { // Xi_bb0 -> Xi_b0
552 ( PDGEncoding > 0 ? PDGEncoding = 5232 : PDGEncoding = -5232 );
553 #ifdef debug_heavyHadrons
554 bottomViolation = 1;
555 #endif
556 } else if ( std::abs( PDGEncoding ) == 5514 ) { // Xi_bb*- -> Xi_b-
557 ( PDGEncoding > 0 ? PDGEncoding = 5132 : PDGEncoding = -5132 );
558 #ifdef debug_heavyHadrons
559 bottomViolation = 1;
560 #endif
561 } else if ( std::abs( PDGEncoding ) == 5524 ) { // Xi_bb*0 -> Xi_b0
562 ( PDGEncoding > 0 ? PDGEncoding = 5232 : PDGEncoding = -5232 );
563 #ifdef debug_heavyHadrons
564 bottomViolation = 1;
565 #endif
566 } else if ( std::abs( PDGEncoding ) == 5532 ) { // Omega_bb- -> Omega_b-
567 ( PDGEncoding > 0 ? PDGEncoding = 5332 : PDGEncoding = -5332 );
568 #ifdef debug_heavyHadrons
569 bottomViolation = 1;
570 #endif
571 } else if ( std::abs( PDGEncoding ) == 5534 ) { // Omega_bb*- -> Omega_b-
572 ( PDGEncoding > 0 ? PDGEncoding = 5332 : PDGEncoding = -5332 );
573 #ifdef debug_heavyHadrons
574 bottomViolation = 1;
575 #endif
576 } else if ( std::abs( PDGEncoding ) == 5542 ) { // Omega_bbc0 -> Xi_b0 (use Xi to conserve charge)
577 ( PDGEncoding > 0 ? PDGEncoding = 5232 : PDGEncoding = -5232 );
578 #ifdef debug_heavyHadrons
579 charmViolation = 1; bottomViolation = 1;
580 #endif
581 } else if ( std::abs( PDGEncoding ) == 5544 ) { // Omega_bbc*0 -> Xi_b0 (use Xi to conserve charge)
582 ( PDGEncoding > 0 ? PDGEncoding = 5232 : PDGEncoding = -5232 );
583 #ifdef debug_heavyHadrons
584 charmViolation = 1; bottomViolation = 1;
585 #endif
586 } else if ( std::abs( PDGEncoding ) == 5554 ) { // Omega_bbb- -> Omega_b-
587 ( PDGEncoding > 0 ? PDGEncoding = 5332 : PDGEncoding = -5332 );
588 #ifdef debug_heavyHadrons
589 bottomViolation = 2;
590 #endif
591 }
592 #ifdef debug_heavyHadrons
593 if ( initialPDGEncoding != PDGEncoding ) {
594 G4cout << "G4HadronBuilder::Barion : forcing (inexisting in G4) heavy baryon with pdgCode="
595 << initialPDGEncoding << " into pdgCode=" << PDGEncoding << G4endl;
596 if ( charmViolation != 0 || bottomViolation != 0 ) {
597 G4cout << "\t --> VIOLATION of " << ( charmViolation != 0 ? " CHARM " : " " )
598 << ( charmViolation != 0 && bottomViolation != 0 ? " and " : " " )
599 << ( bottomViolation != 0 ? " BOTTOM " : " " ) << " quantum number ! " << G4endl;
600 }
601 }
602 #endif
603 // ---------------------------------------------------------------------
604
605 G4ParticleDefinition * BarionDef=
606 G4ParticleTable::GetParticleTable()->FindParticle(PDGEncoding);
607
608 #ifdef debug_Hbuilder
609 if (BarionDef == 0 ) {
610 G4cerr << " G4HadronBuilder - Warning: No particle for PDGcode= "
611 << PDGEncoding << G4endl;
612 } else if ( ( black->GetPDGCharge() + white->GetPDGCharge()
613 - BarionDef->GetPDGCharge() ) > perCent ) {
614 G4cerr << " G4HadronBuilder - Warning: Incorrect Charge : "
615 << " DiQuark/Quark = "
616 << black->GetParticleName() << " / "
617 << white->GetParticleName()
618 << " resulting Hadron " << BarionDef->GetParticleName()
619 << G4endl;
620 }
621 #endif
622
623 return BarionDef;
624}
625
G4double
double G4double
Definition: G4Types.hh:83
G4bool
bool G4bool
Definition: G4Types.hh:86
G4int
int G4int
Definition: G4Types.hh:85
G4cerr
G4GLOB_DLL std::ostream G4cerr
G4endl
#define G4endl
Definition: G4ios.hh:57
G4cout
G4GLOB_DLL std::ostream G4cout
G4UniformRand
#define G4UniformRand()
Definition: Randomize.hh:52
G4HadronBuilder::Build
G4ParticleDefinition * Build(G4ParticleDefinition *black, G4ParticleDefinition *white)
Definition: G4HadronBuilder.cc:59
G4HadronBuilder::BuildLowSpin
G4ParticleDefinition * BuildLowSpin(G4ParticleDefinition *black, G4ParticleDefinition *white)
Definition: G4HadronBuilder.cc:74
G4HadronBuilder::BuildHighSpin
G4ParticleDefinition * BuildHighSpin(G4ParticleDefinition *black, G4ParticleDefinition *white)
Definition: G4HadronBuilder.cc:86
G4HadronBuilder::G4HadronBuilder
G4HadronBuilder(G4double mesonMix, G4double barionMix, std::vector< double > scalarMesonMix, std::vector< double > vectorMesonMix, G4double Eta_cProb, G4double Eta_bProb)
Definition: G4HadronBuilder.cc:45
G4ParticleDefinition::GetPDGCharge
G4double GetPDGCharge() const
Definition: G4ParticleDefinition.hh:112
G4ParticleDefinition::GetParticleName
const G4String & GetParticleName() const
Definition: G4ParticleDefinition.hh:108
G4ParticleDefinition::GetParticleSubType
const G4String & GetParticleSubType() const
Definition: G4ParticleDefinition.hh:130
G4ParticleTable::FindParticle
G4ParticleDefinition * FindParticle(G4int PDGEncoding)
Definition: G4ParticleTable.cc:586
G4ParticleTable::GetParticleTable
static G4ParticleTable * GetParticleTable()
Definition: G4ParticleTable.cc:87