[u/mrichter/AliRoot.git] / ITS / AliITSSortTrkl.cxx

/**************************************************************************
 * Copyright(c) 2009-2010, ALICE Experiment at CERN, All rights reserved. *
 *                                                                        *
 * Author: The ALICE Off-line Project.                                    *
 * Contributors are mentioned in the code where appropriate.              *
 *                                                                        *
 * Permission to use, copy, modify and distribute this software and its   *
 * documentation strictly for non-commercial purposes is hereby granted   *
 * without fee, provided that the above copyright notice appears in all   *
 * copies and that both the copyright notice and this permission notice   *
 * appear in the supporting documentation. The authors make no claims     *
 * about the suitability of this software for any purpose. It is          *
 * provided "as is" without express or implied warranty.                  *
 **************************************************************************/
#include<Riostream.h>
#include<TClonesArray.h>
#include<TMath.h>
#include "AliStrLine.h"
#include "AliITSSortTrkl.h"

/* $Id$ */

////////////////////////////////////////////////////////////////////////
//           Helper class for finding multiple primary vertices       //
//           To be used by AliITSVertexer3D                           //
//           It is based on the association of pairs of tracklets     //
//           obtained by matching reconstructed points onthe first    //
//           2 layers of SPD                                          //
//           Origin M. Masera masera@to.infn.it                       //
////////////////////////////////////////////////////////////////////////


ClassImp(AliITSSortTrkl)

//______________________________________________________________________
AliITSSortTrkl::AliITSSortTrkl():TObject(),
fkSize(0),
fIndex(0),
fPairs(NULL),
fClustersTmp(NULL),
fClusters(NULL),
fNoClus(0),
fSize(NULL),
fMark(),
fCut(0.),
fCoarseMaxRCut(0.)
{
  // Default constructor
}

//______________________________________________________________________
AliITSSortTrkl::AliITSSortTrkl(Int_t n, Double_t cut):TObject(),
fkSize(n),
fIndex(0),
fPairs(),
fClustersTmp(NULL),
fClusters(NULL),
fNoClus(0),
fSize(NULL),
fMark(n),
fCut(cut),
fCoarseMaxRCut(0.) {
  // Standard constructor
  fPairs = new AliITSTracklPairs* [fkSize];
  fClustersTmp = new Int_t* [fkSize-1];
  for(Int_t i=0; i<fkSize-1; i++) fClustersTmp[i]=NULL;
}

//______________________________________________________________________
AliITSSortTrkl::AliITSSortTrkl(TClonesArray &tclo, Int_t n, Double_t cut, Double_t rcut):TObject(),
fkSize(n),
fIndex(0),
fPairs(),
fClustersTmp(NULL),
fClusters(NULL),
fNoClus(0),
fSize(NULL),
fMark(n),
fCut(cut),
fCoarseMaxRCut(rcut){
  // Constructor based on a TClonesArray of AliStrLine
  Int_t numtrack=tclo.GetEntriesFast();
  if(numtrack<2){
    AliFatal(Form("Insufficient number of tracks (%d )",numtrack));
    return;
  }
  TObject *member = tclo[0];
  TString str(member->ClassName());
  if(!(str.Contains("AliStrLine"))){
    AliFatal(Form("Wrong type of class in input TClonesArray (%s )",str.Data()));
    return;
  }
  Int_t siz = numtrack*(numtrack-1)/2;
  if(fkSize < siz){
    AliError(Form("fkSize is too small. It is %d and it should be %d",fkSize,siz)); 
  }
  fPairs = new AliITSTracklPairs* [fkSize];
  fClustersTmp = new Int_t* [fkSize-1];
  for(Int_t i=0; i<fkSize-1; i++) fClustersTmp[i]=NULL;
  for(Int_t i=0;i<numtrack-1;i++){
    AliStrLine *one = (AliStrLine*)tclo[i];
    for(Int_t j=i+1;j<numtrack;j++){
      AliStrLine *two = (AliStrLine*)tclo[j];
      Double_t point[3];
      one->Cross(two,point);
      Double_t dca = one->GetDCA(two);
      if(dca>fCut)continue;
      Double_t rad=TMath::Sqrt(point[0]*point[0]+point[1]*point[1]);
      if(rad>fCoarseMaxRCut)continue;
      AddPairs(i,j,dca,point);
    }
  } 
}

