#include <TRandom.h>
#include <TStopwatch.h>
#include <TSystem.h>
+#include <TMath.h>
/// \cond CLASSIMP
ClassImp(AliMUONCDB)
}
}
+//_____________________________________________________________________________
+Double_t GetRandom(Double_t mean, Double_t sigma, Bool_t mustBePositive)
+{
+ Double_t x(-1);
+ if ( mustBePositive )
+ {
+ while ( x < 0 )
+ {
+ x = gRandom->Gaus(mean,sigma);
+ }
+ }
+ else
+ {
+ x = gRandom->Gaus(mean,sigma);
+ }
+ return x;
+}
+
}
//_____________________________________________________________________________
for ( UInt_t timeStamp = 0; timeStamp < 60*15; timeStamp += 120 )
{
Float_t value = 1500;
- if (!defaultValues) value = gRandom->Gaus(1750,62.5);
+ if (!defaultValues) value = GetRandom(1750,62.5,true);
AliDCSValue* dcsValue = new AliDCSValue(value,timeStamp);
valueSet->Add(dcsValue);
}
}
else
{
- meanPedestal = -1;
- while ( meanPedestal < 0 )
- {
- meanPedestal = gRandom->Gaus(kPedestalMeanMean,kPedestalMeanSigma);
- }
- sigmaPedestal = -1;
- while ( sigmaPedestal < 0 )
- {
- sigmaPedestal = gRandom->Gaus(kPedestalSigmaMean,kPedestalSigmaSigma);
- }
+ Bool_t positive(kTRUE);
+ meanPedestal = GetRandom(kPedestalMeanMean,kPedestalMeanSigma,positive);
+ sigmaPedestal = GetRandom(kPedestalSigmaMean,kPedestalSigmaSigma,positive);
}
ped->SetValueAsFloat(manuChannel,0,meanPedestal);
ped->SetValueAsFloat(manuChannel,1,sigmaPedestal);
}
else
{
- capaValue = -1;
- while ( capaValue < 0 )
- {
- capaValue = gRandom->Gaus(kCapaMean,kCapaSigma);
- }
+ capaValue = GetRandom(kCapaMean,kCapaSigma,kTRUE);
}
capa->SetValueAsFloat(manuChannel,0,capaValue);
}
Int_t
AliMUONCDB::MakeGainStore(AliMUONVStore& gainStore, Bool_t defaultValues)
{
- /// Create a gain store. if defaultValues=true, all gain are 1.0,
- /// otherwise they are from a gaussian with parameters defined in the
- /// kGain* constants below.
+ /// Create a gain store. if defaultValues=true, all gains set so that
+ /// charge = (adc-ped)
+ ///
+ /// otherwise parameters are taken from gaussians with parameters
+ /// defined in the k* constants below.
TIter next(ManuList());
Int_t nchannels(0);
Int_t nmanus(0);
- const Double_t kSaturation(3000);
- const Double_t kGainMean(1.0);
- const Double_t kGainSigma(0.05);
+ const Int_t kSaturation(3000);
+ const Double_t kA0Mean(1.2);
+ const Double_t kA0Sigma(0.1);
+ const Double_t kA1Mean(1E-5);
+ const Double_t kA1Sigma(1E-6);
+ const Double_t kQualMean(0xFF);
+ const Double_t kQualSigma(0x10);
+ const Int_t kThresMean(1600);
+ const Int_t kThresSigma(100);
while ( ( p = (AliMpIntPair*)next() ) )
{
Int_t manuId = p->GetSecond();
AliMUONVCalibParam* gain =
- new AliMUONCalibParamNF(2,AliMpConstants::ManuNofChannels(),
+ new AliMUONCalibParamNF(5,AliMpConstants::ManuNofChannels(),
detElemId,
manuId,
AliMUONVCalibParam::InvalidFloatValue());
++nchannels;
- Float_t meanGain;
- Float_t saturation(kSaturation);
-
if ( defaultValues )
{
- meanGain = 1.0;
+ gain->SetValueAsFloat(manuChannel,0,1.0);
+ gain->SetValueAsFloat(manuChannel,1,0.0);
+ gain->SetValueAsInt(manuChannel,2,4095);
+ gain->SetValueAsInt(manuChannel,3,1);
+ gain->SetValueAsInt(manuChannel,4,kSaturation);
}
else
{
- meanGain = -1;
- while ( meanGain < 0 )
- {
- meanGain = gRandom->Gaus(kGainMean,kGainSigma);
- }
+ Bool_t positive(kTRUE);
+ gain->SetValueAsFloat(manuChannel,0,GetRandom(kA0Mean,kA0Sigma,positive));
+ gain->SetValueAsFloat(manuChannel,1,GetRandom(kA1Mean,kA1Sigma,!positive));
+ gain->SetValueAsInt(manuChannel,2,(Int_t)TMath::Nint(GetRandom(kThresMean,kThresSigma,positive)));
+ gain->SetValueAsInt(manuChannel,3,(Int_t)TMath::Nint(GetRandom(kQualMean,kQualSigma,positive)));
+ gain->SetValueAsInt(manuChannel,4,kSaturation);
}
- gain->SetValueAsFloat(manuChannel,0,meanGain);
- gain->SetValueAsFloat(manuChannel,1,saturation);
}
Bool_t ok = gainStore.Add(gain);
{
/// Calibrate one digit
+ if ( digit.IsCalibrated() )
+ {
+ fLogger->Log("ERROR : trying to calibrate a digit twice");
+ return;
+ }
+
AliMUONVCalibParam* deadmap = static_cast<AliMUONVCalibParam*>
(fStatusMap->FindObject(digit.DetElemId(),digit.ManuId()));
Int_t statusMap = deadmap->ValueAsInt(digit.ManuChannel());
Int_t manuChannel = digit.ManuChannel();
Float_t adc = digit.ADC();
Float_t padc = adc-pedestal->ValueAsFloat(manuChannel,0);
- if ( padc < 3.0*pedestal->ValueAsFloat(manuChannel,1) )
+ Float_t charge(0);
+ if ( padc > 3.0*pedestal->ValueAsFloat(manuChannel,1) )
{
- padc = 0.0;
+ Float_t a0 = gain->ValueAsFloat(manuChannel,0);
+ Float_t a1 = gain->ValueAsFloat(manuChannel,1);
+ Int_t thres = gain->ValueAsInt(manuChannel,2);
+ if ( padc < thres )
+ {
+ charge = a0*padc;
+ }
+ else
+ {
+ charge = a0*thres + a0*padc + a1*padc*padc;
+ }
}
- Float_t charge = padc*gain->ValueAsFloat(manuChannel,0);
digit.SetCharge(charge);
- Int_t saturation = gain->ValueAsInt(manuChannel,1);
+ Int_t saturation = gain->ValueAsInt(manuChannel,4);
if ( charge >= saturation )
{
digit.Saturated(kTRUE);
/// For tracking digits, starting from an ideal digit's charge, we :
///
/// - add some noise (thus leading to a realistic charge), if requested to do so
- /// - divide by a gain (thus decalibrating the digit)
+ /// - "divide" by a gain (thus decalibrating the digit)
/// - add a pedestal (thus decalibrating the digit)
/// - sets the signal to zero if below 3*sigma of the noise
static const Int_t kMaxADC = (1<<12)-1; // We code the charge on a 12 bits ADC.
- Float_t signal = digit.Charge();
+ Float_t charge = digit.Charge();
+
+ // We set the charge to 0, as the only relevant piece of information
+ // after Digitization is the ADC value.
+ digit.SetCharge(0);
if ( !addNoise )
{
- digit.SetADC(TMath::Min(kMaxADC,TMath::Nint(signal)));
+ digit.SetADC(TMath::Min(kMaxADC,TMath::Nint(charge)));
return;
}
fLogger->Log(Form("%s:%d:Could not get pedestal for DE=%4d manuId=%4d. Disabling.",
__FILE__,__LINE__,
detElemId,manuId));
- digit.SetCharge(0);
digit.SetADC(0);
return;
}
fLogger->Log(Form("%s:%d:Could not get gain for DE=%4d manuId=%4d. Disabling.",
__FILE__,__LINE__,
detElemId,manuId));
- digit.SetCharge(0);
digit.SetADC(0);
return;
}
- Float_t gainMean = gain->ValueAsFloat(manuChannel,0);
-
- Float_t adcNoise = gRandom->Gaus(0.0,pedestalSigma);
+
+ Float_t a0 = gain->ValueAsFloat(manuChannel,0);
+ Float_t a1 = gain->ValueAsFloat(manuChannel,1);
+
+ Int_t thres = gain->ValueAsInt(manuChannel,2);
+ Float_t chargeThres = a0*thres;
- Int_t adc;
+ Float_t padc(0); // (adc - ped) value
- if ( gainMean < 1E-6 )
- {
- AliError(Form("Got a too small gain %e for DE=%d manuId=%d manuChannel=%d. "
- "Setting signal to 0.",
- gainMean,detElemId,manuId,manuChannel));
- adc = 0;
- }
- else
+ if ( charge <= chargeThres || TMath::Abs(a1) < 1E-12 )
{
- adc = TMath::Nint( signal / gainMean + pedestalMean + adcNoise);///
+ // linear part only
- if ( adc <= pedestalMean + fgkNSigmas*pedestalSigma )
+ if ( TMath::Abs(a0) > 1E-12 )
{
- adc = 0;
+ padc = charge/a0;
}
}
+ else // charge > chargeThres && a1 not zero
+ {
+ // parabolic part
+ padc = TMath::Sqrt((chargeThres-charge)/a1) + thres;
+ }
+
+ Float_t adcNoise = gRandom->Gaus(0.0,pedestalSigma);
+ Int_t adc(0);
+
+ if ( padc > 0 )
+ {
+ adc = TMath::Nint(padc + pedestalMean + adcNoise);
+ }
// be sure we stick to 12 bits.
if ( adc > kMaxADC )
adc = kMaxADC;
}
- digit.SetCharge(adc);
+ AliDebug(3,Form("DE %4d Manu %4d Ch %2d Charge %e A0 %e A1 %e Thres %d padc %e ADC %4d",
+ detElemId,manuId,manuChannel,charge,
+ gain->ValueAsFloat(manuChannel,0),
+ gain->ValueAsFloat(manuChannel,1),
+ gain->ValueAsInt(manuChannel,2),
+ padc,
+ adc));
+
digit.SetADC(adc);
}
{
ApplyResponseToTriggerDigit(store,*digit);
}
- if ( digit->Charge() > 0 )
+ if ( digit->ADC() > 0 || digit->Charge() > 0 )
{
filteredStore.Add(*digit,AliMUONVDigitStore::kIgnore);
}
d->SetCharge(TMath::Nint(ped+pedestalMean+0.5));
d->NoiseOnly(kTRUE);
ApplyResponseToTrackerDigit(*d,kFALSE);
- if ( d->Charge() > 0 )
+ if ( d->ADC() > 0 )
{
Bool_t ok = digitStore.Add(*d,AliMUONVDigitStore::kDeny);
// this can happen (that we randomly chose a digit that is
///
/// Implementation of AliMUONVSubprocessor class to deal with MUON TRK Gains.
///
-/// Gains are read in from an ascii file, with the format : \n
-///---------------------------------------------------------------------------\n
-///BUS_PATCH MANU CHANNEL Ped. a0 a1 a2 xlim P(chi2) P(chi2)_2 \n
-///---------------------------------------------------------------------------\n
+/// Gains are read in from an ascii file, with the format :
+///
+///---------------------------------------------------------------------------
+///
+///BUS_PATCH MANU CHANNEL a0 a1 thres Qual
+///
+///---------------------------------------------------------------------------
///
/// \author L. Aphecetche
///
Int_t
AliMUONGainSubprocessor::ReadFile(const char* filename)
{
- /// Read the Gains from an ASCII file. \n
- /// Format of that file is one line per channel : \n
- ///-------------------------------------------------------------------------\n
- ///BUS_PATCH MANU CHANNEL Ped. a0 a1 a2 xlim P(chi2) P(chi2)_2 \n
- ///-------------------------------------------------------------------------\n
- /// \n
- /// Return kFALSE if reading was not successfull. \n
+ /// Read the Gains from an ASCII file.
+ /// Format of that file is one line per channel :
+ ///-------------------------------------------------------------------------
+ ///BUS_PATCH MANU CHANNEL a0 a1 thres Qual
+ ///-------------------------------------------------------------------------
+ ///
+ /// Return kFALSE if reading was not successfull.
///
TString sFilename(gSystem->ExpandPathName(filename));
}
char line[256];
Int_t busPatchID, manuID, manuChannel;
- Float_t ped;
- Float_t a0, a1, a2;
- Float_t xlim;
- Float_t chi2, chi22;
+ Float_t a0, a1;
+ Int_t thres;
+ Int_t qual;
+ const Int_t kSaturation(3000); // FIXME: how to get this number ?
static const Int_t kNchannels(AliMpConstants::ManuNofChannels());
Int_t n(0);
{
if ( strlen(line) < 10 ) continue;
if ( line[0] == '/' && line[1] == '/' ) continue;
- std::istringstream sin(line);
+
+ sscanf(line,"%d %d %d %f %f %d %x",&busPatchID,&manuID,&manuChannel,
+ &a0,&a1,&thres,&qual);
AliDebug(3,Form("line=%s",line));
- sin >> busPatchID >> manuID >> manuChannel >> ped >> a0 >> a1 >> a2 >> xlim >> chi2 >> chi22;
Int_t detElemID = AliMpDDLStore::Instance()->GetDEfromBus(busPatchID);
- AliDebug(3,Form("BUSPATCH %3d DETELEMID %4d MANU %3d CH %3d PED %7.2f A0 %7.2f A1 %7.2f A2 %7.2f"
- " XLIM %7.2f CHI2 %7.2f CHI22 %7.2f",
- busPatchID,detElemID,manuID,manuChannel,ped,a0,a1,a2,xlim,chi2,chi22));
- if ( a0==a1 && a1==a2 && a0==-2) continue;
+ AliDebug(3,Form("BUSPATCH %3d DETELEMID %4d MANU %3d CH %3d A0 %7.2f "
+ "A1 %e THRES %5d QUAL %x",
+ busPatchID,detElemID,manuID,manuChannel,a0,a1,thres,qual));
+ if ( qual == 0 ) continue;
AliMUONVCalibParam* gain =
static_cast<AliMUONVCalibParam*>(fGains->FindObject(detElemID,manuID));
if (!gain)
{
- gain = new AliMUONCalibParamNF(6,detElemID,manuID,kNchannels,0);
+ gain = new AliMUONCalibParamNF(5,kNchannels,detElemID,manuID,0);
fGains->Add(gain);
}
gain->SetValueAsFloat(manuChannel,0,a0);
gain->SetValueAsFloat(manuChannel,1,a1);
- gain->SetValueAsFloat(manuChannel,2,a2);
- gain->SetValueAsFloat(manuChannel,3,xlim);
- gain->SetValueAsFloat(manuChannel,4,chi2);
- gain->SetValueAsFloat(manuChannel,5,chi22);
+ gain->SetValueAsInt(manuChannel,2,thres);
+ gain->SetValueAsInt(manuChannel,3,qual);
+ gain->SetValueAsInt(manuChannel,4,kSaturation);
++n;
}
in.close();
git://git.uio.no / u / mrichter / AliRoot.git / commitdiff