clusterizer,reconstructor + many fixes (Magnus,Stefan)
[u/mrichter/AliRoot.git] / EVE / alice-macros / itsU_clusters.C
1 // $Id: its_clusters.C 55060 2012-03-09 18:13:17Z quark $
2 // Main authors: Matevz Tadel & Alja Mrak-Tadel: 2006, 2007
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 #if !defined(__CINT__) || defined(__MAKECINT__)
11 #include <TClonesArray.h>
12 #include <TBranch.h>
13 #include <TTree.h>
14 #include <TEveManager.h>
15 #include <TEveElement.h>
16 #include <TEvePointSet.h>
17
18 #include <AliRunLoader.h>
19 #include <AliCluster.h>
20 #include <AliEveEventManager.h>
21 #include <AliGeomManager.h>
22 #include <AliITSUGeomTGeo.h>
23 #include <AliITSUSegmentationPix.h>
24 #include <TGeoManager.h>
25
26 #else
27 class TEveElement;
28 class TEvePointSet;
29 class TTree;
30 class TBranch;
31 #endif
32
33 void itsU_clusters(TEveElement* cont=0, Float_t maxR=50)
34 {
35   AliEveEventManager::AssertGeometry();
36
37   AliRunLoader* rl = AliEveEventManager::AssertRunLoader();
38   rl->LoadRecPoints("ITS");
39
40   gGeoManager = gEve->GetGeometry("geometry.root");
41   AliITSUGeomTGeo* gm = new AliITSUGeomTGeo(kTRUE);
42   TObjArray segmArr;
43   AliITSUSegmentationPix::LoadSegmentations(&segmArr, AliITSUGeomTGeo::GetITSsegmentationFileName());
44
45   TTree *cTree = rl->GetTreeR("ITS", false);
46   if (cTree == 0)
47     return ;
48
49   TObjArray layerClus;
50   int nlr=0;
51   while(1) {
52     TBranch* br = cTree->GetBranch(Form("ITSRecPoints%d",nlr));
53     if (!br) break;
54     TClonesArray* clr = 0;
55     br->SetAddress(&clr);
56     layerClus.AddLast(clr);
57     nlr++;
58   }
59
60   TEveElementList* evClusters = new TEveElementList("ITS clusters");
61   TEvePointSet* layClusters = 0;
62   //  layClusters->SetOwnIds(kTRUE);
63
64   Int_t layOld =-1;
65   Int_t nentr = (Int_t) cTree->GetEntries();
66   cTree->GetEntry(0);
67   for (int ilr=0;ilr<nlr;ilr++) {
68     TClonesArray* clr = (TClonesArray*)layerClus.At(ilr);
69     int ncl = clr->GetEntries();
70
71     Float_t maxRsqr = maxR*maxR;
72     for (Int_t icl = 0; icl < ncl; ++icl) {
73       AliCluster *c = (AliCluster*) clr->UncheckedAt(icl);
74       Int_t mod = c->GetVolumeId();
75       //      int detType = gm->GetModuleDetTypeID(mod);
76       //      AliITSUSegmentationPix* segm = (AliITSUSegmentationPix*)segmArr.At(detType);
77       int lay,lad,det;
78       gm->GetModuleId(mod, lay,lad,det);
79       
80       if (lay!=layOld) { // assumes order in the digits !
81         layOld=lay;
82         layClusters = new TEvePointSet(Form("ITS clusters - Layer %d",lay),10000);
83         //  layClusters->ApplyVizTag(viz_tag, "Clusters");
84         layClusters->SetMarkerColor(kBlue);
85         layClusters->SetMarkerStyle(21);//kFullDotMedium);
86         layClusters->SetRnrSelf(kTRUE);
87         layClusters->SetOwnIds(kTRUE);
88         
89         evClusters->AddElement(layClusters);   
90       }
91       
92       Double_t l[3]={c->GetX(),c->GetY(),c->GetZ()}; 
93       Double_t g[3]; 
94       gm->LocalToGlobal(mod,l,g);
95       if (g[0]*g[0] + g[1]*g[1] < maxRsqr) {
96         layClusters->SetNextPoint(g[0], g[1], g[2]);
97         AliCluster *atp = new AliCluster(*c);
98         layClusters->SetPointId(atp);
99         //      printf("%d: mod %d: loc(%.4lf,%.4lf,%.4lf); glob(%.4lf,%.4lf,%.4lf); \n",
100         //             icl,atp->GetVolumeId(), l[0],l[1],l[2], g[0],g[1],g[2]);
101       }
102     }
103   }
104   
105   rl->UnloadRecPoints("ITS");
106   
107   //  clusters->SetTitle(Form("N=%d", clusters->Size()));
108   
109   const TString viz_tag("REC Clusters ITS");
110   
111   evClusters->ApplyVizTag(viz_tag, "Clusters");
112   
113   gEve->AddElement(evClusters, cont);
114   
115   gEve->Redraw3D();
116
117 }