include/GaudiCommonSvc/Generic3D.h

0001 /***********************************************************************************\
0002 * (c) Copyright 1998-2019 CERN for the benefit of the LHCb and ATLAS collaborations *
0003 *                                                                                   *
0004 * This software is distributed under the terms of the Apache version 2 licence,     *
0005 * copied verbatim in the file "LICENSE".                                            *
0006 *                                                                                   *
0007 * In applying this licence, CERN does not waive the privileges and immunities       *
0008 * granted to it by virtue of its status as an Intergovernmental Organization        *
0009 * or submit itself to any jurisdiction.                                             *
0010 \***********************************************************************************/
0011 #ifndef GAUDISVC_GENERIC3D_H
0012 #define GAUDISVC_GENERIC3D_H 1
0013
0014 #include "AIDA/IHistogram3D.h"
0015 #include "Annotation.h"
0016 #include "Axis.h"
0017 #include "GaudiKernel/HistogramBase.h"
0018 #include "TFile.h"
0019 #include <memory>
0020 #include <stdexcept>
0021
0022 // Hide warning message:
0023 // warning: 'XYZ' overrides a member function but is not marked 'override'
0024 #ifdef __clang__
0025 #  pragma clang diagnostic push
0026 #  pragma clang diagnostic ignored "-Wsuggest-override"
0027 #  pragma clang diagnostic ignored "-Winconsistent-missing-override"
0028 #elif defined( __GNUC__ ) && __GNUC__ >= 5
0029 #  pragma GCC diagnostic push
0030 #  pragma GCC diagnostic ignored "-Wsuggest-override"
0031 #endif
0032
0033 namespace Gaudi {
0034
0035   /** @class Generic3D Generic3D.h GaudiPI/Generic3D.h
0036    *
0037    * Common AIDA implementation stuff for histograms and profiles
0038    * using ROOT implementations
0039    *
0040    * Credits: This code is the result of some stripdown implementation
0041    * of LCG/PI. Credits to them!
0042    *
0043    *  @author M.Frank
0044    */
0045   template <typename INTERFACE, typename IMPLEMENTATION>
0046   class GAUDI_API Generic3D : virtual public INTERFACE, virtual public HistogramBase {
0047   public:
0048     typedef Generic3D<INTERFACE, IMPLEMENTATION> Base;
0049     /// Default constructor
0050     Generic3D() = default;
0051
0052   protected:
0053     /// constructor
0054     Generic3D( IMPLEMENTATION* p ) : m_rep( p ) {}
0055
0056   public:
0057     /// ROOT object implementation
0058     TObject* representation() const override { return m_rep.get(); }
0059     /// Adopt ROOT histogram representation
0060     void adoptRepresentation( TObject* rep ) override;
0061
0062     /// Get the Histogram's dimension.
0063     int dimension() const override { return 3; }
0064     /// Get the title of the object
0065     std::string title() const override { return m_annotation.value( "Title" ); }
0066     /// Set the title of the object
0067     bool setTitle( const std::string& title ) override;
0068     /// object name
0069     std::string name() const { return m_annotation.value( "Name" ); }
0070     /// Sets the name of the object
0071     bool setName( const std::string& newName );
0072     /// Access annotation object
0073     AIDA::IAnnotation& annotation() override { return m_annotation; }
0074     /// Access annotation object (cons)
0075     const AIDA::IAnnotation& annotation() const override { return m_annotation; }
0076
0077     /// Get the number or all the entries
0078     int entries() const override;
0079     /// Get the number or all the entries, both in range and underflow/overflow bins of the IProfile.
0080     int allEntries() const override;
0081     /// Get the sum of in range bin heights in the IProfile.
0082     double sumBinHeights() const override;
0083     /// Get the sum of all the bins heights (including underflow and overflow bin).
0084     double sumAllBinHeights() const override;
0085     /// Get the sum of the underflow and overflow bin height.
0086     double sumExtraBinHeights() const override { return sumAllBinHeights() - sumBinHeights(); }
0087     /// Get the minimum height of the in-range bins.
0088     double minBinHeight() const override;
0089     /// Get the maximum height of the in-range bins.
0090     double maxBinHeight() const override;
0091
0092     int rIndexX( int index ) const { return m_xAxis.rIndex( index ); }
0093     int rIndexY( int index ) const { return m_yAxis.rIndex( index ); }
0094     int rIndexZ( int index ) const { return m_zAxis.rIndex( index ); }
0095
0096     /// The weighted mean along the x axis of a given bin.
0097     double binMeanX( int indexX, int, int ) const override {
0098       return m_rep->GetXaxis()->GetBinCenter( rIndexX( indexX ) );
0099     }
0100     /// The weighted mean along the y axis of a given bin.
0101     double binMeanY( int, int indexY, int ) const override {
0102       return m_rep->GetYaxis()->GetBinCenter( rIndexY( indexY ) );
0103     }
0104     /// The weighted mean along the z axis of a given bin.
0105     double binMeanZ( int, int, int indexZ ) const override {
0106       return m_rep->GetYaxis()->GetBinCenter( rIndexY( indexZ ) );
0107     }
0108     /// Number of entries in the corresponding bin (ie the number of times fill was calle d for this bin).
0109     int binEntries( int indexX, int indexY, int indexZ ) const override {
0110       if ( binHeight( indexX, indexY, indexZ ) <= 0 ) return 0;
0111       double xx = binHeight( indexX, indexY, indexZ ) / binError( indexX, indexY, indexZ );
0112       return int( xx * xx + 0.5 );
0113     }
0114     /// Sum of all the entries of the bins along a given x bin.
0115     int binEntriesX( int index ) const override {
0116       int n = 0;
0117       for ( int i = -2; i < yAxis().bins(); ++i )
0118         for ( int j = -2; j < zAxis().bins(); ++j ) n += binEntries( index, i, j );
0119       return n;
0120     }
0121     /// Sum of all the entries of the bins along a given y bin.
0122     int binEntriesY( int index ) const override {
0123       int n = 0;
0124       for ( int i = -2; i < xAxis().bins(); ++i )
0125         for ( int j = -2; j < zAxis().bins(); ++j ) n += binEntries( i, index, j );
0126       return n;
0127     }
0128
0129     /// Sum of all the entries of the bins along a given z bin.
0130     int binEntriesZ( int index ) const override {
0131       int n = 0;
0132       for ( int i = -2; i < xAxis().bins(); ++i )
0133         for ( int j = -2; j < yAxis().bins(); ++j ) n += binEntries( i, j, index );
0134       return n;
0135     }
0136
0137     /// Total height of the corresponding bin (ie the sum of the weights in this bin).
0138     double binHeight( int indexX, int indexY, int indexZ ) const {
0139       return m_rep->GetBinContent( rIndexX( indexX ), rIndexY( indexY ), rIndexZ( indexZ ) );
0140     }
0141
0142     /// Sum of all the heights of the bins along a given x bin.
0143     double binHeightX( int index ) const override {
0144       double s = 0;
0145       for ( int i = -2; i < yAxis().bins(); ++i )
0146         for ( int j = -2; j < zAxis().bins(); ++j ) s += binHeight( index, i, j );
0147       return s;
0148     }
0149     /// Sum of all the heights of the bins along a given y bin.
0150     double binHeightY( int index ) const override {
0151       double s = 0;
0152       for ( int i = -2; i < xAxis().bins(); ++i )
0153         for ( int j = -2; j < zAxis().bins(); ++j ) s += binHeight( i, index, j );
0154       return s;
0155     }
0156     /// Sum of all the heights of the bins along a given z bin.
0157     double binHeightZ( int index ) const override {
0158       double s = 0;
0159       for ( int i = -2; i < xAxis().bins(); ++i )
0160         for ( int j = -2; j < yAxis().bins(); ++j ) s += binHeight( i, j, index );
0161       return s;
0162     }
0163     /// The error of a given bin.
0164     double binError( int indexX, int indexY, int indexZ ) const override {
0165       return m_rep->GetBinError( rIndexX( indexX ), rIndexY( indexY ), rIndexZ( indexZ ) );
0166     }
0167     /// The mean of the IHistogram3D along the x axis.
0168     double meanX() const override { return m_rep->GetMean( 1 ); }
0169
0170     /// The mean of the IHistogram3D along the y axis.
0171     double meanY() const override { return m_rep->GetMean( 2 ); }
0172     /// The mean of the IHistogram3D along the z axis.
0173     double meanZ() const override { return m_rep->GetMean( 3 ); }
0174     /// The RMS of the IHistogram3D along the x axis.
0175     double rmsX() const override { return m_rep->GetRMS( 1 ); }
0176     /// The RMS of the IHistogram3D along the y axis.
0177     double rmsY() const override { return m_rep->GetRMS( 2 ); }
0178     /// The RMS of the IHistogram3D along the z axis.
0179     double rmsZ() const override { return m_rep->GetRMS( 3 ); }
0180     /// Get the x axis of the IHistogram3D.
0181     const AIDA::IAxis& xAxis() const override { return m_xAxis; }
0182     /// Get the y axis of the IHistogram3D.
0183     const AIDA::IAxis& yAxis() const override { return m_yAxis; }
0184     /// Get the z axis of the IHistogram3D.
0185     const AIDA::IAxis& zAxis() const override { return m_zAxis; }
0186     /// Get the bin number corresponding to a given coordinate along the x axis.
0187     int coordToIndexX( double coord ) const override { return xAxis().coordToIndex( coord ); }
0188     /// Get the bin number corresponding to a given coordinate along the y axis.
0189     int coordToIndexY( double coord ) const override { return yAxis().coordToIndex( coord ); }
0190     /// Get the bin number corresponding to a given coordinate along the z axis.
0191     int coordToIndexZ( double coord ) const override { return zAxis().coordToIndex( coord ); }
0192
0193     /// Number of equivalent entries, i.e. <tt>SUM[ weight ] ^ 2 / SUM[ weight^2 ]</tt>
0194     double equivalentBinEntries() const override;
0195     /// Scale the weights and the errors of all the IHistogram's bins (in-range and out-of-range ones) by a given scale
0196     /// factor.
0197     bool scale( double scaleFactor ) override;
0198     /// Add to this Histogram3D the contents of another IHistogram3D.
0199     bool add( const INTERFACE& hist ) override {
0200       const Base* p = dynamic_cast<const Base*>( &hist );
0201       if ( !p ) throw std::runtime_error( "Cannot add profile histograms of different implementations." );
0202       m_rep->Add( p->m_rep.get() );
0203       return true;
0204     }
0205
0206     // overwrite extraentries
0207     int extraEntries() const override {
0208       return binEntries( AIDA::IAxis::UNDERFLOW_BIN, AIDA::IAxis::UNDERFLOW_BIN, AIDA::IAxis::UNDERFLOW_BIN ) +
0209              binEntries( AIDA::IAxis::UNDERFLOW_BIN, AIDA::IAxis::UNDERFLOW_BIN, AIDA::IAxis::OVERFLOW_BIN ) +
0210              binEntries( AIDA::IAxis::UNDERFLOW_BIN, AIDA::IAxis::OVERFLOW_BIN, AIDA::IAxis::UNDERFLOW_BIN ) +
0211              binEntries( AIDA::IAxis::OVERFLOW_BIN, AIDA::IAxis::UNDERFLOW_BIN, AIDA::IAxis::UNDERFLOW_BIN ) +
0212              binEntries( AIDA::IAxis::OVERFLOW_BIN, AIDA::IAxis::UNDERFLOW_BIN, AIDA::IAxis::OVERFLOW_BIN ) +
0213              binEntries( AIDA::IAxis::OVERFLOW_BIN, AIDA::IAxis::OVERFLOW_BIN, AIDA::IAxis::UNDERFLOW_BIN ) +
0214              binEntries( AIDA::IAxis::OVERFLOW_BIN, AIDA::IAxis::OVERFLOW_BIN, AIDA::IAxis::OVERFLOW_BIN );
0215     }
0216     /// Print (ASCII) the histogram into the output stream
0217     std::ostream& print( std::ostream& s ) const override;
0218     /// Write (ASCII) the histogram table into the output stream
0219     std::ostream& write( std::ostream& s ) const override;
0220     /// Write (ASCII) the histogram table into a file
0221     int write( const char* file_name ) const override;
0222
0223   protected:
0224     Gaudi::Axis m_xAxis;
0225     Gaudi::Axis m_yAxis;
0226     Gaudi::Axis m_zAxis;
0227     /// Object annotations
0228     mutable AIDA::Annotation m_annotation;
0229     /// Reference to underlying implementation
0230     std::unique_ptr<IMPLEMENTATION> m_rep;
0231     // class type
0232     std::string m_classType;
0233     // cache sumEntries (allEntries)   when setting contents since Root can't compute by himself
0234     int m_sumEntries = 0;
0235   }; // end class IHistogram3D
0236
0237   template <class INTERFACE, class IMPLEMENTATION>
0238   bool Generic3D<INTERFACE, IMPLEMENTATION>::setTitle( const std::string& title ) {
0239     m_rep->SetTitle( title.c_str() );
0240     if ( !annotation().addItem( "Title", title ) ) m_annotation.setValue( "Title", title );
0241     if ( !annotation().addItem( "title", title ) ) annotation().setValue( "title", title );
0242     return true;
0243   }
0244
0245   template <class INTERFACE, class IMPLEMENTATION>
0246   bool Generic3D<INTERFACE, IMPLEMENTATION>::setName( const std::string& newName ) {
0247     m_rep->SetName( newName.c_str() );
0248     m_annotation.setValue( "Name", newName );
0249     return true;
0250   }
0251   template <class INTERFACE, class IMPLEMENTATION>
0252   int Generic3D<INTERFACE, IMPLEMENTATION>::entries() const {
0253     return m_rep->GetEntries();
0254   }
0255
0256   template <class INTERFACE, class IMPLEMENTATION>
0257   int Generic3D<INTERFACE, IMPLEMENTATION>::allEntries() const {
0258     return int( m_rep->GetEntries() );
0259   }
0260
0261   template <class INTERFACE, class IMPLEMENTATION>
0262   double Generic3D<INTERFACE, IMPLEMENTATION>::minBinHeight() const {
0263     return m_rep->GetMinimum();
0264   }
0265
0266   template <class INTERFACE, class IMPLEMENTATION>
0267   double Generic3D<INTERFACE, IMPLEMENTATION>::maxBinHeight() const {
0268     return m_rep->GetMaximum();
0269   }
0270
0271   template <class INTERFACE, class IMPLEMENTATION>
0272   double Generic3D<INTERFACE, IMPLEMENTATION>::sumBinHeights() const {
0273     return m_rep->GetSumOfWeights();
0274   }
0275
0276   template <class INTERFACE, class IMPLEMENTATION>
0277   double Generic3D<INTERFACE, IMPLEMENTATION>::sumAllBinHeights() const {
0278     return m_rep->GetSum();
0279   }
0280
0281   template <class INTERFACE, class IMPLEMENTATION>
0282   double Generic3D<INTERFACE, IMPLEMENTATION>::equivalentBinEntries() const {
0283     if ( sumBinHeights() <= 0 ) return 0;
0284     Stat_t stats[11]; // cover up to 3D...
0285     m_rep->GetStats( stats );
0286     return stats[0] * stats[0] / stats[1];
0287   }
0288
0289   template <class INTERFACE, class IMPLEMENTATION>
0290   bool Generic3D<INTERFACE, IMPLEMENTATION>::scale( double scaleFactor ) {
0291     m_rep->Scale( scaleFactor );
0292     return true;
0293   }
0294
0295   template <class INTERFACE, class IMPLEMENTATION>
0296   std::ostream& Generic3D<INTERFACE, IMPLEMENTATION>::print( std::ostream& s ) const {
0297     /// bin contents and errors are printed for all bins including under and overflows
0298     m_rep->Print( "all" );
0299     return s;
0300   }
0301
0302   /// Write (ASCII) the histogram table into the output stream
0303   template <class INTERFACE, class IMPLEMENTATION>
0304   std::ostream& Generic3D<INTERFACE, IMPLEMENTATION>::write( std::ostream& s ) const {
0305     s << "\n3D Histogram Table: " << std::endl;
0306     s << "BinX, BinY, BinZ, Height, Error " << std::endl;
0307     for ( int i = 0; i < xAxis().bins(); ++i )
0308       for ( int j = 0; j < yAxis().bins(); ++j )
0309         for ( int k = 0; k < zAxis().bins(); ++k )
0310           s << binMeanX( i, j, k ) << ", " << binMeanY( i, j, k ) << ", " << binMeanZ( i, j, k ) << ", "
0311             << binHeight( i, j, k ) << ", " << binError( i, j, k ) << std::endl;
0312     s << std::endl;
0313     return s;
0314   }
0315
0316   /// Write (ASCII) the histogram table into a file
0317   template <class INTERFACE, class IMPLEMENTATION>
0318   int Generic3D<INTERFACE, IMPLEMENTATION>::write( const char* file_name ) const {
0319     TFile* f      = TFile::Open( file_name, "RECREATE" );
0320     Int_t  nbytes = m_rep->Write();
0321     f->Close();
0322     return nbytes;
0323   }
0324 } // namespace Gaudi
0325
0326 #ifdef __clang__
0327 #  pragma clang diagnostic pop
0328 #elif defined( __GNUC__ ) && __GNUC__ >= 5
0329 #  pragma GCC diagnostic pop
0330 #endif
0331
0332 #endif // GAUDIPI_GENERIC3D_H