]>
Commit | Line | Data |
---|---|---|
455b87be | 1 | // Main authors: Matevz Tadel & Alja Mrak-Tadel: 2006, 2007 |
2 | ||
3 | /************************************************************************** | |
4 | * Copyright(c) 1998-2008, ALICE Experiment at CERN, all rights reserved. * | |
5 | * See http://aliceinfo.cern.ch/Offline/AliRoot/License.html for * | |
6 | * full copyright notice. * | |
7 | **************************************************************************/ | |
455b87be | 8 | #ifdef __CINT__ |
9 | ||
a172375c | 10 | class TEveElement; |
11 | class TEvePointSet; | |
455b87be | 12 | |
13 | #else | |
14 | ||
a172375c | 15 | #include <TEveManager.h> |
16 | #include <TEvePointSet.h> | |
17 | #include <EveBase/AliEveEventManager.h> | |
455b87be | 18 | |
19 | #include <AliRunLoader.h> | |
20 | #include <AliCluster.h> | |
21 | ||
22 | #include <TClonesArray.h> | |
23 | ||
24 | #endif | |
25 | ||
26 | TEvePointSet* tof_clusters(TEveElement* cont=0, Float_t maxR=390) | |
27 | { | |
455b87be | 28 | AliEveEventManager::AssertGeometry(); |
29 | ||
30 | AliRunLoader* rl = AliEveEventManager::AssertRunLoader(); | |
31 | rl->LoadRecPoints("TOF"); | |
32 | ||
33 | TTree *cTree = rl->GetTreeR("TOF", false); | |
a172375c | 34 | if (cTree == 0) |
35 | return 0; | |
455b87be | 36 | |
a172375c | 37 | TClonesArray *cl = NULL; |
38 | TBranch *branch = cTree->GetBranch("TOF"); | |
455b87be | 39 | branch->SetAddress(&cl); |
40 | ||
a172375c | 41 | TEvePointSet* clusters = new TEvePointSet(10000); |
42 | clusters->SetOwnIds(kTRUE); | |
43 | ||
455b87be | 44 | Int_t nentr=(Int_t)cTree->GetEntries(); |
45 | for (Int_t i=0; i<nentr; i++) { | |
46 | if (!cTree->GetEvent(i)) continue; | |
47 | ||
48 | Int_t ncl=cl->GetEntriesFast(); | |
49 | cout<<" ncl = "<<ncl<<endl; | |
50 | Float_t maxRsqr = maxR*maxR; | |
51 | while (ncl--) { | |
455b87be | 52 | AliCluster *c=(AliCluster*)cl->UncheckedAt(ncl); |
455b87be | 53 | Float_t g[3]; //global coordinates |
54 | c->GetGlobalXYZ(g); | |
455b87be | 55 | if (g[0]*g[0]+g[1]*g[1] < maxRsqr) |
56 | { | |
57 | clusters->SetNextPoint(g[0], g[1], g[2]); | |
58 | AliCluster *atp = new AliCluster(*c); | |
455b87be | 59 | clusters->SetPointId(atp); |
a172375c | 60 | } |
455b87be | 61 | } |
62 | } | |
63 | ||
a172375c | 64 | if (clusters->Size() == 0 && gEve->GetKeepEmptyCont() == kFALSE) { |
65 | Warning("tof_clusters.C", "No TOF clusters"); | |
455b87be | 66 | delete clusters; |
67 | return 0; | |
68 | } | |
69 | ||
70 | clusters->SetMarkerStyle(2); | |
71 | clusters->SetMarkerSize(0.2); | |
72 | clusters->SetMarkerColor(4); | |
73 | ||
74 | char form[1000]; | |
75 | sprintf(form,"TOF Clusters"); | |
76 | clusters->SetName(form); | |
77 | ||
78 | char tip[1000]; | |
79 | sprintf(tip,"N=%d", clusters->Size()); | |
80 | clusters->SetTitle(tip); | |
455b87be | 81 | gEve->AddElement(clusters, cont); |
82 | gEve->Redraw3D(); | |
83 | ||
84 | return clusters; | |
85 | } |