obsolete/ak/hepmc_writer_tof_pion.C

0001 //

0002 // export LD_LIBRARY_PATH=/home/eic/lib64:${LD_LIBRARY_PATH}

0003 //

0004 // root -l

0005 //

0006 // root [0] gSystem->AddIncludePath("-I/home/eic/include");

0007 // root [1] .x hepmc_writer_tof_pion.C+("out.hepmc", 1000)

0008 //

0009
0010 #include "HepMC3/GenEvent.h"
0011 #include "HepMC3/ReaderAscii.h"
0012 #include "HepMC3/WriterAscii.h"
0013 #include "HepMC3/Print.h"
0014
0015 #include <iostream>
0016 #include <random>
0017 #include <cmath>
0018 #include <math.h>
0019
0020 #include <TMath.h>
0021 #include <TRandom.h>
0022 #include <TDatabasePDG.h>
0023
0024 using namespace HepMC3;
0025
0026 void hepmc_writer_tof_pion(const char* out_fname, int n_events)
0027 {
0028   auto *DatabasePDG = new TDatabasePDG();
0029   auto *pion = DatabasePDG->GetParticle(211);//, *kaon = DatabasePDG->GetParticle(321);

0030
0031   WriterAscii hepmc_output(out_fname);
0032   int events_parsed = 0;
0033   GenEvent evt(Units::GEV, Units::MM);
0034
0035   // Random number generator

0036   TRandom *rdmn_gen = new TRandom(0x12345678);
0037
0038   for (events_parsed = 0; events_parsed < n_events; events_parsed++) {
0039     // type 4 is beam

0040     GenParticlePtr p1 =
0041         std::make_shared<GenParticle>(FourVector(0.0, 0.0, 12.0, 12.0), 11, 4);
0042     GenParticlePtr p2 = std::make_shared<GenParticle>(
0043         FourVector(0.0, 0.0, 100.0, 100.004), 2212, 4);
0044
0045     GenVertexPtr v1 = std::make_shared<GenVertex>();
0046     v1->add_particle_in(p1);
0047     v1->add_particle_in(p2);
0048
0049     for(unsigned iq=0; iq</*2*/1; iq++) {
0050       Double_t eta   = -2.5;//rdmn_gen->Uniform(-2.5, -2.30);

0051       Double_t th    = 2*std::atan(exp(-eta));
0052       Double_t p     = rdmn_gen->Uniform(7.0, 7.1);
0053       Double_t phi   = 90.0*M_PI/180.;//rdmn_gen->Uniform(90.0*M_PI/180., 100.*M_PI/180.);

0054
0055       Double_t px    = p * std::cos(phi) * std::sin(th);
0056       Double_t py    = p * std::sin(phi) * std::sin(th);
0057       Double_t pz    = p * std::cos(th);
0058
0059       //auto particle = pion;//events_parsed%2 ? pion : kaon;

0060       auto particle = pion;//iq%2 ? pion : kaon;

0061       GenParticlePtr pq = std::make_shared<GenParticle>(FourVector(
0062                                    px, py, pz,
0063                                    sqrt(p*p + pow(particle->Mass(), 2))),
0064                             particle->PdgCode(), 1);
0065       v1->add_particle_out(pq);
0066     } //for iq

0067
0068     evt.add_vertex(v1);
0069
0070     if (events_parsed == 0) {
0071       std::cout << "First event: " << std::endl;
0072       Print::listing(evt);
0073     }
0074
0075     hepmc_output.write_event(evt);
0076     if (events_parsed % 10000 == 0) {
0077       std::cout << "Event: " << events_parsed << std::endl;
0078     }
0079     evt.clear();
0080   }
0081   hepmc_output.close();
0082   std::cout << "Events parsed and written: " << events_parsed << std::endl;
0083   exit(0);
0084 }