//______________________________________________________________________
AliITSSortTrkl::~AliITSSortTrkl(){
  // Destructor
  if(fPairs){
    for(Int_t i=0;i<fIndex;i++)delete fPairs[i];
    delete [] fPairs;
  }
  DeleteClustersTmp();
  if(fClusters){
    for(Int_t i=0;i<fNoClus;i++)delete [] fClusters[i];
    delete [] fClusters;
    delete [] fSize;
  } 
}

//______________________________________________________________________
void AliITSSortTrkl::DeleteClustersTmp(){
  // fClustersTmp is deleted
  if(fClustersTmp){
    for(Int_t i=0;i<fIndex-1;i++){
      if(fClustersTmp[i]){
	delete []fClustersTmp[i];
      }
    }
    delete [] fClustersTmp;
    fClustersTmp = NULL;
  }
}

//______________________________________________________________________
Int_t AliITSSortTrkl::AddPairs(Int_t t1, Int_t t2, Double_t dca, Double_t *coo){
  // Add a tracklet pair at current position
  fPairs[fIndex] = new AliITSTracklPairs(t1,t2,dca,coo);
  //  cout<<"Pair "<<fIndex<<" Tracklets "<<t1<<" and "<<t2<<". DCA= "<<dca<<" Crossing "<<coo[0]<<" "<<coo[1]<<" "<<coo[2]<<endl;
  fIndex++;
  return fIndex-1;
}

//______________________________________________________________________
Int_t AliITSSortTrkl::FindClusters(){
  // find clusters
  if(fIndex<2){
    AliWarning(Form("fIndex = %d",fIndex));
    fNoClus = 0;
    return fNoClus;
  }
  fMark.ResetAllBits();
  PrepareClustersTmp();
  //  cout<<"AliITSSortTrkl::FindClusters fkSize "<<fkSize<<"; fIndex "<<fIndex<<endl;
  for(Int_t i=0; i<fIndex-1;i++){
    Int_t *v=fClustersTmp[i];
    AliDebug(25,Form("Starting clustering for pair number %d",i));
    Clustering(i,v);  
    if(v[0]>0){
      AliDebug(25,Form("Clusters starting from pair %d : %d",i,v[0]));
      Int_t dim;
      v[v[0]+1]=i;
      // arr will contain the labels of the tracks associated to
      // the cluster of pairs starting from pair i.
      Int_t *arr=FindLabels(v+1,2*(v[0]+1),dim);

      /*
      cout<<"AliITSSortTrkl::FindClusters: Pairs involved \n";
      for(Int_t j=1;j<=v[0]+1;j++){
	cout<<v[j]<<" ";
	if(j%10==0)cout<<endl;
      }
      cout<<endl;
      cout<<"AliITSSortTrkl::FindClusters: Tracklets involved\n";
      for(Int_t j=0;j<dim;j++){
	cout<<arr[j]<<" ";
	if(j%10==0 && j>0)cout<<endl;
      }
      cout<<endl;
      */

      //In the following loop, all the pairs having 
      // one tracklet already associated to the cluster starting
      // at pair i are marked, in order to be excluded from further
      // searches    
      for(Int_t j=0;j<fIndex;j++){
	if(fMark.TestBitNumber(j))continue;
	for(Int_t k=0;k<dim;k++){
	  if(fPairs[j]->HasTrack(arr[k])){
	    fMark.SetBitNumber(j);
	    // cout<<"Marked pair "<<j<<" because has tracklet "<<arr[k]<<endl;
	    k=dim;
	  }
	}
      }
    
      delete []arr;
    }
  }
  // The following method builds the array fClusters
  Cleanup();
  return fNoClus;
}

