KfMaterialMap.h

Go to the documentation of this file.

/* Copyright (C) 2024 GSI Helmholtzzentrum fuer Schwerionenforschung, Darmstadt
   SPDX-License-Identifier: GPL-3.0-only
   Authors: Igor Kulakov, Sergey Gorbunov, Andrey Lebedev, Sergei Zharko [committer] */
 
#pragma once  // include this header only once per compilation unit
 
// NOTE: No dependency from CaCore is allowed
//#include "CaSimd.h"
//#include "CaUtils.h"
#include "AlgoFairloggerCompat.h"
#include "KfDefs.h"
#include "KfUtils.h"
 
#include <boost/serialization/vector.hpp>
 
#include <iomanip>
#include <sstream>
#include <string>
#include <type_traits>
#include <vector>
 
//TODO: rewrite
namespace cbm::algo::kf
{
  class alignas(VcMemAlign) MaterialMap {
   public:
    MaterialMap() = default;

    MaterialMap(int nBins, float xyMax, float zRef, float zMin, float zMax);

    MaterialMap(const MaterialMap& other) = default;

    MaterialMap& operator=(const MaterialMap& other) = default;

    MaterialMap(MaterialMap&& other) noexcept = default;

    MaterialMap& operator=(MaterialMap&& other) noexcept = default;

    ~MaterialMap() noexcept = default;

    void Add(const MaterialMap& other, float zTarg = defs::Undef<float>);

    int GetBin(float x, float y) const;

    int GetNbins() const { return fNbins; }

    float GetXYmax() const { return fXYmax; }

    float GetZref() const { return fZref; }

    float GetZmin() const { return fZmin; }

    float GetZmax() const { return fZmax; }

    template<typename I>
    I GetThicknessX0(const I& x, const I& y) const
    {
      if constexpr (std::is_same_v<I, fvec>) {
        fvec res;
        for (size_t i = 0; i < utils::simd::Size<I>(); ++i) {
          res[i] = GetThicknessX0(x[i], y[i]);
        }
        return res;
      }
      else {
        I xNew = (x < fXYmax && x >= -fXYmax) ? x : 0;
        I yNew = (y < fXYmax && y >= -fXYmax) ? y : 0;
        int i  = static_cast<int>((xNew + fXYmax) * fFactor);
        int j  = static_cast<int>((yNew + fXYmax) * fFactor);
        i      = (i < fNbins && i >= 0) ? i : fNbins / 2;
        j      = (j < fNbins && j >= 0) ? j : fNbins / 2;
        return fTable[i + j * fNbins];
      }
    }

    template<typename I>
    I GetBinThicknessX0(int iGlob) const
    {
      if constexpr (std::is_same_v<I, fvec>) {
        fvec res;
        for (size_t i = 0; i < utils::simd::Size<I>(); ++i) {
          res[i] = GetBinThicknessX0<fscal>(iGlob);
        }
        return res;
      }
      else {
        return fTable[iGlob];
      }
    }

    template<typename I>
    I GetBinThicknessX0(int iX, int iY) const
    {
      return GetBinThicknessX0<I>(iX + iY * fNbins);
    }

    bool IsUndefined() const
    {
      return utils::IsUndefined(fNbins) || utils::IsUndefined(fXYmax * fFactor * fZref * fZmin * fZmax);
    }

    void Rebin(int factor);

    void SetRadThickBin(int iBinX, int iBinY, float thickness) { fTable[iBinX + fNbins * iBinY] = thickness; }

    std::string ToString(int indentLevel = 0, int verbose = 1) const;

    friend bool operator<(const MaterialMap& lhs, const MaterialMap& rhs) { return lhs.fZref < rhs.fZref; }
 
   private:
    void CheckConsistency() const;
 
 
    int fNbins    = defs::Undef<int>;    
    float fXYmax  = defs::Undef<float>;  
    float fFactor = defs::Undef<float>;  
    float fZref   = defs::Undef<float>;  
    float fZmin   = defs::Undef<float>;  
    float fZmax   = defs::Undef<float>;  
    std::vector<float> fTable{};         

    friend class boost::serialization::access;
    template<class Archive>
    void serialize(Archive& ar, const unsigned int)
    {
      ar& fNbins;
      ar& fXYmax;
      ar& fFactor;
      ar& fZref;
      ar& fZmin;
      ar& fZmax;
      ar& fTable;
    }
  };
}  // namespace cbm::algo::kf