src/framework/Afterburner.cc

0001 /*******************************************************************************
0002  * Copyright (c) The JETSCAPE Collaboration, 2018
0003  *
0004  * Modular, task-based framework for simulating all aspects of heavy-ion collisions
0005  *
0006  * For the list of contributors see AUTHORS.
0007  *
0008  * Report issues at https://github.com/JETSCAPE/JETSCAPE/issues
0009  *
0010  * or via email to bugs.jetscape@gmail.com
0011  *
0012  * Distributed under the GNU General Public License 3.0 (GPLv3 or later).
0013  * See COPYING for details.
0014  ******************************************************************************/
0015 // This is a general basic class for hadronic afterburner
0016
0017 #include "./Afterburner.h"
0018 #include "./JetScapeSignalManager.h"
0019
0020 using namespace std;
0021
0022 namespace Jetscape {
0023 void Afterburner::Init() {
0024   // Makes sure that XML file with options and parameters is loaded
0025   JetScapeModuleBase::Init();
0026   JSINFO << "Initializing Afterburner : " << GetId() << " ...";
0027   // Initialize random number distribution
0028   ZeroOneDistribution = uniform_real_distribution<double>{0.0, 1.0};
0029   InitTask();
0030 }
0031
0032 void Afterburner::Exec() {
0033   VERBOSE(2) << "Afterburner running: " << GetId() << " ...";
0034   ExecuteTask();
0035 }
0036
0037 std::vector<std::vector<std::shared_ptr<Hadron>>> Afterburner::GetSoftParticlizationHadrons() {
0038   auto soft_particlization = JetScapeSignalManager::Instance()->GetSoftParticlizationPointer().lock();
0039   if (!soft_particlization) {
0040     JSWARN << "No soft particlization module found. Check if fragmentation"
0041            << " hadrons are handed to afterburner.";
0042     std::vector<std::shared_ptr<Hadron>> hadrons;
0043     dummy.push_back(hadrons);
0044     return dummy;
0045   } else {
0046     return soft_particlization->Hadron_list_;
0047   }
0048 }
0049
0050 std::vector<shared_ptr<Hadron>> Afterburner::GetFragmentationHadrons() {
0051   JSINFO << "Get fragmentation hadrons in Afterburner";
0052   auto hadronization_mgr = JetScapeSignalManager::Instance()->GetHadronizationManagerPointer().lock();
0053   if (!hadronization_mgr) {
0054     JSWARN << "No hardronization module found. It is necessary to include"
0055           << " fragmentation hadrons to afterburner as requested.";
0056     exit(1);
0057   }
0058   std::vector<shared_ptr<Hadron>> h_list;
0059   hadronization_mgr->GetHadrons(h_list);
0060   JSINFO << "Got " << h_list.size() << " fragmentation hadrons from HadronizationManager.";
0061
0062   std::vector<shared_ptr<Hadron>> h_list_new;
0063   rand_int_ptr_ = (std::make_shared<std::uniform_int_distribution<int>>(0,1));
0064   for (auto h : h_list) {
0065     if (h->has_no_position()) {
0066       JSDEBUG << "Found fragmentation hadron without properly set position in "
0067                 "Afterburner.\nInclusion of fragmentation hadrons only "
0068                 "possible for HybridHadronization.";
0069     }
0070
0071     //move all the fragmentation hadrons a little bit around to avoid having
0072     //multiple hadrons at the same position if they are at the same position
0073     const FourVector r = h->x_in();
0074     const double rand_x = ZeroOneDistribution(*GetMt19937Generator()) * 2e-4 - 1e-4;
0075     const double rand_y = ZeroOneDistribution(*GetMt19937Generator()) * 2e-4 - 1e-4;
0076     const double rand_z = ZeroOneDistribution(*GetMt19937Generator()) * 2e-4 - 1e-4;
0077     double position_smeared[4] = {r.t(), r.x()+rand_x, r.y()+rand_y, r.z()+rand_z};
0078     h->set_x(position_smeared);
0079
0080     if ((std::abs(h->pid())>10) && (h->pid() != 21)) {
0081       if (h->pstat() > 0) {
0082         // convert Kaon-L or Kaon-S into K0 or Anti-K0
0083         if (h->pid() == 310 || h->pid() == 130) {
0084           const int rand_int = (*rand_int_ptr_)(*GetMt19937Generator());
0085           const int id = (rand_int == 0) ? 311 : -311;
0086           h->set_id(id);
0087         }
0088         h_list_new.push_back(h);
0089       } else if(h->pstat() < 0) {
0090         // convert Kaon-L or Kaon-S into K0 or Anti-K0
0091         // change id of negative Kaons to make them consistent with the SMASH output
0092         if (h->pid() == 310 || h->pid() == 130) {
0093           const int rand_int = (*rand_int_ptr_)(*GetMt19937Generator());
0094           const int id = (rand_int == 0) ? 311 : -311;
0095           h->set_id(id);
0096         }
0097       }
0098     } else if((std::abs(h->pid())<10) || (h->pid() == 21)){
0099       JSWARN << "Found a free quark or gluon! This can not be handed over to SMASH.\n"
0100                 "Check confinement in hadronization module!";
0101     }
0102   }
0103   return h_list_new;
0104 }
0105
0106 std::vector<std::vector<std::shared_ptr<Hadron>>> Afterburner::GatherAfterburnerHadrons() {
0107   std::vector<std::vector<shared_ptr<Hadron>>> afterburner_had_events;
0108   afterburner_had_events = GetSoftParticlizationHadrons();
0109
0110   if (GetXMLElementInt({"Afterburner", "output_only_final_state_hadrons"})) {
0111     // clear Hadron_list_ in soft_particlization, otherwise the final hadron
0112     // output of the writer contains also the soft hadrons which were used as
0113     // input for SMASH
0114     auto soft_particlization = JetScapeSignalManager::Instance()->GetSoftParticlizationPointer().lock();
0115     if (soft_particlization) {
0116       soft_particlization->Hadron_list_.clear();
0117     }
0118   }
0119
0120   if (GetXMLElementInt({"Afterburner", "include_fragmentation_hadrons"})) {
0121     if (afterburner_had_events.size() > 1) {
0122       JSWARN << "Fragmentation hadrons in Afterburner are only possible without "
0123                 "repeated sampling from SoftParticlization. Exiting.";
0124       exit(1);
0125     }
0126     std::vector<shared_ptr<Hadron>> frag_hadrons = GetFragmentationHadrons();
0127
0128     if (GetXMLElementInt({"Afterburner", "output_only_final_state_hadrons"})) {
0129       // empty the hadron vector in the hadronization manager to circumvent the
0130       // output of these hadrons if they are implemented in the SMASH afterburner
0131       auto hadronization_mgr = JetScapeSignalManager::Instance()->GetHadronizationManagerPointer().lock();
0132       hadronization_mgr->DeleteRealHadrons();
0133     }
0134
0135     afterburner_had_events[0].insert(afterburner_had_events[0].end(),
0136                                      frag_hadrons.begin(), frag_hadrons.end());
0137     dummy.clear();
0138   }
0139   return afterburner_had_events;
0140 }
0141
0142 } // end namespace Jetscape