From Philippe & Laurent: new variant of MUON visualization.
[u/mrichter/AliRoot.git] / EVE / alice-macros / muon_digits.C
1 // $Id$
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  **************************************************************************/
8
9 // Macro to visualise digits from MUON spectrometer 
10 // (both tracker and trigger).
11 //
12 // Use muon_digits() in order to run it
13 //
14 // Needs that alieve_init() is already called
15
16 #if !defined(__CINT__) || defined(__MAKECINT__)
17
18 #include "AliMUONGeometryTransformer.h"
19 #include "AliMUONVDigit.h"
20 #include "AliMUONVDigitStore.h"
21
22 #include "AliMpPad.h"
23 #include "AliMpSegmentation.h"
24 #include "AliMpVSegmentation.h"
25 #include "AliMpCDB.h"
26
27 #include "AliRunLoader.h"
28
29 #include "EveBase/AliEveEventManager.h"
30
31 #include <TEveManager.h>
32 #include <TEveQuadSet.h>
33
34 #include <TTree.h>
35 #include <TStyle.h>
36 #include <Riostream.h>
37
38 #endif
39
40 AliMUONGeometryTransformer* gMUONGeometryTransformer(0x0);
41
42 //______________________________________________________________________________
43 void add_muon_digits(TIter* next, TEveQuadSet* bending, TEveQuadSet* nonBending, Bool_t fromRaw)
44 {
45   // load mapping
46   AliMpCDB::LoadAll(kFALSE);
47   
48   // load geometry
49   if (gMUONGeometryTransformer == 0) 
50   {
51     AliEveEventManager::AssertGeometry();
52     gMUONGeometryTransformer = new AliMUONGeometryTransformer();
53     gMUONGeometryTransformer->LoadGeometryData();
54   }
55   
56   // loop over digits and produce corresponding graphic objects
57   AliMUONVDigit* digit;
58   while ( ( digit = static_cast<AliMUONVDigit*>((*next)() ) ) )
59   {
60     if (!digit->IsTrigger() && !fromRaw && digit->Charge() < 1.e-3) continue;
61     
62     Int_t detElemId = digit->DetElemId();
63     Int_t manuId = digit->ManuId();
64     
65     const AliMpVSegmentation* vseg =
66       AliMpSegmentation::Instance()->GetMpSegmentation(detElemId, AliMp::GetCathodType(digit->Cathode()));
67     if (!vseg) 
68     {
69       cout << Form("Could not get segmentation for DE %4d MANU %4d",detElemId,manuId) << endl;
70       continue; // should not happen, unless we got a readout error and thus a bad de,manu pair
71     }
72     
73     AliMpPad pad = vseg->PadByLocation(manuId,digit->ManuChannel());
74     
75     Double_t local[] = { pad.GetPositionX(), pad.GetPositionY(), 0.0 };
76     Double_t global[] = { 0.0, 0.0, 0.0 };
77     
78     gMUONGeometryTransformer->Local2Global(detElemId,
79                                            local[0], local[1], local[2],
80                                            global[0], global[1], global[2]);
81     
82     TEveQuadSet* pads = bending;
83     if (vseg->PlaneType()==AliMp::kNonBendingPlane) pads = nonBending;
84     
85     pads->AddQuad(global[0]-pad.GetDimensionX(),global[1]-pad.GetDimensionY(),global[2],
86                   2.*pad.GetDimensionX(),2.*pad.GetDimensionY());
87     
88     if (fromRaw && !digit->IsTrigger()) pads->QuadValue(digit->ADC());
89     else pads->QuadValue(digit->Charge());
90   }
91   
92 }
93
94 //______________________________________________________________________________
95 void muon_digits()
96 {
97   // load digits
98   AliRunLoader* rl =  AliEveEventManager::AssertRunLoader();
99   rl->LoadDigits("MUON");
100   TTree* dt = rl->GetTreeD("MUON", kFALSE);
101   if (!dt) return;
102   AliMUONVDigitStore *digitStore = AliMUONVDigitStore::Create(*dt);
103   digitStore->Clear();
104   digitStore->Connect(*dt,0);
105   dt->GetEvent(0);
106   if (digitStore->GetSize() == 0 && !gEve->GetKeepEmptyCont()) return;
107   
108   // container for graphic representation of digits
109   TEveElementList* cont = new TEveElementList("MUON Digits");
110   
111   TEveQuadSet* bending = new TEveQuadSet(TEveQuadSet::kQT_RectangleXY, kFALSE, 32);
112   bending->SetName("Bending");
113   bending->SetRenderMode(TEveDigitSet::kRM_Fill);
114   bending->SetPickable(kFALSE);
115   cont->AddElement(bending);
116   
117   TEveQuadSet* nonBending = new TEveQuadSet(TEveQuadSet::kQT_RectangleXY, kFALSE, 32);
118   nonBending->SetName("Non bending");
119   nonBending->SetRenderMode(TEveDigitSet::kRM_Line);
120   nonBending->SetPickable(kFALSE);
121   cont->AddElement(nonBending);
122   
123   // add digits to the containers
124   TIter next(digitStore->CreateIterator());
125   add_muon_digits(&next, bending, nonBending, kFALSE);
126   delete digitStore;
127   
128   // set containers' title
129   Int_t nDigitB = bending->GetPlex()->Size();
130   Int_t nDigitNB = nonBending->GetPlex()->Size();
131   cont->SetTitle(Form("N=%d",nDigitB+nDigitNB));
132   bending->SetTitle(Form("N=%d",nDigitB));
133   nonBending->SetTitle(Form("N=%d",nDigitNB));
134   
135   // automatic scaling
136   gStyle->SetPalette(1);
137   bending->AssertPalette();
138   nonBending->AssertPalette();
139   
140   // add graphic containers
141   gEve->DisableRedraw();
142   gEve->AddElement(cont);
143   gEve->EnableRedraw();
144   gEve->Redraw3D();
145 }