#include "AliEMCALRawUtils.h"
#include "AliRun.h"
#include "AliRunLoader.h"
#include "AliAltroBuffer.h"
#include "AliRawReader.h"
#include "AliCaloRawStreamV3.h"
#include "AliDAQ.h"
#include "AliEMCALRecParam.h"
#include "AliEMCALLoader.h"
#include "AliEMCALGeometry.h"
#include "AliEMCALDigit.h"
#include "AliEMCALRawDigit.h"
#include "AliEMCAL.h"
#include "AliCaloCalibPedestal.h"
#include "AliCaloBunchInfo.h"
#include "AliCaloFitResults.h"
#include "AliEMCALTriggerRawDigitMaker.h"
#include "AliEMCALTriggerSTURawStream.h"
#include "AliEMCALTriggerData.h"
#include "AliCaloConstants.h"
#include "AliCaloRawAnalyzer.h"
#include "AliCaloRawAnalyzerFactory.h"
#include "AliEMCALRawResponse.h"
using namespace CALO;
using namespace EMCAL;
using std::vector;
ClassImp(AliEMCALRawUtils)
AliEMCALRawUtils::AliEMCALRawUtils( Algo::fitAlgorithm fitAlgo) : fNoiseThreshold(3),
fNPedSamples(4),
fGeom(0),
fOption(""),
fRemoveBadChannels(kFALSE),
fFittingAlgorithm(fitAlgo),
fTimeMin(-1.),
fTimeMax(1.),
fUseFALTRO(kTRUE),
fRawAnalyzer(0),
fTriggerRawDigitMaker(0x0)
{
SetFittingAlgorithm(fitAlgo);
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);
}
AliRunLoader *rl = AliRunLoader::Instance();
if (rl && rl->GetAliRun())
{
AliEMCAL * emcal = dynamic_cast<AliEMCAL*>(rl->GetAliRun()->GetDetector("EMCAL"));
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());
}
if(!fGeom) AliFatal(Form("Could not get geometry!"));
fTriggerRawDigitMaker = new AliEMCALTriggerRawDigitMaker();
}
AliEMCALRawUtils::~AliEMCALRawUtils()
{
delete fRawAnalyzer;
delete fTriggerRawDigitMaker;
}
void AliEMCALRawUtils::Digits2Raw()
{
AliRunLoader *rl = AliRunLoader::Instance();
AliEMCALLoader *loader = dynamic_cast<AliEMCALLoader*>(rl->GetDetectorLoader("EMCAL"));
loader->LoadDigits("EMCAL");
loader->GetEvent();
TClonesArray* digits = loader->Digits() ;
if (!digits) {
Warning("Digits2Raw", "no digits found !");
return;
}
static const Int_t nDDL = 20*2;
AliAltroBuffer* buffers[nDDL];
for (Int_t i=0; i < nDDL; i++)
buffers[i] = 0;
TArrayI adcValuesLow( TIMEBINS );
TArrayI adcValuesHigh( TIMEBINS );
for (Int_t iDigit = 0; iDigit < digits->GetEntries(); iDigit++) {
AliEMCALDigit* digit = dynamic_cast<AliEMCALDigit *>(digits->At(iDigit)) ;
if(!digit) {
AliFatal("NULL Digit");
} else {
if (digit->GetAmplitude() < AliEMCALRawResponse::GetRawFormatThreshold() ) {
continue;
}
Int_t nSM = 0;
Int_t nIphi = 0;
Int_t nIeta = 0;
Int_t iphi = 0;
Int_t ieta = 0;
Int_t nModule = 0;
fGeom->GetCellIndex(digit->GetId(), nSM, nModule, nIphi, nIeta);
fGeom->GetCellPhiEtaIndexInSModule(nSM, nModule, nIphi, nIeta,iphi, ieta) ;
Int_t iRCU = -111;
if (0<=iphi&&iphi<8) iRCU=0;
else if (8<=iphi&&iphi<16 && 0<=ieta&&ieta<24) iRCU=0;
else if(8<=iphi&&iphi<16 && 24<=ieta&&ieta<48) iRCU=1;
else if(16<=iphi&&iphi<24) iRCU=1;
if (nSM%2==1) iRCU = 1 - iRCU;
if (iRCU<0)
Fatal("Digits2Raw()","Non-existent RCU number: %d", iRCU);
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) {
TString fileName = AliDAQ::DdlFileName("EMCAL",iDDL);
Int_t iRCUside=iRCU+(nSM%2)*2;
buffers[iDDL] = new AliAltroBuffer(fileName.Data(),fMapping[iRCUside]);
buffers[iDDL]->WriteDataHeader(kTRUE, kFALSE);
}
if (digit->GetTimeR() > TIMEBINMAX ) {
AliInfo("Signal is out of time range.\n");
buffers[iDDL]->FillBuffer((Int_t)digit->GetAmplitude());
buffers[iDDL]->FillBuffer( TIMEBINS );
buffers[iDDL]->FillBuffer(3);
buffers[iDDL]->WriteTrailer(3, ieta, iphi, nSM);
} else {
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(), AliEMCALRawResponse::GetRawFormatThreshold() );
}
}
}
}
for (Int_t i=0; i < nDDL; i++) {
if (buffers[i]) {
buffers[i]->Flush();
buffers[i]->WriteDataHeader(kFALSE, kFALSE);
delete buffers[i];
}
}
loader->UnloadDigits();
}
void AliEMCALRawUtils::AddDigit(TClonesArray *digitsArr, Int_t id, Int_t lowGain, Float_t amp, Float_t time, Float_t chi2, Int_t ndf)
{
AliEMCALDigit *digit = 0, *tmpdigit = 0;
TIter nextdigit(digitsArr);
while (digit == 0 && (tmpdigit = (AliEMCALDigit*) nextdigit()))
{
if (tmpdigit->GetId() == id) digit = tmpdigit;
}
if (!digit) {
Int_t type = AliEMCALDigit::kHG;
if (lowGain)
{
amp *= HGLGFACTOR;
type = AliEMCALDigit::kLGnoHG;
}
Int_t idigit = digitsArr->GetEntries();
new((*digitsArr)[idigit]) AliEMCALDigit( -1, -1, id, amp, time, type, idigit, chi2, ndf);
AliDebug(2,Form("Add digit Id %d for the first time, type %d", id, type));
}
else
{
if (lowGain)
{
if (digit->GetAmplitude() > OVERFLOWCUT )
{
digit->SetAmplitude( HGLGFACTOR * amp);
digit->SetTime(time);
digit->SetType(AliEMCALDigit::kLG);
AliDebug(2,Form("Add LG digit ID %d for the second time, type %d", digit->GetId(), digit->GetType()));
}
}
else {
if (amp < OVERFLOWCUT )
{
digit->SetAmplitude(amp);
digit->SetTime(time);
digit->SetType(AliEMCALDigit::kHG);
AliDebug(2,Form("Add HG digit ID %d for the second time, type %d", digit->GetId(), digit->GetType()));
}
else
{
digit->SetType(AliEMCALDigit::kLG);
AliDebug(2,Form("Change LG digit to HG, ID %d, type %d", digit->GetId(), digit->GetType()));
}
}
}
}
void AliEMCALRawUtils::Raw2Digits(AliRawReader* reader,TClonesArray *digitsArr, const AliCaloCalibPedestal* pedbadmap, TClonesArray *digitsTRG, AliEMCALTriggerData* trgData)
{
if ( digitsArr) digitsArr->Clear("C");
if (!digitsArr) { Error("Raw2Digits", "no digits found !");return;}
if (!reader) {Error("Raw2Digits", "no raw reader found !");return;}
AliEMCALTriggerSTURawStream inSTU(reader);
AliCaloRawStreamV3 in(reader,"EMCAL",fMapping);
reader->Select("EMCAL",0,39);
fTriggerRawDigitMaker->Reset();
fTriggerRawDigitMaker->SetIO(reader, in, inSTU, digitsTRG, trgData);
fRawAnalyzer->SetIsZeroSuppressed(true);
Int_t lowGain = 0;
Int_t caloFlag = 0;
Float_t bcTimePhaseCorr = 0;
Int_t bcMod4 = (reader->GetBCID() % 4);
Int_t runNumber = reader->GetRunNumber();
if ((runNumber >130850 ) && (bcMod4==0 || bcMod4==1))
bcTimePhaseCorr = -1e-7;
while (in.NextDDL())
{
while (in.NextChannel())
{
caloFlag = in.GetCaloFlag();
if (caloFlag > 2) continue;
if(caloFlag < 2 && fRemoveBadChannels && pedbadmap->IsBadChannel(in.GetModule(),in.GetColumn(),in.GetRow()))
{
continue;
}
vector<AliCaloBunchInfo> bunchlist;
while (in.NextBunch())
{
bunchlist.push_back( AliCaloBunchInfo(in.GetStartTimeBin(), in.GetBunchLength(), in.GetSignals() ) );
}
if (bunchlist.size() == 0) continue;
if ( caloFlag < 2 )
{
Int_t id = fGeom->GetAbsCellIdFromCellIndexes(in.GetModule(), in.GetRow(), in.GetColumn()) ;
lowGain = in.IsLowGain();
fRawAnalyzer->SetL1Phase( in.GetL1Phase() );
AliCaloFitResults res = fRawAnalyzer->Evaluate( bunchlist, in.GetAltroCFG1(), in.GetAltroCFG2());
if(res.GetAmp() >= fNoiseThreshold )
{
AddDigit(digitsArr, id, lowGain, res.GetAmp(), res.GetTime()+bcTimePhaseCorr, res.GetChi2(), res.GetNdf() );
}
}
else if(fUseFALTRO)
{
fTriggerRawDigitMaker->Add( bunchlist );
}
}
}
fTriggerRawDigitMaker->PostProcess();
TrimDigits(digitsArr);
}
void AliEMCALRawUtils::TrimDigits(TClonesArray *digitsArr)
{
AliEMCALDigit *digit = 0;
Int_t n = 0;
Int_t nDigits = digitsArr->GetEntriesFast();
TIter nextdigit(digitsArr);
while ((digit = (AliEMCALDigit*) nextdigit())) {
if (digit->GetType() == AliEMCALDigit::kLGnoHG) {
AliDebug(1,Form("Remove digit with id %d, LGnoHG",digit->GetId()));
digitsArr->Remove(digit);
}
else if(fTimeMin > digit->GetTime() || fTimeMax < digit->GetTime()) {
digitsArr->Remove(digit);
AliDebug(1,Form("Remove digit with id %d, Bad Time %e",digit->GetId(), digit->GetTime()));
}
else if (0 > digit->GetChi2()) {
digitsArr->Remove(digit);
AliDebug(1,Form("Remove digit with id %d, Bad Chi2 %e",digit->GetId(), digit->GetChi2()));
}
else {
digit->SetIndexInList(n);
n++;
}
}
digitsArr->Compress();
AliDebug(1,Form("N Digits before trimming : %d; after array compression %d",nDigits,digitsArr->GetEntriesFast()));
}
void AliEMCALRawUtils::SetFittingAlgorithm(Int_t fitAlgo)
{
delete fRawAnalyzer;
fRawAnalyzer = AliCaloRawAnalyzerFactory::CreateAnalyzer( fitAlgo );
fRawAnalyzer->SetNsampleCut(5);
fRawAnalyzer->SetOverflowCut ( OVERFLOWCUT );
fRawAnalyzer->SetAmpCut(fNoiseThreshold);
fRawAnalyzer->SetFitArrayCut(fNoiseThreshold);
}