fPhiLumping(4),
fDebug(0),
fCuts(),
- fRecalculateEta(false),
fRecalculatePhi(false),
fMinQuality(10),
fCache(),
fPhiLumping(4),
fDebug(0),
fCuts(),
- fRecalculateEta(false),
fRecalculatePhi(false),
fMinQuality(10),
fCache(),
fPhiLumping(o.fPhiLumping),
fDebug(o.fDebug),
fCuts(o.fCuts),
- fRecalculateEta(o.fRecalculateEta),
fRecalculatePhi(o.fRecalculatePhi),
fMinQuality(o.fMinQuality),
fCache(o.fCache),
fEtaLumping = o.fEtaLumping;
fPhiLumping = o.fPhiLumping;
fCuts = o.fCuts;
- fRecalculateEta = o.fRecalculateEta;
fRecalculatePhi = o.fRecalculatePhi;
fMinQuality = o.fMinQuality;
fCache = o.fCache;
return static_cast<RingHistos*>(fRingHistos.At(idx));
}
+namespace {
+ Double_t Rng2Cut(UShort_t d, Char_t r, Int_t xbin, TH2* h)
+ {
+ Double_t ret = 1024;
+ if (xbin < 1 && xbin > h->GetXaxis()->GetNbins()) return ret;
+ Int_t ybin = 0;
+ switch(d) {
+ case 1: ybin = 1; break;
+ case 2: ybin = (r=='i' || r=='I') ? 2 : 3; break;
+ case 3: ybin = (r=='i' || r=='I') ? 4 : 5; break;
+ default: return ret;
+ }
+ ret = h->GetBinContent(xbin,ybin);
+ return ret;
+ }
+}
+
//____________________________________________________________________
Double_t
AliFMDDensityCalculator::GetMultCut(UShort_t d, Char_t r, Int_t ieta,
- Bool_t errors) const
+ Bool_t /*errors*/) const
{
//
// Get the multiplicity cut. If the user has set fMultCut (via
// Return:
// Lower cut on multiplicity
//
- return fCuts.GetMultCut(d,r,ieta,errors);
+ return Rng2Cut(d, r, ieta, fLowCuts);
+ // return fCuts.GetMultCut(d,r,ieta,errors);
}
//____________________________________________________________________
Double_t
AliFMDDensityCalculator::GetMultCut(UShort_t d, Char_t r, Double_t eta,
- Bool_t errors) const
+ Bool_t /*errors*/) const
{
//
// Get the multiplicity cut. If the user has set fMultCut (via
// Return:
// Lower cut on multiplicity
//
- return fCuts.GetMultCut(d,r,eta,errors);
+ Int_t ieta = fLowCuts->GetXaxis()->FindBin(eta);
+ return Rng2Cut(d, r, ieta, fLowCuts);
+ // return fCuts.GetMultCut(d,r,eta,errors);
}
-
+
+#ifndef NO_TIMING
+# define START_TIMER(T) if (fDoTiming) T.Start(true)
+# define GET_TIMER(T,V) if (fDoTiming) V = T.CpuTime()
+# define ADD_TIMER(T,V) if (fDoTiming) V += T.CpuTime()
+#else
+# define START_TIMER(T) do {} while (false)
+# define GET_TIMER(T,V) do {} while (false)
+# define ADD_TIMER(T,V) do {} while (false)
+#endif
+
//____________________________________________________________________
Bool_t
AliFMDDensityCalculator::Calculate(const AliESDFMD& fmd,
// Copy to cache : fraction of sum 3.9% of total 2.2%
// Poisson calculation : fraction of sum 18.7% of total 10.6%
// Diagnostics : fraction of sum 3.7% of total 2.1%
- Double_t reEtaTime = 0;
Double_t nPartTime = 0;
Double_t corrTime = 0;
Double_t rePhiTime = 0;
Double_t copyTime = 0;
Double_t poissonTime= 0;
Double_t diagTime = 0;
- if (fDoTiming) totalT.Start(true);
+ START_TIMER(totalT);
Double_t etaCache[20*512]; // Same number of strips per ring
Double_t phiCache[20*512]; // whether it is inner our outer.
Double_t oldPhi = phi;
// --- Re-calculate eta - needed for satelittes ------------
- if (fDoTiming) timer.Start(true);
- if (eta == AliESDFMD::kInvalidEta || fRecalculateEta)
- eta = AliForwardUtil::GetEtaFromStrip(d,r,s,t,ip.Z());
- if (fDoTiming) reEtaTime += timer.CpuTime();
+ START_TIMER(timer);
etaCache[s*nt+t] = eta;
// --- Check this strip ------------------------------------
rh->fGood->Fill(eta);
// --- If we asked to re-calculate phi for (x,y) IP --------
- if (fDoTiming) timer.Start(true);
+ START_TIMER(timer);
if (fRecalculatePhi) {
oldPhi = phi;
phi = AliForwardUtil::GetPhiFromStrip(r, t, phi, ip.X(), ip.Y());
}
phiCache[s*nt+t] = phi;
- if (fDoTiming) rePhiTime += timer.CpuTime();
+ ADD_TIMER(timer,rePhiTime);
// --- Apply phi corner correction to eloss ----------------
if (fUsePhiAcceptance == kPhiCorrectELoss)
d, r, s, t, eta);
// --- Now caluculate Nch for this strip using fits --------
- if (fDoTiming) timer.Start(true);
+ START_TIMER(timer);
Double_t n = 0;
if (cut > 0 && mult > cut) n = NParticles(mult,d,r,eta,lowFlux);
rh->fELoss->Fill(mult);
// rh->fEvsN->Fill(mult,n);
// rh->fEtaVsN->Fill(eta, n);
- if (fDoTiming) nPartTime += timer.CpuTime();
+ ADD_TIMER(timer,nPartTime);
// --- Calculate correction if needed ----------------------
- if (fDoTiming) timer.Start(true);
+ START_TIMER(timer);
// Temporary stuff - remove Correction call
Double_t c = 1;
if (fUsePhiAcceptance == kPhiCorrectNch)
c = AcceptanceCorrection(r,t);
// Double_t c = Correction(d,r,t,eta,lowFlux);
- if (fDoTiming) corrTime += timer.CpuTime();
+ ADD_TIMER(timer,corrTime);
fCorrections->Fill(c);
if (c > 0) n /= c;
// rh->fEvsM->Fill(mult,n);
rh->fGood->Divide(rh->fGood, rh->fTotal, 1, 1, "B");
// --- Make a copy and reset as needed -------------------------
- if (fDoTiming) timer.Start(true);
+ START_TIMER(timer);
TH2D* hclone = fCache.Get(d,r);
// hclone->Reset();
// TH2D* hclone = static_cast<TH2D*>(h->Clone("hclone"));
// hclone->Add(h);
h->Reset();
}
- if (fDoTiming) copyTime += timer.CpuTime();
+ ADD_TIMER(timer,copyTime);
// --- Store Poisson result ------------------------------------
- if (fDoTiming) timer.Start(true);
+ START_TIMER(timer);
TH2D* poisson = rh->fPoisson.Result();
for (Int_t t=0; t < poisson->GetNbinsX(); t++) {
for (Int_t s=0; s < poisson->GetNbinsY(); s++) {
Double_t eta = etaCache[s*nt+t];
// Double_t phi = fmd.Phi(d,r,s,t) * TMath::DegToRad();
// Double_t eta = fmd.Eta(d,r,s,t);
- // if (fRecalculateEta)
- // eta = AliForwardUtil::GetEtaFromStrip(d,r,s,t,ip.Z());
- // if (fRecalculatePhi)
- // phi = AliForwardUtil::GetPhiFromStrip(r, t, phi, ip.X(), ip.Y());
if (fUsePoisson) {
h->Fill(eta,phi,poissonV);
rh->fDensity->Fill(eta, phi, poissonV);
hclone->Fill(eta,phi,poissonV);
}
}
- if (fDoTiming) poissonTime += timer.CpuTime();
+ ADD_TIMER(timer,poissonTime);
// --- Make diagnostics - eloss vs poisson ---------------------
- if (fDoTiming) timer.Start(true);
+ START_TIMER(timer);
Int_t nY = h->GetNbinsY();
Int_t nIn = 0; // Count non-outliers
Int_t nOut = 0; // Count outliers
} // if (outlier)
} // for (iphi)
} // for (ieta)
- Double_t outRatio = Double_t(nOut) / (nIn+nOut);
+ Int_t nTotal = (nIn+nOut);
+ Double_t outRatio = (nTotal > 0 ? Double_t(nOut) / nTotal : 0);
rh->fOutliers->Fill(outRatio);
if (outRatio < fMaxOutliers) rh->fPoisson.FillDiagnostics();
else h->SetBit(AliForwardUtil::kSkipRing);
- if (fDoTiming) diagTime += timer.CpuTime();
+ ADD_TIMER(timer,diagTime);
// delete hclone;
} // for q
} // for d
if (fDoTiming) {
- fHTiming->Fill(1,reEtaTime);
+ // fHTiming->Fill(1,reEtaTime);
fHTiming->Fill(2,nPartTime);
fHTiming->Fill(3,corrTime);
fHTiming->Fill(4,rePhiTime);
AliFMDCorrELossFit::ELossFit::fgLeastWeight,
fMaxParticles);
}
+//_____________________________________________________________________
+Int_t
+AliFMDDensityCalculator::FindMaxWeight(const AliFMDCorrELossFit* cor,
+ UShort_t d, Char_t r, Double_t eta) const
+{
+ //
+ // Find the max weight to use for FMD<i>dr</i> in eta bin @a iEta
+ //
+ // Parameters:
+ // cor Correction
+ // d Detector
+ // r Ring
+ // eta Eta
+ //
+ DGUARD(fDebug, 10, "Find maximum weight in FMD density calculator");
+ if(!cor) return -1;
+
+ AliFMDCorrELossFit::ELossFit* fit = cor->FindFit(d,r,eta, -1);
+ if (!fit) {
+ // AliWarning(Form("No energy loss fit for FMD%d%c at eta=%f", d, r, eta));
+ return -1;
+ }
+ return fit->FindMaxWeight(2*AliFMDCorrELossFit::ELossFit::fgMaxRelError,
+ AliFMDCorrELossFit::ELossFit::fgLeastWeight,
+ fMaxParticles);
+}
//_____________________________________________________________________
void
AliForwardCorrectionManager& fcm = AliForwardCorrectionManager::Instance();
const AliFMDCorrELossFit* cor = fcm.GetELossFit();
cor->CacheBins(fMinQuality);
+ cor->Print(fDebug > 5 ? "RCS" : "C");
TAxis eta(axis.GetNbins(),
axis.GetXmin(),
fLowCuts->GetYaxis()->SetBinLabel(3, "FMD2o");
fLowCuts->GetYaxis()->SetBinLabel(4, "FMD3i");
fLowCuts->GetYaxis()->SetBinLabel(5, "FMD3o");
-
+
for (Int_t i = 0; i < nEta; i++) {
+ Double_t leta = fMaxWeights->GetXaxis()->GetBinCenter(i+1);
Double_t w[5];
- w[0] = fFMD1iMax[i] = FindMaxWeight(cor, 1, 'I', i+1);
- w[1] = fFMD2iMax[i] = FindMaxWeight(cor, 2, 'I', i+1);
- w[2] = fFMD2oMax[i] = FindMaxWeight(cor, 2, 'O', i+1);
- w[3] = fFMD3iMax[i] = FindMaxWeight(cor, 3, 'I', i+1);
- w[4] = fFMD3oMax[i] = FindMaxWeight(cor, 3, 'O', i+1);
- Double_t l[5];
- l[0] = GetMultCut(1, 'I', i+1, false);
- l[1] = GetMultCut(2, 'I', i+1, false);
- l[2] = GetMultCut(2, 'O', i+1, false);
- l[3] = GetMultCut(3, 'I', i+1, false);
- l[4] = GetMultCut(3, 'O', i+1, false);
- for (Int_t j = 0; j < 5; j++) {
+ w[0] = fFMD1iMax[i] = FindMaxWeight(cor, 1, 'I', leta);
+ w[1] = fFMD2iMax[i] = FindMaxWeight(cor, 2, 'I', leta);
+ w[2] = fFMD2oMax[i] = FindMaxWeight(cor, 2, 'O', leta);
+ w[3] = fFMD3iMax[i] = FindMaxWeight(cor, 3, 'I', leta);
+ w[4] = fFMD3oMax[i] = FindMaxWeight(cor, 3, 'O', leta);
+ for (Int_t j = 0; j < 5; j++)
if (w[j] > 0) fMaxWeights->SetBinContent(i+1, j+1, w[j]);
- if (l[j] > 0) fLowCuts->SetBinContent(i+1, j+1, l[j]);
- }
}
+
+ // Cache cuts in histogram
+ fCuts.FillHistogram(fLowCuts);
}
//_____________________________________________________________________
// d->Add(sigma);
d->Add(AliForwardUtil::MakeParameter("maxParticle", fMaxParticles));
+ d->Add(AliForwardUtil::MakeParameter("minQuality", fMinQuality));
d->Add(AliForwardUtil::MakeParameter("method", fUsePoisson));
d->Add(AliForwardUtil::MakeParameter("phiAcceptance",fUsePhiAcceptance));
d->Add(AliForwardUtil::MakeParameter("etaLumping", fEtaLumping));
d->Add(AliForwardUtil::MakeParameter("phiLumping", fPhiLumping));
- d->Add(AliForwardUtil::MakeParameter("recalcEta", fRecalculateEta));
d->Add(AliForwardUtil::MakeParameter("recalcPhi", fRecalculatePhi));
d->Add(AliForwardUtil::MakeParameter("maxOutliers", fMaxOutliers));
d->Add(AliForwardUtil::MakeParameter("outlierCut", fOutlierCut));
}
PFV("Max(particles)", fMaxParticles );
- PFV("Poisson method", fUsePoisson );
+ PFV("Min(quality)", fMinQuality );
+ PFB("Poisson method", fUsePoisson );
PFV("Eta lumping", fEtaLumping );
PFV("Phi lumping", fPhiLumping );
- PFV("Recalculate eta", fRecalculateEta );
- PFV("Recalculate phi", fRecalculatePhi );
- PFV("Use phi acceptance", phiM);
+ PFB("Recalculate phi", fRecalculatePhi );
+ PFB("Use phi acceptance", phiM);
PFV("Lower cut", "");
fCuts.Print();
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff