#include <TFolder.h>
#include <TTree.h>
#include <TVirtualMC.h>
+#include <TH1F.h>
+#include <TF1.h>
+#include <TRandom.h>
// --- Standard library ---
AliEMCAL::AliEMCAL():AliDetector()
{
// Default ctor
- fName="EMCAL";
- // fGeom = 0 ;
- fRan = new TRandom(123456) ;
+ fName = "EMCAL" ;
}
//____________________________________________________________________________
AliEMCAL::AliEMCAL(const char* name, const char* title): AliDetector(name,title)
{
// ctor : title is used to identify the layout
- // fGeom = 0;
- fRan = new TRandom(123456) ;
+
+ fTimeMax = 1.28E-5 ;
+ fTimePeak = 2.0E-6 ;
+ fTimeRes = 1.5E-6 ;
+ fHighGainFactor = 40;
+ fHighGainOffset = 0x200 ;
}
//____________________________________________________________________________
// some digitization constants
const Int_t kDDLOffset = 0x800;
- const Double_t kTimeMax = 1.28E-5;
- const Int_t kTimeBins = 256;
- const Double_t kTimePeak = 2.0E-6;
- const Double_t kTimeRes = 1.5E-6;
const Int_t kThreshold = 3;
- const Int_t kHighGainFactor = 40;
- const Int_t kHighGainOffset = 0x200;
- // PHOS has 4 DDL per module; I assume therefore that kChannelsperDDL=896+1 EMCAL channel go to one DDL
const Int_t kChannelsperDDL = 897 ;
AliAltroBuffer* buffer = NULL;
Int_t prevDDL = -1;
- Int_t adcValuesLow[kTimeBins];
- Int_t adcValuesHigh[kTimeBins];
+ Int_t adcValuesLow[fkTimeBins];
+ Int_t adcValuesHigh[fkTimeBins];
// loop over digits (assume ordered digits)
for (Int_t iDigit = 0; iDigit < digits->GetEntries(); iDigit++) {
}
// out of time range signal (?)
- if (digit->GetTime() > kTimeMax) {
+ if (digit->GetTimeR() > fTimeMax) {
buffer->FillBuffer(digit->GetAmp());
- buffer->FillBuffer(kTimeBins); // time bin
+ buffer->FillBuffer(fkTimeBins); // time bin
buffer->FillBuffer(3); // bunch length
buffer->WriteTrailer(3, idDDL, 0, 0); // trailer
} else {
Bool_t highGain = kFALSE;
- // fill time bin values :
- // 1. the signal starts at the time given by the digit
- // 2. the rise is linear and the maximum is reached kTimePeak after start
- // 3. the decay is gaussian with a sigma of kTimeRes
- // 4. the signal is binned into kTimeBins bins
- for (Int_t iTime = 0; iTime < kTimeBins; iTime++) {
- Double_t time = iTime * kTimeMax/kTimeBins;
- Int_t signal = 0;
- if (time < digit->GetTime() + kTimePeak) { // signal is rising
- signal = static_cast<Int_t>((fRan->Rndm() + digit->GetAmp()) *
- (time - digit->GetTime() / kTimePeak) + 0.5);
- } else { // signal is decaying
- signal = static_cast<Int_t>((fRan->Rndm() + digit->GetAmp()) *
- TMath::Gaus(time, digit->GetTime() + kTimePeak, kTimeRes) + 0.5);
- }
- if (signal < 0)
- signal = 0;
- adcValuesLow[iTime] = signal;
- if (signal > 0x3FF) // larger than 10 bits
- adcValuesLow[iTime] = 0x3FF;
- adcValuesHigh[iTime] = static_cast<Int_t>(0.5 + (signal / kHighGainFactor));
- if (adcValuesHigh[iTime] > 0)
- highGain = kTRUE;
- }
-
// write low and eventually high gain channel
- buffer->WriteChannel(idDDL, 0, 0,
- kTimeBins, adcValuesLow, kThreshold);
+ buffer->WriteChannel(iDDL, 0, 0,
+ fkTimeBins, adcValuesLow, kThreshold);
if (highGain) {
- buffer->WriteChannel(idDDL, 0, kHighGainOffset,
- kTimeBins, adcValuesHigh, 1);
+ buffer->WriteChannel(iDDL, 0, fHighGainOffset,
+ fkTimeBins, adcValuesHigh, 1);
}
}
}
-
// write real header and close last file
if (buffer) {
buffer->Flush();
buffer->WriteDataHeader(kFALSE, kFALSE);
delete buffer;
}
-
gime->EmcalLoader()->UnloadDigits();
}
return fLoader;
}
+//__________________________________________________________________
+Double_t AliEMCAL::RawResponseFunction(Double_t *x, Double_t *par)
+{
+ // Shape of the electronics raw reponse:
+ // 1. the signal rises linearly from par[4] to par[1] to reach the maximu par[3]
+ // 2. the signal decays with a gaussian shape for par[4]+par[1] with a sigma of par[2]
+
+ Float_t xx = x[0] ;
+
+ Double_t signal = 0. ;
+
+ if (xx < par[4] + par[1]) // signal is rising
+ signal = (gRandom->Rndm() + par[3]) * (xx - par[4]) / (par[1] - par[4]) ;
+ else // signal is decaying
+ signal = (gRandom->Rndm() + par[3]) * TMath::Gaus(xx, par[4] + par[1], par[2]) ;
+
+ return signal < 0. ? 0. : signal ;
+}
+
+//__________________________________________________________________
+Bool_t AliEMCAL::RawSampledResponse(const Float_t dtime, const Int_t damp, Int_t * adcH, Int_t * adcL) const
+{
+ // for a start time dtime and an amplitude damp given by digit,
+ // calculates the raw sampled response
+
+ const Int_t kRawSignalOverflow = 0x3FF ;
+ Bool_t highGain = kFALSE ;
+
+ TF1 f1("signal", RawResponseFunction, 0, fkTimeBins, 5);
+
+ f1.SetParNames("Time Max", "Peaking time", "Decay width", "Amp max", "Start time") ;
+ f1.SetParameter(0, fTimeMax) ;
+ f1.SetParameter(1, fTimePeak) ;
+ f1.SetParameter(2, fTimeRes) ;
+ f1.SetParameter(3, damp) ;
+ f1.SetParameter(4, dtime) ;
+
+ for (Int_t iTime = 0; iTime < fkTimeBins; iTime++) {
+ Double_t time = iTime * fTimeMax/fkTimeBins;
+ Double_t signal = f1.Eval(time) ;
+ adcL[iTime] = static_cast<Int_t>(signal + 0.5) ;
+ if ( adcL[iTime] > kRawSignalOverflow) // larger than 10 bits
+ adcL[iTime] = kRawSignalOverflow ;
+ adcH[iTime] = static_cast<Int_t>(0.5 + (signal / fHighGainFactor)) ;
+ if (adcH[iTime] > 0)
+ highGain = kTRUE;
+ }
+ return highGain ;
+}
+
//____________________________________________________________________________
void AliEMCAL::SetTreeAddress()
{
{return AliEMCALGeometry::GetInstance(GetTitle(),"") ; }
virtual void Hits2SDigits();
virtual Int_t IsVersion(void) const = 0 ;
+ Int_t GetRawFormatHighGainFactor() const { return fHighGainFactor ; }
+ Int_t GetRawFormatHighGainOffset() const { return fHighGainOffset ; }
+ Int_t GetRawFormatTimeBins() const { return fkTimeBins ; }
+ Double_t GetRawFormatTimeMax() const { return fTimeMax ; }
+ Double_t GetRawFormatTimePeak() const { return fTimePeak ; }
+ Double_t GetRawFormatTimeRes() const { return fTimeRes ; }
virtual AliLoader* MakeLoader(const char* topfoldername);
+ static Double_t RawResponseFunction(Double_t *x, Double_t *par) ;
virtual void SetTreeAddress() ;
virtual const TString Version() const {return TString(" ") ; }
AliEMCAL & operator = (const AliEMCAL & /*rvalue*/) {
Fatal("operator =", "not implemented") ; return *this ; }
-
-
protected:
- //AliEMCALGeometry * fGeom ; // the geometry object
+ Bool_t RawSampledResponse(const Float_t dtime, const Int_t damp, Int_t * adcH, Int_t * adcL) const ;
+
Int_t fBirkC0; // constants for Birk's Law implementation
Double_t fBirkC1; // constants for Birk's Law implementation
Double_t fBirkC2; // constants for Birk's Law implementation
- TRandom * fRan ; //! random number generator
+ Int_t fHighGainFactor ; // High gain attenuation factor of the raw RO signal
+ Int_t fHighGainOffset ; // offset added to the module id to distinguish high and low gain data
+ static const Int_t fkTimeBins = 256 ; // number of sampling bins of the raw RO signal
+ Double_t fTimeMax ; // maximum sampled time of the raw RO signal
+ Double_t fTimePeak ; // peaking time of the raw RO signal
+ Double_t fTimeRes ; // decay rime width of the raw RO signal
ClassDef(AliEMCAL,7) // Electromagnetic calorimeter (base class)
#include <TFolder.h>
#include <TTree.h>
#include <TVirtualMC.h>
+#include <TH1F.h>
+#include <TF1.h>
#include <TRandom.h>
// --- Standard library ---
ClassImp(AliPHOS)
//____________________________________________________________________________
-AliPHOS:: AliPHOS() : AliDetector()
+ AliPHOS:: AliPHOS() : AliDetector()
{
// Default ctor
- fName="PHOS";
+ fName = "PHOS" ;
fQATask = 0;
fTreeQA = 0;
fDebug = 0;
- fRan = new TRandom(123456) ;
}
//____________________________________________________________________________
-AliPHOS::AliPHOS(const char* name, const char* title): AliDetector(name, title)
+AliPHOS::AliPHOS(const char* name, const char* title): AliDetector(name, title)
{
// ctor : title is used to identify the layout
- fQATask = 0;
- fTreeQA = 0;
- fDebug = 0;
- fRan = new TRandom(123456) ;
+ fQATask = 0 ;
+ fTreeQA = 0 ;
+ fDebug = 0 ;
+ fTimeMax = 1.28E-5 ;
+ fTimePeak = 2.0E-6 ;
+ fTimeRes = 1.5E-6 ;
+ fHighGainFactor = 40;
+ fHighGainOffset = 0x200 ;
}
//____________________________________________________________________________
AliPHOS::~AliPHOS()
{
- delete fRan ;
}
//____________________________________________________________________________
}
// some digitization constants
- const Int_t kDDLOffset = 0x600;
- const Double_t kTimeMax = 1.28E-5;
- const Int_t kTimeBins = 256;
- const Double_t kTimePeak = 2.0E-6;
- const Double_t kTimeRes = 1.5E-6;
- const Int_t kThreshold = 3;
- const Int_t kHighGainFactor = 40;
- const Int_t kHighGainOffset = 0x200;
+ const Int_t kDDLOffset = 0x600; // assigned to PHOS
+ const Int_t kThreshold = 3; // skip digits below this threshold
AliAltroBuffer* buffer = NULL;
Int_t prevDDL = -1;
- Int_t adcValuesLow[kTimeBins];
- Int_t adcValuesHigh[kTimeBins];
+ Int_t adcValuesLow[fkTimeBins];
+ Int_t adcValuesHigh[fkTimeBins];
// loop over digits (assume ordered digits)
for (Int_t iDigit = 0; iDigit < digits->GetEntries(); iDigit++) {
}
// out of time range signal (?)
- if (digit->GetTime() > kTimeMax) {
+ if (digit->GetTimeR() > fTimeMax) {
buffer->FillBuffer(digit->GetAmp());
- buffer->FillBuffer(kTimeBins); // time bin
+ buffer->FillBuffer(fkTimeBins); // time bin
buffer->FillBuffer(3); // bunch length
buffer->WriteTrailer(3, relId[3], relId[2], module); // trailer
// simulate linear rise and gaussian decay of signal
} else {
- Bool_t highGain = kFALSE;
-
- // fill time bin values :
- // 1. the signal starts at the time given by the digit
- // 2. the rise is linear and the maximum is reached kTimePeak after start
- // 3. the decay is gaussian with a sigma of kTimeRes
- // 4. the signal is binned into kTimeBins bins
- for (Int_t iTime = 0; iTime < kTimeBins; iTime++) {
- Double_t time = iTime * kTimeMax/kTimeBins;
- Int_t signal = 0;
- if (time < digit->GetTime() + kTimePeak) { // signal is rising
- signal = static_cast<Int_t>(fRan->Rndm() + digit->GetAmp() *
- (time - digit->GetTime()) / kTimePeak);
- } else { // signal is decaying
- signal = static_cast<Int_t>(fRan->Rndm() + digit->GetAmp() *
- TMath::Gaus(time, digit->GetTime() + kTimePeak, kTimeRes));
- }
- if (signal < 0)
- signal = 0;
- adcValuesLow[iTime] = signal;
- if (signal > 0x3FF) // larger than 10 bits
- adcValuesLow[iTime] = 0x3FF;
- adcValuesHigh[iTime] = signal / kHighGainFactor;
- if (adcValuesHigh[iTime] > 0)
- highGain = kTRUE;
- }
-
+ Bool_t highGain = RawSampledResponse(digit->GetTimeR(), digit->GetAmp(), adcValuesHigh, adcValuesLow) ;
+
+
// write low and eventually high gain channel
buffer->WriteChannel(relId[3], relId[2], module,
- kTimeBins, adcValuesLow, kThreshold);
+ fkTimeBins, adcValuesLow, kThreshold);
if (highGain) {
- buffer->WriteChannel(relId[3], relId[2], module + kHighGainOffset,
- kTimeBins, adcValuesHigh, 1);
+ buffer->WriteChannel(relId[3], relId[2], module + fHighGainOffset,
+ fkTimeBins, adcValuesHigh, 1);
}
}
}
-
+
// write real header and close last file
if (buffer) {
buffer->Flush();
buffer->WriteDataHeader(kFALSE, kFALSE);
delete buffer;
}
-
+
gime->PhosLoader()->UnloadDigits();
}
return fLoader;
}
+//__________________________________________________________________
+Double_t AliPHOS::RawResponseFunction(Double_t *x, Double_t *par)
+{
+ // Shape of the electronics raw reponse:
+ // 1. the signal rises linearly from par[4] to par[1] to reach the maximu par[3]
+ // 2. the signal decays with a gaussian shape for par[4]+par[1] with a sigma of par[2]
+
+ Float_t xx = x[0] ;
+
+ Double_t signal = 0. ;
+
+ if (xx < par[4] + par[1]) // signal is rising
+ signal = (gRandom->Rndm() + par[3]) * (xx - par[4]) / (par[1] - par[4]) ;
+ else // signal is decaying
+ signal = (gRandom->Rndm() + par[3]) * TMath::Gaus(xx, par[4] + par[1], par[2]) ;
+
+ return signal < 0. ? 0. : signal ;
+}
+
+//__________________________________________________________________
+Bool_t AliPHOS::RawSampledResponse(const Float_t dtime, const Int_t damp, Int_t * adcH, Int_t * adcL) const
+{
+ // for a start time dtime and an amplitude damp given by digit,
+ // calculates the raw sampled response
+
+ const Int_t kRawSignalOverflow = 0x3FF ;
+ Bool_t highGain = kFALSE ;
+
+ TF1 f1("signal", RawResponseFunction, 0, fkTimeBins, 5);
+
+ f1.SetParNames("Time Max", "Peaking time", "Decay width", "Amp max", "Start time") ;
+ f1.SetParameter(0, fTimeMax) ;
+ f1.SetParameter(1, fTimePeak) ;
+ f1.SetParameter(2, fTimeRes) ;
+ f1.SetParameter(3, damp) ;
+ f1.SetParameter(4, dtime) ;
+
+ for (Int_t iTime = 0; iTime < fkTimeBins; iTime++) {
+ Double_t time = iTime * fTimeMax/fkTimeBins;
+ Double_t signal = f1.Eval(time) ;
+ adcL[iTime] = static_cast<Int_t>(signal + 0.5) ;
+ if ( adcL[iTime] > kRawSignalOverflow) // larger than 10 bits
+ adcL[iTime] = kRawSignalOverflow ;
+ adcH[iTime] = static_cast<Int_t>(0.5 + (signal / fHighGainFactor)) ;
+ if (adcH[iTime] > 0)
+ highGain = kTRUE;
+ }
+ return highGain ;
+}
+
//____________________________________________________________________________
void AliPHOS::SetTreeAddress()
{
{return AliPHOSGeometry::GetInstance(GetTitle(),"") ; }
virtual void Hits2SDigits();
virtual Int_t IsVersion(void) const = 0 ;
+ Int_t GetRawFormatHighGainFactor() const { return fHighGainFactor ; }
+ Int_t GetRawFormatHighGainOffset() const { return fHighGainOffset ; }
+ Int_t GetRawFormatTimeBins() const { return fkTimeBins ; }
+ Double_t GetRawFormatTimeMax() const { return fTimeMax ; }
+ Double_t GetRawFormatTimePeak() const { return fTimePeak ; }
+ Double_t GetRawFormatTimeRes() const { return fTimeRes ; }
virtual AliLoader* MakeLoader(const char* topfoldername);
AliPHOSQAChecker * QAChecker() {return fQATask;}
+ static Double_t RawResponseFunction(Double_t *x, Double_t *par) ;
virtual void SetTreeAddress();
virtual TTree * TreeQA() const {return fTreeQA; }
virtual const TString Version() const {return TString(" ") ; }
AliPHOS & operator = (const AliPHOS & /*rvalue*/) {
Fatal("operator =", "not implemented") ; return *this ; }
+
protected:
-
+
+ Bool_t RawSampledResponse(const Float_t dtime, const Int_t damp, Int_t * adcH, Int_t * adcL) const ;
+
+
AliPHOSQAChecker * fQATask ; //! PHOS checkers container
TTree * fTreeQA ; // the QA tree that contains the alarms
- TRandom * fRan ; //! random number generator
+ Int_t fHighGainFactor ; // High gain attenuation factor of the raw RO signal
+ Int_t fHighGainOffset ; // offset added to the module id to distinguish high and low gain data
+ static const Int_t fkTimeBins = 256 ; // number of sampling bins of the raw RO signal
+ Double_t fTimeMax ; // maximum sampled time of the raw RO signal
+ Double_t fTimePeak ; // peaking time of the raw RO signal
+ Double_t fTimeRes ; // decay rime width of the raw RO signal
- ClassDef(AliPHOS,2) // Photon Spectrometer Detector (base class)
+ ClassDef(AliPHOS,3) // Photon Spectrometer Detector (base class)
} ;
git://git.uio.no / u / mrichter / AliRoot.git / commitdiff