better ESD performance plot
[u/mrichter/AliRoot.git] / TRD / qaRec / macros / AliTrackletsinTRD.C
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
16 ////////////////////////////////////////////////////////////////////////////
17 //                                                                        //
18 //  Reconstruction QA                                                     //
19 //                                                                        //
20 //  Authors: Yvonne C. Pachmayer <pachmay@physi.uni-heidelberg.de>        //
21 //                                                                        //
22 ////////////////////////////////////////////////////////////////////////////
23
24
25 #define AliTrackletsinTRD_cxx
26 #include "AliTrackletsinTRD.h"
27 #include <TH2.h>
28 #include <TStyle.h>
29 #include <TCanvas.h>
30
31 void AliTrackletsinTRD::Loop()
32 {
33     //   In a ROOT session, you can do:
34     //      Root > .L AliTrackletsinTRD.C
35     //      Root > AliTrackletsinTRD t
36     //      Root > t.GetEntry(12); // Fill t data members with entry number 12
37     //      Root > t.Show();       // Show values of entry 12
38     //      Root > t.Show(16);     // Read and show values of entry 16
39     //      Root > t.Loop();       // Loop on all entries
40     //
41
42     //     This is the loop skeleton where:
43     //    jentry is the global entry number in the chain
44     //    ientry is the entry number in the current Tree
45     //  Note that the argument to GetEntry must be:
46     //    jentry for TChain::GetEntry
47     //    ientry for TTree::GetEntry and TBranch::GetEntry
48     //
49     //       To read only selected branches, Insert statements like:
50     // METHOD1:
51     //    fChain->SetBranchStatus("*",0);  // disable all branches
52     //    fChain->SetBranchStatus("branchname",1);  // activate branchname
53     // METHOD2: replace line
54     //    fChain->GetEntry(jentry);       //read all branches
55     //by  b_branchname->GetEntry(ientry); //read only this branch
56
57     gROOT->SetStyle("Plain");
58     gStyle->SetPalette(1);
59
60
61     if (fChain == 0) return;
62
63     // open run loader and load gAlice, kinematics and header
64     AliRunLoader* runLoader = AliRunLoader::Open("galice.root");
65     if (!runLoader)
66     {
67         printf("Error: readKine:\nGetting run loader from file \"%s\" failed", filename);
68         return;
69     }
70
71     runLoader->LoadHeader();
72     runLoader->LoadRecPoints("TRD");
73     TObjArray *module = new TObjArray();
74
75     runLoader->CdGAFile();
76
77     Int_t nEvents = runLoader->GetNumberOfEvents();
78     const Int_t nEventsarray=nEvents;
79     AliTRDcluster *cls = 0;
80
81     Double_t x_clus[540][90];
82     Double_t y_clus[540][90];
83
84     for(Int_t ev = 0; ev < nEvents - 1; ev++)
85     {
86         TTree *tree = runLoader->GetTreeR("TRD", 0);
87         tree->SetBranchAddress("TRDcluster", &module);
88
89         int N = tree->GetEntries(); // number of chamber, max 540
90         // Check number of clusters
91         for(Int_t ind = 0; ind < N; ind++)
92         {
93             tree->GetEntry(ind);
94             Int_t m = module->GetEntries();
95
96             for (Int_t j = 0; j < m; j++)    // loop over clusters of one chamber
97             {
98                 if(j>=90) break;
99                 if (cls != 0) delete cls;
100                 cls = (AliTRDcluster*)module->At(j);
101                 x_clus[ind][j]=cls->GetX();
102                 y_clus[ind][j]=cls->GetY();
103
104             }
105         }
106
107
108
109         // loop over debug file and analysis
110
111         Float_t xbmin = 0., xbmax = 380.;
112         /* // histogram for debugging purpose
113         Int_t nxbins = 3000, nybins = 280;
114         Float_t ybmin = -70., ybmax = 70.;
115         hxy = new TH2F("hxy",";x;y",nxbins,xbmin,xbmax,nybins,ybmin,ybmax);
116         */
117
118         Long64_t nentries = fChain->GetEntriesFast();
119
120         Long64_t nbytes = 0, nb = 0;
121         Float_t xpos[6];
122         Float_t ypos[6];
123         Int_t counter_test;
124         for (Long64_t jentry=0; jentry<nentries;jentry++) {
125             Long64_t ientry = LoadTree(jentry);
126             if (ientry < 0) break;
127             nb = fChain->GetEntry(jentry);   nbytes += nb;
128             xpos[flayer]=fxtrack;
129             ypos[flayer]=fytrack;
130
131
132             if(flayer==5)
133             {
134                 Float_t x[6]= {xpos[0],xpos[1],xpos[2],xpos[3],xpos[4],xpos[5]};
135                 Float_t y[6]= {ypos[0],ypos[1],ypos[2],ypos[3],ypos[4],ypos[5]};
136
137                 gxy = new TGraph(6,x,y);
138
139                 TF1 *f1 = new TF1("f1","pol1",xbmin,xbmax);
140                 gxy->Fit("f1","RQ");
141                 // resulting function: y=ax+b --> x*f1->GetParameter(1)+f1->GetParameter(0)
142
143                 // graph and fit plotting only for checking
144                 /*
145                  c1 = new TCanvas("c1","hough Transform",0,0,800,800);
146                  c1->cd(1);
147                  hxy->Draw();
148                  gxy->SetMarkerStyle(3);
149                  gxy->Draw("p");
150                  */
151
152                 Float_t px;
153                 Float_t py;
154
155                 if(feventcounter==ev)
156                 {
157                     for(Int_t b=0;b<540;b++)
158                     {
159                         if(b==fdettracklet)
160                         {
161                             Int_t counter_distcalc;
162                             Float_t distance_sum=0;
163
164                             for(Int_t c=0;c<90;c++)
165                             {
166                                 px=0;
167                                 py=0;
168
169                                 px=x_clus[b][c];
170                                 py=y_clus[b][c];
171
172                                 if(px!=0 &&  py!=0)
173                                 {
174                                     // Function to calculate the distance of a point to the above fitted function
175                                     Double_t distance = 0;
176                                     if((TMath::Sqrt(f1->GetParameter(1)*f1->GetParameter(1)+1))!=0) distance=TMath::Abs(f1->GetParameter(1)*px-py+f1->GetParameter(0))/(TMath::Sqrt(f1->GetParameter(1)*f1->GetParameter(1)+1));
177                                     // if(distance<10) cout << eventcounter << " " << b << " " << c << " " <<   distance << endl;
178                                     if(distance>0.6 && distance<=2)
179                                     {
180                                         counter_distcalc++;
181                                         distance_sum+=distance;
182                                     }
183                                 }
184
185                             }
186                             if(counter_distcalc!=0)
187                             {
188                                 // these are tracks which have additional clusters close by: 0.6 < distance <= 2
189                             }
190                             else
191                             {
192                                 // these are good tracks no additional clusters close by
193                                cout <<" good tracks " <<  endl;
194                             }
195                             counter_distcalc=0;
196                         }
197                     }
198                 }
199
200                 if(gxy) delete gxy;
201                 if(f1) delete f1;
202             } // end of 5 layer loop
203
204
205            
206
207
208         }  // end loop over debug file and analysis
209
210         runLoader->GetNextEvent();
211     }
212
213 }
214