minor fixes
[u/mrichter/AliRoot.git] / RICH / RichMenu.C
CommitLineData
e30ca504 1AliRun *a; AliStack *s; AliRunLoader *al; TGeoManager *g; //globals for easy manual manipulations
2AliRICH *r; AliLoader *rl; AliRICHParam *rp;
3Bool_t isGeomType=kFALSE;
8493d0aa 4
db910db9 5//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
e30ca504 6void GetParam()
722b211e 7{
e30ca504 8 isGeomType=!isGeomType;
9 if(g) delete g; if(rp) delete rp; //delete current TGeoManager and AliRICHParam
10 if(isGeomType) g=TGeoManager::Import("geometry.root");
11 else g=TGeoManager::Import("misaligned_geometry.root");
12 rp=AliRICHParam::Instance();
722b211e 13}
9a221675 14//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
473cb33a 15void RichMenu()
8493d0aa 16{
9a221675 17 TString status="Status: ";
18 if(gSystem->IsFileInIncludePath("galice.root")){
19 status+="galice.root found";
20 al=AliRunLoader::Open(); //try to open galice.root from current dir
21 if(gAlice) delete gAlice; //in case we execute this in aliroot delete default AliRun object
22 al->LoadgAlice(); a=al->GetAliRun(); //take new AliRun object from galice.root
23 rl=al->GetDetectorLoader("RICH"); r=(AliRICH*)a->GetDetector("RICH"); //get RICH object from galice.root
24
25 status+=Form(" with %i event(s)",al->GetNumberOfEvents()); status+=(r)? " with RICH": " without RICH";
26 }else
27 status+="No galice.root";
28
e30ca504 29 GetParam();
30
9a221675 31 TControlBar *pMenu = new TControlBar("horizontal",status.Data(),0,0);
32 pMenu->AddButton(" ","","");
e30ca504 33 pMenu->AddButton(" General " ,"General()" ,"general items which do not depend on any files");
34 pMenu->AddButton(" ","" ,"");
35 pMenu->AddButton(" Sim data " ,"SimData()" ,"items which expect to have simulated files" );
36 pMenu->AddButton(" ","" ,"");
37 pMenu->AddButton(" Raw data " ,"RawData()" ,"items which expect to have raw files" );
38 pMenu->AddButton(" ","print()" ,"");
39 pMenu->AddButton("Test" ,"Test()" ,"all test utilities");
40 pMenu->AddButton(" ","GetParam()" ,"");
41 pMenu->AddButton("Quit" ,".q" ,"close session" );
9a221675 42 pMenu->Show();
43}
44//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
45void General()
46{
47 TControlBar *pMenu = new TControlBar("vertical","Sim data",60,50);
48 pMenu->AddButton("Debug ON" ,"don();" ,"Switch debug on-off" );
49 pMenu->AddButton("Debug OFF" ,"doff();" ,"Switch debug on-off" );
9a221675 50 pMenu->AddButton("Geo GUI" ,"geo();" ,"Shows geometry" );
51 pMenu->AddButton("Browser" ,"new TBrowser;" ,"Start ROOT TBrowser" );
52 pMenu->Show();
53}//menu()
54//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
55void SimData()
56{
57 TControlBar *pMenu = new TControlBar("vertical","Sim",190,50);
e30ca504 58 pMenu->AddButton("Display ALL chambers" ,"ed();" , "Display Fast");
db910db9 59 pMenu->AddButton("HITS QA" ,"hqa()" ,"QA plots for hits: hqa()");
e30ca504 60
61 pMenu->AddButton("Print hits" ,"hp();" ,"To print hits: hp(evt)");
62 pMenu->AddButton("Print sdigits" ,"sp();" ,"To print sdigits: sp(evt)");
63 pMenu->AddButton("Print digits" ,"dp();" ,"To print digits: dp(evt)");
64 pMenu->AddButton("Print clusters" ,"cp();" ,"To print clusters: cp(evt)");
65
9a221675 66 pMenu->Show();
67}
68//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
69void RawData()
70{
71 TControlBar *pMenu = new TControlBar("vertical","Raw",350,50);
72 pMenu->AddButton("ESD print" ,"ep();" ,"To print ESD info: ep()" );
73 pMenu->AddButton("ESD QA" ,"eq();" ,"To draw ESD hists: eq()" );
74 pMenu->AddButton("Clusters print" ,"cp();" ,"To print clusters: cp()" );
75 pMenu->AddButton("Clusters QA" ,"cq();" ,"To draw clusters hists: cq()" );
76 pMenu->AddButton("Print Matrix" ,"mp();" ,"To print prob matrix: mp()" );
77 pMenu->AddButton("Print occupancy" ,"r->OccupancyPrint(-1);" ,"To print occupancy" );
78 pMenu->AddButton("Print event summary " ,"r->SummaryOfEvent();" ,"To print a summary of the event" );
79 pMenu->Show();
e30ca504 80}//RawData
81//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
82void Test()
83{
84 TControlBar *pMenu = new TControlBar("vertical","Test",400,50);
85 pMenu->AddButton("Hits->Digits" ,"thd();" ,"test hits->sdigits->digits" );
86 pMenu->AddButton("Segmentation" ,"ts()" ,"test segmentation methods" );
87 pMenu->AddButton("Test response" ,"AliRICHParam::TestResp();" ,"Test AliRICHParam response methods" );
88 pMenu->AddButton("Test transformation","AliRICHParam::TestTrans();","Test AliRICHParam transformation methods" );
89 pMenu->AddButton("Test Recon" ,"rec();" ,"Test AliRICHRecon" );
90 pMenu->Show();
91}//Test()
9a221675 92//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
122efc54 93
db910db9 94
9a221675 95void doff(){ Printf("DebugOFF"); AliLog::SetGlobalDebugLevel(0);}
96void don() { Printf("DebugON"); AliLog::SetGlobalDebugLevel(AliLog::kDebug);}
db910db9 97
e30ca504 98void geo ( ) {gGeoManager->GetTopVolume()->Draw("ogl");}
9a221675 99
e30ca504 100void du ( ) {r->Dump ( );} //utility display
101
102void hp (Int_t evt=0 ) {r->HitPrint (evt);} //print hits for requested event
103void hq ( ) {r->HitQA ( );} //hits QA plots for all events
104void sp (Int_t evt=0 ) {r->SdiPrint (evt);} //print sdigits for requested event
105void sq (Int_t evt=0 ) {r->SdiPrint (evt);} //print sdigits for requested event
106void dp (Int_t evt=0 ) {r->DigPrint (evt);} //print digits for requested event
107void dq ( ) {AliRICHReconstructor::DigQA (al );} //digits QA plots for all events
108void cp (Int_t evt=0 ) {r->CluPrint (evt); } //print clusters for requested event
109void cq ( ) {AliRICHReconstructor::CluQA (al );} //clusters QA plots for all events
122efc54 110
e30ca504 111void ep ( ) {AliRICHTracker::EsdQA(1); }
112void eq ( ) {AliRICHTracker::EsdQA(); }
113void mp (Double_t probCut=0.7 ) {AliRICHTracker::MatrixPrint(probCut);}
122efc54 114
115
db910db9 116void t (Int_t evt=0 ) {AliRICHParam::Stack(evt);}
117void tid (Int_t tid,Int_t evt=0) {AliRICHParam::Stack(evt,tid);}
9a221675 118
122efc54 119
db910db9 120Int_t nem (Int_t evt=0) {AliRICHParam::StackCount(kElectron ,evt);} //utility number of electrons
121Int_t nep (Int_t evt=0) {AliRICHParam::StackCount(kPositron ,evt);} //utility number of positrons
122Int_t nmup(Int_t evt=0) {AliRICHParam::StackCount(kMuonPlus ,evt);} //utility number of positive muons
123Int_t nmum(Int_t evt=0) {AliRICHParam::StackCount(kMuonMinus ,evt);} //utility number of negative muons
124Int_t npi0(Int_t evt=0) {AliRICHParam::StackCount(kPi0 ,evt);} //utility number of neutral pions
125Int_t npip(Int_t evt=0) {AliRICHParam::StackCount(kPiPlus ,evt);} //utility number of positive pions
126Int_t npim(Int_t evt=0) {AliRICHParam::StackCount(kPiMinus ,evt);} //utility number of negative pions
127Int_t nk0 (Int_t evt=0) {AliRICHParam::StackCount(kK0 ,evt);} //utility number of neutral kaons
128Int_t nkp (Int_t evt=0) {AliRICHParam::StackCount(kKPlus ,evt);} //utility number of positive kaons
129Int_t nkm (Int_t evt=0) {AliRICHParam::StackCount(kKMinus ,evt);} //utility number of negative kaons
130Int_t npp (Int_t evt=0) {AliRICHParam::StackCount(kProton ,evt);} //utility number of protons
131Int_t npm (Int_t evt=0) {AliRICHParam::StackCount(kProtonBar ,evt);} //utility number of antiprotons
e30ca504 132//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
133
134void tst()
135{
136 Int_t ch=3;
137 TClonesArray hitLst("AliRICHHit");
138 TClonesArray sdiLst("AliRICHDigit");
139 TObjArray digLst; for(Int_t i=0;i<7;i++) digLst.AddAt(new TClonesArray("AliRICHDigit"),i);
140 TClonesArray cluLst("AliRICHCluster");
141
142 Int_t iHitCnt=0; new(hitLst[iHitCnt++]) AliRICHHit(ch,200e-9,kProton,1, 8.40 , 60.14); //c,e,pid,tid,xl,yl
143// new(hitLst[iHitCnt++]) AliRICHHit(ch,e*1e-9,kProton,2,x,y); //c,e,pid,tid,xl,yl,x,y,z
144
145 Printf("HIT------HIT---------HIT--------HIT------HIT------HIT");hitLst.Print();Printf("");
146 AliRICHv1::Hit2Sdi(&hitLst,&sdiLst); Printf("SDI------DIG---------SDI--------SDI------SDI------SDI");sdiLst.Print();Printf("");
147 AliRICHDigitizer::Sdi2Dig(&sdiLst,&digLst); Printf("DIG------DIG---------DIG--------DIG------DIG------DIG");digLst.Print();Printf("");
148 AliRICHReconstructor::Dig2Clu((TClonesArray*)digLst[ch],&cluLst); Printf("CLU------CLU---------CLU--------CLU------CLU------CLU");cluLst.Print();Printf("");
149
150
151}
152
153
154void print()
155{
156
157 Double_t r2d=TMath::RadToDeg();
158
159 Double_t x=AliRICHDigit::SizeAllX(),y=AliRICHDigit::SizeAllY();
160
161 TVector3 c6lt=rp->Lors2Mars(6,0,y); TVector3 c6rt=rp->Lors2Mars(6,x,y);
162 TVector3 c6ce=rp->Lors2Mars(6,x/2,y/2);
163 TVector3 c6lb=rp->Lors2Mars(6,0,0); TVector3 c6rb=rp->Lors2Mars(6,x,0);
164//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
165 TVector3 c5lt=rp->Lors2Mars(5,0,y); TVector3 c5rt=rp->Lors2Mars(5,x,y);
166 TVector3 c5ce=rp->Lors2Mars(5,x/2,y/2);
167 TVector3 c5lb=rp->Lors2Mars(5,0,0); TVector3 c5rb=rp->Lors2Mars(5,x,0);
168//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
169 TVector3 c4lt=rp->Lors2Mars(4,0,y); TVector3 c4rt=rp->Lors2Mars(4,x,y);
170 TVector3 c4ce=rp->Lors2Mars(4,x/2,y/2);
171 TVector3 c4lb=rp->Lors2Mars(4,0,0); TVector3 c4rb=rp->Lors2Mars(4,x,0);
172//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
173 TVector3 c3lt=rp->Lors2Mars(3,0,y); TVector3 c3rt=rp->Lors2Mars(3,x,y);
174 TVector3 c3ce=rp->Lors2Mars(3,x/2,y/2);
175 TVector3 c3lb=rp->Lors2Mars(3,0,0); TVector3 c3rb=rp->Lors2Mars(3,x,0);
176//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
177 TVector3 c2lt=rp->Lors2Mars(2,0,y); TVector3 c2rt=rp->Lors2Mars(2,x,y);
178 TVector3 c2ce=rp->Lors2Mars(2,x/2,y/2);
179 TVector3 c2lb=rp->Lors2Mars(2,0,0); TVector3 c2rb=rp->Lors2Mars(2,x,0);
180//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
181 TVector3 c1lt=rp->Lors2Mars(1,0,y); TVector3 c1rt=rp->Lors2Mars(1,x,y);
182 TVector3 c1ce=rp->Lors2Mars(1,x/2,y/2);
183 TVector3 c1lb=rp->Lors2Mars(1,0,0); TVector3 c1rb=rp->Lors2Mars(1,x,0);
184//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
185 TVector3 c0lt=rp->Lors2Mars(0,0,y); TVector3 c0rt=rp->Lors2Mars(0,x,y);
186 TVector3 c0ce=rp->Lors2Mars(0,x/2,y/2);
187 TVector3 c0lb=rp->Lors2Mars(0,0,0); TVector3 c0rb=rp->Lors2Mars(0,x,0);
188
189
190 Printf("\n\n\n");
191
192 Printf("_______________________________ _______________________________");
193 Printf("|%7.2f %7.2f| |%7.2f %7.2f|",c6lt.Mag() ,c6rt.Mag() ,c5lt.Mag() ,c5rt.Mag() );
194 Printf("|%7.2f %7.2f| |%7.2f %7.2f|",c6lt.Theta()*r2d ,c6rt.Theta()*r2d ,c5lt.Theta()*r2d ,c5rt.Theta()*r2d );
195 Printf("|%7.2f %7.2f| |%7.2f %7.2f|",c6lt.Phi()*r2d ,c6rt.Phi()*r2d ,c5lt.Phi()*r2d ,c5rt.Phi()*r2d );
196 Printf("| | | |");
197 Printf("| %7.2f | | %7.2f | Sensitive area (130.60,126.16)" ,c6ce.Mag() ,c5ce.Mag() );
198 Printf("| %7.2f | | %7.2f | Lors Center ( 65.30, 63.08)" ,c6ce.Theta()*r2d ,c5ce.Theta()*r2d );
199 Printf("| %7.2f | | %7.2f |" ,c6ce.Phi()*r2d ,c5ce.Phi()*r2d );
200 Printf("| | | |");
201 Printf("|%7.2f %7.2f| |%7.2f %7.2f|",c6lb.Mag() ,c6rb.Mag() ,c5lb.Mag() ,c5rb.Mag() );
202 Printf("|%7.2f %7.2f| |%7.2f %7.2f|",c6lb.Theta()*r2d ,c6rb.Theta()*r2d ,c5lb.Theta()*r2d ,c5rb.Theta()*r2d );
203 Printf("|%7.2f %7.2f| |%7.2f %7.2f|",c6lb.Phi()*r2d ,c6rb.Phi()*r2d ,c5lb.Phi()*r2d ,c5rb.Phi()*r2d );
204 Printf("------------------------------- -------------------------------");
205 Printf("");
206 Printf("_______________________________ _______________________________ _______________________________");
207 Printf("|%7.2f %7.2f| |%7.2f %7.2f| |%7.2f %7.2f|",c4lt.Mag() ,c4rt.Mag() ,c3lt.Mag() ,c3rt.Mag() ,c2lt.Mag() ,c2rt.Mag() );
208 Printf("|%7.2f %7.2f| |%7.2f %7.2f| |%7.2f %7.2f|",c4lt.Theta()*r2d,c4rt.Theta()*r2d,c3lt.Theta()*r2d,c3rt.Theta()*r2d,c2lt.Theta()*r2d ,c2rt.Theta()*r2d);
209 Printf("|%7.2f %7.2f| |%7.2f %7.2f| |%7.2f %7.2f|",c4lt.Phi()*r2d ,c4rt.Phi()*r2d ,c3lt.Phi()*r2d ,c3rt.Phi()*r2d ,c2lt.Phi()*r2d ,c2rt.Phi()*r2d );
210 Printf("| | | | | |");
211 Printf("| %7.2f | | %7.2f | | %7.2f |" ,c4ce.Mag() ,c3ce.Mag() ,c2ce.Mag());
212 Printf("| %7.2f | | %7.2f | | %7.2f |" ,c4ce.Theta()*r2d,c3ce.Theta()*r2d,c2ce.Theta()*r2d);
213 Printf("| %7.2f | | %7.2f | | %7.2f |" ,c4ce.Phi()*r2d ,c3ce.Phi()*r2d ,c3ce.Phi()*r2d);
214 Printf("| | | | | |");
215 Printf("|%7.2f %7.2f| |%7.2f %7.2f| |%7.2f %7.2f|",c4lb.Mag() ,c4rb.Mag() ,c3lb.Mag() ,c3rb.Mag() ,c2lt.Mag() ,c2rt.Mag() );
216 Printf("|%7.2f %7.2f| |%7.2f %7.2f| |%7.2f %7.2f|",c4lb.Theta()*r2d,c4rb.Theta()*r2d,c3lb.Theta()*r2d,c3rb.Theta()*r2d,c2lt.Theta()*r2d ,c2rt.Theta()*r2d);
217 Printf("|%7.2f %7.2f| |%7.2f %7.2f| |%7.2f %7.2f|",c4lb.Phi()*r2d ,c4rb.Phi()*r2d ,c3lb.Phi()*r2d ,c3rb.Phi()*r2d ,c2lt.Phi()*r2d ,c2rt.Phi()*r2d );
218 Printf("------------------------------- ------------------------------- -------------------------------");
219 Printf("");
220 Printf(" _______________________________ _______________________________");
221 Printf(" |%7.2f %7.2f| |%7.2f %7.2f|",c1lt.Mag() ,c1rt.Mag() ,c0lt.Mag() ,c0rt.Mag() );
222 Printf(" |%7.2f %7.2f| |%7.2f %7.2f|",c1lt.Theta()*r2d,c1rt.Theta()*r2d,c0lt.Theta()*r2d,c0rt.Theta()*r2d);
223 Printf(" |%7.2f %7.2f| |%7.2f %7.2f|",c1lt.Phi()*r2d ,c1rt.Phi()*r2d ,c0lt.Phi()*r2d ,c0rt.Phi()*r2d );
224 Printf(" | | | |");
225 Printf(" | %7.2f | | %7.2f |" ,c1ce.Mag() ,c0ce.Mag() );
226 Printf(" | %7.2f | | %7.2f |" ,c1ce.Theta()*r2d,c0ce.Theta()*r2d);
227 Printf(" | %7.2f | | %7.2f |" ,c1ce.Phi()*r2d ,c0ce.Phi()*r2d );
228 Printf(" | | | |");
229 Printf(" |%7.2f %7.2f| |%7.2f %7.2f|",c1lb.Mag() ,c1rb.Mag() ,c0lb.Mag() ,c0rb.Mag() );
230 Printf(" |%7.2f %7.2f| |%7.2f %7.2f|",c1lb.Theta()*r2d,c1rb.Theta()*r2d,c0lb.Theta()*r2d,c0rb.Theta()*r2d );
231 Printf(" |%7.2f %7.2f| |%7.2f %7.2f|",c1lb.Phi()*r2d ,c1rb.Phi()*r2d ,c0lb.Phi()*r2d ,c0rb.Phi()*r2d );
232 Printf(" ------------------------------- -------------------------------");
233
234
235 Double_t m[3];
236 for(int i=0;i<1000;i++){
237 Float_t xout=0,xin=gRandom->Rndm()*130.60;
238 Float_t yout=0,yin=gRandom->Rndm()*126.16;
239 Int_t c=gRandom->Rndm()*6;
240 rp->Lors2Mars(c,xin,yin,m);
241 rp->Mars2Lors(c,m,xout,yout);
242 if( (xin-xout) != 0) Printf("Problem in X");
243 if( (yin-yout) != 0) Printf("Problem in Y");
244 }
245}//print()
246
247
248void rec()
249{
250 Double_t th=8*TMath::DegToRad();//gRandom->Rndm()*TMath::PiOver4();
251 Double_t ph=gRandom->Rndm()*TMath::TwoPi(); 172.51*TMath::DegToRad();
252 Double_t x=gRandom->Rndm()*AliRICHDigit::SizeAllX(); //101.59;
253 Double_t y=gRandom->Rndm()*2*AliRICHDigit::SizePcY();//38.06;
254
255
05668f8d 256 Double_t ckovMax=0.75,ckovSim;
e30ca504 257 Int_t nSim=0;
258 while(nSim<3){
259 ckovSim=gRandom->Rndm()*ckovMax;//0.6468;
260 nSim=20*TMath::Power(TMath::Sin(ckovSim)/TMath::Sin(ckovMax),2); //scale number of photons
261 }
262
263
264 TClonesArray *pCluLst=new TClonesArray("AliRICHCluster");
265 TPolyMarker *pMipMap=new TPolyMarker(); pMipMap->SetMarkerStyle(8); pMipMap->SetMarkerColor(kRed);
266 TPolyMarker *pPhoMap=new TPolyMarker(); pPhoMap->SetMarkerStyle(4); pPhoMap->SetMarkerColor(kRed);
267 TPolyMarker *pBkgMap=new TPolyMarker(); pBkgMap->SetMarkerStyle(25); pBkgMap->SetMarkerColor(kRed);
268 TPolyLine *pRing =new TPolyLine; pRing->SetLineColor(kGreen);
269 TPolyLine *pOld =new TPolyLine; pOld->SetLineColor(kBlue);
270
271 Int_t iCluCnt=0; pMipMap->SetPoint(iCluCnt,x,y); new((*pCluLst)[iCluCnt++]) AliRICHCluster(1,x,y,200); //add mip cluster
272
273// for(int j=0;j<30;j++){ //add bkg photons
274// Float_t bkgX=gRandom->Rndm()*AliRICHDigit::SizeAllX();
275// Float_t bkgY=gRandom->Rndm()*AliRICHDigit::SizeAllY();
276// pBkgMap->SetPoint(iCluCnt,bkgX,bkgY); new((*pCluLst)[iCluCnt++]) AliRICHCluster(1,bkgX,bkgY,35);
277// }
278
279
280
281
282 AliRICHRecon rec; rec.SetTrack(th,ph,x,y);
283 AliRICHReconOld old; old.SetTrack(th,ph,x,y);
284
285 TVector2 pos;
286 for(int i=0;i<nSim;i++){
05668f8d 287 rec.TracePhot(ckovSim,gRandom->Rndm()*2*TMath::Pi(),pos); //add photons
e30ca504 288 if(AliRICHDigit::IsInDead(pos.X(),pos.Y())) continue;
289 pPhoMap->SetPoint(iCluCnt,pos.X(),pos.Y()); new((*pCluLst)[iCluCnt++]) AliRICHCluster(1,pos.X(),pos.Y(),35);
290 }
291
292 for(int i=0;i<pCluLst->GetEntries();i++){
293 AliRICHCluster *pClu=(AliRICHCluster*)pCluLst->At(i);
294 Printf("diff phi %f ckov new %f ckov old %f",rec.FindPhotPhi(pClu->X(),pClu->Y())-old.FindPhotPhi(pClu->X(),pClu->Y()),rec.FindPhotCkov(pClu->X(),pClu->Y()),old.FindPhotCkov(pClu->X(),pClu->Y()));
295 }
296
297 Int_t nRec=0,nOld=0;
298 Double_t ckovRec=rec.CkovAngle(pCluLst,nRec); Double_t err=TMath::Sqrt(rec.CkovSigma2());
299 Double_t ckovOld=old.CkovAngle(pCluLst,nOld);
300
301 Printf("---------------- Now reconstructed --------------------");
302
303
05668f8d 304 for(int j=0;j<100;j++){rec.TracePhot(ckovRec,j*0.0628,pos); pRing->SetPoint(j,pos.X(),pos.Y());}
305 for(int j=0;j<100;j++){rec.TracePhot(ckovOld,j*0.0628,pos); pOld->SetPoint(j,pos.X(),pos.Y());}
e30ca504 306
307 new TCanvas; AliRICHDigit::DrawPc(); pMipMap->Draw(); pPhoMap->Draw(); pBkgMap->Draw(); pRing->Draw(); pOld->Draw();
308
309 TLatex txt; txt.SetTextSize(0.03);
310 txt.DrawLatex(65,127,Form("#theta=%.4f#pm%.5f with %2i #check{C}" ,ckovSim, 0.,nSim ));
311 txt.DrawLatex(65,122,Form("#theta=%.4f#pm%.5f with %2i #check{C} Old=%.4f with %i #check{C}" ,ckovRec,err,nRec,ckovOld,nOld));
312 txt.DrawLatex(0 ,127,Form("#theta=%.2f#circ #phi=%.2f#circ @(%.2f,%.2f) ",th*TMath::RadToDeg(),ph*TMath::RadToDeg(),x,y));
313
314// for(int i=0;i<35;i++){
315// Double_t ckov=0.1+i*0.02;
316// Printf("Ckov=%.2f Old=%.3f New=%.3f",ckov,old.FindRingArea(ckov),rec.FindRingArea(ckov));
317// }
318
319
320}//rec()
321
322AliRICHRecon rec; AliRICHReconOld old;
323
324void aaa()
325{
326 rec.SetTrack(10*TMath::DegToRad(),1,30,60); old.SetTrack(10*TMath::DegToRad(),1,30,60);
327 TVector2 pos;
328 Double_t ckovSim=0.6234;
05668f8d 329 rec.TracePhot(ckovSim,1,pos);
e30ca504 330 Printf("ckovSim %f",ckovSim);
331 double ckovRec=rec.FindPhotCkov(pos.X(),pos.Y());
332 double ckovOld=old.FindPhotCkov(pos.X(),pos.Y());
333 Printf("new %f old %f",ckovRec,ckovOld);
334}
335
336
337
338
339
e30ca504 340
341
342void AliRICH::HitQA(Double_t cut,Double_t cutele,Double_t cutR)
343{
344// Provides a set of control plots intended primarily for charged particle flux analisys
345// Arguments: cut (GeV) - cut on momentum of any charged particles but electrons,
346// cetele (GeV) - the same for electrons-positrons
347// cutR (cm) - cut on production vertex radius (cylindrical system)
348 gBenchmark->Start("HitsAna");
349
350 Double_t cutPantiproton =cut;
351 Double_t cutPkaonminus =cut;
352 Double_t cutPpionminus =cut;
353 Double_t cutPmuonminus =cut;
354 Double_t cutPpositron =cutele;
355
356 Double_t cutPelectron =cutele;
357 Double_t cutPmuonplus =cut;
358 Double_t cutPpionplus =cut;
359 Double_t cutPkaonplus =cut;
360 Double_t cutPproton =cut;
361
362
363 TH2F *pEleHitRZ =new TH2F("EleHitRZ" ,Form("e^{+} e^{-} hit %s;z[cm];R[cm]" ,GetName()) , 400,-300,300 ,400,-500,500); //R-z plot 0cm<R<550cm -300cm<z<300cm
364 TH2F *pEleHitRP =new TH2F("EleHitRP" ,Form("e^{+} e^{-} hit %s;p[GeV];R[cm]",GetName()) ,1000,-1 ,1 ,400, 0,550); //R-p plot 0cm<R<550cm -1GeV<p<1GeV
365 TH1F *pEleAllP =new TH1F("EleAllP" , "e^{+} e^{-} all;p[GeV]" ,1000,-1 ,1 );
366 TH1F *pEleHitP =new TH1F("EleHitP" ,Form("e^{+} e^{-} hit %s;p[GeV]" ,GetName()) ,1000,-1 ,1 );
367 TH1F *pMuoHitP =new TH1F("MuoHitP" ,Form("#mu^{-} #mu^{+} hit %s;p[GeV]" ,GetName()) ,1000,-4 ,4 );
368 TH1F *pPioHitP =new TH1F("PioHitP" ,Form("#pi^{-} #pi^{+} hit %s;p[GeV]" ,GetName()) ,1000,-4 ,4 );
369 TH1F *pKaoHitP =new TH1F("KaoHitP" ,Form("K^{-} K^{+} hit %s;p[GeV]" ,GetName()) ,1000,-4 ,4 );
370 TH1F *pProHitP =new TH1F("ProHitP" ,Form("p^{-} p^{+} hit %s;p[GeV]" ,GetName()) ,1000,-4 ,4 );
371 TH2F *pFlux =new TH2F("flux" ,Form("%s flux with Rvertex<%.1fcm" ,GetName(),cutR),10 ,-5 ,5 , 10,0 ,10); //special text hist
372 TH2F *pVertex =new TH2F("vertex" ,Form("%s 2D vertex of RICH hit;x;y" ,GetName()) ,120 ,0 ,600 ,120,0 ,600); //special text hist
373 TH1F *pRho =new TH1F("rho" ,Form("%s r of RICH hit" ,GetName()) ,600 ,0 ,600); //special text hist
374 pFlux->SetStats(0);
375 pFlux->GetXaxis()->SetBinLabel(1 ,Form("p^{-}>%.3fGeV/c" ,cutPantiproton));
376 pFlux->GetXaxis()->SetBinLabel(2 ,Form("K^{-}>%.3fGeV/c" ,cutPkaonminus ));
377 pFlux->GetXaxis()->SetBinLabel(3 ,Form("#pi^{-}>%.3fGeV/c" ,cutPpionminus ));
378 pFlux->GetXaxis()->SetBinLabel(4 ,Form("#mu^{-}>%.3fGeV/c" ,cutPmuonminus ));
379 pFlux->GetXaxis()->SetBinLabel(5 ,Form("e^{+}>%.3fGeV/c" ,cutPpositron ));
380
381 pFlux->GetXaxis()->SetBinLabel(6 ,Form("e^{-}>%.3fGeV/c" ,cutPelectron ));
382 pFlux->GetXaxis()->SetBinLabel(7 ,Form("#mu^{+}>%.3fGeV/c" ,cutPmuonplus ));
383 pFlux->GetXaxis()->SetBinLabel(8 ,Form("#pi^{+}>%.3fGeV/c" ,cutPpionplus ));
384 pFlux->GetXaxis()->SetBinLabel(9 ,Form("K^{+}>%.3fGeV/c" ,cutPkaonplus ));
385 pFlux->GetXaxis()->SetBinLabel(10,Form("p^{+}>%.3fGeV/c" ,cutPproton ));
386
387 pFlux->GetYaxis()->SetBinLabel(1,"sum");
388 pFlux->GetYaxis()->SetBinLabel(2,"ch1");
389 pFlux->GetYaxis()->SetBinLabel(3,"ch2");
390 pFlux->GetYaxis()->SetBinLabel(4,"ch3");
391 pFlux->GetYaxis()->SetBinLabel(5,"ch4");
392 pFlux->GetYaxis()->SetBinLabel(6,"ch5");
393 pFlux->GetYaxis()->SetBinLabel(7,"ch6");
394 pFlux->GetYaxis()->SetBinLabel(8,"ch7");
395 pFlux->GetYaxis()->SetBinLabel(9,"prim");
396 pFlux->GetYaxis()->SetBinLabel(10,"tot");
397
398//end of hists definition
399
400 Int_t iNevents=fLoader->GetRunLoader()->GetAliRun()->GetEventsPerRun(),iCntPrimParts=0,iCntTotParts=0;
401//load all needed trees
402 fLoader->LoadHits();
403 fLoader->GetRunLoader()->LoadHeader();
404 fLoader->GetRunLoader()->LoadKinematics();
405
406 for(Int_t iEvtN=0;iEvtN < iNevents;iEvtN++){//events loop
407 fLoader->GetRunLoader()->GetEvent(iEvtN);
408 AliInfo(Form(" %i event processes",fLoader->GetRunLoader()->GetEventNumber()));
409 AliStack *pStack= fLoader->GetRunLoader()->Stack();
410
411 for(Int_t iParticleN=0;iParticleN<pStack->GetNtrack();iParticleN++){//stack loop
412 TParticle *pPart=pStack->Particle(iParticleN);
413
414 if(iParticleN%10000==0) AliInfo(Form(" %i particles read",iParticleN));
415
416 switch(pPart->GetPdgCode()){
417 case kProtonBar: pFlux->Fill(-4.5,9); if(pPart->Rho()<0.01) pFlux->Fill(-4.5,8); break;
418 case kKMinus: pFlux->Fill(-3.5,9); if(pPart->Rho()<0.01) pFlux->Fill(-3.5,8); break;
419 case kPiMinus: pFlux->Fill(-2.5,9); if(pPart->Rho()<0.01) pFlux->Fill(-2.5,8); break;
420 case kMuonMinus: pFlux->Fill(-1.5,9); if(pPart->Rho()<0.01) pFlux->Fill(-1.5,8); break;
421 case kPositron: pFlux->Fill(-0.5,9); if(pPart->Rho()<0.01) pFlux->Fill(-0.5,8); pEleAllP->Fill(-pPart->P()); break;
422
423 case kElectron: pFlux->Fill( 0.5,9); if(pPart->Rho()<0.01) pFlux->Fill( 0.5,8); pEleAllP->Fill( pPart->P()); break;
424 case kMuonPlus: pFlux->Fill( 1.5,9); if(pPart->Rho()<0.01) pFlux->Fill( 1.5,8); break;
425 case kPiPlus: pFlux->Fill( 2.5,9); if(pPart->Rho()<0.01) pFlux->Fill( 2.5,8); break;
426 case kKPlus: pFlux->Fill( 3.5,9); if(pPart->Rho()<0.01) pFlux->Fill( 3.5,8); break;
427 case kProton: pFlux->Fill( 4.5,9); if(pPart->Rho()<0.01) pFlux->Fill( 4.5,8); break;
428 }//switch
429 }//stack loop
430//now hits analiser
431 for(Int_t iEntryN=0;iEntryN < fLoader->TreeH()->GetEntries();iEntryN++){//TreeH loop
432 fLoader->TreeH()->GetEntry(iEntryN); //get current entry (prim)
433 for(Int_t iHitN=0;iHitN < Hits()->GetEntries();iHitN++){//hits loop
434 AliRICHHit *pHit = (AliRICHHit*)Hits()->At(iHitN); //get current hit
435 TParticle *pPart=pStack->Particle(pHit->GetTrack()); //get stack particle which produced the current hit
436
437 if(pPart->GetPDG()->Charge()!=0&&pPart->Rho()>0.1) pVertex->Fill(pPart->Vx(),pPart->Vy()); //safe margin for sec.
438 if(pPart->GetPDG()->Charge()!=0) pRho->Fill(pPart->Rho()); //safe margin for sec.
439 if(pPart->R()>cutR) continue; //cut on production radius (cylindrical system)
440
441 switch(pPart->GetPdgCode()){
442 case kProtonBar: if(pPart->P()>cutPantiproton) {pProHitP->Fill(-pPart->P()); pFlux->Fill(-4.5,pHit->Ch());}break;
443 case kKMinus : if(pPart->P()>cutPkaonminus) {pKaoHitP->Fill(-pPart->P()); pFlux->Fill(-3.5,pHit->Ch());}break;
444 case kPiMinus : if(pPart->P()>cutPpionminus) {pPioHitP->Fill(-pPart->P()); pFlux->Fill(-2.5,pHit->Ch());}break;
445 case kMuonMinus: if(pPart->P()>cutPmuonminus) {pMuoHitP->Fill(-pPart->P()); pFlux->Fill(-1.5,pHit->Ch());}break;
446 case kPositron : if(pPart->P()>cutPpositron) {pEleHitP->Fill(-pPart->P()); pFlux->Fill(-0.5,pHit->Ch());
447 pEleHitRP->Fill(-pPart->P(),pPart->R()); pEleHitRZ->Fill(pPart->Vz(),pPart->R()); }break;
448
449 case kElectron : if(pPart->P()>cutPelectron) {pEleHitP->Fill( pPart->P()); pFlux->Fill( 0.5,pHit->Ch());
450 pEleHitRP->Fill( pPart->P(),pPart->R()); pEleHitRZ->Fill(pPart->Vz(),pPart->R()); }break;
451 case kMuonPlus : if(pPart->P()>cutPmuonplus) {pMuoHitP->Fill( pPart->P()); pFlux->Fill( 1.5,pHit->Ch());}break;
452 case kPiPlus : if(pPart->P()>cutPpionplus) {pPioHitP->Fill( pPart->P()); pFlux->Fill( 2.5,pHit->Ch());}break;
453 case kKPlus : if(pPart->P()>cutPkaonplus) {pKaoHitP->Fill( pPart->P()); pFlux->Fill( 3.5,pHit->Ch());}break;
454 case kProton : if(pPart->P()>cutPproton) {pProHitP->Fill( pPart->P()); pFlux->Fill( 4.5,pHit->Ch());}break;
455 }
456 }//hits loop
457 }//TreeH loop
458 iCntPrimParts +=pStack->GetNprimary();
459 iCntTotParts +=pStack->GetNtrack();
460 }//events loop
461//unload all loaded staff
462 fLoader->UnloadHits();
463 fLoader->GetRunLoader()->UnloadHeader();
464 fLoader->GetRunLoader()->UnloadKinematics();
465//Calculater some sums
466 Stat_t sum=0;
467//sum row, sum over rows
468 for(Int_t i=1;i<=pFlux->GetNbinsX();i++){
469 sum=0; for(Int_t j=2;j<=8;j++) sum+=pFlux->GetBinContent(i,j);
470 pFlux->SetBinContent(i,1,sum);
471 }
472
473//display everything
474 new TCanvas("canvas1",Form("Events %i Nprims=%i Nparticles=%i",iNevents,iCntPrimParts,iCntTotParts),1000,900); pFlux->Draw("text"); gPad->SetGrid();
475//total prims and particles
476 TLatex latex; latex.SetTextSize(0.02);
477 sum=0; for(Int_t i=1;i<=pFlux->GetNbinsX();i++) sum+=pFlux->GetBinContent(i,10); latex.DrawLatex(5.1,9.5,Form("%.0f",sum));
478 sum=0; for(Int_t i=1;i<=pFlux->GetNbinsX();i++) sum+=pFlux->GetBinContent(i,9); latex.DrawLatex(5.1,8.5,Form("%.0f",sum));
479 for(Int_t iCh=0;iCh<7;iCh++) {
480 sum=0; for(Int_t i=1;i<=pFlux->GetNbinsX();i++) sum+=pFlux->GetBinContent(i,iCh+2);latex.DrawLatex(5.1,iCh+1.5,Form("%.0f",sum));
481 }
482 sum=0; for(Int_t i=1;i<=pFlux->GetNbinsX();i++) sum+=pFlux->GetBinContent(i,1); latex.DrawLatex(5.1,0.5,Form("%.0f",sum));
483
484 new TCanvas("cEleAllP" ,"e" ,200,100); pEleAllP->Draw();
485 new TCanvas("cEleHitRP" ,"e" ,200,100); pEleHitRP->Draw();
486 new TCanvas("cEleHitRZ" ,"e" ,200,100); pEleHitRZ->Draw();
487 new TCanvas("cEleHitP" ,"e" ,200,100); pEleHitP->Draw();
488 new TCanvas("cMuoHitP" ,"mu",200,100); pMuoHitP->Draw();
489 new TCanvas("cPioHitP" ,"pi",200,100); pPioHitP->Draw();
490 new TCanvas("cKaoHitP" ,"K" ,200,100); pKaoHitP->Draw();
491 new TCanvas("cProHitP" ,"p" ,200,100); pProHitP->Draw();
492 new TCanvas("cVertex" ,"2d vertex" ,200,100); pVertex->Draw();
493 new TCanvas("cRho" ,"Rho of sec" ,200,100); pRho->Draw();
494
495 gBenchmark->Show("HitsPlots");
496}//HitQA()
497
498
499
500
501