ThermalParticle.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 THERMALPARTICLE_H
#define THERMALPARTICLE_H


#include <string>
#include <vector>
#include <cmath>

#include "HRGBase/ParticleDecay.h"
#include "HRGBase/ThermalModelParameters.h"
#include "HRGBase/IdealGasFunctions.h"
#include "HRGBase/xMath.h"

namespace thermalfist {

  

  struct ConservedCharge {
    enum Name {
      BaryonCharge = 0,      
      ElectricCharge = 1,    
      StrangenessCharge = 2, 
      CharmCharge = 3        
    };
    static const int NumberOfTypes = 5;
  };

  

  class ThermalParticle
  {
  public:
    typedef std::vector<ParticleDecayChannel> ParticleDecaysVector;

    enum ResonanceWidthShape {
      RelativisticBreitWigner,
      NonRelativisticBreitWigner
    };

    enum ResonanceWidthIntegration {
      ZeroWidth,             
      BWTwoGamma,            
      FullInterval,          
      FullIntervalWeighted,  
      eBW,                   
      eBWconstBR             
    };

    ThermalParticle(bool Stable = true, std::string Name = "hadron", long long PDGID = 0, double Deg = 1., int Stat = 0, double Mass = 0.,
      int Strange = 0, int Baryon = 0, int Charge = 0, double AbsS = 0., double Width = 0., double Threshold = 0., int Charm = 0, double AbsC = 0., int Quark = 0);
    ~ThermalParticle(void);

    void FillCoefficients();

    void FillCoefficientsDynamical();

    double TotalWidtheBW(double M) const;

    std::vector<double> BranchingRatiosM(double M, bool eBW = true) const;

    double ThermalMassDistribution(double M, double T, double Mu, double width);

    double ThermalMassDistribution(double M, double T, double Mu);

    void NormalizeBranchingRatios();

    void RestoreBranchingRatios();

    double Density(const ThermalModelParameters &params, IdealGasFunctions::Quantity type = IdealGasFunctions::ParticleDensity, bool useWidth = 0, double mu = 0.) const;

    double DensityCluster(int n, const ThermalModelParameters &params, IdealGasFunctions::Quantity type = IdealGasFunctions::ParticleDensity, bool useWidth = 0, double mu = 0.) const;

    double chi(int index, const ThermalModelParameters &params, bool useWidth = 0, double mu = 0.) const;

    double ScaledVariance(const ThermalModelParameters &params, bool useWidth = 0, double mu = 0.) const;

    double Skewness(const ThermalModelParameters &params, bool useWidth = 0, double mu = 0.) const;

    double Kurtosis(const ThermalModelParameters &params, bool useWidth = 0, double mu = 0.) const;

    double FD(double k, double T, double mu, double m) const;
    
    double GetAbsQ() const;

    double GetCharge(int index) const;

    double GetAbsCharge(int index) const;

    bool IsNeutral() const;

    bool IsStable() const { return m_Stable; }

    void SetStable(bool stable = true) { m_Stable = stable; }

    bool IsAntiParticle() const { return m_AntiParticle; }

    void SetAntiParticle(bool antpar = true) { m_AntiParticle = antpar; }

    const std::string& Name() const { return m_Name; }

    void SetName(const std::string &name) { m_Name = name; }

    long long  PdgId() const { return m_PDGID; }

    void SetPdgId(long long PdgId) { m_PDGID = PdgId; }

    double Degeneracy() const { return m_Degeneracy; }

    void SetDegeneracy(double deg) { m_Degeneracy = deg; }

    int  Statistics() const { return m_Statistics; }

    void SetStatistics(int stat) { m_Statistics = stat; }

    void UseStatistics(bool enable);

    double Mass() const { return m_Mass; }

    void SetMass(double mass);// { m_Mass = mass; }

    int BaryonCharge() const { return m_Baryon; }

    void SetBaryonCharge(int chg) { m_Baryon = chg; SetAbsoluteQuark(GetAbsQ()); }

    int ElectricCharge() const { return m_ElectricCharge; }

    void SetElectricCharge(int chg) { m_ElectricCharge = chg; }

    int Strangeness() const { return m_Strangeness; }
    void SetStrangenessCharge(int chg) { m_Strangeness = chg; }

    int Charm() const { return m_Charm; }

    void SetCharm(int chg) { m_Charm = chg; }

    int ConservedCharge(ConservedCharge::Name chg) const;

    double ArbitraryCharge() const { return m_ArbitraryCharge; }

    void SetArbitraryCharge(double arbchg) { m_ArbitraryCharge = arbchg; }

    double AbsoluteQuark() const { return m_AbsQuark; }

    void SetAbsoluteQuark(double abschg) { m_AbsQuark = abschg; }

    double AbsoluteStrangeness() const { return m_AbsS; }

    void SetAbsoluteStrangeness(double abschg) { m_AbsS = abschg; SetAbsoluteQuark(GetAbsQ()); }

    double AbsoluteCharm() const { return m_AbsC; }

