EventGeneratorBase.h

Go to the documentation of this file.

/*
 * Thermal-FIST package
 * 
 * Copyright (c) 2014-2019 Volodymyr Vovchenko
 *
 * GNU General Public License (GPLv3 or later)
 */
#ifndef EVENTGENERATORBASE_H
#define EVENTGENERATORBASE_H


#include <sstream>

#include "HRGEventGenerator/SimpleEvent.h"
#include "HRGEventGenerator/Acceptance.h"
#include "HRGEventGenerator/RandomGenerators.h"
#include "HRGBase/xMath.h"
#include "HRGBase/ThermalModelBase.h"

namespace thermalfist {

  template <typename T>
  std::string to_string_fix(T value)
  {
    //create an output string stream
    std::ostringstream os;

    //throw the value into the string stream
    os << value;

    //convert the string stream into a string and return
    return os.str();
  }

  struct EventGeneratorConfiguration {
    enum Ensemble { 
      GCE, 
      CE,  
      SCE, 
      CCE  
    };

    enum ModelType { 
      PointParticle, 
      DiagonalEV,    
      CrosstermsEV,  
      MeanFieldEV,   
      QvdW           
    };
    
    Ensemble fEnsemble;
    
    ModelType fModelType;
    
    ThermalModelParameters CFOParameters;
    
    int B, Q, S, C;

    bool CanonicalB;
    bool CanonicalQ;
    bool CanonicalS;
    bool CanonicalC;

    std::vector< std::vector<double> > bij;

    std::vector< std::vector<double> > aij;

    bool fUsePCE;

    std::vector<double> fPCEChems;

    EventGeneratorConfiguration();
  };

  class EventGeneratorBase
  {
  public:
    EventGeneratorBase() { m_THM = NULL; fCEAccepted = fCETotal = 0; }

    virtual ~EventGeneratorBase();

    void ClearMomentumGenerators();

    void SetCollisionKineticEnergy(double ekin) {
      SetCollisionCMSEnergy(sqrt(2.*xMath::mnucleon()*(ekin + 2. * xMath::mnucleon())));
    }

    void SetCollisionLabEnergy(double elab) {
      SetCollisionCMSEnergy(sqrt(2.*xMath::mnucleon()*(elab + xMath::mnucleon())));
    }

    void SetCollisionCMSEnergy(double ssqrt) {
      m_ssqrt = ssqrt;
      m_ekin = m_ssqrt * m_ssqrt / 2. / xMath::mnucleon() - 2. * xMath::mnucleon();
      m_elab = xMath::mnucleon() + m_ekin;
      double plab = sqrt(m_elab*m_elab - xMath::mnucleon() * xMath::mnucleon());
      m_ycm = 0.5 * log((m_elab + xMath::mnucleon() + plab) / (m_elab + xMath::mnucleon() - plab));
    }

    //std::vector<Acceptance::AcceptanceFunction>& GetAcceptance() { return m_acc; }

    //virtual void ReadAcceptance(std::string accfolder);

    double getYcm() const { return m_ycm; }

    std::pair< std::vector<int>, double > SampleYields() const;

    virtual SimpleEvent SampleMomenta(const std::vector<int>& yields) const;

    virtual SimpleEvent GetEvent(bool PerformDecays = true) const;

    static SimpleEvent PerformDecays(const SimpleEvent& evtin, ThermalParticleSystem* TPS);

    

    static int fCEAccepted, fCETotal;

    void SetVolume(double V);

    void RescaleCEMeans(double Vmod);


  protected:
    void SetConfiguration(ThermalParticleSystem *TPS,
      const EventGeneratorConfiguration& config);

    void PrepareMultinomials();
    
    std::vector<int> GenerateTotals() const;

    std::vector<int> GenerateTotalsGCE() const;

    std::vector<int> GenerateTotalsCE() const;

    std::vector<int> GenerateTotalsSCE() const;

    std::vector<int> GenerateTotalsSCESubVolume(double VolumeSC) const;

    std::vector<int> GenerateTotalsCCE() const;

    std::vector<int> GenerateTotalsCCESubVolume(double VolumeSC) const;

    EventGeneratorConfiguration m_Config;
    ThermalModelBase *m_THM;

    std::vector<double> m_DensitiesIdeal;

    std::vector<RandomGenerators::ParticleMomentumGenerator*>    m_MomentumGens;

    std::vector<RandomGenerators::ThermalBreitWignerGenerator*>  m_BWGens;

  private:

    //static SimpleEvent PerformDecaysAlternativeWay(const SimpleEvent& evtin, ThermalParticleSystem* TPS);

    double m_ekin, m_ycm, m_ssqrt, m_elab;
    // Acceptance discontinued
    //std::vector<Acceptance::AcceptanceFunction> m_acc;

    std::vector< std::pair<double, int> > m_Baryons;
    std::vector< std::pair<double, int> > m_AntiBaryons;
    std::vector< std::pair<double, int> > m_StrangeMesons;
    std::vector< std::pair<double, int> > m_AntiStrangeMesons;
    std::vector< std::pair<double, int> > m_ChargeMesons;
    std::vector< std::pair<double, int> > m_AntiChargeMesons;
    std::vector< std::pair<double, int> > m_CharmMesons;
    std::vector< std::pair<double, int> > m_AntiCharmMesons;
    std::vector< std::pair<double, int> > m_CharmAll;
    std::vector< std::pair<double, int> > m_AntiCharmAll;

    std::vector<double> m_BaryonsProbs;
    std::vector<double> m_AntiBaryonsProbs;
    std::vector<double> m_StrangeMesonsProbs;
    std::vector<double> m_AntiStrangeMesonsProbs;
    std::vector<double> m_ChargeMesonsProbs;
    std::vector<double> m_AntiChargeMesonsProbs;
    std::vector<double> m_CharmMesonsProbs;
    std::vector<double> m_AntiCharmMesonsProbs;
    std::vector<double> m_CharmAllProbs;
    std::vector<double> m_AntiCharmAllProbs;

    double m_MeanB, m_MeanAB;
    double m_MeanSM, m_MeanASM;
    double m_MeanCM, m_MeanACM; 
    double m_MeanCHRMM, m_MeanACHRMM;
    double m_MeanCHRM, m_MeanACHRM;

    static double m_LastWeight;
    static double m_LastLogWeight;
    static double m_LastNormWeight;
  };

} // namespace thermalfist

#endif