[u/mrichter/AliRoot.git] / ITS / AliITSClusterFinderSDD.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.                  *
 **************************************************************************/


#include "AliITSClusterFinderSDD.h"
#include "AliRun.h"


ClassImp(AliITSClusterFinderSDD)

//----------------------------------------------------------
AliITSClusterFinderSDD::AliITSClusterFinderSDD
(AliITSsegmentation *seg, AliITSresponse *response, TClonesArray *digits, TClonesArray *recp)   
{
  // constructor
    fSegmentation=seg;
    fResponse=response;
    fDigits=digits;
    fClusters=recp;
    fNclusters= fClusters->GetEntriesFast();
    printf("SDD: fNclusters %d\n",fNclusters);
    SetCutAmplitude();
    SetDAnode();
    SetDTime();
    SetMap();
    SetMinPeak();
    SetNCells();
}

//_____________________________________________________________________________
AliITSClusterFinderSDD::AliITSClusterFinderSDD()
{
  // constructor
    fSegmentation=0;
    fResponse=0;
    fDigits=0;
    fClusters=0;
    fNclusters=0;
    SetCutAmplitude();
    SetDAnode();
    SetDTime();
    SetMap();
    SetMinPeak();
    SetNCells();

}

//__________________________________________________________________________
AliITSClusterFinderSDD::AliITSClusterFinderSDD(const AliITSClusterFinderSDD &source){
  //     Copy Constructor 
  if(&source == this) return;
  this->fClusters = source.fClusters ;
  this->fNclusters = source.fNclusters ;
  this->fMap = source.fMap ;
  this->fCutAmplitude = source.fCutAmplitude ;
  this->fDAnode = source.fDAnode ;
  this->fDTime = source.fDTime ;
  this->fMinPeak = source.fMinPeak ;
  this->fMinNCells = source.fMinNCells ;
  return;
}

//_________________________________________________________________________
AliITSClusterFinderSDD& 
  AliITSClusterFinderSDD::operator=(const AliITSClusterFinderSDD &source) {
  //    Assignment operator
  if(&source == this) return *this;
  this->fClusters = source.fClusters ;
  this->fNclusters = source.fNclusters ;
  this->fMap = source.fMap ;
  this->fCutAmplitude = source.fCutAmplitude ;
  this->fDAnode = source.fDAnode ;
  this->fDTime = source.fDTime ;
  this->fMinPeak = source.fMinPeak ;
  this->fMinNCells = source.fMinNCells ;
  return *this;
}

//_____________________________________________________________________________
void AliITSClusterFinderSDD::SetMap()
{
  // set map
    if(!fMap) fMap=new AliITSMapA2(fSegmentation);

}
//_____________________________________________________________________________
void AliITSClusterFinderSDD::FillMap()
{
  // fCoord1 = anode #
  // fCoord2 = time sample

  if (!fDigits) return; 
 
  Int_t ndigits = fDigits->GetEntriesFast();
  //printf("FillMap: ndigits %d\n",ndigits);
  if (!ndigits) return;

    AliITSdigitSDD *dig;
    Int_t ndig;
    for(ndig=0; ndig<ndigits; ndig++) {
        dig = (AliITSdigitSDD*)fDigits->UncheckedAt(ndig);
        Double_t signal=dig->fSignal;
	//printf("FillMap: ndig fCoord1 fCoord2 signal %d %d %d %f\n",ndig,dig->fCoord1,dig->fCoord2,signal);
        fMap->SetHit(dig->fCoord1,dig->fCoord2,signal);
    }


}
//_____________________________________________________________________________

