ClusterFactory.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_EMCAL_CLUSTERFACTORY_H_
#define ALICEO2_EMCAL_CLUSTERFACTORY_H_
#include <array>
#include <vector>
#include <optional>
#include <utility>
#include <gsl/span>
#include "Rtypes.h"
#include "fmt/format.h"
#include "DataFormatsEMCAL/Cluster.h"
#include "DataFormatsEMCAL/Digit.h"
#include "DataFormatsEMCAL/Cell.h"
#include "DataFormatsEMCAL/AnalysisCluster.h"
#include "DataFormatsEMCAL/CellLabel.h"
#include "DataFormatsEMCAL/ClusterLabel.h"
#include "EMCALBase/Geometry.h"
#include "MathUtils/Cartesian.h"
 
namespace o2
{
 
namespace emcal
{
 
template <class InputType>
class ClusterFactory
{
 
 public:
  class ClusterRangeException final : public std::exception
  {
   public:
    ClusterRangeException(int clusterIndex, int maxClusters) : std::exception(),
                                                               mClusterID(clusterIndex),
                                                               mMaxClusters(maxClusters),
                                                               mErrorMessage()
    {
      mErrorMessage = fmt::format("Cluster out of range: %d, max %d", mClusterID, mMaxClusters);
    }
 
    ~ClusterRangeException() noexcept final = default;
 
    const char* what() const noexcept final { return mErrorMessage.data(); }
 
    int getClusterID() const { return mClusterID; }
 
    int getMaxNumberOfClusters() const { return mMaxClusters; }
 
   private:
    int mClusterID = 0;        
    int mMaxClusters = 0;      
    std::string mErrorMessage; 
  };
 
  class CellIndexRangeException final : public std::exception
  {
   public:
    CellIndexRangeException(int cellIndex, int maxCellIndex) : std::exception(),
                                                               mCellIndex(cellIndex),
                                                               mMaxCellIndex(maxCellIndex),
                                                               mErrorMessage()
    {
      mErrorMessage = Form("Cell Index out of range: %d, max %d", mCellIndex, mMaxCellIndex);
    }
 
    ~CellIndexRangeException() noexcept final = default;
 
    const char* what() const noexcept final { return mErrorMessage.data(); }
 
    int getCellIndex() const { return mCellIndex; }
 
    int getMaxNumberOfCellIndexs() const { return mMaxCellIndex; }
 
   private:
    int mCellIndex = 0;        
    int mMaxCellIndex = 0;     
    std::string mErrorMessage; 
  };
 
  class GeometryNotSetException final : public std::exception
  {
   public:
    GeometryNotSetException() = default;
    ~GeometryNotSetException() noexcept final = default;
 
    const char* what() const noexcept final { return "Geometry not set"; }
  };
 
  class ClusterIterator
  {
   public:
    ClusterIterator(const ClusterFactory& factory, int clusterIndex, bool forward);
 
    ~ClusterIterator() = default;
 
    bool operator==(const ClusterIterator& rhs) const;
 
    bool operator!=(const ClusterIterator& rhs) const { return !(*this == rhs); }
 
    ClusterIterator& operator++();
 
    ClusterIterator operator++(int);
 
    ClusterIterator& operator--();
 
    ClusterIterator operator--(int);
 
    AnalysisCluster* operator*() { return &mCurrentCluster; }
 
    AnalysisCluster& operator&() { return mCurrentCluster; }
 
    int current_index() const { return mClusterID; }
 
   private:
    const ClusterFactory& mClusterFactory; 
    AnalysisCluster mCurrentCluster;       
    int mClusterID = 0;                    
    bool mForward = true;                  
  };
 
  ClusterFactory() = default;
 
  ClusterFactory(gsl::span<const o2::emcal::Cluster> clustersContainer, gsl::span<const InputType> inputsContainer, gsl::span<const int> cellsIndices);
 
  ClusterFactory(const ClusterFactory& rp) = default;
 
  ClusterFactory& operator=(const ClusterFactory& cf) = default;
 
  ~ClusterFactory() = default;
 
  ClusterIterator begin() const { return ClusterIterator(*this, 0, true); }
 
  ClusterIterator end() const { return ClusterIterator(*this, getNumberOfClusters(), true); }
 
  ClusterIterator rbegin() const { return ClusterIterator(*this, getNumberOfClusters() - 1, false); };
 
  ClusterIterator rend() const { return ClusterIterator(*this, -1, false); };
 
  void reset();
 
  AnalysisCluster buildCluster(int index, o2::emcal::ClusterLabel* clusterLabel = nullptr) const;
 
  void SetECALogWeight(Float_t w) { mLogWeight = w; }
  float GetECALogWeight() const { return mLogWeight; }
 
  void doEvalLocal2tracking(bool justCluster)
  {
    mJustCluster = justCluster;
  }
 
  void evalLocalPosition(gsl::span<const int> inputsIndices, AnalysisCluster& cluster) const;
 
  void evalGlobalPosition(gsl::span<const int> inputsIndices, AnalysisCluster& cluster) const;
 
  void evalLocal2TrackingCSTransform() const;
 
  void evalLocalPositionFit(Double_t deff, Double_t w0, Double_t phiSlope, gsl::span<const int> inputsIndices, AnalysisCluster& cluster) const;
 
  static void getDeffW0(const Double_t esum, Double_t& deff, Double_t& w0);
 
  std::tuple<int, float, float, bool> getMaximalEnergyIndex(gsl::span<const int> inputsIndices) const;
 