//______________________________________________________________________
void AliITSSortTrkl::Cleanup(){
  // empty arrays are eliminated, the others are sorted according
  // to cluster multiplicity
  AliDebug(25,Form("fIndex = %d",fIndex));
  Int_t *siz = new Int_t[fIndex-1];
  Int_t *index = new Int_t[fIndex-1];
  fNoClus=0;
  for(Int_t i=0; i<fIndex-1;i++){
    Int_t *v=fClustersTmp[i];
    if(v[0]>0)fNoClus++;
    siz[i]=v[0];
  }
  if(fNoClus == 0){
    delete []siz;
    delete [] index;
    return;
  }
  AliDebug(25,Form("fNoClus = %d",fNoClus));
  TMath::Sort(fIndex-1,siz,index);
  fClusters = new Int_t* [fNoClus];
  fSize = new Int_t [fNoClus];
  for(Int_t i=0;i<fNoClus;i++){
    //    cout<<"AliITSSortTrkl::Cleanup: Cluster number "<<i<<"; Index= "<<index[i]<<endl;
    Int_t curind = index[i];
    Int_t *v=fClustersTmp[curind];
    fClusters[i] = new Int_t[v[0]+1];
    fSize[i]=v[0]+1;
    //    cout<<"AliITSSortTrkl::Cleanup: Size = "<<fSize[i]<<endl;
    Int_t *vext=fClusters[i];
    vext[0]=curind;
    for(Int_t j=1;j<fSize[i];j++)vext[j]=v[j];
    /*
    for(Int_t j=0;j<fSize[i];j++){
      cout<<vext[j]<<" ";
      if(j%10 == 0 && j!=0)cout<<endl;
    }
    cout<<endl;
    */
  }
  delete []siz;
  delete [] index;
}

//______________________________________________________________________
Int_t* AliITSSortTrkl::GetTrackletsLab(Int_t index, Int_t& dim) const {
  // Returns the tracklet labels corresponding to cluster index
  // Calling code must take care of memory deallocation
  AliDebug(25,Form("called with parameters %d %d",index,dim));
  if(fNoClus <=0){
    dim = 0;
    return NULL;
  }
  Int_t dimmax = 2*GetSizeOfCluster(index);
  AliDebug(25,Form("dimmax = %d",dimmax));
  if(dimmax<=0){
    dim = 0;
    return NULL;
  } 
  Int_t *v = GetClusters(index);
  return FindLabels(v,dimmax,dim);
}

//______________________________________________________________________
Int_t* AliITSSortTrkl::FindLabels(Int_t *v, Int_t dimmax, Int_t& dim) const {
  // Returns the tracklet labels corresponding to te list of pairs 
  // contained in v.
  // Calling code must take care of memory deallocation

  //  cout<<"AliITSSortTrkl::Findlabels parameters "<<v<<" dimmax = "<<dimmax<<" dim "<<dim<<endl;
  Int_t *arr = new Int_t [dimmax];
  Int_t j=0;
  for(Int_t i=0; i<dimmax/2; i++){
    AliITSTracklPairs *pai=GetPairsAt(v[i]);
    arr[j++]=pai->GetTrack1();
    //    cout<<"AliITSSortTrkl::FindLabels - i="<<i<<" arr["<<j-1<<"]= "<<arr[j-1]<<endl;
    arr[j++]=pai->GetTrack2();
    //    cout<<"AliITSSortTrkl::FindLabels arr["<<j-1<<"]= "<<arr[j-1]<<endl;
  }
  SortAndClean(dimmax,arr,dim);
  return arr;
}
//______________________________________________________________________
void AliITSSortTrkl::SortAndClean(Int_t numb, Int_t *arr, Int_t& numb2){
  // Array arr (with numb elements) is sorted in ascending order. 
  // Then possible reoccurrences
  // of elements are eliminated. numb2 is the number of remaining elements
  // after cleanup.
   
  //  cout<<"AliITSSortTrkl::SortAndClean - Parameters: numb= "<<numb<<" numb2= "<<numb2<<endl;
  if(numb<=0)return;
  Int_t* index = new Int_t[numb];
  TMath::Sort(numb,arr,index,kFALSE);
  Int_t* tmp = new Int_t[numb];
  numb2 = 0;
  tmp[0] = arr[index[0]];
  for(Int_t i=1;i<numb;i++){
    if(arr[index[i]] != tmp[numb2]){
      ++numb2;
      tmp[numb2]=arr[index[i]];
    }
  }
  ++numb2;
  /*
  cout<<"AliITSSortTrkl::SortAndClean - numb2 = "<<numb2<<endl;
  for(Int_t i=0;i<numb;i++){
    if(i<numb2){
      arr[i]=tmp[i];
      cout<<"arr["<<i<<"]= "<<arr[i]<<endl;
    }
    else {
      arr[i]=0;
    }
  }
  */
  delete [] index;
  delete [] tmp;
}

