+Double_t AliITSClusterFinderSDD::ChiSqr( Int_t xdim, Int_t zdim, Double_t *spe,
+ Double_t *speFit ) const{
+ // EVALUATES UNNORMALIZED CHI-SQUARED
+ Double_t chi2 = 0.;
+ for( Int_t z=0; z<zdim; z++ ){
+ for( Int_t x=1; x<xdim-1; x++ ){
+ Int_t index = x*zdim+z;
+ Double_t tmp = spe[index] - speFit[index];
+ chi2 += tmp*tmp;
+ } // end for x
+ } // end for z
+ return( chi2 );
+}
+//_______________________________________________________________________
+void AliITSClusterFinderSDD::Minim( Int_t xdim, Int_t zdim, Double_t *param,
+ Double_t *prm0,Double_t *steprm,
+ Double_t *chisqr,Double_t *spe,
+ Double_t *speFit ){
+ //
+ Int_t k, nnn, mmm, i;
+ Double_t p1, delta, d1, chisq1, p2, chisq2, t, p3, chisq3, a, b, p0, chisqt;
+ const Int_t knParam = 5;
+ Int_t npeak = (Int_t)param[0];
+ for( k=1; k<(npeak*knParam+1); k++ ) prm0[k] = param[k];
+ for( k=1; k<(npeak*knParam+1); k++ ){
+ p1 = param[k];
+ delta = steprm[k];
+ d1 = delta;
+ // ENSURE THAT STEP SIZE IS SENSIBLY LARGER THAN MACHINE ROUND OFF
+ if( TMath::Abs( p1 ) > 1.0E-6 )
+ if ( TMath::Abs( delta/p1 ) < 1.0E-4 ) delta = p1/1000;
+ else delta = (Double_t)1.0E-4;
+ // EVALUATE CHI-SQUARED AT FIRST TWO SEARCH POINTS
+ PeakFunc( xdim, zdim, param, speFit );
+ chisq1 = ChiSqr( xdim, zdim, spe, speFit );
+ p2 = p1+delta;
+ param[k] = p2;
+ PeakFunc( xdim, zdim, param, speFit );
+ chisq2 = ChiSqr( xdim, zdim, spe, speFit );
+ if( chisq1 < chisq2 ){
+ // REVERSE DIRECTION OF SEARCH IF CHI-SQUARED IS INCREASING
+ delta = -delta;
+ t = p1;
+ p1 = p2;
+ p2 = t;
+ t = chisq1;
+ chisq1 = chisq2;
+ chisq2 = t;
+ } // end if
+ i = 1; nnn = 0;
+ do { // INCREMENT param(K) UNTIL CHI-SQUARED STARTS TO INCREASE
+ nnn++;
+ p3 = p2 + delta;
+ mmm = nnn - (nnn/5)*5; // multiplo de 5
+ if( mmm == 0 ){
+ d1 = delta;
+ // INCREASE STEP SIZE IF STEPPING TOWARDS MINIMUM IS TOO SLOW
+ delta *= 5;
+ } // end if
+ param[k] = p3;
+ // Constrain paramiters
+ Int_t kpos = (k-1) % knParam;
+ switch( kpos ){
+ case 0 :
+ if( param[k] <= 20 ) param[k] = fMinPeak;
+ break;
+ case 1 :
+ if( TMath::Abs( param[k] - prm0[k] ) > 1.5 ) param[k] = prm0[k];
+ break;
+ case 2 :
+ if( TMath::Abs( param[k] - prm0[k] ) > 1. ) param[k] = prm0[k];
+ break;
+ case 3 :
+ if( param[k] < .5 ) param[k] = .5;
+ break;
+ case 4 :
+ if( param[k] < .288 ) param[k] = .288;// 1/sqrt(12) = 0.288
+ if( param[k] > zdim*.5 ) param[k] = zdim*.5;
+ break;
+ }; // end switch
+ PeakFunc( xdim, zdim, param, speFit );
+ chisq3 = ChiSqr( xdim, zdim, spe, speFit );
+ if( chisq3 < chisq2 && nnn < 50 ){
+ p1 = p2;
+ p2 = p3;
+ chisq1 = chisq2;
+ chisq2 = chisq3;
+ }else i=0;
+ } while( i );
+ // FIND MINIMUM OF PARABOLA DEFINED BY LAST THREE POINTS
+ a = chisq1*(p2-p3)+chisq2*(p3-p1)+chisq3*(p1-p2);
+ b = chisq1*(p2*p2-p3*p3)+chisq2*(p3*p3-p1*p1)+chisq3*(p1*p1-p2*p2);
+ if( a!=0 ) p0 = (Double_t)(0.5*b/a);
+ else p0 = 10000;
+ //--IN CASE OF NEARLY EQUAL CHI-SQUARED AND TOO SMALL STEP SIZE PREVENT
+ // ERRONEOUS EVALUATION OF PARABOLA MINIMUM
+ //---NEXT TWO LINES CAN BE OMITTED FOR HIGHER PRECISION MACHINES
+ //dp = (Double_t) max (TMath::Abs(p3-p2), TMath::Abs(p2-p1));
+ //if( TMath::Abs( p2-p0 ) > dp ) p0 = p2;
+ param[k] = p0;
+ // Constrain paramiters
+ Int_t kpos = (k-1) % knParam;
+ switch( kpos ){
+ case 0 :
+ if( param[k] <= 20 ) param[k] = fMinPeak;
+ break;
+ case 1 :
+ if( TMath::Abs( param[k] - prm0[k] ) > 1.5 ) param[k] = prm0[k];
+ break;
+ case 2 :
+ if( TMath::Abs( param[k] - prm0[k] ) > 1. ) param[k] = prm0[k];
+ break;
+ case 3 :
+ if( param[k] < .5 ) param[k] = .5;
+ break;
+ case 4 :
+ if( param[k] < .288 ) param[k] = .288; // 1/sqrt(12) = 0.288
+ if( param[k] > zdim*.5 ) param[k] = zdim*.5;
+ break;
+ }; // end switch
+ PeakFunc( xdim, zdim, param, speFit );
+ chisqt = ChiSqr( xdim, zdim, spe, speFit );
+ // DO NOT ALLOW ERRONEOUS INTERPOLATION
+ if( chisqt <= *chisqr ) *chisqr = chisqt;
+ else param[k] = prm0[k];
+ // OPTIMIZE SEARCH STEP FOR EVENTUAL NEXT CALL OF MINIM
+ steprm[k] = (param[k]-prm0[k])/5;
+ if( steprm[k] >= d1 ) steprm[k] = d1/5;
+ } // end for k
+ // EVALUATE FIT AND CHI-SQUARED FOR OPTIMIZED PARAMETERS
+ PeakFunc( xdim, zdim, param, speFit );
+ *chisqr = ChiSqr( xdim, zdim, spe, speFit );
+ return;