#include "AliFMDDebug.h"
#include <AliLoader.h>
#include <AliAltroBufferV3.h>
#include "AliFMD.h"
#include "AliFMDParameters.h"
#include "AliFMDDigit.h"
#include "AliFMDRawWriter.h"
#include "AliFMDAltroMapping.h"
#include <TArrayI.h> // ROOT_TArrayI
#include <TArrayF.h> // ROOT_TArrayI
#include <TArrayC.h> // ROOT_TArrayI
#include <TClonesArray.h> // ROOT_TClonesArray
#include <TTree.h>
#include "AliDAQ.h"
ClassImp(AliFMDRawWriter)
#if 0
;
#endif
AliFMDRawWriter::AliFMDRawWriter(AliFMD* fmd)
: TTask("FMDRawWriter", "Writer of Raw ADC values from the FMD"),
fFMD(fmd),
fSampleRate(0),
fChannelsPerAltro(0),
fThreshold(0)
{
AliFMDDebug(5, ("Created AliFMDRawWriter object"));
}
void
AliFMDRawWriter::Exec(Option_t*)
{
AliLoader* loader = fFMD->GetLoader();
loader->LoadDigits("READ");
TTree* digitTree = loader->TreeD();
if (!digitTree) {
AliError("no digit tree");
return;
}
fFMD->SetTreeAddress();
TClonesArray* digits = fFMD->Digits();
Int_t nEvents = Int_t(digitTree->GetEntries());
AliFMDDebug(5, ("Got a total of %5d events from tree", nEvents));
for (Int_t event = 0; event < nEvents; event++) {
fFMD->ResetDigits();
digitTree->GetEvent(event);
WriteDigits(digits);
}
loader->UnloadDigits();
}
#if 1
Long_t
AliFMDRawWriter::WriteDigits(TClonesArray* digits)
{
Int_t nDigits = digits->GetEntries();
if (nDigits < 1) return 0;
AliFMDDebug(5, ("Got a total of %5d digits from tree", nDigits));
AliFMDParameters* pars = AliFMDParameters::Instance();
UShort_t threshold = 0;
UShort_t factor = 0;
UInt_t prevddl = 0xFFFF;
UInt_t prevaddr = 0xFFF;
UShort_t nWords = 0;
UShort_t preSamples = 0;
UShort_t sampleRate = 0;
TArrayI data(pars->GetChannelsPerAltro() * 8);
TArrayF peds(pars->GetChannelsPerAltro() * 8);
TArrayF noise(pars->GetChannelsPerAltro() * 8);
AliAltroBufferV3* altro = 0;
Int_t totalWords = 0;
Int_t nCounts = 0;
Long_t nBits = 0;
UShort_t oldDet = 1000;
for (Int_t i = 0; i < nDigits; i++) {
AliFMDDigit* digit = static_cast<AliFMDDigit*>(digits->At(i));
UShort_t det = digit->Detector();
Char_t ring = digit->Ring();
UShort_t sector = digit->Sector();
UShort_t strip = digit->Strip();
UShort_t ddl, addr, time;
AliFMDDebug(15, ("Processing digit # %5d FMD%d%c[%2d,%3d]",
i, det, ring, sector, strip));
threshold = pars->GetZeroSuppression(det, ring, sector, strip);
sampleRate = pars->GetSampleRate(det, ring, sector, strip);
preSamples = pars->GetPreSamples(det, ring, sector, strip);
factor = UShort_t(pars->GetPedestalFactor());
if (det != oldDet) {
AliFMDDebug(5, ("Got new detector: %d (was %d)", det, oldDet));
oldDet = det;
}
AliFMDDebug(15, ("Sample rate is %d", sampleRate));
for (UShort_t j = 0; j < sampleRate; j++) {
if (!pars->Detector2Hardware(det,ring,sector,strip,j,ddl,addr,time)){
AliError(Form("Failed to get hardware address for FMD%d%c[%2d,%3d]-%d",
det, ring, sector, strip, j));
continue;
}
AliFMDDebug(10, ("FMD%d%c[%2d,%3d]-%d-> 0x%x/0x%x/%04d",
det, ring, sector, strip, j, ddl, addr, time));
if (addr != prevaddr) {
AliFMDDebug(5, ("Now hardware address 0x%x from FMD%d%c[%2d,%3d]-%d"
"(b: 0x%02x, a: 0x%01x, c: 0x%02x, t: %04d), "
"flushing old channel at 0x%x with %d words",
addr, det, ring, sector, strip, j,
(addr >> 7), (addr >> 4) & 0x7, addr & 0xf,
time, prevaddr, nWords));
totalWords += nWords;
ZeroSuppress(data.fArray, nWords, peds.fArray, noise.fArray, threshold);
if (altro)
altro->WriteChannel(prevaddr,nWords,data.fArray,threshold);
data.Reset(-1);
peds.Reset(0);
noise.Reset(0);
nWords = 0;
prevaddr = addr;
}
if (ddl != prevddl) {
AliFMDDebug(10, ("FMD: New DDL, was %d, now %d", prevddl, ddl));
if (altro) {
AliFMDDebug(15, ("Closing output"));
WriteRCUTrailer(altro, prevddl, threshold > 0, factor, sampleRate);
delete altro;
altro = 0;
}
prevddl = ddl;
TString filename(AliDAQ::DdlFileName(fFMD ? fFMD->GetName() : "FMD", ddl));
AliFMDDebug(5, ("New altro buffer with DDL file %s", filename.Data()));
altro = new AliAltroBufferV3(filename.Data());
altro->SetMapping(pars->GetAltroMap());
altro->WriteDataHeader(kTRUE, kFALSE);
}
peds[time] = pars->GetPedestal(det, ring, sector, strip);
noise[time] = pars->GetPedestalWidth(det, ring, sector, strip);
AliFMDDebug(15, ("Storing FMD%d%c[%02d,%03d]-%d in timebin %d (%d)",
det, ring, sector, strip, j, time, preSamples));
data[time] = digit->Count(j);
nWords++;
nCounts++;
if (time == preSamples) {
AliFMDDebug(15, ("Filling in %4d for %d presamples",
data[time], preSamples));
for (int k = 0; k < preSamples; k++) {
peds[k] = peds[time];
noise[k] = noise[time];
data[k] = data[time];
}
nWords += preSamples;
}
}
}
if (altro) {
ZeroSuppress(data.fArray, nWords, peds.fArray, noise.fArray, threshold);
if (nWords > 0)
altro->WriteChannel(prevaddr,nWords,data.fArray,threshold);
WriteRCUTrailer(altro, prevddl, threshold > 0, factor, sampleRate);
delete altro;
}
AliFMDDebug(5, ("Wrote a total of %d words in %ld bytes for %d counts",
nWords, nBits / 8, nCounts));
return nBits;
}
void
AliFMDRawWriter::WriteRCUTrailer(AliAltroBufferV3* altro,
UInt_t ddl,
Bool_t zs,
UShort_t factor,
UShort_t rate) const
{
altro->Flush();
altro->WriteDataHeader(kFALSE, kFALSE);
altro->SetZeroSupp(zs);
altro->SetNNonZSPostsamples(factor);
altro->SetNPretriggerSamples(rate);
altro->SetActiveFECsA((ddl == 0 ? 0x1 : 0x3));
altro->SetActiveFECsB((ddl == 0 ? 0x1 : 0x3));
altro->SetNSamplesPerCh(rate * 128);
AliDebug(5,Form("Writing RCU trailer @ DDL %d w/zs=%d, factor=%d, rate=%d",
ddl, zs > 0, factor, rate));
altro->WriteRCUTrailer(ddl);
}
void
AliFMDRawWriter::ZeroSuppress(Int_t*& data, Int_t nWords,
const Float_t* peds,
const Float_t* noise, UShort_t threshold) const
{
if (threshold <= 0) return;
Bool_t pedSubtract = AliFMDParameters::Instance()->IsZSPedSubtract();
UShort_t pre = AliFMDParameters::Instance()->GetZSPreSamples();
UShort_t post = AliFMDParameters::Instance()->GetZSPostSamples();
Float_t factor = AliFMDParameters::Instance()->GetPedestalFactor();
TArrayC mask(nWords+1);
for (Short_t i = 0; i < nWords; i++) {
Float_t val = data[i] - peds[i] - factor * noise[i];
if (val < 0.5) val = 0;
if (pedSubtract) data[i] = Int_t(val) & 0x3FF;
mask[i] = (val > threshold ? 1 : 0);
AliFMDDebug(10, ("Comparing sample %d %d-%f-%f*%f=%f to %d -> %s",
i, data[i], peds[i], factor, noise[i], val, threshold,
(mask[i] ? "above" : "below")));
}
for (Short_t i = 0; i < nWords; i++) {
if (mask[i]) {
AliFMDDebug(10, ("Sample %d, above", i));
if (i < nWords-1 && !mask[i+1]) {
AliFMDDebug(10, ("At sample %d, next is below, skipping %d to %d",
i, post, i+post));
i += post;
}
continue;
}
Short_t lookahead = TMath::Min(Short_t(nWords), Short_t(i+pre));
AliFMDDebug(10, ("Sample %d, below, look to %d", i, lookahead));
if (mask[lookahead] && pre > 0) {
AliFMDDebug(10, ("Sample %d is a pre-sample to %d", i, lookahead));
i += pre-1;
continue;
}
data[i] = threshold - 1;
}
}
#else
void
AliFMDRawWriter::WriteDigits(TClonesArray* digits)
{
Int_t nDigits = digits->GetEntries();
if (nDigits < 1) return;
AliFMDParameters* pars = AliFMDParameters::Instance();
AliFMDAltroWriter* writer = 0;
Int_t sampleRate = -1;
UShort_t hwaddr = 0;
UShort_t ddl = 0;
std::ofstream* file = 0;
Int_t ret = 0;
for (Int_t i = 0; i < nDigits; i++) {
AliFMDDigit* digit = static_cast<AliFMDDigit*>(digits->At(i));
UInt_t thisDDL, thisHwaddr;
UShort_t det = digit->Detector();
Char_t ring = digit->Ring();
UShort_t sector = digit->Sector();
UShort_t strip = digit->Strip();
if (!pars->Detector2Hardware(det,ring,sector,strip,thisDDL,thisHwaddr)) {
AliError(Form("Failed to get hardware address for FMD%d%c[%2d,%3d]",
det, ring, sector, strip));
continue;
}
AliFMDDebug(40, ("Got DDL=%d and address=%d from FMD%d%c[%2d,%3d]",
thisDDL, thisHwaddr, det, ring, sector, strip));
if (thisHwaddr != hwaddr) {
AliFMDDebug(30, ("Now hardware address 0x%x from FMD%d%c[%2d,%3d] "
"(board 0x%x, chip 0x%x, channel 0x%x)",
thisHwaddr, det, ring, sector, strip,
(thisHwaddr >> 7), (thisHwaddr >> 4) & 0x7,
thisHwaddr & 0xf));
if (writer) writer->AddChannelTrailer(hwaddr);
hwaddr = thisHwaddr;
}
if (ddl != thisDDL) {
if (writer) {
AliFMDDebug(1, ("Closing altro writer %p", writer));
if ((ret = writer->Close()) < 0) {
AliError(Form("Error: %s", writer->ErrorString(ret)));
return;
}
delete writer;
writer = 0;
file->close();
delete file;
file = 0;
}
ddl = thisDDL;
}
if (!writer) {
AliFMDDebug(1, ("Opening new ALTRO writer w/file %s",
AliDAQ::DdlFileName("FMD",ddl)));
file = new std::ofstream(AliDAQ::DdlFileName("FMD",ddl));
if (!file || !*file) {
AliFatal(Form("Failed to open file %s",
AliDAQ::DdlFileName("FMD",ddl)));
return;
}
writer = new AliFMDAltroWriter(*file);
writer->SetThreshold(pars->GetZeroSuppression(det, ring, sector, strip));
}
sampleRate = pars->GetSampleRate(det,ring,sector,strip);
writer->AddSignal(digit->Count1());
if (sampleRate >= 2) writer->AddSignal(digit->Count2());
if (sampleRate >= 3) writer->AddSignal(digit->Count3());
}
if (writer) {
writer->AddChannelTrailer(hwaddr);
writer->Close();
delete writer;
file->close();
delete file;
}
}
#endif