Back to home page

Project CMSSW displayed by LXR

 
 

    


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