X-Git-Url: http://git.uio.no/git/?a=blobdiff_plain;f=TRD%2FAliTRDmcmSim.cxx;h=d6536c4ac4f94ec3172d54befc1dd1f09e04827d;hb=2f65cb9a40e871dedd8fad53a3fa3700b703d04c;hp=2ec8488a533bedecc2627330b3456003d06826de;hpb=6b0948679b4140e9ad2ef54816a6490029605c01;p=u%2Fmrichter%2FAliRoot.git diff --git a/TRD/AliTRDmcmSim.cxx b/TRD/AliTRDmcmSim.cxx index 2ec8488a533..d6536c4ac4f 100644 --- a/TRD/AliTRDmcmSim.cxx +++ b/TRD/AliTRDmcmSim.cxx @@ -35,31 +35,34 @@ #include "TRandom.h" #include "TClonesArray.h" #include "TMath.h" +#include #include "AliLog.h" #include "AliRunLoader.h" #include "AliLoader.h" #include "AliTRDfeeParam.h" +#include "AliTRDcalibDB.h" #include "AliTRDtrapConfig.h" #include "AliTRDdigitsManager.h" #include "AliTRDarrayADC.h" #include "AliTRDarrayDictionary.h" #include "AliTRDtrackletMCM.h" #include "AliTRDmcmSim.h" +#include "TTreeStream.h" ClassImp(AliTRDmcmSim) Bool_t AliTRDmcmSim::fgApplyCut = kTRUE; Int_t AliTRDmcmSim::fgAddBaseline = 0; +Bool_t AliTRDmcmSim::fgStoreClusters = kFALSE; -const Int_t AliTRDmcmSim::fgkFormatIndex = std::ios_base::xalloc(); +const Int_t AliTRDmcmSim::fgkFormatIndex = std::ios_base::xalloc(); -const Int_t AliTRDmcmSim::fgkNADC = AliTRDfeeParam::GetNadcMcm(); -const UShort_t AliTRDmcmSim::fgkFPshifts[4] = {11, 14, 17, 21}; +const UShort_t AliTRDmcmSim::fgkFPshifts[4] = {11, 14, 17, 21}; -AliTRDmcmSim::AliTRDmcmSim() : +AliTRDmcmSim::AliTRDmcmSim() : TObject(), fInitialized(kFALSE), fDetector(-1), @@ -82,22 +85,31 @@ AliTRDmcmSim::AliTRDmcmSim() : fTailAmplLong(NULL), fTailAmplShort(NULL), fNHits(0), - fFitReg(NULL) + fFitReg(NULL), + fDebugStream(0x0) { // // AliTRDmcmSim default constructor // By default, nothing is initialized. // It is necessary to issue Init before use. + + for (Int_t iDict = 0; iDict < 3; iDict++) + fDict[iDict] = 0x0; + + fFitPtr[0] = 0; + fFitPtr[1] = 0; + fFitPtr[2] = 0; + fFitPtr[3] = 0; } -AliTRDmcmSim::~AliTRDmcmSim() +AliTRDmcmSim::~AliTRDmcmSim() { // // AliTRDmcmSim destructor // if(fInitialized) { - for( Int_t iAdc = 0 ; iAdc < fgkNADC; iAdc++ ) { + for( Int_t iAdc = 0 ; iAdc < AliTRDfeeParam::GetNadcMcm(); iAdc++ ) { delete [] fADCR[iAdc]; delete [] fADCF[iAdc]; } @@ -105,57 +117,57 @@ AliTRDmcmSim::~AliTRDmcmSim() delete [] fADCF; delete [] fZSMap; delete [] fMCMT; - + delete [] fPedAcc; delete [] fGainCounterA; delete [] fGainCounterB; delete [] fTailAmplLong; delete [] fTailAmplShort; delete [] fFitReg; - + fTrackletArray->Delete(); delete fTrackletArray; } } -void AliTRDmcmSim::Init( Int_t det, Int_t robPos, Int_t mcmPos, Bool_t /* newEvent */ ) +void AliTRDmcmSim::Init( Int_t det, Int_t robPos, Int_t mcmPos, Bool_t /* newEvent */ ) { // // Initialize the class with new MCM position information // memory is allocated in the first initialization // - + if (!fInitialized) { fFeeParam = AliTRDfeeParam::Instance(); - fTrapConfig = AliTRDtrapConfig::Instance(); + fTrapConfig = AliTRDcalibDB::Instance()->GetTrapConfig(); } fDetector = det; fRobPos = robPos; fMcmPos = mcmPos; - fNTimeBin = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kC13CPUA); fRow = fFeeParam->GetPadRowFromMCM( fRobPos, fMcmPos ); - + if (!fInitialized) { - fADCR = new Int_t *[fgkNADC]; - fADCF = new Int_t *[fgkNADC]; - fZSMap = new Int_t [fgkNADC]; - fGainCounterA = new UInt_t[fgkNADC]; - fGainCounterB = new UInt_t[fgkNADC]; - for( Int_t iAdc = 0 ; iAdc < fgkNADC; iAdc++ ) { + fADCR = new Int_t *[AliTRDfeeParam::GetNadcMcm()]; + fADCF = new Int_t *[AliTRDfeeParam::GetNadcMcm()]; + fZSMap = new Int_t [AliTRDfeeParam::GetNadcMcm()]; + fGainCounterA = new UInt_t[AliTRDfeeParam::GetNadcMcm()]; + fGainCounterB = new UInt_t[AliTRDfeeParam::GetNadcMcm()]; + fNTimeBin = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kC13CPUA, fDetector, fRobPos, fMcmPos); + for( Int_t iAdc = 0 ; iAdc < AliTRDfeeParam::GetNadcMcm(); iAdc++ ) { fADCR[iAdc] = new Int_t[fNTimeBin]; fADCF[iAdc] = new Int_t[fNTimeBin]; } - + // filter registers - fPedAcc = new UInt_t[fgkNADC]; // accumulator for pedestal filter - fTailAmplLong = new UShort_t[fgkNADC]; - fTailAmplShort = new UShort_t[fgkNADC]; - + fPedAcc = new UInt_t[AliTRDfeeParam::GetNadcMcm()]; // accumulator for pedestal filter + fTailAmplLong = new UShort_t[AliTRDfeeParam::GetNadcMcm()]; + fTailAmplShort = new UShort_t[AliTRDfeeParam::GetNadcMcm()]; + // tracklet calculation - fFitReg = new FitReg_t[fgkNADC]; + fFitReg = new FitReg_t[AliTRDfeeParam::GetNadcMcm()]; fTrackletArray = new TClonesArray("AliTRDtrackletMCM", fgkMaxTracklets); - + fMCMT = new UInt_t[fgkMaxTracklets]; } @@ -169,10 +181,10 @@ void AliTRDmcmSim::Reset() // Resets the data values and internal filter registers // by re-initialising them - if( !CheckInitialized() ) + if( !CheckInitialized() ) return; - for( Int_t iAdc = 0 ; iAdc < fgkNADC; iAdc++ ) { + for( Int_t iAdc = 0 ; iAdc < AliTRDfeeParam::GetNadcMcm(); iAdc++ ) { for( Int_t it = 0 ; it < fNTimeBin ; it++ ) { fADCR[iAdc][it] = 0; fADCF[iAdc][it] = 0; @@ -181,41 +193,41 @@ void AliTRDmcmSim::Reset() fGainCounterA[iAdc] = 0; fGainCounterB[iAdc] = 0; } - + for(Int_t i = 0; i < fgkMaxTracklets; i++) { fMCMT[i] = 0; } for (Int_t iDict = 0; iDict < 3; iDict++) fDict[iDict] = 0x0; - + FilterPedestalInit(); FilterGainInit(); FilterTailInit(); } -void AliTRDmcmSim::SetNTimebins(Int_t ntimebins) +void AliTRDmcmSim::SetNTimebins(Int_t ntimebins) { - // Reallocate memory if a change in the number of timebins + // Reallocate memory if a change in the number of timebins // is needed (should not be the case for real data) - if( !CheckInitialized() ) + if( !CheckInitialized() ) return; fNTimeBin = ntimebins; - for( Int_t iAdc = 0 ; iAdc < fgkNADC; iAdc++ ) { - delete fADCR[iAdc]; - delete fADCF[iAdc]; + for( Int_t iAdc = 0 ; iAdc < AliTRDfeeParam::GetNadcMcm(); iAdc++ ) { + delete [] fADCR[iAdc]; + delete [] fADCF[iAdc]; fADCR[iAdc] = new Int_t[fNTimeBin]; fADCF[iAdc] = new Int_t[fNTimeBin]; } } -Bool_t AliTRDmcmSim::LoadMCM(AliRunLoader* const runloader, Int_t det, Int_t rob, Int_t mcm) +Bool_t AliTRDmcmSim::LoadMCM(AliRunLoader* const runloader, Int_t det, Int_t rob, Int_t mcm) { // loads the ADC data as obtained from the digitsManager for the specified MCM. - // This method is meant for rare execution, e.g. in the visualization. When called - // frequently use SetData(...) instead. + // This method is meant for rare execution, e.g. in the visualization. When called + // frequently use SetData(...) instead. Init(det, rob, mcm); @@ -248,35 +260,35 @@ Bool_t AliTRDmcmSim::LoadMCM(AliRunLoader* const runloader, Int_t det, Int_t rob SetData(digits); } - else + else retval = kFALSE; - + delete digMgr; - + return retval; } void AliTRDmcmSim::NoiseTest(Int_t nsamples, Int_t mean, Int_t sigma, Int_t inputGain, Int_t inputTail) { - // This function can be used to test the filters. + // This function can be used to test the filters. // It feeds nsamples of ADC values with a gaussian distribution specified by mean and sigma. // The filter chain implemented here consists of: // Pedestal -> Gain -> Tail - // With inputGain and inputTail the input to the gain and tail filter, respectively, - // can be chosen where + // With inputGain and inputTail the input to the gain and tail filter, respectively, + // can be chosen where // 0: noise input // 1: pedestal output // 2: gain output - // The input has to be chosen from a stage before. - // The filter behaviour is controlled by the TRAP parameters from AliTRDtrapConfig in the + // The input has to be chosen from a stage before. + // The filter behaviour is controlled by the TRAP parameters from AliTRDtrapConfig in the // same way as in normal simulation. // The functions produces four histograms with the values at the different stages. - if( !CheckInitialized() ) + if( !CheckInitialized() ) return; TString nameInputGain; - TString nameInputTail; + TString nameInputTail; switch (inputGain) { case 0: @@ -313,46 +325,46 @@ void AliTRDmcmSim::NoiseTest(Int_t nsamples, Int_t mean, Int_t sigma, Int_t inpu TH1F *h = new TH1F("noise", "Gaussian Noise;sample;ADC count", nsamples, 0, nsamples); TH1F *hfp = new TH1F("ped", "Noise #rightarrow Pedestal filter;sample;ADC count", nsamples, 0, nsamples); - TH1F *hfg = new TH1F("gain", - (nameInputGain + "#rightarrow Gain;sample;ADC count").Data(), + TH1F *hfg = new TH1F("gain", + (nameInputGain + "#rightarrow Gain;sample;ADC count").Data(), nsamples, 0, nsamples); - TH1F *hft = new TH1F("tail", - (nameInputTail + "#rightarrow Tail;sample;ADC count").Data(), + TH1F *hft = new TH1F("tail", + (nameInputTail + "#rightarrow Tail;sample;ADC count").Data(), nsamples, 0, nsamples); h->SetStats(kFALSE); hfp->SetStats(kFALSE); hfg->SetStats(kFALSE); hft->SetStats(kFALSE); - + Int_t value; // ADC count with noise (10 bit) Int_t valuep; // pedestal filter output (12 bit) Int_t valueg; // gain filter output (12 bit) Int_t valuet; // tail filter value (12 bit) - + for (Int_t i = 0; i < nsamples; i++) { - value = (Int_t) gRandom->Gaus(mean, sigma); // generate noise with gaussian distribution + value = (Int_t) gRandom->Gaus(mean, sigma); // generate noise with gaussian distribution h->SetBinContent(i, value); valuep = FilterPedestalNextSample(1, 0, ((Int_t) value) << 2); - + if (inputGain == 0) valueg = FilterGainNextSample(1, ((Int_t) value) << 2); - else - valueg = FilterGainNextSample(1, valuep); - + else + valueg = FilterGainNextSample(1, valuep); + if (inputTail == 0) valuet = FilterTailNextSample(1, ((Int_t) value) << 2); else if (inputTail == 1) - valuet = FilterTailNextSample(1, valuep); + valuet = FilterTailNextSample(1, valuep); else - valuet = FilterTailNextSample(1, valueg); + valuet = FilterTailNextSample(1, valueg); hfp->SetBinContent(i, valuep >> 2); hfg->SetBinContent(i, valueg >> 2); hft->SetBinContent(i, valuet >> 2); } - TCanvas *c = new TCanvas; + TCanvas *c = new TCanvas; c->Divide(2,2); c->cd(1); h->Draw(); @@ -370,7 +382,7 @@ Bool_t AliTRDmcmSim::CheckInitialized() const // Check whether object is initialized // - if( ! fInitialized ) + if( ! fInitialized ) AliError(Form ("AliTRDmcmSim is not initialized but function other than Init() is called.")); return fInitialized; @@ -379,16 +391,16 @@ Bool_t AliTRDmcmSim::CheckInitialized() const void AliTRDmcmSim::Print(Option_t* const option) const { // Prints the data stored and/or calculated for this MCM. - // The output is controlled by option which can be a sequence of any of + // The output is controlled by option which can be a sequence of any of // the following characters: // R - prints raw ADC data - // F - prints filtered data + // F - prints filtered data // H - prints detected hits // T - prints found tracklets - // The later stages are only meaningful after the corresponding calculations + // The later stages are only meaningful after the corresponding calculations // have been performed. - if ( !CheckInitialized() ) + if ( !CheckInitialized() ) return; printf("MCM %i on ROB %i in detector %i\n", fMcmPos, fRobPos, fDetector); @@ -414,39 +426,39 @@ void AliTRDmcmSim::Print(Option_t* const option) const } } -void AliTRDmcmSim::Draw(Option_t* const option) +void AliTRDmcmSim::Draw(Option_t* const option) { // Plots the data stored in a 2-dim. timebin vs. ADC channel plot. - // The option selects what data is plotted and can be a sequence of + // The option selects what data is plotted and can be a sequence of // the following characters: // R - plot raw data (default) // F - plot filtered data (meaningless if R is specified) // In addition to the ADC values: - // H - plot hits + // H - plot hits // T - plot tracklets - if( !CheckInitialized() ) + if( !CheckInitialized() ) return; TString opt = option; TH2F *hist = new TH2F("mcmdata", Form("Data of MCM %i on ROB %i in detector %i", \ fMcmPos, fRobPos, fDetector), \ - fgkNADC, -0.5, fgkNADC-.5, fNTimeBin, -.5, fNTimeBin-.5); + AliTRDfeeParam::GetNadcMcm(), -0.5, AliTRDfeeParam::GetNadcMcm()-.5, fNTimeBin, -.5, fNTimeBin-.5); hist->GetXaxis()->SetTitle("ADC Channel"); hist->GetYaxis()->SetTitle("Timebin"); hist->SetStats(kFALSE); if (opt.Contains("R")) { for (Int_t iTimeBin = 0; iTimeBin < fNTimeBin; iTimeBin++) { - for (Int_t iAdc = 0; iAdc < fgkNADC; iAdc++) { + for (Int_t iAdc = 0; iAdc < AliTRDfeeParam::GetNadcMcm(); iAdc++) { hist->SetBinContent(iAdc+1, iTimeBin+1, fADCR[iAdc][iTimeBin] >> fgkAddDigits); } } } else { for (Int_t iTimeBin = 0; iTimeBin < fNTimeBin; iTimeBin++) { - for (Int_t iAdc = 0; iAdc < fgkNADC; iAdc++) { + for (Int_t iAdc = 0; iAdc < AliTRDfeeParam::GetNadcMcm(); iAdc++) { hist->SetBinContent(iAdc+1, iTimeBin+1, fADCF[iAdc][iTimeBin] >> fgkAddDigits); } } @@ -456,8 +468,8 @@ void AliTRDmcmSim::Draw(Option_t* const option) if (opt.Contains("H")) { TGraph *grHits = new TGraph(); for (Int_t iHit = 0; iHit < fNHits; iHit++) { - grHits->SetPoint(iHit, - fHits[iHit].fChannel + 1 + fHits[iHit].fYpos/256., + grHits->SetPoint(iHit, + fHits[iHit].fChannel + 1 + fHits[iHit].fYpos/256., fHits[iHit].fTimebin); } grHits->Draw("*"); @@ -469,11 +481,13 @@ void AliTRDmcmSim::Draw(Option_t* const option) AliTRDtrackletMCM *trkl = (AliTRDtrackletMCM*) (*fTrackletArray)[iTrkl]; Float_t padWidth = 0.635 + 0.03 * (fDetector % 6); Float_t offset = padWidth/256. * ((((((fRobPos & 0x1) << 2) + (fMcmPos & 0x3)) * 18) << 8) - ((18*4*2 - 18*2 - 3) << 7)); // revert adding offset in FitTracklet - Int_t ndrift = fTrapConfig->GetDmemUnsigned(0xc025, fDetector, fRobPos, fMcmPos) >> 5; - Float_t slope = trkl->GetdY() * 140e-4 / ndrift; + Int_t ndrift = fTrapConfig->GetDmemUnsigned(fgkDmemAddrNdrift, fDetector, fRobPos, fMcmPos) >> 5; + Float_t slope = 0; + if (ndrift) + slope = trkl->GetdY() * 140e-4 / ndrift; - Int_t t0 = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPFS); - Int_t t1 = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPFE); + Int_t t0 = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPFS, fDetector, fRobPos, fMcmPos); + Int_t t1 = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPFE, fDetector, fRobPos, fMcmPos); trklLines[iTrkl].SetX1((offset - (trkl->GetY() - slope * t0)) / padWidth); // ??? sign? trklLines[iTrkl].SetY1(t0); @@ -487,7 +501,7 @@ void AliTRDmcmSim::Draw(Option_t* const option) } } -void AliTRDmcmSim::SetData( Int_t adc, Int_t* const data ) +void AliTRDmcmSim::SetData( Int_t adc, const Int_t* const data ) { // // Store ADC data into array of raw data @@ -495,8 +509,8 @@ void AliTRDmcmSim::SetData( Int_t adc, Int_t* const data ) if( !CheckInitialized() ) return; - if( adc < 0 || adc >= fgkNADC ) { - AliError(Form ("Error: ADC %i is out of range (0 .. %d).", adc, fgkNADC-1)); + if( adc < 0 || adc >= AliTRDfeeParam::GetNadcMcm() ) { + AliError(Form ("Error: ADC %i is out of range (0 .. %d).", adc, AliTRDfeeParam::GetNadcMcm()-1)); return; } @@ -514,8 +528,8 @@ void AliTRDmcmSim::SetData( Int_t adc, Int_t it, Int_t data ) if( !CheckInitialized() ) return; - if( adc < 0 || adc >= fgkNADC ) { - AliError(Form ("Error: ADC %i is out of range (0 .. %d).", adc, fgkNADC-1)); + if( adc < 0 || adc >= AliTRDfeeParam::GetNadcMcm() ) { + AliError(Form ("Error: ADC %i is out of range (0 .. %d).", adc, AliTRDfeeParam::GetNadcMcm()-1)); return; } @@ -527,7 +541,7 @@ void AliTRDmcmSim::SetData(AliTRDarrayADC* const adcArray, AliTRDdigitsManager * { // Set the ADC data from an AliTRDarrayADC - if( !CheckInitialized() ) + if( !CheckInitialized() ) return; fDigitsManager = digitsManager; @@ -535,37 +549,41 @@ void AliTRDmcmSim::SetData(AliTRDarrayADC* const adcArray, AliTRDdigitsManager * for (Int_t iDict = 0; iDict < 3; iDict++) { AliTRDarrayDictionary *newDict = (AliTRDarrayDictionary*) fDigitsManager->GetDictionary(fDetector, iDict); if (fDict[iDict] != 0x0 && newDict != 0x0) { - + if (fDict[iDict] == newDict) continue; fDict[iDict] = newDict; - - if (fDict[iDict]->GetDim() == 0) { - AliError(Form("Dictionary %i of det. %i has dim. 0", fDetector, iDict)); - continue; - } - fDict[iDict]->Expand(); + if(fDict[iDict]->GetDim() != 0) + fDict[iDict]->Expand(); } else { fDict[iDict] = newDict; - if (fDict[iDict]) + if (fDict[iDict] && (fDict[iDict]->GetDim() != 0) ) fDict[iDict]->Expand(); } + + // If there is no data, set dictionary to zero to avoid crashes + if (fDict[iDict]->GetDim() == 0) { + // AliError(Form("Dictionary %i of det. %i has dim. 0", iDict, fDetector)); + fDict[iDict] = 0x0; + } } } if (fNTimeBin != adcArray->GetNtime()) SetNTimebins(adcArray->GetNtime()); - + Int_t offset = (fMcmPos % 4 + 1) * 21 + (fRobPos % 2) * 84 - 1; for (Int_t iTimeBin = 0; iTimeBin < fNTimeBin; iTimeBin++) { - for (Int_t iAdc = 0; iAdc < fgkNADC; iAdc++) { + for (Int_t iAdc = 0; iAdc < AliTRDfeeParam::GetNadcMcm(); iAdc++) { Int_t value = adcArray->GetDataByAdcCol(GetRow(), offset - iAdc, iTimeBin); - if (value < 0 || (offset - iAdc < 1) || (offset - iAdc > 165)) { - fADCR[iAdc][iTimeBin] = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kFPNP) + (fgAddBaseline << fgkAddDigits); - fADCF[iAdc][iTimeBin] = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPFP) + (fgAddBaseline << fgkAddDigits); + // treat 0 as suppressed, + // this is not correct but reported like that from arrayADC + if (value <= 0 || (offset - iAdc < 1) || (offset - iAdc > 165)) { + fADCR[iAdc][iTimeBin] = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kFPNP, fDetector, fRobPos, fMcmPos) + (fgAddBaseline << fgkAddDigits); + fADCF[iAdc][iTimeBin] = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPFP, fDetector, fRobPos, fMcmPos) + (fgAddBaseline << fgkAddDigits); } else { fZSMap[iAdc] = 0; @@ -576,12 +594,12 @@ void AliTRDmcmSim::SetData(AliTRDarrayADC* const adcArray, AliTRDdigitsManager * } } -void AliTRDmcmSim::SetDataByPad(AliTRDarrayADC* const adcArray, AliTRDdigitsManager * const digitsManager) +void AliTRDmcmSim::SetDataByPad(const AliTRDarrayADC* const adcArray, AliTRDdigitsManager * const digitsManager) { - // Set the ADC data from an AliTRDarrayADC + // Set the ADC data from an AliTRDarrayADC // (by pad, to be used during initial reading in simulation) - if( !CheckInitialized() ) + if( !CheckInitialized() ) return; fDigitsManager = digitsManager; @@ -589,41 +607,42 @@ void AliTRDmcmSim::SetDataByPad(AliTRDarrayADC* const adcArray, AliTRDdigitsMana for (Int_t iDict = 0; iDict < 3; iDict++) { AliTRDarrayDictionary *newDict = (AliTRDarrayDictionary*) fDigitsManager->GetDictionary(fDetector, iDict); if (fDict[iDict] != 0x0 && newDict != 0x0) { - + if (fDict[iDict] == newDict) continue; fDict[iDict] = newDict; - - if (fDict[iDict]->GetDim() == 0) { - AliError(Form("Dictionary %i of det. %i has dim. 0", fDetector, iDict)); - continue; - } - fDict[iDict]->Expand(); + fDict[iDict]->Expand(); } else { fDict[iDict] = newDict; if (fDict[iDict]) fDict[iDict]->Expand(); } + + // If there is no data, set dictionary to zero to avoid crashes + if (fDict[iDict]->GetDim() == 0) { + AliError(Form("Dictionary %i of det. %i has dim. 0", iDict, fDetector)); + fDict[iDict] = 0x0; + } } } if (fNTimeBin != adcArray->GetNtime()) SetNTimebins(adcArray->GetNtime()); - + Int_t offset = (fMcmPos % 4 + 1) * 18 + (fRobPos % 2) * 72 + 1; for (Int_t iTimeBin = 0; iTimeBin < fNTimeBin; iTimeBin++) { - for (Int_t iAdc = 0; iAdc < fgkNADC; iAdc++) { + for (Int_t iAdc = 0; iAdc < AliTRDfeeParam::GetNadcMcm(); iAdc++) { Int_t value = -1; Int_t pad = offset - iAdc; - if (pad > -1 && pad < 144) + if (pad > -1 && pad < 144) value = adcArray->GetData(GetRow(), offset - iAdc, iTimeBin); // Int_t value = adcArray->GetDataByAdcCol(GetRow(), offset - iAdc, iTimeBin); if (value < 0 || (offset - iAdc < 1) || (offset - iAdc > 165)) { - fADCR[iAdc][iTimeBin] = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kFPNP) + (fgAddBaseline << fgkAddDigits); - fADCF[iAdc][iTimeBin] = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPFP) + (fgAddBaseline << fgkAddDigits); + fADCR[iAdc][iTimeBin] = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kFPNP, fDetector, fRobPos, fMcmPos) + (fgAddBaseline << fgkAddDigits); + fADCF[iAdc][iTimeBin] = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPFP, fDetector, fRobPos, fMcmPos) + (fgAddBaseline << fgkAddDigits); } else { fZSMap[iAdc] = 0; @@ -640,16 +659,16 @@ void AliTRDmcmSim::SetDataPedestal( Int_t adc ) // Store ADC data into array of raw data // - if( !CheckInitialized() ) + if( !CheckInitialized() ) return; - if( adc < 0 || adc >= fgkNADC ) { + if( adc < 0 || adc >= AliTRDfeeParam::GetNadcMcm() ) { return; } for( Int_t it = 0 ; it < fNTimeBin ; it++ ) { - fADCR[adc][it] = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kFPNP) + (fgAddBaseline << fgkAddDigits); - fADCF[adc][it] = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPFP) + (fgAddBaseline << fgkAddDigits); + fADCR[adc][it] = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kFPNP, fDetector, fRobPos, fMcmPos) + (fgAddBaseline << fgkAddDigits); + fADCF[adc][it] = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPFP, fDetector, fRobPos, fMcmPos) + (fgAddBaseline << fgkAddDigits); } } @@ -659,16 +678,16 @@ Bool_t AliTRDmcmSim::GetHit(Int_t index, Int_t &channel, Int_t &timebin, Int_t & if (index < 0 || index >= fNHits) return kFALSE; - + channel = fHits[index].fChannel; timebin = fHits[index].fTimebin; qtot = fHits[index].fQtot; ypos = fHits[index].fYpos; y = (Float_t) ((((((fRobPos & 0x1) << 2) + (fMcmPos & 0x3)) * 18) << 8) - ((18*4*2 - 18*2 - 1) << 7) - - (channel << 8) - ypos) + (channel << 8) - ypos) * (0.635 + 0.03 * (fDetector % 6)) / 256.0; - label = fHits[index].fLabel; + label = fHits[index].fLabel[0]; return kTRUE; } @@ -679,13 +698,13 @@ Int_t AliTRDmcmSim::GetCol( Int_t adc ) // Return column id of the pad for the given ADC channel // - if( !CheckInitialized() ) + if( !CheckInitialized() ) return -1; Int_t col = fFeeParam->GetPadColFromADC(fRobPos, fMcmPos, adc); - if (col < 0 || col >= fFeeParam->GetNcol()) + if (col < 0 || col >= fFeeParam->GetNcol()) return -1; - else + else return col; } @@ -693,56 +712,62 @@ Int_t AliTRDmcmSim::ProduceRawStream( UInt_t *buf, Int_t bufSize, UInt_t iEv) co { // // Produce raw data stream from this MCM and put in buf - // Returns number of words filled, or negative value + // Returns number of words filled, or negative value // with -1 * number of overflowed words // - if( !CheckInitialized() ) + if( !CheckInitialized() ) return 0; UInt_t x; + UInt_t mcmHeader = 0; + UInt_t adcMask = 0; Int_t nw = 0; // Number of written words Int_t of = 0; // Number of overflowed words Int_t rawVer = fFeeParam->GetRAWversion(); Int_t **adc; Int_t nActiveADC = 0; // number of activated ADC bits in a word - if( !CheckInitialized() ) + if( !CheckInitialized() ) return 0; - if (fTrapConfig->GetTrapReg(AliTRDtrapConfig::kEBSF) != 0) // store unfiltered data + if (fTrapConfig->GetTrapReg(AliTRDtrapConfig::kEBSF, fDetector, fRobPos, fMcmPos) != 0) // store unfiltered data adc = fADCR; - else + else adc = fADCF; - - // Produce MCM header - x = (1<<31) | (fRobPos << 28) | (fMcmPos << 24) | ((iEv % 0x100000) << 4) | 0xC; - - if (nw < bufSize) { - buf[nw++] = x; - } - else { - of++; - } // Produce ADC mask : nncc cccm mmmm mmmm mmmm mmmm mmmm 1100 // n : unused , c : ADC count, m : selected ADCs - if( rawVer >= 3 ) { - x = 0; - for( Int_t iAdc = 0 ; iAdc < fgkNADC ; iAdc++ ) { + if( rawVer >= 3 && + (fTrapConfig->GetTrapReg(AliTRDtrapConfig::kC15CPUA, fDetector, fRobPos, fMcmPos) & (1 << 13))) { // check for zs flag in TRAP configuration + for( Int_t iAdc = 0 ; iAdc < AliTRDfeeParam::GetNadcMcm() ; iAdc++ ) { if( ~fZSMap[iAdc] != 0 ) { // 0 means not suppressed - x = x | (1 << (iAdc+4) ); // last 4 digit reserved for 1100=0xc - nActiveADC++; // number of 1 in mmm....m + adcMask |= (1 << (iAdc+4) ); // last 4 digit reserved for 1100=0xc + nActiveADC++; // number of 1 in mmm....m } } - x = x | (1 << 30) | ( ( 0x3FFFFFFC ) & (~(nActiveADC) << 25) ) | 0xC; // nn = 01, ccccc are inverted, 0xc=1100 - if (nw < bufSize) { - buf[nw++] = x; - } - else { + if ((nActiveADC == 0) && + (fTrapConfig->GetTrapReg(AliTRDtrapConfig::kC15CPUA, fDetector, fRobPos, fMcmPos) & (1 << 8))) // check for DEH flag in TRAP configuration + return 0; + + // assemble adc mask word + adcMask |= (1 << 30) | ( ( 0x3FFFFFFC ) & (~(nActiveADC) << 25) ) | 0xC; // nn = 01, ccccc are inverted, 0xc=1100 + } + + // MCM header + mcmHeader = (1<<31) | (fRobPos << 28) | (fMcmPos << 24) | ((iEv % 0x100000) << 4) | 0xC; + if (nw < bufSize) + buf[nw++] = mcmHeader; + else + of++; + + // ADC mask + if( adcMask != 0 ) { + if (nw < bufSize) + buf[nw++] = adcMask; + else of++; - } } // Produce ADC data. 3 timebins are packed into one 32 bits word @@ -774,26 +799,26 @@ Int_t AliTRDmcmSim::ProduceTrackletStream( UInt_t *buf, Int_t bufSize ) { // // Produce tracklet data stream from this MCM and put in buf - // Returns number of words filled, or negative value + // Returns number of words filled, or negative value // with -1 * number of overflowed words // - if( !CheckInitialized() ) + if( !CheckInitialized() ) return 0; Int_t nw = 0; // Number of written words Int_t of = 0; // Number of overflowed words - - // Produce tracklet data. A maximum of four 32 Bit words will be written per MCM + + // Produce tracklet data. A maximum of four 32 Bit words will be written per MCM // fMCMT is filled continuously until no more tracklet words available for (Int_t iTracklet = 0; iTracklet < fTrackletArray->GetEntriesFast(); iTracklet++) { - if (nw < bufSize) + if (nw < bufSize) buf[nw++] = ((AliTRDtrackletMCM*) (*fTrackletArray)[iTracklet])->GetTrackletWord(); - else + else of++; } - + if( of != 0 ) return -of; else return nw; } @@ -802,20 +827,20 @@ void AliTRDmcmSim::Filter() // // Filter the raw ADC values. The active filter stages and their // parameters are taken from AliTRDtrapConfig. - // The raw data is stored separate from the filtered data. Thus, - // it is possible to run the filters on a set of raw values + // The raw data is stored separate from the filtered data. Thus, + // it is possible to run the filters on a set of raw values // sequentially for parameter tuning. // - if( !CheckInitialized() ) + if( !CheckInitialized() ) return; // Apply filters sequentially. Bypass is handled by filters - // since counters and internal registers may be updated even + // since counters and internal registers may be updated even // if the filter is bypassed. - // The first filter takes the data from fADCR and - // outputs to fADCF. - + // The first filter takes the data from fADCR and + // outputs to fADCF. + // Non-linearity filter not implemented. FilterPedestal(); FilterGain(); @@ -823,31 +848,31 @@ void AliTRDmcmSim::Filter() // Crosstalk filter not implemented. } -void AliTRDmcmSim::FilterPedestalInit(Int_t baseline) +void AliTRDmcmSim::FilterPedestalInit(Int_t baseline) { - // Initializes the pedestal filter assuming that the input has + // Initializes the pedestal filter assuming that the input has // been constant for a long time (compared to the time constant). - UShort_t fptc = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kFPTC); // 0..3, 0 - fastest, 3 - slowest + UShort_t fptc = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kFPTC, fDetector, fRobPos, fMcmPos); // 0..3, 0 - fastest, 3 - slowest - for (Int_t iAdc = 0; iAdc < fgkNADC; iAdc++) - fPedAcc[iAdc] = (baseline << 2) * (1 << fgkFPshifts[fptc]); + for (Int_t iAdc = 0; iAdc < AliTRDfeeParam::GetNadcMcm(); iAdc++) + fPedAcc[iAdc] = (baseline << 2) * (1 << fgkFPshifts[fptc]); } UShort_t AliTRDmcmSim::FilterPedestalNextSample(Int_t adc, Int_t timebin, UShort_t value) { // Returns the output of the pedestal filter given the input value. - // The output depends on the internal registers and, thus, the + // The output depends on the internal registers and, thus, the // history of the filter. - UShort_t fpnp = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kFPNP); // 0..511 -> 0..127.75, pedestal at the output - UShort_t fptc = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kFPTC); // 0..3, 0 - fastest, 3 - slowest - UShort_t fpby = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kFPBY); // 0..1 bypass, active low + UShort_t fpnp = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kFPNP, fDetector, fRobPos, fMcmPos); // 0..511 -> 0..127.75, pedestal at the output + UShort_t fptc = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kFPTC, fDetector, fRobPos, fMcmPos); // 0..3, 0 - fastest, 3 - slowest + UShort_t fpby = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kFPBY, fDetector, fRobPos, fMcmPos); // 0..1 bypass, active low UShort_t accumulatorShifted; Int_t correction; UShort_t inpAdd; - + inpAdd = value + fpnp; accumulatorShifted = (fPedAcc[adc] >> fgkFPshifts[fptc]) & 0x3FF; // 10 bits @@ -856,7 +881,7 @@ UShort_t AliTRDmcmSim::FilterPedestalNextSample(Int_t adc, Int_t timebin, UShort correction = (value & 0x3FF) - accumulatorShifted; fPedAcc[adc] = (fPedAcc[adc] + correction) & 0x7FFFFFFF; // 31 bits } - + if (fpby == 0) return value; @@ -865,9 +890,9 @@ UShort_t AliTRDmcmSim::FilterPedestalNextSample(Int_t adc, Int_t timebin, UShort else { inpAdd = inpAdd - accumulatorShifted; - if (inpAdd > 0xFFF) + if (inpAdd > 0xFFF) return 0xFFF; - else + else return inpAdd; } } @@ -877,14 +902,14 @@ void AliTRDmcmSim::FilterPedestal() // // Apply pedestal filter // - // As the first filter in the chain it reads data from fADCR - // and outputs to fADCF. - // It has only an effect if previous samples have been fed to - // find the pedestal. Currently, the simulation assumes that + // As the first filter in the chain it reads data from fADCR + // and outputs to fADCF. + // It has only an effect if previous samples have been fed to + // find the pedestal. Currently, the simulation assumes that // the input has been stable for a sufficiently long time. for (Int_t iTimeBin = 0; iTimeBin < fNTimeBin; iTimeBin++) { - for (Int_t iAdc = 0; iAdc < fgkNADC; iAdc++) { + for (Int_t iAdc = 0; iAdc < AliTRDfeeParam::GetNadcMcm(); iAdc++) { fADCF[iAdc][iTimeBin] = FilterPedestalNextSample(iAdc, iTimeBin, fADCR[iAdc][iTimeBin]); } } @@ -892,11 +917,11 @@ void AliTRDmcmSim::FilterPedestal() void AliTRDmcmSim::FilterGainInit() { - // Initializes the gain filter. In this case, only threshold + // Initializes the gain filter. In this case, only threshold // counters are reset. - for (Int_t iAdc = 0; iAdc < fgkNADC; iAdc++) { - // these are counters which in hardware continue + for (Int_t iAdc = 0; iAdc < AliTRDfeeParam::GetNadcMcm(); iAdc++) { + // these are counters which in hardware continue // until maximum or reset fGainCounterA[iAdc] = 0; fGainCounterB[iAdc] = 0; @@ -907,35 +932,39 @@ UShort_t AliTRDmcmSim::FilterGainNextSample(Int_t adc, UShort_t value) { // Apply the gain filter to the given value. // BEGIN_LATEX O_{i}(t) = #gamma_{i} * I_{i}(t) + a_{i} END_LATEX - // The output depends on the internal registers and, thus, the + // The output depends on the internal registers and, thus, the // history of the filter. - UShort_t fgby = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kFGBY); // bypass, active low - UShort_t fgf = fTrapConfig->GetTrapReg(AliTRDtrapConfig::TrapReg_t(AliTRDtrapConfig::kFGF0 + adc)); // 0x700 + (0 & 0x1ff); - UShort_t fga = fTrapConfig->GetTrapReg(AliTRDtrapConfig::TrapReg_t(AliTRDtrapConfig::kFGA0 + adc)); // 40; - UShort_t fgta = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kFGTA); // 20; - UShort_t fgtb = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kFGTB); // 2060; + UShort_t fgby = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kFGBY, fDetector, fRobPos, fMcmPos); // bypass, active low + UShort_t fgf = fTrapConfig->GetTrapReg(AliTRDtrapConfig::TrapReg_t(AliTRDtrapConfig::kFGF0 + adc), fDetector, fRobPos, fMcmPos); // 0x700 + (0 & 0x1ff); + UShort_t fga = fTrapConfig->GetTrapReg(AliTRDtrapConfig::TrapReg_t(AliTRDtrapConfig::kFGA0 + adc), fDetector, fRobPos, fMcmPos); // 40; + UShort_t fgta = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kFGTA, fDetector, fRobPos, fMcmPos); // 20; + UShort_t fgtb = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kFGTB, fDetector, fRobPos, fMcmPos); // 2060; + + UInt_t fgfExtended = 0x700 + fgf; // The corr factor which is finally applied has to be extended by 0x700 (hex) or 0.875 (dec) + // because fgf=0 correspons to 0.875 and fgf=511 correspons to 1.125 - 2^(-11) + // (see TRAP User Manual for details) UInt_t corr; // corrected value value &= 0xFFF; - corr = (value * fgf) >> 11; + corr = (value * fgfExtended) >> 11; corr = corr > 0xfff ? 0xfff : corr; corr = AddUintClipping(corr, fga, 12); - // Update threshold counters + // Update threshold counters // not really useful as they are cleared with every new event if (!((fGainCounterA[adc] == 0x3FFFFFF) || (fGainCounterB[adc] == 0x3FFFFFF))) // stop when full { - if (corr >= fgtb) + if (corr >= fgtb) fGainCounterB[adc]++; - else if (corr >= fgta) + else if (corr >= fgta) fGainCounterA[adc]++; } if (fgby == 1) - return corr; + return corr; else return value; } @@ -944,7 +973,7 @@ void AliTRDmcmSim::FilterGain() { // Read data from fADCF and apply gain filter. - for (Int_t iAdc = 0; iAdc < fgkNADC; iAdc++) { + for (Int_t iAdc = 0; iAdc < AliTRDfeeParam::GetNadcMcm(); iAdc++) { for (Int_t iTimeBin = 0; iTimeBin < fNTimeBin; iTimeBin++) { fADCF[iAdc][iTimeBin] = FilterGainNextSample(iAdc, fADCF[iAdc][iTimeBin]); } @@ -953,14 +982,14 @@ void AliTRDmcmSim::FilterGain() void AliTRDmcmSim::FilterTailInit(Int_t baseline) { - // Initializes the tail filter assuming that the input has - // been at the baseline value (configured by FTFP) for a + // Initializes the tail filter assuming that the input has + // been at the baseline value (configured by FTFP) for a // sufficiently long time. // exponents and weight calculated from configuration - UShort_t alphaLong = 0x3ff & fTrapConfig->GetTrapReg(AliTRDtrapConfig::kFTAL); // the weight of the long component - UShort_t lambdaLong = (1 << 10) | (1 << 9) | (fTrapConfig->GetTrapReg(AliTRDtrapConfig::kFTLL) & 0x1FF); // the multiplier - UShort_t lambdaShort = (0 << 10) | (1 << 9) | (fTrapConfig->GetTrapReg(AliTRDtrapConfig::kFTLS) & 0x1FF); // the multiplier + UShort_t alphaLong = 0x3ff & fTrapConfig->GetTrapReg(AliTRDtrapConfig::kFTAL, fDetector, fRobPos, fMcmPos); // the weight of the long component + UShort_t lambdaLong = (1 << 10) | (1 << 9) | (fTrapConfig->GetTrapReg(AliTRDtrapConfig::kFTLL, fDetector, fRobPos, fMcmPos) & 0x1FF); // the multiplier + UShort_t lambdaShort = (0 << 10) | (1 << 9) | (fTrapConfig->GetTrapReg(AliTRDtrapConfig::kFTLS, fDetector, fRobPos, fMcmPos) & 0x1FF); // the multiplier Float_t lambdaL = lambdaLong * 1.0 / (1 << 11); Float_t lambdaS = lambdaShort * 1.0 / (1 << 11); @@ -974,16 +1003,16 @@ void AliTRDmcmSim::FilterTailInit(Int_t baseline) UShort_t aout; if (baseline < 0) - baseline = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kFPNP); - + baseline = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kFPNP, fDetector, fRobPos, fMcmPos); + ql = lambdaL * (1 - lambdaS) * alphaL; qs = lambdaS * (1 - lambdaL) * (1 - alphaL); - for (Int_t iAdc = 0; iAdc < fgkNADC; iAdc++) { + for (Int_t iAdc = 0; iAdc < AliTRDfeeParam::GetNadcMcm(); iAdc++) { Int_t value = baseline & 0xFFF; - Int_t corr = (value * fTrapConfig->GetTrapReg(AliTRDtrapConfig::TrapReg_t(AliTRDtrapConfig::kFGF0 + iAdc))) >> 11; + Int_t corr = (value * fTrapConfig->GetTrapReg(AliTRDtrapConfig::TrapReg_t(AliTRDtrapConfig::kFGF0 + iAdc), fDetector, fRobPos, fMcmPos)) >> 11; corr = corr > 0xfff ? 0xfff : corr; - corr = AddUintClipping(corr, fTrapConfig->GetTrapReg(AliTRDtrapConfig::TrapReg_t(AliTRDtrapConfig::kFGA0 + iAdc)), 12); + corr = AddUintClipping(corr, fTrapConfig->GetTrapReg(AliTRDtrapConfig::TrapReg_t(AliTRDtrapConfig::kFGA0 + iAdc), fDetector, fRobPos, fMcmPos), 12); kt = kdc * baseline; aout = baseline - (UShort_t) kt; @@ -995,35 +1024,35 @@ void AliTRDmcmSim::FilterTailInit(Int_t baseline) UShort_t AliTRDmcmSim::FilterTailNextSample(Int_t adc, UShort_t value) { - // Returns the output of the tail filter for the given input value. - // The output depends on the internal registers and, thus, the + // Returns the output of the tail filter for the given input value. + // The output depends on the internal registers and, thus, the // history of the filter. // exponents and weight calculated from configuration - UShort_t alphaLong = 0x3ff & fTrapConfig->GetTrapReg(AliTRDtrapConfig::kFTAL); // the weight of the long component - UShort_t lambdaLong = (1 << 10) | (1 << 9) | (fTrapConfig->GetTrapReg(AliTRDtrapConfig::kFTLL) & 0x1FF); // the multiplier of the long component - UShort_t lambdaShort = (0 << 10) | (1 << 9) | (fTrapConfig->GetTrapReg(AliTRDtrapConfig::kFTLS) & 0x1FF); // the multiplier of the short component + UShort_t alphaLong = 0x3ff & fTrapConfig->GetTrapReg(AliTRDtrapConfig::kFTAL, fDetector, fRobPos, fMcmPos); // the weight of the long component + UShort_t lambdaLong = (1 << 10) | (1 << 9) | (fTrapConfig->GetTrapReg(AliTRDtrapConfig::kFTLL, fDetector, fRobPos, fMcmPos) & 0x1FF); // the multiplier of the long component + UShort_t lambdaShort = (0 << 10) | (1 << 9) | (fTrapConfig->GetTrapReg(AliTRDtrapConfig::kFTLS, fDetector, fRobPos, fMcmPos) & 0x1FF); // the multiplier of the short component // intermediate signals UInt_t aDiff; UInt_t alInpv; UShort_t aQ; UInt_t tmp; - + UShort_t inpVolt = value & 0xFFF; // 12 bits - + // add the present generator outputs aQ = AddUintClipping(fTailAmplLong[adc], fTailAmplShort[adc], 12); // calculate the difference between the input and the generated signal - if (inpVolt > aQ) + if (inpVolt > aQ) aDiff = inpVolt - aQ; - else + else aDiff = 0; - + // the inputs to the two generators, weighted alInpv = (aDiff * alphaLong) >> 11; - + // the new values of the registers, used next time // long component tmp = AddUintClipping(fTailAmplLong[adc], alInpv, 12); @@ -1033,9 +1062,9 @@ UShort_t AliTRDmcmSim::FilterTailNextSample(Int_t adc, UShort_t value) tmp = AddUintClipping(fTailAmplShort[adc], aDiff - alInpv, 12); tmp = (tmp * lambdaShort) >> 11; fTailAmplShort[adc] = tmp & 0xFFF; - + // the output of the filter - if (fTrapConfig->GetTrapReg(AliTRDtrapConfig::kFTBY) == 0) // bypass mode, active low + if (fTrapConfig->GetTrapReg(AliTRDtrapConfig::kFTBY, fDetector, fRobPos, fMcmPos) == 0) // bypass mode, active low return value; else return aDiff; @@ -1043,10 +1072,10 @@ UShort_t AliTRDmcmSim::FilterTailNextSample(Int_t adc, UShort_t value) void AliTRDmcmSim::FilterTail() { - // Apply tail cancellation filter to all data. + // Apply tail cancellation filter to all data. for (Int_t iTimeBin = 0; iTimeBin < fNTimeBin; iTimeBin++) { - for (Int_t iAdc = 0; iAdc < fgkNADC; iAdc++) { + for (Int_t iAdc = 0; iAdc < AliTRDfeeParam::GetNadcMcm(); iAdc++) { fADCF[iAdc][iTimeBin] = FilterTailNextSample(iAdc, fADCF[iAdc][iTimeBin]); } } @@ -1062,17 +1091,17 @@ void AliTRDmcmSim::ZSMapping() // http://www.kip.uni-heidelberg.de/ti/TRD/doc/trap/TRAP-UserManual.pdf // - if( !CheckInitialized() ) + if( !CheckInitialized() ) return; - Int_t eBIS = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kEBIS); - Int_t eBIT = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kEBIT); - Int_t eBIL = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kEBIL); - Int_t eBIN = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kEBIN); + Int_t eBIS = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kEBIS, fDetector, fRobPos, fMcmPos); + Int_t eBIT = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kEBIT, fDetector, fRobPos, fMcmPos); + Int_t eBIL = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kEBIL, fDetector, fRobPos, fMcmPos); + Int_t eBIN = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kEBIN, fDetector, fRobPos, fMcmPos); Int_t **adc = fADCF; - for (Int_t iAdc = 0; iAdc < fgkNADC; iAdc++) + for (Int_t iAdc = 0; iAdc < AliTRDfeeParam::GetNadcMcm(); iAdc++) fZSMap[iAdc] = -1; for( Int_t it = 0 ; it < fNTimeBin ; it++ ) { @@ -1082,55 +1111,55 @@ void AliTRDmcmSim::ZSMapping() Int_t an; Int_t mask; Int_t supp; // suppression of the current channel (low active) - + // ----- first channel ----- iAdc = 0; - - ap = 0; // previous - ac = adc[iAdc ][it]; // current - an = adc[iAdc+1][it]; // next - + + ap = 0 >> fgkAddDigits; // previous + ac = adc[iAdc ][it] >> fgkAddDigits; // current + an = adc[iAdc+1][it] >> fgkAddDigits; // next + mask = ( ac >= ap && ac >= an ) ? 0 : 0x1; // peak center detection mask += ( ap + ac + an > eBIT ) ? 0 : 0x2; // cluster mask += ( ac > eBIS ) ? 0 : 0x4; // absolute large peak - + supp = (eBIL >> mask) & 1; - + fZSMap[iAdc] &= ~((1-supp) << it); if( eBIN == 0 ) { // neighbour sensitivity fZSMap[iAdc+1] &= ~((1-supp) << it); } - + // ----- last channel ----- - iAdc = fgkNADC - 1; - - ap = adc[iAdc-1][it]; // previous - ac = adc[iAdc ][it]; // current - an = 0; // next - + iAdc = AliTRDfeeParam::GetNadcMcm() - 1; + + ap = adc[iAdc-1][it] >> fgkAddDigits; // previous + ac = adc[iAdc ][it] >> fgkAddDigits; // current + an = 0 >> fgkAddDigits; // next + mask = ( ac >= ap && ac >= an ) ? 0 : 0x1; // peak center detection mask += ( ap + ac + an > eBIT ) ? 0 : 0x2; // cluster mask += ( ac > eBIS ) ? 0 : 0x4; // absolute large peak - + supp = (eBIL >> mask) & 1; - + fZSMap[iAdc] &= ~((1-supp) << it); if( eBIN == 0 ) { // neighbour sensitivity fZSMap[iAdc-1] &= ~((1-supp) << it); } - + // ----- middle channels ----- - for( iAdc = 1 ; iAdc < fgkNADC-1; iAdc++ ) { - ap = adc[iAdc-1][it]; // previous - ac = adc[iAdc ][it]; // current - an = adc[iAdc+1][it]; // next - + for( iAdc = 1 ; iAdc < AliTRDfeeParam::GetNadcMcm()-1; iAdc++ ) { + ap = adc[iAdc-1][it] >> fgkAddDigits; // previous + ac = adc[iAdc ][it] >> fgkAddDigits; // current + an = adc[iAdc+1][it] >> fgkAddDigits; // next + mask = ( ac >= ap && ac >= an ) ? 0 : 0x1; // peak center detection mask += ( ap + ac + an > eBIT ) ? 0 : 0x2; // cluster mask += ( ac > eBIS ) ? 0 : 0x4; // absolute large peak - + supp = (eBIL >> mask) & 1; - + fZSMap[iAdc] &= ~((1-supp) << it); if( eBIN == 0 ) { // neighbour sensitivity fZSMap[iAdc-1] &= ~((1-supp) << it); @@ -1141,23 +1170,23 @@ void AliTRDmcmSim::ZSMapping() } } -void AliTRDmcmSim::AddHitToFitreg(Int_t adc, UShort_t timebin, UShort_t qtot, Short_t ypos, Int_t label) +void AliTRDmcmSim::AddHitToFitreg(Int_t adc, UShort_t timebin, UShort_t qtot, Short_t ypos, Int_t label[]) { - // Add the given hit to the fit register which is lateron used for - // the tracklet calculation. - // In addition to the fit sums in the fit register MC information + // Add the given hit to the fit register which is lateron used for + // the tracklet calculation. + // In addition to the fit sums in the fit register MC information // is stored. - if ((timebin >= fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPQS0)) && - (timebin < fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPQE0))) + if ((timebin >= fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPQS0, fDetector, fRobPos, fMcmPos)) && + (timebin < fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPQE0, fDetector, fRobPos, fMcmPos))) fFitReg[adc].fQ0 += qtot; - - if ((timebin >= fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPQS1)) && - (timebin < fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPQE1))) + + if ((timebin >= fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPQS1, fDetector, fRobPos, fMcmPos)) && + (timebin < fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPQE1, fDetector, fRobPos, fMcmPos))) fFitReg[adc].fQ1 += qtot; - - if ((timebin >= fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPFS) ) && - (timebin < fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPFE))) + + if ((timebin >= fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPFS, fDetector, fRobPos, fMcmPos) ) && + (timebin < fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPFE, fDetector, fRobPos, fMcmPos))) { fFitReg[adc].fSumX += timebin; fFitReg[adc].fSumX2 += timebin*timebin; @@ -1165,6 +1194,9 @@ void AliTRDmcmSim::AddHitToFitreg(Int_t adc, UShort_t timebin, UShort_t qtot, Sh fFitReg[adc].fSumY += ypos; fFitReg[adc].fSumY2 += ypos*ypos; fFitReg[adc].fSumXY += timebin*ypos; + AliDebug(10, Form("fitreg[%2i] in timebin %2i: X=%i, X2=%i, N=%i, Y=%i, Y2=%i, XY=%i, Q0=%i, Q1=%i", + adc, timebin, fFitReg[adc].fSumX, fFitReg[adc].fSumX2, fFitReg[adc].fNhits, + fFitReg[adc].fSumY, fFitReg[adc].fSumY2, fFitReg[adc].fSumXY, fFitReg[adc].fQ0, fFitReg[adc].fQ1)); } // register hits (MC info) @@ -1172,37 +1204,56 @@ void AliTRDmcmSim::AddHitToFitreg(Int_t adc, UShort_t timebin, UShort_t qtot, Sh fHits[fNHits].fQtot = qtot; fHits[fNHits].fYpos = ypos; fHits[fNHits].fTimebin = timebin; - fHits[fNHits].fLabel = label; + fHits[fNHits].fLabel[0] = label[0]; + fHits[fNHits].fLabel[1] = label[1]; + fHits[fNHits].fLabel[2] = label[2]; fNHits++; } -void AliTRDmcmSim::CalcFitreg() +void AliTRDmcmSim::CalcFitreg() { // Preprocessing. // Detect the hits and fill the fit registers. - // Requires 12-bit data from fADCF which means Filter() + // Requires 12-bit data from fADCF which means Filter() // has to be called before even if all filters are bypassed. //??? to be clarified: UInt_t adcMask = 0xffffffff; - - UShort_t timebin, adcch, adcLeft, adcCentral, adcRight, hitQual, timebin1, timebin2, qtotTemp; + + Bool_t hitQual; + Int_t adcLeft, adcCentral, adcRight; + UShort_t timebin, adcch, timebin1, timebin2, qtotTemp; Short_t ypos, fromLeft, fromRight, found; - UShort_t qTotal[19]; // the last is dummy + UShort_t qTotal[19+1]; // the last is dummy UShort_t marked[6], qMarked[6], worse1, worse2; - - timebin1 = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPFS); - if (fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPQS0) - < timebin1) - timebin1 = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPQS0); - timebin2 = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPFE); - if (fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPQE1) - > timebin2) - timebin2 = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPQE1); + + if (fgStoreClusters) { + timebin1 = 0; + timebin2 = fNTimeBin; + } + else { + // find first timebin to be looked at + timebin1 = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPFS, fDetector, fRobPos, fMcmPos); + if (fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPQS0, fDetector, fRobPos, fMcmPos) + < timebin1) + timebin1 = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPQS0, fDetector, fRobPos, fMcmPos); + if (fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPQS1, fDetector, fRobPos, fMcmPos) + < timebin1) + timebin1 = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPQS1, fDetector, fRobPos, fMcmPos); + + // find last timebin to be looked at + timebin2 = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPFE, fDetector, fRobPos, fMcmPos); + if (fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPQE0, fDetector, fRobPos, fMcmPos) + > timebin2) + timebin2 = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPQE0, fDetector, fRobPos, fMcmPos); + if (fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPQE1, fDetector, fRobPos, fMcmPos) + > timebin2) + timebin2 = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPQE1, fDetector, fRobPos, fMcmPos); + } // reset the fit registers - fNHits = 0; - for (adcch = 0; adcch < fgkNADC-2; adcch++) // due to border channels + fNHits = 0; + for (adcch = 0; adcch < AliTRDfeeParam::GetNadcMcm()-2; adcch++) // due to border channels { fFitReg[adcch].fNhits = 0; fFitReg[adcch].fQ0 = 0; @@ -1213,26 +1264,44 @@ void AliTRDmcmSim::CalcFitreg() fFitReg[adcch].fSumY2 = 0; fFitReg[adcch].fSumXY = 0; } - + for (timebin = timebin1; timebin < timebin2; timebin++) { // first find the hit candidates and store the total cluster charge in qTotal array // in case of not hit store 0 there. - for (adcch = 0; adcch < fgkNADC-2; adcch++) { + for (adcch = 0; adcch < AliTRDfeeParam::GetNadcMcm()-2; adcch++) { if ( ( (adcMask >> adcch) & 7) == 7) //??? all 3 channels are present in case of ZS { adcLeft = fADCF[adcch ][timebin]; adcCentral = fADCF[adcch+1][timebin]; adcRight = fADCF[adcch+2][timebin]; - if (fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPVBY) == 1) - hitQual = ( (adcLeft * adcRight) < - (fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPVT) * adcCentral) ); - else - hitQual = 1; + + if (fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPVBY, fDetector, fRobPos, fMcmPos) == 0) { + // bypass the cluster verification + hitQual = kTRUE; + } + else { + hitQual = ( (adcLeft * adcRight) < + ((fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPVT, fDetector, fRobPos, fMcmPos) * adcCentral*adcCentral) >> 10) ); + if (hitQual) + AliDebug(5, Form("cluster quality cut passed with %3i, %3i, %3i - threshold %3i -> %i", + adcLeft, adcCentral, adcRight, + fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPVT, fDetector, fRobPos, fMcmPos), + fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPVT, fDetector, fRobPos, fMcmPos) * adcCentral*adcCentral)); + } + // The accumulated charge is with the pedestal!!! qtotTemp = adcLeft + adcCentral + adcRight; + if ((fDebugStream) && (qtotTemp > 130)) { + (*fDebugStream) << "testtree" + << "qtot=" << qtotTemp + << "qleft=" << adcLeft + << "qcent=" << adcCentral + << "qright=" << adcRight + << "\n"; + } if ( (hitQual) && - (qtotTemp >= fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPHT)) && + (qtotTemp >= fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPHT, fDetector, fRobPos, fMcmPos)) && (adcLeft <= adcCentral) && (adcCentral > adcRight) ) qTotal[adcch] = qtotTemp; @@ -1241,7 +1310,7 @@ void AliTRDmcmSim::CalcFitreg() } else qTotal[adcch] = 0; //jkl - if (qTotal[adcch] != 0) + if (qTotal[adcch] != 0) AliDebug(10,Form("ch %2d qTotal %5d",adcch, qTotal[adcch])); } @@ -1261,7 +1330,7 @@ void AliTRDmcmSim::CalcFitreg() } adcch++; } - + fromRight = -1; adcch = 18; found = 0; @@ -1281,7 +1350,7 @@ void AliTRDmcmSim::CalcFitreg() if ((fromLeft >= 0) && (fromRight >= 0) && (fromLeft < fromRight)) for (adcch = fromLeft+1; adcch < fromRight; adcch++) qTotal[adcch] = 0; - + found = 0; for (adcch = 0; adcch < 19; adcch++) if (qTotal[adcch] > 0) found++; @@ -1295,7 +1364,7 @@ void AliTRDmcmSim::CalcFitreg() qMarked[found] = qTotal[marked[found]] >> 4; AliDebug(10,Form("ch_%d qTotal %d qTotals %d",marked[found],qTotal[marked[found]],qMarked[found])); } - + Sort6To2Worst(marked[0], marked[3], marked[4], marked[1], marked[2], marked[5], qMarked[0], qMarked[3], @@ -1316,7 +1385,7 @@ void AliTRDmcmSim::CalcFitreg() AliDebug(10,Form("Kill ch %d\n",worse2)); } } - + for (adcch = 0; adcch < 19; adcch++) { if (qTotal[adcch] > 0) // the channel is marked for processing { @@ -1325,11 +1394,11 @@ void AliTRDmcmSim::CalcFitreg() adcRight = fADCF[adcch+2][timebin]; // hit detected, in TRAP we have 4 units and a hit-selection, here we proceed all channels! // subtract the pedestal TPFP, clipping instead of wrapping - - Int_t regTPFP = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPFP); + + Int_t regTPFP = fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPFP, fDetector, fRobPos, fMcmPos); AliDebug(10, Form("Hit found, time=%d, adcch=%d/%d/%d, adc values=%d/%d/%d, regTPFP=%d, TPHT=%d\n", - timebin, adcch, adcch+1, adcch+2, adcLeft, adcCentral, adcRight, regTPFP, - fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPHT))); + timebin, adcch, adcch+1, adcch+2, adcLeft, adcCentral, adcRight, regTPFP, + fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPHT, fDetector, fRobPos, fMcmPos))); if (adcLeft < regTPFP) adcLeft = 0; else adcLeft -= regTPFP; if (adcCentral < regTPFP) adcCentral = 0; else adcCentral -= regTPFP; @@ -1339,21 +1408,21 @@ void AliTRDmcmSim::CalcFitreg() // checking for adcCentral != 0 (in case of "bad" configuration) if (adcCentral == 0) continue; - ypos = 128*(adcLeft - adcRight) / adcCentral; + ypos = 128*(adcRight - adcLeft) / adcCentral; if (ypos < 0) ypos = -ypos; // make the correction using the position LUT - ypos = ypos + fTrapConfig->GetTrapReg((AliTRDtrapConfig::TrapReg_t) (AliTRDtrapConfig::kTPL00 + (ypos & 0x7F))); + ypos = ypos + fTrapConfig->GetTrapReg((AliTRDtrapConfig::TrapReg_t) (AliTRDtrapConfig::kTPL00 + (ypos & 0x7F)), + fDetector, fRobPos, fMcmPos); if (adcLeft > adcRight) ypos = -ypos; - // label calculation - Int_t mcLabel = -1; + // label calculation (up to 3) + Int_t mcLabel[] = {-1, -1, -1}; if (fDigitsManager) { - Int_t label[9] = { 0 }; // up to 9 different labels possible - Int_t count[9] = { 0 }; - Int_t maxIdx = -1; - Int_t maxCount = 0; + const Int_t maxLabels = 9; + Int_t label[maxLabels] = { 0 }; // up to 9 different labels possible + Int_t count[maxLabels] = { 0 }; Int_t nLabels = 0; - Int_t padcol[3]; + Int_t padcol[3]; padcol[0] = fFeeParam->GetPadColFromADC(fRobPos, fMcmPos, adcch); padcol[1] = fFeeParam->GetPadColFromADC(fRobPos, fMcmPos, adcch+1); padcol[2] = fFeeParam->GetPadColFromADC(fRobPos, fMcmPos, adcch+2); @@ -1362,37 +1431,40 @@ void AliTRDmcmSim::CalcFitreg() if (!fDict[iDict]) continue; for (Int_t iPad = 0; iPad < 3; iPad++) { - if (padcol[iPad] < 0) + if (padcol[iPad] < 0) continue; - Int_t currLabel = fDict[iDict]->GetData(padrow, padcol[iPad], timebin); //fDigitsManager->GetTrack(iDict, padrow, padcol, timebin, fDetector); + Int_t currLabel = fDict[iDict]->GetData(padrow, padcol[iPad], timebin); AliDebug(10, Form("Read label: %4i for det: %3i, row: %i, col: %i, tb: %i\n", currLabel, fDetector, padrow, padcol[iPad], timebin)); for (Int_t iLabel = 0; iLabel < nLabels; iLabel++) { if (currLabel == label[iLabel]) { count[iLabel]++; - if (count[iLabel] > maxCount) { - maxCount = count[iLabel]; - maxIdx = iLabel; - } currLabel = -1; break; } - } + } if (currLabel >= 0) { - label[nLabels++] = currLabel; + label[nLabels] = currLabel; + count[nLabels] = 1; + nLabels++; } } } - if (maxIdx >= 0) - mcLabel = label[maxIdx]; + Int_t index[2*maxLabels]; + TMath::Sort(maxLabels, count, index); + for (Int_t i = 0; i < 3; i++) { + if (count[index[i]] <= 0) + break; + mcLabel[i] = label[index[i]]; + } } // add the hit to the fitregister - AddHitToFitreg(adcch, timebin, qTotal[adcch], ypos, mcLabel); + AddHitToFitreg(adcch, timebin, qTotal[adcch] >> fgkAddDigits, ypos, mcLabel); } } } - for (Int_t iAdc = 0; iAdc < fgkNADC; iAdc++) { + for (Int_t iAdc = 0; iAdc < AliTRDfeeParam::GetNadcMcm(); iAdc++) { if (fFitReg[iAdc].fNhits != 0) { AliDebug(2, Form("fitreg[%i]: nHits = %i, sumX = %i, sumY = %i, sumX2 = %i, sumY2 = %i, sumXY = %i", iAdc, fFitReg[iAdc].fNhits, @@ -1406,9 +1478,9 @@ void AliTRDmcmSim::CalcFitreg() } } -void AliTRDmcmSim::TrackletSelection() +void AliTRDmcmSim::TrackletSelection() { - // Select up to 4 tracklet candidates from the fit registers + // Select up to 4 tracklet candidates from the fit registers // and assign them to the CPUs. UShort_t adcIdx, i, j, ntracks, tmp; @@ -1416,10 +1488,10 @@ void AliTRDmcmSim::TrackletSelection() ntracks = 0; for (adcIdx = 0; adcIdx < 18; adcIdx++) // ADCs - if ( (fFitReg[adcIdx].fNhits - >= fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPCL)) && + if ( (fFitReg[adcIdx].fNhits + >= fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPCL, fDetector, fRobPos, fMcmPos)) && (fFitReg[adcIdx].fNhits+fFitReg[adcIdx+1].fNhits - >= fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPCT))) + >= fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPCT, fDetector, fRobPos, fMcmPos))) { trackletCand[ntracks][0] = adcIdx; trackletCand[ntracks][1] = fFitReg[adcIdx].fNhits+fFitReg[adcIdx+1].fNhits; @@ -1427,7 +1499,7 @@ void AliTRDmcmSim::TrackletSelection() ntracks++; }; - for (i=0; i 4) @@ -1453,7 +1525,7 @@ void AliTRDmcmSim::TrackletSelection() ntracks = 4; // cut the rest, 4 is the max } // else is not necessary to sort - + // now sort, so that the first tracklet going to CPU0 corresponds to the highest adc channel - as in the TRAP for (j = 0; j < (ntracks-1); j++) { @@ -1482,30 +1554,31 @@ void AliTRDmcmSim::TrackletSelection() void AliTRDmcmSim::FitTracklet() { - // Perform the actual tracklet fit based on the fit sums - // which have been filled in the fit registers. + // Perform the actual tracklet fit based on the fit sums + // which have been filled in the fit registers. // parameters in fitred.asm (fit program) - Int_t decPlaces = 5; Int_t rndAdd = 0; - if (decPlaces > 1) + Int_t decPlaces = 5; // must be larger than 1 or change the following code + // if (decPlaces > 1) rndAdd = (1 << (decPlaces-1)) + 1; - else if (decPlaces == 1) - rndAdd = 1; + // else if (decPlaces == 1) + // rndAdd = 1; + Int_t ndriftDp = 5; // decimal places for drift time Long64_t shift = ((Long64_t) 1 << 32); // calculated in fitred.asm Int_t padrow = ((fRobPos >> 1) << 2) | (fMcmPos >> 2); - Int_t yoffs = (((((fRobPos & 0x1) << 2) + (fMcmPos & 0x3)) * 18) << 8) - + Int_t yoffs = (((((fRobPos & 0x1) << 2) + (fMcmPos & 0x3)) * 18) << 8) - ((18*4*2 - 18*2 - 1) << 7); yoffs = yoffs << decPlaces; // holds position of ADC channel 1 Int_t layer = fDetector % 6; UInt_t scaleY = (UInt_t) ((0.635 + 0.03 * layer)/(256.0 * 160.0e-4) * shift); UInt_t scaleD = (UInt_t) ((0.635 + 0.03 * layer)/(256.0 * 140.0e-4) * shift); - Int_t deflCorr = (Int_t) fTrapConfig->GetDmemUnsigned(AliTRDtrapConfig::fgkDmemAddrDeflCorr, fDetector, fRobPos, fMcmPos); - Int_t ndrift = (Int_t) fTrapConfig->GetDmemUnsigned(AliTRDtrapConfig::fgkDmemAddrNdrift, fDetector, fRobPos, fMcmPos); + Int_t deflCorr = (Int_t) fTrapConfig->GetDmemUnsigned(fgkDmemAddrDeflCorr, fDetector, fRobPos, fMcmPos); + Int_t ndrift = (Int_t) fTrapConfig->GetDmemUnsigned(fgkDmemAddrNdrift, fDetector, fRobPos, fMcmPos); // local variables for calculation Long64_t mult, temp, denom; //??? @@ -1516,7 +1589,7 @@ void AliTRDmcmSim::FitTracklet() Int_t sumY2; // not used in the current TRAP program, now used for error calculation (simulation only) Float_t fitError, fitSlope, fitOffset; FitReg_t *fit0, *fit1; // pointers to relevant fit registers - + // const uint32_t OneDivN[32] = { // 2**31/N : exactly like in the TRAP, the simple division here gives the same result! // 0x00000000, 0x80000000, 0x40000000, 0x2AAAAAA0, 0x20000000, 0x19999990, 0x15555550, 0x12492490, // 0x10000000, 0x0E38E380, 0x0CCCCCC0, 0x0BA2E8B0, 0x0AAAAAA0, 0x09D89D80, 0x09249240, 0x08888880, @@ -1526,7 +1599,7 @@ void AliTRDmcmSim::FitTracklet() for (Int_t cpu = 0; cpu < 4; cpu++) { if (fFitPtr[cpu] == 31) { - fMCMT[cpu] = 0x10001000; //??? AliTRDfeeParam::GetTrackletEndmarker(); + fMCMT[cpu] = 0x10001000; //??? AliTRDfeeParam::GetTrackletEndmarker(); } else { @@ -1541,7 +1614,7 @@ void AliTRDmcmSim::FitTracklet() nHits = fit0->fNhits + fit1->fNhits; // number of hits sumX = fit0->fSumX + fit1->fSumX; sumX2 = fit0->fSumX2 + fit1->fSumX2; - denom = nHits*sumX2 - sumX*sumX; + denom = ((Long64_t) nHits)*((Long64_t) sumX2) - ((Long64_t) sumX)*((Long64_t) sumX); mult = mult / denom; // exactly like in the TRAP program q0 = fit0->fQ0 + fit1->fQ0; @@ -1559,29 +1632,29 @@ void AliTRDmcmSim::FitTracklet() offset = temp >> 32; // take the upper 32 bits offset = offset + yoffs; - AliDebug(10, Form("slope = %i, slope * ndrift = %i, deflCorr: %i", + AliDebug(10, Form("slope = %i, slope * ndrift = %i, deflCorr: %i", slope, slope * ndrift, deflCorr)); slope = ((slope * ndrift) >> ndriftDp) + deflCorr; offset = offset - (fFitPtr[cpu] << (8 + decPlaces)); - + temp = slope; temp = temp * scaleD; slope = (temp >> 32); temp = offset; temp = temp * scaleY; offset = (temp >> 32); - + // rounding, like in the TRAP slope = (slope + rndAdd) >> decPlaces; offset = (offset + rndAdd) >> decPlaces; - AliDebug(5, Form("Det: %3i, ROB: %i, MCM: %2i: deflection: %i, min: %i, max: %i", - fDetector, fRobPos, fMcmPos, slope, - (Int_t) fTrapConfig->GetDmemUnsigned(AliTRDtrapConfig::fgkDmemAddrDeflCutStart + 2*fFitPtr[cpu], fDetector, fRobPos, fMcmPos), - (Int_t) fTrapConfig->GetDmemUnsigned(AliTRDtrapConfig::fgkDmemAddrDeflCutStart + 1 + 2*fFitPtr[cpu], fDetector, fRobPos, fMcmPos))); + AliDebug(5, Form("Det: %3i, ROB: %i, MCM: %2i: deflection: %i, min: %i, max: %i", + fDetector, fRobPos, fMcmPos, slope, + (Int_t) fTrapConfig->GetDmemUnsigned(fgkDmemAddrDeflCutStart + 2*fFitPtr[cpu], fDetector, fRobPos, fMcmPos), + (Int_t) fTrapConfig->GetDmemUnsigned(fgkDmemAddrDeflCutStart + 1 + 2*fFitPtr[cpu], fDetector, fRobPos, fMcmPos))); - AliDebug(5, Form("Fit sums: x = %i, X = %i, y = %i, Y = %i, Z = %i", - sumX, sumX2, sumY, sumY2, sumXY)); + AliDebug(5, Form("Fit sums: x = %i, X = %i, y = %i, Y = %i, Z = %i, q0 = %i, q1 = %i", + sumX, sumX2, sumY, sumY2, sumXY, q0, q1)); fitSlope = (Float_t) (nHits * sumXY - sumX * sumY) / (nHits * sumX2 - sumX*sumX); @@ -1597,8 +1670,8 @@ void AliTRDmcmSim::FitTracklet() Bool_t rejected = kFALSE; // deflection range table from DMEM - if ((slope < ((Int_t) fTrapConfig->GetDmemUnsigned(AliTRDtrapConfig::fgkDmemAddrDeflCutStart + 2*fFitPtr[cpu], fDetector, fRobPos, fMcmPos))) || - (slope > ((Int_t) fTrapConfig->GetDmemUnsigned(AliTRDtrapConfig::fgkDmemAddrDeflCutStart + 1 + 2*fFitPtr[cpu], fDetector, fRobPos, fMcmPos)))) + if ((slope < ((Int_t) fTrapConfig->GetDmemUnsigned(fgkDmemAddrDeflCutStart + 2*fFitPtr[cpu], fDetector, fRobPos, fMcmPos))) || + (slope > ((Int_t) fTrapConfig->GetDmemUnsigned(fgkDmemAddrDeflCutStart + 1 + 2*fFitPtr[cpu], fDetector, fRobPos, fMcmPos)))) rejected = kTRUE; if (rejected && GetApplyCut()) @@ -1608,72 +1681,82 @@ void AliTRDmcmSim::FitTracklet() else { if (slope > 63 || slope < -64) { // wrapping in TRAP! - AliError(Form("Overflow in slope: %i, tracklet discarded!", slope)); + AliDebug(1,Form("Overflow in slope: %i, tracklet discarded!", slope)); fMCMT[cpu] = 0x10001000; continue; } slope = slope & 0x7F; // 7 bit - - if (offset > 0xfff || offset < -0xfff) + + if (offset > 0xfff || offset < -0xfff) AliWarning("Overflow in offset"); offset = offset & 0x1FFF; // 13 bit - pid = GetPID(q0 >> fgkAddDigits, q1 >> fgkAddDigits); // divided by 4 because in simulation there are two additional decimal places + pid = GetPID(q0, q1); if (pid > 0xff) AliWarning("Overflow in PID"); pid = pid & 0xFF; // 8 bit, exactly like in the TRAP program - + // assemble and store the tracklet word fMCMT[cpu] = (pid << 24) | (padrow << 20) | (slope << 13) | offset; - // calculate MC label - Int_t mcLabel = -1; + // calculate number of hits and MC label + Int_t mcLabel[] = { -1, -1, -1}; Int_t nHits0 = 0; Int_t nHits1 = 0; - if (fDigitsManager) { - Int_t label[30] = {0}; // up to 30 different labels possible - Int_t count[30] = {0}; - Int_t maxIdx = -1; - Int_t maxCount = 0; - Int_t nLabels = 0; - for (Int_t iHit = 0; iHit < fNHits; iHit++) { - if ((fHits[iHit].fChannel - fFitPtr[cpu] < 0) || - (fHits[iHit].fChannel - fFitPtr[cpu] > 1)) - continue; - // counting contributing hits - if (fHits[iHit].fTimebin >= fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPQS0) && - fHits[iHit].fTimebin < fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPQE0)) - nHits0++; - if (fHits[iHit].fTimebin >= fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPQS1) && - fHits[iHit].fTimebin < fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPQE1)) - nHits1++; - - Int_t currLabel = fHits[iHit].fLabel; - for (Int_t iLabel = 0; iLabel < nLabels; iLabel++) { - if (currLabel == label[iLabel]) { - count[iLabel]++; - if (count[iLabel] > maxCount) { - maxCount = count[iLabel]; - maxIdx = iLabel; - } - currLabel = -1; - break; - } - } - if (currLabel >= 0) { - label[nLabels++] = currLabel; - } - } - if (maxIdx >= 0) - mcLabel = label[maxIdx]; + const Int_t maxLabels = 30; + Int_t label[maxLabels] = {0}; // up to 30 different labels possible + Int_t count[maxLabels] = {0}; + Int_t nLabels = 0; + + for (Int_t iHit = 0; iHit < fNHits; iHit++) { + if ((fHits[iHit].fChannel - fFitPtr[cpu] < 0) || + (fHits[iHit].fChannel - fFitPtr[cpu] > 1)) + continue; + + // counting contributing hits + if (fHits[iHit].fTimebin >= fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPQS0, fDetector, fRobPos, fMcmPos) && + fHits[iHit].fTimebin < fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPQE0, fDetector, fRobPos, fMcmPos)) + nHits0++; + if (fHits[iHit].fTimebin >= fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPQS1, fDetector, fRobPos, fMcmPos) && + fHits[iHit].fTimebin < fTrapConfig->GetTrapReg(AliTRDtrapConfig::kTPQE1, fDetector, fRobPos, fMcmPos)) + nHits1++; + + // label calculation only if there is a digitsmanager to get the labels from + if (fDigitsManager) { + for (Int_t i = 0; i < 3; i++) { + Int_t currLabel = fHits[iHit].fLabel[i]; + for (Int_t iLabel = 0; iLabel < nLabels; iLabel++) { + if (currLabel == label[iLabel]) { + count[iLabel]++; + currLabel = -1; + break; + } + } + if (currLabel >= 0 && nLabels < maxLabels) { + label[nLabels] = currLabel; + count[nLabels]++; + nLabels++; + } + } + } + + if (fDigitsManager) { + Int_t index[2*maxLabels]; + TMath::Sort(maxLabels, count, index); + for (Int_t i = 0; i < 3; i++) { + if (count[index[i]] <= 0) + break; + mcLabel[i] = label[index[i]]; + } + } } new ((*fTrackletArray)[fTrackletArray->GetEntriesFast()]) AliTRDtrackletMCM((UInt_t) fMCMT[cpu], fDetector*2 + fRobPos%2, fRobPos, fMcmPos); ((AliTRDtrackletMCM*) (*fTrackletArray)[fTrackletArray->GetEntriesFast()-1])->SetLabel(mcLabel); - + ((AliTRDtrackletMCM*) (*fTrackletArray)[fTrackletArray->GetEntriesFast()-1])->SetNHits(fit0->fNhits + fit1->fNhits); ((AliTRDtrackletMCM*) (*fTrackletArray)[fTrackletArray->GetEntriesFast()-1])->SetNHits0(nHits0); ((AliTRDtrackletMCM*) (*fTrackletArray)[fTrackletArray->GetEntriesFast()-1])->SetNHits1(nHits1); @@ -1683,31 +1766,36 @@ void AliTRDmcmSim::FitTracklet() ((AliTRDtrackletMCM*) (*fTrackletArray)[fTrackletArray->GetEntriesFast()-1])->SetOffset(fitOffset); ((AliTRDtrackletMCM*) (*fTrackletArray)[fTrackletArray->GetEntriesFast()-1])->SetError(TMath::Sqrt(TMath::Abs(fitError)/nHits)); -// // cluster information -// Float_t *res = new Float_t[nHits]; -// Float_t *qtot = new Float_t[nHits]; -// Int_t nCls = 0; -// for (Int_t iHit = 0; iHit < fNHits; iHit++) { -// // check if hit contributes -// if (fHits[iHit].fChannel == fFitPtr[cpu]) { -// res[nCls] = fHits[iHit].fYpos - (fitSlope * fHits[iHit].fTimebin + fitOffset); -// qtot[nCls] = fHits[iHit].fQtot; -// nCls++; -// } -// else if (fHits[iHit].fChannel == fFitPtr[cpu] + 1) { -// res[nCls] = fHits[iHit].fYpos + 256 - (fitSlope * fHits[iHit].fTimebin + fitOffset); -// qtot[nCls] = fHits[iHit].fQtot; -// nCls++; -// } -// } -// ((AliTRDtrackletMCM*) (*fTrackletArray)[fTrackletArray->GetEntriesFast()-1])->SetClusters(res, qtot, nCls); -// delete [] res; -// delete [] qtot; + // store cluster information (if requested) + if (fgStoreClusters) { + Float_t *res = new Float_t[fNTimeBin]; + Float_t *qtot = new Float_t[fNTimeBin]; + for (Int_t iTimebin = 0; iTimebin < fNTimeBin; ++iTimebin) { + res[iTimebin] = 0; + qtot[iTimebin] = 0; + } + for (Int_t iHit = 0; iHit < fNHits; iHit++) { + Int_t timebin = fHits[iHit].fTimebin; + + // check if hit contributes + if (fHits[iHit].fChannel == fFitPtr[cpu]) { + res[timebin] = fHits[iHit].fYpos - (fitSlope * timebin + fitOffset); + qtot[timebin] = fHits[iHit].fQtot; + } + else if (fHits[iHit].fChannel == fFitPtr[cpu] + 1) { + res[timebin] = fHits[iHit].fYpos + 256 - (fitSlope * timebin + fitOffset); + qtot[timebin] = fHits[iHit].fQtot; + } + } + ((AliTRDtrackletMCM*) (*fTrackletArray)[fTrackletArray->GetEntriesFast()-1])->SetClusters(res, qtot, fNTimeBin); + delete [] res; + delete [] qtot; + } if (fitError < 0) AliError(Form("Strange fit error: %f from Sx: %i, Sy: %i, Sxy: %i, Sx2: %i, Sy2: %i, nHits: %i", fitError, sumX, sumY, sumXY, sumX2, sumY2, nHits)); - AliDebug(3, Form("fit slope: %f, offset: %f, error: %f", + AliDebug(3, Form("fit slope: %f, offset: %f, error: %f", fitSlope, fitOffset, TMath::Sqrt(TMath::Abs(fitError)/nHits))); } } @@ -1718,7 +1806,7 @@ void AliTRDmcmSim::Tracklet() { // Run the tracklet calculation by calling sequentially: // CalcFitreg(); TrackletSelection(); FitTracklet() - // and store the tracklets + // and store the tracklets if (!fInitialized) { AliError("Called uninitialized! Nothing done!"); @@ -1734,11 +1822,11 @@ void AliTRDmcmSim::Tracklet() FitTracklet(); } -Bool_t AliTRDmcmSim::StoreTracklets() +Bool_t AliTRDmcmSim::StoreTracklets() { // store the found tracklets via the loader - if (fTrackletArray->GetEntriesFast() == 0) + if (fTrackletArray->GetEntriesFast() == 0) return kTRUE; AliRunLoader *rl = AliRunLoader::Instance(); @@ -1755,12 +1843,12 @@ Bool_t AliTRDmcmSim::StoreTracklets() dl->MakeTree(); trackletTree = dl->Tree(); } - + AliTRDtrackletMCM *trkl = 0x0; TBranch *trkbranch = trackletTree->GetBranch(fTrklBranchName.Data()); if (!trkbranch) trkbranch = trackletTree->Branch(fTrklBranchName.Data(), "AliTRDtrackletMCM", &trkl, 32000); - + for (Int_t iTracklet = 0; iTracklet < fTrackletArray->GetEntriesFast(); iTracklet++) { trkl = ((AliTRDtrackletMCM*) (*fTrackletArray)[iTracklet]); trkbranch->SetAddress(&trkl); @@ -1776,14 +1864,14 @@ void AliTRDmcmSim::WriteData(AliTRDarrayADC *digits) // EBSF = 1: unfiltered data; EBSF = 0: filtered data // zero-suppressed valued are written as -1 to digits - if( !CheckInitialized() ) + if( !CheckInitialized() ) return; Int_t offset = (fMcmPos % 4 + 1) * 21 + (fRobPos % 2) * 84 - 1; - if (fTrapConfig->GetTrapReg(AliTRDtrapConfig::kEBSF) != 0) // store unfiltered data + if (fTrapConfig->GetTrapReg(AliTRDtrapConfig::kEBSF, fDetector, fRobPos, fMcmPos) != 0) // store unfiltered data { - for (Int_t iAdc = 0; iAdc < fgkNADC; iAdc++) { + for (Int_t iAdc = 0; iAdc < AliTRDfeeParam::GetNadcMcm(); iAdc++) { if (~fZSMap[iAdc] == 0) { for (Int_t iTimeBin = 0; iTimeBin < fNTimeBin; iTimeBin++) { digits->SetDataByAdcCol(GetRow(), offset - iAdc, iTimeBin, -1); @@ -1797,7 +1885,7 @@ void AliTRDmcmSim::WriteData(AliTRDarrayADC *digits) } } else { - for (Int_t iAdc = 0; iAdc < fgkNADC; iAdc++) { + for (Int_t iAdc = 0; iAdc < AliTRDfeeParam::GetNadcMcm(); iAdc++) { if (~fZSMap[iAdc] != 0) { for (Int_t iTimeBin = 0; iTimeBin < fNTimeBin; iTimeBin++) { digits->SetDataByAdcCol(GetRow(), offset - iAdc, iTimeBin, (fADCF[iAdc][iTimeBin] >> fgkAddDigits) - fgAddBaseline); @@ -1852,133 +1940,60 @@ Int_t AliTRDmcmSim::GetPID(Int_t q0, Int_t q1) ULong64_t addrQ0; ULong64_t addr; - UInt_t nBinsQ0 = fTrapConfig->GetDmemUnsigned(AliTRDtrapConfig::fgkDmemAddrLUTnbins); // number of bins in q0 / 4 !! - UInt_t pidTotalSize = fTrapConfig->GetDmemUnsigned(AliTRDtrapConfig::fgkDmemAddrLUTLength); + UInt_t nBinsQ0 = fTrapConfig->GetDmemUnsigned(fgkDmemAddrLUTnbins, fDetector, fRobPos, fMcmPos); // number of bins in q0 / 4 !! + UInt_t pidTotalSize = fTrapConfig->GetDmemUnsigned(fgkDmemAddrLUTLength, fDetector, fRobPos, fMcmPos); + if(nBinsQ0==0 || pidTotalSize==0) // make sure we don't run into trouble if the value for Q0 is not configured + return 0; // Q1 not configured is ok for 1D LUT - ULong_t corrQ0 = fTrapConfig->GetDmemUnsigned(AliTRDtrapConfig::fgkDmemAddrLUTcor0, fDetector, fRobPos, fMcmPos); - ULong_t corrQ1 = fTrapConfig->GetDmemUnsigned(AliTRDtrapConfig::fgkDmemAddrLUTcor1, fDetector, fRobPos, fMcmPos); + ULong_t corrQ0 = fTrapConfig->GetDmemUnsigned(fgkDmemAddrLUTcor0, fDetector, fRobPos, fMcmPos); + ULong_t corrQ1 = fTrapConfig->GetDmemUnsigned(fgkDmemAddrLUTcor1, fDetector, fRobPos, fMcmPos); + if(corrQ0==0) // make sure we don't run into trouble if one of the values is not configured + return 0; addrQ0 = corrQ0; addrQ0 = (((addrQ0*q0)>>16)>>16); // because addrQ0 = (q0 * corrQ0) >> 32; does not work for unknown reasons - // std::cout << "addrQ0: " << addrQ0 << ", q0: " << q0 << ", corrQ0: " << corrQ0 << std::endl; if(addrQ0 >= nBinsQ0) { // check for overflow - AliDebug(5,Form("Overflow in q0: %i/4 is bigger then %i", addrQ0, nBinsQ0)); + AliDebug(5,Form("Overflow in q0: %llu/4 is bigger then %u", addrQ0, nBinsQ0)); addrQ0 = nBinsQ0 -1; - } + } addr = corrQ1; addr = (((addr*q1)>>16)>>16); addr = addrQ0 + nBinsQ0*addr; // because addr = addrQ0 + nBinsQ0* (((corrQ1*q1)>>32); does not work - // std::cout << "addr: " << addr << ", q1: " << q1 << ", corrQ1: " << corrQ1 << std::endl; if(addr >= pidTotalSize) { - AliDebug(5,Form("Overflow in q1. Address %i/4 is bigger then %i", addr, pidTotalSize)); + AliDebug(5,Form("Overflow in q1. Address %llu/4 is bigger then %u", addr, pidTotalSize)); addr = pidTotalSize -1; - } + } // For a LUT with 11 input and 8 output bits, the first memory address is set to LUT[0] | (LUT[1] << 8) | (LUT[2] << 16) | (LUT[3] << 24) // and so on - UInt_t result = fTrapConfig->GetDmemUnsigned(AliTRDtrapConfig::fgkDmemAddrLUTStart+(addr/4)); + UInt_t result = fTrapConfig->GetDmemUnsigned(fgkDmemAddrLUTStart+(addr/4), fDetector, fRobPos, fMcmPos); return (result>>((addr%4)*8)) & 0xFF; } -void AliTRDmcmSim::SetPIDLut(TH2F * const lut) -{ - // set a user-defined PID LUT from a 2D histogram - - UInt_t nBinsQ0 = lut->GetNbinsX(); - UInt_t nBinsQ1 = lut->GetNbinsY(); - - Double_t scaleQ0 = lut->GetNbinsX() / lut->GetXaxis()->GetXmax(); - Double_t scaleQ1 = lut->GetNbinsY() / lut->GetYaxis()->GetXmax(); - fTrapConfig->SetPIDscale(scaleQ0, scaleQ1); - SetPIDLutScaleDMEM(); - - UInt_t fPIDsizeX = 0; - if(nBinsQ0%4==0) - fPIDsizeX=nBinsQ0/4; - else - fPIDsizeX = (nBinsQ0/4)+1; - - fTrapConfig->SetDmem(AliTRDtrapConfig::fgkDmemAddrLUTnbins, nBinsQ0); // number of bins in q0 - fTrapConfig->SetDmem(AliTRDtrapConfig::fgkDmemAddrLUTLength, nBinsQ0*nBinsQ1); // total size of the table in BYTES (does work because each bin is 8 bit wide) - - UInt_t buffer; - Int_t dmemAddr=0; - if(nBinsQ0*nBinsQ1/4 < AliTRDtrapConfig::fgkDmemAddrLUTEnd - AliTRDtrapConfig::fgkDmemAddrLUTStart) { // /4 because each memory address contains 4 LUT entries - for (UInt_t iy = 0; iy < nBinsQ1; iy++) { - for (UInt_t ix = 0; ix < fPIDsizeX; ix++) { - buffer=0; - for(UInt_t isub=0; isub<4; isub++) { - if(ix*4+isubGetBinContent(ix*4+isub +1, iy +1))) << isub*8 ; - //AliDebug(10, Form("x: %d, y: %d -- %d, %d \n", ix*4+isub, iy, lut->GetBinContent(ix*4+isub, iy), ((Int_t) (255. * lut->GetBinContent(ix*4+isub +1, iy +1))))); - } - else - buffer |= 0<= AliTRDtrapConfig::fgkDmemAddrLUTEnd) { - AliError("LUT table size is too big!"); - } - else { - // AliDebug(8,Form("x: %d, y: %d is memory address %d, setting to %d \n", ix, iy, dmemAddr, buffer)); - fTrapConfig->SetDmem(AliTRDtrapConfig::fgkDmemAddrLUTStart+ix+(iy*fPIDsizeX), buffer); - } - } - } - } - else { - AliError("LUT table is too big!"); - } -} - - -void AliTRDmcmSim::SetPIDLutScaleDMEM() -{ - // set scale factor for PID in DMEM - - Double_t scaleQ[2]; - fTrapConfig->GetPIDscale(scaleQ); - - ULong64_t scale = 1; - scale = scale<<32; - - fTrapConfig->SetDmem(AliTRDtrapConfig::fgkDmemAddrLUTcor0, TMath::Nint(scale*scaleQ[0])); - fTrapConfig->SetDmem(AliTRDtrapConfig::fgkDmemAddrLUTcor1, TMath::Nint(scale*scaleQ[1])); -} - - -void AliTRDmcmSim::SetPIDLut(Int_t* /* lut */, Int_t /* nbinsq0 */, Int_t /* nbinsq1 */) -{ - ; -} - - - - // help functions, to be cleaned up UInt_t AliTRDmcmSim::AddUintClipping(UInt_t a, UInt_t b, UInt_t nbits) const { - // - // This function adds a and b (unsigned) and clips to - // the specified number of bits. - // + // + // This function adds a and b (unsigned) and clips to + // the specified number of bits. + // UInt_t sum = a + b; if (nbits < 32) { UInt_t maxv = (1 << nbits) - 1;; - if (sum > maxv) + if (sum > maxv) sum = maxv; } else { - if ((sum < a) || (sum < b)) + if ((sum < a) || (sum < b)) sum = 0xFFFFFFFF; } return sum; @@ -2150,10 +2165,10 @@ ostream& AliTRDmcmSim::Text(ostream& os) ostream& AliTRDmcmSim::Cfdat(ostream& os) { - // manipulator to activate output in CFDAT format + // manipulator to activate output in CFDAT format // to send to the FEE via SCSN - os.iword(fgkFormatIndex) = 1; + os.iword(fgkFormatIndex) = 1; return os; } @@ -2168,39 +2183,39 @@ ostream& AliTRDmcmSim::Raw(ostream& os) ostream& operator<<(ostream& os, const AliTRDmcmSim& mcm) { // output implementation - + // no output for non-initialized MCM if (!mcm.CheckInitialized()) return os; // ----- human-readable output ----- if (os.iword(AliTRDmcmSim::fgkFormatIndex) == 0) { - - os << "MCM " << mcm.fMcmPos << " on ROB " << mcm.fRobPos << + + os << "MCM " << mcm.fMcmPos << " on ROB " << mcm.fRobPos << " in detector " << mcm.fDetector << std::endl; - + os << "----- Unfiltered ADC data (10 bit) -----" << std::endl; os << "ch "; - for (Int_t iChannel = 0; iChannel < mcm.fgkNADC; iChannel++) + for (Int_t iChannel = 0; iChannel < AliTRDfeeParam::GetNadcMcm(); iChannel++) os << std::setw(5) << iChannel; os << std::endl; for (Int_t iTimeBin = 0; iTimeBin < mcm.fNTimeBin; iTimeBin++) { os << "tb " << std::setw(2) << iTimeBin << ":"; - for (Int_t iChannel = 0; iChannel < mcm.fgkNADC; iChannel++) { + for (Int_t iChannel = 0; iChannel < AliTRDfeeParam::GetNadcMcm(); iChannel++) { os << std::setw(5) << (mcm.fADCR[iChannel][iTimeBin] >> mcm.fgkAddDigits); } os << std::endl; } - + os << "----- Filtered ADC data (10+2 bit) -----" << std::endl; os << "ch "; - for (Int_t iChannel = 0; iChannel < mcm.fgkNADC; iChannel++) + for (Int_t iChannel = 0; iChannel < AliTRDfeeParam::GetNadcMcm(); iChannel++) os << std::setw(4) << iChannel << ((~mcm.fZSMap[iChannel] != 0) ? "!" : " "); os << std::endl; for (Int_t iTimeBin = 0; iTimeBin < mcm.fNTimeBin; iTimeBin++) { os << "tb " << std::setw(2) << iTimeBin << ":"; - for (Int_t iChannel = 0; iChannel < mcm.fgkNADC; iChannel++) { + for (Int_t iChannel = 0; iChannel < AliTRDfeeParam::GetNadcMcm(); iChannel++) { os << std::setw(4) << (mcm.fADCF[iChannel][iTimeBin]) << (((mcm.fZSMap[iChannel] & (1 << iTimeBin)) == 0) ? "!" : " "); } @@ -2213,11 +2228,11 @@ ostream& operator<<(ostream& os, const AliTRDmcmSim& mcm) Int_t dest = 127; Int_t addrOffset = 0x2000; Int_t addrStep = 0x80; - + for (Int_t iTimeBin = 0; iTimeBin < mcm.fNTimeBin; iTimeBin++) { - for (Int_t iChannel = 0; iChannel < mcm.fgkNADC; iChannel++) { - os << std::setw(5) << 10 - << std::setw(5) << addrOffset + iChannel * addrStep + iTimeBin + for (Int_t iChannel = 0; iChannel < AliTRDfeeParam::GetNadcMcm(); iChannel++) { + os << std::setw(5) << 10 + << std::setw(5) << addrOffset + iChannel * addrStep + iTimeBin << std::setw(5) << (mcm.fADCF[iChannel][iTimeBin]) << std::setw(5) << dest << std::endl; } @@ -2229,12 +2244,12 @@ ostream& operator<<(ostream& os, const AliTRDmcmSim& mcm) else if (os.iword(AliTRDmcmSim::fgkFormatIndex) == 2) { Int_t bufSize = 300; UInt_t *buf = new UInt_t[bufSize]; - + Int_t bufLength = mcm.ProduceRawStream(&buf[0], bufSize); - - for (Int_t i = 0; i < bufLength; i++) - std::cout << "0x" << std::hex << buf[i] << std::endl; - + + for (Int_t i = 0; i < bufLength; i++) + std::cout << "0x" << std::hex << buf[i] << std::dec << std::endl; + delete [] buf; } @@ -2264,15 +2279,15 @@ void AliTRDmcmSim::PrintFitRegXml(ostream& os) const os << "" << std::endl; os << " " << std::endl; os << " " << std::endl; - + for(int cpu=0; cpu<4; cpu++) { os << " " << std::endl; if(fFitPtr[cpu] != 31) { for(int adcch=fFitPtr[cpu]; adcch" << adcch << "\">"<< std::endl; + os << " "<< std::endl; os << " " << fFitReg[adcch].fNhits << ""<< std::endl; - os << " " << fFitReg[adcch].fQ0/4 << ""<< std::endl; // divided by 4 because in simulation we have 2 additional decimal places - os << " " << fFitReg[adcch].fQ1/4 << ""<< std::endl; // in the output + os << " " << fFitReg[adcch].fQ0 << ""<< std::endl; + os << " " << fFitReg[adcch].fQ1 << ""<< std::endl; os << " " << fFitReg[adcch].fSumX << ""<< std::endl; os << " " << fFitReg[adcch].fSumX2 << ""<< std::endl; os << " " << fFitReg[adcch].fSumY << ""<< std::endl; @@ -2319,8 +2334,8 @@ void AliTRDmcmSim::PrintTrackletsXml(ostream& os) const offset = (fMCMT[cpu] & 0x1FFF ) ; } - os << " " << pid << "" << " " << padrow << "" - << " " << slope << "" << " " << offset << "" << std::endl; + os << " " << pid << "" << " " << padrow << "" + << " " << slope << "" << " " << offset << "" << "" << std::endl; } os << " " << std::endl; @@ -2332,35 +2347,51 @@ void AliTRDmcmSim::PrintTrackletsXml(ostream& os) const } +void AliTRDmcmSim::PrintAdcDatTxt(ostream& os) const +{ + // print ADC data in text format (suitable as Modelsim stimuli) + + os << "# MCM " << fMcmPos << " on ROB " << fRobPos << + " in detector " << fDetector << std::endl; + + for (Int_t iTimeBin = 0; iTimeBin < fNTimeBin; iTimeBin++) { + for (Int_t iChannel = 0; iChannel < AliTRDfeeParam::GetNadcMcm(); ++iChannel) { + os << std::setw(5) << (fADCR[iChannel][iTimeBin] >> fgkAddDigits); + } + os << std::endl; + } +} + + void AliTRDmcmSim::PrintAdcDatHuman(ostream& os) const { // print ADC data in human-readable format - os << "MCM " << fMcmPos << " on ROB " << fRobPos << + os << "MCM " << fMcmPos << " on ROB " << fRobPos << " in detector " << fDetector << std::endl; - + os << "----- Unfiltered ADC data (10 bit) -----" << std::endl; os << "ch "; - for (Int_t iChannel = 0; iChannel < fgkNADC; iChannel++) + for (Int_t iChannel = 0; iChannel < AliTRDfeeParam::GetNadcMcm(); iChannel++) os << std::setw(5) << iChannel; os << std::endl; for (Int_t iTimeBin = 0; iTimeBin < fNTimeBin; iTimeBin++) { os << "tb " << std::setw(2) << iTimeBin << ":"; - for (Int_t iChannel = 0; iChannel < fgkNADC; iChannel++) { + for (Int_t iChannel = 0; iChannel < AliTRDfeeParam::GetNadcMcm(); iChannel++) { os << std::setw(5) << (fADCR[iChannel][iTimeBin] >> fgkAddDigits); } os << std::endl; } - + os << "----- Filtered ADC data (10+2 bit) -----" << std::endl; os << "ch "; - for (Int_t iChannel = 0; iChannel < fgkNADC; iChannel++) + for (Int_t iChannel = 0; iChannel < AliTRDfeeParam::GetNadcMcm(); iChannel++) os << std::setw(4) << iChannel << ((~fZSMap[iChannel] != 0) ? "!" : " "); os << std::endl; for (Int_t iTimeBin = 0; iTimeBin < fNTimeBin; iTimeBin++) { os << "tb " << std::setw(2) << iTimeBin << ":"; - for (Int_t iChannel = 0; iChannel < fgkNADC; iChannel++) { + for (Int_t iChannel = 0; iChannel < AliTRDfeeParam::GetNadcMcm(); iChannel++) { os << std::setw(4) << (fADCF[iChannel][iTimeBin]) << (((fZSMap[iChannel] & (1 << iTimeBin)) == 0) ? "!" : " "); } @@ -2371,7 +2402,7 @@ void AliTRDmcmSim::PrintAdcDatHuman(ostream& os) const void AliTRDmcmSim::PrintAdcDatXml(ostream& os) const { - // print ADC data in XML format + // print ADC data in XML format os << "" << std::endl; os << "" << std::endl; @@ -2380,7 +2411,7 @@ void AliTRDmcmSim::PrintAdcDatXml(ostream& os) const os << " " << std::endl; os << " " << std::endl; - for(Int_t iChannel = 0; iChannel < fgkNADC; iChannel++) { + for(Int_t iChannel = 0; iChannel < AliTRDfeeParam::GetNadcMcm(); iChannel++) { os << " " << std::endl; for (Int_t iTimeBin = 0; iTimeBin < fNTimeBin; iTimeBin++) { os << "" << fADCF[iChannel][iTimeBin]/4 << ""; @@ -2398,7 +2429,7 @@ void AliTRDmcmSim::PrintAdcDatXml(ostream& os) const -void AliTRDmcmSim::PrintAdcDatDatx(ostream& os, Bool_t broadcast) const +void AliTRDmcmSim::PrintAdcDatDatx(ostream& os, Bool_t broadcast, Int_t timeBinOffset) const { // print ADC data in datx format (to send to FEE) @@ -2408,15 +2439,23 @@ void AliTRDmcmSim::PrintAdcDatDatx(ostream& os, Bool_t broadcast) const Int_t addrOffset = 0x2000; Int_t addrStep = 0x80; Int_t addrOffsetEBSIA = 0x20; - + for (Int_t iTimeBin = 0; iTimeBin < fNTimeBin; iTimeBin++) { - for (Int_t iChannel = 0; iChannel < fgkNADC; iChannel++) { + for (Int_t iChannel = 0; iChannel < AliTRDfeeParam::GetNadcMcm(); iChannel++) { + if ((iTimeBin < timeBinOffset) || (iTimeBin >= fNTimeBin+timeBinOffset)) { if(broadcast==kFALSE) - fTrapConfig->PrintDatx(os, addrOffset+iChannel*addrStep+addrOffsetEBSIA+iTimeBin, (fADCF[iChannel][iTimeBin]/4), GetRobPos(), GetMcmPos()); + fTrapConfig->PrintDatx(os, addrOffset+iChannel*addrStep+addrOffsetEBSIA+iTimeBin, 10, GetRobPos(), GetMcmPos()); else - fTrapConfig->PrintDatx(os, addrOffset+iChannel*addrStep+addrOffsetEBSIA+iTimeBin, (fADCF[iChannel][iTimeBin]/4), 0, 127); - } - os << std::endl; + fTrapConfig->PrintDatx(os, addrOffset+iChannel*addrStep+addrOffsetEBSIA+iTimeBin, 10, 0, 127); + } + else { + if(broadcast==kFALSE) + fTrapConfig->PrintDatx(os, addrOffset+iChannel*addrStep+addrOffsetEBSIA+iTimeBin, (fADCF[iChannel][iTimeBin-timeBinOffset]/4), GetRobPos(), GetMcmPos()); + else + fTrapConfig->PrintDatx(os, addrOffset+iChannel*addrStep+addrOffsetEBSIA+iTimeBin, (fADCF[iChannel][iTimeBin-timeBinOffset]/4), 0, 127); + } + } + os << std::endl; } } @@ -2427,47 +2466,169 @@ void AliTRDmcmSim::PrintPidLutHuman() UInt_t result; - UInt_t addrEnd = AliTRDtrapConfig::fgkDmemAddrLUTStart + fTrapConfig->GetDmemUnsigned(AliTRDtrapConfig::fgkDmemAddrLUTLength)/4; // /4 because each addr contains 4 values - UInt_t nBinsQ0 = fTrapConfig->GetDmemUnsigned(AliTRDtrapConfig::fgkDmemAddrLUTnbins); + UInt_t addrEnd = fgkDmemAddrLUTStart + fTrapConfig->GetDmemUnsigned(fgkDmemAddrLUTLength, fDetector, fRobPos, fMcmPos)/4; // /4 because each addr contains 4 values + UInt_t nBinsQ0 = fTrapConfig->GetDmemUnsigned(fgkDmemAddrLUTnbins, fDetector, fRobPos, fMcmPos); std::cout << "nBinsQ0: " << nBinsQ0 << std::endl; - std::cout << "LUT table length: " << fTrapConfig->GetDmemUnsigned(AliTRDtrapConfig::fgkDmemAddrLUTLength)/4 << std::endl; - - for(UInt_t addr=AliTRDtrapConfig::fgkDmemAddrLUTStart; addr< addrEnd; addr++) { - result = fTrapConfig->GetDmemUnsigned(addr); - std::cout << addr << " # x: " << (addr-AliTRDtrapConfig::fgkDmemAddrLUTStart)%(nBinsQ0/4) << ", y: " <<(addr-AliTRDtrapConfig::fgkDmemAddrLUTStart)/nBinsQ0 << " # " <<((result>>0)&0xFF)/255.0 - << " | " << ((result>>8)&0xFF)/255.0 - << " | " << ((result>>16)&0xFF)/255.0 << " | " << ((result>>24)&0xFF)/255.0 << std::endl; + std::cout << "LUT table length: " << fTrapConfig->GetDmemUnsigned(fgkDmemAddrLUTLength, fDetector, fRobPos, fMcmPos) << std::endl; + + if (nBinsQ0>0) { + for(UInt_t addr=fgkDmemAddrLUTStart; addr< addrEnd; addr++) { + result = fTrapConfig->GetDmemUnsigned(addr, fDetector, fRobPos, fMcmPos); + std::cout << addr << " # x: " << ((addr-fgkDmemAddrLUTStart)%((nBinsQ0)/4))*4 << ", y: " <<(addr-fgkDmemAddrLUTStart)/(nBinsQ0/4) + << " # " <<((result>>0)&0xFF) + << " | " << ((result>>8)&0xFF) + << " | " << ((result>>16)&0xFF) + << " | " << ((result>>24)&0xFF) << std::endl; + } } } -void AliTRDmcmSim::PrintPidLutDatx(ostream& os) const +Bool_t AliTRDmcmSim::ReadPackedConfig(AliTRDtrapConfig *cfg, Int_t hc, UInt_t *data, Int_t size) { - // print PID LUT in datx format (to send to FEE) - - Double_t scaleQ[2]; - fTrapConfig->GetPIDscale(scaleQ); - - ULong64_t scale = 1; - scale = scale<<32; - - os << std::setw(5) << 27; // cmd - os << std::setw(8) << 4; // scaleQ0 - os << std::setw(12) << TMath::Nint(scale*scaleQ[0]); // value - os << std::setw(8) << 1 << std::endl; // destination - - os << std::setw(5) << 27; // cmd - os << std::setw(8) << 5; // scaleQ1 - os << std::setw(12) << TMath::Nint(scale*scaleQ[1]); // value - os << std::setw(8) << 1 << std::endl << std::endl; // destination - - fTrapConfig->PrintMemDatx(os, AliTRDtrapConfig::fgkDmemAddrLUTnbins); - fTrapConfig->PrintMemDatx(os, AliTRDtrapConfig::fgkDmemAddrLUTLength); - - UInt_t addrStart = AliTRDtrapConfig::fgkDmemAddrLUTStart; - UInt_t addrEnd = addrStart + (fTrapConfig->GetDmemUnsigned(AliTRDtrapConfig::fgkDmemAddrLUTLength)/4); // divided by 4 because each addr contains four values - for(UInt_t addr=AliTRDtrapConfig::fgkDmemAddrLUTStart; addr< addrEnd; addr++) { - fTrapConfig->PrintMemDatx(os, addr); - } + // Read the packed configuration from the passed memory block + // + // To be used to retrieve the TRAP configuration from the + // configuration as sent in the raw data. + + AliDebugClass(1, "Reading packed configuration"); + + Int_t det = hc/2; + + Int_t idx = 0; + Int_t err = 0; + Int_t step, bwidth, nwords, exitFlag, bitcnt; + + UShort_t caddr; + UInt_t dat, msk, header, dataHi; + + while (idx < size && *data != 0x00000000) { + + Int_t rob = (*data >> 28) & 0x7; + Int_t mcm = (*data >> 24) & 0xf; + + AliDebugClass(1, Form("Config of det. %3i MCM %i:%02i (0x%08x)", det, rob, mcm, *data)); + data++; + + while (idx < size && *data != 0x00000000) { + + header = *data; + data++; + idx++; + + AliDebugClass(5, Form("read: 0x%08x", header)); + + if (header & 0x01) // single data + { + dat = (header >> 2) & 0xFFFF; // 16 bit data + caddr = (header >> 18) & 0x3FFF; // 14 bit address + + if (caddr != 0x1FFF) // temp!!! because the end marker was wrong + { + if (header & 0x02) // check if > 16 bits + { + dataHi = *data; + AliDebugClass(5, Form("read: 0x%08x", dataHi)); + data++; + idx++; + err += ((dataHi ^ (dat | 1)) & 0xFFFF) != 0; + dat = (dataHi & 0xFFFF0000) | dat; + } + AliDebugClass(5, Form("addr=0x%04x (%s) data=0x%08x\n", caddr, cfg->GetRegName(cfg->GetRegByAddress(caddr)), dat)); + if ( ! cfg->Poke(caddr, dat, det, rob, mcm) ) + AliDebugClass(5, Form("(single-write): non-existing address 0x%04x containing 0x%08x\n", caddr, header)); + if (idx > size) + { + AliDebugClass(5, Form("(single-write): no more data, missing end marker\n")); + return -err; + } + } + else + { + AliDebugClass(5, Form("(single-write): address 0x%04x => old endmarker?\n", caddr)); + return err; + } + } + + else // block of data + { + step = (header >> 1) & 0x0003; + bwidth = ((header >> 3) & 0x001F) + 1; + nwords = (header >> 8) & 0x00FF; + caddr = (header >> 16) & 0xFFFF; + exitFlag = (step == 0) || (step == 3) || (nwords == 0); + + if (exitFlag) + break; + + switch (bwidth) + { + case 15: + case 10: + case 7: + case 6: + case 5: + { + msk = (1 << bwidth) - 1; + bitcnt = 0; + while (nwords > 0) + { + nwords--; + bitcnt -= bwidth; + if (bitcnt < 0) + { + header = *data; + AliDebugClass(5, Form("read 0x%08x", header)); + data++; + idx++; + err += (header & 1); + header = header >> 1; + bitcnt = 31 - bwidth; + } + AliDebugClass(5, Form("addr=0x%04x (%s) data=0x%08x\n", caddr, cfg->GetRegName(cfg->GetRegByAddress(caddr)), header & msk)); + if ( ! cfg->Poke(caddr, header & msk, det, rob, mcm) ) + AliDebugClass(5, Form("(single-write): non-existing address 0x%04x containing 0x%08x\n", caddr, header)); + + caddr += step; + header = header >> bwidth; + if (idx >= size) + { + AliDebugClass(5, Form("(block-write): no end marker! %d words read\n", idx)); + return -err; + } + } + break; + } // end case 5-15 + case 31: + { + while (nwords > 0) + { + header = *data; + AliDebugClass(5, Form("read 0x%08x", header)); + data++; + idx++; + nwords--; + err += (header & 1); + + AliDebugClass(5, Form("addr=0x%04x (%s) data=0x%08x", caddr, cfg->GetRegName(cfg->GetRegByAddress(caddr)), header >> 1)); + if ( ! cfg->Poke(caddr, header >> 1, det, rob, mcm) ) + AliDebugClass(5, Form("(single-write): non-existing address 0x%04x containing 0x%08x\n", caddr, header)); + + caddr += step; + if (idx >= size) + { + AliDebugClass(5, Form("no end marker! %d words read", idx)); + return -err; + } + } + break; + } + default: return err; + } // end switch + } // end block case + } + } // end while + AliDebugClass(5, Form("no end marker! %d words read", idx)); + return -err; // only if the max length of the block reached! }