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

/**************************************************************************
 * Copyright(c) 1998-1999, 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.                  *
 **************************************************************************/

/* $Id$ */

//  Cluster Finder for
//  Silicon
//  Pixel
//  developed by the Dubna group
//
#include <Riostream.h>
#include "AliITSClusterFinderSPDdubna.h"
#include "AliITSMapA1.h"
#include "AliITS.h"
#include "AliITSdigitSPD.h"
#include "AliITSRawClusterSPD.h"
#include "AliITSRecPoint.h"
#include "AliITSsegmentation.h"
#include "AliRun.h"

//#define DEBUG

ClassImp(AliITSClusterFinderSPDdubna)

//----------------------------------------------------------------------
AliITSClusterFinderSPDdubna::AliITSClusterFinderSPDdubna(
           AliITSsegmentation *seg, TClonesArray *digits, TClonesArray *recp){
    // constructor

    fSegmentation = seg;
    fDigits       = digits;
    fClusters     = recp;
    fNclusters    = 0;
    fMap          = new AliITSMapA1(fSegmentation,fDigits);
    SetDx();
    SetDz();
    SetNCells();
}
//_______________________________________________________________________
AliITSClusterFinderSPDdubna::AliITSClusterFinderSPDdubna(){
    // constructor

    fSegmentation = 0;
    fDigits       = 0;
    fClusters     = 0;
    fNclusters    = 0;
    fMap          = 0;
    SetDx();
    SetDz();
    SetNCells();
}
//______________________________________________________________________
AliITSClusterFinderSPDdubna &AliITSClusterFinderSPDdubna::operator=(const
                                          AliITSClusterFinderSPDdubna &s){
    // The = operator for the class AliITSClusterFinderSPDdugna
    // Inputs:
    //    AliITSClusterFinderSPDdubna @s  The AliITSClusterFinderSPDdubna
    //                                    class to be copy constructed.
    // Outputs:
    //    none.
    // Returned:
    //    this, a copy of the class s.

    this->fNclusters = s.fNclusters;
    this->fDz        = s.fDz;
    this->fDx        = s.fDx;
    this->fMinNCells = s.fMinNCells;
    this->fClusters  = new TClonesArray(*(s.fClusters));
    return *this;
}
//______________________________________________________________________
AliITSClusterFinderSPDdubna::AliITSClusterFinderSPDdubna(
      const AliITSClusterFinderSPDdubna &s): AliITSClusterFinder(s){
    // The Copy constructortor the class AliITSClusterFinderSPDdugna
    // It calles the = operator.
    // Inputs:
    //    AliITSClusterFinderSPDdubna @s  The AliITSClusterFinderSPDdubna
    //                                    class to be copy constructed.
    // Outputs:
    //    none.
    // Returned:
    //    none.

    *this = s;
    return;
}
//______________________________________________________________________
AliITSClusterFinderSPDdubna::~AliITSClusterFinderSPDdubna(){
    // destructor

    if (fMap) delete fMap;
}
//______________________________________________________________________
void AliITSClusterFinderSPDdubna::Find1DClusters(Int_t mod){
    // Find one dimensional clusters, i.e.
    // in r*phi(x) direction for each colunm in z direction

    static AliITS *iTS=(AliITS*)gAlice->GetModule("ITS");
  
    // retrieve the parameters 
    Int_t fNofPixels     = fSegmentation->Npz(); 
    Int_t fMaxNofSamples = fSegmentation->Npx();

    // read in digits -> do not apply threshold 
    // signal in fired pixels is always 1
    fMap->FillMap();

    Int_t nofFoundClusters = 0;

    Int_t k,it,m;
    for(k=0;k<fNofPixels;k++) {
	Int_t mmax = 10;  // a size of the window for the cluster finding
	for(it=0;it<fMaxNofSamples;it++) {
	    Int_t lclx   = 0;
	    Int_t xstart = 0;
	    Int_t xstop  = 0;
	    Int_t id     = 0;
	    Int_t ilcl   = 0;
	    for(m=0;m<mmax;m++) {  // find the cluster inside the window
		id = it+m;
		if(id >= fMaxNofSamples) break;    // ! no possible for the fadc
		if(fMap->TestHit(k,id) == kUnused) {   // start of the cluster
		    lclx += 1;
		    if(lclx == 1) xstart = id;
		} // end if fMap
		if(lclx > 0 && fMap->TestHit(k,id) == kEmpty) {  
		    // end of cluster if a gap exists
		    xstop = id-1;
		    ilcl = 1;
		    break;
		} // end if lclx >0 ...
	    }   //  end of m-loop
	    if(lclx == 0 && ilcl == 0) it = id; // no cluster in the window,
	                                        // continue the "it" loop
	    if(id >= fMaxNofSamples && lclx == 0) break; // the x row finished
	    if(id < fMaxNofSamples && ilcl == 0 && lclx > 0) {  
		// cluster end is outside of the window,
		mmax += 5;     // increase mmax and repeat the cluster finding
		it -= 1;
	    } // end if id<
	    if(id >= fMaxNofSamples && lclx > 0) {  // the x row finished but
		xstop = fMaxNofSamples - 1;           // the end cluster exists
		ilcl = 1;
	    }  // end if id>=      
	    // ---  Calculate z and x coordinates for one dimensional clusters
	    if(ilcl == 1) {         // new cluster exists
		it = id;
		mmax = 10;
		nofFoundClusters++;
		Float_t clusterCharge = 0.;
		Float_t zpitch        = fSegmentation->Dpz(k+1);
		Int_t   dummy         = 0;
		Float_t clusterZ, dummyX; 
		fSegmentation->GetPadCxz(dummy,k,dummyX,clusterZ);
		Float_t zstart        = clusterZ - 0.5*zpitch;
		Float_t zstop         = clusterZ + 0.5*zpitch;
		Float_t clusterX      = 0.;
		Int_t   clusterSizeX  = lclx;
		Int_t   clusterSizeZ  = 1;
		Int_t   its;
		for(its=xstart; its<=xstop; its++) {
		    Int_t firedpixel=0;
		    if (fMap->GetHitIndex(k,its)>=0) firedpixel=1; 
		    clusterCharge += firedpixel;
		    clusterX      += its + 0.5;
		} // end for its
		Float_t fRphiPitch = fSegmentation->Dpx(dummy);
		clusterX /= (clusterSizeX/fRphiPitch); // center of gravity for x
		// Write the points (coordinates and some cluster information) 
		// to the AliITSRawClusterSPD object

		AliITSRawClusterSPD clust(clusterZ,clusterX,clusterCharge,
					  clusterSizeZ,clusterSizeX,
					  xstart,xstop,zstart,zstop,k,mod);
		iTS->AddCluster(0,&clust);
	    }    // new cluster (ilcl=1)
	} // X direction loop (it)
    } // Z direction loop (k)
    //fMap->ClearMap();
    return;
}
//______________________________________________________________________
void  AliITSClusterFinderSPDdubna::GroupClusters(){
    // Find two dimensional clusters, i.e. group one dimensional clusters
    // into two dimensional ones (go both in x and z directions).
    // get number of clusters for this module

    Int_t nofClusters = fClusters->GetEntriesFast();
    nofClusters -= fNclusters;
    AliITSRawClusterSPD *clusterI;
    AliITSRawClusterSPD *clusterJ;
    Int_t *label=new Int_t[nofClusters];  
    Int_t i,j;
    for(i=0; i<nofClusters; i++) label[i] = 0;
    for(i=0; i<nofClusters; i++) {
	if(label[i] != 0) continue;
	for(j=i+1; j<nofClusters; j++) { 
	    if(label[j] != 0) continue;
	    clusterI = (AliITSRawClusterSPD*) fClusters->At(i);
	    clusterJ = (AliITSRawClusterSPD*) fClusters->At(j);
	    Bool_t pair = clusterI->Brother(clusterJ,fDz,fDx);
	    if(pair) {
#ifndef DEBUG
		cout << "clusters " << i << "," << j << " before grouping" << endl;
		clusterI->PrintInfo();
		clusterJ->PrintInfo();
#endif
		clusterI->Add(clusterJ);
#ifndef DEBUG
		cout << "remove cluster " << j << endl;
#endif
		label[j] = 1;
		fClusters->RemoveAt(j);
#ifndef DEBUG
		cout << "cluster  " << i << " after grouping" << endl;
		clusterI->PrintInfo();
#endif
	    }  // pair
	} // J clusters  
	label[i] = 1;
    } // I clusters
    fClusters->Compress();
#ifndef DEBUG
    Int_t totalNofClusters = fClusters->GetEntriesFast();
    cout << " Nomber of clusters at the group end ="<< totalNofClusters<<endl;
#endif
    delete [] label;
    return;
}
//______________________________________________________________________
void AliITSClusterFinderSPDdubna::TracksInCluster(){
    // Find tracks creating one cluster
    // get number of clusters for this module
    Int_t nofClusters = fClusters->GetEntriesFast();
    Int_t i, ix, iz, jx, jz, xstart, xstop, zstart, zstop, nclx, nclz;
    const Int_t kTrmax = 100;
    Int_t cltracks[kTrmax], itr, tracki, ii, is, js, ie, ntr, tr0, tr1, tr2;

    nofClusters -= fNclusters;
    for(i=0; i<nofClusters; i++) { 
	ii = 0;
	memset(cltracks,-1,sizeof(int)*kTrmax);
	tr0=tr1=tr2=-1;
	AliITSRawClusterSPD *clusterI = (AliITSRawClusterSPD*) fClusters->At(i);
	nclx = clusterI->NclX();
	nclz = clusterI->NclZ();
	xstart = clusterI->XStartf();
	xstop = clusterI->XStopf();
	zstart = clusterI->Zend()-nclz+1;
	zstop = clusterI->Zend();
	Int_t ind;
	for(iz=0; iz<nclz; iz++) { 
	    jz = zstart + iz;
	    for(ix=0; ix<nclx; ix++) { 
		jx = xstart + ix;
		ind = fMap->GetHitIndex(jz,jx);
		if(ind < 0) continue;
		if(ind == 0 && iz >= 0 && ix > 0) continue;
		if(ind == 0 && iz > 0 && ix >= 0) continue;
		if(ind == 0 && iz == 0 && ix == 0 && i > 0) continue;
		AliITSdigitSPD *dig = (AliITSdigitSPD*)fMap->GetHit(jz,jx);
		for(itr=0; itr<3; itr++) { 
		    tracki = dig->GetTrack(i);
		    if(tracki >= 0) {
			ii += 1;
			if(ii < 99) cltracks[ii-1] = tracki;
		    } // end if tracki >-0
		} // end for itr
	    } // ix pixel
	}  // iz pixel
	for(is=0; is<kTrmax; is++) { 
	    if(cltracks[is]<0) continue;
	    for(js=is+1; js<kTrmax; js++) { 
		if(cltracks[js]<0) continue;
		if(cltracks[js]==cltracks[is]) cltracks[js]=-5;
	    } // end for js
	} // end for is
	ntr = 0;
	for(ie=0; ie<kTrmax; ie++) { 
	    if(cltracks[ie] >= 0) {
		ntr=ntr+1;
		if(ntr==1) tr0=cltracks[ie];
		if(ntr==2) tr1=cltracks[ie];
		if(ntr==3) tr2=cltracks[ie];
	    } // end if cltracks[ie] >-=
	} // end for ie
	// if delta ray only
	if(ntr == 0) ntr = 1;
	clusterI->SetNTracks(ntr);
	clusterI->SetTracks(tr0,tr1,tr2);
    } // I cluster
}
//______________________________________________________________________
void AliITSClusterFinderSPDdubna::GetRecPoints(){
    // get rec points
    static AliITS *iTS = (AliITS*)gAlice->GetModule("ITS");
    // get number of clusters for this module
    Int_t nofClusters   = fClusters->GetEntriesFast();
    const Float_t kconv = 1.0e-4;
    const Float_t kRMSx = 12.0*kconv; // microns -> cm ITS TDR Table 1.3
    const Float_t kRMSz = 70.0*kconv; // microns -> cm ITS TDR Table 1.3
    Float_t spdLength   = fSegmentation->Dz();
    Float_t spdWidth    = fSegmentation->Dx();
    Int_t   dummy       = 0;
    Float_t xpitch      = fSegmentation->Dpx(dummy);
    Int_t i;
    Int_t track0, track1, track2;

    nofClusters -= fNclusters;
    for(i=0; i<nofClusters; i++) {
	AliITSRawClusterSPD *clusterI = (AliITSRawClusterSPD*) fClusters->At(i);
	Int_t clustersizex = clusterI->NclX();
	Int_t clustersizez = clusterI->NclZ();
	Int_t xstart = clusterI->XStartf();
	Float_t clusterx = 0;
	Float_t clusterz = clusterI->Z();
	for(Int_t ii=0; ii<clustersizex; ii++) clusterx += (xstart+0.5+ii)*xpitch;
	clusterx /= clustersizex;
	clusterz /= clustersizez;   
	clusterI->GetTracks(track0, track1, track2); 
	AliITSRecPoint rnew;
	rnew.SetX((clusterx - spdWidth/2)*kconv);
	rnew.SetZ((clusterz - spdLength/2)*kconv);
	rnew.SetQ(1.);
	rnew.SetdEdX(0.);
	rnew.SetSigmaX2(kRMSx*kRMSx);
	rnew.SetSigmaZ2(kRMSz*kRMSz);
	rnew.fTracks[0]=track0;
	rnew.fTracks[1]=track1;
	rnew.fTracks[2]=track2;
	iTS->AddRecPoint(rnew);
    } // I clusters
    fMap->ClearMap();
}
//______________________________________________________________________
void AliITSClusterFinderSPDdubna::FindRawClusters(Int_t mod){
    // find raw clusters

#ifndef DEBUG
    cout<<"SPDdubna: module ="<<mod<<endl;
#endif
    Find1DClusters(mod);
    GroupClusters();
    TracksInCluster();
    GetRecPoints();
}
Commit	Line	Data
409f8c84	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
88cb7938	16	/* $Id$ */
bf3f2830	17
bf3f2830	18	// Cluster Finder for
	19	// Silicon
	20	// Pixel
	21	// developed by the Dubna group
	22	//