  bool isExoticCell(short towerId, float ecell, float const exoticTime) const;
 
  float getECross(short absID, float energy, float const exoticTime) const;
 
  float GetCellWeight(float eCell, float eCluster) const;
 
  int getMultiplicityAtLevel(float level, gsl::span<const int> inputsIndices, AnalysisCluster& clusterAnalysis) const;
 
  int getSuperModuleNumber() const { return mSuperModuleNumber; }
 
  // searches for the local maxima
  // energy above relative level
  // int getNumberOfLocalMax(int nInputMult,
  //                        float locMaxCut, gsl::span<InputType> inputs) const;
 
  // int getNumberOfLocalMax(std::vector<InputType>& maxAt, std::vector<float>& maxAtEnergy,
  //                         float locMaxCut, gsl::span<InputType> inputs) const;
 
  bool sharedCluster() const { return mSharedCluster; }
  void setSharedCluster(bool s) { mSharedCluster = s; }
 
  Double_t tMaxInCm(const Double_t e = 0.0, const int key = 0) const;
 
  bool getLookUpInit() const { return mLookUpInit; }
 
  bool getCoreRadius() const { return mCoreRadius; }
  void setCoreRadius(float radius) { mCoreRadius = radius; }
 
  float getExoticCellFraction() const { return mExoticCellFraction; }
  void setExoticCellFraction(float exoticCellFraction) { mExoticCellFraction = exoticCellFraction; }
 
  float getExoticCellDiffTime() const { return mExoticCellDiffTime; }
  void setExoticCellDiffTime(float exoticCellDiffTime) { mExoticCellDiffTime = exoticCellDiffTime; }
 
  float getExoticCellMinAmplitude() const { return mExoticCellMinAmplitude; }
  void setExoticCellMinAmplitude(float exoticCellMinAmplitude) { mExoticCellMinAmplitude = exoticCellMinAmplitude; }
 
  float getExoticCellInCrossMinAmplitude() const { return mExoticCellInCrossMinAmplitude; }
  void setExoticCellInCrossMinAmplitude(float exoticCellInCrossMinAmplitude) { mExoticCellInCrossMinAmplitude = exoticCellInCrossMinAmplitude; }
 
  bool getUseWeightExotic() const { return mUseWeightExotic; }
  void setUseWeightExotic(float useWeightExotic) { mUseWeightExotic = useWeightExotic; }
 
  void setContainer(gsl::span<const o2::emcal::Cluster> clusterContainer, gsl::span<const InputType> cellContainer, gsl::span<const int> indicesContainer, std::optional<gsl::span<const o2::emcal::CellLabel>> cellLabelContainer = std::nullopt)
  {
    mClustersContainer = clusterContainer;
    mInputsContainer = cellContainer;
    mCellsIndices = indicesContainer;
    if (!getLookUpInit()) {
      setLookUpTable();
    }
    if (cellLabelContainer) {
      mCellLabelContainer = cellLabelContainer.value();
    }
  }
 
  void setLookUpTable(void)
  {
    mLoolUpTowerToIndex.fill(-1);
    for (auto iCellIndex : mCellsIndices) {
      mLoolUpTowerToIndex[mInputsContainer[iCellIndex].getTower()] = iCellIndex;
    }
    mLookUpInit = true;
  }
 
  int getNumberOfClusters() const
  {
    return mClustersContainer.size();
  }
 
  void setGeometry(o2::emcal::Geometry* geometry) { mGeomPtr = geometry; }
 
  class UninitLookUpTableException final : public std::exception
  {
   public:
    UninitLookUpTableException() = default;
 
    ~UninitLookUpTableException() noexcept final = default;
 
    const char* what() const noexcept final { return "Lookup table not initialized, exotics evaluation not possible!"; }
  };
 
 protected:
  void evalCoreEnergy(gsl::span<const int> inputsIndices, AnalysisCluster& clusterAnalysis) const;
 
  void evalDispersion(gsl::span<const int> inputsIndices, AnalysisCluster& clusterAnalysis) const;
 
  void evalElipsAxis(gsl::span<const int> inputsIndices, AnalysisCluster& clusterAnalysis) const;
 
  void evalTime(gsl::span<const int> inputsIndices, AnalysisCluster& clusterAnalysis) const;
 
  float thetaToEta(float arg) const;
 
  float etaToTheta(float arg) const;
 
 private:
  o2::emcal::Geometry* mGeomPtr = nullptr;
 
  float mCoreRadius = 10; 
 
  float mLogWeight = 4.5; 
 
  bool mJustCluster = kFALSE; 
  bool mLookUpInit = false;   
 
  mutable int mSuperModuleNumber = 0;         
  float mDistToBadTower = -1;                 
  bool mSharedCluster = false;                
  float mExoticCellFraction = 0.97;           
  float mExoticCellDiffTime = 1e6;            
  float mExoticCellMinAmplitude = 4.;         
  float mExoticCellInCrossMinAmplitude = 0.1; 
  bool mUseWeightExotic = false;              
 
  gsl::span<const o2::emcal::Cluster> mClustersContainer;    
  gsl::span<const InputType> mInputsContainer;               
  gsl::span<const int> mCellsIndices;                        
  std::array<short, 17664> mLoolUpTowerToIndex;              
  gsl::span<const o2::emcal::CellLabel> mCellLabelContainer; 
 
  ClassDefNV(ClusterFactory, 2);
};
 
} // namespace emcal
} // namespace o2
#endif // ALICEO2_EMCAL_CLUSTERFACTORY_H_