Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/processes/hadronic/models/parton_string/hadronization/src/G4VLongitudinalStringDecay.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 ]

Diff markup

Differences between /processes/hadronic/models/parton_string/hadronization/src/G4VLongitudinalStringDecay.cc (Version 11.3.0) and /processes/hadronic/models/parton_string/hadronization/src/G4VLongitudinalStringDecay.cc (Version 10.4)


  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 // $Id: G4VLongitudinalStringDecay.cc 106967 2017-10-31 08:41:49Z gcosmo $
 27 //                                                 28 //
 28 // -------------------------------------------     29 // -----------------------------------------------------------------------------
 29 //      GEANT 4 class implementation file          30 //      GEANT 4 class implementation file
 30 //                                                 31 //
 31 //      History: first implementation, Maxim K     32 //      History: first implementation, Maxim Komogorov, 1-Jul-1998
 32 //               redesign  Gunter Folger, Augu     33 //               redesign  Gunter Folger, August/September 2001
 33 // -------------------------------------------     34 // -----------------------------------------------------------------------------
 34 #include "G4VLongitudinalStringDecay.hh"           35 #include "G4VLongitudinalStringDecay.hh"
 35 #include "G4PhysicalConstants.hh"                  36 #include "G4PhysicalConstants.hh"
 36 #include "G4SystemOfUnits.hh"                      37 #include "G4SystemOfUnits.hh"
 37 #include "G4ios.hh"                                38 #include "G4ios.hh"
 38 #include "Randomize.hh"                            39 #include "Randomize.hh"
 39 #include "G4FragmentingString.hh"                  40 #include "G4FragmentingString.hh"
 40                                                    41 
 41 #include "G4ParticleDefinition.hh"                 42 #include "G4ParticleDefinition.hh"
 42 #include "G4ParticleTypes.hh"                      43 #include "G4ParticleTypes.hh"
 43 #include "G4ParticleChange.hh"                     44 #include "G4ParticleChange.hh"
 44 #include "G4VShortLivedParticle.hh"                45 #include "G4VShortLivedParticle.hh"
 45 #include "G4ShortLivedConstructor.hh"              46 #include "G4ShortLivedConstructor.hh"
 46 #include "G4ParticleTable.hh"                      47 #include "G4ParticleTable.hh"
 47 #include "G4PhaseSpaceDecayChannel.hh"             48 #include "G4PhaseSpaceDecayChannel.hh"
 48 #include "G4VDecayChannel.hh"                      49 #include "G4VDecayChannel.hh"
 49 #include "G4DecayTable.hh"                         50 #include "G4DecayTable.hh"
 50                                                    51 
 51 #include "G4DiQuarks.hh"                           52 #include "G4DiQuarks.hh"
 52 #include "G4Quarks.hh"                             53 #include "G4Quarks.hh"
 53 #include "G4Gluons.hh"                             54 #include "G4Gluons.hh"
 54                                                    55 
 55 #include "G4Exp.hh"                                56 #include "G4Exp.hh"
 56 #include "G4Log.hh"                                57 #include "G4Log.hh"
 57                                                    58 
 58 #include "G4HadronicException.hh"              << 
 59                                                << 
 60 //------------------------debug switches           59 //------------------------debug switches
 61 //#define debug_VStringDecay                       60 //#define debug_VStringDecay
 62 //#define debug_heavyHadrons                   << 
 63                                                    61 
 64 //******************************************** <<  62 //********************************************************************************
 65 // Constructors                                    63 // Constructors
 66                                                    64 
 67 G4VLongitudinalStringDecay::G4VLongitudinalStr <<  65 G4VLongitudinalStringDecay::G4VLongitudinalStringDecay()
 68   : G4HadronicInteraction(name), ProbCCbar(0.0 << 
 69 {                                                  66 {
 70    MassCut = 210.0*MeV;   // Mpi + Delta       <<  67   MassCut  = 0.35*GeV; 
                                                   >>  68   ClusterMass = 0.15*GeV;
 71                                                    69 
 72    StringLoopInterrupt  = 1000;                <<  70   SmoothParam      = 0.9; 
 73    ClusterLoopInterrupt =  500;                <<  71   StringLoopInterrupt    = 1000;
                                                   >>  72   ClusterLoopInterrupt   =  500;
 74                                                    73 
 75    // Changable Parameters below.              <<  74   // Changable Parameters below.
 76    SigmaQT = 0.5 * GeV;                        <<  75   SigmaQT = 0.5 * GeV;  // 0.5 0.1
 77                                                    76    
 78    StrangeSuppress  = 0.44;    // =0.27/2.27 s <<  77   StrangeSuppress  = 0.44;    //  27 % strange quarks produced, ie. u:d:s=1:1:0.27
 79    DiquarkSuppress  = 0.07;    // Probability  <<  78   DiquarkSuppress  = 0.07;
 80    DiquarkBreakProb = 0.1;     // Probability  <<  79   DiquarkBreakProb = 0.1;
 81                                                    80    
 82    //... pspin_meson is probability to create  <<  81   //... pspin_meson is probability to create pseudo-scalar meson 
 83    pspin_meson.resize(3);                      <<  82   pspin_meson = 0.5;
 84    pspin_meson[0] = 0.5;  // u or d + anti-u o << 
 85    pspin_meson[1] = 0.4;  // one of the quark  << 
 86    pspin_meson[2] = 0.3;  // both of the quark << 
 87                                                << 
 88    //... pspin_barion is probability to create << 
 89    pspin_barion = 0.5;                         << 
 90                                                    83 
 91    //... vectorMesonMix[] is quark mixing para <<  84   //... pspin_barion is probability to create 1/2 barion 
 92    vectorMesonMix.resize(6);                   <<  85   pspin_barion = 0.5;
 93    vectorMesonMix[0] = 0.0;                    << 
 94    vectorMesonMix[1] = 0.5;                    << 
 95    vectorMesonMix[2] = 0.0;                    << 
 96    vectorMesonMix[3] = 0.5;                    << 
 97    vectorMesonMix[4] = 1.0;                    << 
 98    vectorMesonMix[5] = 1.0;                    << 
 99                                                << 
100    //... scalarMesonMix[] is quark mixing para << 
101    scalarMesonMix.resize(6);                   << 
102    scalarMesonMix[0] = 0.5;                    << 
103    scalarMesonMix[1] = 0.25;                   << 
104    scalarMesonMix[2] = 0.5;                    << 
105    scalarMesonMix[3] = 0.25;                   << 
106    scalarMesonMix[4] = 1.0;                    << 
107    scalarMesonMix[5] = 0.5;                    << 
108                                                << 
109    SetProbCCbar(0.0);  // Probability of CCbar << 
110    SetProbEta_c(0.1);  // Mixing of Eta_c and  << 
111    SetProbBBbar(0.0);  // Probability of BBbar << 
112    SetProbEta_b(0.0);  // Mixing of Eta_b and  << 
113                                                << 
114    // Parameters may be changed until the firs << 
115    PastInitPhase=false;                        << 
116    hadronizer = new G4HadronBuilder( pspin_mes << 
117                                      ProbEta_c << 
118                                                    86 
119    MaxMass=-350.0*GeV;  // If there will be a  <<  87   //... vectorMesonMix[] is quark mixing parameters for vector mesons (Variable spin = 3)
                                                   >>  88   vectorMesonMix.resize(6);
                                                   >>  89   vectorMesonMix[0] = 0.0;  //AR-20Oct2014 : it was 0.5
                                                   >>  90   vectorMesonMix[1] = 0.0;
                                                   >>  91   vectorMesonMix[2] = 0.0;  //AR-20Oct2014 : it was 0.5
                                                   >>  92   vectorMesonMix[3] = 0.0;
                                                   >>  93   vectorMesonMix[4] = 1.0;
                                                   >>  94   vectorMesonMix[5] = 1.0; 
                                                   >>  95 
                                                   >>  96   //... scalarMesonMix[] is quark mixing parameters for scalar mesons (Variable spin=1)
                                                   >>  97   scalarMesonMix.resize(6);
                                                   >>  98   scalarMesonMix[0] = 0.5; 
                                                   >>  99   scalarMesonMix[1] = 0.25; 
                                                   >> 100   scalarMesonMix[2] = 0.5; 
                                                   >> 101   scalarMesonMix[3] = 0.25; 
                                                   >> 102   scalarMesonMix[4] = 1.0; 
                                                   >> 103   scalarMesonMix[5] = 0.5; 
                                                   >> 104 
                                                   >> 105   // Parameters may be changed until the first fragmentation starts
                                                   >> 106   PastInitPhase=false;
                                                   >> 107   hadronizer = new G4HadronBuilder(pspin_meson,pspin_barion,scalarMesonMix,vectorMesonMix);
                                                   >> 108   Kappa = 1.0 * GeV/fermi;
                                                   >> 109 }
120                                                   110 
121    SetMinMasses();  // Re-calculation of minim << 
122                                                   111 
123    Kappa = 1.0 * GeV/fermi;                    << 112 G4VLongitudinalStringDecay::~G4VLongitudinalStringDecay()
124    DecayQuark = NewQuark = 0;                  << 113 {
                                                   >> 114   delete hadronizer;
125 }                                                 115 }
126                                                   116 
127 G4VLongitudinalStringDecay::~G4VLongitudinalSt << 117 //=============================================================================
                                                   >> 118 
                                                   >> 119 // Operators
                                                   >> 120 
                                                   >> 121 //-----------------------------------------------------------------------------
                                                   >> 122 
                                                   >> 123 int G4VLongitudinalStringDecay::operator==(const G4VLongitudinalStringDecay &) const
128 {                                                 124 {
129    delete hadronizer;                          << 125   throw G4HadronicException(__FILE__, __LINE__, "G4VLongitudinalStringDecay::operator== forbidden");
                                                   >> 126   return false;
130 }                                                 127 }
131                                                   128 
132 G4HadFinalState*                               << 129 //-------------------------------------------------------------------------------------
133 G4VLongitudinalStringDecay::ApplyYourself(cons << 130 
                                                   >> 131 int G4VLongitudinalStringDecay::operator!=(const G4VLongitudinalStringDecay &) const
134 {                                                 132 {
135   return nullptr;                              << 133   throw G4HadronicException(__FILE__, __LINE__, "G4VLongitudinalStringDecay::operator!= forbidden");
                                                   >> 134   return true;
136 }                                                 135 }
137                                                   136 
138 //============================================ << 137 //***********************************************************************************
139                                                   138 
140 // For changing Mass Cut used for selection of    139 // For changing Mass Cut used for selection of very small mass strings
141 void G4VLongitudinalStringDecay::SetMassCut(G4 << 140 void G4VLongitudinalStringDecay::SetMassCut(G4double aValue){MassCut=aValue;}
142 G4double G4VLongitudinalStringDecay::GetMassCu << 
143                                                   141 
144 //--------------------------------------------    142 //-----------------------------------------------------------------------------
145                                                   143 
146 // For handling a string with very low mass       144 // For handling a string with very low mass
147                                                   145 
148 G4KineticTrackVector* G4VLongitudinalStringDec << 146 G4KineticTrackVector* G4VLongitudinalStringDecay::
                                                   >> 147 LightFragmentationTest(const G4ExcitedString * const string)
149 {                                                 148 {
150         G4KineticTrackVector* result = nullptr << 149   // Check string decay threshold 
151         pDefPair hadrons( nullptr, nullptr );  << 150   G4KineticTrackVector * result=0;  // return 0 when string exceeds the mass cut
152         G4FragmentingString aString( *string ) << 151   
153                                                << 152   pDefPair hadrons((G4ParticleDefinition *)0,(G4ParticleDefinition *)0);
154         #ifdef debug_VStringDecay              << 153 
155         G4cout<<"G4VLongitudinalStringDecay::P << 154   G4FragmentingString aString(*string);
156               <<aString.Mass()<<" MassCut "<<M << 155   if ( sqr(FragmentationMass(&aString,0,&hadrons)+MassCut) < aString.Mass2()) {
157         #endif                                 << 156     return 0;
                                                   >> 157   }
                                                   >> 158 
                                                   >> 159   // The string mass is very low ---------------------------
                                                   >> 160   
                                                   >> 161   result=new G4KineticTrackVector;
158                                                   162         
159         SetMinimalStringMass( &aString );      << 163   if ( hadrons.second ==0 ) {
160         PossibleHadronMass( &aString, 0, &hadr << 164     // Substitute string by light hadron, Note that Energy is not conserved here!
161         result = new G4KineticTrackVector;     << 165 
162         if ( hadrons.first != nullptr ) {      << 166     #ifdef debug_VStringDecay
163            if ( hadrons.second == nullptr ) {  << 167     G4cout << "VlongSF Warning replacing string by single hadron (G4VLongitudinalStringDecay)" <<G4endl;
164                // Substitute string by light h << 168     G4cout << hadrons.first->GetParticleName()<<G4endl
165                                                << 169      << "string .. " << string->Get4Momentum() << " " 
166                #ifdef debug_VStringDecay       << 170      << string->Get4Momentum().m() << G4endl;
167                G4cout << "VlongSD Warning repl << 171     #endif    
168                G4cout << hadrons.first->GetPar << 172 
169                       << "string .. " << strin << 173     G4ThreeVector   Mom3 = string->Get4Momentum().vect();
170                       << string->Get4Momentum( << 174     G4LorentzVector Mom(Mom3, std::sqrt(Mom3.mag2() + sqr(hadrons.first->GetPDGMass())));
171                #endif                          << 175     result->push_back(new G4KineticTrack(hadrons.first, 0, string->GetPosition(), Mom));
172                                                << 176   } else {
173                G4ThreeVector   Mom3 = string-> << 177     //... string was qq--qqbar type: Build two stable hadrons,
174                G4LorentzVector Mom( Mom3, std: << 178 
175                result->push_back( new G4Kineti << 179     #ifdef debug_VStringDecay
176            } else {                            << 180     G4cout << "VlongSF Warning replacing qq-qqbar string by TWO hadrons (G4VLongitudinalStringDecay)" 
177                //... string was qq--qqbar type << 181      << hadrons.first->GetParticleName() << " / " 
178                                                << 182      << hadrons.second->GetParticleName()
179                #ifdef debug_VStringDecay       << 183      << "string .. " << string->Get4Momentum() << " " 
180                G4cout << "VlongSD Warning repl << 184      << string->Get4Momentum().m() << G4endl;
181                       << hadrons.first->GetPar << 185     #endif          
182                       << hadrons.second->GetPa << 186 
183                       << "string .. " << strin << 187     G4LorentzVector  Mom1, Mom2;
184                       << string->Get4Momentum( << 188     Sample4Momentum(&Mom1, hadrons.first->GetPDGMass(), &Mom2,hadrons.second->GetPDGMass(),
185                #endif                          << 189               string->Get4Momentum().mag());
186                                                << 190 
187                G4LorentzVector  Mom1, Mom2;    << 191     result->push_back(new G4KineticTrack(hadrons.first, 0, string->GetPosition(), Mom1));
188                Sample4Momentum( &Mom1, hadrons << 192     result->push_back(new G4KineticTrack(hadrons.second, 0, string->GetPosition(), Mom2));
189                                 &Mom2, hadrons << 193 
190                                 string->Get4Mo << 194     G4ThreeVector Velocity = string->Get4Momentum().boostVector();
191                                                << 195     result->Boost(Velocity);          
192                result->push_back( new G4Kineti << 196   }
193                result->push_back( new G4Kineti << 197 
194                                                << 198   return result;  
195                G4ThreeVector Velocity = string << 
196                result->Boost(Velocity);        << 
197            }                                   << 
198         }                                      << 
199         return result;                         << 
200 }                                                 199 }
201                                                   200 
202 //--------------------------------------------    201 //----------------------------------------------------------------------------------------
203                                                   202 
204 G4double G4VLongitudinalStringDecay::PossibleH << 203 G4double G4VLongitudinalStringDecay::
205                                                << 204 FragmentationMass( const G4FragmentingString * const string, Pcreate build, pDefPair * pdefs )
206 {                                                 205 {
207         G4double mass = 0.0;                   << 206   G4double mass;
                                                   >> 207   static G4ThreadLocal G4bool NeedInit(true);
                                                   >> 208   static G4ThreadLocal std::vector<double> *nomix_G4MT_TLS_ = 0 ; 
                                                   >> 209   if (!nomix_G4MT_TLS_) nomix_G4MT_TLS_ = new  std::vector<double>;
                                                   >> 210   std::vector<double> &nomix = *nomix_G4MT_TLS_;
                                                   >> 211   static G4ThreadLocal G4HadronBuilder * minMassHadronizer;
                                                   >> 212   if ( NeedInit ) 
                                                   >> 213   {
                                                   >> 214     NeedInit = false;
                                                   >> 215     nomix.resize(6);
                                                   >> 216     for ( G4int i=0; i<6 ; i++ ) nomix[i]=0;
                                                   >> 217 
                                                   >> 218     //minMassHadronizer=new G4HadronBuilder(pspin_meson,pspin_barion,nomix,nomix);
                                                   >> 219     minMassHadronizer=hadronizer;
                                                   >> 220   }
208                                                   221 
209   if ( build==0 ) build=&G4HadronBuilder::Buil << 222   if ( build==0 ) build=&G4HadronBuilder::BuildLowSpin;
                                                   >> 223 
                                                   >> 224   G4ParticleDefinition *Hadron1, *Hadron2=0;
                                                   >> 225 
                                                   >> 226   if (!string->FourQuarkString() )
                                                   >> 227   {
                                                   >> 228     // spin 0 meson or spin 1/2 barion will be built
210                                                   229 
211         G4ParticleDefinition* Hadron1 = nullpt << 230     Hadron1 = (minMassHadronizer->*build)(string->GetLeftParton(), string->GetRightParton());
212   G4ParticleDefinition* Hadron2 = nullptr;     << 231 
213                                                << 232     #ifdef debug_VStringDecay
214         if (!string->IsAFourQuarkString() )    << 233     G4cout<<"Quarks at the string ends "<<string->GetLeftParton()->GetParticleName()
215         {                                      << 234           <<" "<<string->GetRightParton()->GetParticleName()<<G4endl;
216            // spin 0 meson or spin 1/2 barion  << 235     G4cout<<"(G4VLongitudinalStringDecay) Hadron "<<Hadron1->GetParticleName()<<" "<<Hadron1->GetPDGMass()<<G4endl;
217                                                << 236     #endif
218            Hadron1 = (hadronizer->*build)(stri << 237 
219            #ifdef debug_VStringDecay           << 238     mass= (Hadron1)->GetPDGMass();
220      G4cout<<"VlongSD PossibleHadronMass"<<G4e << 239   } else {
221            G4cout<<"VlongSD Quarks at the stri << 240     //... string is qq--qqbar: Build two stable hadrons,
222                  <<" "<<string->GetRightParton << 241     //...    with extra uubar or ddbar quark pair
223            if ( Hadron1 != nullptr) {          << 242     G4int iflc = (G4UniformRand() < 0.5)? 1 : 2;
224              G4cout<<"(G4VLongitudinalStringDe << 243     if (string->GetLeftParton()->GetPDGEncoding() < 0) iflc = -iflc;
225                    <<" "<<Hadron1->GetPDGMass( << 244 
226            }                                   << 245     //... theSpin = 4; spin 3/2 baryons will be built
227            #endif                              << 246     Hadron1 = (minMassHadronizer->*build)(string->GetLeftParton(), FindParticle(iflc));
228            if ( Hadron1 != nullptr) { mass = ( << 247     Hadron2 = (minMassHadronizer->*build)(string->GetRightParton(), FindParticle(-iflc));
229            else                  { mass = MaxM << 248     mass = (Hadron1)->GetPDGMass() + (Hadron2)->GetPDGMass();
230         } else                                 << 249   }
231         {                                      << 
232            //... string is qq--qqbar: Build tw << 
233                                                << 
234            #ifdef debug_VStringDecay           << 
235            G4cout<<"VlongSD PossibleHadronMass << 
236            G4cout<<"VlongSD string is qq--qqba << 
237            #endif                              << 
238                                                << 
239            G4double StringMass   = string->Mas << 
240            G4int cClusterInterrupt = 0;        << 
241            do                                  << 
242            {                                   << 
243              if (cClusterInterrupt++ >= Cluste << 
244                                                << 
245              G4int LeftQuark1= string->GetLeft << 
246              G4int LeftQuark2=(string->GetLeft << 
247                                                << 
248              G4int RightQuark1= string->GetRig << 
249              G4int RightQuark2=(string->GetRig << 
250                                                << 
251              if (G4UniformRand()<0.5) {        << 
252                Hadron1 =hadronizer->Build(Find << 
253                Hadron2 =hadronizer->Build(Find << 
254              } else {                          << 
255                Hadron1 =hadronizer->Build(Find << 
256                Hadron2 =hadronizer->Build(Find << 
257              }                                 << 
258              //... repeat procedure, if mass o << 
259              //... ClusterMassCut = 0.15*GeV m << 
260            }                                   << 
261            while ( Hadron1 == nullptr || Hadro << 
262                    ( StringMass <= Hadron1->Ge << 
263                                                << 
264      mass = (Hadron1)->GetPDGMass() + (Hadron2 << 
265         }                                      << 
266                                                << 
267         #ifdef debug_VStringDecay              << 
268         G4cout<<"VlongSD *Hadrons 1 and 2, pro << 
269         #endif                                 << 
270                                                   250   
271   if ( pdefs != 0 )                            << 251   if ( pdefs != 0 ) 
272   { // need to return hadrons as well....      << 252   { // need to return hadrons as well....
273      pdefs->first  = Hadron1;                  << 253     pdefs->first  = Hadron1;
274      pdefs->second = Hadron2;                  << 254     pdefs->second = Hadron2;
275   }                                            << 255   }
276                                                   256      
277         return mass;                           << 257   return mass;
278 }                                                 258 }
279                                                   259 
280 //--------------------------------------------    260 //----------------------------------------------------------------------------
281                                                   261 
282 G4ParticleDefinition* G4VLongitudinalStringDec    262 G4ParticleDefinition* G4VLongitudinalStringDecay::FindParticle(G4int Encoding) 
283 {                                                 263 {
284   /*                                           << 
285   G4cout<<Encoding<<" G4VLongitudinalStringDec << 
286   for (G4int i=4; i<6;i++){                    << 
287     for (G4int j=1;j<6;j++){                   << 
288       G4cout<<i<<" "<<j<<" ";                  << 
289       G4int Code = 1000 * i + 100 * j +1;      << 
290       G4ParticleDefinition* ptr1 = G4ParticleT << 
291       Code +=2;                                << 
292       G4ParticleDefinition* ptr2 = G4ParticleT << 
293       G4cout<<"Code "<<Code - 2<<" ptr "<<ptr1 << 
294     }                                          << 
295     G4cout<<G4endl;                            << 
296   }                                            << 
297   */                                           << 
298                                                << 
299   G4ParticleDefinition* ptr = G4ParticleTable:    264   G4ParticleDefinition* ptr = G4ParticleTable::GetParticleTable()->FindParticle(Encoding);
300                                                << 265   if (ptr == NULL)
301   if (ptr == nullptr)                          << 
302   {                                               266   {
303      for (size_t i=0; i < NewParticles.size(); << 267     G4cout << "Particle with encoding "<<Encoding<<" does not exist!!!"<<G4endl;
304      {                                         << 268     throw G4HadronicException(__FILE__, __LINE__, "Check your particle table");
305        if ( Encoding == NewParticles[i]->GetPD << 
306      }                                         << 
307   }                                               269   }
308                                                << 
309   return ptr;                                     270   return ptr;    
310 }                                                 271 }
311                                                   272 
312 //********************************************    273 //*********************************************************************************
313 //   For decision on continue or stop string f    274 //   For decision on continue or stop string fragmentation
314 //   virtual G4bool StopFragmenting(const G4Fr    275 //   virtual G4bool StopFragmenting(const G4FragmentingString  * const string)=0;
315 //   virtual G4bool IsItFragmentable(const G4F << 276 //   virtual G4bool IsFragmentable(const G4FragmentingString * const string)=0;
316 //                                             << 277 
317 //   If a string can not fragment, make last b    278 //   If a string can not fragment, make last break into 2 hadrons
318 //   virtual G4bool SplitLast(G4FragmentingStr    279 //   virtual G4bool SplitLast(G4FragmentingString * string, 
319 //                            G4KineticTrackVe    280 //                            G4KineticTrackVector * LeftVector,
320 //                            G4KineticTrackVe    281 //                            G4KineticTrackVector * RightVector)=0;
321 //--------------------------------------------    282 //-----------------------------------------------------------------------------
322 //                                                283 //
323 //   If a string can fragment, do the followin << 284 //   If a string fragments, do the following
324 //                                                285 //
325 //   For transver of a string to its CMS frame    286 //   For transver of a string to its CMS frame
326 //--------------------------------------------    287 //-----------------------------------------------------------------------------
327                                                   288 
328 G4ExcitedString *G4VLongitudinalStringDecay::C << 289 G4ExcitedString *G4VLongitudinalStringDecay::CPExcited(const G4ExcitedString & in)
329 {                                                 290 {
330   G4Parton *Left=new G4Parton(*in.GetLeftParto << 291   G4Parton *Left=new G4Parton(*in.GetLeftParton());
331   G4Parton *Right=new G4Parton(*in.GetRightPar << 292   G4Parton *Right=new G4Parton(*in.GetRightParton());
332   return new G4ExcitedString(Left,Right,in.Get << 293   return new G4ExcitedString(Left,Right,in.GetDirection());
333 }                                                 294 }
334                                                   295 
335 //--------------------------------------------    296 //-----------------------------------------------------------------------------
336                                                   297 
337 G4ParticleDefinition * G4VLongitudinalStringDe << 298 G4KineticTrack * G4VLongitudinalStringDecay::
338                                                << 299 Splitup( G4FragmentingString *string, G4FragmentingString *&newString)
339 {                                                 300 {
340    #ifdef debug_VStringDecay                   << 
341    G4cout<<"VlongSD QuarkSplitup: quark ID "<< << 
342    #endif                                      << 
343                                                << 
344    G4int IsParticle=(decay->GetPDGEncoding()>0 << 
345                                                   301 
346    pDefPair QuarkPair = CreatePartonPair(IsPar << 302   #ifdef debug_VStringDecay
347    created = QuarkPair.second;                 << 303   G4cout<<G4endl;
                                                   >> 304   G4cout<<"Start SplitUP (G4VLongitudinalStringDecay) ========================="<<G4endl;
                                                   >> 305   G4cout<<"String partons: " <<string->GetLeftParton()->GetPDGEncoding()<<" "
                                                   >> 306                              <<string->GetRightParton()->GetPDGEncoding()<<" "
                                                   >> 307         <<"Direction "       <<string->GetDecayDirection()<<G4endl;
                                                   >> 308   #endif
                                                   >> 309 
                                                   >> 310   //... random choice of string end to use for creating the hadron (decay)   
                                                   >> 311   G4int SideOfDecay = (G4UniformRand() < 0.5)? 1: -1;
                                                   >> 312   if (SideOfDecay < 0)
                                                   >> 313   {
                                                   >> 314     string->SetLeftPartonStable();
                                                   >> 315   } else {
                                                   >> 316     string->SetRightPartonStable();
                                                   >> 317   }
348                                                   318 
349    DecayQuark = decay->GetPDGEncoding();       << 319   G4ParticleDefinition *newStringEnd;
350    NewQuark   = created->GetPDGEncoding();     << 320   G4ParticleDefinition * HadronDefinition;
                                                   >> 321   if (string->DecayIsQuark())
                                                   >> 322   {
                                                   >> 323     HadronDefinition= QuarkSplitup(string->GetDecayParton(), newStringEnd);
                                                   >> 324   } else {
                                                   >> 325     HadronDefinition= DiQuarkSplitup(string->GetDecayParton(), newStringEnd);
                                                   >> 326   }      
351                                                   327 
352    #ifdef debug_VStringDecay                   << 328   #ifdef debug_VStringDecay
353    G4cout<<"VlongSD QuarkSplitup: "<<decay->Ge << 329   G4cout<<"The parton "<<string->GetDecayParton()->GetPDGEncoding()<<" "
354    G4cout<<"hadronizer->Build(QuarkPair.first, << 330         <<" produces hadron "<<HadronDefinition->GetParticleName()
355    #endif                                      << 331         <<" and is transformed to "<<newStringEnd->GetPDGEncoding()<<G4endl;
                                                   >> 332   G4cout<<"The side of the string decay Left/Right (1/-1) "<<SideOfDecay<<G4endl;
                                                   >> 333   #endif
                                                   >> 334   // create new String from old, ie. keep Left and Right order, but replace decay
                                                   >> 335 
                                                   >> 336   newString=new G4FragmentingString(*string,newStringEnd); // To store possible quark containt of new string
                                                   >> 337 
                                                   >> 338   #ifdef debug_VStringDecay
                                                   >> 339   G4cout<<"An attempt to determine its energy (SplitEandP)"<<G4endl;
                                                   >> 340   #endif
                                                   >> 341 
                                                   >> 342   G4LorentzVector* HadronMomentum=SplitEandP(HadronDefinition, string, newString);
                                                   >> 343 
                                                   >> 344   delete newString; newString=0;
                                                   >> 345   
                                                   >> 346   G4KineticTrack * Hadron =0;
                                                   >> 347   if ( HadronMomentum != 0 ) {
                                                   >> 348 
                                                   >> 349     #ifdef debug_VStringDecay                     
                                                   >> 350     G4cout<<"The attempt was successful"<<G4endl;
                                                   >> 351     #endif
                                                   >> 352 
                                                   >> 353     G4ThreeVector Pos;
                                                   >> 354     Hadron = new G4KineticTrack(HadronDefinition, 0,Pos, *HadronMomentum);
356                                                   355    
357    return hadronizer->Build(QuarkPair.first, d << 356     newString=new G4FragmentingString(*string,newStringEnd,HadronMomentum);
                                                   >> 357 
                                                   >> 358     delete HadronMomentum;
                                                   >> 359   } else {
                                                   >> 360 
                                                   >> 361     #ifdef debug_VStringDecay
                                                   >> 362     G4cout<<"The attempt was not successful !!!"<<G4endl;
                                                   >> 363     #endif
                                                   >> 364   }
                                                   >> 365 
                                                   >> 366   #ifdef debug_VStringDecay
                                                   >> 367   G4cout<<"End SplitUP (G4VLongitudinalStringDecay) ====================="<<G4endl;
                                                   >> 368   #endif
                                                   >> 369 
                                                   >> 370   return Hadron;
358 }                                                 371 }
359                                                   372 
                                                   >> 373 //--------------------------------------------------------------------------------------
                                                   >> 374 
                                                   >> 375 G4ParticleDefinition * G4VLongitudinalStringDecay::
                                                   >> 376 QuarkSplitup(G4ParticleDefinition* decay, G4ParticleDefinition *&created)
                                                   >> 377 {
                                                   >> 378   // if we have a quark, we need antiquark (or diquark)
                                                   >> 379   G4int IsParticle=(decay->GetPDGEncoding()>0) ? -1 : +1;  
                                                   >> 380                                                             
                                                   >> 381   pDefPair QuarkPair = CreatePartonPair(IsParticle);
                                                   >> 382   created = QuarkPair.second;
                                                   >> 383   return hadronizer->Build(QuarkPair.first, decay);
                                                   >> 384 }
                                                   >> 385 
                                                   >> 386 /*   Uzhi June 2014
360 //--------------------------------------------    387 //-----------------------------------------------------------------------------
361                                                   388 
362 G4VLongitudinalStringDecay::pDefPair G4VLongit << 389 G4ParticleDefinition *G4VLongitudinalStringDecay::
363 CreatePartonPair(G4int NeedParticle,G4bool All << 390 DiQuarkSplitup(G4ParticleDefinition* decay, G4ParticleDefinition *&created)
364 {                                                 391 {
365     //  NeedParticle = +1 for Particle, -1 for << 392   //... can Diquark break or not? 
366     if ( AllowDiquarks && G4UniformRand() < Di << 393   if (G4UniformRand() < DiquarkBreakProb ){
                                                   >> 394 
                                                   >> 395     //... Diquark break
                                                   >> 396 
                                                   >> 397     G4int stableQuarkEncoding = decay->GetPDGEncoding()/1000;
                                                   >> 398     G4int decayQuarkEncoding = (decay->GetPDGEncoding()/100)%10;
                                                   >> 399     if (G4UniformRand() < 0.5)
367     {                                             400     {
368       // Create a Diquark - AntiDiquark pair , << 401       G4int Swap = stableQuarkEncoding;
369       #ifdef debug_VStringDecay                << 402       stableQuarkEncoding = decayQuarkEncoding;
370       G4cout<<"VlongSD Create a Diquark - Anti << 403       decayQuarkEncoding = Swap;
371       #endif                                   << 404     }
372       G4int q1(0), q2(0), spin(0), PDGcode(0); << 405 
373                                                << 406     G4int IsParticle=(decayQuarkEncoding>0) ? -1 : +1;  // if we have a quark, we need antiquark
374       q1  = SampleQuarkFlavor();               << 407 
375       q2  = SampleQuarkFlavor();               << 408     pDefPair QuarkPair = CreatePartonPair(IsParticle,false);  // no diquarks wanted
376                                                << 409 
377       spin = (q1 != q2 && G4UniformRand() <= 0 << 410     //... Build new Diquark
378                                      //   conv << 411     G4int QuarkEncoding=QuarkPair.second->GetPDGEncoding();
379       PDGcode = (std::max(q1,q2) * 1000 + std: << 412     G4int i10  = std::max(std::abs(QuarkEncoding), std::abs(stableQuarkEncoding));
380                                                << 413     G4int i20  = std::min(std::abs(QuarkEncoding), std::abs(stableQuarkEncoding));
381       return pDefPair (FindParticle(-PDGcode), << 414     G4int spin = (i10 != i20 && G4UniformRand() <= 0.5)? 1 : 3;
382                                                << 415     G4int NewDecayEncoding = -1*IsParticle*(i10 * 1000 + i20 * 100 + spin);
383     } else {                                   << 416     created = FindParticle(NewDecayEncoding);
384       // Create a Quark - AntiQuark pair, firs << 417     G4ParticleDefinition * decayQuark=FindParticle(decayQuarkEncoding);
385       #ifdef debug_VStringDecay                << 418     G4ParticleDefinition * had=hadronizer->Build(QuarkPair.first, decayQuark);
386       G4cout<<"VlongSD Create a Quark - AntiQu << 419     return had;
387       #endif                                   << 420     //return hadronizer->Build(QuarkPair.first, decayQuark);
388       G4int PDGcode=SampleQuarkFlavor()*NeedPa << 421    
389       return pDefPair (FindParticle(PDGcode),F << 422   } else {
390     }                                          << 423     //... Diquark does not break
                                                   >> 424  
                                                   >> 425     G4int IsParticle=(decay->GetPDGEncoding()>0) ? +1 : -1;  // if we have a diquark, we need quark
                                                   >> 426     pDefPair QuarkPair = CreatePartonPair(IsParticle,false);  // no diquarks wanted
                                                   >> 427     created = QuarkPair.second;
                                                   >> 428 
                                                   >> 429     G4ParticleDefinition * had=hadronizer->Build(QuarkPair.first, decay);
                                                   >> 430     return had;
                                                   >> 431     //return G4ParticleDefinition * had=hadronizer->Build(QuarkPair.first, decay);
                                                   >> 432   }
391 }                                                 433 }
                                                   >> 434 */  // Uzhi June 2014
392                                                   435 
393 //--------------------------------------------    436 //-----------------------------------------------------------------------------
394                                                   437 
395 G4int G4VLongitudinalStringDecay::SampleQuarkF    438 G4int G4VLongitudinalStringDecay::SampleQuarkFlavor(void)
396 {                                                 439 {
397    G4int  quark(1);                            << 440   return (1 + (int)(G4UniformRand()/StrangeSuppress));
398    G4double ksi = G4UniformRand();             << 441 }
399    if ( ksi < ProbCB ) {                       << 442 
400       if ( ksi < ProbCCbar ) {quark = 4;}   // << 443 //-----------------------------------------------------------------------------
401       else                   {quark = 5;}   // << 444 
402       #ifdef debug_heavyHadrons                << 445 G4VLongitudinalStringDecay::pDefPair G4VLongitudinalStringDecay::
403       G4cout << "G4VLongitudinalStringDecay::S << 446 CreatePartonPair(G4int NeedParticle,G4bool AllowDiquarks)
404        << quark << G4endl;                     << 447 {
405       #endif                                   << 448   //  NeedParticle = +1 for Particle, -1 for Antiparticle
406    } else {                                    << 449 
407      quark = 1 + (int)(G4UniformRand()/Strange << 450   if ( AllowDiquarks && G4UniformRand() < DiquarkSuppress )
408    }                                           << 451   {
409    #ifdef debug_VStringDecay                   << 452     // Create a Diquark - AntiDiquark pair , first in pair is anti to IsParticle
410    G4cout<<"VlongSD SampleQuarkFlavor "<<quark << 453     G4int q1  = SampleQuarkFlavor();
411          <<" "<<ProbCCbar<<" "<<ProbBBbar<<" ) << 454     G4int q2  = SampleQuarkFlavor();
412    #endif                                      << 455     G4int spin = (q1 != q2 && G4UniformRand() <= 0.5)? 1 : 3;
413    return quark;                               << 456     // Convention: quark with higher PDG number is first
                                                   >> 457     G4int PDGcode = (std::max(q1,q2) * 1000 + std::min(q1,q2) * 100 + spin) * NeedParticle;
                                                   >> 458     return pDefPair (FindParticle(-PDGcode),FindParticle(PDGcode));
                                                   >> 459   } else {
                                                   >> 460     // Create a Quark - AntiQuark pair, first in pair  IsParticle
                                                   >> 461     G4int PDGcode=SampleQuarkFlavor()*NeedParticle;
                                                   >> 462     return pDefPair (FindParticle(PDGcode),FindParticle(-PDGcode));
                                                   >> 463   }
                                                   >> 464 
414 }                                                 465 }
415                                                   466 
416 //--------------------------------------------    467 //-----------------------------------------------------------------------------
417                                                   468 
418 G4ThreeVector G4VLongitudinalStringDecay::Samp    469 G4ThreeVector G4VLongitudinalStringDecay::SampleQuarkPt(G4double ptMax)
419 {                                                 470 {
420    G4double Pt;                                << 471   G4double Pt;
421    if ( ptMax < 0 ) {                          << 472   if ( ptMax < 0 ) {
422       // sample full gaussian                  << 473     // sample full gaussian
423       Pt = -G4Log(G4UniformRand());            << 474     Pt = -G4Log(G4UniformRand());
424    } else {                                    << 475   } else {
425       // sample in limited range               << 476     // sample in limited range
426       G4double q = ptMax/SigmaQT;              << 477     Pt = -G4Log(G4RandFlat::shoot(G4Exp(-sqr(ptMax)/sqr(SigmaQT)), 1.));
427       G4double ymin = (q > 20.) ? 0.0 : G4Exp( << 478   }
428       Pt = -G4Log(G4RandFlat::shoot(ymin, 1.)) << 479   Pt = SigmaQT * std::sqrt(Pt);
429    }                                           << 480   G4double phi = 2.*pi*G4UniformRand();
430    Pt = SigmaQT * std::sqrt(Pt);               << 481   return G4ThreeVector(Pt * std::cos(phi),Pt * std::sin(phi),0);
431    G4double phi = 2.*pi*G4UniformRand();       << 
432    return G4ThreeVector(Pt * std::cos(phi),Pt  << 
433 }                                                 482 }
434                                                   483 
435 //********************************************    484 //******************************************************************************
436                                                   485 
437 void G4VLongitudinalStringDecay::CalculateHadr << 486 void G4VLongitudinalStringDecay::
438                                                << 487 CalculateHadronTimePosition(G4double theInitialStringMass, G4KineticTrackVector* Hadrons)
439 {                                                 488 {
440    //   `yo-yo` formation time                 << 489   // 'yo-yo' formation time
441    //   const G4double kappa = 1.0 * GeV/fermi << 490   //const G4double kappa = 1.0 * GeV/fermi/4.;      
442    G4double kappa = GetStringTensionParameter( << 491   G4double kappa = GetStringTensionParameter();
443    for (size_t c1 = 0; c1 < Hadrons->size(); c << 492   for(size_t c1 = 0; c1 < Hadrons->size(); c1++)
444    {                                           << 493   {
445       G4double SumPz = 0;                      << 494     G4double SumPz = 0; 
446       G4double SumE  = 0;                      << 495     G4double SumE  = 0;
447       for (size_t c2 = 0; c2 < c1; c2++)       << 496     for(size_t c2 = 0; c2 < c1; c2++)
448       {                                        << 497     {
449          SumPz += Hadrons->operator[](c2)->Get << 498       SumPz += Hadrons->operator[](c2)->Get4Momentum().pz();
450          SumE  += Hadrons->operator[](c2)->Get << 499       SumE  += Hadrons->operator[](c2)->Get4Momentum().e();   
451       }                                        << 500     } 
452       G4double HadronE  = Hadrons->operator[]( << 501     G4double HadronE  = Hadrons->operator[](c1)->Get4Momentum().e();
453       G4double HadronPz = Hadrons->operator[]( << 502     G4double HadronPz = Hadrons->operator[](c1)->Get4Momentum().pz();
454       Hadrons->operator[](c1)->SetFormationTim << 503     Hadrons->operator[](c1)->SetFormationTime( 
455         (theInitialStringMass - 2.*SumPz + Had << 504         (theInitialStringMass - 2.*SumPz + HadronE - HadronPz)/(2.*kappa)/c_light ); 
456       G4ThreeVector aPosition( 0, 0,           << 505 
457         (theInitialStringMass - 2.*SumE  - Had << 506     G4ThreeVector aPosition(0, 0, (theInitialStringMass - 2.*SumE  - HadronE + HadronPz)/(2.*kappa));
458       Hadrons->operator[](c1)->SetPosition(aPo << 507     Hadrons->operator[](c1)->SetPosition(aPosition);
459    }                                           << 508 
                                                   >> 509   } 
460 }                                                 510 }
461                                                   511 
462 //--------------------------------------------    512 //-----------------------------------------------------------------------------
463                                                   513 
464 void G4VLongitudinalStringDecay::SetSigmaTrans    514 void G4VLongitudinalStringDecay::SetSigmaTransverseMomentum(G4double aValue)
465 {                                                 515 {
466    if ( PastInitPhase ) {                      << 516   if ( PastInitPhase ) {
467      throw G4HadronicException(__FILE__, __LIN << 517     throw G4HadronicException(__FILE__, __LINE__, 
468        "G4VLongitudinalStringDecay::SetSigmaTr << 518                               "4VLongitudinalStringDecay::SetSigmaTransverseMomentum after FragmentString() not allowed");
469    } else {                                    << 519   } else {
470      SigmaQT = aValue;                         << 520     SigmaQT = aValue;
471    }                                           << 521   }
472 }                                                 522 }
473                                                   523 
474 //--------------------------------------------    524 //----------------------------------------------------------------------------------------------------------
475                                                   525 
476 void G4VLongitudinalStringDecay::SetStrangenes    526 void G4VLongitudinalStringDecay::SetStrangenessSuppression(G4double aValue)
477 {                                                 527 {
478    StrangeSuppress = aValue;                   << 528   if ( PastInitPhase ) {
                                                   >> 529     throw G4HadronicException(__FILE__, __LINE__, 
                                                   >> 530                               "4VLongitudinalStringDecay::SetStrangenessSuppression after FragmentString() not allowed");
                                                   >> 531   } else {
                                                   >> 532     StrangeSuppress = aValue;
                                                   >> 533   }
479 }                                                 534 }
480                                                   535 
481 //--------------------------------------------    536 //----------------------------------------------------------------------------------------------------------
482                                                   537 
483 void G4VLongitudinalStringDecay::SetDiquarkSup    538 void G4VLongitudinalStringDecay::SetDiquarkSuppression(G4double aValue)
484 {                                                 539 {
485    DiquarkSuppress = aValue;                   << 540   if ( PastInitPhase ) {
                                                   >> 541     throw G4HadronicException(__FILE__, __LINE__, 
                                                   >> 542                               "4VLongitudinalStringDecay::SetDiquarkSuppression after FragmentString() not allowed");
                                                   >> 543   } else {
                                                   >> 544     DiquarkSuppress = aValue;
                                                   >> 545   }
486 }                                                 546 }
487                                                   547 
488 //--------------------------------------------    548 //----------------------------------------------------------------------------------------
489                                                   549 
490 void G4VLongitudinalStringDecay::SetDiquarkBre    550 void G4VLongitudinalStringDecay::SetDiquarkBreakProbability(G4double aValue)
491 {                                                 551 {
492   if ( PastInitPhase ) {                          552   if ( PastInitPhase ) {
493     throw G4HadronicException(__FILE__, __LINE    553     throw G4HadronicException(__FILE__, __LINE__, 
494       "G4VLongitudinalStringDecay::SetDiquarkB << 554                               "4VLongitudinalStringDecay::SetDiquarkBreakProbability after FragmentString() not allowed");
495   } else {                                        555   } else {
496     DiquarkBreakProb = aValue;                    556     DiquarkBreakProb = aValue;
497   }                                               557   }
498 }                                                 558 }
499                                                   559 
500 //--------------------------------------------    560 //----------------------------------------------------------------------------------------------------------
501                                                   561 
                                                   >> 562 void G4VLongitudinalStringDecay::SetVectorMesonProbability(G4double aValue)
                                                   >> 563 {
                                                   >> 564   if ( PastInitPhase ) {
                                                   >> 565     throw G4HadronicException(__FILE__, __LINE__, 
                                                   >> 566                               "G4VLongitudinalStringDecay::SetVectorMesonProbability after FragmentString() not allowed");
                                                   >> 567   } else {
                                                   >> 568     pspin_meson = aValue;
                                                   >> 569     delete hadronizer;
                                                   >> 570     hadronizer = new G4HadronBuilder(pspin_meson,pspin_barion,scalarMesonMix,vectorMesonMix);
                                                   >> 571   }
                                                   >> 572 }
                                                   >> 573 
                                                   >> 574 //----------------------------------------------------------------------------------------------------------
                                                   >> 575 
502 void G4VLongitudinalStringDecay::SetSpinThreeH    576 void G4VLongitudinalStringDecay::SetSpinThreeHalfBarionProbability(G4double aValue)
503 {                                                 577 {
504   if ( PastInitPhase ) {                          578   if ( PastInitPhase ) {
505     throw G4HadronicException(__FILE__, __LINE    579     throw G4HadronicException(__FILE__, __LINE__, 
506       "G4VLongitudinalStringDecay::SetSpinThre << 580                               "G4VLongitudinalStringDecay::SetSpinThreeHalfBarionProbability after FragmentString() not allowed");
507   } else {                                        581   } else {
508     pspin_barion = aValue;                        582     pspin_barion = aValue;
509     delete hadronizer;                            583     delete hadronizer;
510     hadronizer = new G4HadronBuilder( pspin_me << 584     hadronizer = new G4HadronBuilder(pspin_meson,pspin_barion,scalarMesonMix,vectorMesonMix);
511                                       ProbEta_ << 
512   }                                               585   }
513 }                                                 586 }
514                                                   587 
515 //--------------------------------------------    588 //----------------------------------------------------------------------------------------------------------
516                                                   589 
517 void G4VLongitudinalStringDecay::SetScalarMeso    590 void G4VLongitudinalStringDecay::SetScalarMesonMixings(std::vector<G4double> aVector)
518 {                                                 591 {
519   if ( PastInitPhase ) {                          592   if ( PastInitPhase ) {
520     throw G4HadronicException(__FILE__, __LINE    593     throw G4HadronicException(__FILE__, __LINE__, 
521       "G4VLongitudinalStringDecay::SetScalarMe << 594                               "G4VLongitudinalStringDecay::SetScalarMesonMixings after FragmentString() not allowed");
522   } else {                                        595   } else {
523     if ( aVector.size() < 6 )                  << 596     if ( aVector.size() < 6 )  
524       throw G4HadronicException(__FILE__, __LI    597       throw G4HadronicException(__FILE__, __LINE__, 
525         "G4VLongitudinalStringDecay::SetScalar << 598                                 "G4VLongitudinalStringDecay::SetScalarMesonMixings( argument Vector too small");
526     scalarMesonMix[0] = aVector[0];               599     scalarMesonMix[0] = aVector[0];
527     scalarMesonMix[1] = aVector[1];               600     scalarMesonMix[1] = aVector[1];
528     scalarMesonMix[2] = aVector[2];               601     scalarMesonMix[2] = aVector[2];
529     scalarMesonMix[3] = aVector[3];               602     scalarMesonMix[3] = aVector[3];
530     scalarMesonMix[4] = aVector[4];               603     scalarMesonMix[4] = aVector[4];
531     scalarMesonMix[5] = aVector[5];               604     scalarMesonMix[5] = aVector[5];
532     delete hadronizer;                            605     delete hadronizer;
533     hadronizer = new G4HadronBuilder( pspin_me << 606     hadronizer = new G4HadronBuilder(pspin_meson,pspin_barion,scalarMesonMix,vectorMesonMix);
534                                       ProbEta_ << 
535   }                                               607   }
536 }                                                 608 }
537                                                   609 
538 //--------------------------------------------    610 //----------------------------------------------------------------------------------------------------------
539                                                   611 
540 void G4VLongitudinalStringDecay::SetVectorMeso    612 void G4VLongitudinalStringDecay::SetVectorMesonMixings(std::vector<G4double> aVector)
541 {                                                 613 {
542   if ( PastInitPhase ) {                          614   if ( PastInitPhase ) {
543     throw G4HadronicException(__FILE__, __LINE    615     throw G4HadronicException(__FILE__, __LINE__, 
544       "G4VLongitudinalStringDecay::SetVectorMe << 616                               "G4VLongitudinalStringDecay::SetVectorMesonMixings after FragmentString() not allowed");
545   } else {                                        617   } else {
546     if ( aVector.size() < 6 )                     618     if ( aVector.size() < 6 ) 
547       throw G4HadronicException(__FILE__, __LI    619       throw G4HadronicException(__FILE__, __LINE__, 
548         "G4VLongitudinalStringDecay::SetVector << 620                                 "G4VLongitudinalStringDecay::SetVectorMesonMixings( argument Vector too small");
549     vectorMesonMix[0] = aVector[0];               621     vectorMesonMix[0] = aVector[0];
550     vectorMesonMix[1] = aVector[1];               622     vectorMesonMix[1] = aVector[1];
551     vectorMesonMix[2] = aVector[2];               623     vectorMesonMix[2] = aVector[2];
552     vectorMesonMix[3] = aVector[3];               624     vectorMesonMix[3] = aVector[3];
553     vectorMesonMix[4] = aVector[4];               625     vectorMesonMix[4] = aVector[4];
554     vectorMesonMix[5] = aVector[5];               626     vectorMesonMix[5] = aVector[5];
555     delete hadronizer;                            627     delete hadronizer;
556     hadronizer = new G4HadronBuilder( pspin_me << 628     hadronizer = new G4HadronBuilder(pspin_meson,pspin_barion,scalarMesonMix,vectorMesonMix);
557                                       ProbEta_ << 
558   }                                               629   }
559 }                                                 630 }
560                                                   631 
561 //--------------------------------------------    632 //-------------------------------------------------------------------------------------------
562                                                << 
563 void G4VLongitudinalStringDecay::SetProbCCbar( << 
564 {                                              << 
565    ProbCCbar = aValue;                         << 
566    ProbCB = ProbCCbar + ProbBBbar;             << 
567 }                                              << 
568                                                << 
569 //-------------------------------------------- << 
570                                                << 
571 void G4VLongitudinalStringDecay::SetProbEta_c( << 
572 {                                              << 
573    ProbEta_c = aValue;                         << 
574 }                                              << 
575                                                << 
576 //-------------------------------------------- << 
577                                                << 
578 void G4VLongitudinalStringDecay::SetProbBBbar( << 
579 {                                              << 
580    ProbBBbar = aValue;                         << 
581    ProbCB = ProbCCbar + ProbBBbar;             << 
582 }                                              << 
583                                                << 
584 //-------------------------------------------- << 
585                                                << 
586 void G4VLongitudinalStringDecay::SetProbEta_b( << 
587 {                                              << 
588    ProbEta_b = aValue;                         << 
589 }                                              << 
590                                                << 
591 //-------------------------------------------- << 
592                                                << 
593 void G4VLongitudinalStringDecay::SetStringTens    633 void G4VLongitudinalStringDecay::SetStringTensionParameter(G4double aValue)
594 {                                                 634 {
595    Kappa = aValue * GeV/fermi;                 << 635   Kappa = aValue * GeV/fermi;
596 }                                                 636 } 
597                                                << 637 //**************************************************************************************
598 //-------------------------------------------- << 
599                                                << 
600 void G4VLongitudinalStringDecay::SetMinMasses( << 
601 {                                              << 
602     // ------ For estimation of a minimal stri << 
603     Mass_of_light_quark =140.*MeV;             << 
604     Mass_of_s_quark     =500.*MeV;             << 
605     Mass_of_c_quark     =1600.*MeV;            << 
606     Mass_of_b_quark     =4500.*MeV;            << 
607     Mass_of_string_junction=720.*MeV;          << 
608                                                << 
609     // ---------------- Determination of minim << 
610     G4ParticleDefinition * hadron1;    G4int C << 
611     G4ParticleDefinition * hadron2;    G4int C << 
612     for (G4int i=1; i < 6; i++) {              << 
613         Code1 = 100*i + 10*1 + 1;              << 
614         hadron1 = FindParticle(Code1);         << 
615                                                << 
616         if (hadron1 != nullptr) {              << 
617            for (G4int j=1; j < 6; j++) {       << 
618                Code2 = 100*j + 10*1 + 1;       << 
619                hadron2 = FindParticle(Code2);  << 
620                if (hadron2 != nullptr) {       << 
621                  minMassQQbarStr[i-1][j-1] = h << 
622                }                               << 
623            }                                   << 
624         }                                      << 
625     }                                          << 
626                                                << 
627     minMassQQbarStr[1][1] = minMassQQbarStr[0] << 
628                                                << 
629     // ---------------- Determination of minim << 
630     G4ParticleDefinition * hadron3;            << 
631     G4int kfla, kflb;                          << 
632     //  MaxMass = -350.0*GeV;   // If there wi << 
633                                                << 
634     for (G4int i=1; i < 6; i++) {   //i=1      << 
635         Code1 = 100*i + 10*1 + 1;              << 
636         hadron1 = FindParticle(Code1);         << 
637         for (G4int j=1; j < 6; j++) {          << 
638             for (G4int k=1; k < 6; k++) {      << 
639                 kfla = std::max(j,k);          << 
640                 kflb = std::min(j,k);          << 
641                                                << 
642     // Add d-quark                             << 
643                 Code2 = 1000*kfla + 100*kflb + << 
644     if ( (j == 1) && (k==1)) Code2 = 1000*2 +  << 
645                                                << 
646                 hadron2 = G4ParticleTable::Get << 
647                 hadron3 = G4ParticleTable::Get << 
648                                                << 
649                 if ((hadron2 == nullptr) && (h << 
650                                                << 
651                 if ((hadron2 != nullptr) && (h << 
652                    if (hadron2->GetPDGMass() > << 
653                 };                             << 
654                                                << 
655                 if ((hadron2 != nullptr) && (h << 
656                                                << 
657                 if ((hadron2 == nullptr) && (h << 
658                                                << 
659                 minMassQDiQStr[i-1][j-1][k-1]  << 
660             }                                  << 
661         }                                      << 
662     }                                          << 
663                                                << 
664     // ------ An estimated minimal string mass << 
665     MinimalStringMass  = 0.;                   << 
666     MinimalStringMass2 = 0.;                   << 
667     // q charges  d               u            << 
668     Qcharge[0] = -1; Qcharge[1] = 2; Qcharge[2 << 
669                                                << 
670     // For treating of small string decays     << 
671     for (G4int i=0; i<5; i++)                  << 
672     {  for (G4int j=0; j<5; j++)               << 
673        {  for (G4int k=0; k<7; k++)            << 
674           {                                    << 
675             Meson[i][j][k]=0; MesonWeight[i][j << 
676           }                                    << 
677        }                                       << 
678     }                                          << 
679     //--------------------------               << 
680     G4int StrangeQ = 0;                        << 
681     G4int StrangeAQ = 0;                       << 
682     for (G4int i=0; i<5; i++)                  << 
683     {                                          << 
684        if( i >= 2 ) StrangeQ=1;                << 
685        for (G4int j=0; j<5; j++)               << 
686        {                                       << 
687          StrangeAQ = 0;                        << 
688          if( j >= 2 ) StrangeAQ=1;             << 
689          Meson[i][j][0]       = 100 * (std::ma << 
690          MesonWeight[i][j][0] = (   pspin_meso << 
691          Meson[i][j][1]       = 100 * (std::ma << 
692          MesonWeight[i][j][1] = (1.-pspin_meso << 
693        }                                       << 
694     }                                          << 
695                                                << 
696     //qqs                                      << 
697     //dd1 -> scalarMesonMix[0] * 111 + (1-scal << 
698     //dd1 ->                     Pi0           << 
699                                                << 
700     Meson[0][0][0] = 111; MesonWeight[0][0][0] << 
701     Meson[0][0][2] = 221; MesonWeight[0][0][3] << 
702     Meson[0][0][3] = 331; MesonWeight[0][0][4] << 
703                                                << 
704     //dd3 -> (1-vectorMesonMix[1] * 113 + vect << 
705     //dd3 ->                       rho_0       << 
706                                                << 
707     Meson[0][0][1] = 113; MesonWeight[0][0][1] << 
708     Meson[0][0][4] = 223; MesonWeight[0][0][4] << 
709                                                << 
710     //uu1 -> scalarMesonMix[0] * 111 + (1-scal << 
711     //uu1 ->                     Pi0           << 
712                                                << 
713     Meson[1][1][0] = 111; MesonWeight[1][1][0] << 
714     Meson[1][1][2] = 221; MesonWeight[1][1][2] << 
715     Meson[1][1][3] = 331; MesonWeight[1][1][3] << 
716                                                << 
717     //uu3 -> (1-vectorMesonMix[1]) * 113 + vec << 
718     //uu3 ->                        rho_0      << 
719                                                << 
720     Meson[1][1][1] = 113; MesonWeight[1][1][1] << 
721     Meson[1][1][4] = 223; MesonWeight[1][1][4] << 
722                                                << 
723     //ss1     ->                               << 
724     //ss1     ->                               << 
725                                                << 
726     Meson[2][2][0] = 221; MesonWeight[2][2][0] << 
727     Meson[2][2][2] = 331; MesonWeight[2][2][2] << 
728                                                << 
729     //ss3     ->                               << 
730     //ss3     ->                               << 
731                                                << 
732     Meson[2][2][1] = 333; MesonWeight[2][2][1] << 
733                                                << 
734     //cc1     ->    ProbEta_c /(1-pspin_meson) << 
735     //cc3     -> (1-ProbEta_c)/(  pspin_meson) << 
736                                                << 
737     //bb1     ->    ProbEta_b /pspin_meson 551 << 
738     //bb3     -> (1-ProbEta_b)/pspin_meson 553 << 
739                                                << 
740     if ( pspin_meson[2] != 0. ) {              << 
741        Meson[3][3][0] *= (    ProbEta_c)/(   p << 
742        Meson[3][3][1] *= (1.0-ProbEta_c)/(1.-p << 
743                                                << 
744        Meson[4][4][0] *= (    ProbEta_b)/(   p << 
745        Meson[4][4][1] *= (1.0-ProbEta_b)/(1.-p << 
746     }                                          << 
747                                                << 
748     //--------------------------               << 
749                                                << 
750     for (G4int i=0; i<5; i++)                  << 
751     {  for (G4int j=0; j<5; j++)               << 
752        {  for (G4int k=0; k<5; k++)            << 
753           {  for (G4int l=0; l<4; l++)         << 
754              { Baryon[i][j][k][l]=0; BaryonWei << 
755           }                                    << 
756        }                                       << 
757     }                                          << 
758                                                << 
759           kfla =0;  kflb =0;                   << 
760     G4int                   kflc(0), kfld(0),  << 
761     for (G4int i=0; i<5; i++)                  << 
762     {  for (G4int j=0; j<5; j++)               << 
763        {  for (G4int k=0; k<5; k++)            << 
764           {                                    << 
765            kfla = i+1; kflb = j+1; kflc = k+1; << 
766      kfld = std::max(kfla,kflb);               << 
767      kfld = std::max(kfld,kflc);               << 
768                                                << 
769      kflf = std::min(kfla,kflb);               << 
770      kflf = std::min(kflf,kflc);               << 
771                                                << 
772            kfle = kfla + kflb + kflc - kfld -  << 
773                                                << 
774            Baryon[i][j][k][0]       = 1000 * k << 
775            BaryonWeight[i][j][k][0] = (   pspi << 
776            Baryon[i][j][k][1]       = 1000 * k << 
777            BaryonWeight[i][j][k][1] = (1.-pspi << 
778           }                                    << 
779        }                                       << 
780     }                                          << 
781                                                << 
782     // Delta-  ddd - only 1114                 << 
783     Baryon[0][0][0][0] = 1114;    BaryonWeight << 
784     Baryon[0][0][0][1] =    0;    BaryonWeight << 
785                                                << 
786     // Delta++ uuu - only 2224                 << 
787     Baryon[1][1][1][0] = 2224;    BaryonWeight << 
788     Baryon[1][1][1][1] =    0;    BaryonWeight << 
789                                                << 
790     // Omega- sss - only 3334                  << 
791     Baryon[2][2][2][0] = 3334;    BaryonWeight << 
792     Baryon[2][2][2][1] =    0;    BaryonWeight << 
793                                                << 
794     // Omega_cc++ ccc - only 4444              << 
795     Baryon[3][3][3][0] = 4444;    BaryonWeight << 
796     Baryon[3][3][3][1] =    0;    BaryonWeight << 
797                                                << 
798     // Omega_bb-  bbb - only 5554              << 
799     Baryon[4][4][4][0] = 5554;    BaryonWeight << 
800     Baryon[4][4][4][1] =    0;    BaryonWeight << 
801                                                << 
802     // Lambda/Sigma0 sud - 3122/3212           << 
803     Baryon[0][1][2][0] = 3122;    BaryonWeight << 
804     Baryon[0][2][1][0] = 3122;    BaryonWeight << 
805     Baryon[1][0][2][0] = 3122;    BaryonWeight << 
806     Baryon[1][2][0][0] = 3122;    BaryonWeight << 
807     Baryon[2][0][1][0] = 3122;    BaryonWeight << 
808     Baryon[2][1][0][0] = 3122;    BaryonWeight << 
809                                                << 
810     Baryon[0][1][2][2] = 3212;    BaryonWeight << 
811     Baryon[0][2][1][2] = 3212;    BaryonWeight << 
812     Baryon[1][0][2][2] = 3212;    BaryonWeight << 
813     Baryon[1][2][0][2] = 3212;    BaryonWeight << 
814     Baryon[2][0][1][2] = 3212;    BaryonWeight << 
815     Baryon[2][1][0][2] = 3212;    BaryonWeight << 
816                                                << 
817     // Lambda_c+/Sigma_c+ cud - 4122/4212      << 
818     Baryon[0][1][3][0] = 4122;    BaryonWeight << 
819     Baryon[0][3][1][0] = 4122;    BaryonWeight << 
820     Baryon[1][0][3][0] = 4122;    BaryonWeight << 
821     Baryon[1][3][0][0] = 4122;    BaryonWeight << 
822     Baryon[3][0][1][0] = 4122;    BaryonWeight << 
823     Baryon[3][1][0][0] = 4122;    BaryonWeight << 
824                                                << 
825     Baryon[0][1][3][2] = 4212;    BaryonWeight << 
826     Baryon[0][3][1][2] = 4212;    BaryonWeight << 
827     Baryon[1][0][3][2] = 4212;    BaryonWeight << 
828     Baryon[1][3][0][2] = 4212;    BaryonWeight << 
829     Baryon[3][0][1][2] = 4212;    BaryonWeight << 
830     Baryon[3][1][0][2] = 4212;    BaryonWeight << 
831                                                << 
832     // Xi_c+/Xi_c+' cus - 4232/4322            << 
833     Baryon[1][2][3][0] = 4232;    BaryonWeight << 
834     Baryon[1][3][2][0] = 4232;    BaryonWeight << 
835     Baryon[2][1][3][0] = 4232;    BaryonWeight << 
836     Baryon[2][3][1][0] = 4232;    BaryonWeight << 
837     Baryon[3][1][2][0] = 4232;    BaryonWeight << 
838     Baryon[3][2][1][0] = 4232;    BaryonWeight << 
839                                                << 
840     Baryon[1][2][3][2] = 4322;    BaryonWeight << 
841     Baryon[1][3][2][2] = 4322;    BaryonWeight << 
842     Baryon[2][1][3][2] = 4322;    BaryonWeight << 
843     Baryon[2][3][1][2] = 4322;    BaryonWeight << 
844     Baryon[3][1][2][2] = 4322;    BaryonWeight << 
845     Baryon[3][2][1][2] = 4322;    BaryonWeight << 
846                                                << 
847     // Xi_c0/Xi_c0' cus - 4132/4312            << 
848     Baryon[0][2][3][0] = 4132;    BaryonWeight << 
849     Baryon[0][3][2][0] = 4132;    BaryonWeight << 
850     Baryon[2][0][3][0] = 4132;    BaryonWeight << 
851     Baryon[2][3][0][0] = 4132;    BaryonWeight << 
852     Baryon[3][0][2][0] = 4132;    BaryonWeight << 
853     Baryon[3][2][0][0] = 4132;    BaryonWeight << 
854                                                << 
855     Baryon[0][2][3][2] = 4312;    BaryonWeight << 
856     Baryon[0][3][2][2] = 4312;    BaryonWeight << 
857     Baryon[2][0][3][2] = 4312;    BaryonWeight << 
858     Baryon[2][3][0][2] = 4312;    BaryonWeight << 
859     Baryon[3][0][2][2] = 4312;    BaryonWeight << 
860     Baryon[3][2][0][2] = 4312;    BaryonWeight << 
861                                                << 
862     // Lambda_b0/Sigma_b0 bud - 5122/5212      << 
863     Baryon[0][1][4][0] = 5122;    BaryonWeight << 
864     Baryon[0][4][1][0] = 5122;    BaryonWeight << 
865     Baryon[1][0][4][0] = 5122;    BaryonWeight << 
866     Baryon[1][4][0][0] = 5122;    BaryonWeight << 
867     Baryon[4][0][1][0] = 5122;    BaryonWeight << 
868     Baryon[4][1][0][0] = 5122;    BaryonWeight << 
869                                                << 
870     Baryon[0][1][4][2] = 5212;    BaryonWeight << 
871     Baryon[0][4][1][2] = 5212;    BaryonWeight << 
872     Baryon[1][0][4][2] = 5212;    BaryonWeight << 
873     Baryon[1][4][0][2] = 5212;    BaryonWeight << 
874     Baryon[4][0][1][2] = 5212;    BaryonWeight << 
875     Baryon[4][1][0][2] = 5212;    BaryonWeight << 
876                                                << 
877     // Xi_b0/Xi_b0' bus - 5232/5322            << 
878     Baryon[1][2][4][0] = 5232;    BaryonWeight << 
879     Baryon[1][4][2][0] = 5232;    BaryonWeight << 
880     Baryon[2][1][4][0] = 5232;    BaryonWeight << 
881     Baryon[2][4][1][0] = 5232;    BaryonWeight << 
882     Baryon[4][1][2][0] = 5232;    BaryonWeight << 
883     Baryon[4][2][1][0] = 5232;    BaryonWeight << 
884                                                << 
885     Baryon[1][2][4][2] = 5322;    BaryonWeight << 
886     Baryon[1][4][2][2] = 5322;    BaryonWeight << 
887     Baryon[2][1][4][2] = 5322;    BaryonWeight << 
888     Baryon[2][4][1][2] = 5322;    BaryonWeight << 
889     Baryon[4][1][2][2] = 5322;    BaryonWeight << 
890     Baryon[4][2][1][2] = 5322;    BaryonWeight << 
891                                                << 
892     // Xi_b-/Xi_b-' bus - 5132/5312            << 
893     Baryon[0][2][4][0] = 5132;    BaryonWeight << 
894     Baryon[0][4][2][0] = 5132;    BaryonWeight << 
895     Baryon[2][0][4][0] = 5132;    BaryonWeight << 
896     Baryon[2][4][0][0] = 5132;    BaryonWeight << 
897     Baryon[4][0][2][0] = 5132;    BaryonWeight << 
898     Baryon[4][2][0][0] = 5132;    BaryonWeight << 
899                                                << 
900     Baryon[0][2][4][2] = 5312;    BaryonWeight << 
901     Baryon[0][4][2][2] = 5312;    BaryonWeight << 
902     Baryon[2][0][4][2] = 5312;    BaryonWeight << 
903     Baryon[2][4][0][2] = 5312;    BaryonWeight << 
904     Baryon[4][0][2][2] = 5312;    BaryonWeight << 
905     Baryon[4][2][0][2] = 5312;    BaryonWeight << 
906                                                << 
907     for (G4int i=0; i<5; i++)                  << 
908     {  for (G4int j=0; j<5; j++)               << 
909     {  for (G4int k=0; k<5; k++)               << 
910      {  for (G4int l=0; l<4; l++)              << 
911         {                                      << 
912                      G4ParticleDefinition * Te << 
913                        G4ParticleTable::GetPar << 
914                      /*                        << 
915                      G4cout<<i<<" "<<j<<" "<<k << 
916                      if (TestHadron != nullptr << 
917                      if ((TestHadron == nullpt << 
918                      if ((TestHadron == nullpt << 
919                      G4cout<<G4endl;           << 
920                      */                        << 
921                      if ((TestHadron == nullpt << 
922                     }                          << 
923      }                                         << 
924     }                                          << 
925     }                                          << 
926                                                << 
927     // --------- Probabilities of q-qbar pair  << 
928     G4double ProbUUbar = 0.33;                 << 
929     Prob_QQbar[0]=ProbUUbar;         // Probab << 
930     Prob_QQbar[1]=ProbUUbar;         // Probab << 
931     Prob_QQbar[2]=1.0-2.*ProbUUbar;  // Probab << 
932     Prob_QQbar[3]=0.0;               // Probab << 
933     Prob_QQbar[4]=0.0;               // Probab << 
934                                                << 
935     for ( G4int i=0 ; i<350 ; i++ ) { // Must  << 
936       FS_LeftHadron[i] = 0;                    << 
937       FS_RightHadron[i] = 0;                   << 
938       FS_Weight[i] = 0.0;                      << 
939     }                                          << 
940                                                << 
941     NumberOf_FS = 0;                           << 
942 }                                              << 
943                                                << 
944 // ------------------------------------------- << 
945                                                << 
946 void G4VLongitudinalStringDecay::SetMinimalStr << 
947 {                                              << 
948         //MaxMass = -350.0*GeV;                << 
949   G4double EstimatedMass=MaxMass;              << 
950                                                << 
951         G4ParticleDefinition* LeftParton  = st << 
952         G4ParticleDefinition* RightParton = st << 
953         if( LeftParton->GetParticleSubType() = << 
954           if( LeftParton->GetPDGEncoding() * R << 
955             // Not allowed combination of the  << 
956             throw G4HadronicException(__FILE__ << 
957               "G4VLongitudinalStringDecay::Set << 
958           }                                    << 
959         }                                      << 
960         if( LeftParton->GetParticleSubType() ! << 
961           if( LeftParton->GetPDGEncoding() * R << 
962             // Not allowed combination of the  << 
963             throw G4HadronicException(__FILE__ << 
964               "G4VLongitudinalStringDecay::Set << 
965           }                                    << 
966         }                                      << 
967                                                << 
968         G4int Qleft =std::abs(string->GetLeftP << 
969         G4int Qright=std::abs(string->GetRight << 
970                                                << 
971         if ((Qleft < 6) && (Qright < 6)) {   / << 
972           EstimatedMass=minMassQQbarStr[Qleft- << 
973           MinimalStringMass=EstimatedMass;     << 
974           SetMinimalStringMass2(EstimatedMass) << 
975           return;                              << 
976         }                                      << 
977                                                << 
978         if ((Qleft < 6) && (Qright > 1000)) {  << 
979           G4int q1=Qright/1000;                << 
980           G4int q2=(Qright/100)%10;            << 
981           EstimatedMass=minMassQDiQStr[Qleft-1 << 
982           MinimalStringMass=EstimatedMass;     << 
983           SetMinimalStringMass2(EstimatedMass) << 
984           return;                              << 
985         }                                      << 
986                                                << 
987         if ((Qleft > 1000) && (Qright < 6)) {  << 
988           G4int q1=Qleft/1000;                 << 
989           G4int q2=(Qleft/100)%10;             << 
990           EstimatedMass=minMassQDiQStr[Qright- << 
991           MinimalStringMass=EstimatedMass;     << 
992           SetMinimalStringMass2(EstimatedMass) << 
993           return;                              << 
994         }                                      << 
995                                                << 
996         // DiQuark - Anti DiQuark string ----- << 
997                                                << 
998   G4double StringM=string->Get4Momentum().mag( << 
999                                                << 
1000         #ifdef debug_LUNDfragmentation        << 
1001         // G4cout<<"MinStringMass// Input Str << 
1002         #endif                                << 
1003                                               << 
1004         G4int q1= Qleft/1000    ;             << 
1005         G4int q2=(Qleft/100)%10 ;             << 
1006                                               << 
1007         G4int q3= Qright/1000   ;             << 
1008         G4int q4=(Qright/100)%10;             << 
1009                                               << 
1010         // -------------- 2 baryon production << 
1011                                               << 
1012         G4double EstimatedMass1 = minMassQDiQ << 
1013         G4double EstimatedMass2 = minMassQDiQ << 
1014         // Mass is negative if there is no co << 
1015                                               << 
1016         if ( (EstimatedMass1 > 0.) && (Estima << 
1017            EstimatedMass = EstimatedMass1 + E << 
1018            if ( StringM > EstimatedMass ) {   << 
1019               MinimalStringMass=EstimatedMass << 
1020               SetMinimalStringMass2(Estimated << 
1021               return;                         << 
1022           }                                   << 
1023         }                                     << 
1024                                               << 
1025         if ( (EstimatedMass1 < 0.) && (Estima << 
1026            EstimatedMass = MaxMass;           << 
1027            MinimalStringMass=EstimatedMass;   << 
1028            SetMinimalStringMass2(EstimatedMas << 
1029            return;                            << 
1030         }                                     << 
1031                                               << 
1032         if ( (EstimatedMass1 > 0.) && (Estima << 
1033            EstimatedMass = EstimatedMass1;    << 
1034            MinimalStringMass=EstimatedMass;   << 
1035            SetMinimalStringMass2(EstimatedMas << 
1036            return;                            << 
1037         }                                     << 
1038                                               << 
1039         //      if ( EstimatedMass >= StringM << 
1040         // ------------- Re-orangement ------ << 
1041         EstimatedMass=std::min(minMassQQbarSt << 
1042                                minMassQQbarSt << 
1043                                               << 
1044         // In principle, re-arrangement and 2 << 
1045         // More physics consideration is need << 
1046                                               << 
1047         MinimalStringMass=EstimatedMass;      << 
1048         SetMinimalStringMass2(EstimatedMass); << 
1049                                               << 
1050         return;                               << 
1051 }                                             << 
1052                                               << 
1053 //------------------------------------------- << 
1054                                               << 
1055 void G4VLongitudinalStringDecay::SetMinimalSt << 
1056 {                                             << 
1057   MinimalStringMass2=aValue * aValue;         << 
1058 }                                             << 
1059                                                  638 
1060                                                  639