void AliITSClusterFinderSDD::Find1DClusters()
{
  // find 1D clusters

    AliITS *iTS=(AliITS*)gAlice->GetModule("ITS");

   // retrieve the parameters 
   Int_t i;
   Int_t fNofMaps = fSegmentation->Npz();
   Int_t fMaxNofSamples = fSegmentation->Npx();
   Int_t fNofAnodes = fNofMaps/2;
   Int_t dummy=0;
   Float_t fTimeStep = fSegmentation->Dpx(dummy);
   Float_t fSddLength = fSegmentation->Dx();
   Float_t fDriftSpeed = fResponse->DriftSpeed();
   
   Float_t anodePitch = fSegmentation->Dpz(dummy);
   // map the signal
   FillMap();

   // Piergiorgio's code - do not subtract baseline since we start
   // from digits and do not duplicate arrays, i.e. use fMap instead
   // of Float_t fadc[2*fNofAnodes][fMaxNofSamples];
 

  Int_t nofFoundClusters = 0;

  Float_t **dfadc = new Float_t*[fNofMaps];
  for(i=0;i<fNofMaps;i++) dfadc[i] = new Float_t[fMaxNofSamples];
   
  
  Float_t fadc, fadc1, fadc2;
  Int_t j,k,idx,l,m;
  for(j=0;j<2;j++) {
    for(k=0;k<fNofAnodes;k++) {
      idx = j*fNofAnodes+k;
      // signal (fadc) & derivative (dfadc)
      for(l=0; l<fMaxNofSamples; l++) {
	fadc2=(Float_t)fMap->GetSignal(idx,l);
	fadc1=(Float_t)fMap->GetSignal(idx,l-1);
	if(l>0) dfadc[k][l-1] = fadc2-fadc1;
      } // samples
    } // anodes
    
    for(k=0;k<fNofAnodes;k++) {
      idx = j*fNofAnodes+k;
 
      Int_t imax = 0;
      Int_t imaxd = 0;
      Int_t it=0;
      while(it <= fMaxNofSamples-3) {
	//	cout << "sample: " << it << endl;
	
	imax = it;
	imaxd = it;
	// maximum of signal	  
	
	Float_t fadcmax = 0.;
	Float_t dfadcmax = 0.;
	Int_t lthrmina = 1;
	//	if(it >= 60) lthrmina = 2;
        //      if(it >= 100) lthrmina = 3;
	Int_t lthrmint = 2;
	//if(it >= 60) lthrmint = 3;
        //if(it >= 100) lthrmint = 4;

	Int_t lthra = 0;
	Int_t lthrt = 0;
	for(m=0;m<10;m++) {
	  Int_t id = it+m;
	  if(id>=fMaxNofSamples) break;
          fadc=fMap->GetSignal(idx,id);
	  if(fadc > fadcmax) {
	    fadcmax = fadc;
	    if(fadc > fCutAmplitude) { lthra++; lthrt++; }
	    imax = id;
	  }
	  if(dfadc[k][id] > dfadcmax) {
	    dfadcmax = dfadc[k][id];
	    imaxd = id;
	  }
	}
	it = imaxd;
	// skip if no signal over threshold
	if(fMap->GetSignal(idx,imax) < fCutAmplitude) {it++; continue;}

	if(k>0) {
	  if(fMap->GetSignal(idx-1,imax) > fCutAmplitude) lthra++;
	}
	if(k<fNofAnodes-1)
	  if(fMap->GetSignal(idx+1,imax) > fCutAmplitude) lthra++;

	if(imax>0) {
	  if(fMap->GetSignal(idx,imax-1) > fCutAmplitude) lthrt++;
	}
	if(imax<fMaxNofSamples)
	  if(fMap->GetSignal(idx,imax+1) > fCutAmplitude) lthrt++;
	
	// cluster charge
	Int_t tstart = it-1;
	
	Bool_t ilcl = 0;
	if(lthrt >= lthrmint && lthra >= lthrmina) ilcl = 1;
	
	if(ilcl) {
	  nofFoundClusters++;
	  Int_t tstop = tstart;
	  Float_t dfadcmin = 10000.;
          Int_t ij;
	  for(ij=0; ij<10; ij++) {
	    if(dfadc[k][it+ij] < dfadcmin) {
	      tstop = it+ij+1;
	      dfadcmin = dfadc[k][it+ij];
	    }
	  }
	  
	  Float_t clusterCharge = 0.;
	  Float_t clusterAnode = k+0.5;
	  Float_t clusterTime = 0.;
	  Float_t clusterMult = 0.;
	  Float_t clusterPeakAmplitude = 0.;
          Int_t its;
	  for(its=tstart; its<=tstop; its++) {
            fadc=fMap->GetSignal(idx,its);
	    clusterCharge += fadc;
	    if(fadc > clusterPeakAmplitude) clusterPeakAmplitude = fadc;
	    clusterTime += fadc*its;
	    clusterMult++;
	    if(its == tstop) {
	       clusterTime /= (clusterCharge/fTimeStep);   // ns
	       clusterCharge *= (fTimeStep/160.);          // keV
	       if(clusterTime > 58.2) clusterTime -= 58.2;   // ns
	      /*
	      else {
		cout << "Warning: cluster with negative time " << clusterTime << ", peak ampl.: " << clusterPeakAmplitude << ", mult: " << clusterMult << ", charge: " << clusterCharge << endl;
		cout << "Anode: " << k << ", tstart: " << tstart << ", tstop: " << tstop << ", Charge: " << clusterCharge << endl;
	      }
	      */
	      }
	  }
	  // cout << "Anode: " << k << ", tstart: " << tstart << ", tstop: " << tstop << ", Charge: " << clusterCharge << endl;

	  Float_t clusteranodePath = (0.06 + clusterAnode - fNofAnodes/2)*anodePitch;
	  Float_t clusterDriftPath = clusterTime*fDriftSpeed;
	  clusterDriftPath = fSddLength-clusterDriftPath;

          if(clusterCharge < 0.) break;

	  //printf("wing clusterMult clusterAnode clusterTime %d %f %f %f \n",j+1,clusterMult, clusterAnode, clusterTime);

	  AliITSRawClusterSDD *clust = new AliITSRawClusterSDD(j+1,clusterAnode,clusterTime,clusterCharge,clusterPeakAmplitude,0.,0.,clusterDriftPath,clusteranodePath,clusterMult);
	  //fClusters->Add(point);
	  iTS->AddCluster(1,clust);
	  it = tstop;
	} // ilcl
	
	it++;
	
      } // while (samples)
    } // anodes
  } // detectors (2)

  Int_t nofClusters = fClusters->GetEntriesFast();
  nofClusters -= fNclusters;

  //printf("SDD- Find1Dclust: fNclusters nofClusters %d %d \n",fNclusters, nofClusters);

  fMap->ClearMap();
  
  for(i=0;i<fNofMaps;i++) delete[] dfadc[i];
  delete [] dfadc;

  return;

}

//_____________________________________________________________________________
void  AliITSClusterFinderSDD::GroupClusters()
{
  // group clusters
  Int_t dummy=0;
  Float_t fTimeStep = fSegmentation->Dpx(dummy);


  // get number of clusters for this module
  Int_t nofClusters = fClusters->GetEntriesFast();
  nofClusters -= fNclusters;

  //printf("SDD- GroupClusters: fNclusters nofClusters %d %d \n",fNclusters, nofClusters);

  AliITSRawClusterSDD *clusterI;
  AliITSRawClusterSDD *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 = (AliITSRawClusterSDD*) fClusters->At(i);
      clusterJ = (AliITSRawClusterSDD*) fClusters->At(j);
      // 1.3 good
      if(clusterI->T() < fTimeStep*60) fDAnode = 3.2;
      if(clusterI->T() < fTimeStep*10) fDAnode = 1.2;
      Bool_t pair = clusterI->Brother(clusterJ,fDAnode,fDTime);
      if(!pair) continue;
      //      clusterI->Print();
      //      clusterJ->Print();
      clusterI->Add(clusterJ);
      label[j] = 1;
      fClusters->RemoveAt(j);
    } // J clusters  
    label[i] = 1;
  } // I clusters
  fClusters->Compress();
  
  delete [] label;
  return;

}

//_____________________________________________________________________________

void AliITSClusterFinderSDD::SelectClusters()
{
  // get number of clusters for this module
  Int_t nofClusters = fClusters->GetEntriesFast();
  nofClusters -= fNclusters;

  //printf("SDD- SelectClusters: fNclusters nofClusters %d %d \n",fNclusters, nofClusters);

  Int_t i;
  for(i=0; i<nofClusters; i++) { 
    AliITSRawClusterSDD *clusterI = (AliITSRawClusterSDD*) fClusters->At(i);
    Int_t rmflg = 0;
    Float_t wy = 0.;
    if(clusterI->Anodes() != 0.) {
      wy = ((Float_t) clusterI->Samples())/clusterI->Anodes();
    }
    Float_t amp = clusterI->PeakAmpl();
    if(amp < fMinPeak) rmflg = 1;  
    if(wy < fMinNCells) rmflg = 1;
    if(rmflg) fClusters->RemoveAt(i);
  } // I clusters
  fClusters->Compress();
  return;

}

//_____________________________________________________________________________

void AliITSClusterFinderSDD::GetRecPoints()
{
  // get rec points
  //static Int_t counter=0;

  AliITS *iTS=(AliITS*)gAlice->GetModule("ITS");

  // get number of clusters for this module
  Int_t nofClusters = fClusters->GetEntriesFast();
  nofClusters -= fNclusters;

  //printf("SDD- GetRecPoints: fNclusters nofClusters %d %d \n",fNclusters, nofClusters);

//  const Float_t kdEdXtoQ = 2.778e+2; // KeV -> number of e-hole pairs in Si
  const Float_t kconvGeV = 1.e-6; // GeV -> KeV
  const Float_t kconv = 1.0e-4; 
  const Float_t kRMSx = 38.0*kconv; // microns->cm ITS TDR Table 1.3
  const Float_t kRMSz = 28.0*kconv; // microns->cm ITS TDR Table 1.3

  Int_t i;
  for(i=0; i<nofClusters; i++) { 
    AliITSRawClusterSDD *clusterI = (AliITSRawClusterSDD*) fClusters->At(i);
    AliITSRecPoint rnew;
    rnew.SetX(clusterI->X());
    rnew.SetZ(clusterI->Z());
    rnew.SetQ(clusterI->Q());   // in KeV - should be ADC
    //rnew.SetdEdX((clusterI->Q())/kdEdXtoQ);
    rnew.SetdEdX(kconvGeV*clusterI->Q());
    rnew.SetSigmaX2(kRMSx*kRMSx);
    rnew.SetSigmaZ2(kRMSz*kRMSz);
    rnew.SetProbability(1.);
    iTS->AddRecPoint(rnew);
   //counter++;
  } // I clusters
  //printf("counter %d\n",counter);

}

//_____________________________________________________________________________

