TriCubic.h

Go to the documentation of this file.

// Copyright 2019-2020 CERN and copyright holders of ALICE O2.
// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders.
// All rights not expressly granted are reserved.
//
// This software is distributed under the terms of the GNU General Public
// License v3 (GPL Version 3), copied verbatim in the file "COPYING".
//
// In applying this license CERN does not waive the privileges and immunities
// granted to it by virtue of its status as an Intergovernmental Organization
// or submit itself to any jurisdiction.
 
 
#ifndef ALICEO2_TPC_TRICUBIC_H_
#define ALICEO2_TPC_TRICUBIC_H_
 
#include "TPCSpaceCharge/RegularGrid3D.h"
 
// forward declare VC Memory
template <typename DataT, size_t, size_t, bool>
class Memory;
 
namespace o2
{
namespace tpc
{
 
template <typename DataT>
class DataContainer3D;
 
 
template <typename DataT, size_t N>
class Vector;
 
template <typename DataT = double>
class TriCubicInterpolator
{
  using Grid3D = RegularGrid3D<DataT>;
  using DataContainer = DataContainer3D<DataT>;
 
 public:
  TriCubicInterpolator(const DataContainer& gridData, const Grid3D& gridProperties) : mGridData{&gridData}, mGridProperties{&gridProperties} {};
 
  TriCubicInterpolator(TriCubicInterpolator<DataT>&&);
 
  TriCubicInterpolator<DataT>& operator=(TriCubicInterpolator<DataT>&&);
 
  enum class ExtrapolationType {
    Linear = 0,   
    Parabola = 1, 
  };
 
  // interpolate value at given coordinate
  DataT operator()(const DataT z, const DataT r, const DataT phi) const { return interpolateSparse(z, r, phi); }
 
  void setExtrapolationType(const ExtrapolationType extrapolationType) { mExtrapolationType = extrapolationType; }
 
  ExtrapolationType getExtrapolationType() const { return mExtrapolationType; }
 
  void setExtrapolateValues(const bool extraPolate) { mExtraPolateValues = extraPolate; }
 
  bool getExtrapolateValues() const { return mExtraPolateValues; }
 
 private:
  static constexpr unsigned int FDim = Grid3D::getDim();              
  static constexpr unsigned int FZ = Grid3D::getFZ();                 
  static constexpr unsigned int FR = Grid3D::getFR();                 
  static constexpr unsigned int FPHI = Grid3D::getFPhi();             
  const DataContainer* mGridData{};                                   
  const Grid3D* mGridProperties{};                                    
  ExtrapolationType mExtrapolationType = ExtrapolationType::Parabola; 
  bool mExtraPolateValues = true;                                     
 
  //                 DEFINITION OF enum GridPos
  //========================================================
  //              r
  //              |             6------F---7
  //              |           / |        / |
  //              |         K   G YR   L   H
  //              |       /     |    /     |
  //              |      2---B------3      |
  //              |      |      |   |      |
  //              |      |      4---|---E--5
  //              |      C XL  /    D XR  /
  //              |      |   I  YL  |   J
  //              |      | /        | /
  //              |      0---A------1
  //              |------------------------------- z
  //            /
  //          /
  //        /
  //      phi
  //========================================================
 
  enum class GridPos {
    None = 27,
    InnerVolume = 26,
    Edge0 = 0,
    Edge1 = 1,
    Edge2 = 2,
    Edge3 = 3,
    Edge4 = 4,
    Edge5 = 5,
    Edge6 = 6,
    Edge7 = 7,
    LineA = 8,
    LineB = 9,
    LineC = 10,
    LineD = 11,
    LineE = 12,
    LineF = 13,
    LineG = 14,
    LineH = 15,
    LineI = 16,
    LineJ = 17,
    LineK = 18,
    LineL = 19,
    SideXRight = 20,
    SideXLeft = 21,
    SideYRight = 22,
    SideYLeft = 23,
    SideZRight = 24,
    SideZLeft = 25
  };
 
  void setValues(const int iz, const int ir, const int iphi, std::array<Vector<DataT, 4>, 16>& cVals) const;
 
  // interpolate value at given coordinate - this method doesnt compute and stores the coefficients and is faster when quering only a few values per cube
  DataT interpolateSparse(const DataT z, const DataT r, const DataT phi) const;
 
  // for periodic boundary conditions
  void getDataIndexCircularArray(const int index0, const int dim, int arr[]) const;
 
  // this helps to get circular and non circular padding indices
  int getRegulatedDelta(const int index0, const int delta, const unsigned int dim, const int offs) const { return mGridProperties->isIndexInGrid(index0 + delta, dim) ? delta : offs; }
 
  DataT extrapolation(const DataT valk, const DataT valk1, const DataT valk2) const;
 
  DataT linearExtrapolation(const DataT valk, const DataT valk1) const;
 
  DataT parabolExtrapolation(const DataT valk, const DataT valk1, const DataT valk2) const;
 
  GridPos findPos(const int iz, const int ir, const int iphi) const;
 
  bool isInInnerVolume(const int iz, const int ir, const int iphi, GridPos& posType) const;
 
  bool findEdge(const int iz, const int ir, const int iphi, GridPos& posType) const;
 
  bool findLine(const int iz, const int ir, const int iphi, GridPos& posType) const;
 
  bool findSide(const int iz, const int ir, const int iphi, GridPos& posType) const;
 
  bool isSideRight(const int ind, const int dim) const;
 
  bool isSideLeft(const int ind) const;
};
 
} // namespace tpc
} // namespace o2
 
#endif