4ae5bbc4	23	#include <Riostream.h>
409f8c84	24	#include "AliITSClusterFinderSPDdubna.h"
	25	#include "AliITSMapA1.h"
	26	#include "AliITS.h"
e869281d	27	#include "AliITSdigitSPD.h"
41b19549	28	#include "AliITSRawClusterSPD.h"
409f8c84	29	#include "AliITSRecPoint.h"
409f8c84	30	#include "AliITSsegmentation.h"
409f8c84	31	#include "AliRun.h"
409f8c84	32
9e1e0cd7	33	//#define DEBUG
409f8c84	34
	35	ClassImp(AliITSClusterFinderSPDdubna)
	36
9e1e0cd7	37	//----------------------------------------------------------------------
	38	AliITSClusterFinderSPDdubna::AliITSClusterFinderSPDdubna(
	39	AliITSsegmentation seg, TClonesArray digits, TClonesArray *recp){
	40	// constructor
	41
	42	fSegmentation = seg;
	43	fDigits = digits;
	44	fClusters = recp;
	45	fNclusters = 0;
	46	fMap = new AliITSMapA1(fSegmentation,fDigits);
409f8c84	47	SetDx();
	48	SetDz();
	49	SetNCells();
	50	}
9e1e0cd7	51	//_______________________________________________________________________
	52	AliITSClusterFinderSPDdubna::AliITSClusterFinderSPDdubna(){
	53	// constructor
	54
	55	fSegmentation = 0;
	56	fDigits = 0;
	57	fClusters = 0;
	58	fNclusters = 0;
	59	fMap = 0;
	60	SetDx();
	61	SetDz();
	62	SetNCells();
409f8c84	63	}
9e1e0cd7	64	//______________________________________________________________________
579fdb1c	65	AliITSClusterFinderSPDdubna &AliITSClusterFinderSPDdubna::operator=(const
	66	AliITSClusterFinderSPDdubna &s){
	67	// The = operator for the class AliITSClusterFinderSPDdugna
	68	// Inputs:
	69	// AliITSClusterFinderSPDdubna @s The AliITSClusterFinderSPDdubna
	70	// class to be copy constructed.
	71	// Outputs:
	72	// none.
	73	// Returned:
	74	// this, a copy of the class s.
	75
	76	this->fNclusters = s.fNclusters;
	77	this->fDz = s.fDz;
	78	this->fDx = s.fDx;
	79	this->fMinNCells = s.fMinNCells;
	80	this->fClusters = new TClonesArray(*(s.fClusters));
	81	return *this;
	82	}
	83	//______________________________________________________________________
ac74f489	84	AliITSClusterFinderSPDdubna::AliITSClusterFinderSPDdubna(
ac74f489	85	const AliITSClusterFinderSPDdubna &s): AliITSClusterFinder(s){
579fdb1c	86	// The Copy constructortor the class AliITSClusterFinderSPDdugna
	87	// It calles the = operator.
	88	// Inputs:
	89	// AliITSClusterFinderSPDdubna @s The AliITSClusterFinderSPDdubna
	90	// class to be copy constructed.
	91	// Outputs:
	92	// none.
	93	// Returned:
	94	// none.
	95
	96	*this = s;
	97	return;
	98	}
	99	//______________________________________________________________________
9e1e0cd7	100	AliITSClusterFinderSPDdubna::~AliITSClusterFinderSPDdubna(){
9e1e0cd7	101	// destructor
409f8c84	102
9e1e0cd7	103	if (fMap) delete fMap;
409f8c84	104	}
9e1e0cd7	105	//______________________________________________________________________
	106	void AliITSClusterFinderSPDdubna::Find1DClusters(Int_t mod){
	107	// Find one dimensional clusters, i.e.
	108	// in r*phi(x) direction for each colunm in z direction
338e4f06	109
9e1e0cd7	110	static AliITS iTS=(AliITS)gAlice->GetModule("ITS");
409f8c84	111
9e1e0cd7	112	// retrieve the parameters
	113	Int_t fNofPixels = fSegmentation->Npz();
	114	Int_t fMaxNofSamples = fSegmentation->Npx();
	115
	116	// read in digits -> do not apply threshold
	117	// signal in fired pixels is always 1
	118	fMap->FillMap();
	119
	120	Int_t nofFoundClusters = 0;
	121
	122	Int_t k,it,m;
	123	for(k=0;k<fNofPixels;k++) {
	124	Int_t mmax = 10; // a size of the window for the cluster finding
	125	for(it=0;it<fMaxNofSamples;it++) {
	126	Int_t lclx = 0;
	127	Int_t xstart = 0;
	128	Int_t xstop = 0;
	129	Int_t id = 0;
	130	Int_t ilcl = 0;
	131	for(m=0;m<mmax;m++) { // find the cluster inside the window
	132	id = it+m;
	133	if(id >= fMaxNofSamples) break; // ! no possible for the fadc
	134	if(fMap->TestHit(k,id) == kUnused) { // start of the cluster
	135	lclx += 1;
	136	if(lclx == 1) xstart = id;
	137	} // end if fMap
	138	if(lclx > 0 && fMap->TestHit(k,id) == kEmpty) {
	139	// end of cluster if a gap exists
	140	xstop = id-1;
	141	ilcl = 1;
	142	break;
	143	} // end if lclx >0 ...
	144	} // end of m-loop
	145	if(lclx == 0 && ilcl == 0) it = id; // no cluster in the window,
	146	// continue the "it" loop
	147	if(id >= fMaxNofSamples && lclx == 0) break; // the x row finished
	148	if(id < fMaxNofSamples && ilcl == 0 && lclx > 0) {
	149	// cluster end is outside of the window,
	150	mmax += 5; // increase mmax and repeat the cluster finding
	151	it -= 1;
	152	} // end if id<
	153	if(id >= fMaxNofSamples && lclx > 0) { // the x row finished but
	154	xstop = fMaxNofSamples - 1; // the end cluster exists
	155	ilcl = 1;
	156	} // end if id>=
	157	// --- Calculate z and x coordinates for one dimensional clusters
	158	if(ilcl == 1) { // new cluster exists
	159	it = id;
	160	mmax = 10;
	161	nofFoundClusters++;
	162	Float_t clusterCharge = 0.;
	163	Float_t zpitch = fSegmentation->Dpz(k+1);
	164	Int_t dummy = 0;
	165	Float_t clusterZ, dummyX;
	166	fSegmentation->GetPadCxz(dummy,k,dummyX,clusterZ);
	167	Float_t zstart = clusterZ - 0.5*zpitch;
	168	Float_t zstop = clusterZ + 0.5*zpitch;
	169	Float_t clusterX = 0.;
	170	Int_t clusterSizeX = lclx;
	171	Int_t clusterSizeZ = 1;
	172	Int_t its;
	173	for(its=xstart; its<=xstop; its++) {
	174	Int_t firedpixel=0;
	175	if (fMap->GetHitIndex(k,its)>=0) firedpixel=1;
176	clusterCharge += firedpixel;
177	clusterX += its + 0.5;
178	} // end for its
179	Float_t fRphiPitch = fSegmentation->Dpx(dummy);
180	clusterX /= (clusterSizeX/fRphiPitch); // center of gravity for x
181	// Write the points (coordinates and some cluster information)
182	// to the AliITSRawClusterSPD object
183
184	AliITSRawClusterSPD clust(clusterZ,clusterX,clusterCharge,
185	clusterSizeZ,clusterSizeX,
186	xstart,xstop,zstart,zstop,k,mod);
187	iTS->AddCluster(0,&clust);
188	} // new cluster (ilcl=1)
189	} // X direction loop (it)
190	} // Z direction loop (k)
191	//fMap->ClearMap();
192	return;
409f8c84	193	}
9e1e0cd7	194	//______________________________________________________________________
	195	void AliITSClusterFinderSPDdubna::GroupClusters(){
	196	// Find two dimensional clusters, i.e. group one dimensional clusters
	197	// into two dimensional ones (go both in x and z directions).
	198	// get number of clusters for this module
	199
	200	Int_t nofClusters = fClusters->GetEntriesFast();
	201	nofClusters -= fNclusters;
	202	AliITSRawClusterSPD *clusterI;
	203	AliITSRawClusterSPD *clusterJ;
	204	Int_t *label=new Int_t[nofClusters];
	205	Int_t i,j;
	206	for(i=0; i<nofClusters; i++) label[i] = 0;
	207	for(i=0; i<nofClusters; i++) {
	208	if(label[i] != 0) continue;
	209	for(j=i+1; j<nofClusters; j++) {
	210	if(label[j] != 0) continue;
	211	clusterI = (AliITSRawClusterSPD*) fClusters->At(i);
	212	clusterJ = (AliITSRawClusterSPD*) fClusters->At(j);
	213	Bool_t pair = clusterI->Brother(clusterJ,fDz,fDx);
	214	if(pair) {
	215	#ifndef DEBUG
	216	cout << "clusters " << i << "," << j << " before grouping" << endl;
	217	clusterI->PrintInfo();
	218	clusterJ->PrintInfo();
	219	#endif
	220	clusterI->Add(clusterJ);
	221	#ifndef DEBUG
	222	cout << "remove cluster " << j << endl;
	223	#endif
	224	label[j] = 1;
	225	fClusters->RemoveAt(j);
	226	#ifndef DEBUG
	227	cout << "cluster " << i << " after grouping" << endl;
	228	clusterI->PrintInfo();
	229	#endif
	230	} // pair
	231	} // J clusters
	232	label[i] = 1;
	233	} // I clusters
	234	fClusters->Compress();
	235	#ifndef DEBUG
409f8c84	236	Int_t totalNofClusters = fClusters->GetEntriesFast();
409f8c84	237	cout << " Nomber of clusters at the group end ="<< totalNofClusters<<endl;
9e1e0cd7	238	#endif
	239	delete [] label;
	240	return;
409f8c84	241	}
9e1e0cd7	242	//______________________________________________________________________
	243	void AliITSClusterFinderSPDdubna::TracksInCluster(){
	244	// Find tracks creating one cluster
	245	// get number of clusters for this module
	246	Int_t nofClusters = fClusters->GetEntriesFast();
	247	Int_t i, ix, iz, jx, jz, xstart, xstop, zstart, zstop, nclx, nclz;
bf3f2830	248	const Int_t kTrmax = 100;
bf3f2830	249	Int_t cltracks[kTrmax], itr, tracki, ii, is, js, ie, ntr, tr0, tr1, tr2;
9e1e0cd7	250
	251	nofClusters -= fNclusters;
	252	for(i=0; i<nofClusters; i++) {
	253	ii = 0;
bf3f2830	254	memset(cltracks,-1,sizeof(int)*kTrmax);
9e1e0cd7	255	tr0=tr1=tr2=-1;
	256	AliITSRawClusterSPD clusterI = (AliITSRawClusterSPD) fClusters->At(i);
	257	nclx = clusterI->NclX();
	258	nclz = clusterI->NclZ();
	259	xstart = clusterI->XStartf();
	260	xstop = clusterI->XStopf();
	261	zstart = clusterI->Zend()-nclz+1;
	262	zstop = clusterI->Zend();
	263	Int_t ind;
	264	for(iz=0; iz<nclz; iz++) {
	265	jz = zstart + iz;
	266	for(ix=0; ix<nclx; ix++) {
	267	jx = xstart + ix;
	268	ind = fMap->GetHitIndex(jz,jx);
	269	if(ind < 0) continue;
	270	if(ind == 0 && iz >= 0 && ix > 0) continue;
	271	if(ind == 0 && iz > 0 && ix >= 0) continue;
	272	if(ind == 0 && iz == 0 && ix == 0 && i > 0) continue;
	273	AliITSdigitSPD dig = (AliITSdigitSPD)fMap->GetHit(jz,jx);
	274	for(itr=0; itr<3; itr++) {
ecee53fc	275	tracki = dig->GetTrack(i);
9e1e0cd7	276	if(tracki >= 0) {
	277	ii += 1;
	278	if(ii < 99) cltracks[ii-1] = tracki;
	279	} // end if tracki >-0
	280	} // end for itr
	281	} // ix pixel
	282	} // iz pixel
bf3f2830	283	for(is=0; is<kTrmax; is++) {
9e1e0cd7	284	if(cltracks[is]<0) continue;
bf3f2830	285	for(js=is+1; js<kTrmax; js++) {
9e1e0cd7	286	if(cltracks[js]<0) continue;
	287	if(cltracks[js]==cltracks[is]) cltracks[js]=-5;
	288	} // end for js
	289	} // end for is
	290	ntr = 0;
bf3f2830	291	for(ie=0; ie<kTrmax; ie++) {
9e1e0cd7	292	if(cltracks[ie] >= 0) {
	293	ntr=ntr+1;
	294	if(ntr==1) tr0=cltracks[ie];
	295	if(ntr==2) tr1=cltracks[ie];
	296	if(ntr==3) tr2=cltracks[ie];
	297	} // end if cltracks[ie] >-=
	298	} // end for ie
	299	// if delta ray only
	300	if(ntr == 0) ntr = 1;
	301	clusterI->SetNTracks(ntr);
	302	clusterI->SetTracks(tr0,tr1,tr2);
	303	} // I cluster
409f8c84	304	}
9e1e0cd7	305	//______________________________________________________________________
	306	void AliITSClusterFinderSPDdubna::GetRecPoints(){
	307	// get rec points
	308	static AliITS iTS = (AliITS)gAlice->GetModule("ITS");
	309	// get number of clusters for this module
	310	Int_t nofClusters = fClusters->GetEntriesFast();
	311	const Float_t kconv = 1.0e-4;
	312	const Float_t kRMSx = 12.0*kconv; // microns -> cm ITS TDR Table 1.3
	313	const Float_t kRMSz = 70.0*kconv; // microns -> cm ITS TDR Table 1.3
	314	Float_t spdLength = fSegmentation->Dz();
	315	Float_t spdWidth = fSegmentation->Dx();
	316	Int_t dummy = 0;
	317	Float_t xpitch = fSegmentation->Dpx(dummy);
	318	Int_t i;
	319	Int_t track0, track1, track2;
	320
	321	nofClusters -= fNclusters;
	322	for(i=0; i<nofClusters; i++) {
	323	AliITSRawClusterSPD clusterI = (AliITSRawClusterSPD) fClusters->At(i);
	324	Int_t clustersizex = clusterI->NclX();
	325	Int_t clustersizez = clusterI->NclZ();
	326	Int_t xstart = clusterI->XStartf();
	327	Float_t clusterx = 0;
	328	Float_t clusterz = clusterI->Z();
	329	for(Int_t ii=0; ii<clustersizex; ii++) clusterx += (xstart+0.5+ii)*xpitch;
	330	clusterx /= clustersizex;
	331	clusterz /= clustersizez;
	332	clusterI->GetTracks(track0, track1, track2);
	333	AliITSRecPoint rnew;
	334	rnew.SetX((clusterx - spdWidth/2)*kconv);
	335	rnew.SetZ((clusterz - spdLength/2)*kconv);
	336	rnew.SetQ(1.);
	337	rnew.SetdEdX(0.);
	338	rnew.SetSigmaX2(kRMSx*kRMSx);
	339	rnew.SetSigmaZ2(kRMSz*kRMSz);
	340	rnew.fTracks[0]=track0;
	341	rnew.fTracks[1]=track1;
	342	rnew.fTracks[2]=track2;
	343	iTS->AddRecPoint(rnew);
	344	} // I clusters
	345	fMap->ClearMap();
409f8c84	346	}
9e1e0cd7	347	//______________________________________________________________________
	348	void AliITSClusterFinderSPDdubna::FindRawClusters(Int_t mod){
	349	// find raw clusters
	350
	351	#ifndef DEBUG
	352	cout<<"SPDdubna: module ="<<mod<<endl;
	353	#endif
	354	Find1DClusters(mod);
	355	GroupClusters();
	356	TracksInCluster();
	357	GetRecPoints();
409f8c84	358	}