void AliITSClusterFinderSDD::FindRawClusters()
{
  // find raw clusters
    Find1DClusters();
    GroupClusters();
    SelectClusters();
    GetRecPoints();
}
Commit	Line	Data
b0f5e3fc	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	#include "AliITSClusterFinderSDD.h"
	18	#include "AliRun.h"
	19
	20
	21
	22	ClassImp(AliITSClusterFinderSDD)
	23
	24	//----------------------------------------------------------
	25	AliITSClusterFinderSDD::AliITSClusterFinderSDD
	26	(AliITSsegmentation seg, AliITSresponse response, TClonesArray digits, TClonesArray recp)
	27	{
	28	// constructor
	29	fSegmentation=seg;
	30	fResponse=response;
	31	fDigits=digits;
	32	fClusters=recp;
	33	fNclusters= fClusters->GetEntriesFast();
	34	printf("SDD: fNclusters %d\n",fNclusters);
	35	SetCutAmplitude();
	36	SetDAnode();
	37	SetDTime();
	38	SetMap();
	39	SetMinPeak();
	40	SetNCells();
	41	}
	42
	43	//_____________________________________________________________________________
	44	AliITSClusterFinderSDD::AliITSClusterFinderSDD()
	45	{
	46	// constructor
	47	fSegmentation=0;
	48	fResponse=0;
	49	fDigits=0;
	50	fClusters=0;
	51	fNclusters=0;
	52	SetCutAmplitude();
	53	SetDAnode();
	54	SetDTime();
	55	SetMap();
	56	SetMinPeak();
	57	SetNCells();
	58
	59	}
	60
	61	//__________________________________________________________________________
	62	AliITSClusterFinderSDD::AliITSClusterFinderSDD(const AliITSClusterFinderSDD &source){
	63	// Copy Constructor
	64	if(&source == this) return;
65	this->fClusters = source.fClusters ;
66	this->fNclusters = source.fNclusters ;
67	this->fMap = source.fMap ;
68	this->fCutAmplitude = source.fCutAmplitude ;
69	this->fDAnode = source.fDAnode ;
70	this->fDTime = source.fDTime ;
71	this->fMinPeak = source.fMinPeak ;
72	this->fMinNCells = source.fMinNCells ;
73	return;
74	}
75
76	//_________________________________________________________________________
77	AliITSClusterFinderSDD&
78	AliITSClusterFinderSDD::operator=(const AliITSClusterFinderSDD &source) {
79	// Assignment operator
80	if(&source == this) return *this;
81	this->fClusters = source.fClusters ;
82	this->fNclusters = source.fNclusters ;
83	this->fMap = source.fMap ;
84	this->fCutAmplitude = source.fCutAmplitude ;
85	this->fDAnode = source.fDAnode ;
86	this->fDTime = source.fDTime ;
87	this->fMinPeak = source.fMinPeak ;
88	this->fMinNCells = source.fMinNCells ;
89	return *this;
90	}
91
92	//_____________________________________________________________________________
93	void AliITSClusterFinderSDD::SetMap()
94	{
95	// set map
96	if(!fMap) fMap=new AliITSMapA2(fSegmentation);
97
98	}
99	//_____________________________________________________________________________
100	void AliITSClusterFinderSDD::FillMap()
101	{
102	// fCoord1 = anode #
103	// fCoord2 = time sample
104
105	if (!fDigits) return;
106
107	Int_t ndigits = fDigits->GetEntriesFast();
108	//printf("FillMap: ndigits %d\n",ndigits);
109	if (!ndigits) return;
110
111	AliITSdigitSDD *dig;
112	Int_t ndig;
113	for(ndig=0; ndig<ndigits; ndig++) {
114	dig = (AliITSdigitSDD*)fDigits->UncheckedAt(ndig);
115	Double_t signal=dig->fSignal;
116	//printf("FillMap: ndig fCoord1 fCoord2 signal %d %d %d %f\n",ndig,dig->fCoord1,dig->fCoord2,signal);
117	fMap->SetHit(dig->fCoord1,dig->fCoord2,signal);
118	}
119
120
121	}
122	//_____________________________________________________________________________
123
124	void AliITSClusterFinderSDD::Find1DClusters()
125	{
126	// find 1D clusters
127
128	AliITS iTS=(AliITS)gAlice->GetModule("ITS");
129
130	// retrieve the parameters
44b3710f	131	Int_t i;
b0f5e3fc	132	Int_t fNofMaps = fSegmentation->Npz();
	133	Int_t fMaxNofSamples = fSegmentation->Npx();
	134	Int_t fNofAnodes = fNofMaps/2;
	135	Int_t dummy=0;
	136	Float_t fTimeStep = fSegmentation->Dpx(dummy);
	137	Float_t fSddLength = fSegmentation->Dx();
	138	Float_t fDriftSpeed = fResponse->DriftSpeed();
	139
	140	Float_t anodePitch = fSegmentation->Dpz(dummy);
	141	// map the signal
	142	FillMap();
	143
	144	// Piergiorgio's code - do not subtract baseline since we start
	145	// from digits and do not duplicate arrays, i.e. use fMap instead
	146	// of Float_t fadc[2*fNofAnodes][fMaxNofSamples];
	147
	148
	149	Int_t nofFoundClusters = 0;
44b3710f	150
	151	Float_t *dfadc = new Float_t[fNofMaps];
	152	for(i=0;i<fNofMaps;i++) dfadc[i] = new Float_t[fMaxNofSamples];
	153
	154
b0f5e3fc	155	Float_t fadc, fadc1, fadc2;
	156	Int_t j,k,idx,l,m;
	157	for(j=0;j<2;j++) {
	158	for(k=0;k<fNofAnodes;k++) {
	159	idx = j*fNofAnodes+k;
	160	// signal (fadc) & derivative (dfadc)
	161	for(l=0; l<fMaxNofSamples; l++) {
	162	fadc2=(Float_t)fMap->GetSignal(idx,l);
	163	fadc1=(Float_t)fMap->GetSignal(idx,l-1);
	164	if(l>0) dfadc[k][l-1] = fadc2-fadc1;
	165	} // samples
	166	} // anodes
	167
	168	for(k=0;k<fNofAnodes;k++) {
	169	idx = j*fNofAnodes+k;
	170
	171	Int_t imax = 0;
	172	Int_t imaxd = 0;
	173	Int_t it=0;
	174	while(it <= fMaxNofSamples-3) {
	175	// cout << "sample: " << it << endl;
	176
	177	imax = it;
	178	imaxd = it;
	179	// maximum of signal
	180
	181	Float_t fadcmax = 0.;
	182	Float_t dfadcmax = 0.;
	183	Int_t lthrmina = 1;
	184	// if(it >= 60) lthrmina = 2;
	185	// if(it >= 100) lthrmina = 3;
	186	Int_t lthrmint = 2;
	187	//if(it >= 60) lthrmint = 3;
	188	//if(it >= 100) lthrmint = 4;
	189
	190	Int_t lthra = 0;
	191	Int_t lthrt = 0;
	192	for(m=0;m<10;m++) {
	193	Int_t id = it+m;
	194	if(id>=fMaxNofSamples) break;
	195	fadc=fMap->GetSignal(idx,id);
	196	if(fadc > fadcmax) {
	197	fadcmax = fadc;
	198	if(fadc > fCutAmplitude) { lthra++; lthrt++; }
	199	imax = id;
	200	}
	201	if(dfadc[k][id] > dfadcmax) {
	202	dfadcmax = dfadc[k][id];
	203	imaxd = id;
	204	}
	205	}
	206	it = imaxd;
	207	// skip if no signal over threshold
	208	if(fMap->GetSignal(idx,imax) < fCutAmplitude) {it++; continue;}
	209
	210	if(k>0) {
	211	if(fMap->GetSignal(idx-1,imax) > fCutAmplitude) lthra++;
	212	}
	213	if(k<fNofAnodes-1)
	214	if(fMap->GetSignal(idx+1,imax) > fCutAmplitude) lthra++;
	215
	216	if(imax>0) {
	217	if(fMap->GetSignal(idx,imax-1) > fCutAmplitude) lthrt++;
	218	}
219	if(imax<fMaxNofSamples)
220	if(fMap->GetSignal(idx,imax+1) > fCutAmplitude) lthrt++;
221
222	// cluster charge
223	Int_t tstart = it-1;
224
225	Bool_t ilcl = 0;
226	if(lthrt >= lthrmint && lthra >= lthrmina) ilcl = 1;
227
228	if(ilcl) {
229	nofFoundClusters++;
230	Int_t tstop = tstart;
231	Float_t dfadcmin = 10000.;
232	Int_t ij;
233	for(ij=0; ij<10; ij++) {
234	if(dfadc[k][it+ij] < dfadcmin) {
235	tstop = it+ij+1;
236	dfadcmin = dfadc[k][it+ij];
237	}
238	}
239
240	Float_t clusterCharge = 0.;
241	Float_t clusterAnode = k+0.5;
242	Float_t clusterTime = 0.;
243	Float_t clusterMult = 0.;
244	Float_t clusterPeakAmplitude = 0.;
245	Int_t its;
246	for(its=tstart; its<=tstop; its++) {
247	fadc=fMap->GetSignal(idx,its);
248	clusterCharge += fadc;
249	if(fadc > clusterPeakAmplitude) clusterPeakAmplitude = fadc;
250	clusterTime += fadc*its;
251	clusterMult++;
252	if(its == tstop) {
253	clusterTime /= (clusterCharge/fTimeStep); // ns
254	clusterCharge *= (fTimeStep/160.); // keV
255	if(clusterTime > 58.2) clusterTime -= 58.2; // ns
256	/*
257	else {
258	cout << "Warning: cluster with negative time " << clusterTime << ", peak ampl.: " << clusterPeakAmplitude << ", mult: " << clusterMult << ", charge: " << clusterCharge << endl;
259	cout << "Anode: " << k << ", tstart: " << tstart << ", tstop: " << tstop << ", Charge: " << clusterCharge << endl;
260	}
261	*/
262	}
263	}
264	// cout << "Anode: " << k << ", tstart: " << tstart << ", tstop: " << tstop << ", Charge: " << clusterCharge << endl;
265
266	Float_t clusteranodePath = (0.06 + clusterAnode - fNofAnodes/2)*anodePitch;
267	Float_t clusterDriftPath = clusterTime*fDriftSpeed;
268	clusterDriftPath = fSddLength-clusterDriftPath;
269
270	if(clusterCharge < 0.) break;
271
272	//printf("wing clusterMult clusterAnode clusterTime %d %f %f %f \n",j+1,clusterMult, clusterAnode, clusterTime);
273
274	AliITSRawClusterSDD *clust = new AliITSRawClusterSDD(j+1,clusterAnode,clusterTime,clusterCharge,clusterPeakAmplitude,0.,0.,clusterDriftPath,clusteranodePath,clusterMult);
275	//fClusters->Add(point);
276	iTS->AddCluster(1,clust);
277	it = tstop;
278	} // ilcl
279
280	it++;
281
282	} // while (samples)
283	} // anodes
284	} // detectors (2)
285
286	Int_t nofClusters = fClusters->GetEntriesFast();
287	nofClusters -= fNclusters;
288
289	//printf("SDD- Find1Dclust: fNclusters nofClusters %d %d \n",fNclusters, nofClusters);
290
291	fMap->ClearMap();
44b3710f	292
	293	for(i=0;i<fNofMaps;i++) delete[] dfadc[i];
	294	delete [] dfadc;
b0f5e3fc	295
	296	return;
	297
	298	}
	299
	300	//_____________________________________________________________________________
	301	void AliITSClusterFinderSDD::GroupClusters()
	302	{
	303	// group clusters
	304	Int_t dummy=0;
	305	Float_t fTimeStep = fSegmentation->Dpx(dummy);
	306
	307
	308	// get number of clusters for this module
	309	Int_t nofClusters = fClusters->GetEntriesFast();
	310	nofClusters -= fNclusters;
	311
	312	//printf("SDD- GroupClusters: fNclusters nofClusters %d %d \n",fNclusters, nofClusters);
	313
	314	AliITSRawClusterSDD *clusterI;
	315	AliITSRawClusterSDD *clusterJ;
	316
44b3710f	317	Int_t *label = new Int_t [nofClusters];
b0f5e3fc	318	Int_t i,j;
	319	for(i=0; i<nofClusters; i++) label[i] = 0;
	320	for(i=0; i<nofClusters; i++) {
	321	if(label[i] != 0) continue;
	322	for(j=i+1; j<nofClusters; j++) {
	323	if(label[j] != 0) continue;
	324	clusterI = (AliITSRawClusterSDD*) fClusters->At(i);
	325	clusterJ = (AliITSRawClusterSDD*) fClusters->At(j);
	326	// 1.3 good
	327	if(clusterI->T() < fTimeStep*60) fDAnode = 3.2;
	328	if(clusterI->T() < fTimeStep*10) fDAnode = 1.2;
	329	Bool_t pair = clusterI->Brother(clusterJ,fDAnode,fDTime);
	330	if(!pair) continue;
	331	// clusterI->Print();
	332	// clusterJ->Print();
	333	clusterI->Add(clusterJ);
	334	label[j] = 1;
	335	fClusters->RemoveAt(j);
	336	} // J clusters
	337	label[i] = 1;
	338	} // I clusters
	339	fClusters->Compress();
44b3710f	340
44b3710f	341	delete [] label;
b0f5e3fc	342	return;
	343
	344	}
	345
	346	//_____________________________________________________________________________
	347
	348	void AliITSClusterFinderSDD::SelectClusters()
	349	{
	350	// get number of clusters for this module
	351	Int_t nofClusters = fClusters->GetEntriesFast();
	352	nofClusters -= fNclusters;
	353
	354	//printf("SDD- SelectClusters: fNclusters nofClusters %d %d \n",fNclusters, nofClusters);
	355
	356	Int_t i;
	357	for(i=0; i<nofClusters; i++) {
	358	AliITSRawClusterSDD clusterI = (AliITSRawClusterSDD) fClusters->At(i);
	359	Int_t rmflg = 0;
	360	Float_t wy = 0.;
	361	if(clusterI->Anodes() != 0.) {
	362	wy = ((Float_t) clusterI->Samples())/clusterI->Anodes();
	363	}
	364	Float_t amp = clusterI->PeakAmpl();
	365	if(amp < fMinPeak) rmflg = 1;
	366	if(wy < fMinNCells) rmflg = 1;
	367	if(rmflg) fClusters->RemoveAt(i);
	368	} // I clusters
	369	fClusters->Compress();
	370	return;
	371
	372	}
	373
	374	//_____________________________________________________________________________
	375
	376	void AliITSClusterFinderSDD::GetRecPoints()
	377	{
	378	// get rec points
	379	//static Int_t counter=0;
	380
	381	AliITS iTS=(AliITS)gAlice->GetModule("ITS");
	382
	383	// get number of clusters for this module
	384	Int_t nofClusters = fClusters->GetEntriesFast();
	385	nofClusters -= fNclusters;
	386
	387	//printf("SDD- GetRecPoints: fNclusters nofClusters %d %d \n",fNclusters, nofClusters);
	388
	389	// const Float_t kdEdXtoQ = 2.778e+2; // KeV -> number of e-hole pairs in Si
	390	const Float_t kconvGeV = 1.e-6; // GeV -> KeV
	391	const Float_t kconv = 1.0e-4;
	392	const Float_t kRMSx = 38.0*kconv; // microns->cm ITS TDR Table 1.3
	393	const Float_t kRMSz = 28.0*kconv; // microns->cm ITS TDR Table 1.3
	394
	395	Int_t i;
	396	for(i=0; i<nofClusters; i++) {
	397	AliITSRawClusterSDD clusterI = (AliITSRawClusterSDD) fClusters->At(i);
	398	AliITSRecPoint rnew;
	399	rnew.SetX(clusterI->X());
	400	rnew.SetZ(clusterI->Z());
	401	rnew.SetQ(clusterI->Q()); // in KeV - should be ADC
	402	//rnew.SetdEdX((clusterI->Q())/kdEdXtoQ);
	403	rnew.SetdEdX(kconvGeV*clusterI->Q());
	404	rnew.SetSigmaX2(kRMSx*kRMSx);
	405	rnew.SetSigmaZ2(kRMSz*kRMSz);
406	rnew.SetProbability(1.);
407	iTS->AddRecPoint(rnew);
408	//counter++;
409	} // I clusters
410	//printf("counter %d\n",counter);
411
412	}
413
414	//_____________________________________________________________________________
415
416	void AliITSClusterFinderSDD::FindRawClusters()
417	{
418	// find raw clusters
419	Find1DClusters();
420	GroupClusters();
421	SelectClusters();
422	GetRecPoints();
423	}
424
425
426
427