123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124 |
- /**
- * MpdHistoGenerator.cxx
- * *@author A.Zinchenko <Alexander.Zinchenko@jinr.ru>
- * The MpdHistoGenerator produces particles according to their Eta-Pt histogram
- **/
- #include "MpdHistoGenerator.h"
- #include "FairMCEventHeader.h"
- #include "FairPrimaryGenerator.h"
- #include "TDatabasePDG.h"
- #include "TFile.h"
- #include "TH2D.h"
- #include "TLorentzVector.h"
- #include "TMath.h"
- #include "TRandom.h"
- #include "TVirtualMC.h"
- #include <fstream>
- #include <iostream>
- #include <iomanip>
- // #define debug_hsd
- using namespace std;
- // ------------------------------------------------------------------------
- MpdHistoGenerator::MpdHistoGenerator() :
- FairGenerator(), fPdgCode(-211), fMult(1), fHist(NULL), fYield(-1.0),
- fFileName("")
- {
- // Constructor
- //gRandom->SetSeed(0); // change random seed
- // Get histograms Pt vs Eta (He3)
- //TFile f("he3_eta_pt.root");
- //f.ReadAll();
- //fHist = (TH2D*) f.FindObject("hPtEta");
- //fHist->SetDirectory(0);
- }
- // ------------------------------------------------------------------------
- MpdHistoGenerator::MpdHistoGenerator(Int_t pdgid, Int_t mult, Double_t yield) :
- FairGenerator(), fPdgCode(pdgid), fMult(mult), fHist(NULL), fYield(yield),
- fFileName("")
- {
- // Constructor
- //gRandom->SetSeed(0); // change random seed
- // Get histograms Pt vs Eta (He3)
- //TFile f("he3_eta_pt.root");
- //f.ReadAll();
- //fHist = (TH2D*) f.FindObject("hPtVsEta");
- //fHist->SetDirectory(0);
- }
- // ------------------------------------------------------------------------
- Bool_t MpdHistoGenerator::Init()
- {
- // Get histograms Pt vs Eta (He3)
- if (fFileName == "") fFileName = "genPtVsEta.root";
- TFile *f = TFile::Open(fFileName);
- f->ReadAll();
- fHist = (TH2D*) f->FindObject("hPtVsEta");
- fHist->SetDirectory(0);
- return kTRUE;
- }
- // ------------------------------------------------------------------------
- MpdHistoGenerator::~MpdHistoGenerator()
- {
- // Destructor
- }
- //-------------------------------------------------------------------------
- Bool_t MpdHistoGenerator::ReadEvent(FairPrimaryGenerator* primGen)
- {
- // Generate event
- static Int_t first = 1;
- if (first) {
- first = 0;
- //cout << " histo " << fHist << " " << fHist->GetNbinsX() << " " << fHist->GetNbinsY() << endl;
- if (TDatabasePDG::Instance()->GetParticle(fPdgCode) == NULL) {
- // Define particle with fPdgCode (take it from DCM-QGSM input file)
- // To be done ...
- }
- }
- if (fYield > 0 && gRandom->Rndm() > fYield) return kTRUE;
- Double_t pt, eta, phi;
- TVector3 part;
- for (Int_t i = 0; i < fMult; ++i) {
- phi = gRandom->Uniform(TMath::TwoPi());
- fHist->GetRandom2(eta, pt);
-
- part.SetPtEtaPhi(pt, eta, phi);
- //FairMCEventHeader* eventHeader = primGen->GetEvent();
- /* Set event impact parameter in MCEvent if not yet done */
- /*
- if (eventHeader && (!eventHeader->IsSet()) ) {
- eventHeader->SetB(fB);
- eventHeader->MarkSet(kTRUE);
- }
- */
- //primGen->AddTrack(1000020030, He3.Px(), He3.Py(), He3.Pz(), 0., 0., 0., -1);
- primGen->AddTrack(fPdgCode, part.Px(), part.Py(), part.Pz(), 0., 0., 0., -1);
- }
- return kTRUE;
- }
- // ------------------------------------------------------------------------
|