* about the suitability of this software for any purpose. It is *
* provided "as is" without express or implied warranty. *
**************************************************************************/
-#include <iostream.h>
-#include <TFile.h>
+/*
+ $Id$
+ $Log$
+ Revision 1.29 2002/10/25 18:54:22 barbera
+ Various improvements and updates from B.S.Nilsen and T. Virgili
+
+ Revision 1.28 2002/10/22 14:45:29 alibrary
+ Introducing Riostream.h
+
+ Revision 1.27 2002/10/14 14:57:00 hristov
+ Merging the VirtualMC branch to the main development branch (HEAD)
+
+ Revision 1.23.4.2 2002/10/14 13:14:07 hristov
+ Updating VirtualMC to v3-09-02
+
+ Revision 1.26 2002/09/09 17:23:28 nilsen
+ Minor changes in support of changes to AliITSdigitS?D class'.
+
+ Revision 1.25 2002/05/10 22:29:40 nilsen
+ Change my Massimo Masera in the default constructor to bring things into
+ compliance.
+
+ Revision 1.24 2002/04/24 22:02:31 nilsen
+ New SDigits and Digits routines, and related changes, (including new
+ noise values).
+
+ */
+//
+// Cluster finder
+// for Silicon
+// Drift Detector
+//
+#include <Riostream.h>
+
#include <TMath.h>
#include <math.h>
SetCutAmplitude();
SetDAnode();
SetDTime();
- SetMinPeak();
+ SetMinPeak((Int_t)(((AliITSresponseSDD*)fResponse)->GetNoiseAfterElectronics()*5));
+ // SetMinPeak();
SetMinNCells();
SetMaxNCells();
SetTimeCorr();
fNclusters = 0;
fMap = 0;
fCutAmplitude = 0;
+ fDAnode = 0;
+ fDTime = 0;
+ fMinPeak = 0;
+ fMinNCells = 0;
+ fMaxNCells = 0;
+ fTimeCorr = 0;
+ fMinCharge = 0;
+ /*
SetDAnode();
SetDTime();
- SetMinPeak();
+ SetMinPeak((Int_t)(((AliITSresponseSDD*)fResponse)->GetNoiseAfterElectronics()*5));
SetMinNCells();
SetMaxNCells();
SetTimeCorr();
SetMinCharge();
+ */
}
//____________________________________________________________________________
AliITSClusterFinderSDD::~AliITSClusterFinderSDD(){
//______________________________________________________________________
void AliITSClusterFinderSDD::SetCutAmplitude(Float_t nsigma){
// set the signal threshold for cluster finder
- Float_t baseline,noise,noise_after_el;
+ Float_t baseline,noise,noiseAfterEl;
fResponse->GetNoiseParam(noise,baseline);
- noise_after_el = ((AliITSresponseSDD*)fResponse)->GetNoiseAfterElectronics();
- fCutAmplitude = (Int_t)((baseline + nsigma*noise_after_el));
+ noiseAfterEl = ((AliITSresponseSDD*)fResponse)->GetNoiseAfterElectronics();
+ fCutAmplitude = (Int_t)((baseline + nsigma*noiseAfterEl));
}
//______________________________________________________________________
void AliITSClusterFinderSDD::Find1DClusters(){
Int_t i,j;
// search peaks
for( Int_t z=1; z<zdim-1; z++ ){
- for( Int_t x=2; x<xdim-3; x++ ){
+ for( Int_t x=1; x<xdim-2; x++ ){
Float_t sxz = spect[x*zdim+z];
Float_t sxz1 = spect[(x+1)*zdim+z];
Float_t sxz2 = spect[(x-1)*zdim+z];
for( Int_t i=0; i<npeak; i++ ){
if( integral != 0 ) integral[i] = 0.;
Float_t sigmaA2 = par[k+4]*par[k+4]*2.;
- Float_t T2 = par[k+3]; // PASCAL
- if( electronics == 2 ) { T2 *= T2; T2 *= 2; } // OLA
+ Float_t t2 = par[k+3]; // PASCAL
+ if( electronics == 2 ) { t2 *= t2; t2 *= 2; } // OLA
for( Int_t z=0; z<zdim; z++ ){
for( Int_t x=0; x<xdim; x++ ){
Float_t z2 = (z-par[k+2])*(z-par[k+2])/sigmaA2;
Float_t x2 = 0.;
Float_t signal = 0.;
if( electronics == 1 ){ // PASCAL
- x2 = (x-par[k+1]+T2)/T2;
+ x2 = (x-par[k+1]+t2)/t2;
signal = (x2>0.) ? par[k]*x2*exp(-x2+1.-z2) :0.0; // RCCR2
// signal =(x2>0.) ? par[k]*x2*x2*exp(-2*x2+2.-z2 ):0.0;//RCCR
}else if( electronics == 2 ) { // OLA
- x2 = (x-par[k+1])*(x-par[k+1])/T2;
+ x2 = (x-par[k+1])*(x-par[k+1])/t2;
signal = par[k] * exp( -x2 - z2 );
} else {
- cout << "Wrong SDD Electronics =" << electronics << endl;
+ Warning("PeakFunc","Wrong SDD Electronics = %d",electronics);
// exit( 1 );
} // end if electronicx
spe[x*zdim+z] += signal;
}
//__________________________________________________________________________
Float_t AliITSClusterFinderSDD::ChiSqr( Int_t xdim, Int_t zdim, Float_t *spe,
- Float_t *speFit ){
+ Float_t *speFit ) const{
// EVALUATES UNNORMALIZED CHI-SQUARED
Float_t chi2 = 0.;
for( Int_t z=0; z<zdim; z++ ){
} // end if
Int_t xdim = tstop-tstart+3;
Int_t zdim = astop-astart+3;
- if(xdim > 50 || zdim > 30) { cout << "Warning: xdim: " << xdim << ", zdim: " << zdim << endl; continue; }
+ if(xdim > 50 || zdim > 30) {
+ Warning("ResolveClustersE","xdim: %d , zdim: %d ",xdim,zdim);
+ continue;
+ }
Float_t *sp = new Float_t[ xdim*zdim+1 ];
memset( sp, 0, sizeof(Float_t)*(xdim*zdim+1) );
// clusterJ->PrintInfo();
Float_t *par = new Float_t[npeak*5+1];
par[0] = (Float_t)npeak;
- // Initial paramiters in cell dimentions
+ // Initial parameters in cell dimentions
Int_t k1 = 1;
for( i=0; i<npeak; i++ ){
par[k1] = peakAmp1[i];
newiTimef *= fTimeStep;
if( newiTimef > fTimeCorr ) newiTimef -= fTimeCorr;
if( electronics == 1 ){
- newiTimef *= 0.999438; // PASCAL
- newiTimef += (6./fDriftSpeed - newiTimef/3000.);
+ // newiTimef *= 0.999438; // PASCAL
+ // newiTimef += (6./fDriftSpeed - newiTimef/3000.);
}else if( electronics == 2 )
newiTimef *= 0.99714; // OLA
Float_t driftPath = fSddLength - newiTimef * fDriftSpeed;
clusterI.SetAsigma( sigma[i]*anodePitch );
clusterI.SetTsigma( tau[i]*fTimeStep );
clusterI.SetQ( integral[i] );
- // clusterI.PrintInfo();
+ // clusterI.PrintInfo();
iTS->AddCluster( 1, &clusterI );
} // end for i
fClusters->RemoveAt( j );
delete [] par;
- } else {
- clusterJ->PrintInfo();
- cout <<" --- Peak not found!!!! minpeak=" << fMinPeak<<
- " cluster peak=" << clusterJ->PeakAmpl() <<
- " npeak=" << npeak << endl << "xdim=" << xdim-2 << " zdim="
- << zdim-2 << endl << endl;
+ } else { // something odd
+ Warning("ResolveClustersE","--- Peak not found!!!! minpeak=%d ,cluster peak= %f , module= %d",fMinPeak,clusterJ->PeakAmpl(),fModule);
+ clusterJ->PrintInfo();
+ Warning("ResolveClustersE"," xdim= %d zdim= %d",xdim-2,zdim-2);
}
delete [] sp;
} // cluster loop
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;
+ Int_t i,j;
Int_t ix, iz, idx=-1;
AliITSdigitSDD *dig=0;
Int_t ndigits=fDigits->GetEntriesFast();
rnew.SetdEdX(kconvGeV*clusterI->Q());
rnew.SetSigmaX2(kRMSx*kRMSx);
rnew.SetSigmaZ2(kRMSz*kRMSz);
- if(dig) rnew.fTracks[0]=dig->fTracks[0];
- if(dig) rnew.fTracks[1]=dig->fTracks[1];
- if(dig) rnew.fTracks[2]=dig->fTracks[2];
+ if(dig){
+ rnew.fTracks[0] = dig->fTracks[0];
+ rnew.fTracks[1] = -3;
+ rnew.fTracks[2] = -3;
+ j=1;
+ while(rnew.fTracks[0]==dig->fTracks[j] &&
+ j<dig->GetNTracks()) j++;
+ if(j<dig->GetNTracks()){
+ rnew.fTracks[1] = dig->fTracks[j];
+ while((rnew.fTracks[0]==dig->fTracks[j] ||
+ rnew.fTracks[1]==dig->fTracks[j] )&&
+ j<dig->GetNTracks()) j++;
+ if(j<dig->GetNTracks()) rnew.fTracks[2] = dig->fTracks[j];
+ } // end if
+ } // end if
//printf("SDD: i %d track1 track2 track3 %d %d %d x y %f %f\n",
// i,rnew.fTracks[0],rnew.fTracks[1],rnew.fTracks[2],c
// lusterI->X(),clusterI->Z());
GetRecPoints();
}
//_______________________________________________________________________
-void AliITSClusterFinderSDD::Print(){
+void AliITSClusterFinderSDD::Print() const{
// Print SDD cluster finder Parameters
cout << "**************************************************" << endl;