| Geant4 Cross Reference |
1 // 1 //
2 // ******************************************* 2 // ********************************************************************
3 // * License and Disclaimer 3 // * License and Disclaimer *
4 // * 4 // * *
5 // * The Geant4 software is copyright of th 5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided 6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License 7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/ 8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. 9 // * include a list of copyright holders. *
10 // * 10 // * *
11 // * Neither the authors of this software syst 11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing fin 12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warran 13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assum 14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file 15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitatio 16 // * for the full disclaimer and the limitation of liability. *
17 // * 17 // * *
18 // * This code implementation is the result 18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboratio 19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distri 20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you ag 21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publicati 22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Sof 23 // * acceptance of all terms of the Geant4 Software license. *
24 // ******************************************* 24 // ********************************************************************
25 // 25 //
26 // 26 //
27 // 27 //
28 28
29 #include <utility> 29 #include <utility>
30 30
31 #include "G4FTFParameters.hh" 31 #include "G4FTFParameters.hh"
32 32
33 #include "G4ios.hh" 33 #include "G4ios.hh"
34 #include "G4PhysicalConstants.hh" 34 #include "G4PhysicalConstants.hh"
35 #include "G4SystemOfUnits.hh" 35 #include "G4SystemOfUnits.hh"
36 36
37 #include "G4ParticleDefinition.hh" 37 #include "G4ParticleDefinition.hh"
38 #include "G4Proton.hh" 38 #include "G4Proton.hh"
39 #include "G4Neutron.hh" 39 #include "G4Neutron.hh"
40 #include "G4PionPlus.hh" 40 #include "G4PionPlus.hh"
41 #include "G4PionMinus.hh" 41 #include "G4PionMinus.hh"
42 #include "G4KaonPlus.hh" 42 #include "G4KaonPlus.hh"
43 #include "G4KaonMinus.hh" 43 #include "G4KaonMinus.hh"
44 44
45 #include "G4CrossSectionDataSetRegistry.hh" 45 #include "G4CrossSectionDataSetRegistry.hh"
46 #include "G4VComponentCrossSection.hh" 46 #include "G4VComponentCrossSection.hh"
47 #include "G4ComponentGGHadronNucleusXsc.hh" 47 #include "G4ComponentGGHadronNucleusXsc.hh"
48 #include "G4LundStringFragmentation.hh" 48 #include "G4LundStringFragmentation.hh"
49 49
50 #include "G4Exp.hh" 50 #include "G4Exp.hh"
51 #include "G4Log.hh" 51 #include "G4Log.hh"
52 #include "G4Pow.hh" 52 #include "G4Pow.hh"
53 53
54 #include "G4HadronicDeveloperParameters.hh" 54 #include "G4HadronicDeveloperParameters.hh"
55 #include "G4HadronicParameters.hh" 55 #include "G4HadronicParameters.hh"
56 56
57 //============================================ 57 //============================================================================
58 58
59 //#define debugFTFparams 59 //#define debugFTFparams
60 60
61 //============================================ 61 //============================================================================
62 62
63 G4FTFParameters::G4FTFParameters() 63 G4FTFParameters::G4FTFParameters()
64 { 64 {
65 // Set-up alternative sets of FTF parameters 65 // Set-up alternative sets of FTF parameters (called "tunes").
66 // Note that the very first tune (with index 66 // Note that the very first tune (with indexTune == 0) corresponds to the default
67 // set of parameters, which does not need to 67 // set of parameters, which does not need to be set-up explicitly: that's why
68 // the for loop below starts from 1 and not 68 // the for loop below starts from 1 and not from 0.
69 // The check whether an alternative tune has 69 // The check whether an alternative tune has been switched on is done at the
70 // level of the G4FTFParamCollection::SetTun 70 // level of the G4FTFParamCollection::SetTune method.
71 for ( G4int indexTune = 1; indexTune < G4FTF 71 for ( G4int indexTune = 1; indexTune < G4FTFTunings::sNumberOfTunes; ++indexTune ) {
72 fArrayParCollBaryonProj[indexTune].SetTune 72 fArrayParCollBaryonProj[indexTune].SetTune(indexTune);
73 fArrayParCollMesonProj[indexTune].SetTune( 73 fArrayParCollMesonProj[indexTune].SetTune(indexTune);
74 fArrayParCollPionProj[indexTune].SetTune(i 74 fArrayParCollPionProj[indexTune].SetTune(indexTune);
75 } 75 }
76 76
77 StringMass = new G4LundStringFragmentation; 77 StringMass = new G4LundStringFragmentation; // for estimation of min. mass of diffr. states
78 Reset(); 78 Reset();
79 csGGinstance = 79 csGGinstance =
80 G4CrossSectionDataSetRegistry::Instance()- 80 G4CrossSectionDataSetRegistry::Instance()->GetComponentCrossSection("Glauber-Gribov");
81 if (!csGGinstance) { 81 if (!csGGinstance) {
82 csGGinstance = new G4ComponentGGHadronNucl 82 csGGinstance = new G4ComponentGGHadronNucleusXsc();
83 } 83 }
84 84
85 EnableDiffDissociationForBGreater10 = G4Hadr 85 EnableDiffDissociationForBGreater10 = G4HadronicParameters::Instance()->EnableDiffDissociationForBGreater10();
86 86
87 // Set parameters of a string kink 87 // Set parameters of a string kink
88 SetPt2Kink( 0.0*GeV*GeV ); // To switch off 88 SetPt2Kink( 0.0*GeV*GeV ); // To switch off kinky strings (bad results obtained with 6.0*GeV*GeV)
89 G4double Puubar( 1.0/3.0 ), Pddbar( 1.0/3.0 89 G4double Puubar( 1.0/3.0 ), Pddbar( 1.0/3.0 ), Pssbar( 1.0/3.0 ); // SU(3) symmetry
90 //G4double Puubar( 0.41 ), Pddbar( 0.41 ), P 90 //G4double Puubar( 0.41 ), Pddbar( 0.41 ), Pssbar( 0.18 ); // Broken SU(3) symmetry
91 SetQuarkProbabilitiesAtGluonSplitUp( Puubar, 91 SetQuarkProbabilitiesAtGluonSplitUp( Puubar, Pddbar, Pssbar );
92 } 92 }
93 93
94 //============================================ 94 //============================================================================
95 95
96 void G4FTFParameters::InitForInteraction( cons 96 void G4FTFParameters::InitForInteraction( const G4ParticleDefinition* particle,
97 G4in 97 G4int theA, G4int theZ, G4double PlabPerParticle )
98 { 98 {
99 Reset(); 99 Reset();
100 100
101 G4int ProjectilePDGcode = particle->Ge 101 G4int ProjectilePDGcode = particle->GetPDGEncoding();
102 G4int ProjectileabsPDGcode = std::abs( Pr 102 G4int ProjectileabsPDGcode = std::abs( ProjectilePDGcode );
103 G4double ProjectileMass = particle->Ge 103 G4double ProjectileMass = particle->GetPDGMass();
104 G4double ProjectileMass2 = ProjectileMa 104 G4double ProjectileMass2 = ProjectileMass * ProjectileMass;
105 105
106 G4int ProjectileBaryonNumber( 0 ), AbsProjec 106 G4int ProjectileBaryonNumber( 0 ), AbsProjectileBaryonNumber( 0 ), AbsProjectileCharge( 0 );
107 G4bool ProjectileIsNucleus = false; 107 G4bool ProjectileIsNucleus = false;
108 108
109 if ( std::abs( particle->GetBaryonNumber() ) 109 if ( std::abs( particle->GetBaryonNumber() ) > 1 ) { // The projectile is a nucleus
110 ProjectileIsNucleus = true; 110 ProjectileIsNucleus = true;
111 ProjectileBaryonNumber = particle->GetB 111 ProjectileBaryonNumber = particle->GetBaryonNumber();
112 AbsProjectileBaryonNumber = std::abs( Proj 112 AbsProjectileBaryonNumber = std::abs( ProjectileBaryonNumber );
113 AbsProjectileCharge = std::abs( G4in 113 AbsProjectileCharge = std::abs( G4int( particle->GetPDGCharge() ) );
114 if ( ProjectileBaryonNumber > 1 ) { 114 if ( ProjectileBaryonNumber > 1 ) {
115 ProjectilePDGcode = 2212; ProjectileabsP 115 ProjectilePDGcode = 2212; ProjectileabsPDGcode = 2212; // Proton
116 } else { 116 } else {
117 ProjectilePDGcode = -2212; Projectileabs 117 ProjectilePDGcode = -2212; ProjectileabsPDGcode = 2212; // Anti-Proton
118 } 118 }
119 ProjectileMass = G4Proton::Proton()->GetP 119 ProjectileMass = G4Proton::Proton()->GetPDGMass();
120 ProjectileMass2 = sqr( ProjectileMass ); 120 ProjectileMass2 = sqr( ProjectileMass );
121 } 121 }
122 122
123 G4double TargetMass = G4Proton::Proton()->G 123 G4double TargetMass = G4Proton::Proton()->GetPDGMass();
124 G4double TargetMass2 = TargetMass * TargetMa 124 G4double TargetMass2 = TargetMass * TargetMass;
125 125
126 G4double Plab = PlabPerParticle; 126 G4double Plab = PlabPerParticle;
127 G4double Elab = std::sqrt( Plab*Plab + Proje 127 G4double Elab = std::sqrt( Plab*Plab + ProjectileMass2 );
128 G4double KineticEnergy = Elab - ProjectileMa 128 G4double KineticEnergy = Elab - ProjectileMass;
129 129
130 G4double S = ProjectileMass2 + TargetMass2 + 130 G4double S = ProjectileMass2 + TargetMass2 + 2.0*TargetMass*Elab;
131 131
132 #ifdef debugFTFparams 132 #ifdef debugFTFparams
133 G4cout << "--------- FTF Parameters -------- 133 G4cout << "--------- FTF Parameters --------------" << G4endl << "Proj Plab "
134 << ProjectilePDGcode << " " << Plab < 134 << ProjectilePDGcode << " " << Plab << G4endl << "Mass KinE " << ProjectileMass
135 << " " << KineticEnergy << G4endl << 135 << " " << KineticEnergy << G4endl << " A Z " << theA << " " << theZ << G4endl;
136 #endif 136 #endif
137 137
138 G4double Ylab, Xtotal( 0.0 ), Xelastic( 0.0 138 G4double Ylab, Xtotal( 0.0 ), Xelastic( 0.0 ), Xannihilation( 0.0 );
139 G4int NumberOfTargetNucleons; 139 G4int NumberOfTargetNucleons;
140 140
141 Ylab = 0.5 * G4Log( (Elab + Plab)/(Elab - Pl 141 Ylab = 0.5 * G4Log( (Elab + Plab)/(Elab - Plab) );
142 142
143 G4double ECMSsqr = S/GeV/GeV; 143 G4double ECMSsqr = S/GeV/GeV;
144 G4double SqrtS = std::sqrt( S )/GeV; 144 G4double SqrtS = std::sqrt( S )/GeV;
145 145
146 #ifdef debugFTFparams 146 #ifdef debugFTFparams
147 G4cout << "Sqrt(s) " << SqrtS << G4endl; 147 G4cout << "Sqrt(s) " << SqrtS << G4endl;
148 #endif 148 #endif
149 149
150 TargetMass /= GeV; TargetMass2 /= (G 150 TargetMass /= GeV; TargetMass2 /= (GeV*GeV);
151 ProjectileMass /= GeV; ProjectileMass2 /= (G 151 ProjectileMass /= GeV; ProjectileMass2 /= (GeV*GeV);
152 152
153 Plab /= GeV; 153 Plab /= GeV;
154 G4double Xftf = 0.0; 154 G4double Xftf = 0.0;
155 155
156 G4int NumberOfTargetProtons = theZ; 156 G4int NumberOfTargetProtons = theZ;
157 G4int NumberOfTargetNeutrons = theA - theZ; 157 G4int NumberOfTargetNeutrons = theA - theZ;
158 NumberOfTargetNucleons = NumberOfTargetProto 158 NumberOfTargetNucleons = NumberOfTargetProtons + NumberOfTargetNeutrons;
159 159
160 // ---------- hadron projectile ------------ 160 // ---------- hadron projectile ----------------
161 if ( AbsProjectileBaryonNumber <= 1 ) { // 161 if ( AbsProjectileBaryonNumber <= 1 ) { // Projectile is hadron or baryon
162 162
163 // Interaction on P 163 // Interaction on P
164 G4double xTtP = csGGinstance->GetTotalIsot 164 G4double xTtP = csGGinstance->GetTotalIsotopeCrossSection( particle, KineticEnergy, 1, 1);
165 G4double xElP = csGGinstance->GetElasticIs 165 G4double xElP = csGGinstance->GetElasticIsotopeCrossSection(particle, KineticEnergy, 1, 1);
166 166
167 // Interaction on N 167 // Interaction on N
168 G4double xTtN = csGGinstance->GetTotalIsot 168 G4double xTtN = csGGinstance->GetTotalIsotopeCrossSection( particle, KineticEnergy, 0, 1);
169 G4double xElN = csGGinstance->GetElasticIs 169 G4double xElN = csGGinstance->GetElasticIsotopeCrossSection(particle, KineticEnergy, 0, 1);
170 170
171 // Average properties of h+N interactions 171 // Average properties of h+N interactions
172 Xtotal = ( NumberOfTargetProtons * xTtP 172 Xtotal = ( NumberOfTargetProtons * xTtP + NumberOfTargetNeutrons * xTtN ) / NumberOfTargetNucleons;
173 Xelastic = ( NumberOfTargetProtons * xElP 173 Xelastic = ( NumberOfTargetProtons * xElP + NumberOfTargetNeutrons * xElN ) / NumberOfTargetNucleons;
174 Xannihilation = 0.0; 174 Xannihilation = 0.0;
175 175
176 Xtotal /= millibarn; 176 Xtotal /= millibarn;
177 Xelastic /= millibarn; 177 Xelastic /= millibarn;
178 178
179 #ifdef debugFTFparams 179 #ifdef debugFTFparams
180 G4cout<<"Estimated cross sections (total a 180 G4cout<<"Estimated cross sections (total and elastic) of h+N interactions "<<Xtotal<<" "<<Xelastic<<" (mb)"<<G4endl;
181 #endif 181 #endif
182 } 182 }
183 183
184 // ---------- nucleus projectile ----------- 184 // ---------- nucleus projectile ----------------
185 if ( ProjectileIsNucleus && ProjectileBary 185 if ( ProjectileIsNucleus && ProjectileBaryonNumber > 1 ) {
186 186
187 #ifdef debugFTFparams 187 #ifdef debugFTFparams
188 G4cout<<"Projectile is a nucleus: A and Z 188 G4cout<<"Projectile is a nucleus: A and Z - "<<ProjectileBaryonNumber<<" "<<ProjectileCharge<<G4endl;
189 #endif 189 #endif
190 190
191 const G4ParticleDefinition* Proton = G4Pro 191 const G4ParticleDefinition* Proton = G4Proton::Proton();
192 // Interaction on P 192 // Interaction on P
193 G4double XtotPP = csGGinstance->GetTotalIs 193 G4double XtotPP = csGGinstance->GetTotalIsotopeCrossSection(Proton, KineticEnergy, 1, 1);
194 G4double XelPP = csGGinstance->GetElastic 194 G4double XelPP = csGGinstance->GetElasticIsotopeCrossSection(Proton, KineticEnergy, 1, 1);
195 195
196 const G4ParticleDefinition* Neutron = G4Ne 196 const G4ParticleDefinition* Neutron = G4Neutron::Neutron();
197 // Interaction on N 197 // Interaction on N
198 G4double XtotPN = csGGinstance->GetTotalIs 198 G4double XtotPN = csGGinstance->GetTotalIsotopeCrossSection(Neutron, KineticEnergy, 0, 1);
199 G4double XelPN = csGGinstance->GetElastic 199 G4double XelPN = csGGinstance->GetElasticIsotopeCrossSection(Neutron, KineticEnergy, 0, 1);
200 200
201 #ifdef debugFTFparams 201 #ifdef debugFTFparams
202 G4cout << "XsPP (total and elastic) " << X 202 G4cout << "XsPP (total and elastic) " << XtotPP/millibarn << " " << XelPP/millibarn <<" (mb)"<< G4endl
203 << "XsPN (total and elastic) " << X 203 << "XsPN (total and elastic) " << XtotPN/millibarn << " " << XelPN/millibarn <<" (mb)"<< G4endl;
204 #endif 204 #endif
205 205
206 Xtotal = ( 206 Xtotal = (
207 AbsProjectileCharge * Number 207 AbsProjectileCharge * NumberOfTargetProtons * XtotPP +
208 ( AbsProjectileBaryonNumber 208 ( AbsProjectileBaryonNumber - AbsProjectileCharge ) *
209 NumberOfTargetNeutrons * 209 NumberOfTargetNeutrons * XtotPP
210 + 210 +
211 ( AbsProjectileCharge * Numb 211 ( AbsProjectileCharge * NumberOfTargetNeutrons +
212 ( AbsProjectileBaryonNumbe 212 ( AbsProjectileBaryonNumber - AbsProjectileCharge ) *
213 NumberOfTargetProtons 213 NumberOfTargetProtons ) * XtotPN
214 ) / ( AbsProjectileBaryonNumbe 214 ) / ( AbsProjectileBaryonNumber * NumberOfTargetNucleons );
215 Xelastic= ( 215 Xelastic= (
216 AbsProjectileCharge * Number 216 AbsProjectileCharge * NumberOfTargetProtons * XelPP +
217 ( AbsProjectileBaryonNumber 217 ( AbsProjectileBaryonNumber - AbsProjectileCharge ) *
218 NumberOfTargetNeutrons * 218 NumberOfTargetNeutrons * XelPP
219 + 219 +
220 ( AbsProjectileCharge * Numb 220 ( AbsProjectileCharge * NumberOfTargetNeutrons +
221 ( AbsProjectileBaryonNumbe 221 ( AbsProjectileBaryonNumber - AbsProjectileCharge ) *
222 NumberOfTargetProtons 222 NumberOfTargetProtons ) * XelPN
223 ) / ( AbsProjectileBaryonNumbe 223 ) / ( AbsProjectileBaryonNumber * NumberOfTargetNucleons );
224 224
225 Xannihilation = 0.0; 225 Xannihilation = 0.0;
226 Xtotal /= millibarn; 226 Xtotal /= millibarn;
227 Xelastic /= millibarn; 227 Xelastic /= millibarn;
228 } 228 }
229 229
230 // ---------- The projectile is anti-baryon 230 // ---------- The projectile is anti-baryon or anti-nucleus ----------------
231 // anti Sigma^0_c 231 // anti Sigma^0_c anti Delta^-
232 if ( ProjectilePDGcode >= -4112 && Project 232 if ( ProjectilePDGcode >= -4112 && ProjectilePDGcode <= -1114 ) {
233 // Only non-strange and strange baryons ar 233 // Only non-strange and strange baryons are considered
234 234
235 #ifdef debugFTFparams 235 #ifdef debugFTFparams
236 G4cout<<"Projectile is a anti-baryon or an 236 G4cout<<"Projectile is a anti-baryon or anti-nucleus - "<<ProjectileBaryonNumber<<" "<<ProjectileCharge<<G4endl;
237 G4cout<<"(Only non-strange and strange bar 237 G4cout<<"(Only non-strange and strange baryons are considered)"<<G4endl;
238 #endif 238 #endif
239 239
240 G4double X_a( 0.0 ), X_b( 0.0 ), X_c( 0.0 240 G4double X_a( 0.0 ), X_b( 0.0 ), X_c( 0.0 ), X_d( 0.0 );
241 G4double MesonProdThreshold = ProjectileMa 241 G4double MesonProdThreshold = ProjectileMass + TargetMass +
242 ( 2.0 * 0.14 242 ( 2.0 * 0.14 + 0.016 ); // 2 Mpi + DeltaE;
243 243
244 if ( PlabPerParticle < 40.0*MeV ) { // Low 244 if ( PlabPerParticle < 40.0*MeV ) { // Low energy limits. Projectile at rest.
245 Xtotal = 1512.9; // mb 245 Xtotal = 1512.9; // mb
246 Xelastic = 473.2; // mb 246 Xelastic = 473.2; // mb
247 X_a = 625.1; // mb 247 X_a = 625.1; // mb
248 X_b = 9.780; // mb 248 X_b = 9.780; // mb
249 X_c = 49.989; // mb 249 X_c = 49.989; // mb
250 X_d = 6.614; // mb 250 X_d = 6.614; // mb
251 } else { // Total and elastic cross sectio 251 } else { // Total and elastic cross section of PbarP interactions a'la Arkhipov
252 G4double LogS = G4Log( ECMSsqr / 33.0625 252 G4double LogS = G4Log( ECMSsqr / 33.0625 );
253 G4double Xasmpt = 36.04 + 0.304*LogS*Log 253 G4double Xasmpt = 36.04 + 0.304*LogS*LogS; // mb
254 LogS = G4Log( SqrtS / 20.74 ); 254 LogS = G4Log( SqrtS / 20.74 );
255 G4double Basmpt = 11.92 + 0.3036*LogS*Lo 255 G4double Basmpt = 11.92 + 0.3036*LogS*LogS; // GeV^(-2)
256 G4double R0 = std::sqrt( 0.40874044*Xasm 256 G4double R0 = std::sqrt( 0.40874044*Xasmpt - Basmpt ); // GeV^(-1)
257 257
258 G4double FlowF = SqrtS / std::sqrt( ECMS 258 G4double FlowF = SqrtS / std::sqrt( ECMSsqr*ECMSsqr + ProjectileMass2*ProjectileMass2 +
259 Targ 259 TargetMass2*TargetMass2 - 2.0*ECMSsqr*ProjectileMass2
260 - 2. 260 - 2.0*ECMSsqr*TargetMass2
261 - 2. 261 - 2.0*ProjectileMass2*TargetMass2 );
262 262
263 Xtotal = Xasmpt * ( 1.0 + 13.55*FlowF/R0 263 Xtotal = Xasmpt * ( 1.0 + 13.55*FlowF/R0/R0/R0*
264 (1.0 - 4.47/Sq 264 (1.0 - 4.47/SqrtS + 12.38/ECMSsqr - 12.43/SqrtS/ECMSsqr) ); // mb
265 265
266 Xasmpt = 4.4 + 0.101*LogS*LogS; // mb 266 Xasmpt = 4.4 + 0.101*LogS*LogS; // mb
267 Xelastic = Xasmpt * ( 1.0 + 59.27*FlowF/ 267 Xelastic = Xasmpt * ( 1.0 + 59.27*FlowF/R0/R0/R0*
268 (1.0 - 6.95/ 268 (1.0 - 6.95/SqrtS + 23.54/ECMSsqr - 25.34/SqrtS/ECMSsqr ) ); // mb
269 269
270 X_a = 25.0*FlowF; // mb, 3-shirts diagr 270 X_a = 25.0*FlowF; // mb, 3-shirts diagram
271 271
272 if ( SqrtS < MesonProdThreshold ) { 272 if ( SqrtS < MesonProdThreshold ) {
273 X_b = 3.13 + 140.0*G4Pow::GetInstance( 273 X_b = 3.13 + 140.0*G4Pow::GetInstance()->powA( MesonProdThreshold - SqrtS, 2.5 ); // mb anti-quark-quark annihilation
274 Xelastic -= 3.0*X_b; // Xel-X(PbarP-> 274 Xelastic -= 3.0*X_b; // Xel-X(PbarP->NNbar)
275 } else { 275 } else {
276 X_b = 6.8/SqrtS; // mb anti-quark-qua 276 X_b = 6.8/SqrtS; // mb anti-quark-quark annihilation
277 Xelastic -= 3.0*X_b; // Xel-X(PbarP-> 277 Xelastic -= 3.0*X_b; // Xel-X(PbarP->NNbar)
278 } 278 }
279 279
280 X_c = 2.0*FlowF*sqr( ProjectileMass + Ta 280 X_c = 2.0*FlowF*sqr( ProjectileMass + TargetMass )/ECMSsqr; // mb rearrangement
281 281
282 X_d = 23.3/ECMSsqr; // mb anti-quark-qu 282 X_d = 23.3/ECMSsqr; // mb anti-quark-quark string creation
283 } 283 }
284 284
285 G4double Xann_on_P( 0.0), Xann_on_N( 0.0 ) 285 G4double Xann_on_P( 0.0), Xann_on_N( 0.0 );
286 286
287 if ( ProjectilePDGcode == -2212 ) { 287 if ( ProjectilePDGcode == -2212 ) { // Pbar+P/N
288 Xann_on_P = X_a + X_b*5.0 + X_c*5.0 + X_ 288 Xann_on_P = X_a + X_b*5.0 + X_c*5.0 + X_d*6.0;
289 Xann_on_N = X_a + X_b*4.0 + X_c*4.0 + X_ 289 Xann_on_N = X_a + X_b*4.0 + X_c*4.0 + X_d*4.0;
290 } else if ( ProjectilePDGcode == -2112 ) { 290 } else if ( ProjectilePDGcode == -2112 ) { // NeutrBar+P/N
291 Xann_on_P = X_a + X_b*4.0 + X_c*4.0 + X_ 291 Xann_on_P = X_a + X_b*4.0 + X_c*4.0 + X_d*4.0;
292 Xann_on_N = X_a + X_b*5.0 + X_c*5.0 + X_ 292 Xann_on_N = X_a + X_b*5.0 + X_c*5.0 + X_d*6.0;
293 } else if ( ProjectilePDGcode == -3122 ) { 293 } else if ( ProjectilePDGcode == -3122 ) { // LambdaBar+P/N
294 Xann_on_P = X_a + X_b*3.0 + X_c*3.0 + X_ 294 Xann_on_P = X_a + X_b*3.0 + X_c*3.0 + X_d*2.0;
295 Xann_on_N = X_a + X_b*3.0 + X_c*3.0 + X_ 295 Xann_on_N = X_a + X_b*3.0 + X_c*3.0 + X_d*2.0;
296 } else if ( ProjectilePDGcode == -3112 ) { 296 } else if ( ProjectilePDGcode == -3112 ) { // Sigma-Bar+P/N
297 Xann_on_P = X_a + X_b*2.0 + X_c*2.0 + X_ 297 Xann_on_P = X_a + X_b*2.0 + X_c*2.0 + X_d*0.0;
298 Xann_on_N = X_a + X_b*4.0 + X_c*4.0 + X_ 298 Xann_on_N = X_a + X_b*4.0 + X_c*4.0 + X_d*2.0;
299 } else if ( ProjectilePDGcode == -3212 ) { 299 } else if ( ProjectilePDGcode == -3212 ) { // Sigma0Bar+P/N
300 Xann_on_P = X_a + X_b*3.0 + X_c*3.0 + X_ 300 Xann_on_P = X_a + X_b*3.0 + X_c*3.0 + X_d*2.0;
301 Xann_on_N = X_a + X_b*3.0 + X_c*3.0 + X_ 301 Xann_on_N = X_a + X_b*3.0 + X_c*3.0 + X_d*2.0;
302 } else if ( ProjectilePDGcode == -3222 ) { 302 } else if ( ProjectilePDGcode == -3222 ) { // Sigma+Bar+P/N
303 Xann_on_P = X_a + X_b*4.0 + X_c*4.0 + X_ 303 Xann_on_P = X_a + X_b*4.0 + X_c*4.0 + X_d*2.0;
304 Xann_on_N = X_a + X_b*2.0 + X_c*2.0 + X_ 304 Xann_on_N = X_a + X_b*2.0 + X_c*2.0 + X_d*0.0;
305 } else if ( ProjectilePDGcode == -3312 ) { 305 } else if ( ProjectilePDGcode == -3312 ) { // Xi-Bar+P/N
306 Xann_on_P = X_a + X_b*1.0 + X_c*1.0 + X_ 306 Xann_on_P = X_a + X_b*1.0 + X_c*1.0 + X_d*0.0;
307 Xann_on_N = X_a + X_b*2.0 + X_c*2.0 + X_ 307 Xann_on_N = X_a + X_b*2.0 + X_c*2.0 + X_d*0.0;
308 } else if ( ProjectilePDGcode == -3322 ) { 308 } else if ( ProjectilePDGcode == -3322 ) { // Xi0Bar+P/N
309 Xann_on_P = X_a + X_b*2.0 + X_c*2.0 + X_ 309 Xann_on_P = X_a + X_b*2.0 + X_c*2.0 + X_d*0.0;
310 Xann_on_N = X_a + X_b*1.0 + X_c*1.0 + X_ 310 Xann_on_N = X_a + X_b*1.0 + X_c*1.0 + X_d*0.0;
311 } else if ( ProjectilePDGcode == -3334 ) { 311 } else if ( ProjectilePDGcode == -3334 ) { // Omega-Bar+P/N
312 Xann_on_P = X_a + X_b*0.0 + X_c*0.0 + X_ 312 Xann_on_P = X_a + X_b*0.0 + X_c*0.0 + X_d*0.0;
313 Xann_on_N = X_a + X_b*0.0 + X_c*0.0 + X_ 313 Xann_on_N = X_a + X_b*0.0 + X_c*0.0 + X_d*0.0;
314 } else { 314 } else {
315 G4cout << "Unknown anti-baryon for FTF a 315 G4cout << "Unknown anti-baryon for FTF annihilation" << G4endl;
316 } 316 }
317 317
318 //G4cout << "Sum " << Xann_on_P < 318 //G4cout << "Sum " << Xann_on_P << G4endl;
319 319
320 if ( ! ProjectileIsNucleus ) { // Project 320 if ( ! ProjectileIsNucleus ) { // Projectile is anti-baryon
321 Xannihilation = ( NumberOfTargetProtons 321 Xannihilation = ( NumberOfTargetProtons * Xann_on_P + NumberOfTargetNeutrons * Xann_on_N )
322 / NumberOfTargetNucleons 322 / NumberOfTargetNucleons;
323 } else { // Projectile is a nucleus 323 } else { // Projectile is a nucleus
324 Xannihilation = ( 324 Xannihilation = (
325 ( AbsProjectileCharge 325 ( AbsProjectileCharge * NumberOfTargetProtons +
326 ( AbsProjectileBaryo 326 ( AbsProjectileBaryonNumber - AbsProjectileCharge ) *
327 NumberOfTargetNeutro 327 NumberOfTargetNeutrons ) * Xann_on_P
328 + 328 +
329 ( AbsProjectileCharge 329 ( AbsProjectileCharge * NumberOfTargetNeutrons +
330 ( AbsProjectileBaryo 330 ( AbsProjectileBaryonNumber - AbsProjectileCharge ) *
331 NumberOfTargetProton 331 NumberOfTargetProtons ) * Xann_on_N
332 ) / ( AbsProjectileBaryo 332 ) / ( AbsProjectileBaryonNumber * NumberOfTargetNucleons );
333 } 333 }
334 334
335 //G4double Xftf = 0.0; 335 //G4double Xftf = 0.0;
336 MesonProdThreshold = ProjectileMass + Targ 336 MesonProdThreshold = ProjectileMass + TargetMass + (0.14 + 0.08); // Mpi + DeltaE
337 if ( SqrtS > MesonProdThreshold ) { 337 if ( SqrtS > MesonProdThreshold ) {
338 Xftf = 36.0 * ( 1.0 - MesonProdThreshold 338 Xftf = 36.0 * ( 1.0 - MesonProdThreshold/SqrtS );
339 } 339 }
340 340
341 Xtotal = Xelastic + Xannihilation + Xftf; 341 Xtotal = Xelastic + Xannihilation + Xftf;
342 342
343 #ifdef debugFTFparams 343 #ifdef debugFTFparams
344 G4cout << "Plab Xtotal, Xelastic Xinel Xf 344 G4cout << "Plab Xtotal, Xelastic Xinel Xftf " << Plab << " " << Xtotal << " " << Xelastic
345 << " " << Xtotal - Xelastic << " " 345 << " " << Xtotal - Xelastic << " " << Xtotal - Xelastic - Xannihilation << " (mb)"<< G4endl
346 << "Plab Xelastic/Xtotal, Xann/Xin 346 << "Plab Xelastic/Xtotal, Xann/Xin " << Plab << " " << Xelastic/Xtotal << " "
347 << Xannihilation/(Xtotal - Xelastic 347 << Xannihilation/(Xtotal - Xelastic) << G4endl;
348 #endif 348 #endif
349 349
350 } 350 }
351 351
352 if ( Xtotal == 0.0 ) { // Projectile is und 352 if ( Xtotal == 0.0 ) { // Projectile is undefined, Nucleon assumed
353 353
354 const G4ParticleDefinition* Proton = G4Pro 354 const G4ParticleDefinition* Proton = G4Proton::Proton();
355 // Interaction on P 355 // Interaction on P
356 G4double XtotPP = csGGinstance->GetTotalIs 356 G4double XtotPP = csGGinstance->GetTotalIsotopeCrossSection(Proton, KineticEnergy, 1, 1);
357 G4double XelPP = csGGinstance->GetElastic 357 G4double XelPP = csGGinstance->GetElasticIsotopeCrossSection(Proton, KineticEnergy, 1, 1);
358 358
359 // Interaction on N 359 // Interaction on N
360 G4double XtotPN = csGGinstance->GetTotalIs 360 G4double XtotPN = csGGinstance->GetTotalIsotopeCrossSection(Proton, KineticEnergy, 0, 1);
361 G4double XelPN = csGGinstance->GetElastic 361 G4double XelPN = csGGinstance->GetElasticIsotopeCrossSection(Proton, KineticEnergy, 0, 1);
362 362
363 Xtotal = ( NumberOfTargetProtons * Xtot 363 Xtotal = ( NumberOfTargetProtons * XtotPP + NumberOfTargetNeutrons * XtotPN )
364 / NumberOfTargetNucleons; 364 / NumberOfTargetNucleons;
365 Xelastic = ( NumberOfTargetProtons * XelP 365 Xelastic = ( NumberOfTargetProtons * XelPP + NumberOfTargetNeutrons * XelPN )
366 / NumberOfTargetNucleons; 366 / NumberOfTargetNucleons;
367 Xannihilation = 0.0; 367 Xannihilation = 0.0;
368 Xtotal /= millibarn; 368 Xtotal /= millibarn;
369 Xelastic /= millibarn; 369 Xelastic /= millibarn;
370 }; 370 };
371 371
372 // Geometrical parameters 372 // Geometrical parameters
373 SetTotalCrossSection( Xtotal ); 373 SetTotalCrossSection( Xtotal );
374 SetElasticCrossSection( Xelastic ); << 374 SetElastisCrossSection( Xelastic );
375 SetInelasticCrossSection( Xtotal - Xelastic 375 SetInelasticCrossSection( Xtotal - Xelastic );
376 376
377 // Interactions with elastic and inelastic c 377 // Interactions with elastic and inelastic collisions
378 SetProbabilityOfElasticScatt( Xtotal, Xelast 378 SetProbabilityOfElasticScatt( Xtotal, Xelastic );
379 379
380 SetRadiusOfHNinteractions2( Xtotal/pi/10.0 ) 380 SetRadiusOfHNinteractions2( Xtotal/pi/10.0 );
381 381
382 if ( ( Xtotal - Xelastic ) == 0.0 ) { 382 if ( ( Xtotal - Xelastic ) == 0.0 ) {
383 SetProbabilityOfAnnihilation( 0.0 ); 383 SetProbabilityOfAnnihilation( 0.0 );
384 } else { 384 } else {
385 SetProbabilityOfAnnihilation( Xannihilatio 385 SetProbabilityOfAnnihilation( Xannihilation / (Xtotal - Xelastic) );
386 } 386 }
387 387
388 if(Xelastic > 0.0) { 388 if(Xelastic > 0.0) {
389 SetSlope( Xtotal*Xtotal/16.0/pi/Xelastic/0 389 SetSlope( Xtotal*Xtotal/16.0/pi/Xelastic/0.3894 );// Slope parameter of elastic scattering
390 390 // (GeV/c)^(-2))
391 // Parameters of elastic scattering 391 // Parameters of elastic scattering
392 // Gaussian parametrization of elastic sca 392 // Gaussian parametrization of elastic scattering amplitude assumed
393 SetAvaragePt2ofElasticScattering( 1.0/( Xt 393 SetAvaragePt2ofElasticScattering( 1.0/( Xtotal*Xtotal/16.0/pi/Xelastic/0.3894 )*GeV*GeV );
394 } else { 394 } else {
395 SetSlope(1.0); 395 SetSlope(1.0);
396 SetAvaragePt2ofElasticScattering( 0.0); 396 SetAvaragePt2ofElasticScattering( 0.0);
397 } 397 }
398 SetGamma0( GetSlope()*Xtotal/10.0/2.0/pi ); 398 SetGamma0( GetSlope()*Xtotal/10.0/2.0/pi );
399 399
400 G4double Xinel = Xtotal - Xelastic; 400 G4double Xinel = Xtotal - Xelastic;
401 401
402 #ifdef debugFTFparams 402 #ifdef debugFTFparams
403 G4cout<< "Slope of hN elastic scattering" << 403 G4cout<< "Slope of hN elastic scattering" << GetSlope() << G4endl;
404 G4cout << "AvaragePt2ofElasticScattering " < 404 G4cout << "AvaragePt2ofElasticScattering " << GetAvaragePt2ofElasticScattering() << G4endl;
405 G4cout<<"Parameters of excitation for projec 405 G4cout<<"Parameters of excitation for projectile "<<ProjectilePDGcode<< G4endl;
406 #endif 406 #endif
407 407
408 if ( (ProjectilePDGcode == 2212) || (Project 408 if ( (ProjectilePDGcode == 2212) || (ProjectilePDGcode == 2112) ) { // Projectile is proton or neutron
409 409
410 const G4int indexTune = G4FTFTunings::Inst 410 const G4int indexTune = G4FTFTunings::Instance()->GetIndexTune( particle, KineticEnergy );
411 411
412 // A process probability is parameterized 412 // A process probability is parameterized as Prob = A_1*exp(-A_2*y) + A_3*exp(-A_4*y) + A_top
413 // y is a rapidity of a partcle in the tar 413 // y is a rapidity of a partcle in the target nucleus. Ymin is a minimal rapidity below it X=0
414 414
415 // Proc# A1 B1 A2 415 // Proc# A1 B1 A2 B2 A3 Atop Ymin
416 /* original hadr-string-diff-V10-03-07 (si 416 /* original hadr-string-diff-V10-03-07 (similar to 10.3.x)
417 SetParams( 0, 13.71, 1.75, -2 417 SetParams( 0, 13.71, 1.75, -214.5, 4.25, 0.0, 0.5 , 1.1 ); // Qexchange without Exc.
418 SetParams( 1, 25.0, 1.0, -5 418 SetParams( 1, 25.0, 1.0, -50.34, 1.5 , 0.0, 0.0 , 1.4 ); // Qexchange with Exc.
419 SetParams( 2, 6.0/Xinel, 0.0 ,-6.0/Xinel*1 419 SetParams( 2, 6.0/Xinel, 0.0 ,-6.0/Xinel*16.28, 3.0 , 0.0, 0.0 , 0.93); // Projectile diffraction
420 SetParams( 3, 6.0/Xinel, 0.0 ,-6.0/Xinel*1 420 SetParams( 3, 6.0/Xinel, 0.0 ,-6.0/Xinel*16.28, 3.0 , 0.0, 0.0 , 0.93); // Target diffraction
421 SetParams( 4, 1.0, 0.0 , -2 421 SetParams( 4, 1.0, 0.0 , -2.01 , 0.5 , 0.0, 0.0 , 1.4 ); // Qexchange with Exc. Additional multiplier
422 */ 422 */
423 423
424 // Proc# 424 // Proc#
425 SetParams( 0, fArrayParCollBaryonProj[inde 425 SetParams( 0, fArrayParCollBaryonProj[indexTune].GetProc0A1(),
426 fArrayParCollBaryonProj[indexTune] 426 fArrayParCollBaryonProj[indexTune].GetProc0B1(),
427 fArrayParCollBaryonProj[indexTune].GetPr 427 fArrayParCollBaryonProj[indexTune].GetProc0A2(),
428 fArrayParCollBaryonProj[indexTune] 428 fArrayParCollBaryonProj[indexTune].GetProc0B2(),
429 fArrayParCollBaryonProj[indexTune].GetPr 429 fArrayParCollBaryonProj[indexTune].GetProc0A3(),
430 fArrayParCollBaryonProj[indexTune] 430 fArrayParCollBaryonProj[indexTune].GetProc0Atop(),
431 fArrayParCollBaryonProj[indexTune].GetPr 431 fArrayParCollBaryonProj[indexTune].GetProc0Ymin() ); // Qexchange without Exc.
432 SetParams( 1, fArrayParCollBaryonProj[inde 432 SetParams( 1, fArrayParCollBaryonProj[indexTune].GetProc1A1(),
433 fArrayParCollBaryonProj[indexTune] 433 fArrayParCollBaryonProj[indexTune].GetProc1B1(),
434 fArrayParCollBaryonProj[indexTune].GetPr 434 fArrayParCollBaryonProj[indexTune].GetProc1A2(),
435 fArrayParCollBaryonProj[indexTune] 435 fArrayParCollBaryonProj[indexTune].GetProc1B2(),
436 fArrayParCollBaryonProj[indexTune].GetPr 436 fArrayParCollBaryonProj[indexTune].GetProc1A3(),
437 fArrayParCollBaryonProj[indexTune] 437 fArrayParCollBaryonProj[indexTune].GetProc1Atop(),
438 fArrayParCollBaryonProj[indexTune].GetPr 438 fArrayParCollBaryonProj[indexTune].GetProc1Ymin() ); // Qexchange with Exc.
439 if ( Xinel > 0.0 ) { 439 if ( Xinel > 0.0 ) {
440 SetParams( 2, 6.0/Xinel, 0.0, -6.0/Xinel 440 SetParams( 2, 6.0/Xinel, 0.0, -6.0/Xinel*16.28, 3.0, 0.0, 0.0, 0.93 ); // Projectile diffraction
441 SetParams( 3, 6.0/Xinel, 0.0, -6.0/Xinel 441 SetParams( 3, 6.0/Xinel, 0.0, -6.0/Xinel*16.28, 3.0, 0.0, 0.0, 0.93 ); // Target diffraction
442 442
443 SetParams( 4, fArrayParCollBaryonProj[in 443 SetParams( 4, fArrayParCollBaryonProj[indexTune].GetProc4A1(),
444 fArrayParCollBaryonProj[indexTune].Get 444 fArrayParCollBaryonProj[indexTune].GetProc4B1(),
445 fArrayParCollBaryonProj[indexTune].Get 445 fArrayParCollBaryonProj[indexTune].GetProc4A2(),
446 fArrayParCollBaryonProj[indexTune].Get 446 fArrayParCollBaryonProj[indexTune].GetProc4B2(),
447 fArrayParCollBaryonProj[indexTune].Get 447 fArrayParCollBaryonProj[indexTune].GetProc4A3(),
448 fArrayParCollBaryonProj[indexTune].Get 448 fArrayParCollBaryonProj[indexTune].GetProc4Atop(),
449 fArrayParCollBaryonProj[indexTune].Get 449 fArrayParCollBaryonProj[indexTune].GetProc4Ymin() ); // Qexchange with Exc. Additional multiplier
450 } else { // if Xinel=0., zero everything 450 } else { // if Xinel=0., zero everything out (obviously)
451 SetParams( 2, 0.0, 0.0, 0.0, 0.0, 0.0, 0 451 SetParams( 2, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 );
452 SetParams( 3, 0.0, 0.0, 0.0, 0.0, 0.0, 0 452 SetParams( 3, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 );
453 SetParams( 4, 0.0, 0.0, 0.0, 0.0, 0.0, 0 453 SetParams( 4, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 );
454 } 454 }
455 455
456 if ( (AbsProjectileBaryonNumber > 10 || Nu 456 if ( (AbsProjectileBaryonNumber > 10 || NumberOfTargetNucleons > 10) && !EnableDiffDissociationForBGreater10 ) {
457 // It is not decided what to do with dif 457 // It is not decided what to do with diffraction dissociation in Had-Nucl and Nucl-Nucl interactions
458 // For the moment both ProjDiffDisso & T 458 // For the moment both ProjDiffDisso & TgtDiffDisso for A > 10 are set to false,
459 // so both projectile and target diffrac 459 // so both projectile and target diffraction are turned OFF
460 if ( ! fArrayParCollBaryonProj[indexTune 460 if ( ! fArrayParCollBaryonProj[indexTune].IsProjDiffDissociation() )
461 SetParams( 2, 0.0, 0.0, 0.0, 0.0, 0.0 461 SetParams( 2, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, -100.0 ); // Projectile diffraction
462 if ( ! fArrayParCollBaryonProj[indexTune 462 if ( ! fArrayParCollBaryonProj[indexTune].IsTgtDiffDissociation() )
463 SetParams( 3, 0.0, 0.0, 0.0, 0.0, 0.0 463 SetParams( 3, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, -100.0 ); // Target diffraction
464 } 464 }
465 465
466 SetDeltaProbAtQuarkExchange( fArrayParColl 466 SetDeltaProbAtQuarkExchange( fArrayParCollBaryonProj[indexTune].GetDeltaProbAtQuarkExchange() );
467 467
468 if ( NumberOfTargetNucleons > 26 ) { 468 if ( NumberOfTargetNucleons > 26 ) {
469 SetProbOfSameQuarkExchange( 1.0 ); 469 SetProbOfSameQuarkExchange( 1.0 );
470 } else { 470 } else {
471 SetProbOfSameQuarkExchange( fArrayParCol 471 SetProbOfSameQuarkExchange( fArrayParCollBaryonProj[indexTune].GetProbOfSameQuarkExchange() );
472 } 472 }
473 473
474 SetProjMinDiffMass( fArrayParCollBaryon 474 SetProjMinDiffMass( fArrayParCollBaryonProj[indexTune].GetProjMinDiffMass() ); // GeV
475 SetProjMinNonDiffMass( fArrayParCollBaryon 475 SetProjMinNonDiffMass( fArrayParCollBaryonProj[indexTune].GetProjMinNonDiffMass() ); // GeV
476 476
477 SetTarMinDiffMass( fArrayParCollBaryon 477 SetTarMinDiffMass( fArrayParCollBaryonProj[indexTune].GetTgtMinDiffMass() ); // GeV
478 SetTarMinNonDiffMass( fArrayParCollBaryon 478 SetTarMinNonDiffMass( fArrayParCollBaryonProj[indexTune].GetTgtMinNonDiffMass() ); // GeV
479 479
480 SetAveragePt2( fArrayParCollBaryon 480 SetAveragePt2( fArrayParCollBaryonProj[indexTune].GetAveragePt2() ); // GeV^2
481 SetProbLogDistrPrD( fArrayParCollBaryon 481 SetProbLogDistrPrD( fArrayParCollBaryonProj[indexTune].GetProbLogDistrPrD() );
482 SetProbLogDistr( fArrayParCollBaryon 482 SetProbLogDistr( fArrayParCollBaryonProj[indexTune].GetProbLogDistr() );
483 483
484 } else if ( ProjectilePDGcode == -2212 || 484 } else if ( ProjectilePDGcode == -2212 || ProjectilePDGcode == -2112 ) { // Projectile is anti_proton or anti_neutron
485 485
486 // Below, in the call to the G4FTFTunings: 486 // Below, in the call to the G4FTFTunings::GetIndexTune method, we pass the proton
487 // as projectile, instead of the real one, 487 // as projectile, instead of the real one, because for switching on/off diffraction
488 // we assume the same treatment for anti_p 488 // we assume the same treatment for anti_proton/anti_neutron as for proton/neutron,
489 // whereas all other parameters for anti_p 489 // whereas all other parameters for anti_proton/anti_neutron are hardwired.
490 const G4int indexTune = G4FTFTunings::Inst 490 const G4int indexTune = G4FTFTunings::Instance()->GetIndexTune( G4Proton::Definition(), KineticEnergy );
491 491
492 // Proc# A1 B1 A2 492 // Proc# A1 B1 A2 B2 A3 Atop Ymin
493 SetParams( 0, 0.0 , 0.0 , 0 493 SetParams( 0, 0.0 , 0.0 , 0.0 , 0.0 , 0.0, 0.0 , 1000.0 ); // Qexchange without Exc.
494 SetParams( 1, 0.0 , 0.0 , 0 494 SetParams( 1, 0.0 , 0.0 , 0.0 , 0.0 , 0.0, 0.0 , 1000.0 ); // Qexchange with Exc.
495 if ( Xinel > 0.) { 495 if ( Xinel > 0.) {
496 SetParams( 2, 6.0/Xinel, 0.0 ,-6.0/Xinel 496 SetParams( 2, 6.0/Xinel, 0.0 ,-6.0/Xinel*16.28, 3.0 , 0.0, 0.0 , 0.93 ); // Projectile diffraction
497 SetParams( 3, 6.0/Xinel, 0.0 ,-6.0/Xinel 497 SetParams( 3, 6.0/Xinel, 0.0 ,-6.0/Xinel*16.28, 3.0 , 0.0, 0.0 , 0.93 ); // Target diffraction
498 SetParams( 4, 1.0, 0.0 , 498 SetParams( 4, 1.0, 0.0 , 0.0, 0.0 , 0.0, 0.0 , 0.93 ); // Qexchange with Exc. Additional multiply
499 } else { 499 } else {
500 SetParams( 2, 0.0, 0.0, 0.0, 0.0, 0.0, 0 500 SetParams( 2, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 );
501 SetParams( 3, 0.0, 0.0, 0.0, 0.0, 0.0, 0 501 SetParams( 3, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 );
502 SetParams( 4, 0.0, 0.0, 0.0, 0.0, 0.0, 0 502 SetParams( 4, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 );
503 } 503 }
504 504
505 if ( AbsProjectileBaryonNumber > 10 || N 505 if ( AbsProjectileBaryonNumber > 10 || NumberOfTargetNucleons > 10 ) {
506 // It is not decided what to do with dif 506 // It is not decided what to do with diffraction dissociation in Had-Nucl and Nucl-Nucl interactions
507 // For the moment both ProjDiffDisso & T 507 // For the moment both ProjDiffDisso & TgtDiffDisso are set to false,
508 // so both projectile and target diffrac 508 // so both projectile and target diffraction are turned OFF
509 if ( ! fArrayParCollBaryonProj[indexTune 509 if ( ! fArrayParCollBaryonProj[indexTune].IsProjDiffDissociation() )
510 SetParams( 2, 0.0, 0.0 , 510 SetParams( 2, 0.0, 0.0 , 0.0 , 0.0 , 0.0, 0.0 , -100.0 ); // Projectile diffraction
511 if ( ! fArrayParCollBaryonProj[indexTune 511 if ( ! fArrayParCollBaryonProj[indexTune].IsTgtDiffDissociation() )
512 SetParams( 3, 0.0, 0.0 , 512 SetParams( 3, 0.0, 0.0 , 0.0 , 0.0 , 0.0, 0.0 , -100.0 ); // Target diffraction
513 } 513 }
514 514
515 SetDeltaProbAtQuarkExchange( 0.0 ); 515 SetDeltaProbAtQuarkExchange( 0.0 );
516 SetProbOfSameQuarkExchange( 0.0 ); 516 SetProbOfSameQuarkExchange( 0.0 );
517 SetProjMinDiffMass( ProjectileMass + 0.22 517 SetProjMinDiffMass( ProjectileMass + 0.22 ); // GeV
518 SetProjMinNonDiffMass( ProjectileMass + 0. 518 SetProjMinNonDiffMass( ProjectileMass + 0.22 ); // GeV
519 SetTarMinDiffMass( TargetMass + 0.22 ); 519 SetTarMinDiffMass( TargetMass + 0.22 ); // GeV
520 SetTarMinNonDiffMass( TargetMass + 0.22 ); 520 SetTarMinNonDiffMass( TargetMass + 0.22 ); // GeV
521 SetAveragePt2( 0.3 ); 521 SetAveragePt2( 0.3 ); // GeV^2
522 SetProbLogDistrPrD( 0.55 ); 522 SetProbLogDistrPrD( 0.55 );
523 SetProbLogDistr( 0.55 ); 523 SetProbLogDistr( 0.55 );
524 524
525 } else if ( ProjectileabsPDGcode == 211 || 525 } else if ( ProjectileabsPDGcode == 211 || ProjectilePDGcode == 111 ) { // Projectile is Pion
526 526
527 const G4int indexTune = G4FTFTunings::Inst 527 const G4int indexTune = G4FTFTunings::Instance()->GetIndexTune( particle, KineticEnergy );
528 528
529 // Proc# A1 B1 A2 529 // Proc# A1 B1 A2 B2 A3 Atop Ymin
530 /* --> original code 530 /* --> original code
531 SetParams( 0, 150.0, 1.8 , -247. 531 SetParams( 0, 150.0, 1.8 , -247.3, 2.3, 0., 1. , 2.3 );
532 SetParams( 1, 5.77, 0.6 , -5.7 532 SetParams( 1, 5.77, 0.6 , -5.77, 0.8, 0., 0. , 0.0 );
533 SetParams( 2, 2.27, 0.5 , -98052. 533 SetParams( 2, 2.27, 0.5 , -98052.0, 4.0, 0., 0. , 3.0 );
534 SetParams( 3, 7.0, 0.9, -85.2 534 SetParams( 3, 7.0, 0.9, -85.28, 1.9, 0.08, 0. , 2.2 );
535 SetParams( 4, 1.0, 0.0 , -11.0 535 SetParams( 4, 1.0, 0.0 , -11.02, 1.0, 0.0, 0. , 2.4 ); // Qexchange with Exc. Additional multiply
536 */ 536 */
537 // Proc# 537 // Proc#
538 SetParams( 0, fArrayParCollPionProj[indexT 538 SetParams( 0, fArrayParCollPionProj[indexTune].GetProc0A1(),
539 fArrayParCollPionProj[indexTune].G 539 fArrayParCollPionProj[indexTune].GetProc0B1(),
540 fArrayParCollPionProj[indexTune].GetProc 540 fArrayParCollPionProj[indexTune].GetProc0A2(),
541 fArrayParCollPionProj[indexTune].G 541 fArrayParCollPionProj[indexTune].GetProc0B2(),
542 fArrayParCollPionProj[indexTune].GetProc 542 fArrayParCollPionProj[indexTune].GetProc0A3(),
543 fArrayParCollPionProj[indexTune].G 543 fArrayParCollPionProj[indexTune].GetProc0Atop(),
544 fArrayParCollPionProj[indexTune].GetProc 544 fArrayParCollPionProj[indexTune].GetProc0Ymin() ); // Qexchange without Exc.
545 SetParams( 1, fArrayParCollPionProj[indexT 545 SetParams( 1, fArrayParCollPionProj[indexTune].GetProc1A1(),
546 fArrayParCollPionProj[indexTune].G 546 fArrayParCollPionProj[indexTune].GetProc1B1(),
547 fArrayParCollPionProj[indexTune].GetProc 547 fArrayParCollPionProj[indexTune].GetProc1A2(),
548 fArrayParCollPionProj[indexTune].G 548 fArrayParCollPionProj[indexTune].GetProc1B2(),
549 fArrayParCollPionProj[indexTune].GetProc 549 fArrayParCollPionProj[indexTune].GetProc1A3(),
550 fArrayParCollPionProj[indexTune].G 550 fArrayParCollPionProj[indexTune].GetProc1Atop(),
551 fArrayParCollPionProj[indexTune].GetProc 551 fArrayParCollPionProj[indexTune].GetProc1Ymin() ); // Qexchange with Exc.
552 SetParams( 2, fArrayParCollPionProj[indexT 552 SetParams( 2, fArrayParCollPionProj[indexTune].GetProc2A1(),
553 fArrayParCollPionProj[indexTune].G 553 fArrayParCollPionProj[indexTune].GetProc2B1(),
554 fArrayParCollPionProj[indexTune].GetProc 554 fArrayParCollPionProj[indexTune].GetProc2A2(),
555 fArrayParCollPionProj[indexTune].G 555 fArrayParCollPionProj[indexTune].GetProc2B2(),
556 fArrayParCollPionProj[indexTune].GetProc 556 fArrayParCollPionProj[indexTune].GetProc2A3(),
557 fArrayParCollPionProj[indexTune].G 557 fArrayParCollPionProj[indexTune].GetProc2Atop(),
558 fArrayParCollPionProj[indexTune].GetProc 558 fArrayParCollPionProj[indexTune].GetProc2Ymin() ); // Projectile diffraction
559 SetParams( 3, fArrayParCollPionProj[indexT 559 SetParams( 3, fArrayParCollPionProj[indexTune].GetProc3A1(),
560 fArrayParCollPionProj[indexTune].G 560 fArrayParCollPionProj[indexTune].GetProc3B1(),
561 fArrayParCollPionProj[indexTune].GetProc 561 fArrayParCollPionProj[indexTune].GetProc3A2(),
562 fArrayParCollPionProj[indexTune].G 562 fArrayParCollPionProj[indexTune].GetProc3B2(),
563 fArrayParCollPionProj[indexTune].GetProc 563 fArrayParCollPionProj[indexTune].GetProc3A3(),
564 fArrayParCollPionProj[indexTune].G 564 fArrayParCollPionProj[indexTune].GetProc3Atop(),
565 fArrayParCollPionProj[indexTune].GetProc 565 fArrayParCollPionProj[indexTune].GetProc3Ymin() ); // Target diffraction
566 SetParams( 4, fArrayParCollPionProj[indexT 566 SetParams( 4, fArrayParCollPionProj[indexTune].GetProc4A1(),
567 fArrayParCollPionProj[indexTune].G 567 fArrayParCollPionProj[indexTune].GetProc4B1(),
568 fArrayParCollPionProj[indexTune].GetProc 568 fArrayParCollPionProj[indexTune].GetProc4A2(),
569 fArrayParCollPionProj[indexTune].G 569 fArrayParCollPionProj[indexTune].GetProc4B2(),
570 fArrayParCollPionProj[indexTune].GetProc 570 fArrayParCollPionProj[indexTune].GetProc4A3(),
571 fArrayParCollPionProj[indexTune].G 571 fArrayParCollPionProj[indexTune].GetProc4Atop(),
572 fArrayParCollPionProj[indexTune].GetProc 572 fArrayParCollPionProj[indexTune].GetProc4Ymin() ); // Qexchange with Exc. Additional multiply
573 573
574 // NOTE: how can it be |ProjectileBaryonNu 574 // NOTE: how can it be |ProjectileBaryonNumber| > 10 if projectile is a pion ???
575 // 575 //
576 if ( AbsProjectileBaryonNumber > 10 || N 576 if ( AbsProjectileBaryonNumber > 10 || NumberOfTargetNucleons > 10 ) {
577 if ( ! fArrayParCollPionProj[indexTune] 577 if ( ! fArrayParCollPionProj[indexTune].IsProjDiffDissociation() )
578 SetParams( 2, 0.0 , 0.0 , 578 SetParams( 2, 0.0 , 0.0 , 0.0 , 0.0 , 0.0, 0.0 , -100.0 ); // Projectile diffraction
579 if ( ! fArrayParCollPionProj[indexTune] 579 if ( ! fArrayParCollPionProj[indexTune].IsTgtDiffDissociation() )
580 SetParams( 3, 0.0 , 0.0 , 580 SetParams( 3, 0.0 , 0.0 , 0.0 , 0.0 , 0.0, 0.0 , -100.0 ); // Target diffraction
581 } 581 }
582 582
583 /* original code --> 583 /* original code -->
584 SetDeltaProbAtQuarkExchange( 0.56 ); 584 SetDeltaProbAtQuarkExchange( 0.56 );
585 SetProjMinDiffMass( 1.0 ); // G 585 SetProjMinDiffMass( 1.0 ); // GeV
586 SetProjMinNonDiffMass( 1.0 ); // G 586 SetProjMinNonDiffMass( 1.0 ); // GeV
587 SetTarMinDiffMass( 1.16 ); // G 587 SetTarMinDiffMass( 1.16 ); // GeV
588 SetTarMinNonDiffMass( 1.16 ); // G 588 SetTarMinNonDiffMass( 1.16 ); // GeV
589 SetAveragePt2( 0.3 ); // G 589 SetAveragePt2( 0.3 ); // GeV^2
590 SetProbLogDistrPrD( 0.55 ); 590 SetProbLogDistrPrD( 0.55 );
591 SetProbLogDistr( 0.55 ); 591 SetProbLogDistr( 0.55 );
592 */ 592 */
593 593
594 // JVY update, Aug.8, 2018 --> Feb.14, 201 594 // JVY update, Aug.8, 2018 --> Feb.14, 2019
595 // 595 //
596 SetDeltaProbAtQuarkExchange( fArrayParColl 596 SetDeltaProbAtQuarkExchange( fArrayParCollPionProj[indexTune].GetDeltaProbAtQuarkExchange() );
597 SetProjMinDiffMass( fArrayParCollPionProj[ 597 SetProjMinDiffMass( fArrayParCollPionProj[indexTune].GetProjMinDiffMass() ); // GeV
598 SetProjMinNonDiffMass( fArrayParCollPionPr 598 SetProjMinNonDiffMass( fArrayParCollPionProj[indexTune].GetProjMinNonDiffMass() ); // GeV
599 SetTarMinDiffMass( fArrayParCollPionProj[i 599 SetTarMinDiffMass( fArrayParCollPionProj[indexTune].GetTgtMinDiffMass() ); // GeV
600 SetTarMinNonDiffMass( fArrayParCollPionPro 600 SetTarMinNonDiffMass( fArrayParCollPionProj[indexTune].GetTgtMinNonDiffMass() ); // GeV
601 SetAveragePt2( fArrayParCollPionProj[index 601 SetAveragePt2( fArrayParCollPionProj[indexTune].GetAveragePt2() ); // GeV^2
602 SetProbLogDistrPrD( fArrayParCollPionProj[ 602 SetProbLogDistrPrD( fArrayParCollPionProj[indexTune].GetProbLogDistrPrD() );
603 SetProbLogDistr( fArrayParCollPionProj[ind 603 SetProbLogDistr( fArrayParCollPionProj[indexTune].GetProbLogDistr() );
604 604
605 // ---> end update 605 // ---> end update
606 606
607 } else if ( ProjectileabsPDGcode == 321 || 607 } else if ( ProjectileabsPDGcode == 321 || ProjectileabsPDGcode == 311 ||
608 ProjectilePDGcode == 130 || 608 ProjectilePDGcode == 130 || ProjectilePDGcode == 310 ) { // Projectile is Kaon
609 609
610 // Proc# A1 B1 A2 610 // Proc# A1 B1 A2 B2 A3 Atop Ymin
611 SetParams( 0, 60.0 , 2.5 , 0 611 SetParams( 0, 60.0 , 2.5 , 0.0 , 0.0 , 0.0, 0.0 , -100.0 ); // Qexchange without Exc.
612 SetParams( 1, 6.0 , 1.0 , -24 612 SetParams( 1, 6.0 , 1.0 , -24.33 , 2.0 , 0.0, 0.0 , 1.40 ); // Qexchange with Exc.
613 SetParams( 2, 2.76, 1.2 , -22 613 SetParams( 2, 2.76, 1.2 , -22.5 , 2.7 ,0.04, 0.0 , 1.40 ); // Projectile diffraction
614 SetParams( 3, 1.09, 0.5 , -8 614 SetParams( 3, 1.09, 0.5 , -8.88 , 2. ,0.05, 0.0 , 1.40 ); // Target diffraction
615 SetParams( 4, 1.0, 0.0 , 0 615 SetParams( 4, 1.0, 0.0 , 0.0 , 0.0 , 0.0, 0.0 , 0.93 ); // Qexchange with Exc. Additional multiply
616 if ( AbsProjectileBaryonNumber > 10 || N 616 if ( AbsProjectileBaryonNumber > 10 || NumberOfTargetNucleons > 10 ) {
617 SetParams( 2, 0.0 , 0.0 , 617 SetParams( 2, 0.0 , 0.0 , 0.0 , 0.0 , 0.0, 0.0 , -100.0 ); // Projectile diffraction
618 SetParams( 3, 0.0 , 0.0 , 618 SetParams( 3, 0.0 , 0.0 , 0.0 , 0.0 , 0.0, 0.0 , -100.0 ); // Target diffraction
619 } 619 }
620 620
621 SetDeltaProbAtQuarkExchange( 0.6 ); 621 SetDeltaProbAtQuarkExchange( 0.6 );
622 SetProjMinDiffMass( 0.7 ); // GeV 622 SetProjMinDiffMass( 0.7 ); // GeV
623 SetProjMinNonDiffMass( 0.7 ); // GeV 623 SetProjMinNonDiffMass( 0.7 ); // GeV
624 SetTarMinDiffMass( 1.16 ); // GeV 624 SetTarMinDiffMass( 1.16 ); // GeV
625 SetTarMinNonDiffMass( 1.16 ); // GeV 625 SetTarMinNonDiffMass( 1.16 ); // GeV
626 SetAveragePt2( 0.3 ); // GeV^2 626 SetAveragePt2( 0.3 ); // GeV^2
627 SetProbLogDistrPrD( 0.55 ); 627 SetProbLogDistrPrD( 0.55 );
628 SetProbLogDistr( 0.55 ); 628 SetProbLogDistr( 0.55 );
629 629
630 } else { // Projectile is not p, n, Pi0, Pi 630 } else { // Projectile is not p, n, Pi0, Pi+, Pi-, K+, K-, K0 or their anti-particles
631 631
632 if ( ProjectileabsPDGcode > 1000 ) { // T 632 if ( ProjectileabsPDGcode > 1000 ) { // The projectile is a baryon as P or N
633 // Proc# A1 B1 633 // Proc# A1 B1 A2 B2 A3 Atop Ymin
634 SetParams( 0, 13.71, 1.75, 634 SetParams( 0, 13.71, 1.75, -30.69, 3.0 , 0.0, 1.0 , 0.93 ); // Qexchange without Exc.
635 SetParams( 1, 25.0, 1.0, - 635 SetParams( 1, 25.0, 1.0, -50.34, 1.5 , 0.0, 0.0 , 1.4 ); // Qexchange with Exc.
636 if ( Xinel > 0.) { 636 if ( Xinel > 0.) {
637 SetParams( 2, 6.0/Xinel, 0.0 ,-6.0/Xin 637 SetParams( 2, 6.0/Xinel, 0.0 ,-6.0/Xinel*16.28, 3.0 , 0.0, 0.0 , 0.93); // Projectile diffraction
638 SetParams( 3, 6.0/Xinel, 0.0 ,-6.0/Xin 638 SetParams( 3, 6.0/Xinel, 0.0 ,-6.0/Xinel*16.28, 3.0 , 0.0, 0.0 , 0.93); // Target diffraction
639 SetParams( 4, 1.0, 0.0 , 639 SetParams( 4, 1.0, 0.0 , -2.01 , 0.5 , 0.0, 0.0 , 1.4 ); // Qexchange with Exc. Additional multiply
640 } else { 640 } else {
641 SetParams( 2, 0.0, 0.0 ,0.0, 0.0 , 0.0 641 SetParams( 2, 0.0, 0.0 ,0.0, 0.0 , 0.0, 0.0 , 0.0);
642 SetParams( 3, 0.0, 0.0 ,0.0, 0.0 , 0.0 642 SetParams( 3, 0.0, 0.0 ,0.0, 0.0 , 0.0, 0.0 , 0.0);
643 SetParams( 4, 0.0, 0.0 ,0.0, 0.0 , 0.0 643 SetParams( 4, 0.0, 0.0 ,0.0, 0.0 , 0.0, 0.0 , 0.0);
644 } 644 }
645 645
646 } else { // The projectile is a meson as 646 } else { // The projectile is a meson as K+-0
647 // Proc# A1 B1 647 // Proc# A1 B1 A2 B2 A3 Atop Ymin
648 SetParams( 0, 60.0 , 2.5 , 648 SetParams( 0, 60.0 , 2.5 , 0.0 , 0.0 , 0.0, 0.0 , -100.0 ); // Qexchange without Exc.
649 SetParams( 1, 6.0 , 1.0 , - 649 SetParams( 1, 6.0 , 1.0 , -24.33 , 2.0 , 0.0, 0.0 , 1.40 ); // Qexchange with Exc.
650 SetParams( 2, 2.76, 1.2 , - 650 SetParams( 2, 2.76, 1.2 , -22.5 , 2.7 ,0.04, 0.0 , 1.40 ); // Projectile diffraction
651 SetParams( 3, 1.09, 0.5 , 651 SetParams( 3, 1.09, 0.5 , -8.88 , 2. ,0.05, 0.0 , 1.40 ); // Target diffraction
652 SetParams( 4, 1.0, 0.0 , 652 SetParams( 4, 1.0, 0.0 , 0.0 , 0.0 , 0.0, 0.0 , 0.93 ); // Qexchange with Exc. Additional multiply
653 } 653 }
654 654
655 if ( AbsProjectileBaryonNumber > 10 || N 655 if ( AbsProjectileBaryonNumber > 10 || NumberOfTargetNucleons > 10 ) {
656 SetParams( 2, 0.0 , 0.0 , 656 SetParams( 2, 0.0 , 0.0 , 0.0 , 0.0 , 0.0, 0.0 , -100.0 ); // Projectile diffraction
657 SetParams( 3, 0.0 , 0.0 , 657 SetParams( 3, 0.0 , 0.0 , 0.0 , 0.0 , 0.0, 0.0 , -100.0 ); // Target diffraction
658 } 658 }
659 659
660 SetDeltaProbAtQuarkExchange( 0.0 ); 660 SetDeltaProbAtQuarkExchange( 0.0 );
661 SetProbOfSameQuarkExchange( 0.0 ); 661 SetProbOfSameQuarkExchange( 0.0 );
662 662
663 SetProjMinDiffMass( GetMinMass(particle 663 SetProjMinDiffMass( GetMinMass(particle)/GeV );
664 SetProjMinNonDiffMass( GetMinMass(particle 664 SetProjMinNonDiffMass( GetMinMass(particle)/GeV );
665 665
666 const G4ParticleDefinition* Neutron = G4Ne 666 const G4ParticleDefinition* Neutron = G4Neutron::Neutron();
667 SetTarMinDiffMass( GetMinMass(Neutron)/ 667 SetTarMinDiffMass( GetMinMass(Neutron)/GeV );
668 SetTarMinNonDiffMass( GetMinMass(Neutron)/ 668 SetTarMinNonDiffMass( GetMinMass(Neutron)/GeV );
669 669
670 SetAveragePt2( 0.3 ); // GeV^2 670 SetAveragePt2( 0.3 ); // GeV^2
671 SetProbLogDistrPrD( 0.55 ); 671 SetProbLogDistrPrD( 0.55 );
672 SetProbLogDistr( 0.55 ); 672 SetProbLogDistr( 0.55 );
673 673
674 } 674 }
675 675
676 #ifdef debugFTFparams 676 #ifdef debugFTFparams
677 G4cout<<"DeltaProbAtQuarkExchange "<< GetDel 677 G4cout<<"DeltaProbAtQuarkExchange "<< GetDeltaProbAtQuarkExchange() << G4endl;
678 G4cout<<"ProbOfSameQuarkExchange "<< GetPro 678 G4cout<<"ProbOfSameQuarkExchange "<< GetProbOfSameQuarkExchange() << G4endl;
679 G4cout<<"ProjMinDiffMass "<< GetPro 679 G4cout<<"ProjMinDiffMass "<< GetProjMinDiffMass()/GeV <<" GeV"<< G4endl;
680 G4cout<<"ProjMinNonDiffMass "<< GetPro 680 G4cout<<"ProjMinNonDiffMass "<< GetProjMinNonDiffMass() <<" GeV"<< G4endl;
681 G4cout<<"TarMinDiffMass "<< GetTar 681 G4cout<<"TarMinDiffMass "<< GetTarMinDiffMass() <<" GeV"<< G4endl;
682 G4cout<<"TarMinNonDiffMass "<< GetTar 682 G4cout<<"TarMinNonDiffMass "<< GetTarMinNonDiffMass() <<" GeV"<< G4endl;
683 G4cout<<"AveragePt2 "<< GetAve 683 G4cout<<"AveragePt2 "<< GetAveragePt2() <<" GeV^2"<< G4endl;
684 G4cout<<"ProbLogDistrPrD "<< GetPro 684 G4cout<<"ProbLogDistrPrD "<< GetProbLogDistrPrD() << G4endl;
685 G4cout<<"ProbLogDistrTrD "<< GetPro 685 G4cout<<"ProbLogDistrTrD "<< GetProbLogDistr() << G4endl;
686 #endif 686 #endif
687 687
688 // Set parameters of nuclear destruction 688 // Set parameters of nuclear destruction
689 689
690 if ( ProjectileabsPDGcode < 1000 ) { // Mes 690 if ( ProjectileabsPDGcode < 1000 ) { // Meson projectile
691 691
692 const G4int indexTune = G4FTFTunings::Inst 692 const G4int indexTune = G4FTFTunings::Instance()->GetIndexTune( particle, KineticEnergy );
693 693
694 SetMaxNumberOfCollisions( Plab, 2.0 ); // 694 SetMaxNumberOfCollisions( Plab, 2.0 ); // 3.0 )
695 // 695 //
696 // target destruction 696 // target destruction
697 // 697 //
698 /* original code ---> 698 /* original code --->
699 SetCofNuclearDestruction( 0.00481*G4double 699 SetCofNuclearDestruction( 0.00481*G4double(NumberOfTargetNucleons)*
700 G4Exp( 4.0*(Ylab - 2.1) )/( 1.0 + 700 G4Exp( 4.0*(Ylab - 2.1) )/( 1.0 + G4Exp( 4.0*(Ylab - 2.1) ) ) );
701 701
702 SetR2ofNuclearDestruction( 1.5*fermi*fermi 702 SetR2ofNuclearDestruction( 1.5*fermi*fermi );
703 SetDofNuclearDestruction( 0.3 ); 703 SetDofNuclearDestruction( 0.3 );
704 SetPt2ofNuclearDestruction( ( 0.035 + 0.04 704 SetPt2ofNuclearDestruction( ( 0.035 + 0.04*G4Exp( 4.0*(Ylab - 2.5) )/
705 ( 1.0 705 ( 1.0 + G4Exp( 4.0*(Ylab - 2.5) ) ) )*GeV*GeV );
706 SetMaxPt2ofNuclearDestruction( 1.0*GeV*GeV 706 SetMaxPt2ofNuclearDestruction( 1.0*GeV*GeV );
707 SetExcitationEnergyPerWoundedNucleon( 40.0 707 SetExcitationEnergyPerWoundedNucleon( 40.0*MeV );
708 */ 708 */
709 double coeff = fArrayParCollMesonProj[inde 709 double coeff = fArrayParCollMesonProj[indexTune].GetNuclearTgtDestructP1();
710 // 710 //
711 // NOTE (JVY): Set this switch to false/tr 711 // NOTE (JVY): Set this switch to false/true on line 138
712 // 712 //
713 if ( fArrayParCollMesonProj[indexTune].IsN 713 if ( fArrayParCollMesonProj[indexTune].IsNuclearTgtDestructP1_ADEP() )
714 { 714 {
715 coeff *= G4double(NumberOfTargetNucleons 715 coeff *= G4double(NumberOfTargetNucleons);
716 } 716 }
717 double exfactor = G4Exp( fArrayParCollMeso 717 double exfactor = G4Exp( fArrayParCollMesonProj[indexTune].GetNuclearTgtDestructP2()
718 * (Ylab-fArrayParCo 718 * (Ylab-fArrayParCollMesonProj[indexTune].GetNuclearTgtDestructP3()) );
719 coeff *= exfactor; 719 coeff *= exfactor;
720 coeff /= ( 1.+ exfactor ); 720 coeff /= ( 1.+ exfactor );
721 721
722 SetCofNuclearDestruction( coeff ); 722 SetCofNuclearDestruction( coeff );
723 723
724 SetR2ofNuclearDestruction( fArrayParCollMe 724 SetR2ofNuclearDestruction( fArrayParCollMesonProj[indexTune].GetR2ofNuclearDestruct() );
725 SetDofNuclearDestruction( fArrayParCollMes 725 SetDofNuclearDestruction( fArrayParCollMesonProj[indexTune].GetDofNuclearDestruct() );
726 coeff = fArrayParCollMesonProj[indexTune]. 726 coeff = fArrayParCollMesonProj[indexTune].GetPt2NuclearDestructP2();
727 exfactor = G4Exp( fArrayParCollMesonProj[i 727 exfactor = G4Exp( fArrayParCollMesonProj[indexTune].GetPt2NuclearDestructP3()
728 * (Ylab-fArrayParCollMeson 728 * (Ylab-fArrayParCollMesonProj[indexTune].GetPt2NuclearDestructP4()) );
729 coeff *= exfactor; 729 coeff *= exfactor;
730 coeff /= ( 1. + exfactor ); 730 coeff /= ( 1. + exfactor );
731 SetPt2ofNuclearDestruction( (fArrayParColl 731 SetPt2ofNuclearDestruction( (fArrayParCollMesonProj[indexTune].GetPt2NuclearDestructP1()+coeff)*CLHEP::GeV*CLHEP::GeV );
732 732
733 SetMaxPt2ofNuclearDestruction( fArrayParCo 733 SetMaxPt2ofNuclearDestruction( fArrayParCollMesonProj[indexTune].GetMaxPt2ofNuclearDestruct() );
734 SetExcitationEnergyPerWoundedNucleon( fArr 734 SetExcitationEnergyPerWoundedNucleon( fArrayParCollMesonProj[indexTune].GetExciEnergyPerWoundedNucleon() );
735 735
736 } else if ( ProjectilePDGcode == -2212 || 736 } else if ( ProjectilePDGcode == -2212 || ProjectilePDGcode == -2112 ) { // for anti-baryon projectile
737 737
738 SetMaxNumberOfCollisions( Plab, 2.0 ); 738 SetMaxNumberOfCollisions( Plab, 2.0 );
739 739
740 SetCofNuclearDestruction( 0.00481*G4double 740 SetCofNuclearDestruction( 0.00481*G4double(NumberOfTargetNucleons)*
741 G4Exp( 4.0*(Ylab - 2.1) )/( 1.0 + G 741 G4Exp( 4.0*(Ylab - 2.1) )/( 1.0 + G4Exp( 4.0*(Ylab - 2.1) ) ) );
742 SetR2ofNuclearDestruction( 1.5*fermi*fermi 742 SetR2ofNuclearDestruction( 1.5*fermi*fermi );
743 SetDofNuclearDestruction( 0.3 ); 743 SetDofNuclearDestruction( 0.3 );
744 SetPt2ofNuclearDestruction( ( 0.035 + 0.04 744 SetPt2ofNuclearDestruction( ( 0.035 + 0.04*G4Exp( 4.0*(Ylab - 2.5) )/
745 ( 1.0 745 ( 1.0 + G4Exp( 4.0*(Ylab - 2.5) ) ) )*GeV*GeV );
746 SetMaxPt2ofNuclearDestruction( 1.0*GeV*GeV 746 SetMaxPt2ofNuclearDestruction( 1.0*GeV*GeV );
747 SetExcitationEnergyPerWoundedNucleon( 40.0 747 SetExcitationEnergyPerWoundedNucleon( 40.0*MeV );
748 if ( Plab < 2.0 ) { // 2 GeV/c 748 if ( Plab < 2.0 ) { // 2 GeV/c
749 // For slow anti-baryon we have to garan 749 // For slow anti-baryon we have to garanty putting on mass-shell
750 SetCofNuclearDestruction( 0.0 ); 750 SetCofNuclearDestruction( 0.0 );
751 SetR2ofNuclearDestruction( 1.5*fermi*fer 751 SetR2ofNuclearDestruction( 1.5*fermi*fermi ); // this is equivalent to setting a few line above
752 752 // is it even necessary ?
753 SetDofNuclearDestruction( 0.01 ); 753 SetDofNuclearDestruction( 0.01 );
754 SetPt2ofNuclearDestruction( 0.035*GeV*Ge 754 SetPt2ofNuclearDestruction( 0.035*GeV*GeV );
755 SetMaxPt2ofNuclearDestruction( 0.04*GeV* 755 SetMaxPt2ofNuclearDestruction( 0.04*GeV*GeV );
756 } 756 }
757 757
758 } else { // Projectile baryon assumed 758 } else { // Projectile baryon assumed
759 759
760 // Below, in the call to the G4FTFTunings: 760 // Below, in the call to the G4FTFTunings::GetIndexTune method, we pass the proton
761 // as projectile, instead of the real one, 761 // as projectile, instead of the real one, because for the treatment of nuclear
762 // destruction, we assume for this categor 762 // destruction, we assume for this category of hadron projectiles the same treatment
763 // as for "baryon". 763 // as for "baryon".
764 const G4int indexTune = G4FTFTunings::Inst 764 const G4int indexTune = G4FTFTunings::Instance()->GetIndexTune( G4Proton::Definition(), KineticEnergy );
765 765
766 // NOTE (JVY) FIXME !!! Will decide later 766 // NOTE (JVY) FIXME !!! Will decide later how/if to make this one configurable...
767 // 767 //
768 SetMaxNumberOfCollisions( Plab, 2.0 ); 768 SetMaxNumberOfCollisions( Plab, 2.0 );
769 769
770 // projectile destruction - does NOT reall 770 // projectile destruction - does NOT really matter for particle projectile, only for a nucleus projectile
771 // 771 //
772 double coeff = 0.; 772 double coeff = 0.;
773 coeff = fArrayParCollBaryonProj[indexTune] 773 coeff = fArrayParCollBaryonProj[indexTune].GetNuclearProjDestructP1();
774 // 774 //
775 // NOTE (JVY): Set this switch to false/tr 775 // NOTE (JVY): Set this switch to false/true on line 136
776 // 776 //
777 if ( fArrayParCollBaryonProj[indexTune].Is 777 if ( fArrayParCollBaryonProj[indexTune].IsNuclearProjDestructP1_NBRNDEP() )
778 { 778 {
779 coeff *= G4double(AbsProjectileBaryonNum 779 coeff *= G4double(AbsProjectileBaryonNumber);
780 } 780 }
781 double exfactor = G4Exp( fArrayParCollBary 781 double exfactor = G4Exp( fArrayParCollBaryonProj[indexTune].GetNuclearProjDestructP2()*
782 (Ylab-fArrayParCollBaryonProj[index 782 (Ylab-fArrayParCollBaryonProj[indexTune].GetNuclearProjDestructP3()) );
783 coeff *= exfactor; 783 coeff *= exfactor;
784 coeff /= ( 1.+ exfactor ); 784 coeff /= ( 1.+ exfactor );
785 SetCofNuclearDestructionPr( coeff ); 785 SetCofNuclearDestructionPr( coeff );
786 786
787 // target desctruction 787 // target desctruction
788 // 788 //
789 coeff = fArrayParCollBaryonProj[indexTune] 789 coeff = fArrayParCollBaryonProj[indexTune].GetNuclearTgtDestructP1();
790 // 790 //
791 // NOTE (JVY): Set this switch to false/tr 791 // NOTE (JVY): Set this switch to false/true on line 138
792 // 792 //
793 if ( fArrayParCollBaryonProj[indexTune].Is 793 if ( fArrayParCollBaryonProj[indexTune].IsNuclearTgtDestructP1_ADEP() )
794 { 794 {
795 coeff *= G4double(NumberOfTargetNucleons 795 coeff *= G4double(NumberOfTargetNucleons);
796 } 796 }
797 exfactor = G4Exp( fArrayParCollBaryonProj[ 797 exfactor = G4Exp( fArrayParCollBaryonProj[indexTune].GetNuclearTgtDestructP2()*
798 (Ylab-fArrayParCollBaryonProj[indexT 798 (Ylab-fArrayParCollBaryonProj[indexTune].GetNuclearTgtDestructP3()) );
799 coeff *= exfactor; 799 coeff *= exfactor;
800 coeff /= ( 1.+ exfactor ); 800 coeff /= ( 1.+ exfactor );
801 SetCofNuclearDestruction( coeff ); 801 SetCofNuclearDestruction( coeff );
802 802
803 SetR2ofNuclearDestruction( fArrayParCollBa 803 SetR2ofNuclearDestruction( fArrayParCollBaryonProj[indexTune].GetR2ofNuclearDestruct() );
804 SetDofNuclearDestruction( fArrayParCollBar 804 SetDofNuclearDestruction( fArrayParCollBaryonProj[indexTune].GetDofNuclearDestruct() );
805 805
806 coeff = fArrayParCollBaryonProj[indexTune] 806 coeff = fArrayParCollBaryonProj[indexTune].GetPt2NuclearDestructP2();
807 exfactor = G4Exp( fArrayParCollBaryonProj[ 807 exfactor = G4Exp( fArrayParCollBaryonProj[indexTune].GetPt2NuclearDestructP3()*
808 (Ylab-fArrayParCollBaryonProj[indexT 808 (Ylab-fArrayParCollBaryonProj[indexTune].GetPt2NuclearDestructP4()) );
809 coeff *= exfactor; 809 coeff *= exfactor;
810 coeff /= ( 1. + exfactor ); 810 coeff /= ( 1. + exfactor );
811 SetPt2ofNuclearDestruction( (fArrayParColl 811 SetPt2ofNuclearDestruction( (fArrayParCollBaryonProj[indexTune].GetPt2NuclearDestructP1()+coeff)*CLHEP::GeV*CLHEP::GeV );
812 812
813 SetMaxPt2ofNuclearDestruction( fArrayParCo 813 SetMaxPt2ofNuclearDestruction( fArrayParCollBaryonProj[indexTune].GetMaxPt2ofNuclearDestruct() );
814 SetExcitationEnergyPerWoundedNucleon( fArr 814 SetExcitationEnergyPerWoundedNucleon( fArrayParCollBaryonProj[indexTune].GetExciEnergyPerWoundedNucleon() );
815 815
816 } 816 }
817 817
818 #ifdef debugFTFparams 818 #ifdef debugFTFparams
819 G4cout<<"CofNuclearDestructionPr " 819 G4cout<<"CofNuclearDestructionPr "<< GetCofNuclearDestructionPr() << G4endl;
820 G4cout<<"CofNuclearDestructionTr " 820 G4cout<<"CofNuclearDestructionTr "<< GetCofNuclearDestruction() << G4endl;
821 G4cout<<"R2ofNuclearDestruction " 821 G4cout<<"R2ofNuclearDestruction "<< GetR2ofNuclearDestruction()/fermi/fermi <<" fermi^2"<< G4endl;
822 G4cout<<"DofNuclearDestruction " 822 G4cout<<"DofNuclearDestruction "<< GetDofNuclearDestruction() << G4endl;
823 G4cout<<"Pt2ofNuclearDestruction " 823 G4cout<<"Pt2ofNuclearDestruction "<< GetPt2ofNuclearDestruction()/GeV/GeV <<" GeV^2"<< G4endl;
824 G4cout<<"ExcitationEnergyPerWoundedNucleon " 824 G4cout<<"ExcitationEnergyPerWoundedNucleon "<< GetExcitationEnergyPerWoundedNucleon() <<" MeV"<< G4endl;
825 #endif 825 #endif
826 826
827 //SetCofNuclearDestruction( 0.47*G4Exp( 2.0* 827 //SetCofNuclearDestruction( 0.47*G4Exp( 2.0*(Ylab - 2.5) )/( 1.0 + G4Exp( 2.0*(Ylab - 2.5) ) ) );
828 //SetPt2ofNuclearDestruction( ( 0.035 + 0.1* 828 //SetPt2ofNuclearDestruction( ( 0.035 + 0.1*G4Exp( 4.0*(Ylab - 3.0) )/( 1.0 + G4Exp( 4.0*(Ylab - 3.0) ) ) )*GeV*GeV );
829 829
830 //SetMagQuarkExchange( 120.0 ); // 210 830 //SetMagQuarkExchange( 120.0 ); // 210.0 PipP
831 //SetSlopeQuarkExchange( 2.0 ); 831 //SetSlopeQuarkExchange( 2.0 );
832 //SetDeltaProbAtQuarkExchange( 0.6 ); 832 //SetDeltaProbAtQuarkExchange( 0.6 );
833 //SetProjMinDiffMass( 0.7 ); // GeV 833 //SetProjMinDiffMass( 0.7 ); // GeV 1.1
834 //SetProjMinNonDiffMass( 0.7 ); // GeV 834 //SetProjMinNonDiffMass( 0.7 ); // GeV
835 //SetProbabilityOfProjDiff( 0.0); // 0.8 835 //SetProbabilityOfProjDiff( 0.0); // 0.85*G4Pow::GetInstance()->powA( s/GeV/GeV, -0.5 ) ); // 40/32 X-dif/X-inel
836 //SetTarMinDiffMass( 1.1 ); // GeV 836 //SetTarMinDiffMass( 1.1 ); // GeV
837 //SetTarMinNonDiffMass( 1.1 ); // GeV 837 //SetTarMinNonDiffMass( 1.1 ); // GeV
838 //SetProbabilityOfTarDiff( 0.0 ); // 0.8 838 //SetProbabilityOfTarDiff( 0.0 ); // 0.85*G4Pow::GetInstance()->powA( s/GeV/GeV, -0.5 ) ); // 40/32 X-dif/X-inel
839 839
840 //SetAveragePt2( 0.0 ); // GeV 840 //SetAveragePt2( 0.0 ); // GeV^2 0.3
841 //------------------------------------ 841 //------------------------------------
842 //SetProbabilityOfElasticScatt( 1.0, 1.0); 842 //SetProbabilityOfElasticScatt( 1.0, 1.0); //(Xtotal, Xelastic);
843 //SetProbabilityOfProjDiff( 1.0*0.62*G4Pow:: 843 //SetProbabilityOfProjDiff( 1.0*0.62*G4Pow::GetInstance()->powA( s/GeV/GeV, -0.51 ) ); // 0->1
844 //SetProbabilityOfTarDiff( 4.0*0.62*G4Pow::G 844 //SetProbabilityOfTarDiff( 4.0*0.62*G4Pow::GetInstance()->powA( s/GeV/GeV, -0.51 ) ); // 2->4
845 //SetAveragePt2( 0.3 ); 845 //SetAveragePt2( 0.3 ); // (0.15)
846 //SetAvaragePt2ofElasticScattering( 0.0 ); 846 //SetAvaragePt2ofElasticScattering( 0.0 );
847 847
848 //SetMaxNumberOfCollisions( Plab, 6.0 ); //( 848 //SetMaxNumberOfCollisions( Plab, 6.0 ); //(4.0*(Plab + 0.01), Plab); // 6.0 );
849 //SetAveragePt2( 0.15 ); 849 //SetAveragePt2( 0.15 );
850 //SetCofNuclearDestruction(-1.);//( 0.75 ); 850 //SetCofNuclearDestruction(-1.);//( 0.75 ); // (0.25)
851 //SetExcitationEnergyPerWoundedNucleon(0.);/ 851 //SetExcitationEnergyPerWoundedNucleon(0.);//( 30.0*MeV ); // (75.0*MeV)
852 //SetDofNuclearDestruction(0.);//( 0.2 ); // 852 //SetDofNuclearDestruction(0.);//( 0.2 ); //0.4 // 0.3 0.5
853 853
854 /* 854 /*
855 SetAveragePt2(0.3); 855 SetAveragePt2(0.3);
856 SetCofNuclearDestructionPr(0.); 856 SetCofNuclearDestructionPr(0.);
857 SetCofNuclearDestruction(0.); // 857 SetCofNuclearDestruction(0.); //( 0.5 ); (0.25)
858 SetExcitationEnergyPerWoundedNucleon(0.); // 858 SetExcitationEnergyPerWoundedNucleon(0.); // 30.0*MeV; (75.0*MeV)
859 SetDofNuclearDestruction(0.); // 859 SetDofNuclearDestruction(0.); // 0.2; 0.4; 0.3; 0.5
860 SetPt2ofNuclearDestruction(0.); // 860 SetPt2ofNuclearDestruction(0.); //(2.*0.075*GeV*GeV); ( 0.3*GeV*GeV ); (0.168*GeV*GeV)
861 */ 861 */
862 862
863 //SetExcitationEnergyPerWoundedNucleon(0.001 863 //SetExcitationEnergyPerWoundedNucleon(0.001);
864 //SetPt2Kink( 0.0*GeV*GeV ); 864 //SetPt2Kink( 0.0*GeV*GeV );
865 865
866 //SetRadiusOfHNinteractions2( Xtotal/pi/10.0 866 //SetRadiusOfHNinteractions2( Xtotal/pi/10.0 /2.);
867 //SetRadiusOfHNinteractions2( (Xtotal - Xela 867 //SetRadiusOfHNinteractions2( (Xtotal - Xelastic)/pi/10.0 );
868 //SetProbabilityOfElasticScatt( 1.0, 0.0); 868 //SetProbabilityOfElasticScatt( 1.0, 0.0);
869 /* 869 /*
870 G4cout << "Pt2 " << GetAveragePt2()<<" "<<Ge 870 G4cout << "Pt2 " << GetAveragePt2()<<" "<<GetAveragePt2()/GeV/GeV<<G4endl;
871 G4cout << "Cnd " << GetCofNuclearDestruction 871 G4cout << "Cnd " << GetCofNuclearDestruction() << G4endl;
872 G4cout << "Dnd " << GetDofNuclearDestruction 872 G4cout << "Dnd " << GetDofNuclearDestruction() << G4endl;
873 G4cout << "Pt2 " << GetPt2ofNuclearDestructi 873 G4cout << "Pt2 " << GetPt2ofNuclearDestruction()/GeV/GeV << G4endl;
874 */ 874 */
875 875
876 } 876 }
877 877
878 //============================================ 878 //============================================================================
879 879
880 G4double G4FTFParameters::GetMinMass( const G4 880 G4double G4FTFParameters::GetMinMass( const G4ParticleDefinition* aParticle ) {
881 // The code is used for estimating the mini 881 // The code is used for estimating the minimal string mass produced in diffraction dissociation.
882 // The indices used for minMassQDiQStr must 882 // The indices used for minMassQDiQStr must be between 1 and 5, corresponding to the 5 considered
883 // quarks: d, u, s, c and b; enforcing this 883 // quarks: d, u, s, c and b; enforcing this explicitly avoids compilation errors.
884 G4double EstimatedMass = 0.0; 884 G4double EstimatedMass = 0.0;
885 G4int partID = std::abs(aParticle->GetPDGEn 885 G4int partID = std::abs(aParticle->GetPDGEncoding());
886 G4int Qleft = std::max( partID/100, 1 886 G4int Qleft = std::max( partID/100, 1 );
887 G4int Qright = std::max( (partID/ 10)%10, 1 887 G4int Qright = std::max( (partID/ 10)%10, 1 );
888 if ( Qleft < 6 && Qright < 6 ) { 888 if ( Qleft < 6 && Qright < 6 ) { // Q-Qbar string
889 EstimatedMass = StringMass->minMassQQbarS 889 EstimatedMass = StringMass->minMassQQbarStr[Qleft-1][Qright-1];
890 } else if ( Qleft < 6 && Qright > 6 ) { 890 } else if ( Qleft < 6 && Qright > 6 ) { // Q - DiQ string
891 G4int q1 = std::max( std::min( Qright/10, 891 G4int q1 = std::max( std::min( Qright/10, 5 ), 1 );
892 G4int q2 = std::max( std::min( Qright%10, 892 G4int q2 = std::max( std::min( Qright%10, 5 ), 1 );
893 EstimatedMass = StringMass->minMassQDiQSt 893 EstimatedMass = StringMass->minMassQDiQStr[Qleft-1][q1-1][q2-1];
894 } else if ( Qleft > 6 && Qright < 6 ) { 894 } else if ( Qleft > 6 && Qright < 6 ) { // DiQ - Q string
895 G4int q1 = std::max( std::min( Qleft/10, 895 G4int q1 = std::max( std::min( Qleft/10, 5 ), 1 );
896 G4int q2 = std::max( std::min( Qleft%10, 896 G4int q2 = std::max( std::min( Qleft%10, 5 ), 1 );
897 EstimatedMass = StringMass->minMassQDiQSt 897 EstimatedMass = StringMass->minMassQDiQStr[Qright-1][q1-1][q2-1];
898 } 898 }
899 return EstimatedMass; 899 return EstimatedMass;
900 } 900 }
901 901
902 //============================================ 902 //============================================================================
903 903
904 G4double G4FTFParameters::GetProcProb( const G 904 G4double G4FTFParameters::GetProcProb( const G4int ProcN, const G4double y ) {
905 G4double Prob( 0.0 ); 905 G4double Prob( 0.0 );
906 if ( y < ProcParams[ProcN][6] ) { 906 if ( y < ProcParams[ProcN][6] ) {
907 Prob = ProcParams[ProcN][5]; 907 Prob = ProcParams[ProcN][5];
908 if (Prob < 0.) Prob=0.; 908 if (Prob < 0.) Prob=0.;
909 return Prob; 909 return Prob;
910 } 910 }
911 Prob = ProcParams[ProcN][0] * G4Exp( -ProcPa 911 Prob = ProcParams[ProcN][0] * G4Exp( -ProcParams[ProcN][1]*y ) +
912 ProcParams[ProcN][2] * G4Exp( -ProcPa 912 ProcParams[ProcN][2] * G4Exp( -ProcParams[ProcN][3]*y ) +
913 ProcParams[ProcN][4]; 913 ProcParams[ProcN][4];
914 if (Prob < 0.) Prob=0.; 914 if (Prob < 0.) Prob=0.;
915 return Prob; 915 return Prob;
916 } 916 }
917 917
918 //============================================ 918 //============================================================================
919 919
920 G4FTFParameters::~G4FTFParameters() { 920 G4FTFParameters::~G4FTFParameters() {
921 if ( StringMass ) delete StringMass; 921 if ( StringMass ) delete StringMass;
922 } 922 }
923 923
924 //============================================ 924 //============================================================================
925 925
926 void G4FTFParameters::Reset() 926 void G4FTFParameters::Reset()
927 { 927 {
928 FTFhNcmsEnergy = 0.0; 928 FTFhNcmsEnergy = 0.0;
929 FTFXtotal = 0.0; 929 FTFXtotal = 0.0;
930 FTFXelastic = 0.0; 930 FTFXelastic = 0.0;
931 FTFXinelastic = 0.0; 931 FTFXinelastic = 0.0;
932 FTFXannihilation = 0.0; 932 FTFXannihilation = 0.0;
933 ProbabilityOfAnnihilation = 0.0; 933 ProbabilityOfAnnihilation = 0.0;
934 ProbabilityOfElasticScatt = 0.0; 934 ProbabilityOfElasticScatt = 0.0;
935 RadiusOfHNinteractions2 = 0.0; 935 RadiusOfHNinteractions2 = 0.0;
936 FTFSlope = 0.0; 936 FTFSlope = 0.0;
937 AvaragePt2ofElasticScattering = 0.0; 937 AvaragePt2ofElasticScattering = 0.0;
938 FTFGamma0 = 0.0; 938 FTFGamma0 = 0.0;
939 DeltaProbAtQuarkExchange = 0.0; 939 DeltaProbAtQuarkExchange = 0.0;
940 ProbOfSameQuarkExchange = 0.0; 940 ProbOfSameQuarkExchange = 0.0;
941 ProjMinDiffMass = 0.0; 941 ProjMinDiffMass = 0.0;
942 ProjMinNonDiffMass = 0.0; 942 ProjMinNonDiffMass = 0.0;
943 ProbLogDistrPrD = 0.0; 943 ProbLogDistrPrD = 0.0;
944 TarMinDiffMass = 0.0; 944 TarMinDiffMass = 0.0;
945 TarMinNonDiffMass = 0.0; 945 TarMinNonDiffMass = 0.0;
946 AveragePt2 = 0.0; 946 AveragePt2 = 0.0;
947 ProbLogDistr = 0.0; 947 ProbLogDistr = 0.0;
948 Pt2kink = 0.0; 948 Pt2kink = 0.0;
949 MaxNumberOfCollisions = 0.0; 949 MaxNumberOfCollisions = 0.0;
950 ProbOfInelInteraction = 0.0; 950 ProbOfInelInteraction = 0.0;
951 CofNuclearDestructionPr = 0.0; 951 CofNuclearDestructionPr = 0.0;
952 CofNuclearDestruction = 0.0; 952 CofNuclearDestruction = 0.0;
953 R2ofNuclearDestruction = 0.0; 953 R2ofNuclearDestruction = 0.0;
954 ExcitationEnergyPerWoundedNucleon = 0.0; 954 ExcitationEnergyPerWoundedNucleon = 0.0;
955 DofNuclearDestruction = 0.0; 955 DofNuclearDestruction = 0.0;
956 Pt2ofNuclearDestruction = 0.0; 956 Pt2ofNuclearDestruction = 0.0;
957 MaxPt2ofNuclearDestruction = 0.0; 957 MaxPt2ofNuclearDestruction = 0.0;
958 958
959 for ( G4int i = 0; i < 4; i++ ) { 959 for ( G4int i = 0; i < 4; i++ ) {
960 for ( G4int j = 0; j < 7; j++ ) { 960 for ( G4int j = 0; j < 7; j++ ) {
961 ProcParams[i][j] = 0.0; 961 ProcParams[i][j] = 0.0;
962 } 962 }
963 } 963 }
964 964
965 return; 965 return;
966 } 966 }
967 967
968 968