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