LitMaterialGrid.h

Go to the documentation of this file.

/* Copyright (C) 2013 GSI/JINR-LIT, Darmstadt/Dubna
   SPDX-License-Identifier: GPL-3.0-only
   Authors: Andrey Lebedev [committer] */
 
#ifndef LITMATERIALGRID_H_
#define LITMATERIALGRID_H_
 
#include "LitTypes.h"
 
#include <iomanip>
#include <sstream>
#include <string>
#include <vector>
 
using std::ostream;
using std::right;
using std::setfill;
using std::setw;
using std::string;
using std::stringstream;
using std::vector;
 
namespace lit
{
  namespace parallel
  {
 
    class LitMaterialGrid {
     public:
      //  fscal GetZ() const { return fZ; }
 
      //  void SetZ(fscal Z) { fZ = Z; }
 
      void SetMaterial(const vector<vector<fscal>>& material, fscal xmin, fscal xmax, fscal ymin, fscal ymax,
                       int nofBinsX, int nofBinsY)
      {
        fMaterial = material;
        fXMin     = xmin;
        fXMax     = xmax;
        fYMin     = ymin;
        fYMax     = ymax;
        fNofBinsX = nofBinsX;
        fNofBinsY = nofBinsY;
        fBinSizeX = ((xmax - xmin) / nofBinsX);
        fBinSizeY = ((ymax - ymin) / nofBinsY);
      }
 
      fscal GetMaterial(fscal x, fscal y) const
      {
        // Check bound conditions and if out of bounds return zero.
        // Can be removed considering performance!
        if (x < fXMin || x > fXMax || y < fYMin || y > fXMax) {
          //   std::cout << "LitMaterialGrid::GetMaterial: out of bounds x=" << x << " y=" << y << std::endl;
          return 0.;
        }
        // Calculate bin indices for X and Y
        short ix = short((x - fXMin) / fBinSizeX);
        short iy = short((y - fYMin) / fBinSizeY);
 
        // Check bound conditions and if out of bounds return zero.
        // Can be removed considering performance!
        if (ix < 0 || iy < 0 || ix >= fNofBinsX - 1 || iy >= fNofBinsY - 1) {
          //  std::cout << "LitMaterialGrid::GetMaterial: out of bounds ix=" << ix << " iy=" << iy << std::endl;
          return 0.;
        }
 
        return fMaterial[ix][iy];
 
 
        // Implementation based on the weighted mean of material thicknesses.
 
 
        //      // Material thickness on the bin nodes
        //      const fscal v1 = fMaterial[ix  ][iy];
        //      const fscal v2 = fMaterial[ix+1][iy];
        //      const fscal v3 = fMaterial[ix  ][iy+1];
        //      const fscal v4 = fMaterial[ix+1][iy+1];
        //      // Calculate weights depending on the distance to the bin nodes
        //      fscal dx1 = (x - ix * fBinSizeX - fXMin);
        //      fscal dx2 = (x - (ix + 1) * fBinSizeX - fXMin);
        //      fscal dy1 = (y - iy * fBinSizeY - fYMin);
        //      fscal dy2 = (y - (iy + 1) * fBinSizeY - fYMin);
        //      fscal w1 = 1./(dx1 * dx1 + dy1 * dy1);
        //      fscal w2 = 1./(dx2 * dx2 + dy1 * dy1);
        //      fscal w3 = 1./(dx1 * dx1 + dy2 * dy2);
        //      fscal w4 = 1./(dx2 * dx2 + dy2 * dy2);
        //      fscal wsum = w1 + w2 + w3 + w4;
        //      if (wsum == 0.) { // Can be removed considering performance!
        //         return 0.;
        //      }
        //      // Calculate output weighted mean material thickness
        //      return (w1 * v1 + w2 * v2 + w3 * v3 + w4 * v4) / wsum;
      }
 
      fvec GetMaterialValue(fvec x, fvec y) const
      {
        fvec v;
        // Get material thickness for each packed value
        for (size_t i = 0; i < fvec::size(); i++) {
          v[i] = GetMaterial(x[i], y[i]);
        }
        return v;
      }
 
      bool IsEmpty() const { return fMaterial.empty(); }
 
      string ToString() const
      {
        stringstream ss;
        ss << "LitMaterialGrid:" /*Z=" << fZ*/ << " Xmin=" << fXMin << " Xmax=" << fXMax << " Ymin=" << fYMin
           << " Ymax=" << fYMax << " nofBinsX=" << fNofBinsX << " nofBinsY=" << fNofBinsY << " binSizeX=" << fBinSizeX
           << " binSizeY=" << fBinSizeY << " material.size=" << fMaterial.size();
        if (fNofBinsX > 0 && fNofBinsY > 0) {
          unsigned int stepX = fNofBinsX / 10;
          unsigned int stepY = fNofBinsY / 10;
          ss << "\nGrid: stepX=" << stepX << " stepY=" << stepY << "\n";
          for (unsigned int i = 0; i < 11; i++) {
            for (unsigned int j = 0; j < 11; j++) {
              if (i < 10 && j < 10) ss << right << setfill(' ') << setw(10) << fMaterial[i * stepX][j * stepY] << " ";
              if (j == 10 && i != 10)
                ss << right << setfill(' ') << setw(10) << fMaterial[i * stepX][fNofBinsY - 1] << " ";
              if (i == 10 && j != 10)
                ss << right << setfill(' ') << setw(10) << fMaterial[fNofBinsX - 1][j * stepY] << " ";
              if (i == 10 && j == 10)
                ss << right << setfill(' ') << setw(10) << fMaterial[fNofBinsX - 1][fNofBinsY - 1] << " ";
            }
            ss << "\n";
          }
        }
        return ss.str();
      }
 
      friend ostream& operator<<(ostream& strm, const LitMaterialGrid& grid)
      {
        strm << grid.ToString();
        return strm;
      }
 
     private:
      fscal fXMin, fXMax;        // Maximum and minimum grid size in X [cm]
      fscal fYMin, fYMax;        // Maximum and minimum grid size in Y [cm]
                                 //   fscal fZ; // Z position of grid slice
      unsigned short fNofBinsX;  // Number of bins along X
      unsigned short fNofBinsY;  // Number of bins along Y
      fscal fBinSizeX;           // Bin size along X [cm]
      fscal fBinSizeY;           // Bin size along Y [cm]
      // Material thickness in bin nodes.
      // Total number of values is
      // (fNofBinsX + 1) * (fNofBinsY + 1)
      vector<vector<fscal>> fMaterial;
    } _fvecalignment;
 
  }  // namespace parallel
}  // namespace lit
 
 
#endif /* LITMATERIALGRID_H_ */