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