2 #include "HepMC3/ReaderFactory.h"
3 #include "HepMC3ViewerFrame.h"
11 #include <graphviz/gvc.h>
12 #define CONSERVATION_TOLERANCE 1e-5
14 static char* create_image_from_dot(
char* m_buffer)
16 GVC_t * gvc=gvContext();
17 Agraph_t *g= agmemread(m_buffer);
18 gvLayout(gvc,g,
"dot");
26 err = gvRenderData(gvc, g,
"png", &data, &length);
29 data = (
char*)realloc(data, length + 1);
35 static bool show_as_parton(HepMC3::ConstGenParticlePtr p )
40 if (pd==81||pd==82||pd<25) parton=
true;
42 (pd/1000==1||pd/1000==2||pd/1000==3||pd/1000==4||pd/1000==5)
43 &&(pd%1000/100==1||pd%1000/100==2||pd%1000/100==3||pd%1000/100==4)
44 &&(pd%100==1||pd%100==3)
47 if (p->status()==4) parton=
true;
51 static char* write_event_to_dot(
char* used_cursor,
const HepMC3::GenEvent &evt,
int used_style=1)
53 used_cursor += sprintf(used_cursor,
"digraph graphname%d {\n",evt.
event_number());
54 used_cursor += sprintf(used_cursor,
"v0[label=\"Machine\"];\n");
61 if (v->status()==2) used_cursor += sprintf(used_cursor,
"node [color=\"green\"];\n");
62 else used_cursor += sprintf(used_cursor,
"node [color=\"black\"];\n");
68 for(
auto p1: v->particles_in() ) {
70 energy+=
std::abs(p1->momentum().e());
72 for(
auto p2: v->particles_out() ) {
74 energy+=
std::abs(p2->momentum().e());
79 double energyviolation=std::sqrt(momviolation.length2() +momviolation.e()*momviolation.e() );
81 if (energyviolation>CONSERVATION_TOLERANCE*energy) violation=
true;
85 used_cursor += sprintf(used_cursor,
"node [shape=rectangle];\n");
86 used_cursor += sprintf(used_cursor,
"v%d [label=\"%d\nd=%4.2f\"];\n", -v->id(),v->id(),energyviolation);
90 used_cursor += sprintf(used_cursor,
"node [shape=ellipse];\n");
91 used_cursor += sprintf(used_cursor,
"v%d[label=\"%d\"];\n", -v->id(),v->id());
94 used_cursor += sprintf(used_cursor,
"node [shape=ellipse];\n");
96 for(
auto p: evt.
beams() )
98 if (!p->end_vertex())
continue;
99 used_cursor += sprintf(used_cursor,
"node [shape=point];\n");
100 used_cursor += sprintf(used_cursor,
"v0 -> v%d [label=\"%d(%d)\"];\n", -p->end_vertex()->id(),p->id(),p->pid());
106 for(
auto p: v->particles_out() )
113 if (show_as_parton(p)&&p->status()!=1) used_cursor += sprintf(used_cursor,
"edge [color=\"red\"];\n");
114 else used_cursor +=sprintf(used_cursor,
"edge [color=\"black\"];\n");
117 if (!p->end_vertex())
119 used_cursor += sprintf(used_cursor,
"node [shape=point];\n");
120 used_cursor += sprintf(used_cursor,
"v%d -> o%d [label=\"%d(%d)\"];\n", -v->id(),p->id(),p->id(),p->pid());
124 used_cursor += sprintf(used_cursor,
"v%d -> v%d [label=\"%d(%d)\"];\n", -v->id(),-p->end_vertex()->id(),p->id(),p->pid());
128 used_cursor += sprintf(used_cursor,
"labelloc=\"t\";\nlabel=\"Event %d; Vertices %lu; Particles %lu;\";\n", evt.
event_number(), evt.
vertices().size(), evt.
particles().size());
129 used_cursor += sprintf(used_cursor,
"}\n\n");
138 char* m_cursor=m_buffer;
140 char *buf=create_image_from_dot(m_buffer);
155 TPad *p1 =
new TPad(
"i1",
"i1", 0.05, 0.05, 0.05+d*
fGraphImage->GetWidth()/
fGraphImage->GetHeight(), 0.95);
174 TH1S* particles1=
new TH1S();
175 fAnalysisH[
"particles1"]=particles1;
176 particles1->SetTitle(
"Flavour: all particles; PDG ID; Number of particles");
177 particles1->SetFillColor(kBlue);
179 particles1->Fill((std::to_string(p->pid())).c_str(),1.0);
180 particles1->LabelsOption(
">",
"X");
182 TH1S* particles2=
new TH1S();
183 fAnalysisH[
"particles2"]=particles2;
184 particles2->SetTitle(
"Flavour: particles with status 1; PDG ID; Number of particles");
185 particles2->SetFillColor(kBlue);
187 if(p->status()==1) particles2->Fill((std::to_string(p->pid())).c_str(),1.0);
188 particles2->LabelsOption(
">",
"X");
190 std::vector<double> masses;
192 if(show_as_parton(p)) masses.push_back(p->momentum().m());
193 TH1D* particles3=
new TH1D(
"particles3",
"Mass: parton particles; Mass, GeV; Number of particles",masses.size(),
194 0,*std::max_element(masses.begin(),masses.end()));
195 fAnalysisH[
"particles3"]=particles3;
196 particles3->SetFillColor(kBlue);
197 for(
auto m: masses) particles3->Fill(m);
201 TPad *p1 =
new TPad(
"i1",
"i1", 0.00, 0.75, 1.0, 1.0);
206 TPad *p2 =
new TPad(
"i2",
"i2", 0.00, 0.50, 1.0, 0.75);
211 TPad *p3 =
new TPad(
"i3",
"i3", 0.00, 0.25, 1.0, 0.50);
231 if (pos==
fEventsCache.end()) printf(
"This event was not found in the cache. Cache size is %zu\n",
fEventsCache.size());
244 bool ok=
fReader->read_event(*(evt1));
263 static const char *FileType[] = {
"All",
"*.*",
"HepMC",
"*.hepmc*",
"LHEF",
"*.lhe*",
"ROOT",
"*.root", 0, 0 };
264 static TString dir(
"./");
266 fi.fFileTypes = FileType;
267 fi.fIniDir = StrDup(dir);
268 new TGFileDialog(gClient->GetRoot(),
this, kFDOpen, &fi);
276 fMainFrame =
new TGCompositeFrame(
this, 1350, 500, kHorizontalFrame|kFixedWidth);
290 fChooseInput->Connect(
"Clicked()",
"HepMC3ViewerFrame",
this,
"ChooseInput()");
297 fNextEvent->Connect(
"Clicked()",
"HepMC3ViewerFrame",
this,
"NextEvent()");
298 fNextEvent->SetToolTipText(
"Click to display next event");
303 fPreviousEvent->Connect(
"Clicked()",
"HepMC3ViewerFrame",
this,
"PreviousEvent()");
304 fPreviousEvent->SetToolTipText(
"Click to display previous event");
309 fClearEventCache->Connect(
"Clicked()",
"HepMC3ViewerFrame",
this,
"ClearEventCache()");
314 fExit->SetToolTipText(
"Click to exit");
315 fButtonFrame->AddFrame(
fExit,
new TGLayoutHints( kLHintsExpandX|kLHintsLeft,1,1,2,2));
317 AddFrame(
fMainFrame,
new TGLayoutHints(kLHintsTop |kLHintsExpandX| kLHintsExpandY, 1, 1, 2, 2));
319 SetWindowName(
"Event viewer");
321 Resize(GetDefaultSize());
const std::vector< ConstGenVertexPtr > & vertices() const
Get list of vertices (const)
TGCompositeFrame * fButtonFrame
Button frame.
std::map< std::string, TH1 * > fAnalysisH
Analysis histograms.
TGTextButton * fClearEventCache
Button.
Definition of class GenParticle.
TCanvas * fEventImageCanvas
Event canvas.
void DrawEvent()
Draw evemt.
Definition of class GenVertex.
Definition of class ReaderRootTree.
int event_number() const
Get event number.
std::shared_ptr< Reader > deduce_reader(std::istream &stream)
This function will deduce the type of input stream based on its content and will return appropriate R...
TGTextButton * fNextEvent
Button.
TCanvas * fAnalysisCanvas
Analysis canvas.
HepMC3::GenEvent * fCurrentEvent
Event.
TRootEmbeddedCanvas * fEmbEventImageCanvas
Event canvas.
std::vector< ConstGenParticlePtr > beams() const
Vector of beam particles.
std::shared_ptr< HepMC3::Reader > fReader
Reader.
Stores event-related information.
std::vector< HepMC3::GenEvent * > fEventsCache
Cache of events.
void ClearEventCache()
slot
TGTextButton * fExit
Button.
virtual ~HepMC3ViewerFrame()
Destructor.
TGTextButton * fChooseInput
Button.
static const size_t m_char_buffer_size
Size of writer buffer.
TImage * fGraphImage
Image passed from graphviz.
TGCompositeFrame * fMainFrame
Main frame.
TRootEmbeddedCanvas * fEmbAnalysisCanvas
Analysis canvas.
Definition of class GenEvent.
HepMC3ViewerFrame(const TGWindow *p, UInt_t w, UInt_t h)
Constructor.
void DoAnalysis()
Do analysis.
void ReadFile(const char *a)
Open file.
const std::vector< ConstGenParticlePtr > & particles() const
Get list of particles (const)
TGTextButton * fPreviousEvent
Button.
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...