//______________________________________________________________________
void AliITSSortTrkl::PrepareClustersTmp(){
  // prepare arrays of clusters
  for(Int_t i=0; i<fIndex-1;i++){
    fClustersTmp[i] = new Int_t [fIndex+1];
    Int_t *v = fClustersTmp[i];
    v[0]=0;
    for(Int_t j=1;j<fIndex+1;j++)v[j]=-1;
  }
}


//______________________________________________________________________
void AliITSSortTrkl::Clustering(Int_t i,Int_t *v){
  // recursive method to build up clusters starting from point i
  AliDebug(25,Form("Clustering called for pair %d",i));
  if(fMark.TestBitNumber(i)){
    AliDebug(25,Form("Leaving Clustering for pair %d - nothing done",i));
  return;
  }
  AliITSTracklPairs *p1 = fPairs[i];
  for(Int_t j=0;j<fIndex;j++){
    if(j == i) continue;
    if(fMark.TestBitNumber(j))continue;
    AliITSTracklPairs *p2 = fPairs[j];
    Double_t dist = p1->GetDistance(*p2);
    //    AliInfo(Form("  ******* i %d , j %d . Distance %g ",i,j,dist));
    if(dist<=fCut){
      Int_t dimclu=v[0];
      Bool_t already = kFALSE;
      for(Int_t k=1;k<=dimclu;k++){
	if(v[k]==j)already=kTRUE;
      }
      if(!already){
	++dimclu;
	v[0]=dimclu;
	fMark.SetBitNumber(i);
	AliDebug(25,Form("Marked pair %d - and call Clustering for pair %d",i,j));
	v[dimclu]=j;
	Clustering(j,v);
	fMark.SetBitNumber(j);
	AliDebug(25,Form("Marked pair %d",j));
      }
    }
  }
  AliDebug(25,Form("Leaving Clustering for pair %d ",i));
}
Commit	Line	Data
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
	33	ClassImp(AliITSSortTrkl)
	34
	35	//______________________________________________________________________
	36	AliITSSortTrkl::AliITSSortTrkl():TObject(),
	37	fkSize(0),
	38	fIndex(0),
	39	fPairs(NULL),
	40	fClustersTmp(NULL),
	41	fClusters(NULL),
	42	fNoClus(0),
	43	fSize(NULL),
	44	fMark(),
	45	fCut(0.),
	46	fCoarseMaxRCut(0.)
	47	{
	48	// Default constructor
	49	}
	50
	51	//______________________________________________________________________
	52	AliITSSortTrkl::AliITSSortTrkl(Int_t n, Double_t cut):TObject(),
	53	fkSize(n),
	54	fIndex(0),
	55	fPairs(),
