New method Invert() for changing alpha by pi (forbiden operation via Rotate())
[u/mrichter/AliRoot.git] / EVE / alice-macros / hmpid_digits.C
CommitLineData
ba978640 1#if !defined(__CINT__) || defined(__MAKECINT__)
2#include <TClonesArray.h>
3#include <TBranch.h>
4#include <TGeoMatrix.h>
5#include <TTree.h>
6#include <TStyle.h>
7#include <TEveElement.h>
8#include <TEveFrameBox.h>
9#include <TEveManager.h>
10#include <TEvePointSet.h>
11#include <TEveRGBAPalette.h>
12#include <TEveTrans.h>
13#include <TEveQuadSet.h>
14
15#include <HMPID/AliHMPIDDigit.h>
16#include <HMPID/AliHMPIDv3.h>
17#include <STEER/STEER/AliCluster3D.h>
18#include <STEER/STEER/AliRunLoader.h>
19#include <EVE/EveBase/AliEveEventManager.h>
20#endif
21
391fa967 22void hmpid_digits()
23{
24 const Char_t *name[]={ "HMPID0", "HMPID1", "HMPID2", "HMPID3",
25 "HMPID4", "HMPID5", "HMPID6" };
26
27 AliRunLoader* rl = AliEveEventManager::AssertRunLoader();
28 rl->LoadDigits("HMPID");
29
30 TTree *dTree = rl->GetTreeD("HMPID", kFALSE);
ba978640 31 if (!dTree) return;
391fa967 32
33 TEveElementList* list = new TEveElementList("HMPID Digits");
34 gEve->AddElement(list);
35
36 gStyle->SetPalette(1, 0);
37
38 TEveRGBAPalette *pal = new TEveRGBAPalette(0, 3000);
39 pal->SetMax(1000);
40 TEveFrameBox *box = new TEveFrameBox();
41 box->SetAAQuadXY(0, 0, 0, 144, 121);
42 box->SetFrameColor(kGray);
43
44 TClonesArray* digits = new TClonesArray("AliHMPIDDigit");
45 for (Int_t iCh = 0; iCh < 7; ++iCh)
46 {
47 TBranch *br = dTree->GetBranch(name[iCh]);
48 br->SetAddress(&digits);
49 br->GetEntry(0);
50
51 TEveQuadSet* q = new TEveQuadSet(Form("Chamber %d", iCh));
52 q->SetOwnIds(kTRUE);
53 q->SetPalette(pal);
54 q->SetFrame(box);
55 q->SetAntiFlick(kTRUE);
56 q->SetPickable(kTRUE);
57
58 q->Reset(TEveQuadSet::kQT_RectangleXYFixedDimZ, kFALSE, 64);
59 q->SetDefCoord(0);
60 q->SetDefHeight(0.84f);
61 q->SetDefWidth(0.8f);
62
63 for(Int_t iDig = 0; iDig < digits->GetEntriesFast(); ++iDig)
64 {
65 AliHMPIDDigit *pDig = (AliHMPIDDigit*) digits->At(iDig);
66
67 q->AddQuad(pDig->PadChX()*0.8f, pDig->PadChY()*0.84f);
68 q->QuadValue(TMath::Nint(pDig->Q()));
69 q->QuadId(new AliHMPIDDigit(*pDig));
70 }
71
72 q->RefitPlex();
73
74 TGeoHMatrix mat;
75 AliHMPIDv3::IdealPosition(iCh, &mat);
76 q->RefMainTrans().SetFrom(mat);
77 q->RefMainTrans().Move3LF(-0.5f*144, -0.5f*121, 0);
78
79 list->AddElement(q);
80 }
81
82 delete digits;
83 rl->UnloadDigits("HMPID");
84
85 gEve->Redraw3D();
86}