Few minor parameter fixes

Implements JETSCAPE PR 177 in X-SCAPE.

config/jetscape_main.xml

@@ -145,12 +145,14 @@
     <tStart> 0.6 </tStart> <!-- Start time of jet quenching, proper time, fm/c   -->
     <mutex>ON</mutex>
     <AddLiquefier> false </AddLiquefier>
+    <lambdaQCD>0.2</lambdaQCD> <!-- 0.4 is the value chosen in JETSET -->
 
     <Matter>
       <name>Matter</name>
       <matter_on> 1 </matter_on>
       <Q0> 2.0 </Q0>
       <T0> 0.16 </T0>
+      <!--if negative vir_factor: parton 3-momentum of hard parton is used to compute initial virtuality in shower; if positive: transverse momentum of parton is used-->
       <vir_factor> 0.25 </vir_factor>
       <in_vac> 1 </in_vac>
       <recoil_on> 0 </recoil_on>
@@ -351,20 +353,30 @@
    <JetHadronization>
        <name>colored/colorless/hybrid</name>
        <take_recoil>1</take_recoil>
-       <!--eCMforHadronization only for pp collisions-->
-       <!--in hybrid put the full eCM here, this is distributed to the beam partons-->
+       <!--use eCMforHadronization only for pp collisions to introduce beam partons-->
+       <!--recommended value for e+e- and AA = 0-->
+       <!--in all hadronization modules put the full eCM here, this is distributed to the beam partons (1/6 for each one)-->
        <eCMforHadronization>5020</eCMforHadronization>
+       <!--maximum energy level for mesons-->
        <reco_Mlevelmax>1</reco_Mlevelmax>
+        <!--maximum energy level for baryons-->
        <reco_Blevelmax>1</reco_Blevelmax>
+       <!--switch for recombination of Goldstone bosons-->
        <reco_goldstone>0</reco_goldstone>
+       <!--squared maximum distance for recombination with thermal partons-->
        <thermreco_distmax>5.0</thermreco_distmax>
+       <!--factor for recombination in the shower-->
        <shower_recofactor>1.0</shower_recofactor>
+       <!--factor for recombination with thermal partons-->
        <thermal_recofactor>1.0</thermal_recofactor>
+       <!--switch for recombination biased by the time of the partons-->
        <recobias_t>0</recobias_t>
        <hydro_Tc>0.16</hydro_Tc>
        <eta_max_boost_inv>2.0</eta_max_boost_inv> <!--only relevant for boost invariant hydro, max eta for thermal parton sampler-->
+       <!--free-straming propagation of hadrons after hadronization-->
        <had_postprop>0.0</had_postprop>
-       <part_prop>1.0</part_prop>
+       <!--free-streaming of partons before hadronization-->
+       <part_prop>0.0</part_prop>
        <pythia_decays>on</pythia_decays> <!-- lets the particles given to pythia decay-->
        <tau0Max>10.0</tau0Max> <!-- only particles with tau0 < tau0Max (given in mm/c) can decay, increase to include weak decays-->
        <weak_decays>depracted</weak_decays> <!-- use the parameters pythia_decays and tau0Max, this parameter is only a dummy to make old xml files work-->

src/framework/JetScapeConstants.h

@@ -31,9 +31,6 @@ static double Ca = 3.0;
 
 static double Nc = 3.0;
 
-static double Lambda_QCD = 0.2;
-// 0.4 is the value chosen in JETSET
-
 static const double hbarC = 0.197327053;
 
 static const double fmToGeVinv = 1.0 / hbarC;

src/hadronization/ColoredHadronization.cc