1ca5c3f5	56	fClustersTmp(NULL),
108cb8b9	57	fClusters(NULL),
	58	fNoClus(0),
	59	fSize(NULL),
	60	fMark(n),
	61	fCut(cut),
	62	fCoarseMaxRCut(0.) {
	63	// Standard constructor
	64	fPairs = new AliITSTracklPairs* [fkSize];
1ca5c3f5	65	fClustersTmp = new Int_t* [fkSize-1];
1ca5c3f5	66	for(Int_t i=0; i<fkSize-1; i++) fClustersTmp[i]=NULL;
108cb8b9	67	}
	68
	69	//______________________________________________________________________
	70	AliITSSortTrkl::AliITSSortTrkl(TClonesArray &tclo, Int_t n, Double_t cut, Double_t rcut):TObject(),
	71	fkSize(n),
	72	fIndex(0),
	73	fPairs(),
1ca5c3f5	74	fClustersTmp(NULL),
108cb8b9	75	fClusters(NULL),
	76	fNoClus(0),
	77	fSize(NULL),
	78	fMark(n),
	79	fCut(cut),
	80	fCoarseMaxRCut(rcut){
	81	// Constructor based on a TClonesArray of AliStrLine
	82	Int_t numtrack=tclo.GetEntriesFast();
	83	if(numtrack<2){
	84	AliFatal(Form("Insufficient number of tracks (%d )",numtrack));
	85	return;
	86	}
	87	TObject *member = tclo[0];
	88	TString str(member->ClassName());
	89	if(!(str.Contains("AliStrLine"))){
	90	AliFatal(Form("Wrong type of class in input TClonesArray (%s )",str.Data()));
	91	return;
	92	}
	93	Int_t siz = numtrack*(numtrack-1)/2;
	94	if(fkSize < siz){
	95	AliError(Form("fkSize is too small. It is %d and it should be %d",fkSize,siz));
	96	}
	97	fPairs = new AliITSTracklPairs* [fkSize];
	98	fClustersTmp = new Int_t* [fkSize-1];
1ca5c3f5	99	for(Int_t i=0; i<fkSize-1; i++) fClustersTmp[i]=NULL;
108cb8b9	100	for(Int_t i=0;i<numtrack-1;i++){
	101	AliStrLine one = (AliStrLine)tclo[i];
	102	for(Int_t j=i+1;j<numtrack;j++){
	103	AliStrLine two = (AliStrLine)tclo[j];
	104	Double_t point[3];
	105	one->Cross(two,point);
	106	Double_t dca = one->GetDCA(two);
	107	if(dca>fCut)continue;
	108	Double_t rad=TMath::Sqrt(point[0]point[0]+point[1]point[1]);
	109	if(rad>fCoarseMaxRCut)continue;
	110	AddPairs(i,j,dca,point);
	111	}
	112	}
	113	}
	114
	115	//______________________________________________________________________
	116	AliITSSortTrkl::~AliITSSortTrkl(){
	117	// Destructor
	118	if(fPairs){
	119	for(Int_t i=0;i<fIndex;i++)delete fPairs[i];
	120	delete [] fPairs;
	121	}
	122	DeleteClustersTmp();
	123	if(fClusters){
6fefa5ce	124	for(Int_t i=0;i<fNoClus;i++)delete [] fClusters[i];
	125	delete [] fClusters;
	126	delete [] fSize;
108cb8b9	127	}
	128	}
	129
	130	//______________________________________________________________________
	131	void AliITSSortTrkl::DeleteClustersTmp(){
	132	// fClustersTmp is deleted
	133	if(fClustersTmp){
1ca5c3f5	134	for(Int_t i=0;i<fIndex-1;i++){
	135	if(fClustersTmp[i]){
	136	delete []fClustersTmp[i];
	137	}
	138	}
6fefa5ce	139	delete [] fClustersTmp;
108cb8b9	140	fClustersTmp = NULL;
	141	}
	142	}
	143
	144	//______________________________________________________________________
	145	Int_t AliITSSortTrkl::AddPairs(Int_t t1, Int_t t2, Double_t dca, Double_t *coo){
	146	// Add a tracklet pair at current position
	147	fPairs[fIndex] = new AliITSTracklPairs(t1,t2,dca,coo);
	148	// cout<<"Pair "<<fIndex<<" Tracklets "<<t1<<" and "<<t2<<". DCA= "<<dca<<" Crossing "<<coo[0]<<" "<<coo[1]<<" "<<coo[2]<<endl;
	149	fIndex++;
	150	return fIndex-1;
	151	}
	152
	153	//______________________________________________________________________
	154	Int_t AliITSSortTrkl::FindClusters(){
	155	// find clusters
1ca5c3f5	156	if(fIndex<2){
	157	AliWarning(Form("fIndex = %d",fIndex));
	158	fNoClus = 0;
	159	return fNoClus;
	160	}
108cb8b9	161	fMark.ResetAllBits();
108cb8b9	162	PrepareClustersTmp();
1ca5c3f5	163	// cout<<"AliITSSortTrkl::FindClusters fkSize "<<fkSize<<"; fIndex "<<fIndex<<endl;
108cb8b9	164	for(Int_t i=0; i<fIndex-1;i++){
108cb8b9	165	Int_t *v=fClustersTmp[i];
1ca5c3f5	166	AliDebug(25,Form("Starting clustering for pair number %d",i));
1ca5c3f5	167	Clustering(i,v);
108cb8b9	168	if(v[0]>0){
1ca5c3f5	169	AliDebug(25,Form("Clusters starting from pair %d : %d",i,v[0]));
108cb8b9	170	Int_t dim;
108cb8b9	171	v[v[0]+1]=i;
1ca5c3f5	172	// arr will contain the labels of the tracks associated to
1ca5c3f5	173	// the cluster of pairs starting from pair i.
108cb8b9	174	Int_t arr=FindLabels(v+1,2(v[0]+1),dim);
1ca5c3f5	175
	176	/*
	177	cout<<"AliITSSortTrkl::FindClusters: Pairs involved \n";
	178	for(Int_t j=1;j<=v[0]+1;j++){
108cb8b9	179	cout<<v[j]<<" ";
	180	if(j%10==0)cout<<endl;
	181	}
	182	cout<<endl;
1ca5c3f5	183	cout<<"AliITSSortTrkl::FindClusters: Tracklets involved\n";
108cb8b9	184	for(Int_t j=0;j<dim;j++){
	185	cout<<arr[j]<<" ";
	186	if(j%10==0 && j>0)cout<<endl;
	187	}
	188	cout<<endl;
1ca5c3f5	189	*/
	190
	191	//In the following loop, all the pairs having
	192	// one tracklet already associated to the cluster starting
	193	// at pair i are marked, in order to be excluded from further
	194	// searches
108cb8b9	195	for(Int_t j=0;j<fIndex;j++){
	196	if(fMark.TestBitNumber(j))continue;
	197	for(Int_t k=0;k<dim;k++){
	198	if(fPairs[j]->HasTrack(arr[k])){
	199	fMark.SetBitNumber(j);
1ca5c3f5	200	// cout<<"Marked pair "<<j<<" because has tracklet "<<arr[k]<<endl;
108cb8b9	201	k=dim;
	202	}
	203	}
	204	}
1ca5c3f5	205
108cb8b9	206	delete []arr;
	207	}
	208	}
1ca5c3f5	209	// The following method builds the array fClusters
108cb8b9	210	Cleanup();
	211	return fNoClus;
	212	}
	213
	214	//______________________________________________________________________
	215	void AliITSSortTrkl::Cleanup(){
	216	// empty arrays are eliminated, the others are sorted according
	217	// to cluster multiplicity
1ca5c3f5	218	AliDebug(25,Form("fIndex = %d",fIndex));
9e817851	219	Int_t *siz = new Int_t[fIndex-1];
9e817851	220	Int_t *index = new Int_t[fIndex-1];
108cb8b9	221	fNoClus=0;
	222	for(Int_t i=0; i<fIndex-1;i++){
	223	Int_t *v=fClustersTmp[i];
	224	if(v[0]>0)fNoClus++;
	225	siz[i]=v[0];
	226	}
9e817851	227	if(fNoClus == 0){
	228	delete []siz;
	229	delete [] index;
	230	return;
	231	}
1ca5c3f5	232	AliDebug(25,Form("fNoClus = %d",fNoClus));
108cb8b9	233	TMath::Sort(fIndex-1,siz,index);
	234	fClusters = new Int_t* [fNoClus];
	235	fSize = new Int_t [fNoClus];
	236	for(Int_t i=0;i<fNoClus;i++){
1ca5c3f5	237	// cout<<"AliITSSortTrkl::Cleanup: Cluster number "<<i<<"; Index= "<<index[i]<<endl;
108cb8b9	238	Int_t curind = index[i];
	239	Int_t *v=fClustersTmp[curind];
	240	fClusters[i] = new Int_t[v[0]+1];
	241	fSize[i]=v[0]+1;
1ca5c3f5	242	// cout<<"AliITSSortTrkl::Cleanup: Size = "<<fSize[i]<<endl;
108cb8b9	243	Int_t *vext=fClusters[i];
	244	vext[0]=curind;
	245	for(Int_t j=1;j<fSize[i];j++)vext[j]=v[j];
1ca5c3f5	246	/*
	247	for(Int_t j=0;j<fSize[i];j++){
	248	cout<<vext[j]<<" ";
	249	if(j%10 == 0 && j!=0)cout<<endl;
	250	}
	251	cout<<endl;
	252	*/
108cb8b9	253	}
9e817851	254	delete []siz;
9e817851	255	delete [] index;
108cb8b9	256	}
	257
	258	//______________________________________________________________________
	259	Int_t* AliITSSortTrkl::GetTrackletsLab(Int_t index, Int_t& dim) const {
	260	// Returns the tracklet labels corresponding to cluster index
	261	// Calling code must take care of memory deallocation
1ca5c3f5	262	AliDebug(25,Form("called with parameters %d %d",index,dim));
108cb8b9	263	if(fNoClus <=0){
	264	dim = 0;
	265	return NULL;
	266	}
	267	Int_t dimmax = 2*GetSizeOfCluster(index);
1ca5c3f5	268	AliDebug(25,Form("dimmax = %d",dimmax));
108cb8b9	269	if(dimmax<=0){
	270	dim = 0;
	271	return NULL;
	272	}
	273	Int_t *v = GetClusters(index);
	274	return FindLabels(v,dimmax,dim);
	275	}
	276
	277	//______________________________________________________________________
	278	Int_t* AliITSSortTrkl::FindLabels(Int_t *v, Int_t dimmax, Int_t& dim) const {
	279	// Returns the tracklet labels corresponding to te list of pairs
	280	// contained in v.
	281	// Calling code must take care of memory deallocation
1ca5c3f5	282
1ca5c3f5	283	// cout<<"AliITSSortTrkl::Findlabels parameters "<<v<<" dimmax = "<<dimmax<<" dim "<<dim<<endl;
108cb8b9	284	Int_t *arr = new Int_t [dimmax];
	285	Int_t j=0;
	286	for(Int_t i=0; i<dimmax/2; i++){
	287	AliITSTracklPairs *pai=GetPairsAt(v[i]);
	288	arr[j++]=pai->GetTrack1();
1ca5c3f5	289	// cout<<"AliITSSortTrkl::FindLabels - i="<<i<<" arr["<<j-1<<"]= "<<arr[j-1]<<endl;
108cb8b9	290	arr[j++]=pai->GetTrack2();
1ca5c3f5	291	// cout<<"AliITSSortTrkl::FindLabels arr["<<j-1<<"]= "<<arr[j-1]<<endl;
108cb8b9	292	}
	293	SortAndClean(dimmax,arr,dim);
	294	return arr;
	295	}
	296	//______________________________________________________________________
	297	void AliITSSortTrkl::SortAndClean(Int_t numb, Int_t *arr, Int_t& numb2){
	298	// Array arr (with numb elements) is sorted in ascending order.
	299	// Then possible reoccurrences
	300	// of elements are eliminated. numb2 is the number of remaining elements
	301	// after cleanup.
1ca5c3f5	302
1ca5c3f5	303	// cout<<"AliITSSortTrkl::SortAndClean - Parameters: numb= "<<numb<<" numb2= "<<numb2<<endl;
108cb8b9	304	if(numb<=0)return;
9e817851	305	Int_t* index = new Int_t[numb];
108cb8b9	306	TMath::Sort(numb,arr,index,kFALSE);
9e817851	307	Int_t* tmp = new Int_t[numb];
108cb8b9	308	numb2 = 0;
	309	tmp[0] = arr[index[0]];
	310	for(Int_t i=1;i<numb;i++){
	311	if(arr[index[i]] != tmp[numb2]){
	312	++numb2;
	313	tmp[numb2]=arr[index[i]];
	314	}
	315	}
	316	++numb2;
1ca5c3f5	317	/*
1ca5c3f5	318	cout<<"AliITSSortTrkl::SortAndClean - numb2 = "<<numb2<<endl;
108cb8b9	319	for(Int_t i=0;i<numb;i++){
	320	if(i<numb2){
	321	arr[i]=tmp[i];
1ca5c3f5	322	cout<<"arr["<<i<<"]= "<<arr[i]<<endl;
108cb8b9	323	}
	324	else {
	325	arr[i]=0;
	326	}
	327	}
1ca5c3f5	328	*/
9e817851	329	delete [] index;
9e817851	330	delete [] tmp;
108cb8b9	331	}
	332
	333	//______________________________________________________________________
	334	void AliITSSortTrkl::PrepareClustersTmp(){
	335	// prepare arrays of clusters
	336	for(Int_t i=0; i<fIndex-1;i++){
1ca5c3f5	337	fClustersTmp[i] = new Int_t [fIndex+1];
108cb8b9	338	Int_t *v = fClustersTmp[i];
108cb8b9	339	v[0]=0;
1ca5c3f5	340	for(Int_t j=1;j<fIndex+1;j++)v[j]=-1;
108cb8b9	341	}
	342	}
	343
