// Need to add concurrent high and low-gain info in the future
// No pedestal is added to the raw signal.
//*-- Author: Marco van Leeuwen (LBL)
-
+//*-- Major refactoring by Per Thomas Hille
#include "AliEMCALRawUtils.h"
-#include "TF1.h"
-#include "TGraph.h"
-#include <TRandom.h>
#include "AliRun.h"
#include "AliRunLoader.h"
#include "AliAltroBuffer.h"
#include "AliCaloConstants.h"
#include "AliCaloRawAnalyzer.h"
#include "AliCaloRawAnalyzerFactory.h"
+#include "AliEMCALRawResponse.h"
using namespace CALO;
using namespace EMCAL;
-Double_t AliEMCALRawUtils::fgTimeTrigger = 600E-9 ; // the time of the trigger as approximately seen in the data
-Int_t AliEMCALRawUtils::fgThreshold = 1;
-Int_t AliEMCALRawUtils::fgPedestalValue = 0; // pedestal value for digits2raw, default generate ZS data
-Double_t AliEMCALRawUtils::fgFEENoise = 3.; // 3 ADC channels of noise (sampled)
-
ClassImp(AliEMCALRawUtils)
fGeom(0),
fOption(""),
fRemoveBadChannels(kFALSE),
- fFittingAlgorithm(0),
+ fFittingAlgorithm(fitAlgo),
fTimeMin(-1.),
fTimeMax(1.),
fUseFALTRO(kTRUE),
fTriggerRawDigitMaker(0x0)
{
SetFittingAlgorithm(fitAlgo);
- // SetFittingAlgorithm( Algo::kLMSOffline);
-
- //Get Mapping RCU files from the AliEMCALRecParam
const TObjArray* maps = AliEMCALRecParam::GetMappings();
if(!maps) AliFatal("Cannot retrieve ALTRO mappings!!");
-
- for(Int_t i = 0; i < 4; i++) {
- fMapping[i] = (AliAltroMapping*)maps->At(i);
- }
-
- //To make sure we match with the geometry in a simulation file,
- //let's try to get it first. If not, take the default geometry
+ for(Int_t i = 0; i < 4; i++)
+ {
+ fMapping[i] = (AliAltroMapping*)maps->At(i);
+ }
+
AliRunLoader *rl = AliRunLoader::Instance();
- if (rl && rl->GetAliRun()) {
+ if (rl && rl->GetAliRun())
+ {
AliEMCAL * emcal = dynamic_cast<AliEMCAL*>(rl->GetAliRun()->GetDetector("EMCAL"));
- if(emcal)fGeom = emcal->GetGeometry();
- else {
+ if(emcal)
+ {
+ fGeom = emcal->GetGeometry();
+ }
+ else
+ {
+ AliDebug(1, Form("Using default geometry in raw reco"));
+ fGeom = AliEMCALGeometry::GetInstance(AliEMCALGeometry::GetDefaultGeometryName());
+ }
+ }
+ else
+ {
AliDebug(1, Form("Using default geometry in raw reco"));
fGeom = AliEMCALGeometry::GetInstance(AliEMCALGeometry::GetDefaultGeometryName());
}
-
- } else {
- AliDebug(1, Form("Using default geometry in raw reco"));
- fGeom = AliEMCALGeometry::GetInstance(AliEMCALGeometry::GetDefaultGeometryName());
- }
-
- if(!fGeom) AliFatal(Form("Could not get geometry!"));
-
- fTriggerRawDigitMaker = new AliEMCALTriggerRawDigitMaker();
-
-}
-
-
-//____________________________________________________________________________
-AliEMCALRawUtils::AliEMCALRawUtils(AliEMCALGeometry *pGeometry, Algo::fitAlgorithm fitAlgo) : //fHighLowGainFactor(16.),
- //fOrder(2),
- // fTau(2.35),
- fNoiseThreshold(3),
- fNPedSamples(4),
- fGeom(pGeometry),
- fOption(""),
- fRemoveBadChannels(kFALSE),fFittingAlgorithm(0),
- fTimeMin(-1.),fTimeMax(1.),
- fUseFALTRO(kTRUE),fRawAnalyzer(0),
- fTriggerRawDigitMaker(0x0)
-{
-
-
- // Initialize with the given geometry - constructor required by HLT
- // HLT does not use/support AliRunLoader(s) instances
- // This is a minimum intervention solution
- // Comment by MPloskon@lbl.gov
- SetFittingAlgorithm(fitAlgo);
- // SetFittingAlgorithm( Algo::kLMSOffline);
- //Get Mapping RCU files from the AliEMCALRecParam
- const TObjArray* maps = AliEMCALRecParam::GetMappings();
- if(!maps) AliFatal("Cannot retrieve ALTRO mappings!!");
- for(Int_t i = 0; i < 4; i++)
- {
- fMapping[i] = (AliAltroMapping*)maps->At(i);
- }
-
if(!fGeom) AliFatal(Form("Could not get geometry!"));
- fTriggerRawDigitMaker = new AliEMCALTriggerRawDigitMaker();
+ fTriggerRawDigitMaker = new AliEMCALTriggerRawDigitMaker();
}
-//____________________________________________________________________________
AliEMCALRawUtils::~AliEMCALRawUtils()
{
//dtor
}
-//____________________________________________________________________________
void AliEMCALRawUtils::Digits2Raw()
{
// convert digits of the current event to raw data
AliRunLoader *rl = AliRunLoader::Instance();
AliEMCALLoader *loader = dynamic_cast<AliEMCALLoader*>(rl->GetDetectorLoader("EMCAL"));
-
- // get the digits
loader->LoadDigits("EMCAL");
loader->GetEvent();
TClonesArray* digits = loader->Digits() ;
}
else
{
- if (digit->GetAmplitude() < fgThreshold)
+ if (digit->GetAmplitude() < AliEMCALRawResponse::GetRawFormatThreshold() )
{
continue;
}
//Check which is the RCU, 0 or 1, of the cell.
Int_t iRCU = -111;
- //RCU0
- if (0<=iphi&&iphi<8) iRCU=0; // first cable row
- else if (8<=iphi&&iphi<16 && 0<=ieta&&ieta<24) iRCU=0; // first half;
- //second cable row
- //RCU1
- else if(8<=iphi&&iphi<16 && 24<=ieta&&ieta<48) iRCU=1; // second half;
- //second cable row
- else if(16<=iphi&&iphi<24) iRCU=1; // third cable row
-
- if (nSM%2==1) iRCU = 1 - iRCU; // swap for odd=C side, to allow us to cable both sides the same
-
- if (iRCU<0)
- Fatal("Digits2Raw()","Non-existent RCU number: %d", iRCU);
+ if (0<=iphi&&iphi<8) iRCU=0; // first cable row
+ else if (8<=iphi&&iphi<16 && 0<=ieta&&ieta<24) iRCU=0; // first half;
+ else if(8<=iphi&&iphi<16 && 24<=ieta&&ieta<48) iRCU=1; // second half;
+ //second cable row
+ else if(16<=iphi&&iphi<24) iRCU=1; // third cable row
+
+ if (nSM%2==1) iRCU = 1 - iRCU; // swap for odd=C side, to allow us to cable both sides the same
+
+ if (iRCU<0)
+ Fatal("Digits2Raw()","Non-existent RCU number: %d", iRCU);
- //Which DDL?
- Int_t iDDL = NRCUSPERMODULE*nSM + iRCU;
- if (iDDL < 0 || iDDL >= nDDL){
- Fatal("Digits2Raw()","Non-existent DDL board number: %d", iDDL);
- }
- else{
- if (buffers[iDDL] == 0) {
- // open new file and write dummy header
- TString fileName = AliDAQ::DdlFileName("EMCAL",iDDL);
- //Select mapping file RCU0A, RCU0C, RCU1A, RCU1C
+ //Which DDL?
+ Int_t iDDL = NRCUSPERMODULE*nSM + iRCU;
+ if (iDDL < 0 || iDDL >= nDDL){
+ Fatal("Digits2Raw()","Non-existent DDL board number: %d", iDDL);
+ }
+ else{
+ if (buffers[iDDL] == 0)
+ {
+ // open new file and write dummy header
+ TString fileName = AliDAQ::DdlFileName("EMCAL",iDDL);
+ //Select mapping file RCU0A, RCU0C, RCU1A, RCU1C
Int_t iRCUside=iRCU+(nSM%2)*2;
//iRCU=0 and even (0) SM -> RCU0A.data 0
//iRCU=1 and even (0) SM -> RCU1A.data 1
//iRCU=0 and odd (1) SM -> RCU0C.data 2
//iRCU=1 and odd (1) SM -> RCU1C.data 3
- //cout<<" nSM "<<nSM<<"; iRCU "<<iRCU<<"; iRCUside "<<iRCUside<<endl;
- buffers[iDDL] = new AliAltroBuffer(fileName.Data(),fMapping[iRCUside]);
+ buffers[iDDL] = new AliAltroBuffer(fileName.Data(),fMapping[iRCUside]);
buffers[iDDL]->WriteDataHeader(kTRUE, kFALSE); //Dummy;
}
buffers[iDDL]->WriteTrailer(3, ieta, iphi, nSM); // trailer
// calculate the time response function
} else {
- Bool_t lowgain = RawSampledResponse(digit->GetTimeR(), digit->GetAmplitude(), adcValuesHigh.GetArray(), adcValuesLow.GetArray()) ;
- if (lowgain)
- buffers[iDDL]->WriteChannel(ieta, iphi, 0, TIMEBINS, adcValuesLow.GetArray(), fgThreshold);
+ Bool_t lowgain = AliEMCALRawResponse::RawSampledResponse(digit->GetTimeR(), digit->GetAmplitude(),
+ adcValuesHigh.GetArray(), adcValuesLow.GetArray()) ;
+
+ if (lowgain)
+ buffers[iDDL]->WriteChannel(ieta, iphi, 0, TIMEBINS, adcValuesLow.GetArray(), AliEMCALRawResponse::GetRawFormatThreshold() );
else
- buffers[iDDL]->WriteChannel(ieta,iphi, 1, TIMEBINS, adcValuesHigh.GetArray(), fgThreshold);
+ buffers[iDDL]->WriteChannel(ieta,iphi, 1, TIMEBINS, adcValuesHigh.GetArray(), AliEMCALRawResponse::GetRawFormatThreshold() );
}
}// iDDL under the limits
}//digit exists
}
-//____________________________________________________________________________
+
void AliEMCALRawUtils::AddDigit(TClonesArray *digitsArr, Int_t id, Int_t lowGain, Float_t amp, Float_t time, Float_t chi2, Int_t ndf)
{
+ // comment
AliEMCALDigit *digit = 0, *tmpdigit = 0;
TIter nextdigit(digitsArr);
}
-//____________________________________________________________________________
void AliEMCALRawUtils::Raw2Digits(AliRawReader* reader,TClonesArray *digitsArr, const AliCaloCalibPedestal* pedbadmap, TClonesArray *digitsTRG, AliEMCALTriggerData* trgData)
{
-
+ //conversion of raw data to digits
if(digitsArr) digitsArr->Clear("C");
if (!digitsArr) { Error("Raw2Digits", "no digits found !");return;}
if (!reader) {Error("Raw2Digits", "no raw reader found !");return;}
while (in.NextDDL())
{
- // fprintf(fp," TP1\n");
while (in.NextChannel())
{
- // fprintf(fp," TP2\n");
caloFlag = in.GetCaloFlag();
if (caloFlag > 2) continue; // Work with ALTRO and FALTRO
if(caloFlag < 2 && fRemoveBadChannels && pedbadmap->IsBadChannel(in.GetModule(),in.GetColumn(),in.GetRow()))
}
-Double_t
-AliEMCALRawUtils::RawResponseFunction(Double_t *x, Double_t *par)
-{
- Double_t signal = 0.;
- Double_t tau = par[2];
- Double_t n = par[3];
- Double_t ped = par[4];
- Double_t xx = ( x[0] - par[1] + tau ) / tau ;
- if (xx <= 0)
- signal = ped ;
- else
- {
- signal = ped + par[0] * TMath::Power(xx , n) * TMath::Exp(n * (1 - xx )) ;
- }
- return signal ;
-}
-
-
-Bool_t AliEMCALRawUtils::RawSampledResponse(const Double_t dtime, const Double_t damp, Int_t * adcH,
- Int_t * adcL, const Int_t keyErr) const
-{
- Bool_t lowGain = kFALSE ;
- TF1 signalF("signal", RawResponseFunction, 0, TIMEBINS, 5);
- signalF.SetParameter(0, damp) ;
- signalF.SetParameter(1, (dtime + fgTimeTrigger)/ TIMEBINWITH) ;
- signalF.SetParameter(2, TAU) ;
- signalF.SetParameter(3, ORDER);
- signalF.SetParameter(4, fgPedestalValue);
-
- Double_t signal=0.0, noise=0.0;
- for (Int_t iTime = 0; iTime < TIMEBINS; iTime++) {
- signal = signalF.Eval(iTime) ;
-
- if(keyErr>0) {
- noise = gRandom->Gaus(0.,fgFEENoise);
- signal += noise;
- }
-
- adcH[iTime] = static_cast<Int_t>(signal + 0.5) ;
- if ( adcH[iTime] > MAXBINVALUE ){ // larger than 10 bits
- adcH[iTime] = MAXBINVALUE ;
- lowGain = kTRUE ;
- }
- signal /= HGLGFACTOR;
- adcL[iTime] = static_cast<Int_t>(signal + 0.5) ;
- if ( adcL[iTime] > MAXBINVALUE ) // larger than 10 bits
- adcL[iTime] = MAXBINVALUE ;
- }
- return lowGain ;
-}
-
-
-//__________________________________________________________________
-Double_t AliEMCALRawUtils::RawResponseFunctionLog(Double_t *x, Double_t *par)
-{
- Double_t signal = 0. ;
- Double_t tau = par[2];
- Double_t n = par[3];
- Double_t xx = ( x[0] - par[1] + tau ) / tau ;
-
- if (xx < 0)
- {
- signal = par[0] - n*TMath::Log(TMath::Abs(xx)) + n * (1 - xx ) ;
- }
- else
- {
- signal = par[0] + n*TMath::Log(xx) + n * (1 - xx ) ;
- }
- return signal ;
-}
-
-
-
-//__________________________________________________________________
void AliEMCALRawUtils::SetFittingAlgorithm(Int_t fitAlgo)
{
- // fRawAnalyzer = AliCaloRawAnalyzerFactory::CreateAnalyzer( Algo::kStandard );
fRawAnalyzer = AliCaloRawAnalyzerFactory::CreateAnalyzer( fitAlgo );
-
- //fRawAnalyzer = AliCaloRawAnalyzerFactory::CreateAnalyzer( kStandard );
-
fRawAnalyzer->SetNsampleCut(5); // requirement for fits to be done, for the new methods
fRawAnalyzer->SetOverflowCut ( OVERFLOWCUT );
fRawAnalyzer->SetAmpCut(fNoiseThreshold);
fRawAnalyzer->SetFitArrayCut(fNoiseThreshold);
- // return;
}