@@ -41,7 +41,7 @@ void ColoredHadronization::InitTask() {
 
   double p_read_xml =
       GetXMLElementDouble({"JetHadronization", "eCMforHadronization"});
-  p_fake = p_read_xml;
+  p_fake = p_read_xml / 6.;
 
   /*std::string weak_decays =
       GetXMLElementText({"JetHadronization", "weak_decays"});*/

src/hadronization/ColorlessHadronization.cc

@@ -54,7 +54,7 @@ void ColorlessHadronization::InitTask() {
   // Read sqrts to know remnants energies
   double p_read_xml =
       GetXMLElementDouble({"JetHadronization", "eCMforHadronization"});
-  p_fake = p_read_xml;
+  p_fake = p_read_xml / 6.;
 
   take_recoil = GetXMLElementInt({"JetHadronization", "take_recoil"});
 
@@ -110,6 +110,8 @@ void ColorlessHadronization::InitTask() {
     JSINFO << "Also reading in: " << s;
     pythia.readString(s);
   }
+
+  Lambda_QCD = GetXMLElementDouble({"Eloss","lambdaQCD"});
 
   // Initialize random number distribution
   ZeroOneDistribution = std::uniform_real_distribution<double> { 0.0, 1.0 };

src/hadronization/ColorlessHadronization.h

@@ -37,6 +37,7 @@ class ColorlessHadronization
 private:
   double p_fake;
   bool take_recoil;
+  double Lambda_QCD;
 
   // Allows the registration of the module so that it is available to be used by the Jetscape framework.
   static RegisterJetScapeModule<ColorlessHadronization> reg;

src/hadronization/ThermPtnSampler.cc

@@ -1106,6 +1106,8 @@ void ThermalPartonSampler::sample_2p1d(double eta_max){
 				// Sample rest frame momentum given T and mass of light quark
 				MCSampler(TRead, 1); // NewP=P, NewX=Px, ...
 
+				Vel[1] /= cosh(eta_slice);
+				Vel[2] /= cosh(eta_slice);
 				Vel[3] = tanh(eta_slice);
 				double vsquare = Vel[1]*Vel[1] + Vel[2]*Vel[2] + Vel[3]*Vel[3];
 				if(vsquare < 10e-16){
@@ -1206,6 +1208,8 @@ void ThermalPartonSampler::sample_2p1d(double eta_max){
 				// Sample rest frame momentum given T and mass of s quark
 				MCSampler(TRead, 2); // NewP=P, NewX=Px, ...
 
+				Vel[1] /= cosh(eta_slice);
+				Vel[2] /= cosh(eta_slice);
 				Vel[3] = tanh(eta_slice);
 				double vsquare = Vel[1]*Vel[1] + Vel[2]*Vel[2] + Vel[3]*Vel[3];
 				if(vsquare < 10e-16){

src/initialstate/epemGun.cc

@@ -140,6 +140,8 @@ void epemGun::InitTask() {
     throw std::runtime_error("Pythia init() failed.");
   }
 
+  Lambda_QCD = GetXMLElementDouble({"Eloss","lambdaQCD"});
+
   // Initialize random number distribution
   ZeroOneDistribution = uniform_real_distribution<double>{0.0, 1.0};
 

src/initialstate/epemGun.h

@@ -31,6 +31,7 @@ class epemGun : public HardProcess, public Pythia8::Pythia {
   //double pTHatMax;
   double eCM;
   //bool FSR_on;
+  double Lambda_QCD;
 
   // Allows the registration of the module so that it is available to be used by the Jetscape framework.
   static RegisterJetScapeModule<epemGun> reg;

src/jet/Matter.cc

@@ -142,6 +142,7 @@ void Matter::InitTask() {
   hydro_Tc = GetXMLElementDouble({"Eloss", "Matter", "hydro_Tc"});
   brick_length = GetXMLElementDouble({"Eloss", "Matter", "brick_length"});
   vir_factor = GetXMLElementDouble({"Eloss", "Matter", "vir_factor"});
+  Lambda_QCD = GetXMLElementDouble({"Eloss","lambdaQCD"});
 
   if (vir_factor < 0.0) {
     cout << "Reminder: negative vir_factor is set, initial energy will be used "

src/jet/Matter.h

@@ -107,6 +107,7 @@ class Matter : public JetEnergyLossModule<Matter> //, public std::enable_shared_
   double initR0, initRx, initRy, initRz, initVx, initVy, initVz, initRdotV,
       initVdotV, initEner;
   double Q00, Q0, T0;
+  double Lambda_QCD;
 
   static const int dimQhatTab = 151;
   double qhatTab1D[dimQhatTab] = {0.0};