Generic3D.h

Go to the documentation of this file.

/***********************************************************************************\
* (c) Copyright 1998-2019 CERN for the benefit of the LHCb and ATLAS collaborations *
*                                                                                   *
* This software is distributed under the terms of the Apache version 2 licence,     *
* copied verbatim in the file "LICENSE".                                            *
*                                                                                   *
* In applying this licence, 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 GAUDISVC_GENERIC3D_H
#define GAUDISVC_GENERIC3D_H 1

#include "AIDA/IHistogram3D.h"
#include "Annotation.h"
#include "Axis.h"
#include "GaudiKernel/HistogramBase.h"
#include "TFile.h"
#include <memory>
#include <stdexcept>

// Hide warning message:
// warning: 'XYZ' overrides a member function but is not marked 'override'
#ifdef __clang__
#  pragma clang diagnostic push
#  pragma clang diagnostic ignored "-Winconsistent-missing-override"
#elif defined( __GNUC__ ) && __GNUC__ >= 5
#  pragma GCC diagnostic push
#  pragma GCC diagnostic ignored "-Wsuggest-override"
#endif

namespace Gaudi {

  template <typename INTERFACE, typename IMPLEMENTATION>
  class GAUDI_API Generic3D : virtual public INTERFACE, virtual public HistogramBase {
  public:
    typedef Generic3D<INTERFACE, IMPLEMENTATION> Base;
    Generic3D() = default;

  protected:
    Generic3D( IMPLEMENTATION* p ) : m_rep( p ) {}

  public:
    TObject* representation() const override { return m_rep.get(); }
    void adoptRepresentation( TObject* rep ) override;

    int dimension() const override { return 3; }
    std::string title() const override { return m_annotation.value( "Title" ); }
    bool setTitle( const std::string& title ) override;
    std::string name() const { return m_annotation.value( "Name" ); }
    bool setName( const std::string& newName );
    AIDA::IAnnotation& annotation() override { return m_annotation; }
    const AIDA::IAnnotation& annotation() const override { return m_annotation; }

    int entries() const override;
    int allEntries() const override;
    double sumBinHeights() const override;
    double sumAllBinHeights() const override;
    double sumExtraBinHeights() const override { return sumAllBinHeights() - sumBinHeights(); }
    double minBinHeight() const override;
    double maxBinHeight() const override;

    int rIndexX( int index ) const { return m_xAxis.rIndex( index ); }
    int rIndexY( int index ) const { return m_yAxis.rIndex( index ); }
    int rIndexZ( int index ) const { return m_zAxis.rIndex( index ); }

    double binMeanX( int indexX, int, int ) const override {
      return m_rep->GetXaxis()->GetBinCenter( rIndexX( indexX ) );
    }
    double binMeanY( int, int indexY, int ) const override {
      return m_rep->GetYaxis()->GetBinCenter( rIndexY( indexY ) );
    }
    double binMeanZ( int, int, int indexZ ) const override {
      return m_rep->GetYaxis()->GetBinCenter( rIndexY( indexZ ) );
    }
    int binEntries( int indexX, int indexY, int indexZ ) const override {
      if ( binHeight( indexX, indexY, indexZ ) <= 0 ) return 0;
      double xx = binHeight( indexX, indexY, indexZ ) / binError( indexX, indexY, indexZ );
      return int( xx * xx + 0.5 );
    }
    int binEntriesX( int index ) const override {
      int n = 0;
      for ( int i = -2; i < yAxis().bins(); ++i )
        for ( int j = -2; j < zAxis().bins(); ++j ) n += binEntries( index, i, j );
      return n;
    }
    int binEntriesY( int index ) const override {
      int n = 0;
      for ( int i = -2; i < xAxis().bins(); ++i )
        for ( int j = -2; j < zAxis().bins(); ++j ) n += binEntries( i, index, j );
      return n;
    }

    int binEntriesZ( int index ) const override {
      int n = 0;
      for ( int i = -2; i < xAxis().bins(); ++i )
        for ( int j = -2; j < yAxis().bins(); ++j ) n += binEntries( i, j, index );
      return n;
    }

    double binHeight( int indexX, int indexY, int indexZ ) const {
      return m_rep->GetBinContent( rIndexX( indexX ), rIndexY( indexY ), rIndexZ( indexZ ) );
    }

    double binHeightX( int index ) const override {
      double s = 0;
      for ( int i = -2; i < yAxis().bins(); ++i )
        for ( int j = -2; j < zAxis().bins(); ++j ) s += binHeight( index, i, j );
      return s;
    }
    double binHeightY( int index ) const override {
      double s = 0;
      for ( int i = -2; i < xAxis().bins(); ++i )
        for ( int j = -2; j < zAxis().bins(); ++j ) s += binHeight( i, index, j );
      return s;
    }
    double binHeightZ( int index ) const override {
      double s = 0;
      for ( int i = -2; i < xAxis().bins(); ++i )
        for ( int j = -2; j < yAxis().bins(); ++j ) s += binHeight( i, j, index );
      return s;
    }
    double binError( int indexX, int indexY, int indexZ ) const override {
      return m_rep->GetBinError( rIndexX( indexX ), rIndexY( indexY ), rIndexZ( indexZ ) );
    }
    double meanX() const override { return m_rep->GetMean( 1 ); }

    double meanY() const override { return m_rep->GetMean( 2 ); }
    double meanZ() const override { return m_rep->GetMean( 3 ); }
    double rmsX() const override { return m_rep->GetRMS( 1 ); }
    double rmsY() const override { return m_rep->GetRMS( 2 ); }
    double rmsZ() const override { return m_rep->GetRMS( 3 ); }
    const AIDA::IAxis& xAxis() const override { return m_xAxis; }
    const AIDA::IAxis& yAxis() const override { return m_yAxis; }
    const AIDA::IAxis& zAxis() const override { return m_zAxis; }
    int coordToIndexX( double coord ) const override { return xAxis().coordToIndex( coord ); }
    int coordToIndexY( double coord ) const override { return yAxis().coordToIndex( coord ); }
    int coordToIndexZ( double coord ) const override { return zAxis().coordToIndex( coord ); }

    double equivalentBinEntries() const override;
    bool scale( double scaleFactor ) override;
    bool add( const INTERFACE& hist ) override {
      const Base* p = dynamic_cast<const Base*>( &hist );
      if ( !p ) throw std::runtime_error( "Cannot add profile histograms of different implementations." );
      m_rep->Add( p->m_rep.get() );
      return true;
    }

    // overwrite extraentries
    int extraEntries() const override {
      return binEntries( AIDA::IAxis::UNDERFLOW_BIN, AIDA::IAxis::UNDERFLOW_BIN, AIDA::IAxis::UNDERFLOW_BIN ) +
             binEntries( AIDA::IAxis::UNDERFLOW_BIN, AIDA::IAxis::UNDERFLOW_BIN, AIDA::IAxis::OVERFLOW_BIN ) +
             binEntries( AIDA::IAxis::UNDERFLOW_BIN, AIDA::IAxis::OVERFLOW_BIN, AIDA::IAxis::UNDERFLOW_BIN ) +
             binEntries( AIDA::IAxis::OVERFLOW_BIN, AIDA::IAxis::UNDERFLOW_BIN, AIDA::IAxis::UNDERFLOW_BIN ) +
             binEntries( AIDA::IAxis::OVERFLOW_BIN, AIDA::IAxis::UNDERFLOW_BIN, AIDA::IAxis::OVERFLOW_BIN ) +
             binEntries( AIDA::IAxis::OVERFLOW_BIN, AIDA::IAxis::OVERFLOW_BIN, AIDA::IAxis::UNDERFLOW_BIN ) +
             binEntries( AIDA::IAxis::OVERFLOW_BIN, AIDA::IAxis::OVERFLOW_BIN, AIDA::IAxis::OVERFLOW_BIN );
    }
    std::ostream& print( std::ostream& s ) const override;
    std::ostream& write( std::ostream& s ) const override;
    int write( const char* file_name ) const override;

  protected:
    Gaudi::Axis m_xAxis;
    Gaudi::Axis m_yAxis;
    Gaudi::Axis m_zAxis;
    mutable AIDA::Annotation m_annotation;
    std::unique_ptr<IMPLEMENTATION> m_rep;
    // class type
    std::string m_classType;
    // cache sumEntries (allEntries)   when setting contents since Root can't compute by himself
    int m_sumEntries = 0;
  }; // end class IHistogram3D

  template <class INTERFACE, class IMPLEMENTATION>
  bool Generic3D<INTERFACE, IMPLEMENTATION>::setTitle( const std::string& title ) {
    m_rep->SetTitle( title.c_str() );
    if ( !annotation().addItem( "Title", title ) ) m_annotation.setValue( "Title", title );
    if ( !annotation().addItem( "title", title ) ) annotation().setValue( "title", title );
    return true;
  }

  template <class INTERFACE, class IMPLEMENTATION>
  bool Generic3D<INTERFACE, IMPLEMENTATION>::setName( const std::string& newName ) {
    m_rep->SetName( newName.c_str() );
    m_annotation.setValue( "Name", newName );
    return true;
  }
  template <class INTERFACE, class IMPLEMENTATION>
  int Generic3D<INTERFACE, IMPLEMENTATION>::entries() const {
    return m_rep->GetEntries();
  }

  template <class INTERFACE, class IMPLEMENTATION>
  int Generic3D<INTERFACE, IMPLEMENTATION>::allEntries() const {
    return int( m_rep->GetEntries() );
  }

  template <class INTERFACE, class IMPLEMENTATION>
  double Generic3D<INTERFACE, IMPLEMENTATION>::minBinHeight() const {
    return m_rep->GetMinimum();
  }

  template <class INTERFACE, class IMPLEMENTATION>
  double Generic3D<INTERFACE, IMPLEMENTATION>::maxBinHeight() const {
    return m_rep->GetMaximum();
  }

  template <class INTERFACE, class IMPLEMENTATION>
  double Generic3D<INTERFACE, IMPLEMENTATION>::sumBinHeights() const {
    return m_rep->GetSumOfWeights();
  }

  template <class INTERFACE, class IMPLEMENTATION>
  double Generic3D<INTERFACE, IMPLEMENTATION>::sumAllBinHeights() const {
    return m_rep->GetSum();
  }

  template <class INTERFACE, class IMPLEMENTATION>
  double Generic3D<INTERFACE, IMPLEMENTATION>::equivalentBinEntries() const {
    if ( sumBinHeights() <= 0 ) return 0;
    Stat_t stats[11]; // cover up to 3D...
    m_rep->GetStats( stats );
    return stats[0] * stats[0] / stats[1];
  }

  template <class INTERFACE, class IMPLEMENTATION>
  bool Generic3D<INTERFACE, IMPLEMENTATION>::scale( double scaleFactor ) {
    m_rep->Scale( scaleFactor );
    return true;
  }

  template <class INTERFACE, class IMPLEMENTATION>
  std::ostream& Generic3D<INTERFACE, IMPLEMENTATION>::print( std::ostream& s ) const {
    m_rep->Print( "all" );
    return s;
  }

  template <class INTERFACE, class IMPLEMENTATION>
  std::ostream& Generic3D<INTERFACE, IMPLEMENTATION>::write( std::ostream& s ) const {
    s << "\n3D Histogram Table: " << std::endl;
    s << "BinX, BinY, BinZ, Height, Error " << std::endl;
    for ( int i = 0; i < xAxis().bins(); ++i )
      for ( int j = 0; j < yAxis().bins(); ++j )
        for ( int k = 0; k < zAxis().bins(); ++k )
          s << binMeanX( i, j, k ) << ", " << binMeanY( i, j, k ) << ", " << binMeanZ( i, j, k ) << ", "
            << binHeight( i, j, k ) << ", " << binError( i, j, k ) << std::endl;
    s << std::endl;
    return s;
  }

  template <class INTERFACE, class IMPLEMENTATION>
  int Generic3D<INTERFACE, IMPLEMENTATION>::write( const char* file_name ) const {
    TFile* f      = TFile::Open( file_name, "RECREATE" );
    Int_t  nbytes = m_rep->Write();
    f->Close();
    return nbytes;
  }
} // namespace Gaudi

#ifdef __clang__
#  pragma clang diagnostic pop
#elif defined( __GNUC__ ) && __GNUC__ >= 5
#  pragma GCC diagnostic pop
#endif

#endif // GAUDIPI_GENERIC3D_H