1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15 #include "AliRICHDisplFast.h"
17 #include "AliRICHChamber.h"
18 #include "AliRICHParam.h"
19 #include <AliLoader.h>
21 #include <TPolyLine.h>
22 #include <TParticle.h>
27 ClassImp(AliRICHDisplFast)
29 //__________________________________________________________________________________________________
30 void AliRICHDisplFast::Exec()
32 TH2F *pHitsH2 = new TH2F("pHitsH2" , "Event Display",165,-AliRICHParam::PcSizeX()/2,AliRICHParam::PcSizeX()/2,
33 144,-AliRICHParam::PcSizeY()/2,AliRICHParam::PcSizeY()/2);
34 pHitsH2->SetStats(kFALSE);
35 TH2F *pDigitsH2 = new TH2F("pDigitsH2" ,"Event Display",165,-AliRICHParam::PcSizeX()/2,AliRICHParam::PcSizeX()/2,
36 144,-AliRICHParam::PcSizeY()/2,AliRICHParam::PcSizeY()/2);
37 TH2F *pClustersH2 = new TH2F("pClustersH2","Event Display",165,-AliRICHParam::PcSizeX()/2,AliRICHParam::PcSizeX()/2,
38 144,-AliRICHParam::PcSizeY()/2,AliRICHParam::PcSizeY()/2);
40 TCanvas *Display = new TCanvas("Display","RICH Display",0,0,600,600);
42 gStyle->SetPalette(1);
44 AliRICH *pRich = (AliRICH*)gAlice->GetDetector("RICH");
45 pRich->GetLoader()->LoadRecPoints();
46 pRich->GetLoader()->LoadDigits();
47 pRich->GetLoader()->LoadHits();
48 for(Int_t iEventN=0;iEventN<gAlice->GetEventsPerRun();iEventN++){//events Loop
49 pRich->GetLoader()->GetRunLoader()->GetEvent(iEventN);
50 pRich->GetLoader()->TreeD()->GetEntry(0);
51 pRich->GetLoader()->TreeR()->GetEntry(0);
54 Int_t nPrimaries = (Int_t)pRich->GetLoader()->TreeH()->GetEntries();
55 TObjArray * Hits = new TObjArray[nPrimaries];
57 for(Int_t i=0;i<nPrimaries;i++) {
58 pRich->GetLoader()->TreeH()->GetEntry(i);
59 Int_t nHits = pRich->Hits()->GetEntries();
60 for(Int_t k=0;k<nHits;k++) Hits[i].Add(pRich->Hits()->At(k));
64 for(Int_t iChamber=1;iChamber<=7;iChamber++){//modules loop
66 Int_t nDigits = pRich->Digits(iChamber)->GetEntries();
67 Int_t nClusters = pRich->Clusters(iChamber)->GetEntries();
69 pHitsH2->Reset(); pDigitsH2->Reset(); pClustersH2->Reset();
73 for(Int_t i=0;i<nPrimaries;i++){//prims loop
74 pRich->GetLoader()->TreeH()->GetEntry(i);
75 Int_t nHits = pRich->Hits()->GetEntries();
76 for(Int_t j=0;j<nHits;j++){//hits loop
77 AliRICHhit *pHit = (AliRICHhit*)Hits[i].At(j);
78 if(pHit->C()==iChamber){
79 TVector3 xyzhit(pHit->X(),pHit->Y(),pHit->Z());
80 TVector3 hitlocal = pRich->C(iChamber)->Glob2Loc(xyzhit);
81 pHitsH2->Fill(hitlocal.X(),hitlocal.Y(),200);
86 for(Int_t j=0;j<nDigits;j++){//digits loop
87 AliRICHdigit *pDigit = (AliRICHdigit*)pRich->Digits(iChamber)->At(j);
88 AliRICHParam::Pad2Loc(pDigit->X(),pDigit->Y(),xpad,ypad);
89 pDigitsH2->Fill(xpad,ypad,100);
92 for(Int_t j=0;j<nClusters;j++){//clusters loop
93 AliRICHcluster *pCluster = (AliRICHcluster*)pRich->Clusters(iChamber)->At(j);
94 pClustersH2->Fill(pCluster->X(),pCluster->Y(),50);
97 pHitsH2->SetTitle(Form("event %i module %2i",iEventN,iChamber));
98 pHitsH2->SetMarkerColor(kRed); pHitsH2->SetMarkerStyle(29); pHitsH2->SetMarkerSize(0.4);
105 pDigitsH2->SetMarkerColor(kGreen); pDigitsH2->SetMarkerStyle(29); pDigitsH2->SetMarkerSize(0.4);
106 pDigitsH2->Draw("same");
111 pClustersH2->SetMarkerColor(kBlue); pClustersH2->SetMarkerStyle(29); pClustersH2->SetMarkerSize(0.4);
112 pClustersH2->Draw("same");
119 pRich->GetLoader()->UnloadRecPoints();
120 pRich->GetLoader()->UnloadDigits();
121 pRich->GetLoader()->UnloadHits();
123 //__________________________________________________________________________________________________
124 void AliRICHDisplFast::DrawSectors()
126 Double_t x1[5] = {-AliRICHParam::PcSizeX()/2,
127 -AliRICHParam::SectorSizeX()/2-AliRICHParam::DeadZone(),
128 -AliRICHParam::SectorSizeX()/2-AliRICHParam::DeadZone(),
129 -AliRICHParam::PcSizeX()/2,
130 -AliRICHParam::PcSizeX()/2};
131 Double_t y1[5] = {AliRICHParam::DeadZone()/2,
132 AliRICHParam::DeadZone()/2,
133 AliRICHParam::PcSizeY()/2,
134 AliRICHParam::PcSizeY()/2,
135 AliRICHParam::DeadZone()/2};
136 Double_t x2[5] = {-AliRICHParam::SectorSizeX()/2,
137 AliRICHParam::SectorSizeX()/2,
138 AliRICHParam::SectorSizeX()/2,
139 -AliRICHParam::SectorSizeX()/2,
140 -AliRICHParam::SectorSizeX()/2};
141 Double_t y2[5] = {AliRICHParam::DeadZone()/2,
142 AliRICHParam::DeadZone()/2,
143 AliRICHParam::PcSizeY()/2,
144 AliRICHParam::PcSizeY()/2,
145 AliRICHParam::DeadZone()/2};
146 Double_t x3[5] = { AliRICHParam::SectorSizeX()/2+AliRICHParam::DeadZone(),
147 AliRICHParam::PcSizeX()/2,
148 AliRICHParam::PcSizeX()/2,
149 AliRICHParam::SectorSizeX()/2+AliRICHParam::DeadZone(),
150 AliRICHParam::SectorSizeX()/2+AliRICHParam::DeadZone()};
151 Double_t y3[5] = {AliRICHParam::DeadZone()/2,
152 AliRICHParam::DeadZone()/2,
153 AliRICHParam::PcSizeY()/2,
154 AliRICHParam::PcSizeY()/2,
155 AliRICHParam::DeadZone()/2};
156 Double_t x4[5] = {-AliRICHParam::PcSizeX()/2,
157 -AliRICHParam::SectorSizeX()/2-AliRICHParam::DeadZone(),
158 -AliRICHParam::SectorSizeX()/2-AliRICHParam::DeadZone(),
159 -AliRICHParam::PcSizeX()/2,
160 -AliRICHParam::PcSizeX()/2};
161 Double_t y4[5] = {-AliRICHParam::PcSizeY()/2,
162 -AliRICHParam::PcSizeY()/2,
163 -AliRICHParam::DeadZone()/2,
164 -AliRICHParam::DeadZone()/2,
165 -AliRICHParam::PcSizeY()/2};
166 Double_t x5[5] = {-AliRICHParam::SectorSizeX()/2,
167 AliRICHParam::SectorSizeX()/2,
168 AliRICHParam::SectorSizeX()/2,
169 -AliRICHParam::SectorSizeX()/2,
170 -AliRICHParam::SectorSizeX()/2};
171 Double_t y5[5] = {-AliRICHParam::PcSizeY()/2,
172 -AliRICHParam::PcSizeY()/2,
173 -AliRICHParam::DeadZone()/2,
174 -AliRICHParam::DeadZone()/2,
175 -AliRICHParam::PcSizeY()/2};
176 Double_t x6[5] = { AliRICHParam::SectorSizeX()/2+AliRICHParam::DeadZone(),
177 AliRICHParam::PcSizeX()/2,
178 AliRICHParam::PcSizeX()/2,
179 AliRICHParam::SectorSizeX()/2+AliRICHParam::DeadZone(),
180 AliRICHParam::SectorSizeX()/2+AliRICHParam::DeadZone()};
181 Double_t y6[5] = {-AliRICHParam::PcSizeY()/2,
182 -AliRICHParam::PcSizeY()/2,
183 -AliRICHParam::DeadZone()/2,
184 -AliRICHParam::DeadZone()/2,
185 -AliRICHParam::PcSizeY()/2};
186 TPolyLine *sector1 = new TPolyLine(5,x1,y1);
187 TPolyLine *sector2 = new TPolyLine(5,x2,y2);
188 TPolyLine *sector3 = new TPolyLine(5,x3,y3);
189 TPolyLine *sector4 = new TPolyLine(5,x4,y4);
190 TPolyLine *sector5 = new TPolyLine(5,x5,y5);
191 TPolyLine *sector6 = new TPolyLine(5,x6,y6);
192 sector1->SetLineColor(21);
193 sector2->SetLineColor(21);
194 sector3->SetLineColor(21);
195 sector4->SetLineColor(21);
196 sector5->SetLineColor(21);
197 sector6->SetLineColor(21);