Splines.h

/*
 * HalSplines.h
 *
 *  Created on: 30-04-2022
 *      Author: Daniel Wielanek
 *      E-mail: daniel.wielanek@gmail.com
 *      Warsaw University of Technology, Faculty of Physics
 */
#ifndef HALSPLINES_H_
#define HALSPLINES_H_
 
#include "Array.h"
 
#include <TString.h>
 
 
class TAxis;
class TH1;
class TH1D;
class TH2;
class TH2D;
class TH3;
class TH3D;
 
namespace Hal {
  class Spline1D : public TObject {
    TAxis* fXaxis;
    Int_t fNbins;
    Double_t* fA;   //[fNbins]
    Double_t* fB;   //[fNbins]
    Double_t* fC;   //[fNbins]
    Double_t* fAe;  //[fNbins]
    Double_t* fBe;  //[fNbins]
    void CalcParams(Double_t x1,
                    Double_t y1,
                    Double_t x2,
                    Double_t y2,
                    Double_t x3,
                    Double_t y3,
                    Double_t& A,
                    Double_t& B,
                    Double_t& C);
    void SetFirstBin(TH1* h, Double_t begval, Bool_t use_begval = kFALSE);
    void SetLastBin(TH1* h, Double_t enval, Bool_t use_endval = kFALSE);
 
  public:
    Spline1D(TH1* h = NULL, Double_t begval = 0, Double_t endval = 0);
    void FastOverwrite(TH1D* h, Double_t begval = 0, Double_t endval = 0);
    void Refit();
    Double_t Eval(Double_t x) const;
    Double_t GetError(Double_t x) const;
    Double_t EvalBin(Int_t bin) const;
    Double_t ErrorBin(Int_t bin) const;
    void Eval(Double_t x, Double_t& f, Double_t& error) const;
    virtual TH1D* GetTHD(TString name = "spline_thd") const;
    TAxis* GetXaxis() const { return fXaxis; };
    virtual ~Spline1D();
    ClassDef(Spline1D, 1)
  };
 
  class Spline2D : public TObject {
    TAxis* fXaxis;
    TAxis* fYaxis;
    Int_t fNbinsX;
    Int_t fNbinsY;
    Array_3<Double_t> fA;
    Array_2<Double_t>* fAe;
    void CalcParams(Double_t x[3], Double_t y[3], Double_t z[3][3], Double_t params[9]);
    void Extrapolate(TH2* h, Option_t* extraopolation_opt) const;
    Double_t Extrapolate(TH2* h, Int_t ix1, Int_t ix2, Int_t x, Int_t iy1, Int_t iy2, Int_t y, Int_t opt) const;
    // TODO default constructor & copy constructor
  public:
    Spline2D(TH2* h = NULL, Option_t* interpolation = "");
    void Refit();
    Double_t GetExtrapolationParam(Double_t x, Double_t y, Int_t param) const;
    Double_t Eval(Double_t x, Double_t y) const;
    Double_t GetError(Double_t x, Double_t y) const;
    Double_t EvalBin(Int_t bin_x, Int_t bin_y) const;
    Double_t ErrorBin(Int_t bin_x, Int_t bin_y) const;
    void SetPoint(Int_t xbin, Int_t ybin, Double_t val, Double_t err = 0);
    void PrintParams(Double_t x, Double_t y);
    void Eval(Double_t x, Double_t y, Double_t& f, Double_t& error) const;
    TAxis* GetXaxis() const { return fXaxis; };
    TAxis* GetYaxis() const { return fYaxis; };
    virtual TH2D* GetTHD(TString name = "spline_thd") const;
    virtual ~Spline2D();
    ClassDef(Spline2D, 1)
  };
  class Spline3D : public TObject {
    TAxis* fXaxis;
    TAxis* fYaxis;
    TAxis* fZaxis;
    Int_t fNbinsX;
    Int_t fNbinsY;
    Int_t fNbinsZ;
    Array_4<Double_t> fA;
    Array_3<Double_t> fAe;
    Double_t EvalBin(Int_t ix, Int_t iy, Int_t iz, Double_t x, Double_t y, Double_t z) const;
    void CalcParams(Double_t x[2], Double_t y[2], Double_t z[2], Double_t v[2][2][2], Double_t params[8]);
 
  public:
    Spline3D(TH3* h = NULL, Option_t* interpolation = "");
    void FastOverwrite(TH3* h, Option_t* interpolation = "");
    Double_t Eval(Double_t x, Double_t y, Double_t z) const;
    Double_t GetError(Double_t x, Double_t y, Double_t z) const;
    Double_t EvalBin(Int_t x, Int_t y, Int_t z) const;
    Double_t ErrorBin(Int_t x, Int_t y, Int_t z) const;
    void Refit();
    void PrintParams(Double_t x, Double_t y, Double_t z) const;
    void Eval(Double_t x, Double_t y, Double_t z, Double_t& f, Double_t& error) const;
    TAxis* GetXaxis() const { return fXaxis; };
    TAxis* GetYaxis() const { return fYaxis; };
    TAxis* GetZaxis() const { return fZaxis; };
    virtual TH3D* GetTHD(TString name = "spline_thd") const;
    virtual ~Spline3D();
    ClassDef(Spline3D, 1)
  };
}  // namespace Hal
#endif /* HALSPLINES_H_ */