]> git.uio.no Git - u/mrichter/AliRoot.git/blob - TPC/fastSimul/testEdgeFit.C
Added speed test macro
[u/mrichter/AliRoot.git] / TPC / fastSimul / testEdgeFit.C
1 /*
2   .L $ALICE_ROOT/TPC/fastSimul/testEdgeFit.C+
3   testEdge(50000);
4 */
5
6 #include "TVectorD.h"
7 #include "TMath.h"
8 #include "TLinearFitter.h"
9 #include "TRandom.h"
10 #include "TH2F.h"
11 #include "TFile.h"
12 #include "TTree.h"
13 #include "TLegend.h"
14 #include "TPad.h"
15
16 #include "TTreeStream.h"
17 //#include "TVectorD.h"
18
19
20 void testEdge(Int_t ntest){
21
22   TTreeSRedirector *cstream = new TTreeSRedirector("testEdge.root");
23   Double_t yth=1;
24   Double_t ymc0=0;
25   Double_t kx=0;
26   Double_t sy=0.25;
27   Int_t npoints=100;
28   TVectorD vecX;
29   TVectorD vecY0;
30   TVectorD vecYC;
31   TVectorD vecYT;
32
33   TVectorD vecU;
34   TVectorD vecYTa;
35   TVectorD vecYTb;
36   TVectorD vecYEa;
37   TVectorD vecYEb;
38   //
39   TVectorD vecYTM5;
40   TVectorD vecYTk5;
41   //
42   TVectorD vecE;
43   //
44   //
45   for (Int_t it=0; it<ntest;it++){
46     //
47     ymc0 = gRandom->Rndm()*1;
48     kx = gRandom->Rndm()*1.;
49     npoints = TMath::Nint(gRandom->Rndm()*50);
50     vecX.ResizeTo(npoints);
51     vecY0.ResizeTo(npoints);
52     vecYC.ResizeTo(npoints);
53     vecYT.ResizeTo(npoints);
54     //
55     vecYTa.ResizeTo(npoints);
56     vecYTb.ResizeTo(npoints);
57     vecYEa.ResizeTo(npoints);
58     vecYEb.ResizeTo(npoints);
59     vecU.ResizeTo(npoints);
60     //
61     TLinearFitter fitter(2,"pol1");
62     for (Int_t ipoint=npoints-1;ipoint>0;ipoint--){
63       Double_t x[10];
64       //
65       x[0]=(Double_t)(ipoint)*0.2;
66       Double_t y0   = ymc0+Double_t(x[0])*kx; 
67       Double_t yc   = y0+gRandom->Gaus()*sy;
68       //
69       if (yc<yth)   yc=yth;
70       vecX[ipoint]   = x[0];;
71       vecY0[ipoint]  = y0;
72       vecYC[ipoint]  = yc;
73       vecU[ipoint]   = 0;
74       vecYT[ipoint]  = -10;
75       vecYTa[ipoint] = -10;
76       vecYTb[ipoint] = -10;
77       fitter.AddPoint(x,yc,sy);
78       if (fitter.GetNpoints()>4){
79         fitter.Eval();
80         fitter.GetErrors(vecE);
81         vecYTa[ipoint]=fitter.GetParameter(0);
82         vecYTb[ipoint]=fitter.GetParameter(1);
83         vecYEa[ipoint]=vecE[0];
84         vecYEb[ipoint]=vecE[1];
85         vecU[ipoint]=fitter.GetNpoints();
86         vecYT[ipoint]=fitter.GetParameter(0)+fitter.GetParameter(1)*x[0];
87       }
88     }
89     (*cstream)<<"Dump"<<
90       "n="<<npoints<<
91       "ymc="<<ymc0<<
92       "kx="<<kx<<
93       "vX.="<<&vecX<<
94       "vU.="<<&vecU<<
95       "vY0.="<<&vecY0<<
96       "vYC.="<<&vecYC<<
97       "vYT.="<<&vecYT<<
98       "vYTax.="<<&vecYTa<<
99       "vYTb.="<<&vecYTb<<
100       "vYEa.="<<&vecYEa<<
101       "vYEb.="<<&vecYEb<<
102       "\n";
103   }
104   delete cstream;
105 }
106
107
108 void MakePicDy(){
109   TFile f("testEdge.root");
110   TTree * tree = (TTree*)f.Get("Dump");
111   TObjArray arrfdy0,arrfdyT;
112   TH2F *phisdY0 = new TH2F("hisdY0","hisdY0",25,1,3,200,-3,3);
113   TH2F *phisdYT = new TH2F("hisdYT","hisdYT",25,1,3,200,-3,3);
114   tree->Draw("10*(vYT.fElements-vY0.fElements):vY0.fElements>>hisdY0","vU.fElements>20&&kx>0.1","colz");
115   phisdY0->SetDirectory(0);
116   tree->Draw("10*(vYT.fElements-vY0.fElements):vYT.fElements>>hisdYT","vU.fElements>20&&kx>0.1","colz");
117   phisdYT->SetDirectory(0);
118
119   phisdY0->FitSlicesY(0,0,-1,0,"QNR",&arrfdy0); 
120   phisdYT->FitSlicesY(0,0,-1,0,"QNR",&arrfdyT);
121   //
122   TH1 *hisdy0 = (TH1*)arrfdy0.At(1);
123   TH1 *hisdyT = (TH1*)arrfdyT.At(1);
124   hisdy0->SetLineColor(1);
125   hisdyT->SetLineColor(2);
126   hisdy0->SetMarkerColor(1);
127   hisdyT->SetMarkerColor(2);
128   hisdy0->SetMarkerStyle(22);
129   hisdyT->SetMarkerStyle(24);
130   hisdyT->SetXTitle("Y (cm)");
131   hisdyT->SetYTitle("Y_{t}-Y_{mc} (mm)");
132   hisdyT->SetMinimum(-hisdyT->GetMaximum());
133   hisdyT->Draw();
134   hisdy0->Draw("same");
135   TLegend *legend = new TLegend(0.25,0.60,0.85,0.85, "Track residuals close to edge (Edge= 1 cm, #sigma_{cl}=0.25 cm)");
136   legend->AddEntry(hisdy0,"Y=Y_{mc}");
137   legend->AddEntry(hisdyT,"Y=Y_{t}");
138   legend->Draw();
139   gPad->SaveAs("picEdge/dYedgeMC.eps");
140   gPad->SaveAs("picEdge/dYedgeMC.gif");
141 }
142
143 void MakePickDy(){
144   TFile f("testEdge.root");
145   TTree * tree = (TTree*)f.Get("Dump");
146   TObjArray arrfdky0,arrfdkyT;
147   TH2F *phisdKY0 = new TH2F("hisdKY0","hisdKY0",25,1,3,200,-150,150);
148   TH2F *phisdKYT = new TH2F("hisdKYT","hisdKYT",25,1,3,200,-150,150);
149   tree->Draw("1000*(vYTb.fElements-kx):vY0.fElements>>hisdKY0","vU.fElements>20","colz");
150   phisdKY0->SetDirectory(0);
151   tree->Draw("1000*(vYTb.fElements-kx):vYT.fElements>>hisdKYT","vU.fElements>20","colz");
152   phisdKYT->SetDirectory(0);
153
154   phisdKY0->FitSlicesY(0,0,-1,0,"QNR",&arrfdky0); 
155   phisdKYT->FitSlicesY(0,0,-1,0,"QNR",&arrfdkyT);
156   //
157   TH1 *hisdky0 = (TH1*)arrfdky0.At(1);
158   TH1 *hisdkyT = (TH1*)arrfdkyT.At(1);
159   hisdky0->SetLineColor(1);
160   hisdkyT->SetLineColor(2);
161   hisdky0->SetMarkerColor(1);
162   hisdkyT->SetMarkerColor(2);
163   hisdky0->SetMarkerStyle(22);
164   hisdkyT->SetMarkerStyle(24);
165   hisdkyT->SetMinimum(-20);
166   hisdkyT->SetMaximum(20);
167   hisdkyT->SetXTitle("Y (cm)");
168   hisdkyT->SetYTitle("k_{yt}-k_{ymc} (mrad)");
169   hisdkyT->Draw();
170   hisdky0->Draw("same");
171   TLegend *legend = new TLegend(0.25,0.60,0.85,0.85, "Track residuals close to edge (Edge= 1 cm, #sigma_{cl}=0.25 cm)");
172   legend->AddEntry(hisdky0,"Y=Y_{mc}");
173   legend->AddEntry(hisdkyT,"Y=Y_{t}");
174   legend->Draw();
175   gPad->SaveAs("picEdge/dKYedgeMC.eps");
176   gPad->SaveAs("picEdge/dKYedgeMC.gif");
177 }
178
179