HepMC3 event record library
HepMC3ViewerFrame.cc
2 #include "HepMC3/ReaderFactory.h"
3 #include "HepMC3ViewerFrame.h"
4 
5 
6 
7 #include "HepMC3/GenEvent.h"
8 #include "HepMC3/GenParticle.h"
9 #include "HepMC3/GenVertex.h"
10 
11 #include <graphviz/gvc.h>
12 #define CONSERVATION_TOLERANCE 1e-5
13 
14 static char* create_image_from_dot(char* m_buffer)
15 {
16  GVC_t * gvc=gvContext();
17  Agraph_t *g= agmemread(m_buffer);
18  gvLayout(gvc,g,"dot");
19 
20  int err;
21  char *data;
22  unsigned int length;
23 
24  if (!g)
25  return NULL;
26  err = gvRenderData(gvc, g, "png", &data, &length);
27  if (err)
28  return NULL;
29  data = (char*)realloc(data, length + 1);
30  delete g;
31  delete gvc;
32  return data;
33 }
34 
35 static bool show_as_parton(HepMC3::ConstGenParticlePtr p )
36 {
37  const int pd=std::abs(p->pid());
38  bool parton=false;
39 
40  if (pd==81||pd==82||pd<25) parton=true;
41  if (
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)
45  )
46  parton=true;
47  if (p->status()==4) parton=true;
48  return parton;
49 }
50 
51 static char* write_event_to_dot(char* used_cursor,const HepMC3::GenEvent &evt,int used_style=1)
52 {
53  used_cursor += sprintf(used_cursor, "digraph graphname%d {\n",evt.event_number());
54  used_cursor += sprintf(used_cursor, "v0[label=\"Machine\"];\n");
55  for(auto v: evt.vertices() )
56  {
57  if (used_style!=0)
58  {
59  if (used_style==1) //paint decay and fragmentation vertices in green
60  {
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");
63  }
64  }
67  double energy=0;
68  for(auto p1: v->particles_in() ) {
69  in+=p1->momentum();
70  energy+=std::abs(p1->momentum().e());
71  }
72  for(auto p2: v->particles_out() ) {
73  out+=p2->momentum();
74  energy+=std::abs(p2->momentum().e());
75  }
76  HepMC3::FourVector momviolation(0,0,0,0);
77  momviolation+=in;
78  momviolation-=out;
79  double energyviolation=std::sqrt(momviolation.length2() +momviolation.e()*momviolation.e() );
80  bool violation=false;
81  if (energyviolation>CONSERVATION_TOLERANCE*energy) violation=true;
82 
83  if(violation)
84  {
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);
87  }
88  else
89  {
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());
92  }
93 
94  used_cursor += sprintf(used_cursor, "node [shape=ellipse];\n");
95  }
96  for(auto p: evt.beams() )
97  {
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());
101  }
102 
103  for(auto v: evt.vertices() )
104  {
105 
106  for(auto p: v->particles_out() )
107  {
108  {
109  if (used_style!=0)
110  {
111  if (used_style==1) //paint suspected partons and 81/82 in red
112  {
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");
115  }
116  }
117  if (!p->end_vertex())
118  {
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());
121  continue;
122  }
123  else
124  used_cursor += sprintf(used_cursor, "v%d -> v%d [label=\"%d(%d)\"];\n", -v->id(),-p->end_vertex()->id(),p->id(),p->pid());
125  }
126  }
127  }
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");
130 
131  return used_cursor;
132 }
133 
134 
136 {
137  char* m_buffer = new char[m_char_buffer_size]();
138  char* m_cursor=m_buffer;
139  m_cursor=write_event_to_dot(m_cursor,*(fCurrentEvent));
140  char *buf=create_image_from_dot(m_buffer);
141  fEmbEventImageCanvas->MapSubwindows();
142 
143  if(!fEventImageCanvas) fEventImageCanvas=new TCanvas("fEmbEventImageCanvas","fEmbEventImageCanvas",1024,768);
144 
145  fEventImageCanvas->cd();
146  fEventImageCanvas->Clear();
147  double d=0.60;
148 
149  fGraphImage = TImage::Create();
150  fGraphImage->SetName("Event");
151  fGraphImage->SetImageBuffer(&buf, TImage::kPng);
152 
153  fGraphImage->SetConstRatio(kFALSE);
154 
155  TPad *p1 = new TPad("i1", "i1", 0.05, 0.05, 0.05+d*fGraphImage->GetWidth()/fGraphImage->GetHeight(), 0.95);
156  p1->Draw();
157  p1->cd();
158 
159  fGraphImage->Draw("xxx");
160  delete [] m_buffer;
161  gPad->Update();
162  DoAnalysis();
163 }
164 
166 {
167  fEmbAnalysisCanvas->MapSubwindows();
168  fAnalysisCanvas->cd();
169  fAnalysisCanvas->Clear();
170  for (auto h: fAnalysisH) h.second->Delete();
171  fAnalysisH.clear();
172 
173  /* */
174  TH1S* particles1= new TH1S();
175  fAnalysisH["particles1"]=particles1;
176  particles1->SetTitle("Flavour: all particles; PDG ID; Number of particles");
177  particles1->SetFillColor(kBlue);
178  for(auto p: fCurrentEvent->particles() )
179  particles1->Fill((std::to_string(p->pid())).c_str(),1.0);
180  particles1->LabelsOption(">","X");
181  /* */
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);
186  for(auto p: fCurrentEvent->particles() )
187  if(p->status()==1) particles2->Fill((std::to_string(p->pid())).c_str(),1.0);
188  particles2->LabelsOption(">","X");
189  /* */
190  std::vector<double> masses;
191  for(auto p: fCurrentEvent->particles() )
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);
198 
199 
200  fAnalysisCanvas->cd();
201  TPad *p1 = new TPad("i1", "i1", 0.00, 0.75, 1.0, 1.0);
202  p1->Draw();
203  p1->cd();
204  particles1->Draw();
205  fAnalysisCanvas->cd();
206  TPad *p2 = new TPad("i2", "i2", 0.00, 0.50, 1.0, 0.75);
207  p2->Draw();
208  p2->cd();
209  particles2->Draw();
210  fAnalysisCanvas->cd();
211  TPad *p3 = new TPad("i3", "i3", 0.00, 0.25, 1.0, 0.50);
212  p3->Draw();
213  p3->cd();
214  particles3->Draw();
215 
216  gPad->Update();
217 }
218 
220 {
221  fEventsCache.clear();
222  fCurrentEvent=nullptr;
223 }
224 
226 {
227  auto pos=find(fEventsCache.begin(),fEventsCache.end(),fCurrentEvent);
228  if (pos==fEventsCache.begin()) return;
229  pos--;
230  fCurrentEvent=*(pos);
231  if (pos==fEventsCache.end()) printf("This event was not found in the cache. Cache size is %zu\n",fEventsCache.size());
232  DrawEvent();
233 }
234 
235 void HepMC3ViewerFrame::ReadFile(const char* a) {
237 }
238 
240 {
241  if (fCurrentEvent==nullptr||fEventsCache.back()==fCurrentEvent)
242  {
243  HepMC3::GenEvent* evt1=new HepMC3::GenEvent(HepMC3::Units::GEV,HepMC3::Units::MM);
244  bool ok=fReader->read_event(*(evt1));
245  ok=(ok&&!fReader->failed());
246  if (ok)
247  {
248  fEventsCache.push_back(evt1);
249  fCurrentEvent=evt1;
250  }
251  else return;
252  }
253  else
254  {
255  auto pos=find(fEventsCache.begin(),fEventsCache.end(),fCurrentEvent);
256  pos++;
257  fCurrentEvent=*(pos);
258  }
259  DrawEvent();
260 }
262 {
263  static const char *FileType[] = {"All", "*.*","HepMC", "*.hepmc*","LHEF", "*.lhe*","ROOT", "*.root", 0, 0 };
264  static TString dir("./");
265  TGFileInfo fi;
266  fi.fFileTypes = FileType;
267  fi.fIniDir = StrDup(dir);
268  new TGFileDialog(gClient->GetRoot(), this, kFDOpen, &fi);
269  if (fReader) fReader->close();
270  fReader=HepMC3::deduce_reader(fi.fFilename);
271 }
272 
273 HepMC3ViewerFrame::HepMC3ViewerFrame(const TGWindow *p, UInt_t w, UInt_t h) :
274  TGMainFrame(p, w, h)
275 {
276  fMainFrame = new TGCompositeFrame(this, 1350, 500, kHorizontalFrame|kFixedWidth);
277  fButtonFrame = new TGCompositeFrame(fMainFrame, 150, 200, kFixedWidth);
278 
279  fEmbEventImageCanvas =new TRootEmbeddedCanvas("MainCanvaslegent", fMainFrame, 850, 500);
280 
281  fEmbAnalysisCanvas =new TRootEmbeddedCanvas("EmbAnalysisCanvaslegend", fMainFrame, 350, 500);
282 
283 
284  fMainFrame->AddFrame(fEmbEventImageCanvas,new TGLayoutHints(kLHintsTop | kLHintsExpandX| kLHintsExpandY, 1, 1, 2, 2));
285  fMainFrame->AddFrame(fEmbAnalysisCanvas,new TGLayoutHints(kLHintsTop | kFixedWidth| kLHintsExpandY, 1, 1, 2, 2));
286  fMainFrame->AddFrame(fButtonFrame,new TGLayoutHints(kLHintsTop, 1, 1, 2, 2));
287 
288 
289  fChooseInput = new TGTextButton(fButtonFrame, "&Choose input");
290  fChooseInput->Connect("Clicked()", "HepMC3ViewerFrame", this, "ChooseInput()");
291  fChooseInput->SetToolTipText("Click to choose file");
292  fButtonFrame->AddFrame(fChooseInput, new TGLayoutHints(kLHintsTop | kLHintsExpandX, 1, 1, 2, 2));
293 
294 
295 
296  fNextEvent = new TGTextButton(fButtonFrame, "&Next event");
297  fNextEvent->Connect("Clicked()", "HepMC3ViewerFrame", this, "NextEvent()");
298  fNextEvent->SetToolTipText("Click to display next event");
299  fButtonFrame->AddFrame(fNextEvent, new TGLayoutHints(kLHintsExpandX|kLHintsLeft, 1, 1, 2, 2));
300 
301 
302  fPreviousEvent = new TGTextButton(fButtonFrame, "&Previous event");
303  fPreviousEvent->Connect("Clicked()", "HepMC3ViewerFrame", this, "PreviousEvent()");
304  fPreviousEvent->SetToolTipText("Click to display previous event");
305  fButtonFrame->AddFrame(fPreviousEvent, new TGLayoutHints( kLHintsExpandX|kLHintsLeft, 1, 1, 2, 2));
306 
307 
308  fClearEventCache = new TGTextButton(fButtonFrame, "&Clear event cache");
309  fClearEventCache->Connect("Clicked()", "HepMC3ViewerFrame", this, "ClearEventCache()");
310  fClearEventCache->SetToolTipText("Click to clear event cache ");
311  fButtonFrame->AddFrame(fClearEventCache, new TGLayoutHints( kLHintsExpandX|kLHintsLeft, 1, 1, 2, 2));
312 
313  fExit = new TGTextButton(fButtonFrame, "&Exit ","gApplication->Terminate(0)");
314  fExit->SetToolTipText("Click to exit");
315  fButtonFrame->AddFrame(fExit, new TGLayoutHints( kLHintsExpandX|kLHintsLeft,1,1,2,2));
316 
317  AddFrame(fMainFrame, new TGLayoutHints(kLHintsTop |kLHintsExpandX| kLHintsExpandY, 1, 1, 2, 2));
318 
319  SetWindowName("Event viewer");
320  MapSubwindows();
321  Resize(GetDefaultSize());
322  MapWindow();
323 
324  fReader=nullptr;
326  fAnalysisCanvas=fEmbAnalysisCanvas->GetCanvas();
327  fCurrentEvent=nullptr;
328  fGraphImage = TImage::Create();
329 }
331 {
332  fMainFrame->Cleanup();
333  fReader->close();
334  Cleanup();
335 }
const std::vector< ConstGenVertexPtr > & vertices() const
Get list of vertices (const)
Definition: GenEvent.cc:43
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.
Definition: GenEvent.h:135
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.
Definition: GenEvent.cc:420
std::shared_ptr< HepMC3::Reader > fReader
Reader.
Stores event-related information.
Definition: GenEvent.h:41
Generic 4-vector.
Definition: FourVector.h:35
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)
Definition: GenEvent.cc:39
TGTextButton * fPreviousEvent
Button.
Feature< Feature_type > abs(const Feature< Feature_type > &input)
Obtain the absolute value of a Feature. This works as you&#39;d expect. If foo is a valid Feature...
Definition: Feature.h:316