Fork me on GitHub

source: git/display/DelphesEventDisplay.cc@ 369744d

ImprovedOutputFile Timing dual_readout llp
Last change on this file since 369744d was 369744d, checked in by Christophe Delaere <christophe.delaere@…>, 10 years ago

Legal notice
  • Property mode set to 100644
File size: 18.6 KB
Line 
1/*
2 * Delphes: a framework for fast simulation of a generic collider experiment
3 * Copyright (C) 2012-2014 Universite catholique de Louvain (UCL), Belgium
4 *
5 * This program is free software: you can redistribute it and/or modify
6 * it under the terms of the GNU General Public License as published by
7 * the Free Software Foundation, either version 3 of the License, or
8 * (at your option) any later version.
9 *
10 * This program is distributed in the hope that it will be useful,
11 * but WITHOUT ANY WARRANTY; without even the implied warranty of
12 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
13 * GNU General Public License for more details.
14 *
15 * You should have received a copy of the GNU General Public License
16 * along with this program. If not, see <http://www.gnu.org/licenses/>.
17 */
18
19#include <cassert>
20#include <iostream>
21#include <utility>
22#include <algorithm>
23#include "TGeoManager.h"
24#include "TGeoVolume.h"
25#include "external/ExRootAnalysis/ExRootConfReader.h"
26#include "external/ExRootAnalysis/ExRootTreeReader.h"
27#include "display/DelphesCaloData.h"
28#include "display/DelphesBranchElement.h"
29#include "display/Delphes3DGeometry.h"
30#include "display/DelphesEventDisplay.h"
31#include "classes/DelphesClasses.h"
32#include "TEveElement.h"
33#include "TEveJetCone.h"
34#include "TEveTrack.h"
35#include "TEveTrackPropagator.h"
36#include "TEveCalo.h"
37#include "TEveManager.h"
38#include "TEveGeoNode.h"
39#include "TEveTrans.h"
40#include "TEveViewer.h"
41#include "TEveBrowser.h"
42#include "TEveArrow.h"
43#include "TMath.h"
44#include "TSystem.h"
45#include "TRootBrowser.h"
46#include "TGButton.h"
47#include "TClonesArray.h"
48
49DelphesEventDisplay::DelphesEventDisplay()
50{
51 event_id_ = 0;
52 tkRadius_ = 1.29;
53 totRadius_ = 2.0;
54 tkHalfLength_ = 3.0;
55 bz_ = 3.8;
56 chain_ = new TChain("Delphes");
57 treeReader_ = 0;
58 delphesDisplay_ = 0;
59}
60
61DelphesEventDisplay::~DelphesEventDisplay()
62{
63 delete chain_;
64}
65
66DelphesEventDisplay::DelphesEventDisplay(const char *configFile, const char *inputFile, Delphes3DGeometry& det3D)
67{
68 event_id_ = 0;
69 tkRadius_ = 1.29;
70 totRadius_ = 2.0;
71 tkHalfLength_ = 3.0;
72 bz_ = 3.8;
73 chain_ = new TChain("Delphes");
74 treeReader_ = 0;
75 delphesDisplay_ = 0;
76
77 // initialize the application
78 TEveManager::Create(kTRUE, "IV");
79 TGeoManager* geom = gGeoManager;
80
81 // build the detector
82 tkRadius_ = det3D.getTrackerRadius();
83 totRadius_ = det3D.getDetectorRadius();
84 tkHalfLength_ = det3D.getTrackerHalfLength();
85 bz_ = det3D.getBField();
86 TGeoVolume* top = det3D.getDetector(false);
87 geom->SetTopVolume(top);
88 TEveElementList *geometry = new TEveElementList("Geometry");
89 TObjArray* nodes = top->GetNodes();
90 TIter itNodes(nodes);
91 TGeoNode* nodeobj;
92 TEveGeoTopNode* node;
93 while((nodeobj = (TGeoNode*)itNodes.Next())) {
94 node = new TEveGeoTopNode(gGeoManager,nodeobj);
95 node->UseNodeTrans();
96 geometry->AddElement(node);
97 }
98
99 // Create chain of root trees
100 chain_->Add(inputFile);
101
102 // Create object of class ExRootTreeReader
103 printf("*** Opening Delphes data file ***\n");
104 treeReader_ = new ExRootTreeReader(chain_);
105
106 // prepare data collections
107 readConfig(configFile, det3D, elements_, arrays_);
108 for(std::vector<DelphesBranchBase*>::iterator element = elements_.begin(); element<elements_.end(); ++element) {
109 DelphesBranchElement<TEveTrackList>* item_v1 = dynamic_cast<DelphesBranchElement<TEveTrackList>*>(*element);
110 DelphesBranchElement<TEveElementList>* item_v2 = dynamic_cast<DelphesBranchElement<TEveElementList>*>(*element);
111 if(item_v1) gEve->AddElement(item_v1->GetContainer());
112 if(item_v2) gEve->AddElement(item_v2->GetContainer());
113 }
114
115 // viewers and scenes
116 delphesDisplay_ = new DelphesDisplay;
117 gEve->AddGlobalElement(geometry);
118 delphesDisplay_->ImportGeomRPhi(geometry);
119 delphesDisplay_->ImportGeomRhoZ(geometry);
120 // find the first calo data and use that to initialize the calo display
121 for(std::vector<DelphesBranchBase*>::iterator data=elements_.begin();data<elements_.end();++data) {
122 if(TString((*data)->GetType())=="tower") { // we could also use GetClassName()=="DelphesCaloData"
123 DelphesCaloData* container = dynamic_cast<DelphesBranchElement<DelphesCaloData>*>((*data))->GetContainer();
124 assert(container);
125 TEveCalo3D *calo3d = new TEveCalo3D(container);
126 calo3d->SetBarrelRadius(tkRadius_);
127 calo3d->SetEndCapPos(tkHalfLength_);
128 gEve->AddGlobalElement(calo3d);
129 delphesDisplay_->ImportCaloRPhi(calo3d);
130 delphesDisplay_->ImportCaloRhoZ(calo3d);
131 TEveCaloLego *lego = new TEveCaloLego(container);
132 lego->InitMainTrans();
133 lego->RefMainTrans().SetScale(TMath::TwoPi(), TMath::TwoPi(), TMath::Pi());
134 lego->SetAutoRebin(kFALSE);
135 lego->Set2DMode(TEveCaloLego::kValSizeOutline);
136 delphesDisplay_->ImportCaloLego(lego);
137 break;
138 }
139 }
140
141 make_gui();
142 load_event();
143 gEve->Redraw3D(kTRUE);
144
145}
146
147// function that parses the config to extract the branches of interest and prepare containers
148void DelphesEventDisplay::readConfig(const char *configFile, Delphes3DGeometry& det3D, std::vector<DelphesBranchBase*>& elements, std::vector<TClonesArray*>& arrays) {
149 ExRootConfReader *confReader = new ExRootConfReader;
150 confReader->ReadFile(configFile);
151 Double_t tk_radius = det3D.getTrackerRadius();
152 Double_t tk_length = det3D.getTrackerHalfLength();
153 Double_t tk_Bz = det3D.getBField();
154 Double_t mu_radius = det3D.getDetectorRadius();
155 Double_t mu_length = det3D.getDetectorHalfLength();
156 TAxis* etaAxis = det3D.getCaloAxes().first;
157 TAxis* phiAxis = det3D.getCaloAxes().second;
158 ExRootConfParam branches = confReader->GetParam("TreeWriter::Branch");
159 Int_t nBranches = branches.GetSize()/3;
160 DelphesBranchElement<TEveTrackList>* tlist;
161 DelphesBranchElement<DelphesCaloData>* clist;
162 DelphesBranchElement<TEveElementList>* elist;
163 // first loop with all but tracks
164 for(Int_t b = 0; b<nBranches; ++b) {
165 TString input = branches[b*3].GetString();
166 TString name = branches[b*3+1].GetString();
167 TString className = branches[b*3+2].GetString();
168 if(className=="Tower") {
169 if(input.Contains("eflow",TString::kIgnoreCase) || name.Contains("eflow",TString::kIgnoreCase)) continue; //no eflow
170 clist = new DelphesBranchElement<DelphesCaloData>(name,"tower",kBlack);
171 clist->GetContainer()->SetEtaBins(etaAxis);
172 clist->GetContainer()->SetPhiBins(phiAxis);
173 elements.push_back(clist);
174 } else if(className=="Jet") {
175 if(input.Contains("GenJetFinder")) {
176 elist = new DelphesBranchElement<TEveElementList>(name,"jet",kCyan);
177 elist->GetContainer()->SetRnrSelf(false);
178 elist->GetContainer()->SetRnrChildren(false);
179 elements.push_back(elist);
180 } else {
181 elements.push_back(new DelphesBranchElement<TEveElementList>(name,"jet",kYellow));
182 }
183 } else if(className=="Electron") {
184 tlist = new DelphesBranchElement<TEveTrackList>(name,"electron",kRed);
185 elements.push_back(tlist);
186 TEveTrackPropagator *trkProp = tlist->GetContainer()->GetPropagator();
187 trkProp->SetMagField(0., 0., -tk_Bz);
188 trkProp->SetMaxR(tk_radius);
189 trkProp->SetMaxZ(tk_length);
190 } else if(className=="Photon") {
191 tlist = new DelphesBranchElement<TEveTrackList>(name,"photon",kYellow);
192 elements.push_back(tlist);
193 TEveTrackPropagator *trkProp = tlist->GetContainer()->GetPropagator();
194 trkProp->SetMagField(0., 0., 0.);
195 trkProp->SetMaxR(tk_radius);
196 trkProp->SetMaxZ(tk_length);
197 } else if(className=="Muon") {
198 tlist = new DelphesBranchElement<TEveTrackList>(name,"muon",kGreen);
199 elements.push_back(tlist);
200 TEveTrackPropagator *trkProp = tlist->GetContainer()->GetPropagator();
201 trkProp->SetMagField(0., 0., -tk_Bz);
202 trkProp->SetMaxR(mu_radius);
203 trkProp->SetMaxZ(mu_length);
204 } else if(className=="MissingET") {
205 elements.push_back(new DelphesBranchElement<TEveElementList>(name,"vector",kViolet));
206 } else if(className=="GenParticle") {
207 tlist = new DelphesBranchElement<TEveTrackList>(name,"genparticle",kCyan);
208 elements.push_back(tlist);
209 tlist->GetContainer()->SetRnrSelf(false);
210 tlist->GetContainer()->SetRnrChildren(false);
211 TEveTrackPropagator *trkProp = tlist->GetContainer()->GetPropagator();
212 trkProp->SetMagField(0., 0., -tk_Bz);
213 trkProp->SetMaxR(tk_radius);
214 trkProp->SetMaxZ(tk_length);
215 } else {
216 continue;
217 }
218 arrays.push_back(treeReader_->UseBranch(name));
219 }
220 // second loop for tracks
221 for(Int_t b = 0; b<nBranches; ++b) {
222 TString input = branches[b*3].GetString();
223 TString name = branches[b*3+1].GetString();
224 TString className = branches[b*3+2].GetString();
225 if(className=="Track") {
226 if(input.Contains("eflow",TString::kIgnoreCase) || name.Contains("eflow",TString::kIgnoreCase)) continue; //no eflow
227 tlist = new DelphesBranchElement<TEveTrackList>(name,"track",kBlue);
228 elements.push_back(tlist);
229 TEveTrackPropagator *trkProp = tlist->GetContainer()->GetPropagator();
230 trkProp->SetMagField(0., 0., -tk_Bz);
231 trkProp->SetMaxR(tk_radius);
232 trkProp->SetMaxZ(tk_length);
233 arrays.push_back(treeReader_->UseBranch(name));
234 }
235 }
236}
237
238void DelphesEventDisplay::load_event()
239{
240 // Load event specified in global event_id_.
241 // The contents of previous event are removed.
242
243 // safety
244 if(event_id_ >= treeReader_->GetEntries() || event_id_<0 ) return;
245
246 //TODO move this to the status bar ???
247 printf("Loading event %d.\n", event_id_);
248
249 // clear the previous event
250 gEve->GetViewers()->DeleteAnnotations();
251 for(std::vector<DelphesBranchBase*>::iterator data=elements_.begin();data<elements_.end();++data) {
252 (*data)->Reset();
253 }
254
255 // Load selected branches with data from specified event
256 treeReader_->ReadEntry(event_id_);
257
258 // loop over selected branches, and apply the proper recipe to fill the collections.
259 // this is basically to loop on arrays_ to fill elements_.
260 std::vector<TClonesArray*>::iterator data = arrays_.begin();
261 std::vector<DelphesBranchBase*>::iterator element = elements_.begin();
262 for(; data<arrays_.end() && element<elements_.end(); ++data, ++element) {
263 TString type = (*element)->GetType();
264 // branch on the element type
265 if(type=="tower") delphes_read_towers(*data,*element);
266 else if(type=="track" || type=="photon" || type=="electron" || type=="muon" || type=="genparticle") delphes_read_tracks(*data,*element);
267 else if(type=="jet") delphes_read_jets(*data,*element);
268 else if(type=="vector") delphes_read_vectors(*data,*element);
269 }
270
271 // update display
272 TEveElement* top = (TEveElement*)gEve->GetCurrentEvent();
273 delphesDisplay_->DestroyEventRPhi();
274 delphesDisplay_->ImportEventRPhi(top);
275 delphesDisplay_->DestroyEventRhoZ();
276 delphesDisplay_->ImportEventRhoZ(top);
277 //update_html_summary();
278
279 gEve->Redraw3D(kFALSE, kTRUE);
280}
281
282void DelphesEventDisplay::delphes_read_towers(TClonesArray* data, DelphesBranchBase* element) {
283 DelphesCaloData* container = dynamic_cast<DelphesBranchElement<DelphesCaloData>*>(element)->GetContainer();
284 assert(container);
285 // Loop over all towers
286 TIter itTower(data);
287 Tower *tower;
288 while((tower = (Tower *) itTower.Next()))
289 {
290 container->AddTower(tower->Edges[0], tower->Edges[1], tower->Edges[2], tower->Edges[3]);
291 container->FillSlice(0, tower->Eem);
292 container->FillSlice(1, tower->Ehad);
293 }
294 container->DataChanged();
295}
296
297void DelphesEventDisplay::delphes_read_tracks(TClonesArray* data, DelphesBranchBase* element) {
298 TEveTrackList* container = dynamic_cast<DelphesBranchElement<TEveTrackList>*>(element)->GetContainer();
299 assert(container);
300 TString type = element->GetType();
301 TIter itTrack(data);
302 Int_t counter = 0;
303 TEveTrack *eveTrack;
304 TEveTrackPropagator *trkProp = container->GetPropagator();
305 if(type=="track") {
306 // Loop over all tracks
307 Track *track;
308 while((track = (Track *) itTrack.Next())) {
309 TParticle pb(track->PID, 1, 0, 0, 0, 0,
310 track->P4().Px(), track->P4().Py(),
311 track->P4().Pz(), track->P4().E(),
312 track->X, track->Y, track->Z, 0.0);
313 eveTrack = new TEveTrack(&pb, counter, trkProp);
314 eveTrack->SetName(Form("%s [%d]", pb.GetName(), counter++));
315 eveTrack->SetStdTitle();
316 eveTrack->SetAttLineAttMarker(container);
317 container->AddElement(eveTrack);
318 eveTrack->SetLineColor(element->GetColor());
319 eveTrack->MakeTrack();
320 }
321 } else if(type=="electron") {
322 // Loop over all electrons
323 Electron *electron;
324 while((electron = (Electron *) itTrack.Next())) {
325 TParticle pb(electron->Charge<0?11:-11, 1, 0, 0, 0, 0,
326 electron->P4().Px(), electron->P4().Py(),
327 electron->P4().Pz(), electron->P4().E(),
328 0., 0., 0., 0.);
329
330 eveTrack = new TEveTrack(&pb, counter, trkProp);
331 eveTrack->SetName(Form("%s [%d]", pb.GetName(), counter++));
332 eveTrack->SetStdTitle();
333 eveTrack->SetAttLineAttMarker(container);
334 container->AddElement(eveTrack);
335 eveTrack->SetLineColor(element->GetColor());
336 eveTrack->MakeTrack();
337 }
338 } else if(type=="muon") {
339 // Loop over all muons
340 Muon *muon;
341 while((muon = (Muon *) itTrack.Next())) {
342 TParticle pb(muon->Charge<0?13:-13, 1, 0, 0, 0, 0,
343 muon->P4().Px(), muon->P4().Py(),
344 muon->P4().Pz(), muon->P4().E(),
345 0., 0., 0., 0.);
346
347 eveTrack = new TEveTrack(&pb, counter, trkProp);
348 eveTrack->SetName(Form("%s [%d]", pb.GetName(), counter++));
349 eveTrack->SetStdTitle();
350 eveTrack->SetAttLineAttMarker(container);
351 container->AddElement(eveTrack);
352 eveTrack->SetLineColor(element->GetColor());
353 eveTrack->MakeTrack();
354 }
355 } else if(type=="photon") {
356 // Loop over all photons
357 Photon *photon;
358 while((photon = (Photon *) itTrack.Next())) {
359 TParticle pb(22, 1, 0, 0, 0, 0,
360 photon->P4().Px(), photon->P4().Py(),
361 photon->P4().Pz(), photon->P4().E(),
362 0., 0., 0., 0.);
363 eveTrack = new TEveTrack(&pb, counter, trkProp);
364 eveTrack->SetName(Form("%s [%d]", pb.GetName(), counter++));
365 eveTrack->SetStdTitle();
366 eveTrack->SetAttLineAttMarker(container);
367 eveTrack->SetLineStyle(7);
368 container->AddElement(eveTrack);
369 eveTrack->SetLineColor(element->GetColor());
370 eveTrack->MakeTrack();
371 }
372 } else if(type=="genparticle") {
373 // Loop over all particles
374 GenParticle *particle;
375 while((particle = (GenParticle *) itTrack.Next())) {
376 TParticle pb(particle->PID, particle->Status, particle->M1, particle->M2, particle->D1, particle->D2,
377 particle->P4().Px(), particle->P4().Py(),
378 particle->P4().Pz(), particle->P4().E(),
379 particle->X, particle->Y, particle->Z, particle->T);
380 eveTrack = new TEveTrack(&pb, counter, trkProp);
381 eveTrack->SetName(Form("%s [%d]", pb.GetName(), counter++));
382 eveTrack->SetStdTitle();
383 eveTrack->SetAttLineAttMarker(container);
384 container->AddElement(eveTrack);
385 eveTrack->SetLineColor(element->GetColor());
386 if(particle->Charge==0) eveTrack->SetLineStyle(7);
387 eveTrack->MakeTrack();
388 }
389 }
390}
391
392void DelphesEventDisplay::delphes_read_jets(TClonesArray* data, DelphesBranchBase* element) {
393 TEveElementList* container = dynamic_cast<DelphesBranchElement<TEveElementList>*>(element)->GetContainer();
394 assert(container);
395 TIter itJet(data);
396 Jet *jet;
397 TEveJetCone *eveJetCone;
398 // Loop over all jets
399 Int_t counter = 0;
400 while((jet = (Jet *) itJet.Next()))
401 {
402 eveJetCone = new TEveJetCone();
403 eveJetCone->SetTitle(Form("jet [%d]: Pt=%f, Eta=%f, \nPhi=%f, M=%f",counter,jet->PT, jet->Eta, jet->Phi, jet->Mass));
404 eveJetCone->SetName(Form("jet [%d]", counter++));
405 eveJetCone->SetMainTransparency(60);
406 eveJetCone->SetLineColor(element->GetColor());
407 eveJetCone->SetFillColor(element->GetColor());
408 eveJetCone->SetCylinder(tkRadius_ - 10, tkHalfLength_ - 10);
409 eveJetCone->SetPickable(kTRUE);
410 eveJetCone->AddEllipticCone(jet->Eta, jet->Phi, jet->DeltaEta, jet->DeltaPhi);
411 container->AddElement(eveJetCone);
412 }
413}
414
415void DelphesEventDisplay::delphes_read_vectors(TClonesArray* data, DelphesBranchBase* element) {
416 TEveElementList* container = dynamic_cast<DelphesBranchElement<TEveElementList>*>(element)->GetContainer();
417 assert(container);
418 TIter itMet(data);
419 MissingET *MET;
420 TEveArrow *eveMet;
421 // Missing Et
422 Double_t maxPt = 50.;
423 // TODO to be changed as we don't have access to maxPt anymore. MET scale could be a general parameter set in GUI
424 while((MET = (MissingET*) itMet.Next())) {
425 eveMet = new TEveArrow((tkRadius_ * MET->MET/maxPt)*cos(MET->Phi), (tkRadius_ * MET->MET/maxPt)*sin(MET->Phi), 0., 0., 0., 0.);
426 eveMet->SetMainColor(element->GetColor());
427 eveMet->SetTubeR(0.04);
428 eveMet->SetConeR(0.08);
429 eveMet->SetConeL(0.10);
430 eveMet->SetPickable(kTRUE);
431 eveMet->SetName("Missing Et");
432 eveMet->SetTitle(Form("Missing Et (%.1f GeV)",MET->MET));
433 container->AddElement(eveMet);
434 }
435}
436
437/******************************************************************************/
438// GUI
439/******************************************************************************/
440
441void DelphesEventDisplay::make_gui()
442{
443 // Create minimal GUI for event navigation.
444 // TODO: better GUI could be made based on the ch15 of the manual (Writing a GUI)
445
446 // add a tab on the left
447 TEveBrowser* browser = gEve->GetBrowser();
448 browser->StartEmbedding(TRootBrowser::kLeft);
449
450 // set the main title
451 TGMainFrame* frmMain = new TGMainFrame(gClient->GetRoot(), 1000, 600);
452 frmMain->SetWindowName("Delphes Event Display");
453 frmMain->SetCleanup(kDeepCleanup);
454
455 // build the navigation menu
456 TGHorizontalFrame* hf = new TGHorizontalFrame(frmMain);
457 {
458 TString icondir;
459 if(gSystem->Getenv("ROOTSYS"))
460 icondir = Form("%s/icons/", gSystem->Getenv("ROOTSYS"));
461 if(!gSystem->OpenDirectory(icondir))
462 icondir = Form("%s/icons/", (const char*)gSystem->GetFromPipe("root-config --etcdir") );
463 TGPictureButton* b = 0;
464
465 b = new TGPictureButton(hf, gClient->GetPicture(icondir+"GoBack.gif"));
466 hf->AddFrame(b);
467 b->Connect("Clicked()", "DelphesEventDisplay", this, "Bck()");
468
469 b = new TGPictureButton(hf, gClient->GetPicture(icondir+"GoForward.gif"));
470 hf->AddFrame(b);
471 b->Connect("Clicked()", "DelphesEventDisplay", this, "Fwd()");
472 }
473 frmMain->AddFrame(hf);
474 frmMain->MapSubwindows();
475 frmMain->Resize();
476 frmMain->MapWindow();
477 browser->StopEmbedding();
478 browser->SetTabTitle("Event Control", 0);
479}
480
Note: See TracBrowser for help on using the repository browser.