8 #include "HepMC3TestUtils.h"
13 using namespace HepMC3;
14 const Int_t kMaxparticles = 2000;
15 const Int_t kMaxvertices = 2000;
16 #ifndef DOXYGEN_SHOULD_SKIP_THIS
27 Int_t particles_pid[kMaxparticles];
28 Int_t particles_status[kMaxparticles];
29 Bool_t particles_is_mass_set[kMaxparticles];
30 Double_t particles_mass[kMaxparticles];
31 Double_t particles_momentum_m_v1[kMaxparticles];
32 Double_t particles_momentum_m_v2[kMaxparticles];
33 Double_t particles_momentum_m_v3[kMaxparticles];
34 Double_t particles_momentum_m_v4[kMaxparticles];
36 Int_t vertices_status[kMaxvertices];
37 Double_t vertices_position_m_v1[kMaxvertices];
38 Double_t vertices_position_m_v2[kMaxvertices];
39 Double_t vertices_position_m_v3[kMaxvertices];
40 Double_t vertices_position_m_v4[kMaxvertices];
41 vector<double> weights;
42 Double_t event_pos_m_v1;
43 Double_t event_pos_m_v2;
44 Double_t event_pos_m_v3;
45 Double_t event_pos_m_v4;
48 vector<int> attribute_id;
49 vector<string> attribute_name;
50 vector<string> attribute_string;
53 TBranch *b_hepmc3_event_event_number;
54 TBranch *b_hepmc3_event_momentum_unit;
55 TBranch *b_hepmc3_event_length_unit;
56 TBranch *b_hepmc3_event_particles_;
57 TBranch *b_particles_pid;
58 TBranch *b_particles_status;
59 TBranch *b_particles_is_mass_set;
60 TBranch *b_particles_mass;
61 TBranch *b_particles_momentum_m_v1;
62 TBranch *b_particles_momentum_m_v2;
63 TBranch *b_particles_momentum_m_v3;
64 TBranch *b_particles_momentum_m_v4;
65 TBranch *b_hepmc3_event_vertices_;
66 TBranch *b_vertices_status;
67 TBranch *b_vertices_position_m_v1;
68 TBranch *b_vertices_position_m_v2;
69 TBranch *b_vertices_position_m_v3;
70 TBranch *b_vertices_position_m_v4;
71 TBranch *b_hepmc3_event_weights;
72 TBranch *b_hepmc3_event_event_pos_m_v1;
73 TBranch *b_hepmc3_event_event_pos_m_v2;
74 TBranch *b_hepmc3_event_event_pos_m_v3;
75 TBranch *b_hepmc3_event_event_pos_m_v4;
76 TBranch *b_hepmc3_event_links1;
77 TBranch *b_hepmc3_event_links2;
78 TBranch *b_hepmc3_event_attribute_id;
79 TBranch *b_hepmc3_event_attribute_name;
80 TBranch *b_hepmc3_event_attribute_string;
83 void Init(TChain *tree)
87 fChain->SetMakeClass(1);
89 fChain->SetBranchAddress(
"event_number", &event_number, &b_hepmc3_event_event_number);
90 fChain->SetBranchAddress(
"momentum_unit", &momentum_unit, &b_hepmc3_event_momentum_unit);
91 fChain->SetBranchAddress(
"length_unit", &length_unit, &b_hepmc3_event_length_unit);
92 fChain->SetBranchAddress(
"particles", &particles_, &b_hepmc3_event_particles_);
93 fChain->SetBranchAddress(
"particles.pid", particles_pid, &b_particles_pid);
94 fChain->SetBranchAddress(
"particles.status", particles_status, &b_particles_status);
95 fChain->SetBranchAddress(
"particles.is_mass_set", particles_is_mass_set, &b_particles_is_mass_set);
96 fChain->SetBranchAddress(
"particles.mass", particles_mass, &b_particles_mass);
97 fChain->SetBranchAddress(
"particles.momentum.m_v1", particles_momentum_m_v1, &b_particles_momentum_m_v1);
98 fChain->SetBranchAddress(
"particles.momentum.m_v2", particles_momentum_m_v2, &b_particles_momentum_m_v2);
99 fChain->SetBranchAddress(
"particles.momentum.m_v3", particles_momentum_m_v3, &b_particles_momentum_m_v3);
100 fChain->SetBranchAddress(
"particles.momentum.m_v4", particles_momentum_m_v4, &b_particles_momentum_m_v4);
101 fChain->SetBranchAddress(
"vertices", &vertices_, &b_hepmc3_event_vertices_);
102 fChain->SetBranchAddress(
"vertices.status", vertices_status, &b_vertices_status);
103 fChain->SetBranchAddress(
"vertices.position.m_v1", vertices_position_m_v1, &b_vertices_position_m_v1);
104 fChain->SetBranchAddress(
"vertices.position.m_v2", vertices_position_m_v2, &b_vertices_position_m_v2);
105 fChain->SetBranchAddress(
"vertices.position.m_v3", vertices_position_m_v3, &b_vertices_position_m_v3);
106 fChain->SetBranchAddress(
"vertices.position.m_v4", vertices_position_m_v4, &b_vertices_position_m_v4);
107 fChain->SetBranchAddress(
"weights", &weights, &b_hepmc3_event_weights);
108 fChain->SetBranchAddress(
"event_pos.m_v1", &event_pos_m_v1, &b_hepmc3_event_event_pos_m_v1);
109 fChain->SetBranchAddress(
"event_pos.m_v2", &event_pos_m_v2, &b_hepmc3_event_event_pos_m_v2);
110 fChain->SetBranchAddress(
"event_pos.m_v3", &event_pos_m_v3, &b_hepmc3_event_event_pos_m_v3);
111 fChain->SetBranchAddress(
"event_pos.m_v4", &event_pos_m_v4, &b_hepmc3_event_event_pos_m_v4);
112 fChain->SetBranchAddress(
"links1", &links1, &b_hepmc3_event_links1);
113 fChain->SetBranchAddress(
"links2", &links2, &b_hepmc3_event_links2);
114 fChain->SetBranchAddress(
"attribute_id", &attribute_id, &b_hepmc3_event_attribute_id);
115 fChain->SetBranchAddress(
"attribute_name", &attribute_name, &b_hepmc3_event_attribute_name);
116 fChain->SetBranchAddress(
"attribute_string", &attribute_string, &b_hepmc3_event_attribute_string);
118 fChain->SetBranchStatus(
"*",0);
119 fChain->SetBranchStatus(
"particles",1);
120 fChain->SetBranchStatus(
"particles.pid",1);
121 fChain->SetBranchStatus(
"particles.status",1);
122 fChain->SetBranchStatus(
"particles.momentum.m_v1",1);
123 fChain->SetBranchStatus(
"particles.momentum.m_v2",1);
126 SomeAnalysis(
const std::string& file)
128 TChain* TempChain=
new TChain(
"hepmc3_tree");
129 TempChain->Add(file.c_str());
138 TH1D H1(
"H1",
"Pt of pions;Events/100MeV;P_{T},GeV",1000,0,100);
139 SomeAnalysis* A=
new SomeAnalysis(
"inputIO4.root");
140 if (!A->fChain->GetEntries())
return 10001;
141 for (
int entry=0; entry<A->fChain->GetEntries(); entry++)
143 A->fChain->GetEntry(entry);
144 for (
int i=0; i<A->particles_; i++)
145 if (A->particles_status[i]==1&&(
std::abs(A->particles_pid[i])==211||
std::abs(A->particles_pid[i])==11))
146 H1.Fill(std::sqrt(A->particles_momentum_m_v1[i]*A->particles_momentum_m_v1[i]+A->particles_momentum_m_v2[i]*A->particles_momentum_m_v2[i]) );
150 TH1D H2(
"H2",
"Pt of pions;Events/100MeV;P_{T},GeV",1000,0,100);
152 if(inputA.failed())
return 10002;
153 while( !inputA.failed() )
156 inputA.read_event(evt);
157 if( inputA.failed() ) {printf(
"End of file reached. Exit.\n");
break;}
158 for (ConstGenParticlePtr p: evt.particles())
160 H2.Fill( p->momentum().perp());
166 for (
int i=0; i<H1.GetNbinsX(); i++)
168 double eps=
std::abs(H1.GetBinContent(i)-H2.GetBinContent(i));
169 if (eps<1e-5)
continue;
170 std::cout<<
"Bin: "<<i<<
" "<<H1.GetBinContent(i)<<
" "<<H2.GetBinContent(i)<<std::endl;
Definition of class GenParticle.
Definition of class WriterRootTree.
Definition of class ReaderRootTree.
Definition of class ReaderAsciiHepMC2.
Stores event-related information.
Definition of class WriterAsciiHepMC2.
GenEvent I/O parsing and serialization for root files based on root TTree.
int main(int argc, char **argv)
Definition of class GenEvent.
Feature< Feature_type > abs(const Feature< Feature_type > &input)
Obtain the absolute value of a Feature. This works as you'd expect. If foo is a valid Feature...