AliEveMacro
[u/mrichter/AliRoot.git] / EVE / alice-macros / vplot_tpc.C
1 // $Id$
2 // Main authors: Adam Jacholkowski & Matevz Tadel: 2009
3
4 /**************************************************************************
5  * Copyright(c) 1998-2008, ALICE Experiment at CERN, all rights reserved. *
6  * See http://aliceinfo.cern.ch/Offline/AliRoot/License.html for          *
7  * full copyright notice.                                                 *
8  **************************************************************************/
9
10 #ifndef __CINT__
11
12 #include <TGLViewer.h>
13 #include <TGLCameraOverlay.h>
14
15 #include <TEveManager.h>
16 #include <TEveBrowser.h>
17 #include <TEveViewer.h>
18 #include <TEveScene.h>
19 #include <TEvePointSet.h>
20
21 #include <EveBase/AliEveEventManager.h>
22
23 #include <AliRunLoader.h>
24 #include <AliCluster.h>
25 #include <TPC/AliTPCClustersRow.h>
26
27 #endif
28
29 TEveViewer *gVPTPCView   = 0;
30 TEveScene  *gVPTPCScene  = 0;
31
32 TEvePointSet* vplot_tpc(TEveElement* cont=0, Float_t maxR=270)
33 {
34   if (gVPTPCView == 0)
35   {
36     TEveWindowSlot *slot    = 0;
37     TEveBrowser    *browser = gEve->GetBrowser();
38
39     slot = TEveWindow::CreateWindowInTab(browser->GetTabRight());
40     slot->MakeCurrent();
41     gVPTPCView  = gEve->SpawnNewViewer("V-Plot", "");
42     gVPTPCScene = gEve->SpawnNewScene("V-Plot", "Scene holding elements for the V-Plot TPC.");
43     gVPTPCView->AddScene(gVPTPCScene);
44
45     TGLViewer *glv = gVPTPCView->GetGLViewer();
46     glv->SetCurrentCamera(TGLViewer::kCameraOrthoXOY);
47     glv->ResetCamerasAfterNextUpdate();
48
49     TGLCameraOverlay* co = glv->GetCameraOverlay();
50     co->SetShowOrthographic(true); //(false);
51     co->SetOrthographicMode(TGLCameraOverlay::kAxis); // ::kPlaneIntersect or ::kBar
52   }
53
54   const Int_t kMaxCl=100*160;
55
56   AliEveEventManager::AssertGeometry();
57
58   AliRunLoader* rl = AliEveEventManager::AssertRunLoader();
59   rl->LoadRecPoints("TPC");
60
61   TTree *cTree = rl->GetTreeR("TPC", false);
62   if (cTree == 0)
63     return 0;
64
65   AliTPCClustersRow *clrow = new AliTPCClustersRow();
66   clrow->SetClass("AliTPCclusterMI");
67   clrow->SetArray(kMaxCl);
68   cTree->SetBranchAddress("Segment", &clrow);
69
70   TEvePointSet* vplot = new TEvePointSet(kMaxCl);
71   vplot->SetOwnIds(kTRUE);
72
73   const Float_t phimin = -3.15;
74   const Float_t phimax =  3.15;
75   const Float_t etamin = -1.2;
76   const Float_t etamax =  1.2;
77   const Float_t vconst =  0.0003;
78
79   Float_t rhomax = 246.6;
80   Float_t rholim, rhoup;
81   Float_t zmax = 250.0;
82   Float_t rho, eta, phi, theta;
83   Float_t r3d, r3dmax, r3d1, r3d2;
84   //
85   Float_t maxRsqr = maxR*maxR;
86
87   Int_t nentr = (Int_t) cTree->GetEntries();
88   for (Int_t i = 0; i < nentr; ++i)
89   {
90     if (!cTree->GetEvent(i)) continue;
91
92     TClonesArray *cl = clrow->GetArray();
93     Int_t ncl = cl->GetEntriesFast();
94
95     while (ncl--)
96     {
97       AliCluster *c = (AliCluster*) cl->UncheckedAt(ncl);
98       Float_t g[3]; //global coordinates
99       c->GetGlobalXYZ(g);
100       if (g[0]*g[0] + g[1]*g[1] < maxRsqr)
101       {
102         phi    =  TMath::ATan2(g[1], g[0]);
103         rho    =  TMath::Sqrt(g[0]*g[0] + g[1]*g[1]);
104         theta  =  TMath::ATan2(rho, g[2]);
105         eta    = -1.0*TMath::Log(TMath::Tan(0.5*theta));
106         rhoup  =  zmax*rho / TMath::Abs(g[2]);
107         rholim =  TMath::Min(rhoup,rhomax);
108         // a version using rather r3d
109         r3d    =  TMath::Sqrt(rho*rho + g[2]*g[2]);
110         r3d1   =  rhomax / TMath::Sin(theta);
111         r3d2   =  TMath::Abs(zmax/TMath::Cos(theta));
112         r3dmax =  TMath::Min(r3d1, r3d2);
113
114         if (eta>etamin && eta<etamax && phi>phimin && phi<phimax)
115         {
116           Float_t deta = vconst*(r3dmax - r3d);
117
118           vplot->SetNextPoint(eta + deta, phi, 0);
119           vplot->SetNextPoint(eta - deta, phi, 0);
120         }
121       }
122     }
123     cl->Clear();
124   }
125
126   delete clrow;
127
128   rl->UnloadRecPoints("TPC");
129
130   if (vplot->Size() == 0 && gEve->GetKeepEmptyCont() == kFALSE)
131   {
132     Warning("vplot_tpc.C", "No TPC clusters were found.");
133     delete vplot;
134     return 0;
135   }
136
137   vplot->SetName("V=Plot TPC");
138   vplot->SetTitle(Form("N=%d", vplot->Size() / 2));
139
140   vplot->SetMainColor (kOrange);
141   vplot->SetMarkerSize(0.2);
142   vplot->SetMarkerStyle(1);
143   // vplot->ApplyVizTag("V-Plot TPC", "V-Plot");
144
145   if (cont)
146   {
147     cont->AddElement(vplot);
148   }
149   else
150   {
151     gVPTPCScene->AddElement(vplot);
152   }
153   AliEveEventManager::RegisterTransient(vplot);
154
155   gEve->Redraw3D();
156
157   return vplot;
158 }