* @ingroup pwglf_forward_eloss
*/
AliAnalysisTask*
-AddTaskFMDELoss(Bool_t mc, Bool_t useCent, Int_t debug=0)
+AddTaskFMDELoss(Bool_t mc, Bool_t useCent, Int_t debug=0,
+ const Char_t* residuals="")
{
// --- Load libraries ----------------------------------------------
gROOT->LoadClass("AliAODForwardMult", "libPWGLFforward2");
// Debug
task->SetDebug(debug);
+ TString resi(residuals);
+ resi.ToUpper();
+ AliFMDEnergyFitter::EResidualMethod rm = AliFMDEnergyFitter::kNoResiduals;
+ if (!resi.IsNull() && !resi.BeginsWith("no")) {
+ if (resi.BeginsWith("square"))
+ rm = AliFMDEnergyFitter::kResidualSquareDifference;
+ else if (resi.BeginsWith("scale"))
+ rm = AliFMDEnergyFitter::kResidualScaledDifference;
+ else // Anything else gives plain difference and errors in errors
+ rm = AliFMDEnergyFitter::kResidualDifference;
+ }
+ Printf("Got residual: \"%s\" -> %d", resi.Data(), rm);
+ task->GetEnergyFitter().SetStoreResiduals(rm);
+
// --- Set limits on fits the energy -------------------------------
// DO NOT CHANGE THESE UNLESS YOU KNOW WHAT YOU ARE DOING
// Maximum relative error on parameters
if (!opt.Contains("SAME")) {
TH1* frame = new TH1F(GetName(),
- Form("FMD%d%c, eta bin %d",fDet,fRing,fBin),
- 100, 0, 10);
+ Form("FMD%d%c, eta bin %d",fDet,fRing,fBin),
+ 100, 0, 10);
frame->SetMinimum(max/10000);
frame->SetMaximum(max*1.4);
frame->SetXTitle("#Delta / #Delta_{mip}");
frame->Draw();
opt.Append(" SAME");
}
- tot->DrawCopy(opt.Data());
+ TF1* cpy = tot->DrawCopy(opt.Data());
cleanup.Add(tot);
if (vals) {
cleanup.Add(f);
}
min /= 100;
- tot->GetHistogram()->SetMaximum(max);
- tot->GetHistogram()->SetMinimum(min);
- tot->GetHistogram()->GetYaxis()->SetRangeUser(min, max);
+ if (max <= 0) max = 0.1;
+ if (min <= 0) min = 1e-4;
+ cpy->SetMaximum(max);
+ cpy->SetMinimum(min);
+ cpy->GetHistogram()->SetMaximum(max);
+ cpy->GetHistogram()->SetMinimum(min);
+ cpy->GetHistogram()->GetYaxis()->SetRangeUser(min, max);
if (l) l->Draw();
gPad->cd();
if (CHECKPAR(fXi, fEXi, maxRelError)) qual++;
if (CHECKPAR(fSigma, fESigma, maxRelError)) qual++;
if (CHECKPAR(fSigmaN, fESigmaN, maxRelError)) qual++;
+ // Large penalty for large sigma - often a bad fit.
+ if (fSigma > 10*fXi) qual -= 4;
qual += FindMaxWeight(2*maxRelError, leastWeight, fN);
fQuality = qual;
}
fMaxRelParError(AliFMDCorrELossFit::ELossFit::fgMaxRelError),
fMaxChi2PerNDF(AliFMDCorrELossFit::ELossFit::fgMaxChi2nu),
fMinWeight(AliFMDCorrELossFit::ELossFit::fgLeastWeight),
- fDebug(0)
+ fDebug(0),
+ fResidualMethod(kNoResiduals),
+ fSkips(0),
+ fRegularizationCut(3e6)
{
//
// Default Constructor - do not use
fMaxRelParError(AliFMDCorrELossFit::ELossFit::fgMaxRelError),
fMaxChi2PerNDF(AliFMDCorrELossFit::ELossFit::fgMaxChi2nu),
fMinWeight(AliFMDCorrELossFit::ELossFit::fgLeastWeight),
- fDebug(3)
+ fDebug(3),
+ fResidualMethod(kNoResiduals),
+ fSkips(0),
+ fRegularizationCut(3e6)
{
//
// Constructor
fMaxRelParError(o.fMaxRelParError),
fMaxChi2PerNDF(o.fMaxChi2PerNDF),
fMinWeight(o.fMinWeight),
- fDebug(o.fDebug)
+ fDebug(o.fDebug),
+ fResidualMethod(o.fResidualMethod),
+ fSkips(o.fSkips),
+ fRegularizationCut(o.fRegularizationCut)
{
//
// Copy constructor
fCentralityAxis.Set(o.fCentralityAxis.GetNbins(),
o.fCentralityAxis.GetXmin(),
o.fCentralityAxis.GetXmax());
- fDebug = o.fDebug;
- fMaxRelParError= o.fMaxRelParError;
- fMaxChi2PerNDF = o.fMaxChi2PerNDF;
- fMinWeight = o.fMinWeight;
+ fDebug = o.fDebug;
+ fMaxRelParError = o.fMaxRelParError;
+ fMaxChi2PerNDF = o.fMaxChi2PerNDF;
+ fMinWeight = o.fMinWeight;
+ fResidualMethod = o.fResidualMethod;
+ fSkips = o.fSkips;
+ fRegularizationCut = o.fRegularizationCut;
fRingHistos.Clear();
TIter next(&o.fRingHistos);
TIter next(&fRingHistos);
RingHistos* o = 0;
while ((o = static_cast<RingHistos*>(next()))) {
+ if (CheckSkip(o->fDet, o->fRing, fSkips)) {
+ AliWarningF("Skipping FMD%d%c for fitting", o->fDet, o->fRing);
+ continue;
+ }
+
TObjArray* l = o->Fit(d, fEtaAxis, fLowCut, fNParticles,
fMinEntries, fFitRangeBinWidth,
fMaxRelParError, fMaxChi2PerNDF,
- fMinWeight);
+ fMinWeight, fRegularizationCut,
+ fResidualMethod);
if (!l) continue;
for (Int_t i = 0; i < l->GetEntriesFast()-1; i++) { // Last is status
stack[i % nStack]->Add(static_cast<TH1*>(l->At(i)));
TIter next(&fRingHistos);
RingHistos* o = 0;
while ((o = static_cast<RingHistos*>(next()))) {
+ if (CheckSkip(o->fDet, o->fRing, fSkips)) {
+ AliWarningF("Skipping FMD%d%c for correction object", o->fDet, o->fRing);
+ continue;
+ }
+
o->FindBestFits(d, *obj, fEtaAxis, fMaxRelParError,
fMaxChi2PerNDF, fMinWeight);
}
d->Add(AliForwardUtil::MakeParameter("maxRelPerError",fMaxRelParError));
d->Add(AliForwardUtil::MakeParameter("maxChi2PerNDF", fMaxChi2PerNDF));
d->Add(AliForwardUtil::MakeParameter("minWeight", fMinWeight));
+ d->Add(AliForwardUtil::MakeParameter("regCut", fRegularizationCut));
TIter next(&fRingHistos);
RingHistos* o = 0;
o->fDebug = dbg;
}
+//____________________________________________________________________
+namespace {
+ template <typename T>
+ void GetParam(Bool_t& ret, const TCollection* col,
+ const TString& name, T& var)
+ {
+ TObject* o = col->FindObject(name);
+ if (o) AliForwardUtil::GetParameter(o,var);
+ else ret = false;
+ }
+}
+
+//____________________________________________________________________
+Bool_t
+AliFMDEnergyFitter::ReadParameters(const TCollection* col)
+{
+ // Read parameters of this object from a collection
+ //
+ // Parameters:
+ // col Collection to read parameters from
+ //
+ // Return value:
+ // true on success, false otherwise
+ //
+ if (!col) return false;
+ Bool_t ret = true;
+ TAxis* axis = static_cast<TAxis*>(col->FindObject("etaAxis"));
+ if (!axis) ret = false;
+ else {
+ if (axis->GetXbins()->GetArray())
+ fEtaAxis.Set(axis->GetNbins(), axis->GetXbins()->GetArray());
+ else
+ fEtaAxis.Set(axis->GetNbins(), axis->GetXmin(), axis->GetXmax());
+ }
+ GetParam(ret,col,"lowCut", fLowCut);
+ GetParam(ret,col,"nParticles", fNParticles);
+ GetParam(ret,col,"minEntries", fMinEntries);
+ GetParam(ret,col,"subtractBins", fFitRangeBinWidth);
+ GetParam(ret,col,"doFits", fDoFits);
+ GetParam(ret,col,"doObject", fDoMakeObject);
+ GetParam(ret,col,"maxE", fMaxE);
+ GetParam(ret,col,"nEbins", fNEbins);
+ GetParam(ret,col,"increasingBins",fUseIncreasingBins);
+ GetParam(ret,col,"maxRelPerError",fMaxRelParError);
+ GetParam(ret,col,"maxChi2PerNDF", fMaxChi2PerNDF);
+ GetParam(ret,col,"minWeight", fMinWeight);
+ Bool_t dummy;
+ GetParam(dummy,col,"regCut", fRegularizationCut);
+
+ return ret;
+}
+
+//____________________________________________________________________
+Bool_t
+AliFMDEnergyFitter::CheckSkip(UShort_t d, Char_t r, UShort_t skips)
+{
+ UShort_t q = (r == 'I' || r == 'i' ? 0 : 1);
+ UShort_t c = 1 << (d-1);
+ UShort_t t = 1 << (c+q-1);
+
+ return (t & skips) == t;
+}
+
//____________________________________________________________________
void
AliFMDEnergyFitter::Print(Option_t*) const
<< (fUseIncreasingBins ?"yes\n" : "no\n")
<< ind << " max(delta p/p): " << fMaxRelParError << '\n'
<< ind << " max(chi^2/nu): " << fMaxChi2PerNDF << '\n'
- << ind << " min(a_i): " << fMinWeight << std::endl;
+ << ind << " min(a_i): " << fMinWeight << '\n'
+ << ind << " Residuals: ";
+ switch (fResidualMethod) {
+ case kNoResiduals: std::cout << "None"; break;
+ case kResidualDifference: std::cout << "Difference"; break;
+ case kResidualScaledDifference: std::cout << "Scaled difference"; break;
+ case kResidualSquareDifference: std::cout << "Square difference"; break;
+ }
+ std::cout << std::endl;
}
//====================================================================
// Default CTOR
//
DGUARD(fDebug, 3, "Default CTOR AliFMDEnergyFitter::RingHistos");
+ fBest.Expand(0);
}
//____________________________________________________________________
//
DGUARD(fDebug, 3, "Named CTOR AliFMDEnergyFitter::RingHistos: FMD%d%c",
d, r);
+ fBest.Expand(0);
}
//____________________________________________________________________
AliFMDEnergyFitter::RingHistos::RingHistos(const RingHistos& o)
UShort_t minusBins,
Double_t relErrorCut,
Double_t chi2nuCut,
- Double_t minWeight) const
+ Double_t minWeight,
+ Double_t regCut,
+ EResidualMethod residuals) const
{
//
// Fit each histogram to up to @a nParticles particle responses.
return 0;
}
+ // Optional container for residuals
+ TList* resi = 0;
+ if (residuals != kNoResiduals) {
+ resi = new TList();
+ resi->SetName("residuals");
+ resi->SetOwner();
+ l->Add(resi);
+ }
+
// Container of the fit results histograms
TObjArray* pars = new TObjArray(3+nParticles+1);
pars->SetName("FitResults");
}
// Now fit
- AliFMDCorrELossFit::ELossFit* res = FitHist(dist,i+1,lowCut,
- nParticles,minusBins,
- relErrorCut,chi2nuCut,
- minWeight);
+ ELossFit_t* res = FitHist(dist,i+1,lowCut, nParticles,minusBins,
+ relErrorCut,chi2nuCut,minWeight,regCut);
if (!res) continue;
nFitted++;
fBest.AddAt(res, i+1);
// dist->GetListOfFunctions()->Add(res);
+
+ if (residuals != kNoResiduals && resi)
+ CalculateResiduals(residuals, lowCut, dist, res, resi);
// Store eta limits
low = TMath::Min(low,i+1);
Double_t max = total->GetMaximum();
if (max > 0) total->Scale(1/max);
- AliFMDCorrELossFit::ELossFit* resT =
- FitHist(total,nDists+1,lowCut, nParticles,minusBins,
- relErrorCut,chi2nuCut,minWeight);
+ ELossFit_t* resT = FitHist(total,nDists+1,lowCut, nParticles,minusBins,
+ relErrorCut,chi2nuCut,minWeight,regCut);
if (resT) {
// Make histograms for the result of this fit
Double_t chi2 = resT->GetChi2();
}
//____________________________________________________________________
-AliFMDCorrELossFit::ELossFit*
+AliFMDEnergyFitter::RingHistos::ELossFit_t*
AliFMDEnergyFitter::RingHistos::FitHist(TH1* dist,
UShort_t bin,
Double_t lowCut,
UShort_t minusBins,
Double_t relErrorCut,
Double_t chi2nuCut,
- Double_t minWeight) const
+ Double_t minWeight,
+ Double_t regCut) const
{
//
// Fit a signal histogram. First, the bin @f$ b_{min}@f$ with
// Create a fitter object
AliForwardUtil::ELossFitter f(lowCut, maxRange, minusBins);
f.Clear();
- f.SetDebug(fDebug > 2);
-
+ f.SetDebug(fDebug > 2);
+
+ // regularization cut - should be a parameter of the class
+ if (dist->GetEntries() > regCut) {
+ // We should rescale the errors
+ Double_t s = TMath::Sqrt(dist->GetEntries() / regCut);
+ if (fDebug > 0) printf("Error scale: %f ", s);
+ for (Int_t i = 1; i <= dist->GetNbinsX(); i++) {
+ Double_t e = dist->GetBinError(i);
+ dist->SetBinError(i, e * s);
+ }
+ }
// If we are only asked to fit a single particle, return this fit,
// no matter what.
if (nParticles == 1) {
if (!r) return 0;
TF1* ff = new TF1(*r);
dist->GetListOfFunctions()->Add(ff);
- AliFMDCorrELossFit::ELossFit* ret =
- new AliFMDCorrELossFit::ELossFit(0, *ff);
+ ELossFit_t* ret = new ELossFit_t(0, *ff);
ret->CalculateQuality(chi2nuCut, relErrorCut, minWeight);
return ret;
}
// Here, we use the real quality assesor instead of the old
// `CheckResult' to ensure consitency in all output.
- AliFMDCorrELossFit::ELossFit* ret = FindBestFit(bin,
- dist,
- relErrorCut,
- chi2nuCut,
- minWeight);
+ ELossFit_t* ret = FindBestFit(bin, dist, relErrorCut, chi2nuCut, minWeight);
return ret;
}
-#if 0
-//____________________________________________________________________
-TF1*
-AliFMDEnergyFitter::RingHistos::FitHist(TH1* dist,
- Double_t lowCut,
- UShort_t nParticles,
- UShort_t minusBins,
- Double_t relErrorCut,
- Double_t chi2nuCut,
- Double_t minWeight) const
-{
- //
- // Fit a signal histogram. First, the bin @f$ b_{min}@f$ with
- // maximum bin content in the range @f$ [E_{min},\infty]@f$ is
- // found. Then the fit range is set to the bin range
- // @f$ [b_{min}-\Delta b,b_{min}+2\Delta b]@f$, and a 1
- // particle signal is fitted to that. The parameters of that fit
- // is then used as seeds for a fit of the @f$ N@f$ particle response
- // to the data in the range
- // @f$ [b_{min}-\Delta b,N(\Delta_1+\xi_1\log(N))+2N\xi@f$
- //
- // Parameters:
- // dist Histogram to fit
- // lowCut Lower cut @f$ E_{min}@f$ on signal
- // nParticles Max number @f$ N@f$ of convolved landaus to fit
- // minusBins Number of bins @f$ \Delta b@f$ from peak to
- // subtract to get the fit range
- // relErrorCut Cut applied to relative error of parameter.
- // Note, for multi-particle weights, the cut
- // is loosend by a factor of 2
- // chi2nuCut Cut on @f$ \chi^2/\nu@f$ -
- // the reduced @f$\chi^2@f$
- //
- // Return:
- // The best fit function
- //
- DGUARD(fDebug, 3, "Fit histogram in AliFMDEnergyFitter::RingHistos: %s",
- dist->GetName());
- Double_t maxRange = 10;
-
- // Create a fitter object
- AliForwardUtil::ELossFitter f(lowCut, maxRange, minusBins);
- f.Clear();
- f.SetDebug(fDebug > 2);
-
- // If we are only asked to fit a single particle, return this fit,
- // no matter what.
- if (nParticles == 1) {
- TF1* r = f.Fit1Particle(dist, 0);
- if (!r) return 0;
- TF1* ret = new TF1(*r);
- dist->GetListOfFunctions()->Add(ret);
- return ret;
- }
-
- // Fit from 2 upto n particles
- for (Int_t i = 2; i <= nParticles; i++) f.FitNParticle(dist, i, 0);
-
-
-
- // Now, we need to select the best fit
- Int_t nFits = f.GetFitResults().GetEntriesFast();
- TF1* good[nFits];
- for (Int_t i = nFits-1; i >= 0; i--) {
- good[i] = 0;
- TF1* ff = static_cast<TF1*>(f.GetFunctions().At(i));
- // ff->SetLineColor(Color());
- // ff->SetRange(0, maxRange);
- dist->GetListOfFunctions()->Add(new TF1(*ff));
- if (CheckResult(static_cast<TFitResult*>(f.GetFitResults().At(i)),
- relErrorCut, chi2nuCut, minWeight)) {
- good[i] = ff;
- ff->SetLineWidth(2);
- if (fDebug > 1) {
- AliInfoF("Candiate fit: %s", ff->GetName());
- f.GetFitResults().At(i)->Print("V");
- }
- }
- }
- // If all else fails, use the 1 particle fit
- TF1* ret = static_cast<TF1*>(f.GetFunctions().At(0));
-
- // Find the fit with the most valid particles
- for (Int_t i = nFits-1; i >= 0; i--) {
- if (!good[i]) continue;
- if (fDebug > 0)
- Printf("%30s: Choosing fit with n=%d %s",
- dist->GetName(), i+1, good[i]->GetName());
- if (fDebug >= 3)
- f.GetFitResults().At(i)->Print();
- ret = good[i];
- break;
- }
- // Give a warning if we're using fall-back
- if (ret == f.GetFunctions().At(0)) {
- Printf("%30s: Choosing fall-back 1 particle fit", dist->GetName());
- }
- // Copy our result and return (the functions are owned by the fitter)
- TF1* fret = new TF1(*ret);
- return fret;
-}
-
-//____________________________________________________________________
-Bool_t
-AliFMDEnergyFitter::RingHistos::CheckResult(TFitResult* r,
- Double_t relErrorCut,
- Double_t chi2nuCut,
- Double_t minWeight) const
-{
- //
- // Check the result of the fit. Returns true if
- // - @f$ \chi^2/\nu < \max{\chi^2/\nu}@f$
- // - @f$ \Delta p_i/p_i < \delta_e@f$ for all parameters. Note,
- // for multi-particle fits, this requirement is relaxed by a
- // factor of 2
- // - @f$ a_{n} > 10^{-7}@f$ when fitting to an @f$ n@f$
- // particle response
- //
- // Parameters:
- // r Result to check
- // relErrorCut Cut @f$ \delta_e@f$ applied to relative error
- // of parameter.
- // chi2nuCut Cut @f$ \max{\chi^2/\nu}@f$
- //
- // Return:
- // true if fit is good.
- //
- if (fDebug > 10) r->Print();
- TString n = r->GetName();
- n.ReplaceAll("TFitResult-", "");
- Double_t chi2 = r->Chi2();
- Int_t ndf = r->Ndf();
- Double_t red = (ndf >= 0 ? chi2/ndf : 999);
- // Double_t prob = r.Prob();
- Bool_t ret = kTRUE;
-
- // Check that the reduced chi square isn't larger than cut
- if (red > chi2nuCut) {
- if (fDebug > 2) {
- AliWarning(Form("%s: chi^2/ndf=%12.5f/%3d=%12.5f>%12.5f",
- n.Data(), chi2, ndf, (ndf<0 ? 0 : chi2/ndf),
- chi2nuCut)); }
- ret = kFALSE;
- }
-
- // Check each parameter
- UShort_t nPar = r->NPar();
- for (UShort_t i = 0; i < nPar; i++) {
- if (i == kN) continue; // Skip the number parameter
-
- // Get value and error and check value
- Double_t v = r->Parameter(i);
- Double_t e = r->ParError(i);
- if (v == 0) continue;
-
- // Calculate the relative error and check it
- Double_t re = e / v;
- Double_t cut = relErrorCut * (i < kN ? 1 : 2);
- if (re > cut) {
- if (fDebug > 2) {
- AliWarning(Form("%s: Delta %s/%s=%9.5f/%9.5f=%5.1f%%>%5.1f%%",
- n.Data(), r->ParName(i).c_str(),
- r->ParName(i).c_str(), e, v,
- 100*(v == 0 ? 0 : e/v),
- 100*(cut))); }
- ret = kFALSE;
- }
- }
-
- // Check if we have scale parameters
- if (nPar > kN) {
-
- // Check that the last particle has a significant contribution
- Double_t lastScale = r->Parameter(nPar-1);
- if (lastScale <= minWeight) {
- if (fDebug > 2) {
- AliWarningF("%s: %s=%9.6f<%g",
- n.Data(), r->ParName(nPar-1).c_str(),
- lastScale, minWeight);
- }
- ret = kFALSE;
- }
- }
- return ret;
-}
-#endif
-
//__________________________________________________________________
void
return;
}
Int_t nBin = eta.GetNbins();
-
+ if (fBest.GetEntriesFast() <= 0) {
+ AliWarningF("No fits found for %s", GetName());
+ return;
+ }
+
for (Int_t b = 1; b <= nBin; b++) {
TString n(Form(fgkEDistFormat, fName.Data(), b));
TH1D* dist = static_cast<TH1D*>(dists->FindObject(n));
// Ignore empty histograms altoghether
if (!dist || dist->GetEntries() <= 0) continue;
- AliFMDCorrELossFit::ELossFit* best =
- static_cast<AliFMDCorrELossFit::ELossFit*>(fBest.At(b));
- // FindBestFit(b, dist, relErrorCut, chi2nuCut, minWeightCut);
+ ELossFit_t* best = static_cast<ELossFit_t*>(fBest.At(b));
+ if (!best) {
+ AliErrorF("No best fit found @ %d for %s", b, GetName());
+ continue;
+ }
+ // FindBestFit(b, dist, relErrorCut, chi2nuCut, minWeightCut);
best->fDet = fDet;
best->fRing = fRing;
best->fBin = b; //
}
// Double_t eta = fAxis->GetBinCenter(b);
obj.SetFit(fDet, fRing, b, best);
- // new AliFMDCorrELossFit::ELossFit(*best));
+ // new ELossFit_t(*best));
}
}
//__________________________________________________________________
-AliFMDCorrELossFit::ELossFit*
+AliFMDEnergyFitter::RingHistos::ELossFit_t*
AliFMDEnergyFitter::RingHistos::FindBestFit(UShort_t b,
const TH1* dist,
Double_t relErrorCut,
printf("Find best fit for %s ... ", dist->GetName());
// Info("FindBestFit", "%s", dist->GetName());
while ((func = static_cast<TF1*>(next()))) {
- AliFMDCorrELossFit::ELossFit* fit =
- new(fFits[i++]) AliFMDCorrELossFit::ELossFit(0,*func);
+ ELossFit_t* fit = new(fFits[i++]) ELossFit_t(0,*func);
fit->fDet = fDet;
fit->fRing = fRing;
fit->fBin = b;
}
fFits.Sort();
if (fDebug > 1) fFits.Print("s");
- AliFMDCorrELossFit::ELossFit* ret =
- static_cast<AliFMDCorrELossFit::ELossFit*>(fFits.At(i-1));
+ ELossFit_t* ret = static_cast<ELossFit_t*>(fFits.At(i-1));
if (!ret) {
AliWarningF("No fit found for %s, bin %3d", GetName(), b);
return 0;
}
// We have to make a copy here, because other wise the clones array
// will overwrite the address
- return new AliFMDCorrELossFit::ELossFit(*ret);
+ return new ELossFit_t(*ret);
}
+//____________________________________________________________________
+void
+AliFMDEnergyFitter::RingHistos::CalculateResiduals(EResidualMethod mode,
+ Double_t lowCut,
+ TH1* dist,
+ ELossFit_t* fit,
+ TCollection* out) const
+{
+
+ // Clone the input, and reset
+ TH1* resi = static_cast<TH1*>(dist->Clone());
+ TString tit(resi->GetTitle());
+ tit.ReplaceAll("#DeltaE/#DeltaE_{mip}", "Residuals");
+ resi->SetTitle(tit);
+ resi->SetDirectory(0);
+
+ // Set title on Y axis
+ switch (mode) {
+ case kResidualDifference:
+ resi->SetYTitle("h_{i}-f(#Delta_{i}) #pm #delta_{i}");
+ break;
+ case kResidualScaledDifference:
+ resi->SetYTitle("[h_{i}-f(#Delta_{i})]/#delta_{i}"); break;
+ case kResidualSquareDifference:
+ resi->SetYTitle("#chi_{i}^{2}=[h_{i}-f(#Delta_{i})]^{2}/#delta^{2}_{i}");
+ break;
+ default:
+ resi->SetYTitle("Unknown");
+ break;
+ }
+ out->Add(resi);
+
+ // Try to find the function
+ Double_t highCut = dist->GetXaxis()->GetXmax();
+ TString funcName("landau1");
+ if (fit->GetN() > 1)
+ funcName = Form("nlandau%d", fit->GetN());
+ TF1* func = dist->GetFunction(funcName);
+ if (func) func->GetRange(lowCut, highCut);
+ resi->Reset("ICES");
+ resi->GetListOfFunctions()->Clear();
+ resi->SetUniqueID(mode);
+
+ // Reset histogram
+ Int_t nX = resi->GetNbinsX();
+ for (Int_t i = 1; i <= nX; i++) {
+ Double_t x = dist->GetBinCenter(i);
+ if (x < lowCut) continue;
+ if (x > highCut) break;
+
+ Double_t h = dist->GetBinContent(i);
+ Double_t e = dist->GetBinError(i);
+ Double_t r = 0;
+ Double_t er = 0;
+ if (h > 0 && e > 0) {
+ Double_t f = fit->GetC() * fit->Evaluate(x);
+ if (f > 0) {
+ r = h-f;
+ switch (mode) {
+ case kResidualDifference: er = e; break;
+ case kResidualScaledDifference: r /= e; break;
+ case kResidualSquareDifference: r *= r; r /= (e*e); break;
+ default: r = 0; break;
+ }
+ }
+ }
+ resi->SetBinContent(i, r);
+ resi->SetBinError(i, er);
+ }
+}
//____________________________________________________________________
void
class AliFMDEnergyFitter : public TNamed
{
public:
- enum {
- kC = AliForwardUtil::ELossFitter::kC,
- kDelta = AliForwardUtil::ELossFitter::kDelta,
- kXi = AliForwardUtil::ELossFitter::kXi,
- kSigma = AliForwardUtil::ELossFitter::kSigma,
- kSigmaN = AliForwardUtil::ELossFitter::kSigmaN,
- kN = AliForwardUtil::ELossFitter::kN,
- kA = AliForwardUtil::ELossFitter::kA
- };
+ /**
+ * Enumeration of parameters
+ */
+ enum {
+ /** Index of pre-constant @f$ C@f$ */
+ kC = AliForwardUtil::ELossFitter::kC,
+ /** Index of most probable value @f$ \Delta_p@f$ */
+ kDelta = AliForwardUtil::ELossFitter::kDelta,
+ /** Index of Landau width @f$ \xi@f$ */
+ kXi = AliForwardUtil::ELossFitter::kXi,
+ /** Index of Gaussian width @f$ \sigma@f$ */
+ kSigma = AliForwardUtil::ELossFitter::kSigma,
+ /** Index of Gaussian additional width @f$ \sigma_n@f$ */
+ kSigmaN = AliForwardUtil::ELossFitter::kSigmaN,
+ /** Index of Number of particles @f$ N@f$ */
+ kN = AliForwardUtil::ELossFitter::kN,
+ /** Base index of particle strengths @f$ a_i@f$ for
+ @f$i=2,\ldots,N@f$ */
+ kA = AliForwardUtil::ELossFitter::kA
+ };
+ /**
+ * Enumeration of residual methods
+ */
+ enum EResidualMethod {
+ /** Do not calculate residuals */
+ kNoResiduals = 0,
+ /** The residuals stored are the difference, and the errors are
+ stored in the error bars of the histogram. */
+ kResidualDifference,
+ /** The residuals stored are the differences scaled to the error
+ on the data */
+ kResidualScaledDifference,
+ /** The residuals stored are the square difference scale to the
+ square error on the data. */
+ kResidualSquareDifference
+ };
+ /**
+ * FMD ring bits for skipping
+ */
+ enum FMDRingBits {
+ /** FMD1i */
+ kFMD1I=0x01,
+ /** All of FMD1 */
+ kFMD1 =kFMD1I,
+ /** FMD2i */
+ kFMD2I=0x02,
+ /** FMD2o */
+ kFMD2O=0x04,
+ /** All of FMD2 */
+ kFMD2 =kFMD2I|kFMD2O,
+ /** FMD3i */
+ kFMD3I=0x08,
+ /** FMD3o */
+ kFMD3O=0x10,
+ /** All of FMD3 */
+ kFMD3 =kFMD3I|kFMD3O
+ };
/**
* Destructor
*/
* has already been set (using SetEtaAxis), then this parameter will be
* ignored
*/
- void SetupForData(const TAxis& etaAxis);
+ void SetupForData(const TAxis& etaAxis);
/**
* Set the eta axis to use. This will force the code to use this
* eta axis definition - irrespective of whatever axis is passed to
* @param x Wheter to use increasing bin sizes
*/
void SetUseIncreasingBins(Bool_t x) { fUseIncreasingBins = x; }
+ /**
+ * Set whether to make residuals, and in that case how.
+ *
+ * - Square difference: @f$chi_i^2=(h_i - f(x_i))^2/\delta_i^2@f$
+ * - Scaled difference: @f$(h_i - f(x_i))/\delta_i@f$
+ * - Difference: @f$(h_i - f(x_i)) \pm\delta_i@f$
+ *
+ * where @f$h_i, x_i, \delta_i@f$ is the bin content, bin center,
+ * and bin error for bin @f$i@f$ respectively, and @f$ f@f$ is the
+ * fitted function.
+ *
+ * @param x Residual method
+ */
+ void SetStoreResiduals(EResidualMethod x=kResidualDifference)
+ {
+ fResidualMethod = x;
+ }
+ /**
+ * Set the regularization cut @f$c_{R}@f$. If a @f$\Delta@f$
+ * distribution has more entries @f$ N_{dist}@f$ than @f$c_{R}@f$,
+ * then we modify the errors of the the distribution with the factor
+ *
+ * @f[
+ * \sqrt{N_{dist}/c_{R}}
+ * @f]
+ *
+ * to keep the @f$\chi^2/\nu@f$ within resonable limits.
+ *
+ * The large residuals @f$chi_i^2=(h_i - f(x_i))^2/\delta_i^2@f$
+ * (see also SetStoreResiduals) comes about on the boundary between
+ * the @f$N@f$ and @f$N+1@f$ particle contributions, and seems to
+ * fall off for larger @f$N@f$. This may indicate that there's a
+ * component in the distributions that the function
+ *
+ * @f[
+ * f(\Delta;\Delta_p,\xi,\sigma,\mathbf{a}) = \sum_i=1^{n} a_i\int
+ * d\Delta' L(\Delta;\Delta',\xi) G(\Delta';\Delta_p,\sigma)
+ * @f]
+ *
+ * does not capture.
+ *
+ * @param cut
+ */
+ void SetRegularizationCut(Double_t cut=3e6)
+ {
+ fRegularizationCut = cut;
+ }
/**
* Fitter the input AliESDFMD object
*
* @param option Not used
*/
void Print(Option_t* option="") const;
+ /**
+ * Read the parameters from a list - used when re-running the code
+ *
+ * @param list Input list
+ *
+ * @return true if the parameter where read
+ */
+ Bool_t ReadParameters(const TCollection* list);
+ void SetSkips(UShort_t skip) { fSkips = skip; }
protected:
/**
* Internal data structure to keep track of the histograms
*/
struct RingHistos : public AliForwardUtil::RingHistos
{
+ typedef AliFMDCorrELossFit::ELossFit ELossFit_t;
/**
* Default CTOR
*/
* is loosend by a factor of 2
* @param chi2nuCut Cut on @f$ \chi^2/\nu@f$ -
* the reduced @f$\chi^2@f$
+ * @param regCut Regularization cut-off
+ * @param residuals Mode for residual plots
*
* @return List of fits
*/
- TObjArray* Fit(TList* dir,
- const TAxis& eta,
- Double_t lowCut,
- UShort_t nParticles,
- UShort_t minEntries,
- UShort_t minusBins,
- Double_t relErrorCut,
- Double_t chi2nuCut,
- Double_t minWeight) const;
+ TObjArray* Fit(TList* dir,
+ const TAxis& eta,
+ Double_t lowCut,
+ UShort_t nParticles,
+ UShort_t minEntries,
+ UShort_t minusBins,
+ Double_t relErrorCut,
+ Double_t chi2nuCut,
+ Double_t minWeight,
+ Double_t regCut,
+ EResidualMethod residuals) const;
/**
* Fit a signal histogram. First, the bin @f$ b_{min}@f$ with
* maximum bin content in the range @f$ [E_{min},\infty]@f$ is
* is loosend by a factor of 2
* @param chi2nuCut Cut on @f$ \chi^2/\nu@f$ -
* the reduced @f$\chi^2@f$
+ * @param regCut Regularization cut-off
*
* @return The best fit function
*/
- AliFMDCorrELossFit::ELossFit* FitHist(TH1* dist,
- UShort_t bin,
- Double_t lowCut,
- UShort_t nParticles,
- UShort_t minusBins,
- Double_t relErrorCut,
- Double_t chi2nuCut,
- Double_t minWeight) const;
-#if 0
+ ELossFit_t* FitHist(TH1* dist,
+ UShort_t bin,
+ Double_t lowCut,
+ UShort_t nParticles,
+ UShort_t minusBins,
+ Double_t relErrorCut,
+ Double_t chi2nuCut,
+ Double_t minWeight,
+ Double_t regCut) const;
/**
- * Fit a signal histogram. First, the bin @f$ b_{min}@f$ with
- * maximum bin content in the range @f$ [E_{min},\infty]@f$ is
- * found. Then the fit range is set to the bin range
- * @f$ [b_{min}-\Delta b,b_{min}+2\Delta b]@f$, and a 1
- * particle signal is fitted to that. The parameters of that fit
- * is then used as seeds for a fit of the @f$ N@f$ particle response
- * to the data in the range
- * @f$ [b_{min}-\Delta b,N(\Delta_1+\xi_1\log(N))+2N\xi@f$
- *
- * @param dist Histogram to fit
- * @param lowCut Lower cut @f$ E_{min}@f$ on signal
- * @param nParticles Max number @f$ N@f$ of convolved landaus to fit
- * @param minusBins Number of bins @f$ \Delta b@f$ from peak to
- * subtract to get the fit range
- * @param relErrorCut Cut applied to relative error of parameter.
- * Note, for multi-particle weights, the cut
- * is loosend by a factor of 2
- * @param chi2nuCut Cut on @f$ \chi^2/\nu@f$ -
- * the reduced @f$\chi^2@f$
+ * Calculate residuals of the fit
*
- * @return The best fit function
+ * @param mode How to calculate
+ * @param lowCut Lower cut
+ * @param dist Distribution
+ * @param fit Function fitted to distribution
+ * @param out Output list to store residual histogram in
*/
- TF1* FitHist(TH1* dist,
- Double_t lowCut,
- UShort_t nParticles,
- UShort_t minusBins,
- Double_t relErrorCut,
- Double_t chi2nuCut,
- Double_t minWeight) const;
-#endif
+ void CalculateResiduals(EResidualMethod mode,
+ Double_t lowCut,
+ TH1* dist,
+ ELossFit_t* fit,
+ TCollection* out) const;
/**
* Find the best fits
*
*
* @return Best fit
*/
- AliFMDCorrELossFit::ELossFit* FindBestFit(UShort_t b,
- const TH1* dist,
- Double_t relErrorCut,
- Double_t chi2nuCut,
- Double_t minWeightCut) const;
+ ELossFit_t* FindBestFit(UShort_t b,
+ const TH1* dist,
+ Double_t relErrorCut,
+ Double_t chi2nuCut,
+ Double_t minWeightCut) const;
#if 0
/**
* Check the result of the fit. Returns true if
* @return Ring histogram container
*/
RingHistos* GetRingHistos(UShort_t d, Char_t r) const;
+ /**
+ * Check if the detector @a d, ring @a r is listed <i>in</i> the @a
+ * skips bit mask. If the detector/ring is in the mask, return true.
+ *
+ * That is, use case is
+ * @code
+ * for (UShort_t d=1. d<=3, d++) {
+ * UShort_t nr = (d == 1 ? 1 : 2);
+ * for (UShort_t q = 0; q < nr; q++) {
+ * Char_t r = (q == 0 ? 'I' : 'O');
+ * if (CheckSkips(d, r, skips)) continue;
+ * // Process detector/ring
+ * }
+ * }
+ * @endcode
+ *
+ * @param d Detector
+ * @param r Ring
+ * @param skips Mask of detector/rings to skip
+ *
+ * @return True if detector @a d, ring @a r is in the mask @a skips
+ */
+ static Bool_t CheckSkip(UShort_t d, Char_t r, UShort_t skips);
TList fRingHistos; // List of histogram containers
Double_t fLowCut; // Low cut on energy
Double_t fMaxChi2PerNDF; // chi^2/nu cit
Double_t fMinWeight; // Minimum weight value
Int_t fDebug; // Debug level
-
+ EResidualMethod fResidualMethod;// Whether to store residuals (debugging)
+ UShort_t fSkips; // Rings to skip when fitting
+ Double_t fRegularizationCut; // When to regularize the chi^2
- ClassDef(AliFMDEnergyFitter,4); //
+ ClassDef(AliFMDEnergyFitter,5); //
};
#endif
#include <TDirectory.h>
#include <TTree.h>
#include <TFile.h>
+#include <TROOT.h>
#include "AliMCEvent.h"
#include "AliGenHijingEventHeader.h"
#include "AliHeader.h"
TAxis vAxis(10,-10,10); // Default only
fEnergyFitter.SetupForData(eAxis);
fEventInspector.SetupForData(vAxis);
-
+ Print();
}
//____________________________________________________________________
// cnt++;
// Get the input data
DGUARD(fDebug,3,"Analyse event of AliFMDEnergyFitterTask");
-
- AliMCEvent* mcevent = MCEvent();
- if(mcevent) {
- AliHeader* header = mcevent->Header();
- AliGenHijingEventHeader* hijingHeader =
- dynamic_cast<AliGenHijingEventHeader*>(header->GenEventHeader());
- if(hijingHeader) {
- Float_t b = hijingHeader->ImpactParameter();
- if(b<fbLow || b>fbHigh) return;
- else
- std::cout<<"Selecting event with impact parameter "<<b<<std::endl;
- }
- }
-
AliESDEvent* esd = dynamic_cast<AliESDEvent*>(InputEvent());
// AliInfo(Form("Event # %6d (esd=%p)", cnt, esd));
if (!esd) {
// On the first event, initialize the parameters
if (fFirstEvent && esd->GetESDRun()) {
+ fEventInspector.SetMC(MCEvent());
fEventInspector.ReadRunDetails(esd);
AliInfo(Form("Initializing with parameters from the ESD:\n"
//____________________________________________________________________
void
-AliFMDEnergyFitterTask::Print(Option_t*) const
+AliFMDEnergyFitterTask::Print(Option_t* option) const
{
//
// Print information
// Parameters:
// option Not used
//
+ std::cout << ClassName() << ": " << GetName() << std::endl;
+ gROOT->IncreaseDirLevel();
+
+ fEventInspector.Print(option);
+ fEnergyFitter.Print(option);
+
+ gROOT->DecreaseDirLevel();
}
//
fHCentVsQual(0),
fHStatus(0),
fHVtxStatus(0),
- fHTrgStatus(0),
+ fHTrgStatus(0),
fLowFluxCut(1000),
fMaxVzErr(0.2),
fList(0),
fVtxAxis(10,-10,10),
fUseFirstPhysicsVertex(false),
fUseV0AND(false),
- fMinPileupContrib(3),
+ fMinPileupContrib(3),
fMinPileupDistance(0.8),
- fUseDisplacedVertices(false),
+ fUseDisplacedVertices(false),
fDisplacedVertex(),
fCollWords(),
fBgWords(),
fMinCent(-1.0),
fMaxCent(-1.0),
fUsepA2012Vertex(false),
- fRunNumber(0)
+ fRunNumber(0),
+ fMC(false)
{
//
// Constructor
fMinCent(-1.0),
fMaxCent(-1.0),
fUsepA2012Vertex(false),
- fRunNumber(0)
+ fRunNumber(0),
+ fMC(false)
{
//
// Constructor
fMinCent(o.fMinCent),
fMaxCent(o.fMaxCent),
fUsepA2012Vertex(o.fUsepA2012Vertex),
- fRunNumber(o.fRunNumber)
+ fRunNumber(o.fRunNumber),
+ fMC(o.fMC)
{
//
fMaxCent = o.fMaxCent;
fUsepA2012Vertex = o.fUsepA2012Vertex;
fRunNumber = o.fRunNumber;
+ fMC = o.fMC;
if (fList) {
fList->SetName(GetName());
AliForwardUtil::AliROOTRevision()));
fList->Add(AliForwardUtil::MakeParameter("alirootBranch",
AliForwardUtil::AliROOTBranch()));
+ fList->Add(AliForwardUtil::MakeParameter("mc", fMC));
}
// triggers are fast.
if (TMath::Abs(fEnergy - 2750.) > 20) return false;
if (fCollisionSystem != AliForwardUtil::kPP) return false;
+ if (fMC) return false; // All MC events for pp @ 2.76TeV are `fast'
if (fastonly)
DMSG(fDebug,1,"Fast trigger in pp @ sqrt(s)=2.76TeV removed");
<< ind << " CMS energy per nucleon: " << sNN << '\n'
<< ind << " Field: " << field << '\n'
<< ind << " Satellite events: " << fUseDisplacedVertices<<'\n'
+ << ind << " Use 2012 pA vertex: " << fUsepA2012Vertex << '\n'
+ << ind << " Use PWG-UD vertex: " <<fUseFirstPhysicsVertex<<'\n'
+ << ind << " Simulation input: " << fMC << '\n'
<< ind << " Centrality method: " << fCentMethod << '\n'
<< std::noboolalpha;
if (!fCentAxis) { std::cout << std::flush; return; }
* @param dbg Debug level
*/
void SetDebug(Int_t dbg=1) { fDebug = dbg; }
+ /**
+ * Set whether this is MC or not. Needed by energy loss fitter task
+ * that never instantices AliFMDMCEventInspector. In particular, we
+ * need this to make sure we ignore the FAST partition flag in MC
+ * for 2.76TeV pp.
+ *
+ * @param isMC If true, assume MC input
+ */
+ void SetMC(Bool_t isMC=true) { fMC = isMC; }
/**
* Fetch our histograms from the passed list
*
Double_t fMaxCent; //max centrailty
Bool_t fUsepA2012Vertex;// flag to use pA2012 Veretx selection
ULong_t fRunNumber; // Current run number
-
- ClassDef(AliFMDEventInspector,10); // Inspect the event
+ Bool_t fMC; // Is this MC input
+ ClassDef(AliFMDEventInspector,11); // Inspect the event
};
#endif
// Parameters:
// name Name of object
//
+ fMC = true;
}
//____________________________________________________________________
if (fUseDisplacedVertices) fDisplacedVertex.SetupForData(fList, "", true);
}
-//____________________________________________________________________
-void AliFMDMCEventInspector::StoreInformation()
-{
- // Store information about running conditions in the output list
- if (!fList) return;
- AliFMDEventInspector::StoreInformation();
- Bool_t mc = true;
- fList->Add(AliForwardUtil::MakeParameter("mc", mc));
- // , fProduction.Data());
-}
-
namespace
{
TString& AppendField(TString& s, bool test, const char* f)
Double_t cent, Double_t mcC,
Double_t b,
Int_t npart, Int_t nbin);
- /**
- * Store information about running conditions in output list
- *
- * The 3 TNamed objects from AliFMDEventInspector::StoreInformation
- * are defined. In addition, a fourth TNamed object is defined.
- * The presence of this indicate MC data.
- *
- * - mc Nothing special, and unique id is 1
- */
- virtual void StoreInformation();
/**
* Read the production details
*
#include <TMath.h>
#include <TError.h>
#include <TROOT.h>
+#define FIT_OPTIONS "RNQS"
//====================================================================
ULong_t AliForwardUtil::AliROOTRevision()
{
if (!o) return;
TParameter<int>* p = static_cast<TParameter<int>*>(o);
- if (p->TestBit(BIT(17)))
+ if (p->TestBit(BIT(19)))
value = p->GetVal();
else
value = o->GetUniqueID();
{
if (!o) return;
TParameter<int>* p = static_cast<TParameter<int>*>(o);
- if (p->TestBit(BIT(17)))
+ if (p->TestBit(BIT(19)))
value = p->GetVal();
else
value = o->GetUniqueID();
{
if (!o) return;
TParameter<Long_t>* p = static_cast<TParameter<Long_t>*>(o);
- if (p->TestBit(BIT(17)))
+ if (p->TestBit(BIT(19)))
value = p->GetVal();
else
value = o->GetUniqueID();
{
if (!o) return;
TParameter<double>* p = static_cast<TParameter<double>*>(o);
- if (p->TestBit(BIT(17)))
+ if (p->TestBit(BIT(19)))
value = p->GetVal(); // o->GetUniqueID();
else {
UInt_t i = o->GetUniqueID();
{
if (!o) return;
TParameter<bool>* p = static_cast<TParameter<bool>*>(o);
- if (p->TestBit(BIT(17)))
+ if (p->TestBit(BIT(19)))
value = p->GetVal(); // o->GetUniqueID();
else
value = o->GetUniqueID();
// Do the fit, getting the result object
if (fDebug)
::Info("Fit1Particle", "Fitting in the range %f,%f", minE, maxE);
- TFitResultPtr r = dist->Fit(landau1, "RNQS", "", minE, maxE);
+ TFitResultPtr r = dist->Fit(landau1, FIT_OPTIONS, "", minE, maxE);
// landau1->SetRange(minE, fMaxRange);
fFitResults.AddAtAndExpand(new TFitResult(*r), 0);
fFunctions.AddAtAndExpand(landau1, 0);
// Do the fit
if (fDebug)
::Info("FitNParticle", "Fitting in the range %f,%f (%d)", minE, maxEi, n);
- TFitResultPtr tr = dist->Fit(landaun, "RSQN", "", minE, maxEi);
+ TFitResultPtr tr = dist->Fit(landaun, FIT_OPTIONS, "", minE, maxEi);
// landaun->SetRange(minE, fMaxRange);
fFitResults.AddAtAndExpand(new TFitResult(*tr), n-1);
// Do the fit, getting the result object
if (fDebug)
::Info("FitComposite", "Fitting seed in the range %f,%f", minE, maxE);
- /* TFitResultPtr r = */ dist->Fit(seed, "RNQS", "", minE, maxE);
+ /* TFitResultPtr r = */ dist->Fit(seed, FIT_OPTIONS, "", minE, maxE);
maxE = dist->GetXaxis()->GetXmax();
TF1* comp = new TF1("composite", landauGausComposite,
// Do the fit, getting the result object
if (fDebug)
::Info("FitComposite", "Fitting composite in the range %f,%f", minE, maxE);
- /* TFitResultPtr r = */ dist->Fit(comp, "RNQS", "", minE, maxE);
+ /* TFitResultPtr r = */ dist->Fit(comp, FIT_OPTIONS, "", minE, maxE);
#if 0
TF1* part1 = static_cast<TF1*>(seed->Clone("part1"));
# --- Extract and upload ---------------------------------------------
extract_upload()
{
- dummy extract_upload
+ echo "=== Download, extract, and uploade in `basename $PWD` ==="
+ _extract
+ upload
}
# --- Draw -----------------------------------------------------------
_draw()
{
script $@
}
+# --- Generic draw ---------------------------------------------------
+draw()
+{
+ _draw $@
+}
# --- Draw dN/deta results -------------------------------------------
_dndeta_draw()
{
(cd $d && \
draw ${fwd_dir}/DrawdNdetaSummary.C && \
draw ${scr})
+ echo "Back to `pwd`"
+}
+# --- Draw dN/deta ---------------------------------------------------
+dndeta_draw()
+{
+ _dndeta_draw $@
}
-
# === Task functions =================================================
# --- Common options help --------------------------------------------
usage()
fi
done
+ check setup
+
# Set the date/time string
now=`date '+%Y%m%d_%H%M'`
ln -sf ${name}_acc_${now} last_${name}_acc
ln -sf ${name}_corrs_${now} last_${name}_corrs
cat <<-EOF > ${corrdir}/Browse.C
- TObject* Browse()
+ TObject* Browse()
{
const char* fwd = "$ALICE_ROOT/PWGLF/FORWARD/analysis2";
if (!gROOT->GetClass("AliOADBForward"))
gROOT->Macro(Form("%s/scripts/LoadLibs.C", fwd));
- gROOT->LoadMacro(Form("%s/corrs/ForwardOADBGui.C++g", fwd));
+ gROOT->LoadMacro(Form("%s/corrs/ForwardOADBGui.C", fwd));
AliOADBForward* db = new AliOADBForward;
db->Open("fmd_corrections.root", "*");
EOF
echo "=== Make acceptance corrections"
(cd ${name}_acc_${now} && \
- accGen `$run_for_acc` && \
+ accGen `run_for_acc` && \
upload )
fi
for i in fmd_corrections.root spd_corrections.root deadstrips.C ; do
echo "Linking ${corrdir}/$i here"
ln -fs ${corrdir}/$i .
done
+ print
+
}
# --- Default implementation -----------------------------------------
check_token
local run=$1
script ${fwd_dir}/corrs/ExtractAcceptance.C "${run}"
-
+ if test -f deadstrips.C ; then
+ cp ${here}/${name}_corrs_${now}/
+ fi
}
# === Check function =================================================
case $type in
mc*) mc=1 ;;
esac
-
+
+ opts=
+ url=
train_opts $mc $type $trig
url_opts $mc $type $trig
local files=
case $d in
corr) files="forward_mccorr.pdf" ;;
- eloss) files="corrs*.pdf" ;;
+ eloss) files="forward_elossfits.pdf" ;;
aod) files="forward.pdf" ;;
dndeta*) files="forward_dndeta.pdf dNdeta*.pdf" ;;
multdists) files="forward_multdists.pdf" ;;
pdfnup -q --nup 2x1 -o ${name}_dndeta_${now}.pdf tmp.pdf && \
rm -f tmp.pdf)
echo "Made ${name}_summary_${now}.pdf and ${name}_dndeta_${now}.pdf"
+ rm -f last_${name}_pdfs
+ ln -s ${out} last_${name}_pdfs
}
_collect_files()
*/forward.pdf) tgt=summary_${d}_${M} ;;
*/forward_dndeta.pdf) tgt=summary_${d}_${M} ;;
*/forward_multdists.pdf) tgt=summary_${d}_${M} ;;
- */eloss*.pdf) tgt=summary_${d}_${M} ;;
+ */forward_elossfits.pdf) tgt=summary_${d}_${M} ;;
*/dNdeta*.pdf) tgt=${d}_${M} ;;
*) echo "Don't know how to deal with $ff" >/dev/stderr
esac
# include <TString.h>
# include <TError.h>
#else
-// class SummaryDrawer;
+class SummaryDrawer;
+class TObject;
// class TAxis;
-// class AliFMDCorrAcceptance;
-// class AliFMDCorrSecondaryMap;
-// class AliFMDCorrELossFit;
+class AliFMDCorrAcceptance;
+class AliFMDCorrSecondaryMap;
+class AliFMDCorrELossFit;
+#include <TString.h>
#endif
class CorrDrawer : public SummaryDrawer
Bool_t mc=false,
Bool_t sat=false)
{
+ out = TString::Format("forward_%s.pdf", prefix.Data());
+#if 0
out = TString::Format("%s_run%09lu_%s_%04dGeV_%c%dkG_%s_%s.pdf",
prefix.Data(), runNo,
(sys == 1 ? "pp" :
(field >= 0 ? 'p' : 'm'), TMath::Abs(field),
(mc ? "MC" : "real"),
(sat ? "satellite" : "nominal"));
+#endif
}
/**
* Draw corrections using the correction manager to get them
*
* @param o Object to draw
*/
- void Draw(const TObject* o)
+ virtual void Draw(const TObject* o)
{
if (!o) return;
Warning("CorrDrawer", "Don't know how to draw a %s object",
*
* @param acc Acceptance correction
*/
- void Draw(const AliFMDCorrAcceptance* acc) { Summarize(acc, false); }
+ virtual void Draw(const AliFMDCorrAcceptance* acc) { Summarize(acc, false); }
/**
* Draw a single plot of the mean secondary correction
*
* @param sec Secondary correction
*/
- void Draw(const AliFMDCorrSecondaryMap* sec) { Summarize(sec, false); }
+ virtual void Draw(const AliFMDCorrSecondaryMap* sec) { Summarize(sec, false);}
/**
* Draw a single plot summarizing the energy loss fits
*
* @param sec Energy loss fits
*/
- void Draw(const AliFMDCorrELossFit* fits) { Summarize(fits, false); }
+ virtual void Draw(const AliFMDCorrELossFit* fits) { Summarize(fits, false); }
/**
* A generalized entry to the summarization functions
*
* @param options Options
* @param local Local storage
*/
- void Summarize(const TString& what,
- ULong_t runNo,
- const Char_t* sys,
- UShort_t sNN,
- Short_t field,
- Bool_t mc=false,
- Bool_t sat=false,
- Option_t* options="",
- const char* local="")
+ virtual void Summarize(const TString& what,
+ ULong_t runNo,
+ const Char_t* sys,
+ UShort_t sNN,
+ Short_t field,
+ Bool_t mc=false,
+ Bool_t sat=false,
+ Option_t* options="",
+ const char* local="")
{
Summarize(AliForwardCorrectionManager::ParseFields(what),
runNo, AliForwardUtil::ParseCollisionSystem(sys),
* @param options Options
* @param local Local storage
*/
- void Summarize(UShort_t what,
- ULong_t runNo,
- UShort_t sys,
- UShort_t sNN,
- Short_t field,
- Bool_t mc=false,
- Bool_t sat=false,
- Option_t* options="",
- const char* local="")
+ virtual void Summarize(UShort_t what,
+ ULong_t runNo,
+ UShort_t sys,
+ UShort_t sNN,
+ Short_t field,
+ Bool_t mc=false,
+ Bool_t sat=false,
+ Option_t* options="",
+ const char* local="")
{
AliForwardCorrectionManager& mgr = AliForwardCorrectionManager::Instance();
mgr.SetDebug(true);
* @param o Object to draw
* @param pdf Not used
*/
- void Summarize(const TObject* o, Bool_t pdf=true)
+ virtual void Summarize(const TObject* o, Bool_t pdf=true)
{
if (!o) return;
Warning("CorrDrawer", "Don't know how to draw a %s object",
* @param acc Acceptance correction
* @param pdf If true, do multiple plots. Otherwise a single summary plot
*/
- void Summarize(const AliFMDCorrAcceptance* acc, Bool_t pdf=true)
+ virtual void Summarize(const AliFMDCorrAcceptance* acc, Bool_t pdf=true)
{
CreateCanvas("acceptance.pdf", false, pdf);
DrawIt(acc, pdf);
* @param sec Secondary correction
* @param pdf If true, do multiple plots. Otherwise a single summary plot
*/
- void Summarize(const AliFMDCorrSecondaryMap* sec, Bool_t pdf=true)
+ virtual void Summarize(const AliFMDCorrSecondaryMap* sec, Bool_t pdf)
{
CreateCanvas("scondarymap.pdf", false, pdf);
DrawIt(sec, pdf);
* @param sec Energy loss fits
* @param pdf If true, do multiple plots. Otherwise a single summary plot
*/
- void Summarize(const AliFMDCorrELossFit* fits, Bool_t pdf=true)
+ virtual void Summarize(const AliFMDCorrELossFit* fits, Bool_t pdf)
{
CreateCanvas("elossfits.pdf", false, pdf);
DrawIt(fits, pdf);
if (pdf) CloseCanvas();
}
- static void Summarize(const TString& what = "",
+ static void Summarize(const TString& what = TString(""),
Bool_t mc = false,
- const TString& output = "forward_eloss.root",
- const TString& local = "fmd_corrections.root",
+ const TString& output = TString("forward_eloss.root"),
+ const TString& local = TString("fmd_corrections.root"),
Option_t* options= "")
{
Summarize(AliForwardCorrectionManager::ParseFields(what), mc,
*
* @param o Object to summarize
*/
- void DrawIt(const TObject* o)
+ virtual void DrawIt(const TObject* o)
{
if (!o) return;
Warning("CorrDrawer", "Don't know how to summarize a %s object",
* @param corr Correction
* @param details If true, make a multipage PDF, otherwise plot the mean.
*/
- void DrawIt(const AliFMDCorrAcceptance* corr, Bool_t details=true)
+ virtual void DrawIt(const AliFMDCorrAcceptance* corr, Bool_t details=true)
{
if (!corr || !fCanvas) return;
* @param corr Correction
* @param details If true, make a multipage PDF, otherwise plot the mean.
*/
- void DrawIt(const AliFMDCorrSecondaryMap* corr, bool details=true)
+ virtual void DrawIt(const AliFMDCorrSecondaryMap* corr, bool details)
{
if (!corr || !fCanvas) return;
* @param details If true, make a multipage PDF,
* otherwise plot the parameters.
*/
- void DrawIt(const AliFMDCorrELossFit* corr, bool details=true)
+ virtual void DrawIt(const AliFMDCorrELossFit* corr, bool details)
{
if (!corr || !fCanvas) return;
savDir->cd();
}
fBody->cd();
- TLatex* ll = new TLatex(.5,.8, fCanvas->GetTitle());
+ TLatex* ll = new TLatex(.5,.8, "ESD #rightarrow #Delta-fits"
+ /* fCanvas->GetTitle() */);
ll->SetTextAlign(22);
ll->SetTextSize(0.05);
ll->SetNDC();
for (UShort_t q = 0; q < nQ; q++) {
Char_t r = (q == 0 ? 'I' : 'O');
TList* dists = 0;
+ TList* resis = 0;
if (fitter) {
// Info("", "Fitter: %s", fitter->GetName());
TList* dl =
// Info("", "Detector list: %s", dl->GetName());
dists = static_cast<TList*>(dl->FindObject("EDists"));
// Info("", "Got distributions -> %p", dists);
+ resis = static_cast<TList*>(dl->FindObject("residuals"));
}
}
ClearCanvas();
TObjArray* ra = fits->GetRingArray(d, r);
if (!ra) continue;
- DrawELossFits(d, r, ra, dists);
+ DrawELossFits(d, r, ra, dists, resis);
}
}
}
* @param r
* @param ra
*/
- void DrawELossFits(UShort_t d, Char_t r, TObjArray* ra, TList* dists)
+ void DrawELossFits(UShort_t d, Char_t r, TObjArray* ra,
+ TList* dists, TList* resis)
{
Int_t nPad = 6;
AliFMDCorrELossFit::ELossFit* fit = 0;
Bool_t last = j == nPad-1;
if (j == 0) DivideForRings(true, true);
- Bool_t same = false;
+ Bool_t same = false;
+ TVirtualPad* drawPad = fBody->GetPad(j+1);
+ Int_t subPad = 0;
if (dists) {
// Info("", "Distributions: %s", dists->GetName());
- TH1* dist =
- static_cast<TH1*>(dists->FindObject(Form("FMD%d%c_etabin%03d",
- d,r,fit->GetBin())));
+ TString hName(Form("FMD%d%c_etabin%03d", d,r,fit->GetBin()));
+ TH1* dist = static_cast<TH1*>(dists->FindObject(hName));
+ TH1* resi = 0;
+ if (resis) resi = static_cast<TH1*>(resis->FindObject(hName));
// Info("", "Got histogram -> %p", dist);
+ if (resi) {
+ Bool_t err = resi->GetUniqueID() <= 1;
+ drawPad->SetGridx();
+ if (err) {
+ resi->SetYTitle("#chi^{2}_{bin}=(h-f)^{2}/#delta^{2}h");
+ for (Int_t k=1; k<=resi->GetNbinsX(); k++) {
+ Double_t c = resi->GetBinContent(k);
+ Double_t e = resi->GetBinError(k);
+ if (e <= 0) continue;
+ c *= c;
+ c /= (e*e);
+ resi->SetBinContent(k, c);
+ resi->SetBinError(k, 0);
+ }
+ }
+ drawPad->Divide(1,2,0,0);
+ DrawInPad(drawPad, 2, resi, "HIST", 0x100);
+ subPad = 1;
+ Double_t red = fit->GetNu() > 0 ? fit->GetChi2() / fit->GetNu() : 0;
+ if (red > 0) {
+ drawPad->cd(2);
+ TLine* l = new TLine(resi->GetXaxis()->GetXmin(), red,
+ resi->GetXaxis()->GetXmax(), red);
+ l->SetLineWidth(2);
+ l->SetLineStyle(2);
+ l->Draw();
+ TLatex* cltx = new TLatex(0.5, 0.5,
+ Form("#chi^{2}/#nu=%6.2f", red));
+ cltx->SetNDC();
+ cltx->SetTextAlign(22);
+ cltx->SetTextFont(42);
+ cltx->SetTextSize(0.07);
+ cltx->Draw();
+ cltx->DrawLatex(0.5,0.4,Form("%g", dist->GetEntries()));
+ }
+ }
if (dist) {
// Info("", "Histogram: %s", dist->GetName());
- DrawInPad(fBody, j+1, dist, "HIST", 0x2);
+ DrawInPad(drawPad, subPad, dist, "HIST E", (subPad * 0x100) + 0x2);
same = true;
}
}
// if (same)
- DrawInPad(fBody, j+1, fit,
+ DrawInPad(drawPad, subPad, fit,
Form("comp good values legend %s", (same ? "same" : "")),
0x2);
if (fit->GetQuality() < fMinQuality) {
*
* @ingroup pwglf_forward_scripts_corr
*/
+void Setup(Bool_t compile)
+{
+ const char* post = (compile ? "++g" : "");
+ const char* fwd = "$ALICE_ROOT/PWGLF/FORWARD/analysis2";
+ if (!gROOT->GetClass("AliFMDCorrELossFit"))
+ gROOT->Macro(Form("%s/scripts/LoadLibs.C", fwd));
+ gSystem->AddIncludePath(Form("-I%s/scripts -I%s/corrs -I%s "
+ "-I$ALICE_ROOT/include", fwd, fwd, fwd));
+ gROOT->LoadMacro(Form("%s/scripts/SummaryDrawer.C%s", fwd, post));
+ gROOT->LoadMacro(Form("%s/corrs/CorrDrawer.C%s", fwd, post));
+
+}
/**
* Draw energy loss fits to a multi-page PDF.
*
//__________________________________________________________________
// Load libraries and object
// const char* fwd = "$ALICE_ROOT/PWGLF/FORWARD/analysis2";
- const char* fwd = "$ALICE_ROOT/PWGLF/FORWARD/analysis2";
- gROOT->Macro(Form("%s/scripts/LoadLibs.C", fwd));
- gROOT->LoadMacro(Form("%s/scripts/SummaryDrawer.C", fwd));
- gROOT->LoadMacro(Form("%s/corrs/CorrDrawer.C", fwd));
+ Setup(false);
CorrDrawer d;
d.Summarize(AliForwardCorrectionManager::kELossFits, runNo, sys, sNN, field,
const char* file="forward_eloss.root",
const char* local="fmd_corrections.root")
{
- const char* fwd = "$ALICE_ROOT/PWGLF/FORWARD/analysis2";
- gROOT->Macro(Form("%s/scripts/LoadLibs.C", fwd));
- gROOT->LoadMacro(Form("%s/scripts/SummaryDrawer.C", fwd));
- gROOT->LoadMacro(Form("%s/corrs/CorrDrawer.C", fwd));
+ Setup(true);
CorrDrawer::Summarize(AliForwardCorrectionManager::kELossFits,
mc, file, local);
DrawCorrELoss(Bool_t mc, Bool_t dummy,
const char* file="forward_eloss_rerun.root")
{
- const char* fwd = "$ALICE_ROOT/PWGLF/FORWARD/analysis2";
- if (!gROOT->GetClass("AliAODForwardMult"))
- gROOT->Macro(Form("%s/scripts/LoadLibs.C", fwd));
- gROOT->LoadMacro(Form("%s/scripts/SummaryDrawer.C", fwd));
- gROOT->LoadMacro(Form("%s/corrs/CorrDrawer.C", fwd));
+ Setup(true);
+ Printf("Drawing fit results from %s", file);
TFile* hist = TFile::Open(file, "READ");
if (!hist) {
Error("DrawCorrELoss", "Failed to open %s", file);
CorrDrawer* cd = new CorrDrawer;
cd->fELossExtra = file;
cd->fMinQuality = 8;
- cd->Summarize(fits);
+ cd->Summarize(const_cast<const AliFMDCorrELossFit*>(fits), true);
}
//
#include <TDatime.h>
#include <TParameter.h>
#include <TPaveText.h>
+#include <TBrowser.h>
namespace {
void
#else
class AliOADBForward;
class AliOADBForward::Entry;
+#if 0
class TGFrame;
-class TGLVEtry;
+class TGLVEntry;
class TGHorizontalFrame;
class TGTextButton;
class TGTextEntry;
class TGLayoutHints;
class TGNumberEntry;
#endif
+#endif
struct ForwardOADBGUI
{
void HandleOpen()
{
if (fDB) HandleClose();
- fDB = new AliOADBForward();
+ fDB = new AliOADBForward;
Info("HandleOpen", "Opening DB file %s for tables %s",
fFileText.GetText(), fTablesText.GetText());
- if (!fDB->Open(fFileText.GetText(), fTablesText.GetText(),
- false, true, true)) {
+ TString fn(fFileText.GetText());
+ TString tn(fTablesText.GetText());
+ if (!fDB->Open(fn, tn, false, true, true)) {
Error("HandleOpen", "Failed to open database");
delete fDB;
fDB = 0;
TGMainFrame* ForwardOADBGui(AliOADBForward* db=0)
{
const char* fwd = "$ALICE_ROOT/PWGLF/FORWARD/analysis2";
- // if (!gROOT->GetClass("AliOADBForward"))
- // gSystem->Load("libGui");
- gROOT->Macro(Form("%s/scripts/LoadLibs.C", fwd));
+ if (!gROOT->GetClass("AliOADBForward"))
+ // gSystem->Load("libGui");
+ gROOT->Macro(Form("%s/scripts/LoadLibs.C", fwd));
// gSystem->AddIncludePath(Form("-I%s", fwd));
// gROOT->LoadMacro(Form("%s/corrs/ForwardOADBGUI.C", fwd));
* to. If this is not specified, it defaults to the name of the input
* file with "_rerun" attached to the base name
*/
-void RerunELossFits(const TString& input="forward_eloss.root",
+void RerunELossFits(Bool_t forceSet=false,
+ const TString& input="forward_eloss.root",
const TString& output="")
{
const char* fwd = "$ALICE_ROOT/PWGLF/FORWARD/analysis2";
gROOT->Macro(Form("%s/scripts/LoadLibs.C", fwd));
- TFile* inFile = 0;
- TFile* outFile = 0;
+ TFile* inFile = 0;
+ TFile* outFile = 0;
+ TString outName(output);
+ if (outName.IsNull()) {
+ outName = input;
+ outName.ReplaceAll(".root", "_rerun.root");
+ }
try {
// --- Open input file ---------------------------------------------
TCollection* inEFRes = GetCollection(inFwdRes, "fmdEnergyFitter");
if (!inEFRes) throw TString("Cannot proceed without previous results");
- const TAxis* etaAxis = static_cast<TAxis*>(GetObject(inEFSum, "etaAxis"));
- if (!etaAxis) throw TString("Cannot proceed without eta axis");
-
// --- Open output file --------------------------------------------
- TString outName(output);
- if (outName.IsNull()) {
- outName = input;
- outName.ReplaceAll(".root", "_rerun.root");
- }
outFile = TFile::Open(outName, "RECREATE");
if (!outFile)
- throw TString::Format("Failed to open %s", input.Data());
+ throw TString::Format("Failed to open %s", outName.Data());
// --- Write copy of sum collection to output --------------------
TCollection* outFwdSum = static_cast<TCollection*>(inFwdSum->Clone());
// --- Make our fitter object ------------------------------------
AliFMDEnergyFitter* fitter = new AliFMDEnergyFitter("energy");
- fitter->SetEtaAxis(*etaAxis);
- // Here, we should get the parameters from the file, but since
- // that isn't implemented yet, we do it by hand
+ if (forceSet || !fitter->ReadParameters(inEFSum)) {
- // Set maximum energy loss to consider
- fitter->SetMaxE(15);
- // Set number of energy loss bins
- fitter->SetNEbins(500);
- // Set whether to use increasing bin sizes
- // fitter->SetUseIncreasingBins(true);
- // Set whether to do fit the energy distributions
- fitter->SetDoFits(kTRUE);
- // Set whether to make the correction object
- fitter->SetDoMakeObject(kTRUE);
- // Set the low cut used for energy
- fitter->SetLowCut(0.4);
+ const TAxis* etaAxis = static_cast<TAxis*>(GetObject(inEFSum,"etaAxis"));
+ if (!etaAxis) throw TString("Cannot proceed without eta axis");
+ fitter->SetEtaAxis(*etaAxis);
+
+ // Set maximum energy loss to consider
+ fitter->SetMaxE(15);
+ // Set number of energy loss bins
+ fitter->SetNEbins(500);
+ // Set whether to use increasing bin sizes
+ // fitter->SetUseIncreasingBins(true);
+ // Set whether to do fit the energy distributions
+ fitter->SetDoFits(kTRUE);
+ // Set whether to make the correction object
+ fitter->SetDoMakeObject(kTRUE);
+ // Set the low cut used for energy
+ fitter->SetLowCut(0.4);
+ // Set the number of bins to subtract from maximum of distributions
+ // to get the lower bound of the fit range
+ fitter->SetFitRangeBinWidth(4);
+ // Set the maximum number of landaus to try to fit (max 5)
+ fitter->SetNParticles(5);
+ // Set the minimum number of entries in the distribution before
+ // trying to fit to the data - 10k seems the least we can do
+ fitter->SetMinEntries(10000);
+ // fitter->SetMaxChi2PerNDF(10);
+ // Enable debug
+ }
+ fitter->SetDebug(1);
+ fitter->SetStoreResiduals(AliFMDEnergyFitter::kResidualSquareDifference);
+ // fitter->SetRegularizationCut(3e6);
// Set the number of bins to subtract from maximum of distributions
// to get the lower bound of the fit range
- fitter->SetFitRangeBinWidth(4);
- // Set the maximum number of landaus to try to fit (max 5)
- fitter->SetNParticles(5);
- // Set the minimum number of entries in the distribution before
- // trying to fit to the data - 10k seems the least we can do
- fitter->SetMinEntries(10000);
- // fitter->SetMaxChi2PerNDF(10);
- // Enable debug
- fitter->SetDebug(1);
+ // fitter->SetFitRangeBinWidth(2);
+ // Skip all of FMD2 and 3
+ // fitter->SetSkips(AliFMDEnergyFitter::kFMD2|AliFMDEnergyFitter::kFMD3);
// --- Now do the fits -------------------------------------------
fitter->Print();
// --- Write out new results folder ------------------------------
outFile->cd();
outFwdRes->Write("ForwardResults", TObject::kSingleKey);
+ Printf("Wrote results to \"%s\" (%s)", outName.Data(), outFile->GetName());
}
catch (const TString& e) {
Error("RerunELossFits", e);
}
if (inFile) inFile->Close();
- if (outFile) outFile->Close();
+ if (outFile) {
+ Printf("Wrote new output to \"%s\"", outName.Data());
+ outFile->Close();
+ }
- gROOT->Macro(Form("%s/corrs/DrawCorrELoss.C(false,false)", fwd));
+ gROOT->Macro(Form("%s/corrs/DrawCorrELoss.C(false,false,\"%s\")",
+ fwd, outName.Data()));
}
unzip ../$i > /dev/null 2>&1
touch .zip
fi
- extract ../Extract.C
+ _extract ../Extract.C
upload
cd ..
done
# --- Draw -----------------------------------------------------------
draw()
{
- local scr=$1 ; shift
+ local dd=`pwd`
+ dd=`basename $dd`
+ echo "=== Drawing in $dd using $1"
+ local scr=$1 ; shift
+ case x$scr in
+ x/*) ;;
+ x*) scr=../$scr ;;
+ esac
# local args=$1 ; shift
download
for i in *.zip ; do
if test "X$i" = "X*.zip" ; then continue ; fi
- echo "Will extract $i"
+ echo "--- Will extract $i in $dd"
d=`basename $i .zip`
if test ! -d $d ; then
mkdir -p $d
}
dndeta_draw()
{
+ echo "=== Drawing dN/deta in $1"
+ local top=$1 ; shift
+ cd $top
download
for i in *.zip ; do
if test "X$i" = "X*.zip" ; then continue ; fi
- echo "Will extract $i"
+ echo "--- Will extract $i"
d=`basename $i .zip`
if test ! -d $d ; then
mkdir -p $d
unzip $i -d $d
fi
- _dndeta_draw $d ../Draw.C $@
+ (cd $d && \
+ script ${fwd_dir}/DrawdNdetaSummary.C && \
+ script ../Draw.C)
done
+ cd ..
}
# === Script specific functions ======================================
# --- Extract corrections --------------------------------------------
download()
{
+ echo "=== Executing download in `pwd`"
if test -f .download ; then
- echo "Already downloaded in `basename $PWD`"
+ echo "--- Already downloaded in `basename $PWD`"
return 0
fi
script Download.C
run_for_acc()
{
- local r=`grep -v ^# $runs | awk '{FS=" \n\t"}{printf "%d\n", $1}' | head -n 1`
+ local r=`grep -v ^# ../$runs | awk '{FS=" \n\t"}{printf "%d\n", $1}' | head -n 1`
if test x$r = "x" || test $r -lt 1; then
echo "No run for acceptance correction specified" > /dev/stderr
exit 1
{
local type=$1 ; shift
local trig=$1 ; shift
- opts="--batch"
- _allAboard "$type" "$trig" $@
+ _allAboard "$type" "$trig" --batch $@
if test $watch -lt 1 ; then
cat <<-EOF
// --- Make our canvas -------------------------------------------
TString pdfName(filename);
pdfName.ReplaceAll(".root", ".pdf");
- CreateCanvas(pdfName, what & 0x100);
+ CreateCanvas(pdfName, what & kLandscape);
+ // --- Make a Title page -------------------------------------------
+ DrawTitlePage(file);
+
// --- Possibly make a chapter here ------------------------------
if (what & kCentral && GetCollection(file, "CentralCorrSums"))
MakeChapter("Forward");
CloseCanvas();
}
protected:
+ //____________________________________________________________________
+ void DrawTitlePage(TFile* file)
+ {
+ TCollection* c = GetCollection(file, "ForwardCorrResults");
+
+ fBody->cd();
+
+ Double_t y = .9;
+ TLatex* ltx = new TLatex(.5, y, "ESD+MC #rightarrow Corrections");
+ ltx->SetTextSize(0.07);
+ ltx->SetTextFont(42);
+ ltx->SetTextAlign(22);
+ ltx->SetNDC();
+ ltx->Draw();
+ y -= .075;
+
+ TCollection* ei = GetCollection(c, "fmdEventInspector");
+ if (ei) {
+ Int_t sys=0, sNN=0, field=0, runNo=0;
+
+ if (GetParameter(ei, "sys", sys))
+ DrawParameter(y, "System", (sys == 1 ? "pp" : sys == 2 ? "PbPb" :
+ sys == 3 ? "pPb" : "unknown"));
+ if (GetParameter(ei, "sNN", sNN)) {
+ TString tsNN = TString::Format("%dGeV", sNN);
+ if (sNN >= 10000)
+ tsNN = TString::Format("%5.2f", float(sNN)/1000);
+ else if (sNN >= 1000)
+ tsNN = TString::Format("%4.2f", float(sNN)/1000);
+ DrawParameter(y, "#sqrt{s_{NN}}", tsNN);
+ }
+
+ if (GetParameter(ei, "field", field))
+ DrawParameter(y, "L3 B field", Form("%+2dkG", field));
+
+ if (GetParameter(ei, "runNo", runNo))
+ DrawParameter(y, "Run #", Form("%6d", runNo));
+ }
+
+ PrintCanvas("MC Corrections");
+ }
//____________________________________________________________________
TCollection* GetVertexList(TCollection* parent, const TAxis& axis, Int_t bin)
{
: TrainSetup(name)
{
fOptions.Add("cent", "Use centrality");
+ fOptions.Add("residuals", "MODE", "Optional calculation of residuals", "");
+ fOptions.Set("type", "ESD");
}
protected:
//------------------------------------------------------------------
gROOT->GetMacroPath()));
// --- Check if this is MC ---------------------------------------
- Bool_t mc = mgr->GetMCtruthEventHandler() != 0;
- Bool_t cent = fOptions.Has("cent");
- Int_t verb = fOptions.AsInt("verbose");
+ Bool_t mc = mgr->GetMCtruthEventHandler() != 0;
+ Bool_t cent = fOptions.Has("cent");
+ Int_t verb = fOptions.AsInt("verbose");
+ TString resi = "";
+ if (fOptions.Has("residuals")) resi = fOptions.Get("residuals");
// --- Add the task ----------------------------------------------
- gROOT->Macro(Form("AddTaskFMDELoss.C(%d,%d,%d)", mc, cent, verb));
+ gROOT->Macro(Form("AddTaskFMDELoss.C(%d,%d,%d,\"%s\")",
+ mc, cent, verb, resi.Data()));
}
/**
* Create entrality selection if enabled
{
fOptions.Add("cent", "Use centrality");
fOptions.Set("type", "ESD");
+ fOptions.Add("corr", "DIR", "Corrections dir", "");
}
protected:
//__________________________________________________________________
Form("%d,%d", mc, fOptions.Has("cent"))))
Fatal("CreateTasks", "Failed to add ForwardQA task");
#endif
+
+ TString cor = "";
+ if (fOptions.Has("corr")) cor = fOptions.Get("corr");
+ if (!cor.IsNull()) {
+ fHelper->LoadAux(Form("%s/fmd_corrections.root",cor.Data()), true);
+ }
}
/**
* Create entrality selection if enabled
git://git.uio.no / u / mrichter / AliRoot.git / commitdiff