HepMC3 event record library
HepMC3ViewerFrame.cc
3#include "HepMC3/GenEvent.h"
5#include "HepMC3/GenVertex.h"
6#include "HepMC3ViewerFrame.h"
7
8/* Older graphviz versions can have conflictiong declarations of memcmp/strcmp function
9 * This can break compilation with -pedantic. Uncomenting line below can fix it.
10 */
11// #define _PACKAGE_ast 1
12
13#include <graphviz/gvc.h>
14#define CONSERVATION_TOLERANCE 1e-5
15
16static char* create_image_from_dot(char* m_buffer)
17{
18 GVC_t *gvc = gvContext();
19 Agraph_t *g = agmemread(m_buffer);
20 gvLayout(gvc, g, "dot");
21
22 int err;
23 char *data;
24 unsigned int length;
25
26 if (!g)
27 return nullptr;
28 err = gvRenderData(gvc, g, "png", &data, &length);
29 if (err)
30 return nullptr;
31 data = static_cast<char*>(realloc(data, length + 1));
32 gvFreeLayout(gvc, g);
33 agclose(g);
34 gvFreeContext(gvc);
35 return data;
36}
37
38static bool show_as_parton(HepMC3::ConstGenParticlePtr p )
39{
40 const int pd=std::abs(p->pid());
41 bool parton=false;
42
43 if (pd==81||pd==82||pd<25) parton=true;
44 if (
45 (pd/1000==1||pd/1000==2||pd/1000==3||pd/1000==4||pd/1000==5)
46 &&(pd%1000/100==1||pd%1000/100==2||pd%1000/100==3||pd%1000/100==4)
47 &&(pd%100==1||pd%100==3)
48 )
49 parton=true;
50 if (p->status()==4) parton=true;
51 return parton;
52}
53
54static char* write_event_to_dot(char* used_cursor,const HepMC3::GenEvent &evt,int used_style=1)
55{
56 used_cursor += sprintf(used_cursor, "digraph graphname%d {\n",evt.event_number());
57 used_cursor += sprintf(used_cursor, "v0[label=\"Machine\"];\n");
58 for(auto v: evt.vertices() )
59 {
60 if (used_style!=0)
61 {
62 if (used_style==1) //paint decay and fragmentation vertices in green
63 {
64 if (v->status()==2) used_cursor += sprintf(used_cursor, "node [color=\"green\"];\n");
65 else used_cursor += sprintf(used_cursor, "node [color=\"black\"];\n");
66 }
67 }
68 HepMC3::FourVector in(0,0,0,0);
69 HepMC3::FourVector out(0,0,0,0);
70 double energy=0;
71 for(auto p1: v->particles_in() ) {
72 in+=p1->momentum();
73 energy+=std::abs(p1->momentum().e());
74 }
75 for(auto p2: v->particles_out() ) {
76 out+=p2->momentum();
77 energy+=std::abs(p2->momentum().e());
78 }
79 HepMC3::FourVector momviolation(0,0,0,0);
80 momviolation+=in;
81 momviolation-=out;
82 double energyviolation=std::sqrt(momviolation.length2() +momviolation.e()*momviolation.e() );
83 bool violation=false;
84 if (energyviolation>CONSERVATION_TOLERANCE*energy) violation=true;
85
86 if(violation)
87 {
88 used_cursor += sprintf(used_cursor, "node [shape=rectangle];\n");
89 used_cursor += sprintf(used_cursor, "v%d [label=\"%d\nd=%4.2f\"];\n", -v->id(),v->id(),energyviolation);
90 }
91 else
92 {
93 used_cursor += sprintf(used_cursor, "node [shape=ellipse];\n");
94 used_cursor += sprintf(used_cursor, "v%d[label=\"%d\"];\n", -v->id(),v->id());
95 }
96
97 used_cursor += sprintf(used_cursor, "node [shape=ellipse];\n");
98 }
99 for(auto p: evt.beams() )
100 {
101 if (!p->end_vertex()) continue;
102 used_cursor += sprintf(used_cursor, "node [shape=point];\n");
103 used_cursor += sprintf(used_cursor, "v0 -> v%d [label=\"%d(%d)\"];\n", -p->end_vertex()->id(),p->id(),p->pid());
104 }
105
106 for(auto v: evt.vertices() )
107 {
108
109 for(auto p: v->particles_out() )
110 {
111 {
112 if (used_style!=0)
113 {
114 if (used_style==1) //paint suspected partons and 81/82 in red
115 {
116 if (show_as_parton(p)&&p->status()!=1) used_cursor += sprintf(used_cursor, "edge [color=\"red\"];\n");
117 else used_cursor +=sprintf(used_cursor, "edge [color=\"black\"];\n");
118 }
119 }
120 if (!p->end_vertex())
121 {
122 used_cursor += sprintf(used_cursor, "node [shape=point];\n");
123 used_cursor += sprintf(used_cursor, "v%d -> o%d [label=\"%d(%d)\"];\n", -v->id(),p->id(),p->id(),p->pid());
124 continue;
125 }
126 else
127 used_cursor += sprintf(used_cursor, "v%d -> v%d [label=\"%d(%d)\"];\n", -v->id(),-p->end_vertex()->id(),p->id(),p->pid());
128 }
129 }
130 }
131 used_cursor += sprintf(used_cursor, "labelloc=\"t\";\nlabel=\"Event %d; Vertices %lu; Particles %lu;\";\n", evt.event_number(), evt.vertices().size(), evt.particles().size());
132 used_cursor += sprintf(used_cursor,"}\n\n");
133
134 return used_cursor;
135}
136
137
139{
140 char* m_buffer = new char[m_char_buffer_size]();
141 char* m_cursor=m_buffer;
142 m_cursor=write_event_to_dot(m_cursor,*(fCurrentEvent));
143 char *buf=create_image_from_dot(m_buffer);
144 fEmbEventImageCanvas->MapSubwindows();
145
146 if(!fEventImageCanvas) fEventImageCanvas=new TCanvas("fEmbEventImageCanvas","fEmbEventImageCanvas",1024,768);
147
148 fEventImageCanvas->cd();
149 fEventImageCanvas->Clear();
150 double d=0.60;
151
152 fGraphImage = TImage::Create();
153 fGraphImage->SetName("Event");
154 fGraphImage->SetImageBuffer(&buf, TImage::kPng);
155
156 fGraphImage->SetConstRatio(kFALSE);
157
158 TPad *p1 = new TPad("i1", "i1", 0.05, 0.05, 0.05+d*fGraphImage->GetWidth()/fGraphImage->GetHeight(), 0.95);
159 p1->Draw();
160 p1->cd();
161
162 fGraphImage->Draw("xxx");
163 delete [] m_buffer;
164 gPad->Update();
165 DoAnalysis();
166}
167
169{
170 fEmbAnalysisCanvas->MapSubwindows();
171 fAnalysisCanvas->cd();
172 fAnalysisCanvas->Clear();
173 for (auto h: fAnalysisH) h.second->Delete();
174 fAnalysisH.clear();
175
176 /* */
177 TH1S* particles1= new TH1S();
178 fAnalysisH["particles1"]=particles1;
179 particles1->SetTitle("Flavour: all particles; PDG ID; Number of particles");
180 particles1->SetFillColor(kBlue);
181 for(auto p: fCurrentEvent->particles() )
182 particles1->Fill((std::to_string(p->pid())).c_str(),1.0);
183 particles1->LabelsOption(">","X");
184 /* */
185 TH1S* particles2= new TH1S();
186 fAnalysisH["particles2"]=particles2;
187 particles2->SetTitle("Flavour: particles with status 1; PDG ID; Number of particles");
188 particles2->SetFillColor(kBlue);
189 for(auto p: fCurrentEvent->particles() )
190 if(p->status()==1) particles2->Fill((std::to_string(p->pid())).c_str(),1.0);
191 particles2->LabelsOption(">","X");
192 /* */
193 std::vector<double> masses;
194 for(auto p: fCurrentEvent->particles() )
195 if(show_as_parton(p)) masses.push_back(p->momentum().m());
196 TH1D* particles3= new TH1D("particles3","Mass: parton particles; Mass, GeV; Number of particles",masses.size(),
197 0,*std::max_element(masses.begin(),masses.end()));
198 fAnalysisH["particles3"]=particles3;
199 particles3->SetFillColor(kBlue);
200 for(auto m: masses) particles3->Fill(m);
201
202
203 fAnalysisCanvas->cd();
204 TPad *p1 = new TPad("i1", "i1", 0.00, 0.75, 1.0, 1.0);
205 p1->Draw();
206 p1->cd();
207 particles1->Draw();
208 fAnalysisCanvas->cd();
209 TPad *p2 = new TPad("i2", "i2", 0.00, 0.50, 1.0, 0.75);
210 p2->Draw();
211 p2->cd();
212 particles2->Draw();
213 fAnalysisCanvas->cd();
214 TPad *p3 = new TPad("i3", "i3", 0.00, 0.25, 1.0, 0.50);
215 p3->Draw();
216 p3->cd();
217 particles3->Draw();
218
219 gPad->Update();
220}
221
223{
224 fEventsCache.clear();
225 fCurrentEvent=nullptr;
226}
227
229{
230 auto pos=find(fEventsCache.begin(),fEventsCache.end(),fCurrentEvent);
231 if (pos==fEventsCache.begin()) return;
232 pos--;
233 fCurrentEvent=*(pos);
234 if (pos==fEventsCache.end()) printf("This event was not found in the cache. Cache size is %zu\n",fEventsCache.size());
235 DrawEvent();
236}
237
241
243{
244 if (fCurrentEvent==nullptr||fEventsCache.back()==fCurrentEvent)
245 {
246 HepMC3::GenEvent* evt1=new HepMC3::GenEvent(HepMC3::Units::GEV,HepMC3::Units::MM);
247 bool ok=fReader->read_event(*(evt1));
248 ok=(ok&&!fReader->failed());
249 if (ok)
250 {
251 fEventsCache.push_back(evt1);
252 fCurrentEvent=evt1;
253 }
254 else return;
255 }
256 else
257 {
258 auto pos=find(fEventsCache.begin(),fEventsCache.end(),fCurrentEvent);
259 pos++;
260 fCurrentEvent=*(pos);
261 }
262 DrawEvent();
263}
265{
266 static const char *FileType[] = {"All", "*.*","HepMC", "*.hepmc*","LHEF", "*.lhe*","ROOT", "*.root", 0, 0 };
267 static TString dir("./");
268 TGFileInfo fi;
269 fi.fFileTypes = FileType;
270 fi.fIniDir = StrDup(dir);
271 new TGFileDialog(gClient->GetRoot(), this, kFDOpen, &fi);
272 if (fReader) fReader->close();
273 fReader=HepMC3::deduce_reader(fi.fFilename);
274}
275
276HepMC3ViewerFrame::HepMC3ViewerFrame(const TGWindow *p, UInt_t w, UInt_t h) :
277 TGMainFrame(p, w, h)
278{
279 fMainFrame = new TGCompositeFrame(this, 1350, 500, kHorizontalFrame|kFixedWidth);
280 fButtonFrame = new TGCompositeFrame(fMainFrame, 150, 200, kFixedWidth);
281
282 fEmbEventImageCanvas =new TRootEmbeddedCanvas("MainCanvaslegent", fMainFrame, 850, 500);
283
284 fEmbAnalysisCanvas =new TRootEmbeddedCanvas("EmbAnalysisCanvaslegend", fMainFrame, 350, 500);
285
286
287 fMainFrame->AddFrame(fEmbEventImageCanvas,new TGLayoutHints(kLHintsTop | kLHintsExpandX| kLHintsExpandY, 1, 1, 2, 2));
288 fMainFrame->AddFrame(fEmbAnalysisCanvas,new TGLayoutHints(kLHintsTop | kFixedWidth| kLHintsExpandY, 1, 1, 2, 2));
289 fMainFrame->AddFrame(fButtonFrame,new TGLayoutHints(kLHintsTop, 1, 1, 2, 2));
290
291
292 fChooseInput = new TGTextButton(fButtonFrame, "&Choose input");
293 fChooseInput->Connect("Clicked()", "HepMC3ViewerFrame", this, "ChooseInput()");
294 fChooseInput->SetToolTipText("Click to choose file");
295 fButtonFrame->AddFrame(fChooseInput, new TGLayoutHints(kLHintsTop | kLHintsExpandX, 1, 1, 2, 2));
296
297
298
299 fNextEvent = new TGTextButton(fButtonFrame, "&Next event");
300 fNextEvent->Connect("Clicked()", "HepMC3ViewerFrame", this, "NextEvent()");
301 fNextEvent->SetToolTipText("Click to display next event");
302 fButtonFrame->AddFrame(fNextEvent, new TGLayoutHints(kLHintsExpandX|kLHintsLeft, 1, 1, 2, 2));
303
304
305 fPreviousEvent = new TGTextButton(fButtonFrame, "&Previous event");
306 fPreviousEvent->Connect("Clicked()", "HepMC3ViewerFrame", this, "PreviousEvent()");
307 fPreviousEvent->SetToolTipText("Click to display previous event");
308 fButtonFrame->AddFrame(fPreviousEvent, new TGLayoutHints( kLHintsExpandX|kLHintsLeft, 1, 1, 2, 2));
309
310
311 fClearEventCache = new TGTextButton(fButtonFrame, "&Clear event cache");
312 fClearEventCache->Connect("Clicked()", "HepMC3ViewerFrame", this, "ClearEventCache()");
313 fClearEventCache->SetToolTipText("Click to clear event cache ");
314 fButtonFrame->AddFrame(fClearEventCache, new TGLayoutHints( kLHintsExpandX|kLHintsLeft, 1, 1, 2, 2));
315
316 fExit = new TGTextButton(fButtonFrame, "&Exit ","gApplication->Terminate(0)");
317 fExit->SetToolTipText("Click to exit");
318 fButtonFrame->AddFrame(fExit, new TGLayoutHints( kLHintsExpandX|kLHintsLeft,1,1,2,2));
319
320 AddFrame(fMainFrame, new TGLayoutHints(kLHintsTop |kLHintsExpandX| kLHintsExpandY, 1, 1, 2, 2));
321
322 SetWindowName("Event viewer");
323 MapSubwindows();
324 Resize(GetDefaultSize());
325 MapWindow();
326
327 fReader=nullptr;
330 fCurrentEvent=nullptr;
331 fGraphImage = TImage::Create();
332}
334{
335 fMainFrame->Cleanup();
336 fReader->close();
337 Cleanup();
338}
Definition of class GenEvent.
Definition of class GenParticle.
Definition of class GenVertex.
Declaration of deduce_reader and related functions.
Definition of class ReaderRootTree.
void DrawEvent()
Draw evemt.
virtual ~HepMC3ViewerFrame()
Destructor.
TCanvas * fAnalysisCanvas
Analysis canvas.
void ReadFile(const char *a)
Open file.
std::map< std::string, TH1 * > fAnalysisH
Analysis histograms.
HepMC3ViewerFrame(const TGWindow *p, UInt_t w, UInt_t h)
Constructor.
static const size_t m_char_buffer_size
Size of writer buffer.
TGCompositeFrame * fButtonFrame
Button frame.
std::shared_ptr< HepMC3::Reader > fReader
Reader.
TRootEmbeddedCanvas * fEmbAnalysisCanvas
Analysis canvas.
TCanvas * fEventImageCanvas
Event canvas.
TGCompositeFrame * fMainFrame
Main frame.
TGTextButton * fChooseInput
Button.
std::vector< HepMC3::GenEvent * > fEventsCache
Cache of events.
TGTextButton * fExit
Button.
TImage * fGraphImage
Image passed from graphviz.
TGTextButton * fPreviousEvent
Button.
TGTextButton * fClearEventCache
Button.
TGTextButton * fNextEvent
Button.
HepMC3::GenEvent * fCurrentEvent
Event.
void DoAnalysis()
Do analysis.
TRootEmbeddedCanvas * fEmbEventImageCanvas
Event canvas.
Generic 4-vector.
Definition FourVector.h:36
Stores event-related information.
Definition GenEvent.h:47
int event_number() const
Get event number.
Definition GenEvent.h:155
std::vector< ConstGenParticlePtr > beams() const
Vector of beam particles.
Definition GenEvent.cc:429
const std::vector< ConstGenVertexPtr > & vertices() const
Get list of vertices (const)
Definition GenEvent.cc:43
const std::vector< ConstGenParticlePtr > & particles() const
Get list of particles (const)
Definition GenEvent.cc:39
std::shared_ptr< Reader > deduce_reader(const std::string &filename)
This function deduces the type of input file based on the name/URL and its content,...