Monitoring the eta-phi plots
[u/mrichter/AliRoot.git] / ITS / AliITSSortTrkl.cxx
CommitLineData
108cb8b9 1/**************************************************************************
2 * Copyright(c) 2009-2010, 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<Riostream.h>
16#include<TClonesArray.h>
17#include<TMath.h>
18#include "AliStrLine.h"
19#include "AliITSSortTrkl.h"
20
21/* $Id$ */
22
23////////////////////////////////////////////////////////////////////////
24// Helper class for finding multiple primary vertices //
25// To be used by AliITSVertexer3D //
26// It is based on the association of pairs of tracklets //
27// obtained by matching reconstructed points onthe first //
28// 2 layers of SPD //
29// Origin M. Masera masera@to.infn.it //
30////////////////////////////////////////////////////////////////////////
31
32
33ClassImp(AliITSSortTrkl)
34
35//______________________________________________________________________
36AliITSSortTrkl::AliITSSortTrkl():TObject(),
37fkSize(0),
38fIndex(0),
39fPairs(NULL),
40fClustersTmp(NULL),
41fClusters(NULL),
42fNoClus(0),
43fSize(NULL),
44fMark(),
45fCut(0.),
46fCoarseMaxRCut(0.)
47{
48 // Default constructor
49}
50
51//______________________________________________________________________
52AliITSSortTrkl::AliITSSortTrkl(Int_t n, Double_t cut):TObject(),
53fkSize(n),
54fIndex(0),
55fPairs(),
56fClustersTmp(),
57fClusters(NULL),
58fNoClus(0),
59fSize(NULL),
60fMark(n),
61fCut(cut),
62fCoarseMaxRCut(0.) {
63 // Standard constructor
64 fPairs = new AliITSTracklPairs* [fkSize];
65 fClustersTmp = new Int_t* [fkSize-1];
66}
67
68//______________________________________________________________________
69AliITSSortTrkl::AliITSSortTrkl(TClonesArray &tclo, Int_t n, Double_t cut, Double_t rcut):TObject(),
70fkSize(n),
71fIndex(0),
72fPairs(),
73fClustersTmp(),
74fClusters(NULL),
75fNoClus(0),
76fSize(NULL),
77fMark(n),
78fCut(cut),
79fCoarseMaxRCut(rcut){
80 // Constructor based on a TClonesArray of AliStrLine
81 Int_t numtrack=tclo.GetEntriesFast();
82 if(numtrack<2){
83 AliFatal(Form("Insufficient number of tracks (%d )",numtrack));
84 return;
85 }
86 TObject *member = tclo[0];
87 TString str(member->ClassName());
88 if(!(str.Contains("AliStrLine"))){
89 AliFatal(Form("Wrong type of class in input TClonesArray (%s )",str.Data()));
90 return;
91 }
92 Int_t siz = numtrack*(numtrack-1)/2;
93 if(fkSize < siz){
94 AliError(Form("fkSize is too small. It is %d and it should be %d",fkSize,siz));
95 }
96 fPairs = new AliITSTracklPairs* [fkSize];
97 fClustersTmp = new Int_t* [fkSize-1];
98 for(Int_t i=0;i<numtrack-1;i++){
99 AliStrLine *one = (AliStrLine*)tclo[i];
100 for(Int_t j=i+1;j<numtrack;j++){
101 AliStrLine *two = (AliStrLine*)tclo[j];
102 Double_t point[3];
103 one->Cross(two,point);
104 Double_t dca = one->GetDCA(two);
105 if(dca>fCut)continue;
106 Double_t rad=TMath::Sqrt(point[0]*point[0]+point[1]*point[1]);
107 if(rad>fCoarseMaxRCut)continue;
108 AddPairs(i,j,dca,point);
109 }
110 }
111}
112
113//______________________________________________________________________
114AliITSSortTrkl::~AliITSSortTrkl(){
115 // Destructor
116 if(fPairs){
117 for(Int_t i=0;i<fIndex;i++)delete fPairs[i];
118 delete [] fPairs;
119 }
120 DeleteClustersTmp();
121 if(fClusters){
122 for(Int_t i=0;i<fNoClus;i++)delete []fClusters[i];
123 delete fClusters;
124 delete []fSize;
125 }
126}
127
128//______________________________________________________________________
129void AliITSSortTrkl::DeleteClustersTmp(){
130 // fClustersTmp is deleted
131 if(fClustersTmp){
132 for(Int_t i=0;i<fIndex-1;i++)delete []fClustersTmp[i];
133 delete fClustersTmp;
134 fClustersTmp = NULL;
135 }
136}
137
138//______________________________________________________________________
139Int_t AliITSSortTrkl::AddPairs(Int_t t1, Int_t t2, Double_t dca, Double_t *coo){
140 // Add a tracklet pair at current position
141 fPairs[fIndex] = new AliITSTracklPairs(t1,t2,dca,coo);
142 // cout<<"Pair "<<fIndex<<" Tracklets "<<t1<<" and "<<t2<<". DCA= "<<dca<<" Crossing "<<coo[0]<<" "<<coo[1]<<" "<<coo[2]<<endl;
143 fIndex++;
144 return fIndex-1;
145}
146
147//______________________________________________________________________
148Int_t AliITSSortTrkl::FindClusters(){
149 // find clusters
150 fMark.ResetAllBits();
151 PrepareClustersTmp();
152 for(Int_t i=0; i<fIndex-1;i++){
153 Int_t *v=fClustersTmp[i];
154 Clustering(i,v);
155 if(v[0]>0){
156 AliInfo(Form("Clusters starting from pair %d : %d",i,v[0]));
157 Int_t dim;
158 v[v[0]+1]=i;
159 Int_t *arr=FindLabels(v+1,2*(v[0]+1),dim);
160 cout<<"Pairs involved \n";
161 for(Int_t j=1;j<=v[0];j++){
162 cout<<v[j]<<" ";
163 if(j%10==0)cout<<endl;
164 }
165 cout<<endl;
166 cout<<"Tracklets involved\n";
167 for(Int_t j=0;j<dim;j++){
168 cout<<arr[j]<<" ";
169 if(j%10==0 && j>0)cout<<endl;
170 }
171 cout<<endl;
172 for(Int_t j=0;j<fIndex;j++){
173 if(fMark.TestBitNumber(j))continue;
174 for(Int_t k=0;k<dim;k++){
175 if(fPairs[j]->HasTrack(arr[k])){
176 fMark.SetBitNumber(j);
177 // cout<<"Marked pair "<<j<<" because has tracklet "<<arr[k]<<endl;
178 k=dim;
179 }
180 }
181 }
182 delete []arr;
183 }
184 }
185 Cleanup();
186 return fNoClus;
187}
188
189//______________________________________________________________________
190void AliITSSortTrkl::Cleanup(){
191 // empty arrays are eliminated, the others are sorted according
192 // to cluster multiplicity
9e817851 193 Int_t *siz = new Int_t[fIndex-1];
194 Int_t *index = new Int_t[fIndex-1];
108cb8b9 195 fNoClus=0;
196 for(Int_t i=0; i<fIndex-1;i++){
197 Int_t *v=fClustersTmp[i];
198 if(v[0]>0)fNoClus++;
199 siz[i]=v[0];
200 }
9e817851 201 if(fNoClus == 0){
202 delete []siz;
203 delete [] index;
204 return;
205 }
108cb8b9 206 TMath::Sort(fIndex-1,siz,index);
207 fClusters = new Int_t* [fNoClus];
208 fSize = new Int_t [fNoClus];
209 for(Int_t i=0;i<fNoClus;i++){
210 Int_t curind = index[i];
211 Int_t *v=fClustersTmp[curind];
212 fClusters[i] = new Int_t[v[0]+1];
213 fSize[i]=v[0]+1;
214 Int_t *vext=fClusters[i];
215 vext[0]=curind;
216 for(Int_t j=1;j<fSize[i];j++)vext[j]=v[j];
217 }
9e817851 218 delete []siz;
219 delete [] index;
108cb8b9 220}
221
222//______________________________________________________________________
223Int_t* AliITSSortTrkl::GetTrackletsLab(Int_t index, Int_t& dim) const {
224 // Returns the tracklet labels corresponding to cluster index
225 // Calling code must take care of memory deallocation
226 if(fNoClus <=0){
227 dim = 0;
228 return NULL;
229 }
230 Int_t dimmax = 2*GetSizeOfCluster(index);
231 if(dimmax<=0){
232 dim = 0;
233 return NULL;
234 }
235 Int_t *v = GetClusters(index);
236 return FindLabels(v,dimmax,dim);
237}
238
239//______________________________________________________________________
240Int_t* AliITSSortTrkl::FindLabels(Int_t *v, Int_t dimmax, Int_t& dim) const {
241 // Returns the tracklet labels corresponding to te list of pairs
242 // contained in v.
243 // Calling code must take care of memory deallocation
244 Int_t *arr = new Int_t [dimmax];
245 Int_t j=0;
246 for(Int_t i=0; i<dimmax/2; i++){
247 AliITSTracklPairs *pai=GetPairsAt(v[i]);
248 arr[j++]=pai->GetTrack1();
249 arr[j++]=pai->GetTrack2();
250 }
251 SortAndClean(dimmax,arr,dim);
252 return arr;
253}
254//______________________________________________________________________
255void AliITSSortTrkl::SortAndClean(Int_t numb, Int_t *arr, Int_t& numb2){
256 // Array arr (with numb elements) is sorted in ascending order.
257 // Then possible reoccurrences
258 // of elements are eliminated. numb2 is the number of remaining elements
259 // after cleanup.
260 if(numb<=0)return;
9e817851 261 Int_t* index = new Int_t[numb];
108cb8b9 262 TMath::Sort(numb,arr,index,kFALSE);
9e817851 263 Int_t* tmp = new Int_t[numb];
108cb8b9 264 numb2 = 0;
265 tmp[0] = arr[index[0]];
266 for(Int_t i=1;i<numb;i++){
267 if(arr[index[i]] != tmp[numb2]){
268 ++numb2;
269 tmp[numb2]=arr[index[i]];
270 }
271 }
272 ++numb2;
273 for(Int_t i=0;i<numb;i++){
274 if(i<numb2){
275 arr[i]=tmp[i];
276 }
277 else {
278 arr[i]=0;
279 }
280 }
9e817851 281 delete [] index;
282 delete [] tmp;
108cb8b9 283}
284
285//______________________________________________________________________
286void AliITSSortTrkl::PrepareClustersTmp(){
287 // prepare arrays of clusters
288 for(Int_t i=0; i<fIndex-1;i++){
289 fClustersTmp[i] = new Int_t [fIndex-i];
290 Int_t *v = fClustersTmp[i];
291 v[0]=0;
292 for(Int_t j=1;j<fIndex-i;j++)v[j]=-1;
293 }
294}
295
296//______________________________________________________________________
297void AliITSSortTrkl::Clustering(Int_t i,Int_t *v){
298 // recursive method to build up clusters starting from point i
299 if(fMark.TestBitNumber(i))return;
300 AliITSTracklPairs *p1 = fPairs[i];
301 for(Int_t j=i+1;j<fIndex;j++){
302 if(fMark.TestBitNumber(j))continue;
303 AliITSTracklPairs *p2 = fPairs[j];
304 Double_t dist = p1->GetDistance(*p2);
305 // AliInfo(Form(" ******* i %d , j %d . Distance %g ",i,j,dist));
306 if(dist<=fCut){
307 Int_t dimclu=v[0];
308 Bool_t already = kFALSE;
309 for(Int_t k=1;k<=dimclu;k++){
310 if(v[k]==j)already=kTRUE;
311 }
312 if(!already){
313 ++dimclu;
314 v[0]=dimclu;
315 fMark.SetBitNumber(j);
316 fMark.SetBitNumber(i);
317 v[dimclu]=j;
318 Clustering(j,v);
319 }
320 }
321 }
322}
323
324
325