1ca5c3f5	344
108cb8b9	345	//______________________________________________________________________
	346	void AliITSSortTrkl::Clustering(Int_t i,Int_t *v){
	347	// recursive method to build up clusters starting from point i
1ca5c3f5	348	AliDebug(25,Form("Clustering called for pair %d",i));
	349	if(fMark.TestBitNumber(i)){
	350	AliDebug(25,Form("Leaving Clustering for pair %d - nothing done",i));
	351	return;
	352	}
108cb8b9	353	AliITSTracklPairs *p1 = fPairs[i];
1ca5c3f5	354	for(Int_t j=0;j<fIndex;j++){
1ca5c3f5	355	if(j == i) continue;
108cb8b9	356	if(fMark.TestBitNumber(j))continue;
	357	AliITSTracklPairs *p2 = fPairs[j];
	358	Double_t dist = p1->GetDistance(*p2);
	359	// AliInfo(Form(" ******* i %d , j %d . Distance %g ",i,j,dist));
	360	if(dist<=fCut){
	361	Int_t dimclu=v[0];
	362	Bool_t already = kFALSE;
	363	for(Int_t k=1;k<=dimclu;k++){
	364	if(v[k]==j)already=kTRUE;
	365	}
	366	if(!already){
	367	++dimclu;
	368	v[0]=dimclu;
108cb8b9	369	fMark.SetBitNumber(i);
1ca5c3f5	370	AliDebug(25,Form("Marked pair %d - and call Clustering for pair %d",i,j));
108cb8b9	371	v[dimclu]=j;
108cb8b9	372	Clustering(j,v);
1ca5c3f5	373	fMark.SetBitNumber(j);
1ca5c3f5	374	AliDebug(25,Form("Marked pair %d",j));
108cb8b9	375	}
	376	}
	377	}
1ca5c3f5	378	AliDebug(25,Form("Leaving Clustering for pair %d ",i));
108cb8b9	379	}
	380
	381
	382
1ca5c3f5	383
1ca5c3f5	384