X-Git-Url: http://git.uio.no/git/?a=blobdiff_plain;f=EMCAL%2FAliCaloRawAnalyzer.cxx;h=50bd8ff2b3caf7b5f39d01bab56d5e84b9d743b2;hb=4df4ca8c8cac78a926f02beba5d596afb7986cda;hp=6f5ea4e2f0e9c3a1749ff79b49abce0a5f27ce61;hpb=e7acbc379ffd4fea420f28f4029615233624a18f;p=u%2Fmrichter%2FAliRoot.git diff --git a/EMCAL/AliCaloRawAnalyzer.cxx b/EMCAL/AliCaloRawAnalyzer.cxx index 6f5ea4e2f0e..50bd8ff2b3c 100644 --- a/EMCAL/AliCaloRawAnalyzer.cxx +++ b/EMCAL/AliCaloRawAnalyzer.cxx @@ -1,3 +1,4 @@ +// -*- mode: c++ -*- /************************************************************************** * This file is property of and copyright by * * the Relativistic Heavy Ion Group (RHIG), Yale University, US, 2009 * @@ -32,17 +33,19 @@ #include using namespace std; -//ClassImp(AliCaloRawAnalyzer) +ClassImp(AliCaloRawAnalyzer) AliCaloRawAnalyzer::AliCaloRawAnalyzer(const char *name, const char *nameshort) : TObject(), - fMinTimeIndex(-1), - fMaxTimeIndex(-1), - fFitArrayCut(5), - fAmpCut(4), - fNsampleCut(5), - fNsamplePed(3), - fIsZerosupressed( false ), - fVerbose( false ) + fMinTimeIndex(-1), + fMaxTimeIndex(-1), + fFitArrayCut(5), + fAmpCut(4), + fNsampleCut(5), + fOverflowCut(950), + fNsamplePed(3), + fIsZerosupressed( false ), + fVerbose( false ), + fAlgo(Algo::kNONE) { //Comment sprintf(fName, "%s", name); @@ -101,25 +104,42 @@ AliCaloRawAnalyzer::SelectSubarray( const Double_t *fData, const int length, con //Until the ADC value is less that fFitArrayCut, or derivative changes sign (data jump) int tmpfirst = maxindex; int tmplast = maxindex; - Double_t valFirst = fData[maxindex]; - Double_t valLast = fData[maxindex]; - - while( (tmpfirst >= 0) && (fData[tmpfirst] >= fFitArrayCut) && - (fData[tmpfirst]= 0) && (fData[tmpfirst] >= fFitArrayCut) && (!firstJump) ) { - valFirst = fData[tmpfirst]; + // jump check: + if (tmpfirst != maxindex) { // neighbor to maxindex can share peak with maxindex + if (fData[tmpfirst] >= prevFirst) { + firstJump = true; + } + } + prevFirst = fData[tmpfirst]; tmpfirst -- ; } - while( (tmplast < length) && (fData[tmplast] >= fFitArrayCut) && - (fData[tmplast]= fFitArrayCut) && (!lastJump) ) { - valLast = fData[tmplast]; + // jump check: + if (tmplast != maxindex) { // neighbor to maxindex can share peak with maxindex + if (fData[tmplast] >= prevLast) { + lastJump = true; + } + } + prevLast = fData[tmplast]; tmplast ++; } - *first = tmpfirst +1; - *last = tmplast -1; + // we keep one pre- or post- sample if we can (as in online) + // check though if we ended up on a 'jump', or out of bounds: if so, back up + if (firstJump || tmpfirst<0) tmpfirst ++; + if (lastJump || tmplast>=length) tmplast --; + + *first = tmpfirst; + *last = tmplast; return; } @@ -296,14 +316,106 @@ AliCaloRawAnalyzer::PrintBunch( const AliCaloBunchInfo &bunch ) const } -AliCaloFitResults -AliCaloRawAnalyzer::Evaluate( const vector &/*bunchvector*/, const UInt_t /*altrocfg1*/, const UInt_t /*altrocfg2*/) -{ // method to do the selection of what should possibly be fitted - // not implemented for base class - return AliCaloFitResults( 0, 0 ); +Double_t +AliCaloRawAnalyzer::CalculateChi2(const Double_t amp, const Double_t time, + const Int_t first, const Int_t last, + const Double_t adcErr, + const Double_t tau) +{ + // Input: + // amp - max amplitude; + // time - time of max amplitude; + // first, last - sample array indices to be used + // adcErr - nominal error of amplitude measurement (one value for all channels) + // if adcErr<0 that mean adcErr=1. + // tau - filter time response (in timebin units) + // Output: + // chi2 - chi2 + + if (first == last || first<0 ) { // signal consists of single sample, chi2 estimate (0) not too well defined.. + // or, first is negative, the indices are not valid + return Ret::kDummy; + } + + int nsamples = last - first + 1; + // printf(" AliCaloRawAnalyzer::CalculateChi2 : first %i last %i : nsamples %i : amp %3.2f time %3.2f \n", first, last, nsamples, amp, time); + + Int_t x = 0; + Double_t chi2 = 0; + Double_t dy = 0.0, xx = 0.0, f=0.0; + + for (Int_t i=0; i 0) { + f = amp * xx*xx * TMath::Exp(2 * (1 - xx )) ; + } + dy = fReversed[x] - f; + chi2 += dy*dy; + // printf(" AliCaloRawAnalyzer::CalculateChi2 : %i : y %f -> f %f : dy %f \n", i, fReversed[first+i], f, dy); + } + + if (adcErr>0.0) { // weight chi2 + chi2 /= (adcErr*adcErr); + } + return chi2; +} + + +void +AliCaloRawAnalyzer::CalculateMeanAndRMS(const Int_t first, const Int_t last, + Double_t & mean, Double_t & rms) +{ + // Input: + // first, last - sample array indices to be used + // Output: + // mean and RMS of samples + // + // To possibly be used to differentiate good signals from bad before fitting + // + mean = Ret::kDummy; + rms = Ret::kDummy; + + if (first == last || first<0 ) { // signal consists of single sample, chi2 estimate (0) not too well defined.. + // or, first is negative, the indices are not valid + return; + } + + int nsamples = last - first + 1; + // printf(" AliCaloRawAnalyzer::CalculateMeanAndRMS : first %i last %i : nsamples %i \n", first, last, nsamples); + + int x = 0; + Double_t sampleSum = 0; // sum of samples + Double_t squaredSampleSum = 0; // sum of samples squared + + for (Int_t i=0; i &/*bunchvector*/, const UInt_t /*altrocfg1*/, const UInt_t /*altrocfg2*/) +// { // method to do the selection of what should possibly be fitted +// // not implemented for base class +// cout << __FILE__ << ":" << __LINE__ << " " << endl; + +// return AliCaloFitResults( 0, 0 ); +// } + + + int AliCaloRawAnalyzer::PreFitEvaluateSamples( const vector &bunchvector, const UInt_t altrocfg1, const UInt_t altrocfg2, Int_t & index, Float_t & maxf, short & maxamp, short & maxrev, Float_t & ped, int & first, int & last) { // method to do the selection of what should possibly be fitted @@ -312,13 +424,13 @@ AliCaloRawAnalyzer::PreFitEvaluateSamples( const vector &bunc index = SelectBunch( bunchvector, &maxampindex, &maxamp ); // select the bunch with the highest amplitude unless any time constraints is set - if( index >= 0 && maxamp > fAmpCut) // something valid was found, and non-zero amplitude + if( index >= 0 && maxamp >= fAmpCut) // something valid was found, and non-zero amplitude { // use more convenient numbering and possibly subtract pedestal ped = ReverseAndSubtractPed( &(bunchvector.at(index)), altrocfg1, altrocfg2, fReversed ); maxf = TMath::MaxElement( bunchvector.at(index).GetLength(), fReversed ); - if ( maxf > fAmpCut ) // possibly significant signal + if ( maxf >= fAmpCut ) // possibly significant signal { // select array around max to possibly be used in fit maxrev = maxampindex - bunchvector.at(index).GetStartBin();