b171449cb530a0e99f675173ce96c334b0d064c1
[u/mrichter/AliRoot.git] / RICH / AliRICHDisplFast.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3  *                                                                        *
4  * Author: The ALICE Off-line Project.                                    *
5  * Contributors are mentioned in the code where appropriate.              *
6  *                                                                        *
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"
16 #include "AliRICH.h"
17 #include "AliRICHChamber.h"
18 #include "AliRICHParam.h"
19 #include <AliLoader.h>
20 #include <TCanvas.h>
21 #include <TPolyLine.h>
22 #include <TParticle.h>
23 #include <TStyle.h>
24 #include <TH2.h>
25 #include <TMath.h>
26
27 ClassImp(AliRICHDisplFast)
28
29 //__________________________________________________________________________________________________
30 void AliRICHDisplFast::Exec()
31 {
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);
39
40   TCanvas *Display = new TCanvas("Display","RICH Display",0,0,600,600);
41     
42   gStyle->SetPalette(1);
43
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);
52
53   
54     Int_t nPrimaries = (Int_t)pRich->GetLoader()->TreeH()->GetEntries();
55     TObjArray * Hits = new TObjArray[nPrimaries];
56   
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));
61     
62     }
63            
64     for(Int_t iChamber=1;iChamber<=7;iChamber++){//modules loop
65     
66      Int_t nDigits   = pRich->Digits(iChamber)->GetEntries();
67      Int_t nClusters = pRich->Clusters(iChamber)->GetEntries();
68    
69      pHitsH2->Reset();     pDigitsH2->Reset();     pClustersH2->Reset();
70
71
72       for(Int_t i=0;i<nPrimaries;i++){//prims loop
73         pRich->GetLoader()->TreeH()->GetEntry(i);
74         Int_t nHits = pRich->Hits()->GetEntries();
75         for(Int_t j=0;j<nHits;j++){//hits loop
76           AliRICHhit *pHit = (AliRICHhit*)Hits[i].At(j);
77           if(pHit->C()==iChamber){
78             TVector3 xyzhit(pHit->X(),pHit->Y(),pHit->Z());
79             TVector2 hitlocal = pRich->C(iChamber)->Glob2Loc(xyzhit);
80             pHitsH2->Fill(hitlocal.X(),hitlocal.Y(),200);
81           }//if
82         }//hits loop         
83       }//prims loop
84      
85       for(Int_t j=0;j<nDigits;j++){//digits loop
86         AliRICHdigit *pDigit = (AliRICHdigit*)pRich->Digits(iChamber)->At(j);
87         TVector2 x2=AliRICHParam::Pad2Loc(pDigit->X(),pDigit->Y());
88         pDigitsH2->Fill(x2.X(),x2.Y(),100);
89       }//digits loop
90         
91       for(Int_t j=0;j<nClusters;j++){//clusters loop
92         AliRICHcluster *pCluster = (AliRICHcluster*)pRich->Clusters(iChamber)->At(j);
93         pClustersH2->Fill(pCluster->X(),pCluster->Y(),50);
94       }//clusters loop
95
96       pHitsH2->SetTitle(Form("event %i module %2i",iEventN,iChamber));
97       pHitsH2->SetMarkerColor(kRed); pHitsH2->SetMarkerStyle(29); pHitsH2->SetMarkerSize(0.4);
98       pHitsH2->Draw();
99       DrawSectors();
100       Display->Update();
101       Display->Modified();
102       getchar();
103              
104       pDigitsH2->SetMarkerColor(kGreen); pDigitsH2->SetMarkerStyle(29); pDigitsH2->SetMarkerSize(0.4);
105       pDigitsH2->Draw("same");
106       Display->Update();
107       Display->Modified();       
108       getchar();
109       
110       pClustersH2->SetMarkerColor(kBlue); pClustersH2->SetMarkerStyle(29);  pClustersH2->SetMarkerSize(0.4);
111       pClustersH2->Draw("same");
112       Display->Update();
113       Display->Modified();
114       getchar();
115      }//modules loop
116     delete [] Hits;
117   }////events Loop
118   pRich->GetLoader()->UnloadRecPoints();
119   pRich->GetLoader()->UnloadDigits();
120   pRich->GetLoader()->UnloadHits();
121 }//Exec()
122 //__________________________________________________________________________________________________
123 void AliRICHDisplFast::DrawSectors() 
124
125   Double_t x1[5] = {-AliRICHParam::PcSizeX()/2,
126                     -AliRICHParam::SectorSizeX()/2-AliRICHParam::DeadZone(),
127                     -AliRICHParam::SectorSizeX()/2-AliRICHParam::DeadZone(),
128                     -AliRICHParam::PcSizeX()/2,
129                     -AliRICHParam::PcSizeX()/2};
130   Double_t y1[5] = {AliRICHParam::DeadZone()/2,
131                     AliRICHParam::DeadZone()/2,
132                     AliRICHParam::PcSizeY()/2,
133                     AliRICHParam::PcSizeY()/2,
134                     AliRICHParam::DeadZone()/2};  
135   Double_t x2[5] = {-AliRICHParam::SectorSizeX()/2,
136                      AliRICHParam::SectorSizeX()/2,
137                      AliRICHParam::SectorSizeX()/2,
138                     -AliRICHParam::SectorSizeX()/2,
139                     -AliRICHParam::SectorSizeX()/2};
140   Double_t y2[5] = {AliRICHParam::DeadZone()/2,
141                     AliRICHParam::DeadZone()/2,
142                     AliRICHParam::PcSizeY()/2,
143                     AliRICHParam::PcSizeY()/2,
144                     AliRICHParam::DeadZone()/2};
145   Double_t x3[5] = { AliRICHParam::SectorSizeX()/2+AliRICHParam::DeadZone(),
146                      AliRICHParam::PcSizeX()/2,
147                      AliRICHParam::PcSizeX()/2,
148                      AliRICHParam::SectorSizeX()/2+AliRICHParam::DeadZone(),
149                      AliRICHParam::SectorSizeX()/2+AliRICHParam::DeadZone()};
150   Double_t y3[5] = {AliRICHParam::DeadZone()/2,
151                     AliRICHParam::DeadZone()/2,
152                     AliRICHParam::PcSizeY()/2,
153                     AliRICHParam::PcSizeY()/2,
154                     AliRICHParam::DeadZone()/2};
155   Double_t x4[5] = {-AliRICHParam::PcSizeX()/2,
156                     -AliRICHParam::SectorSizeX()/2-AliRICHParam::DeadZone(),
157                     -AliRICHParam::SectorSizeX()/2-AliRICHParam::DeadZone(),
158                     -AliRICHParam::PcSizeX()/2,
159                     -AliRICHParam::PcSizeX()/2};
160   Double_t y4[5] = {-AliRICHParam::PcSizeY()/2,
161                     -AliRICHParam::PcSizeY()/2,
162                     -AliRICHParam::DeadZone()/2,
163                     -AliRICHParam::DeadZone()/2,
164                     -AliRICHParam::PcSizeY()/2};
165    Double_t x5[5] = {-AliRICHParam::SectorSizeX()/2,
166                      AliRICHParam::SectorSizeX()/2,
167                      AliRICHParam::SectorSizeX()/2,
168                     -AliRICHParam::SectorSizeX()/2,
169                     -AliRICHParam::SectorSizeX()/2};
170   Double_t y5[5] = {-AliRICHParam::PcSizeY()/2,
171                     -AliRICHParam::PcSizeY()/2,
172                     -AliRICHParam::DeadZone()/2,
173                     -AliRICHParam::DeadZone()/2,
174                     -AliRICHParam::PcSizeY()/2};
175   Double_t x6[5] = { AliRICHParam::SectorSizeX()/2+AliRICHParam::DeadZone(),
176                      AliRICHParam::PcSizeX()/2,
177                      AliRICHParam::PcSizeX()/2,
178                      AliRICHParam::SectorSizeX()/2+AliRICHParam::DeadZone(),
179                      AliRICHParam::SectorSizeX()/2+AliRICHParam::DeadZone()};
180   Double_t y6[5] = {-AliRICHParam::PcSizeY()/2,
181                     -AliRICHParam::PcSizeY()/2,
182                     -AliRICHParam::DeadZone()/2,
183                     -AliRICHParam::DeadZone()/2,
184                     -AliRICHParam::PcSizeY()/2};
185   TPolyLine *sector1 = new TPolyLine(5,x1,y1);
186   TPolyLine *sector2 = new TPolyLine(5,x2,y2);
187   TPolyLine *sector3 = new TPolyLine(5,x3,y3);
188   TPolyLine *sector4 = new TPolyLine(5,x4,y4);
189   TPolyLine *sector5 = new TPolyLine(5,x5,y5);
190   TPolyLine *sector6 = new TPolyLine(5,x6,y6);
191   sector1->SetLineColor(21);
192   sector2->SetLineColor(21);
193   sector3->SetLineColor(21);
194   sector4->SetLineColor(21);
195   sector5->SetLineColor(21);
196   sector6->SetLineColor(21);
197   sector1->Draw();
198   sector2->Draw();
199   sector3->Draw();
200   sector4->Draw();
201   sector5->Draw();
202   sector6->Draw();
203 }