2 // Author: Stefano Carrazza 2010, CERN, stefano.carrazza@cern.ch
4 /**************************************************************************
5 * Copyright(c) 1998-2009, ALICE Experiment at CERN, all rights reserved. *
6 * See http://aliceinfo.cern.ch/Offline/AliRoot/License.html for *
7 * full copyright notice. *
8 **************************************************************************/
10 #include "AliESDEvent.h"
11 #include "AliEveLego.h"
12 #include "AliEveEventManager.h"
13 #include "AliEveMultiView.h"
14 #include "AliPhysicsSelection.h"
18 #include "TGLViewer.h"
19 #include "TEveWindow.h"
20 #include "TEveManager.h"
21 #include "TEveBrowser.h"
22 #include "TEveViewer.h"
23 #include "TEveScene.h"
24 #include "TEveCaloLegoOverlay.h"
26 #include "TEveCaloData.h"
27 #include "TEveLegoEventHandler.h"
28 #include "TEveTrans.h"
29 #include "TEveProjectionManager.h"
30 #include "TEveProjectionAxes.h"
31 #include "TGLWidget.h"
32 #include "TStopwatch.h"
34 //______________________________________________________________________________
35 // This class provides the following features:
37 // 1) 2D and 3D calorimeter like histograms of tracks distribution,
38 // providing particle selection by charge (positive and negative),
39 // and by specie: electrons, muons, pions, kaons and protons.
41 // 2) Histograms are plotted around detectors, allowing track selection: all tracks,
42 // only primary tracks. It allows pT maximum and minimum threashold.
44 // 3) 2D and 3D calorimeter like histograms of particles distribution obtained
47 // 4) Possibility to use AliPhysicsSelection during the all events histograms creation.
48 // It is also possible to switch between real data and simulation data (MC).
52 Double_t kPi = TMath::Pi();
54 //______________________________________________________________________________
55 AliEveLego::AliEveLego(const char* name) :
56 TEveElementList(name),
58 fCollisionCandidatesOnly(kFALSE),
68 fHistoposAllEvents(0),
70 fHistonegAllEvents(0),
72 fHistoElectronsAllEvents(0),
74 fHistoMuonsAllEvents(0),
76 fHistoPionsAllEvents(0),
78 fHistoKaonsAllEvents(0),
80 fHistoProtonsAllEvents(0),
91 fHisto2dAllEventsv0(0),
92 fHisto2dAllEventsv1(0),
93 fHisto2dAllEventsv2(0),
94 fHisto2dAllEventsv3(0),
95 fHisto2dAllEventss0(0),
96 fHisto2dAllEventss1(0),
97 fHisto2dAllEventss2(0),
98 fHisto2dAllEventss3(0),
100 fHisto2dLegoOverlay(0),
101 fHisto2dAllEventsLegoOverlay(0),
102 fHisto2dAllEventsSlot(0)
105 gEve->AddToListTree(this,0);
107 // Get Current ESD event
108 fEsd = AliEveEventManager::AssertESD();
110 // Particles types per default
111 fParticleTypeId = new Bool_t[7];
112 fParticleTypeIdAE = new Bool_t[7];
114 for (Int_t s = 0; s < 7; s++)
118 fParticleTypeId[s] = kFALSE;
119 fParticleTypeIdAE[s] = kFALSE;
121 fParticleTypeId[s] = kTRUE;
122 fParticleTypeIdAE[s] = kTRUE;
126 // Loading Physics Selection to determine the collision candidates
127 fPhysicsSelection = new AliPhysicsSelection();
128 fPhysicsSelection->Initialize(fEsd);
130 fHistopos = new TH2F("histopos","Histo 2d positive", 100, -1.5, 1.5, 80, -kPi, kPi);
131 fHistoneg = new TH2F("histoneg","Histo 2d negative", 100, -1.5, 1.5, 80, -kPi, kPi);
132 fHistoElectrons = new TH2F("histoele","Histo 2d electron", 100, -1.5, 1.5, 80, -kPi, kPi);
133 fHistoMuons = new TH2F("histomuo","Histo 2d muons ", 100, -1.5, 1.5, 80, -kPi, kPi);
134 fHistoPions = new TH2F("histopio","Histo 2d pions ", 100, -1.5, 1.5, 80, -kPi, kPi);
135 fHistoKaons = new TH2F("histokao","Histo 2d kaons ", 100, -1.5, 1.5, 80, -kPi, kPi);
136 fHistoProtons = new TH2F("histopro","Histo 2d protons ", 100, -1.5, 1.5, 80, -kPi, kPi);
138 fHistopos->SetDirectory(0);
139 fHistoneg->SetDirectory(0);
140 fHistoElectrons->SetDirectory(0);
141 fHistoMuons->SetDirectory(0);
142 fHistoPions->SetDirectory(0);
143 fHistoKaons->SetDirectory(0);
144 fHistoProtons->SetDirectory(0);
146 // Colors from get_pdg_color() in /alice-macros/kine_tracks.C
147 static Color_t fElectro = 5;
148 static Color_t fMuon = 30;
149 static Color_t fPion = 3;
150 static Color_t fKaon = 38;
151 static Color_t fProton = 10;
153 // Adding data to TEveCaloDataHist
154 fData = new TEveCaloDataHist();
155 fData->AddHistogram(fHistoneg);
156 fData->RefSliceInfo(0).Setup("NegCg:", 0, kBlue);
157 fData->AddHistogram(fHistopos);
158 fData->RefSliceInfo(1).Setup("PosCg:", 0, kRed);
159 fData->AddHistogram(fHistoElectrons);
160 fData->RefSliceInfo(2).Setup("Elect.:", 0, fElectro);
161 fData->AddHistogram(fHistoMuons);
162 fData->RefSliceInfo(3).Setup("Muons:", 0, fMuon);
163 fData->AddHistogram(fHistoPions);
164 fData->RefSliceInfo(4).Setup("Pions:", 0, fPion);
165 fData->AddHistogram(fHistoKaons);
166 fData->RefSliceInfo(5).Setup("Kaons:", 0, fKaon);
167 fData->AddHistogram(fHistoProtons);
168 fData->RefSliceInfo(6).Setup("Proto.:", 0, fProton);
170 fData->GetEtaBins()->SetTitleFont(120);
171 fData->GetEtaBins()->SetTitle("h");
172 fData->GetPhiBins()->SetTitleFont(120);
173 fData->GetPhiBins()->SetTitle("f");
174 fData->IncDenyDestroy();
176 // Setting up the position of the 3D calorimeter histogram view.
177 fCalo3d = new TEveCalo3D(fData);
178 fCalo3d->SetBarrelRadius(550);
179 fCalo3d->SetEndCapPos(550);
181 // Adding data to the Lego object
182 fLego = new TEveCaloLego(fData);
184 // Creating projections
185 fAl = AliEveMultiView::Instance();
186 fAl->ImportEventRPhi(fCalo3d);
187 fAl->ImportEventRhoZ(fCalo3d);
189 // Update viewers and scenes
193 //______________________________________________________________________________
194 AliEveLego::~AliEveLego()
196 // Deleting variables
198 delete fPhysicsSelection;
200 delete fHistoposAllEvents;
202 delete fHistonegAllEvents;
203 delete fHistoElectrons;
204 delete fHistoElectronsAllEvents;
206 delete fHistoMuonsAllEvents;
208 delete fHistoPionsAllEvents;
210 delete fHistoKaonsAllEvents;
211 delete fHistoProtons;
212 delete fHistoProtonsAllEvents;
213 delete[] fParticleTypeId;
214 delete[] fParticleTypeIdAE;
216 delete fDataAllEvents;
218 delete fLegoAllEvents;
220 delete fCalo3dAllEvents;
225 delete fHisto2dAllEventsv0;
226 delete fHisto2dAllEventsv1;
227 delete fHisto2dAllEventsv2;
228 delete fHisto2dAllEventsv3;
229 delete fHisto2dAllEventss0;
230 delete fHisto2dAllEventss1;
231 delete fHisto2dAllEventss2;
232 delete fHisto2dAllEventss3;
234 delete fHisto2dLegoOverlay;
235 delete fHisto2dAllEventsLegoOverlay;
236 delete fHisto2dAllEventsSlot;
241 //____________________________________________________________________________
242 Double_t getphi(Double_t phi)
244 // Phi correction for alice
246 if (phi > TMath::Pi()) {
247 phi -= TMath::TwoPi();
253 //______________________________________________________________________________
254 TEveCaloDataHist* AliEveLego::LoadData()
256 // Load data from ESD tree
259 fHistoElectrons->Reset();
260 fHistoMuons->Reset();
261 fHistoPions->Reset();
262 fHistoKaons->Reset();
263 fHistoProtons->Reset();
265 // Getting current tracks, filling histograms
266 for (int n = 0; n < fEsd->GetNumberOfTracks(); ++n) {
268 AliESDtrack *track = fEsd->GetTrack(n);
269 const Double_t sign = fEsd->GetTrack(n)->GetSign();
270 const Int_t prob = GetParticleType(track);
272 // Filling histograms
274 fHistopos->Fill(track->Eta(), getphi(track->Phi()), fabs(track->Pt()));
277 fHistoneg->Fill(track->Eta(), getphi(track->Phi()), fabs(track->Pt()));
280 fHistoElectrons->Fill(track->Eta(), getphi(track->Phi()), fabs(track->Pt()));
283 fHistoMuons->Fill(track->Eta(), getphi(track->Phi()), fabs(track->Pt()));
286 fHistoPions->Fill(track->Eta(), getphi(track->Phi()), fabs(track->Pt()));
289 fHistoKaons->Fill(track->Eta(), getphi(track->Phi()), fabs(track->Pt()));
292 fHistoProtons->Fill(track->Eta(), getphi(track->Phi()), fabs(track->Pt()));
295 fData->DataChanged();
302 //______________________________________________________________________________
303 TEveCaloDataHist* AliEveLego::LoadAllData()
305 // Load data from all events ESD
306 fHistoposAllEvents->Reset();
307 fHistonegAllEvents->Reset();
308 fHistoElectronsAllEvents->Reset();
309 fHistoMuonsAllEvents->Reset();
310 fHistoPionsAllEvents->Reset();
311 fHistoKaonsAllEvents->Reset();
312 fHistoProtonsAllEvents->Reset();
314 TTree* t = AliEveEventManager::GetMaster()->GetESDTree();
316 // Getting current tracks for each event, filling histograms
317 Int_t fAcceptedEvents = 0;
318 for (int event = 0; event < t->GetEntries(); event++) {
321 if (fCollisionCandidatesOnly == kTRUE)
322 if (fPhysicsSelection->IsCollisionCandidate(fEsd) == kFALSE) continue;
326 for (int n = 0; n < fEsd->GetNumberOfTracks(); ++n) {
328 AliESDtrack *track = fEsd->GetTrack(n);
329 const Double_t sign = track->GetSign();
330 const Int_t prob = GetParticleType(track);
333 fHistoposAllEvents->Fill(track->Eta(), getphi(track->Phi()), fabs(track->Pt()));
336 fHistonegAllEvents->Fill(track->Eta(), getphi(track->Phi()), fabs(track->Pt()));
339 fHistoElectronsAllEvents->Fill(track->Eta(), getphi(track->Phi()), fabs(track->Pt()));
342 fHistoMuonsAllEvents->Fill(track->Eta(), getphi(track->Phi()), fabs(track->Pt()));
345 fHistoPionsAllEvents->Fill(track->Eta(), getphi(track->Phi()), fabs(track->Pt()));
348 fHistoKaonsAllEvents->Fill(track->Eta(), getphi(track->Phi()), fabs(track->Pt()));
351 fHistoProtonsAllEvents->Fill(track->Eta(), getphi(track->Phi()), fabs(track->Pt()));
355 // Setting the current view to the first event
358 // Usefull information,
359 // with this we can estimate the event efficiency
360 printf("Number of events loaded: %i, with AliPhysicsSelection: %i\n",fAcceptedEvents,fCollisionCandidatesOnly);
362 fDataAllEvents->DataChanged();
364 return fDataAllEvents;
367 //______________________________________________________________________________
368 TEveCaloDataHist* AliEveLego::FilterData()
371 if ( fTracksId == 2 )
375 fHistoElectrons->Reset();
376 fHistoMuons->Reset();
377 fHistoPions->Reset();
378 fHistoKaons->Reset();
379 fHistoProtons->Reset();
381 const AliESDVertex *pv = fEsd->GetPrimaryVertex();
383 for (Int_t n = 0; n < pv->GetNIndices(); n++ )
385 AliESDtrack *at = fEsd->GetTrack(pv->GetIndices()[n]);
386 const Double_t sign = at->GetSign();
387 const Int_t prob = GetParticleType(at);
390 fHistopos->Fill(at->Eta(), getphi(at->Phi()), fabs(at->Pt()));
393 fHistoneg->Fill(at->Eta(), getphi(at->Phi()), fabs(at->Pt()));
396 fHistoElectrons->Fill(at->Eta(), getphi(at->Phi()), fabs(at->Pt()));
399 fHistoMuons->Fill(at->Eta(), getphi(at->Phi()), fabs(at->Pt()));
402 fHistoPions->Fill(at->Eta(), getphi(at->Phi()), fabs(at->Pt()));
405 fHistoKaons->Fill(at->Eta(), getphi(at->Phi()), fabs(at->Pt()));
408 fHistoProtons->Fill(at->Eta(), getphi(at->Phi()), fabs(at->Pt()));
412 fData->DataChanged();
415 if (GetPtMax() >= fMaxPt){
416 for (Int_t binx = 1; binx <= 100; binx++) {
417 for (Int_t biny = 1; biny <= 80; biny++) {
419 if (fHistopos->GetBinContent(binx, biny) >= fMaxPt)
420 fHistopos->SetBinContent(binx, biny, fMaxPt);
422 if (fHistoneg->GetBinContent(binx, biny) >= fMaxPt)
423 fHistoneg->SetBinContent(binx, biny, fMaxPt);
425 if (fHistoElectrons->GetBinContent(binx, biny) >= fMaxPt)
426 fHistoElectrons->SetBinContent(binx, biny, fMaxPt);
428 if (fHistoMuons->GetBinContent(binx, biny) >= fMaxPt)
429 fHistoMuons->SetBinContent(binx, biny, fMaxPt);
431 if (fHistoPions->GetBinContent(binx, biny) >= fMaxPt)
432 fHistoPions->SetBinContent(binx, biny, fMaxPt);
434 if (fHistoKaons->GetBinContent(binx, biny) >= fMaxPt)
435 fHistoKaons->SetBinContent(binx, biny, fMaxPt);
437 if (fHistoProtons->GetBinContent(binx, biny) >= fMaxPt)
438 fHistoProtons->SetBinContent(binx, biny, fMaxPt);
443 // Particle type filter
444 if ( fParticleTypeId[0] == kFALSE) fHistopos->Reset();
445 if ( fParticleTypeId[1] == kFALSE) fHistoneg->Reset();
446 if ( fParticleTypeId[2] == kFALSE) fHistoElectrons->Reset();
447 if ( fParticleTypeId[3] == kFALSE) fHistoMuons->Reset();
448 if ( fParticleTypeId[4] == kFALSE) fHistoPions->Reset();
449 if ( fParticleTypeId[5] == kFALSE) fHistoKaons->Reset();
450 if ( fParticleTypeId[6] == kFALSE) fHistoProtons->Reset();
452 fData->DataChanged();
457 //______________________________________________________________________________
458 TEveCaloDataHist* AliEveLego::FilterAllData()
461 if ( fTracksIdAE == 2 )
463 fHistoposAllEvents->Reset();
464 fHistonegAllEvents->Reset();
465 fHistoElectronsAllEvents->Reset();
466 fHistoMuonsAllEvents->Reset();
467 fHistoPionsAllEvents->Reset();
468 fHistoKaonsAllEvents->Reset();
469 fHistoProtonsAllEvents->Reset();
471 TTree* t = AliEveEventManager::GetMaster()->GetESDTree();
473 // Getting current tracks for each event, filling histograms
474 Int_t fAcceptedEvents = 0;
475 for (int event = 0; event < t->GetEntries(); event++) {
478 const AliESDVertex *pv = fEsd->GetPrimaryVertex();
480 if (fCollisionCandidatesOnly == kTRUE)
481 if (fPhysicsSelection->IsCollisionCandidate(fEsd) == kFALSE) continue;
485 for (Int_t n = 0; n < pv->GetNIndices(); n++ )
487 AliESDtrack *at = fEsd->GetTrack(pv->GetIndices()[n]);
488 const Double_t sign = at->GetSign();
489 const Int_t prob = GetParticleType(at);
492 fHistoposAllEvents->Fill(at->Eta(), getphi(at->Phi()), fabs(at->Pt()));
495 fHistonegAllEvents->Fill(at->Eta(), getphi(at->Phi()), fabs(at->Pt()));
498 fHistoElectronsAllEvents->Fill(at->Eta(), getphi(at->Phi()), fabs(at->Pt()));
501 fHistoMuonsAllEvents->Fill(at->Eta(), getphi(at->Phi()), fabs(at->Pt()));
504 fHistoPionsAllEvents->Fill(at->Eta(), getphi(at->Phi()), fabs(at->Pt()));
507 fHistoKaonsAllEvents->Fill(at->Eta(), getphi(at->Phi()), fabs(at->Pt()));
510 fHistoProtonsAllEvents->Fill(at->Eta(), getphi(at->Phi()), fabs(at->Pt()));
514 printf("Number of events loaded: %i, with AliPhysicsSelection: %i\n",fAcceptedEvents,fCollisionCandidatesOnly);
522 fDataAllEvents->DataChanged();
525 if (GetPtMaxAE() >= fMaxPtAE){
526 for (Int_t binx = 1; binx <= 100; binx++) {
527 for (Int_t biny = 1; biny <= 80; biny++) {
529 if (fHistoposAllEvents->GetBinContent(binx, biny) >= fMaxPtAE)
530 fHistoposAllEvents->SetBinContent(binx, biny, fMaxPtAE);
532 if (fHistonegAllEvents->GetBinContent(binx, biny) >= fMaxPtAE)
533 fHistonegAllEvents->SetBinContent(binx, biny, fMaxPtAE);
535 if (fHistoElectronsAllEvents->GetBinContent(binx, biny) >= fMaxPt)
536 fHistoElectronsAllEvents->SetBinContent(binx, biny, fMaxPt);
538 if (fHistoMuonsAllEvents->GetBinContent(binx, biny) >= fMaxPt)
539 fHistoMuonsAllEvents->SetBinContent(binx, biny, fMaxPt);
541 if (fHistoPionsAllEvents->GetBinContent(binx, biny) >= fMaxPt)
542 fHistoPionsAllEvents->SetBinContent(binx, biny, fMaxPt);
544 if (fHistoKaonsAllEvents->GetBinContent(binx, biny) >= fMaxPt)
545 fHistoKaonsAllEvents->SetBinContent(binx, biny, fMaxPt);
547 if (fHistoProtonsAllEvents->GetBinContent(binx, biny) >= fMaxPt)
548 fHistoProtonsAllEvents->SetBinContent(binx, biny, fMaxPt);
553 // Particles species and charges filter
554 if ( fParticleTypeIdAE[0] == kFALSE ) fHistoposAllEvents->Reset();
555 if ( fParticleTypeIdAE[1] == kFALSE ) fHistonegAllEvents->Reset();
556 if ( fParticleTypeIdAE[2] == kFALSE ) fHistoElectronsAllEvents->Reset();
557 if ( fParticleTypeIdAE[3] == kFALSE ) fHistoMuonsAllEvents->Reset();
558 if ( fParticleTypeIdAE[4] == kFALSE ) fHistoPionsAllEvents->Reset();
559 if ( fParticleTypeIdAE[5] == kFALSE ) fHistoKaonsAllEvents->Reset();
560 if ( fParticleTypeIdAE[6] == kFALSE ) fHistoProtonsAllEvents->Reset();
562 fDataAllEvents->DataChanged();
564 gEve->Redraw3D(kTRUE);
566 return fDataAllEvents;
569 //______________________________________________________________________________
570 void AliEveLego::Update()
575 // Create new histo2d
581 // Update the viewers
582 gEve->Redraw3D(kTRUE);
585 //______________________________________________________________________________
586 TEveCaloLego* AliEveLego::CreateHistoLego()
588 // Viewer initialization, tab creation
589 if (fHisto2dv == 0) {
590 TEveWindowSlot *fslot = 0;
591 TEveBrowser *fbrowser = gEve->GetBrowser();
593 fslot = TEveWindow::CreateWindowInTab(fbrowser->GetTabRight());
594 fslot->MakeCurrent();
595 fHisto2dv = gEve->SpawnNewViewer("2D Lego Histogram", "2D Lego Histogram");
596 fHisto2ds = gEve->SpawnNewScene("2D Lego Histogram", "2D Lego Histogram");
597 fHisto2dv->AddScene(fHisto2ds);
598 fHisto2dv->SetElementName("2D Lego Viewer");
599 fHisto2ds->SetElementName("2D Lego Scene");
601 fGlv = fHisto2dv->GetGLViewer();
602 fHisto2dLegoOverlay = new TEveCaloLegoOverlay();
603 fGlv->AddOverlayElement(fHisto2dLegoOverlay);
604 fGlv->SetCurrentCamera(TGLViewer::kCameraPerspXOY);
606 fHisto2ds->AddElement(fLego);
608 // move to real world coordinates
609 fLego->InitMainTrans();
610 Float_t sc = TMath::Min(fLego->GetEtaRng(), fLego->GetPhiRng());
611 fLego->RefMainTrans().SetScale(sc, sc, sc);
613 // set event handler to move from perspective to orthographic view.
614 fGlv->SetEventHandler(new TEveLegoEventHandler(fGlv->GetGLWidget(), fGlv, fLego));
616 fHisto2dLegoOverlay->SetCaloLego(fLego);
622 //______________________________________________________________________________
623 TEveCaloLego* AliEveLego::CreateHistoLego(TEveWindowSlot *slot)
625 // Viewer initialization, tab creation
626 if (fHisto2dAllEventsv0 == 0) {
629 fHisto2dAllEventsv0 = gEve->SpawnNewViewer("2D Lego Histogram", "2D Lego Histogram");
630 fHisto2dAllEventss0 = gEve->SpawnNewScene("2D Lego Histogram", "2D Lego Histogram");
631 fHisto2dAllEventsv0->AddScene(fHisto2dAllEventss0);
632 fHisto2dAllEventsv0->SetElementName("2D Lego Viewer");
633 fHisto2dAllEventss0->SetElementName("2D Lego Scene");
635 TGLViewer* glv = fHisto2dAllEventsv0->GetGLViewer();
636 fHisto2dAllEventsLegoOverlay = new TEveCaloLegoOverlay();
637 glv->AddOverlayElement(fHisto2dAllEventsLegoOverlay);
638 glv->SetCurrentCamera(TGLViewer::kCameraPerspXOY);
640 // Plotting histogram lego
641 fLegoAllEvents = new TEveCaloLego(fDataAllEvents);
642 fHisto2dAllEventss0->AddElement(fLegoAllEvents);
644 // Move to real world coordinates
645 fLegoAllEvents->InitMainTrans();
646 Float_t sc = TMath::Min(fLegoAllEvents->GetEtaRng(), fLegoAllEvents->GetPhiRng());
647 fLegoAllEvents->RefMainTrans().SetScale(sc, sc, sc);
649 // Set event handler to move from perspective to orthographic view.
650 glv->SetEventHandler(new TEveLegoEventHandler(glv->GetGLWidget(), glv, fLegoAllEvents));
652 fHisto2dAllEventsLegoOverlay->SetCaloLego(fLegoAllEvents);
655 return fLegoAllEvents;
658 //______________________________________________________________________________
659 TEveCalo3D* AliEveLego::Create3DView()
662 if (fHisto2ds2 == 0) {
663 fHisto2ds2 = gEve->SpawnNewScene("3D Histogram", "3D Histogram");
664 gEve->GetDefaultViewer()->AddScene(fHisto2ds2);
665 fHisto2ds2->SetElementName("3D Histogram Scene");
666 fHisto2ds2->AddElement(fCalo3d);
672 //______________________________________________________________________________
673 TEveCalo3D* AliEveLego::Create3DView(TEveWindowSlot *slot)
675 // Creates a 3d view for the 3d histogram
676 if ( fHisto2dAllEventsv1 == 0 ) {
679 fHisto2dAllEventsv1 = gEve->SpawnNewViewer("3D Histogram", "3D Histogram");
680 fHisto2dAllEventss1 = gEve->SpawnNewScene("3D Histogram", "3D Histogram");
681 fHisto2dAllEventsv1->AddScene(fHisto2dAllEventss1);
682 fHisto2dAllEventsv1->SetElementName("3D Histogram Viewer");
683 fHisto2dAllEventss1->SetElementName("3D Histogram Scene");
685 fCalo3dAllEvents = new TEveCalo3D(fDataAllEvents);
687 fCalo3dAllEvents->SetBarrelRadius(550);
688 fCalo3dAllEvents->SetEndCapPos(550);
689 fHisto2dAllEventss1->AddElement(fCalo3dAllEvents);
692 return fCalo3dAllEvents;
695 //______________________________________________________________________________
696 void AliEveLego::CreateProjections(TEveWindowSlot* slot1, TEveWindowSlot* slot2)
698 // Create projections
699 if (fHisto2dAllEventsv2 == 0) {
701 slot1->MakeCurrent();
702 fHisto2dAllEventsv2 = gEve->SpawnNewViewer("RPhi projection", "RPhi projection");
703 fHisto2dAllEventss2 = gEve->SpawnNewScene("RPhi projection", "RPhi projection");
704 fHisto2dAllEventsv2->AddScene(fHisto2dAllEventss2);
705 fHisto2dAllEventsv2->SetElementName("RPhi Projection Viewer");
706 fHisto2dAllEventss2->SetElementName("RPhi Projection Scene");
708 TEveProjectionManager* mng1 = new TEveProjectionManager();
709 mng1->SetProjection(TEveProjection::kPT_RPhi);
711 TEveProjectionAxes* axeghisto2dAllEventss1 = new TEveProjectionAxes(mng1);
712 fHisto2dAllEventss2->AddElement(axeghisto2dAllEventss1);
713 TEveCalo2D* fcalo2d1 = (TEveCalo2D*) mng1->ImportElements(fCalo3dAllEvents);
714 fHisto2dAllEventss2->AddElement(fcalo2d1);
716 fHisto2dAllEventsv2->GetGLViewer()->SetCurrentCamera(TGLViewer::kCameraOrthoXOY);
719 if (fHisto2dAllEventsv3 == 0)
721 slot2->MakeCurrent();
722 fHisto2dAllEventsv3 = gEve->SpawnNewViewer("RhoZ projection", "RhoZ projection");
723 fHisto2dAllEventss3 = gEve->SpawnNewScene("RhoZ projection", "RhoZ projection");
724 fHisto2dAllEventsv3->AddScene(fHisto2dAllEventss3);
725 fHisto2dAllEventsv3->SetElementName("RhoZ Projection Viewer");
726 fHisto2dAllEventss3->SetElementName("RhoZ Projection Viewer");
728 TEveProjectionManager* mng2 = new TEveProjectionManager();
729 mng2->SetProjection(TEveProjection::kPT_RhoZ);
731 TEveProjectionAxes* axeghisto2dAllEventss2 = new TEveProjectionAxes(mng2);
732 fHisto2dAllEventss3->AddElement(axeghisto2dAllEventss2);
733 TEveCalo2D* fcalo2d2 = (TEveCalo2D*) mng2->ImportElements(fCalo3dAllEvents);
734 fHisto2dAllEventss3->AddElement(fcalo2d2);
736 fHisto2dAllEventsv3->GetGLViewer()->SetCurrentCamera(TGLViewer::kCameraOrthoXOY);
741 //______________________________________________________________________________
742 TEveCaloDataHist* AliEveLego::LoadAllEvents()
744 // load all events data from ESD
745 if ( fHisto2dAllEventsSlot == 0 ) {
747 printf("Filling histogram...\n");
751 // Creating 2D histograms
752 fHistoposAllEvents = new TH2F("fHistoposAllEvents","Histo 2d positive",
753 100,-1.5,1.5,80,-kPi,kPi);
754 fHistonegAllEvents = new TH2F("fHistonegAllEvents","Histo 2d negative",
755 100,-1.5,1.5,80,-kPi,kPi);
756 fHistoElectronsAllEvents = new TH2F("fHistoElectronsAllEvents","Histo 2d electrons",
757 100,-1.5,1.5,80,-kPi,kPi);
758 fHistoMuonsAllEvents = new TH2F("fHistoMuonsAllEvents","Histo 2d muons",
759 100,-1.5,1.5,80,-kPi,kPi);
760 fHistoPionsAllEvents = new TH2F("fHistoPionsAllEvents","Histo 2d pions",
761 100,-1.5,1.5,80,-kPi,kPi);
762 fHistoKaonsAllEvents = new TH2F("fHistoKaonsAllEvents","Histo 2d kaons",
763 100,-1.5,1.5,80,-kPi,kPi);
764 fHistoProtonsAllEvents = new TH2F("fHistoProtonsAllEvents","Histo 2d protons",
765 100,-1.5,1.5,80,-kPi,kPi);
767 fHistoposAllEvents->SetDirectory(0);
768 fHistonegAllEvents->SetDirectory(0);
769 fHistoElectronsAllEvents->SetDirectory(0);
770 fHistoMuonsAllEvents->SetDirectory(0);
771 fHistoPionsAllEvents->SetDirectory(0);
772 fHistoKaonsAllEvents->SetDirectory(0);
773 fHistoProtonsAllEvents->SetDirectory(0);
775 // colors from get_pdg_color() in /alice-macros/kine_tracks.C
776 static Color_t fElectro = 5;
777 static Color_t fMuon = 30;
778 static Color_t fPion = 3;
779 static Color_t fKaon = 38;
780 static Color_t fProton = 10;
782 fDataAllEvents = new TEveCaloDataHist();
783 fDataAllEvents->AddHistogram(fHistonegAllEvents);
784 fDataAllEvents->RefSliceInfo(0).Setup("NegCg:", 0, kBlue);
785 fDataAllEvents->AddHistogram(fHistoposAllEvents);
786 fDataAllEvents->RefSliceInfo(1).Setup("PosCg:", 0, kRed);
787 fDataAllEvents->AddHistogram(fHistoElectronsAllEvents);
788 fDataAllEvents->RefSliceInfo(2).Setup("Electrons:", 0, fElectro);
789 fDataAllEvents->AddHistogram(fHistoMuonsAllEvents);
790 fDataAllEvents->RefSliceInfo(3).Setup("Muons:", 0, fMuon);
791 fDataAllEvents->AddHistogram(fHistoPionsAllEvents);
792 fDataAllEvents->RefSliceInfo(4).Setup("Pions:", 0, fPion);
793 fDataAllEvents->AddHistogram(fHistoKaonsAllEvents);
794 fDataAllEvents->RefSliceInfo(5).Setup("Kaons:", 0, fKaon);
795 fDataAllEvents->AddHistogram(fHistoProtonsAllEvents);
796 fDataAllEvents->RefSliceInfo(6).Setup("Protons:", 0, fProton);
798 fDataAllEvents->GetEtaBins()->SetTitleFont(120);
799 fDataAllEvents->GetEtaBins()->SetTitle("h");
800 fDataAllEvents->GetPhiBins()->SetTitleFont(120);
801 fDataAllEvents->GetPhiBins()->SetTitle("f");
802 fDataAllEvents->IncDenyDestroy();
805 fHisto2dAllEventsSlot = TEveWindow::CreateWindowInTab(gEve->GetBrowser()->GetTabRight());
806 TEveWindowPack* packH = fHisto2dAllEventsSlot->MakePack();
807 packH->SetElementName("Projections");
808 packH->SetHorizontal();
809 packH->SetShowTitleBar(kFALSE);
811 fHisto2dAllEventsSlot = packH->NewSlot();
812 TEveWindowPack* pack0 = fHisto2dAllEventsSlot->MakePack();
813 pack0->SetShowTitleBar(kFALSE);
814 TEveWindowSlot* slotLeftTop = pack0->NewSlot();
815 TEveWindowSlot* slotLeftBottom = pack0->NewSlot();
817 fHisto2dAllEventsSlot = packH->NewSlot();
818 TEveWindowPack* pack1 = fHisto2dAllEventsSlot->MakePack();
819 pack1->SetShowTitleBar(kFALSE);
820 TEveWindowSlot* slotRightTop = pack1->NewSlot();
821 TEveWindowSlot* slotRightBottom = pack1->NewSlot();
823 // Creating viewers and scenes
824 Create3DView(slotLeftTop);
825 CreateHistoLego(slotLeftBottom);
826 CreateProjections(slotRightTop, slotRightBottom);
830 gEve->Redraw3D(kTRUE);
832 printf("Filling histogram... Finished\n");
836 return fDataAllEvents;
839 //______________________________________________________________________________
840 void AliEveLego::ApplyParticleTypeSelectionAE()
842 // Reload all events applying particle type selection
846 //______________________________________________________________________________
847 Float_t AliEveLego::GetPtMax()
850 return fData->GetMaxVal(fLego->GetPlotEt());
853 //______________________________________________________________________________
854 Float_t AliEveLego::GetPtMaxAE()
856 // Return pT max from all events
857 return fDataAllEvents->GetMaxVal(fLegoAllEvents->GetPlotEt());
860 //______________________________________________________________________________
861 Int_t AliEveLego::GetParticleType(AliESDtrack *track)
863 // Determine the particle type
864 Double_t prob[5]={0.};
865 track->GetESDpid(prob);
867 Double_t max = prob[0];
870 for (Int_t i = 1 ; i < 5; i++)
880 //______________________________________________________________________________
881 void AliEveLego::SetParticleType(Int_t id, Bool_t status)
883 // Activate/deactivate particles types
884 fParticleTypeId[id] = status;
889 //______________________________________________________________________________
890 void AliEveLego::SetParticleTypeAE(Int_t id, Bool_t status)
892 // Activate/deactivate particles types
893 fParticleTypeIdAE[id] = status;
896 //______________________________________________________________________________
897 void AliEveLego::SetMaxPt(Double_t val)
904 //______________________________________________________________________________
905 void AliEveLego::SetMaxPtAE(Double_t val)
912 //______________________________________________________________________________
913 void AliEveLego::SetThreshold(Double_t val)
915 // Setting up the new threshold for all histograms
916 fData->SetSliceThreshold(0,val);
917 fData->SetSliceThreshold(1,val);
918 fData->SetSliceThreshold(2,val);
919 fData->SetSliceThreshold(3,val);
920 fData->SetSliceThreshold(4,val);
921 fData->SetSliceThreshold(5,val);
922 fData->SetSliceThreshold(6,val);
923 fData->DataChanged();
925 gEve->Redraw3D(kTRUE);
928 //______________________________________________________________________________
929 void AliEveLego::SetThresholdAE(Double_t val)
931 // Setting up the new threshold for all histograms
932 fDataAllEvents->SetSliceThreshold(0,val);
933 fDataAllEvents->SetSliceThreshold(1,val);
934 fDataAllEvents->SetSliceThreshold(2,val);
935 fDataAllEvents->SetSliceThreshold(3,val);
936 fDataAllEvents->SetSliceThreshold(4,val);
937 fDataAllEvents->SetSliceThreshold(5,val);
938 fDataAllEvents->SetSliceThreshold(6,val);
939 fDataAllEvents->DataChanged();
941 gEve->Redraw3D(kTRUE);
944 //______________________________________________________________________________
945 void AliEveLego::SwitchDataType(Bool_t status)
947 // Activate/deactivate MC / real data type
950 // Removing defaul physics selection
951 delete fPhysicsSelection;
952 fPhysicsSelection = NULL;
954 // Re-initialization of physics selection
955 fPhysicsSelection = new AliPhysicsSelection();
956 fPhysicsSelection->SetAnalyzeMC(fIsMC);
957 fPhysicsSelection->Initialize(fEsd);
961 //______________________________________________________________________________
962 void AliEveLego::SetCollisionCandidatesOnly()
964 // Activate/deactivate MC / real data type
965 if (fCollisionCandidatesOnly == 0)
967 fCollisionCandidatesOnly = 1;
969 fCollisionCandidatesOnly = 0;
974 /******************************************************************************/