    void SetAbsoluteCharm(double abschg) { m_AbsC = abschg; SetAbsoluteQuark(GetAbsQ()); }

    bool ZeroWidthEnforced() const;

    double ResonanceWidth() const { return m_Width; }

    void SetResonanceWidth(double width);

    double DecayThresholdMass() const { return m_Threshold; }
    void SetDecayThresholdMass(double threshold);

    

    double DecayThresholdMassDynamical() const { return m_ThresholdDynamical; }

    void SetDecayThresholdMassDynamical(double threshold);

    void CalculateAndSetDynamicalThreshold();

    ResonanceWidthShape GetResonanceWidthShape() const { return m_ResonanceWidthShape; }

    void SetResonanceWidthShape(ResonanceWidthShape shape);

    ResonanceWidthIntegration GetResonanceWidthIntegrationType() const { return m_ResonanceWidthIntegrationType; }
    
    void SetResonanceWidthIntegrationType(ResonanceWidthIntegration type);

    double MassDistribution(double m) const;

    // Resonance mass distribution with manually input width
    double MassDistribution(double m, double width) const;

    double Weight() const { return m_Weight; }

    void SetWeight(double weight) { m_Weight = weight; }

    ParticleDecayType::DecayType DecayType() const { return m_DecayType; }

    void SetDecayType(ParticleDecayType::DecayType type) { m_DecayType = type; }

    const ParticleDecaysVector& Decays() const { return m_Decays; }

    ParticleDecaysVector& Decays() { return m_Decays; }

    void SetDecays(const ParticleDecaysVector &Decays) { m_Decays = Decays; }

    void ClearDecays() { m_Decays.resize(0); }

    const ParticleDecaysVector& DecaysOriginal() const { return m_DecaysOrig; }
    ParticleDecaysVector& DecaysOriginal() { return m_DecaysOrig; }
    void SetDecaysOriginal(const ParticleDecaysVector &DecaysOrig) { m_DecaysOrig = DecaysOrig; }

    void ReadDecays(std::string filename = "");

    
    void CalculateThermalBranchingRatios(const ThermalModelParameters &params, bool useWidth = 0, double mu = 0.);

    void SetCalculationType(IdealGasFunctions::QStatsCalculationType type) { m_QuantumStatisticsCalculationType = type; }
    
    IdealGasFunctions::QStatsCalculationType CalculationType()       const { return m_QuantumStatisticsCalculationType; }

    int ClusterExpansionOrder() const { return m_ClusterExpansionOrder; }

    void SetClusterExpansionOrder(int order) { m_ClusterExpansionOrder = order; }

    std::vector<double> BranchingRatioWeights(const std::vector<double> & ms) const;

    const std::vector<double>& Nch() const { return m_Nch; }
    std::vector<double>&  Nch() { return m_Nch; }

    const std::vector<double>& DeltaNch() const { return m_DeltaNch; }
    std::vector<double>&  DeltaNch() { return m_DeltaNch; }

    ThermalParticle GenerateAntiParticle(/*ThermalParticleSystem *TPS = NULL*/) const;

    bool operator==(const ThermalParticle &rhs) const; // TODO: improve
    bool operator!=(const ThermalParticle &rhs) const { return !(*this == rhs); }

  private:
    std::vector<double> m_xlag32, m_wlag32;
    std::vector<double> m_xleg, m_wleg;
    std::vector<double> m_xleg32, m_wleg32;
    std::vector<double> m_brweight;


    std::vector<double> m_xlegdyn, m_wlegdyn, m_vallegdyn;
    std::vector<double> m_xlegpdyn, m_wlegpdyn, m_vallegpdyn;
    std::vector<double> m_xlagdyn, m_wlagdyn, m_vallagdyn;

    std::vector<double> m_xalldyn, m_walldyn, m_densalldyn;


    bool m_Stable;                
    ParticleDecayType::DecayType m_DecayType;        
    bool m_AntiParticle;          
    std::string m_Name;           
    long long m_PDGID;            
    double m_Degeneracy;          
    int m_Statistics;             
    int m_StatisticsOrig;         
    double m_Mass;                
    IdealGasFunctions::QStatsCalculationType m_QuantumStatisticsCalculationType;

    int m_ClusterExpansionOrder;

    int m_Baryon;                 
    int m_ElectricCharge;         
    int m_Strangeness;            
    int m_Charm;                  
    int m_Quark;                  
    double m_ArbitraryCharge;     
    double m_AbsQuark;            
    double m_AbsS;                
    double m_AbsC;                
    double m_Width;               
    double m_Threshold;           
    double m_ThresholdDynamical;  
    ResonanceWidthShape m_ResonanceWidthShape;                  
    ResonanceWidthIntegration m_ResonanceWidthIntegrationType;  
    double m_Radius;              
    double m_Vo;                  
    double m_Weight;              
    ParticleDecaysVector m_Decays;      
    ParticleDecaysVector m_DecaysOrig;


    std::vector<double> m_Nch;
    std::vector<double> m_DeltaNch;

    bool m_LastDensityOk;
  };

} // namespace thermalfist

#endif