+// -*- mode: c++ -*-
/**************************************************************************
- * This file is property of and copyright by the Experimental Nuclear *
- * Physics Group, Dep. of Physics *
- * University of Oslo, Norway, 2007 *
+ * This file is property of and copyright by *
+ * the Relativistic Heavy Ion Group (RHIG), Yale University, US, 2009 *
+ * *
+ * Primary Author: Per Thomas Hille <perthomas.hille@yale.edu> *
* *
- * Author: Per Thomas Hille <perthi@fys.uio.no> for the ALICE HLT Project.*
* Contributors are mentioned in the code where appropriate. *
- * Please report bugs to perthi@fys.uio.no *
+ * Please report bugs to perthomas.hille@yale.edu *
* *
* Permission to use, copy, modify and distribute this software and its *
* documentation strictly for non-commercial purposes is hereby granted *
#include "AliCaloRawAnalyzerPeakFinder.h"
#include "AliCaloBunchInfo.h"
#include "AliCaloFitResults.h"
-#include <iostream>
-#include "unistd.h"
#include "TMath.h"
#include "AliLog.h"
+#include "AliCDBEntry.h"
+#include "AliCDBManager.h"
+#include "TFile.h"
+#include "AliCaloPeakFinderVectors.h"
+#include <iostream>
using namespace std;
-AliCaloRawAnalyzerPeakFinder::AliCaloRawAnalyzerPeakFinder() :AliCaloRawAnalyzer(),
- fTof(0),
- fAmp(0)
-{
- //comment
- fNsampleCut = 5;
- for(int i=0; i < MAXSTART; i++)
+ClassImp( AliCaloRawAnalyzerPeakFinder )
+
+
+AliCaloRawAnalyzerPeakFinder::AliCaloRawAnalyzerPeakFinder() :AliCaloRawAnalyzer("Peak-Finder", "PF"),
+ fPeakFinderVectors(0),
+ fRunOnAlien(false),
+ fIsInitialized(false)
+{
+ //Comment
+ fAlgo= Algo::kPeakFinder;
+ fPeakFinderVectors = new AliCaloPeakFinderVectors() ;
+ ResetVectors();
+ LoadVectorsOCDB();
+}
+
+void
+AliCaloRawAnalyzerPeakFinder::ResetVectors()
+{
+ //As name implies
+ for(int i=0; i < PF::MAXSTART; i++)
{
- for(int j=0; j < SAMPLERANGE; j++ )
+ for(int j=0; j < PF::SAMPLERANGE; j++ )
{
- fPFAmpVectors[i][j] = new double[100];
- fPFTofVectors[i][j] = new double[100];
-
for(int k=0; k < 100; k++ )
{
fPFAmpVectors[i][j][k] = 0;
fPFTofVectors[i][j][k] = 0;
+ fPFAmpVectorsCoarse[i][j][k] = 0;
+ fPFTofVectorsCoarse[i][j][k] = 0;
}
}
}
- LoadVectors();
}
AliCaloRawAnalyzerPeakFinder::~AliCaloRawAnalyzerPeakFinder()
{
//comment
- for(int i=0; i < MAXSTART; i++)
+}
+
+
+Double_t
+AliCaloRawAnalyzerPeakFinder::ScanCoarse(const Double_t *const array, const int length ) const
+{
+ // Fisrt (coarce) estimate of Amplitude using the Peak-Finder.
+ // The output of the first iteration is sued to select vectors
+ // for the second iteration.
+
+ Double_t tmpTof = 0;
+ Double_t tmpAmp= 0;
+
+ for(int i=0; i < length; i++)
{
- for(int j=0; j < SAMPLERANGE; j++ )
- {
- delete[] fPFAmpVectors[i][j];
- delete[] fPFTofVectors[i][j];
- }
+ tmpTof += fPFTofVectorsCoarse[0][length][i]*array[i];
+ tmpAmp += fPFAmpVectorsCoarse[0][length][i]*array[i];
}
+
+ tmpTof = tmpTof / tmpAmp ;
+ return tmpTof;
}
AliCaloFitResults
AliCaloRawAnalyzerPeakFinder::Evaluate( const vector<AliCaloBunchInfo> &bunchvector, const UInt_t altrocfg1, const UInt_t altrocfg2 )
{
+ // Evaluation of amplitude and TOF
+ if( fIsInitialized == false )
+ {
+ cout << __FILE__ << ":" << __LINE__ << "ERROR, peakfinder vectors not loaded" << endl;
+ return AliCaloFitResults(kInvalid, kInvalid);
+ }
+
// Extracting the amplitude using the Peak-Finder algorithm
// The amplitude is a weighted sum of the samples using
// optimum weights.
short maxampindex; //index of maximum amplitude
short maxamp; //Maximum amplitude
fAmp = 0;
- fAmpA[0] = 0;
- fAmpA[1] = 0;
- fAmpA[2] = 0;
-
int index = SelectBunch( bunchvector, &maxampindex, &maxamp );
-
+
if( index >= 0)
{
Float_t ped = ReverseAndSubtractPed( &(bunchvector.at(index)) , altrocfg1, altrocfg2, fReversed );
Float_t maxf = TMath::MaxElement( bunchvector.at(index).GetLength(), fReversed );
-
- if( maxf < fAmpCut || ( maxamp - ped) > 900 ) // (maxamp - ped) > 900 = Close to saturation (use low gain then)
+ short timebinOffset = maxampindex - (bunchvector.at( index ).GetLength()-1);
+
+ if( maxf < fAmpCut || ( maxamp - ped) > fOverflowCut ) // (maxamp - ped) > fOverflowCut = Close to saturation (use low gain then)
{
- // cout << __FILE__ << __LINE__ <<":, maxamp = " << maxamp << ", ped = "<< ped << ",. maxf = "<< maxf << ", maxampindex = "<< maxampindex << endl;
- return AliCaloFitResults( maxamp, ped, -1, maxf, maxampindex, -1, -1 );
- }
-
- int first;
- int last;
-
- if ( maxf > fAmpCut )
+ return AliCaloFitResults( maxamp, ped, Ret::kCrude, maxf, timebinOffset);
+ }
+ else if ( maxf >= fAmpCut )
{
- SelectSubarray( fReversed, bunchvector.at(index).GetLength(), maxampindex - bunchvector.at(index).GetStartBin(), &first, &last);
+ int first = 0;
+ int last = 0;
+ short maxrev = maxampindex - bunchvector.at(index).GetStartBin();
+ SelectSubarray( fReversed, bunchvector.at(index).GetLength(), maxrev, &first, &last, fFitArrayCut);
int nsamples = last - first;
- if( ( nsamples ) >= fNsampleCut )
+
+ if( ( nsamples ) >= fNsampleCut ) // no if statement needed really; keep for readability
{
int startbin = bunchvector.at(index).GetStartBin();
- int n = last -first;
+ int n = last - first;
int pfindex = n - fNsampleCut;
- pfindex = pfindex > SAMPLERANGE ? SAMPLERANGE : pfindex;
-
+ pfindex = pfindex > PF::SAMPLERANGE ? PF::SAMPLERANGE : pfindex;
int dt = maxampindex - startbin -2;
- for(int i=0; i < SAMPLERANGE; i++ )
+ int tmpindex = 0;
+ Float_t tmptof = ScanCoarse( &fReversed[dt] , n );
+
+ if( tmptof < -1 )
{
- // int dt = maxampindex - startbin -2;
-
- // double tmp[3];
-
- // tmp[0] = fReversed[ dt +i -1];
- // tmp[1] = fReversed[ dt +i];
- // tmp[2] = fReversed[ dt +i +1];
-
- for(int j = 0; j < 3; j++)
- {
- // fAmpA[j] += fPFAmpVectors[0][pfindex][i]*tmp[j];
- fAmpA[j] += fPFAmpVectors[0][pfindex][i]*fReversed[ dt +i +j -1];
-
- }
+ tmpindex = 0;
}
-
- double diff = 9999;
-
- int tmpindex = 0;
+ else
+ if( tmptof > -1 && tmptof < 100 )
+ {
+ tmpindex = 1;
+ }
+ else
+ {
+ tmpindex = 2;
+ }
- for(int k=0; k < 3; k ++)
+ double tof = 0;
+ for(int k=0; k < PF::SAMPLERANGE; k++ )
{
- // cout << __FILE__ << __LINE__ << "amp[="<< k <<"] = " << fAmpA[k] << endl;
- if( TMath::Abs(fAmpA[k] - ( maxamp - ped) ) < diff)
- {
- diff = TMath::Abs(fAmpA[k] - ( maxamp - ped));
- tmpindex = k;
- }
+ tof += fPFTofVectors[0][pfindex][k]*fReversed[ dt +k + tmpindex -1 ];
+ }
+ for( int i=0; i < PF::SAMPLERANGE; i++ )
+ {
+ {
+
+ fAmp += fPFAmpVectors[0][pfindex][i]*fReversed[ dt +i +tmpindex -1 ];
+ }
}
-
- double tof = 0;
- for(int k=0; k < SAMPLERANGE; k++ )
+ if( TMath::Abs( (maxf - fAmp )/maxf ) > 0.1 )
{
- tof += fPFTofVectors[0][pfindex][k]*fReversed[ dt +k + tmpindex -1 ];
+ fAmp = maxf;
}
- tof = tof / fAmpA[tmpindex];
-
- // return AliCaloFitResults( maxamp, ped , -1, fAmp, -1, -1, -1 );
- return AliCaloFitResults( maxamp, ped , -1, fAmpA[tmpindex], tof, -2, -3 );
-
+ tof = timebinOffset - 0.01*tof/fAmp - fL1Phase/TIMEBINWITH; // clock
+ Float_t chi2 = CalculateChi2(fAmp, tof-timebinOffset+maxrev, first, last);
+ Int_t ndf = last - first - 1; // nsamples - 2
+ return AliCaloFitResults( maxamp, ped , Ret::kFitPar, fAmp, tof,
+ timebinOffset, chi2, ndf,
+ Ret::kDummy, AliCaloFitSubarray(index, maxrev, first, last) );
}
-
else
{
- return AliCaloFitResults( maxamp, ped , -5, maxf, -6, -7, -8 );
+ Float_t chi2 = CalculateChi2(maxf, maxrev, first, last);
+ Int_t ndf = last - first - 1; // nsamples - 2
+ return AliCaloFitResults( maxamp, ped , Ret::kCrude, maxf, timebinOffset,
+ timebinOffset, chi2, ndf, Ret::kDummy, AliCaloFitSubarray(index, maxrev, first, last) );
}
- }
+ } // ampcut
}
-
- // cout << __FILE__ << __LINE__ << "WARNING, returning amp = -1 " << endl;
-
- return AliCaloFitResults(-1, -1);
+ return AliCaloFitResults(kInvalid, kInvalid);
}
-void
-AliCaloRawAnalyzerPeakFinder::LoadVectors()
+void
+AliCaloRawAnalyzerPeakFinder::CopyVectors( const AliCaloPeakFinderVectors *const pfv )
{
- //Read in the Peak finder vecors from file
- for(int i = 0; i < MAXSTART ; i++)
+ // As name implies
+ if ( pfv != 0)
{
- for( int j=0; j < SAMPLERANGE; j++)
+ for(int i = 0; i < PF::MAXSTART ; i++)
{
- char filename[256];
- int n = j+fNsampleCut;
- double start = (double)i+0.5;
- sprintf(filename, "%s/EMCAL/vectors-emcal/start%.1fN%dtau0.235fs10dt1.0.txt", getenv("ALICE_ROOT"), start, n);
- FILE *fp = fopen(filename, "r");
-
- if(fp == 0)
+ for( int j=0; j < PF::SAMPLERANGE; j++)
{
- AliFatal( Form( "could not open file: %s", filename ) );
+ pfv->GetVector( i, j, fPFAmpVectors[i][j] , fPFTofVectors[i][j],
+ fPFAmpVectorsCoarse[i][j] , fPFTofVectorsCoarse[i][j] );
+
+ fPeakFinderVectors->SetVector( i, j, fPFAmpVectors[i][j], fPFTofVectors[i][j],
+ fPFAmpVectorsCoarse[i][j], fPFTofVectorsCoarse[i][j] );
}
- else
+ }
+ }
+ else
+ {
+ AliFatal( "pfv = ZERO !!!!!!!");
+ }
+}
+
+
+void
+AliCaloRawAnalyzerPeakFinder::LoadVectorsOCDB()
+{
+ //Loading of Peak-Finder vectors from the
+ //Offline Condition Database (OCDB)
+ AliCDBEntry* entry = AliCDBManager::Instance()->Get("EMCAL/Calib/PeakFinder/");
+
+ if( entry != 0 )
+ {
+ //cout << __FILE__ << ":" << __LINE__ << ": Printing metadata !! " << endl;
+ //entry->PrintMetaData();
+ AliCaloPeakFinderVectors *pfv = (AliCaloPeakFinderVectors *)entry->GetObject();
+ if( pfv == 0 )
+ {
+ cout << __FILE__ << ":" << __LINE__ << "_ ERRROR " << endl;
+ }
+ CopyVectors( pfv );
+
+ if( pfv != 0 )
+ {
+ fIsInitialized = true;
+ }
+ }
+}
+
+
+void
+AliCaloRawAnalyzerPeakFinder::WriteRootFile() const
+{ // Utility function to write Peak-Finder vectors to an root file
+ // The output is used to create an OCDB entry.
+ fPeakFinderVectors->PrintVectors();
+ TFile *f = new TFile("peakfindervectors2.root", "recreate" );
+ fPeakFinderVectors->Write();
+ f->Close();
+ delete f;
+}
+
+
+void
+AliCaloRawAnalyzerPeakFinder::PrintVectors()
+{ // Utility function to write Peak-Finder vectors
+ for(int i=0; i < 20; i++)
+ {
+ for( int j = 0; j < PF::MAXSTART; j ++ )
+ {
+ for( int k=0; k < PF::SAMPLERANGE; k++ )
{
- for(int m = 0; m < n ; m++ )
- {
- fscanf(fp, "%lf\t", &fPFAmpVectors[i][j][m] );
- }
- fscanf(fp, "\n");
-
- for(int m = 0; m < n ; m++ )
- {
- fPFTofVectors[i][j][m] = 1;
- // fscanf(fp, "%lf\t", &fPFTofVectors[i][j][m] );
- }
-
- fclose (fp);
+ cout << fPFAmpVectors[j][k][i] << "\t" ;
}
}
+ cout << endl;
}
+ cout << __FILE__ << ":" << __LINE__ << ":.... DONE !!" << endl;
}