Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-04-06 12:13:48

0001 #ifndef HydjetHadronizer_h
0002 #define HydjetHadronizer_h
0003 
0004 /**
0005    \class HydjetHadronizer
0006    \brief Interface to the HYDJET generator (since core v. 1.9.1), produces HepMC events
0007    \version 2.0
0008    \authors Camelia Mironov, Yetkin Yilmaz, Andrey Belyaev
0009 */
0010 
0011 #define PYCOMP pycomp_
0012 
0013 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0014 #include "GeneratorInterface/Core/interface/BaseHadronizer.h"
0015 
0016 #include "FWCore/Utilities/interface/EDGetToken.h"
0017 #include "FWCore/Framework/interface/ConsumesCollector.h"
0018 #include "SimDataFormats/CrossingFrame/interface/MixCollection.h"
0019 #include "SimDataFormats/PileupSummaryInfo/interface/PileupMixingContent.h"
0020 
0021 #include <map>
0022 #include <string>
0023 #include <vector>
0024 #include <cmath>
0025 
0026 namespace CLHEP {
0027   class HepRandomEngine;
0028 }
0029 
0030 namespace HepMC {
0031   class GenEvent;
0032   class GenParticle;
0033   class GenVertex;
0034   class FourVector;
0035 }  // namespace HepMC
0036 
0037 namespace gen {
0038   class Pythia6Service;
0039 
0040   class HydjetHadronizer : public BaseHadronizer {
0041   public:
0042     HydjetHadronizer(const edm::ParameterSet&, edm::ConsumesCollector&&);
0043     ~HydjetHadronizer() override;
0044 
0045     bool generatePartonsAndHadronize();
0046     bool hadronize();
0047     bool decay();
0048     bool residualDecay();
0049     bool readSettings(int);
0050     bool initializeForExternalPartons();
0051     bool initializeForInternalPartons();
0052     bool declareStableParticles(const std::vector<int>&);
0053     bool declareSpecialSettings(const std::vector<std::string>&) { return true; }
0054 
0055     void finalizeEvent();
0056     void statistics();
0057     const char* classname() const;
0058 
0059   private:
0060     void doSetRandomEngine(CLHEP::HepRandomEngine* v) override;
0061     std::vector<std::string> const& doSharedResources() const override { return theSharedResources; }
0062 
0063     static const std::vector<std::string> theSharedResources;
0064 
0065     void add_heavy_ion_rec(HepMC::GenEvent* evt);
0066     HepMC::GenParticle* build_hyjet(int index, int barcode);
0067     HepMC::GenVertex* build_hyjet_vertex(int i, int id);
0068     bool get_particles(HepMC::GenEvent* evt);
0069     bool call_hyinit(double energy, double a, int ifb, double bmin, double bmax, double bfix, int nh);
0070     bool hydjet_init(const edm::ParameterSet& pset);
0071     inline double nuclear_radius() const;
0072     void rotateEvtPlane();
0073     int convertStatus(int st);
0074 
0075     HepMC::GenEvent* evt;
0076     edm::ParameterSet pset_;
0077     double abeamtarget_;                  ///< beam/target atomic mass number
0078     int angularspecselector_;             ///< angular emitted gluon spectrum selection
0079     double bfixed_;                       ///< fixed impact param (fm); valid only if cflag_=0
0080     double bmax_;                         ///< max impact param;
0081                                           ///< units of nucl radius
0082     double bmin_;                         ///< min impact param;
0083                                           ///< units of nucl radius
0084     int cflag_;                           ///< centrality flag
0085                                           ///< =  0 fixed impact param,
0086                                           ///< <> 0 between bmin and bmax
0087     bool embedding_;                      ///< Switch for embedding mode
0088     double comenergy;                     ///< collision energy
0089     bool doradiativeenloss_;              ///< DEFAULT = true
0090     bool docollisionalenloss_;            ///< DEFAULT = true
0091     double fracsoftmult_;                 ///< fraction of soft hydro induced hadronic multiplicity
0092                                           ///< proportional to no of nucleon participants
0093                                           ///< (1-fracsoftmult_)--- fraction of soft
0094                                           ///< multiplicity proportional to the numebr
0095                                           ///< of nucleon-nucleon binary collisions
0096                                           ///< DEFAULT=1., allowed range [0.01,1]
0097     double hadfreeztemp_;                 ///< hadron freez-out temperature
0098                                           ///< DEFAULT=0.14MeV, allowed ranges [0.08,0.2]MeV
0099     std::string hymode_;                  ///< Hydjet running mode
0100     unsigned int maxEventsToPrint_;       ///< Events to print if verbosity
0101     double maxlongy_;                     ///< max longitudinal collective rapidity:
0102                                           ///< controls width of eta-spectra
0103                                           ///< DEFAULT=4, allowed range [0.01,7.0]
0104     double maxtrany_;                     ///< max transverse collective rapidity:
0105                                           ///< controls slope of low-pt spectra
0106                                           ///< DEFAULT=1.5, allowed range [0.01,3.0]
0107     int nsub_;                            ///< number of sub-events
0108     int nhard_;                           ///< multiplicity of PYTHIA(+PYQUEN)-induced particles in event
0109     int nmultiplicity_;                   ///< mean soft multiplicity in central PbPb
0110                                           ///< automatically calculated for other centralitie and beams
0111     int nsoft_;                           ///< multiplicity of HYDRO-induced particles in event
0112     unsigned int nquarkflavor_;           ///< number of active quark flavors in qgp
0113                                           ///< DEFAULT=0; allowed values: 0,1,2,3.
0114     unsigned int pythiaPylistVerbosity_;  ///< pythia verbosity; def=1
0115     double qgpt0_;                        ///< initial temperature of QGP
0116                                           ///< DEFAULT = 1GeV; allowed range [0.2,2.0]GeV;
0117     double qgptau0_;                      ///< proper time of QGP formation
0118                                           ///< DEFAULT = 0.1 fm/c; allowed range [0.01,10.0]fm/
0119 
0120     double phi0_;  ///< Event plane angle
0121     double sinphi0_;
0122     double cosphi0_;
0123     bool rotate_;  ///< Switch to rotate event plane
0124 
0125     unsigned int shadowingswitch_;   ///< shadowing switcher
0126                                      ///< 1-ON, 0-OFF
0127     double signn_;                   ///< inelastic nucleon nucleon cross section [mb]
0128                                      ///< DEFAULT= 58 mb
0129     HepMC::FourVector* fVertex_;     ///< Event signal vertex
0130     std::vector<double> signalVtx_;  ///< Pset double vector to set event signal vertex
0131 
0132     Pythia6Service* pythia6Service_;
0133     edm::EDGetTokenT<CrossingFrame<edm::HepMCProduct> > src_;
0134   };
0135 
0136   /**
0137     Return the nuclear radius derived from the
0138     beam/target atomic mass number.
0139  */
0140   double HydjetHadronizer::nuclear_radius() const { return 1.15 * pow((double)abeamtarget_, 1. / 3.); }
0141 }  // namespace gen
0142 #endif