13#include <TAttMarker.h>
27#include "CorrFit3DCF.h"
28#include "CorrFit3DCFPainter.h"
29#include "CorrFitHDFunc3D.h"
30#include "CorrFitMask.h"
31#include "CorrFitMask3D.h"
33#include "DividedHisto.h"
41 CorrFit3DCF::CorrFit3DCF(e3DMode mode, Int_t parameters) : CorrFitFunc3D(mode, parameters, 3) {
42 for (
int i = 0; i < 3; i++) {
44 fRange[2 * i + 1] = 1.0;
47 Cout::PrintInfo(Form(
"%s must have at least 3 parameters", this->ClassName()), EInfo::kWarning);
53 SetParameterName(NormID(),
"N");
56 CorrFit3DCF::CorrFit3DCF(Int_t parameters) :
CorrFit3DCF(e3DMode::kNormal3R, parameters) {
69 Double_t CorrFit3DCF::GetFunX(Double_t* x, Double_t* params)
const {
80 return GetScaledValue(
EvalCF(X, params), params);
86 Int_t binY = (fYbins.
GetSize() - 1) / 2;
87 Int_t binZ = (fZbins.
GetSize() - 1) / 2;
88 for (
int i = 0; i < fYbins.
GetSize(); i++) {
89 Double_t tempX[3] = {x[0], fYbins.
Get(i), fZbins.
Get(binZ)};
92 for (
int j = 0; j < fZbins.
GetSize(); j++) {
93 Double_t tempX[3] = {x[0], fYbins.
Get(binY), fZbins.
Get(j)};
96 return GetScaledValue(sum, params) / ctn;
99 Double_t CorrFit3DCF::GetFunY(Double_t* x, Double_t* params)
const {
109 return GetScaledValue(
EvalCF(X, params), params);
113 Int_t binZ = (fZbins.
GetSize() - 1) / 2;
114 Int_t binX = (fXbins.
GetSize() - 1) / 2;
115 for (
int i = 0; i < fXbins.
GetSize(); i++) {
116 Double_t tempX[3] = {fXbins.
Get(i), x[0], fZbins.
Get(binZ)};
117 sum += GetScaledValue(
CalculateCF(tempX, params), params);
119 for (
int j = 0; j < fZbins.
GetSize(); j++) {
120 Double_t tempX[3] = {fXbins.
Get(binX), x[0], fZbins.
Get(j)};
121 sum += GetScaledValue(
CalculateCF(tempX, params), params);
126 Double_t CorrFit3DCF::GetFunZ(Double_t* x, Double_t* params)
const {
136 return GetScaledValue(
EvalCF(X, params), params);
140 Int_t binX = (fXbins.
GetSize() - 1) / 2;
141 Int_t binY = (fYbins.
GetSize() - 1) / 2;
142 for (
int i = 0; i < fXbins.
GetSize(); i++) {
143 Double_t tempX[3] = {fXbins.
Get(i), fYbins.
Get(binY), x[0]};
144 sum += GetScaledValue(
CalculateCF(tempX, params), params);
146 for (
int j = 0; j < fYbins.
GetSize(); j++) {
147 Double_t tempX[3] = {fXbins.
Get(binX), fYbins.
Get(j), x[0]};
148 sum += GetScaledValue(
CalculateCF(tempX, params), params);
162 void CorrFit3DCF::SetParametersToTF1(TF1* f)
const {
169 for (
int i = 0; i <= num->GetXaxis()->GetNbins(); i++) {
170 for (
int j = 0; j <= num->GetYaxis()->GetNbins(); j++) {
171 for (
int k = 0; k <= num->GetZaxis()->GetNbins(); k++) {
172 Double_t ea = num->GetBinError(i, j, k);
173 Double_t eb = den->GetBinError(i, j, k);
174 num->SetBinError(i, j, k, TMath::Sqrt(ea * ea + eb * eb));
191 double CorrFit3DCF::GetChiTFD(
const double* )
const {
194 Double_t e, e1, e2, chi;
195 Bool_t useHD = kFALSE;
198 for (
int i = 0; i < cf->GetNBins(); i++) {
211 if (e == 0)
continue;
214#ifdef CF_FIT_TRACKING
216 std::cout << Form(
"%4.3f\t", par[i]);
218 std::cout <<
"->" << f / z << std::endl;
223 double CorrFit3DCF::GetChiTF(
const double* )
const {
227 Bool_t useHD = kFALSE;
229 CorrFitHDFunc3D* cf =
static_cast<CorrFitHDFunc3D*
>(
fHDMaps);
230 for (
int i = 0; i < cf->GetNBins(); i++) {
231 fBinX = cf->GetXBin(i);
232 fBinY = cf->GetYBin(i);
233 fBinZ = cf->GetZBin(i);
238 C = cf->GetBinCFVal(
fBinX, fBinY, fBinZ, useHD);
242 if (e == 0)
continue;
245#ifdef CF_FIT_TRACKING
247 std::cout << Form(
"%4.3f\t", par[i]);
249 std::cout <<
"->" << f << std::endl;
254 double CorrFit3DCF::GetLogTFD(
const double* )
const {
257 Bool_t useHD = kFALSE;
259 CorrFitHDFunc3D* cf =
static_cast<CorrFitHDFunc3D*
>(
fHDMaps);
260 for (
int i = 0; i < cf->GetNBins(); i++) {
261 fBinX = cf->GetXBin(i);
262 fBinY = cf->GetYBin(i);
263 fBinZ = cf->GetZBin(i);
266 C = cf->GetBinCFVal(
fBinX, fBinY, fBinZ, useHD);
267 Double_t logA = (C * (A + B)) / (A * (C + 1.0));
268 Double_t logB = (A + B) / (B * (C + 1.0));
269 Double_t step = -(A * TMath::Log(logA) + B * TMath::Log(logB));
271 if (TMath::IsNaN(step))
continue;
274#ifdef CF_FIT_TRACKING
276 std::cout << Form(
"%4.3f\t", par[i]);
278 std::cout <<
"->" << f << std::endl;
285 if (
fMask ==
nullptr) {
299 GetMask()->Reset(
true);
308 Bool_t useHD = kFALSE;
312 Double_t free_parameters = 0;
325 CorrFit3DCF::~CorrFit3DCF() {}
328 Double_t q[3] = {x, y, z};
339 x[0] = x[1] - fXBinf;
340 x[2] = x[1] + fXBinf;
343 y[0] = y[1] - fYBinf;
344 y[2] = y[1] + fYBinf;
347 z[0] = z[1] - fZBinf;
348 z[2] = z[1] + fZBinf;
349 Int_t xbin0 =
fBinX * 2 - 2;
350 Int_t ybin0 = fBinY * 2 - 2;
351 Int_t zbin0 = fBinZ * 2 - 2;
355 for (
int i = 0; i < 3; i++) {
357 for (
int j = 0; j < 3; j++) {
359 for (
int k = 0; k < 3; k++) {
361 Double_t n_temp = cf->GetDenominatorHD()[xbin0 + i][ybin0 + j][zbin0 + k];
362 if (n_temp < 0)
continue;
368 return val * cf->GetDenominatorSum()[
fBinX][fBinY][fBinZ];
373 Double_t CorrFit3DCF::GetFunXYpp(Double_t* x, Double_t* params)
const {
383 return GetScaledValue(
EvalCF(X, params), params);
388 Double_t tempX[3] = {x[0], x[0], fZbins.
Get(0)};
389 return GetScaledValue(
CalculateCF(tempX, params), params);
392 Double_t CorrFit3DCF::GetFunXYpm(Double_t* x, Double_t* params)
const {
403 return GetScaledValue(
EvalCF(X, params), params);
408 Double_t sY = 2.0 * fYAxisf - x[0];
409 Double_t tempX[3] = {x[0], sY, fZbins.
Get(0)};
410 return GetScaledValue(
CalculateCF(tempX, params), params);
413 Double_t CorrFit3DCF::GetFunXZpp(Double_t* x, Double_t* params)
const {
423 return GetScaledValue(
EvalCF(X, params), params);
428 Double_t tempX[3] = {x[0], fYbins.
Get(0), x[0]};
429 return GetScaledValue(
CalculateCF(tempX, params), params);
432 Double_t CorrFit3DCF::GetFunXZpm(Double_t* x, Double_t* params)
const {
443 return GetScaledValue(
EvalCF(X, params), params);
448 Double_t sZ = 2.0 * fZAxisf - x[0];
449 Double_t tempX[3] = {x[0], fYbins.
Get(0), sZ};
450 return GetScaledValue(
CalculateCF(tempX, params), params);
453 Double_t CorrFit3DCF::GetFunYZpp(Double_t* x, Double_t* params)
const {
463 return GetScaledValue(
EvalCF(X, params), params);
468 Double_t tempX[3] = {fXbins.
Get(0), x[0], x[0]};
469 return GetScaledValue(
CalculateCF(tempX, params), params);
472 Double_t CorrFit3DCF::GetFunYZpm(Double_t* x, Double_t* params)
const {
483 return GetScaledValue(
EvalCF(X, params), params);
488 Double_t sZ = 2.0 * fZAxisf - x[0];
489 Double_t tempX[3] = {fXbins.
Get(0), x[0], sZ};
490 return GetScaledValue(
CalculateCF(tempX, params), params);
493 Double_t CorrFit3DCF::GetFunXYZppp(Double_t* x, Double_t* params)
const {
503 return GetScaledValue(
EvalCF(X, params), params);
508 Double_t tempX[3] = {x[0], x[0], x[0]};
509 return GetScaledValue(
CalculateCF(tempX, params), params);
512 Double_t CorrFit3DCF::GetFunXYZpmp(Double_t* x, Double_t* params)
const {
523 return GetScaledValue(
EvalCF(X, params), params);
528 Double_t sY = 2.0 * fYAxisf - x[0];
529 Double_t tempX[3] = {x[0], sY, x[0]};
530 return GetScaledValue(
CalculateCF(tempX, params), params);
533 Double_t CorrFit3DCF::GetFunXYZpmm(Double_t* x, Double_t* params)
const {
545 return GetScaledValue(
EvalCF(X, params), params);
550 Double_t sY = 2.0 * fYAxisf - x[0];
551 Double_t sZ = 2.0 * fZAxisf - x[0];
552 Double_t tempX[3] = {x[0], sY, sZ};
553 return GetScaledValue(
CalculateCF(tempX, params), params);
556 Double_t CorrFit3DCF::GetFunXYZppm(Double_t* x, Double_t* params)
const {
567 return GetScaledValue(
EvalCF(X, params), params);
572 Double_t sZ = 2.0 * fZAxisf - x[0];
573 Double_t tempX[3] = {x[0], x[0], sZ};
574 return GetScaledValue(
CalculateCF(tempX, params), params);
577 Double_t CorrFit3DCF::GetFunXY2d(Double_t* x, Double_t* params)
const {
587 return GetScaledValue(
EvalCF(X, params), params);
592 Double_t tempX[3] = {x[1], x[0], fZbins.
Get(0)};
593 return GetScaledValue(
CalculateCF(tempX, params), params);
596 Double_t CorrFit3DCF::GetFunXZ2d(Double_t* x, Double_t* params)
const {
606 return GetScaledValue(
EvalCF(X, params), params);
611 Double_t tempX[3] = {x[1], fYbins.
Get(0), x[0]};
612 return GetScaledValue(
CalculateCF(tempX, params), params);
615 Double_t CorrFit3DCF::GetFunYZ2d(Double_t* x, Double_t* params)
const {
625 return GetScaledValue(
EvalCF(X, params), params);
630 Double_t tempX[3] = {fXbins.
Get(0), x[1], x[0]};
631 return GetScaledValue(
CalculateCF(tempX, params), params);
634 Double_t CorrFit3DCF::GetScaledValue(Double_t x, Double_t* params)
const {
640 void CorrFit3DCF::Calculatef(Double_t width) {
641 TH3* h = (TH3*) ((DividedHisto1D*)
fCF)->GetNum();
645 Hal::Std::GetAxisPar(*h, nbins, min, max,
"x");
646 fXAxisf = 0.5 * (min + max);
647 Hal::Std::GetAxisPar(*h, nbins, min, max,
"y");
648 fYAxisf = 0.5 * (min + max);
649 Hal::Std::GetAxisPar(*h, nbins, min, max,
"z");
650 fZAxisf = 0.5 * (min + max);
651 const Femto3DCF* CF =
static_cast<Femto3DCF*
>(
fCF);
652 Int_t middle_x = CF->GetNum()->GetXaxis()->FindBin(0.0);
653 Int_t middle_y = CF->GetNum()->GetYaxis()->FindBin(0.0);
654 Int_t middle_z = CF->GetNum()->GetZaxis()->FindBin(0.0);
658 fXbins[0] = CF->GetNum()->GetXaxis()->GetBinCenter(middle_x);
659 fYbins[0] = CF->GetNum()->GetYaxis()->GetBinCenter(middle_y);
660 fZbins[0] = CF->GetNum()->GetZaxis()->GetBinCenter(middle_z);
665 fXbins[0] = CF->GetNum()->GetXaxis()->GetBinLowEdge(middle_x);
666 fYbins[0] = CF->GetNum()->GetYaxis()->GetBinLowEdge(middle_y);
667 fZbins[0] = CF->GetNum()->GetZaxis()->GetBinLowEdge(middle_z);
668 fXbins[1] = CF->GetNum()->GetXaxis()->GetBinCenter(middle_x);
669 fYbins[1] = CF->GetNum()->GetYaxis()->GetBinCenter(middle_y);
670 fZbins[1] = CF->GetNum()->GetZaxis()->GetBinCenter(middle_z);
671 fXbins[2] = CF->GetNum()->GetXaxis()->GetBinUpEdge(middle_x);
672 fYbins[2] = CF->GetNum()->GetYaxis()->GetBinUpEdge(middle_y);
673 fZbins[2] = CF->GetNum()->GetZaxis()->GetBinUpEdge(middle_z);
675 }
else if (width != 0) {
676 fXbins.
Resize(2 * width + 1);
677 fYbins.
Resize(2 * width + 1);
678 fZbins.
Resize(2 * width + 1);
679 for (
int i = 0; i <= width * 2; i++) {
680 fXbins[i] = CF->GetNum()->GetXaxis()->GetBinCenter(middle_x - width + i);
681 fYbins[i] = CF->GetNum()->GetYaxis()->GetBinCenter(middle_y - width + i);
682 fZbins[i] = CF->GetNum()->GetZaxis()->GetBinCenter(middle_z - width + i);
690 for (
int a = 0; a < cf->GetBinsHDX().GetSize(); a++) {
691 Int_t i = cf->GetBinsHDX()[a];
692 Int_t j = cf->GetBinsHDY()[a];
693 Int_t k = cf->GetBinsHDZ()[a];
697 X[0] = cf->EvalHDX(i);
698 X[1] = cf->EvalHDY(j);
699 X[2] = cf->EvalHDZ(k);
701 cf->CFMapHD()[i][j][k] = CF;
707 if (mask ==
nullptr)
return;
708 if (!mask->AreCompatible(
fCF))
return;
714 TF1* CorrFit3DCF::GetDrawableFunc(TString option) {
716 TString className = this->ClassName();
718 if (option.EqualTo(
"x"))
719 func =
new TF1(
"funcX",
this, &CorrFit3DCF::GetFunX,
fRange[0],
fRange[1], nPar, className,
"GetFunX");
720 if (option.EqualTo(
"y"))
721 func =
new TF1(
"funcY",
this, &CorrFit3DCF::GetFunY,
fRange[2],
fRange[3], nPar, className,
"GetFunY");
722 if (option.EqualTo(
"z"))
723 func =
new TF1(
"funcZ",
this, &CorrFit3DCF::GetFunZ,
fRange[4],
fRange[5], nPar, className,
"GetFunZ");
725 if (option.EqualTo(
"xy++"))
726 func =
new TF1(
"funcXY++",
this, &CorrFit3DCF::GetFunXYpp,
fRange[0],
fRange[1], nPar, className,
"GetFunXYpp");
727 if (option.EqualTo(
"xy+-"))
728 func =
new TF1(
"funcXY+-",
this, &CorrFit3DCF::GetFunXYpm,
fRange[0],
fRange[1], nPar, className,
"GetFunXYpm");
729 if (option.EqualTo(
"yz++"))
730 func =
new TF1(
"funcYZ++",
this, &CorrFit3DCF::GetFunYZpp,
fRange[2],
fRange[3], nPar, className,
"GetFunXZpp");
731 if (option.EqualTo(
"yz+-"))
732 func =
new TF1(
"funcYZ+-",
this, &CorrFit3DCF::GetFunYZpm,
fRange[2],
fRange[3], nPar, className,
"GetFunYZpm");
733 if (option.EqualTo(
"xz++"))
734 func =
new TF1(
"funcXZ++",
this, &CorrFit3DCF::GetFunXZpp,
fRange[4],
fRange[5], nPar, className,
"GetFunYZpp");
735 if (option.EqualTo(
"xz+-"))
736 func =
new TF1(
"funcXZ+-",
this, &CorrFit3DCF::GetFunXZpm,
fRange[2],
fRange[3], nPar, className,
"GetFunXZpm");
739 if (option.EqualTo(
"xyz+++"))
740 func =
new TF1(
"funcXYZ+++",
this, &CorrFit3DCF::GetFunXYZppp,
fRange[0],
fRange[1], nPar, className,
"GetFunXYZppp");
741 if (option.EqualTo(
"xyz+-+"))
742 func =
new TF1(
"funcXYZ+-+",
this, &CorrFit3DCF::GetFunXYZpmp,
fRange[0],
fRange[1], nPar, className,
"GetFunXYZpmp");
743 if (option.EqualTo(
"xyz+--"))
744 func =
new TF1(
"funcXYZ+--",
this, &CorrFit3DCF::GetFunXYZpmm,
fRange[0],
fRange[1], nPar, className,
"GetFunXYZpmm");
745 if (option.EqualTo(
"xyz++-"))
746 func =
new TF1(
"funcXYZ++-",
this, &CorrFit3DCF::GetFunXYZppm,
fRange[0],
fRange[1], nPar, className,
"GetFunXYZppm");
748 if (option.EqualTo(
"xy"))
750 "fun2dxy",
this, &CorrFit3DCF::GetFunXY2d,
fRange[0],
fRange[1],
fRange[2],
fRange[3], nPar, className,
"GetFunXY2d");
751 if (option.EqualTo(
"xz"))
753 "fun2dxz",
this, &CorrFit3DCF::GetFunXZ2d,
fRange[0],
fRange[1],
fRange[4],
fRange[5], nPar, className,
"GetFunXZ2d");
754 if (option.EqualTo(
"yz"))
756 "fun2dyz",
this, &CorrFit3DCF::GetFunYZ2d,
fRange[2],
fRange[3],
fRange[4],
fRange[5], nPar, className,
"GetFunYZ2d");
758 if (!func) { std::cout << __LINE__ << __FILE__ <<
" wrong option " << option << std::endl; }
759 SetParametersToTF1(func);
void Resize(Int_t new_dim)
void SetRadiusLimits(Double_t min, Double_t max)
virtual Double_t GetNumericalError(Int_t, Int_t, Int_t) const
virtual Double_t CalculateCF(const Double_t *x, const Double_t *params) const =0
Double_t EvalDenominator(Double_t x, Double_t y, Double_t z) const
virtual void MakePainter(TString options)
virtual Double_t EvalCF(const Double_t *x, const Double_t *params) const
virtual void RecalculateSmoothFunction() const
virtual void EstimateActiveBins()
Double_t Eval(Double_t x, Double_t y, Double_t z)
void SetErrors(TH1 *num, const TH1 *den) const
void SetFittingMask(const CorrFitMask &map)
void SetFuncRange(Double_t x_min, Double_t x_max, Double_t y_min, Double_t y_max, Double_t z_min, Double_t z_max)
void SetRoutLimits(Double_t min, Double_t max)
void SetRsideLimits(Double_t min, Double_t max)
void SetRlongLimits(Double_t min, Double_t max)
TH1 * fDenominatorHistogram
CorrFitPainter * fPainter
TH1 * fCorrelationFunctionHistogram
Femto::EKinematics fKinematics
TH1 * fNumeratorHistogram
Array_1< Double_t > fRange
Double_t GetBinCFVal(Int_t binX, Int_t binY, Int_t binZ, Bool_t extrapolated) const
Int_t GetYBin(Int_t i) const
Int_t GetXBin(Int_t i) const
Int_t GetZBin(Int_t i) const
virtual void SetMask(const CorrFitMask &mask, TH1 *denominator, Bool_t hd)=0
Int_t HDBinToBin(Int_t hd_bin) const
Width_t GetLineWidth() const
Int_t GetParametersNo() const
Double_t GetParameter(Int_t par) const
Double_t * fTempParamsEval
Style_t GetLineStyle() const
Color_t GetLineColor() const
void SetParameterName(Int_t par, TString name)
std::vector< FitParam > fParameters
virtual void SetOption(TString option)