Skip to content
Snippets Groups Projects
Commit e0c779c9 authored by Sebastien Ponce's avatar Sebastien Ponce
Browse files

Fixed float comparisons in histograms

parent a6483977
No related branches found
No related tags found
1 merge request!1490Fixed unsafe floating point comparisons
...@@ -125,22 +125,17 @@ bool Gaudi::Histogram1D::setBinContents( int i, int entries, double height, doub ...@@ -125,22 +125,17 @@ bool Gaudi::Histogram1D::setBinContents( int i, int entries, double height, doub
return true; return true;
} }
#ifdef __ICC
// disable icc remark #1572: floating-point equality and inequality comparisons are unreliable
// The comparison are meant
# pragma warning( push )
# pragma warning( disable : 1572 )
#endif
bool Gaudi::Histogram1D::setRms( double rms ) { bool Gaudi::Histogram1D::setRms( double rms ) {
m_rep->SetEntries( m_sumEntries ); m_rep->SetEntries( m_sumEntries );
std::vector<double> stat( 11 ); std::vector<double> stat( 11 );
// sum weights // sum weights
stat[0] = sumBinHeights(); stat[0] = sumBinHeights();
stat[1] = 0; stat[1] = 0;
if ( equivalentBinEntries() != 0 ) stat[1] = ( sumBinHeights() * sumBinHeights() ) / equivalentBinEntries(); if ( abs( equivalentBinEntries() ) > std::numeric_limits<double>::epsilon() )
stat[1] = ( sumBinHeights() * sumBinHeights() ) / equivalentBinEntries();
stat[2] = m_sumwx; stat[2] = m_sumwx;
double mean = 0; double mean = 0;
if ( sumBinHeights() != 0 ) mean = m_sumwx / sumBinHeights(); if ( abs( sumBinHeights() ) > std::numeric_limits<double>::epsilon() ) mean = m_sumwx / sumBinHeights();
stat[3] = ( mean * mean + rms * rms ) * sumBinHeights(); stat[3] = ( mean * mean + rms * rms ) * sumBinHeights();
m_rep->PutStats( &stat.front() ); m_rep->PutStats( &stat.front() );
return true; return true;
...@@ -155,7 +150,8 @@ bool Gaudi::Histogram1D::setStatistics( int allEntries, double eqBinEntries, dou ...@@ -155,7 +150,8 @@ bool Gaudi::Histogram1D::setStatistics( int allEntries, double eqBinEntries, dou
stat[0] = sumBinHeights(); stat[0] = sumBinHeights();
// sum weights **2 // sum weights **2
stat[1] = 0; stat[1] = 0;
if ( eqBinEntries != 0 ) stat[1] = ( sumBinHeights() * sumBinHeights() ) / eqBinEntries; if ( abs( eqBinEntries ) > std::numeric_limits<double>::epsilon() )
stat[1] = ( sumBinHeights() * sumBinHeights() ) / eqBinEntries;
// sum weights * x // sum weights * x
stat[2] = mean * sumBinHeights(); stat[2] = mean * sumBinHeights();
// sum weight * x **2 // sum weight * x **2
...@@ -167,7 +163,7 @@ bool Gaudi::Histogram1D::setStatistics( int allEntries, double eqBinEntries, dou ...@@ -167,7 +163,7 @@ bool Gaudi::Histogram1D::setStatistics( int allEntries, double eqBinEntries, dou
bool Gaudi::Histogram1D::fill( double x, double weight ) { bool Gaudi::Histogram1D::fill( double x, double weight ) {
// avoid race conditions when filling the histogram // avoid race conditions when filling the histogram
auto guard = std::scoped_lock{ m_fillSerialization }; auto guard = std::scoped_lock{ m_fillSerialization };
( weight == 1. ) ? m_rep->Fill( x ) : m_rep->Fill( x, weight ); m_rep->Fill( x, weight );
return true; return true;
} }
...@@ -192,7 +188,8 @@ void Gaudi::Histogram1D::copyFromAida( const AIDA::IHistogram1D& h ) { ...@@ -192,7 +188,8 @@ void Gaudi::Histogram1D::copyFromAida( const AIDA::IHistogram1D& h ) {
double sumw = h.sumBinHeights(); double sumw = h.sumBinHeights();
// sumw2 // sumw2
double sumw2 = 0; double sumw2 = 0;
if ( h.equivalentBinEntries() != 0 ) sumw2 = ( sumw * sumw ) / h.equivalentBinEntries(); if ( abs( h.equivalentBinEntries() ) > std::numeric_limits<double>::epsilon() )
sumw2 = ( sumw * sumw ) / h.equivalentBinEntries();
double sumwx = h.mean() * h.sumBinHeights(); double sumwx = h.mean() * h.sumBinHeights();
double sumwx2 = ( h.mean() * h.mean() + h.rms() * h.rms() ) * h.sumBinHeights(); double sumwx2 = ( h.mean() * h.mean() + h.rms() * h.rms() ) * h.sumBinHeights();
...@@ -215,11 +212,6 @@ void Gaudi::Histogram1D::copyFromAida( const AIDA::IHistogram1D& h ) { ...@@ -215,11 +212,6 @@ void Gaudi::Histogram1D::copyFromAida( const AIDA::IHistogram1D& h ) {
m_rep->PutStats( &stat.front() ); m_rep->PutStats( &stat.front() );
} }
#ifdef __ICC
// re-enable icc remark #1572
# pragma warning( pop )
#endif
StreamBuffer& Gaudi::Histogram1D::serialize( StreamBuffer& s ) { StreamBuffer& Gaudi::Histogram1D::serialize( StreamBuffer& s ) {
// DataObject::serialize(s); // DataObject::serialize(s);
std::string title; std::string title;
......
...@@ -173,16 +173,10 @@ bool Gaudi::Histogram2D::reset() { ...@@ -173,16 +173,10 @@ bool Gaudi::Histogram2D::reset() {
return Base::reset(); return Base::reset();
} }
#ifdef __ICC
// disable icc remark #1572: floating-point equality and inequality comparisons are unreliable
// The comparison are meant
# pragma warning( push )
# pragma warning( disable : 1572 )
#endif
bool Gaudi::Histogram2D::fill( double x, double y, double weight ) { bool Gaudi::Histogram2D::fill( double x, double y, double weight ) {
// avoid race conditiosn when filling the histogram // avoid race conditiosn when filling the histogram
auto guard = std::scoped_lock{ m_fillSerialization }; auto guard = std::scoped_lock{ m_fillSerialization };
( weight == 1. ) ? m_rep->Fill( x, y ) : m_rep->Fill( x, y, weight ); m_rep->Fill( x, y, weight );
return true; return true;
} }
...@@ -191,14 +185,17 @@ bool Gaudi::Histogram2D::setRms( double rmsX, double rmsY ) { ...@@ -191,14 +185,17 @@ bool Gaudi::Histogram2D::setRms( double rmsX, double rmsY ) {
std::vector<double> stat( 11 ); std::vector<double> stat( 11 );
stat[0] = sumBinHeights(); stat[0] = sumBinHeights();
stat[1] = 0; stat[1] = 0;
if ( equivalentBinEntries() != 0 ) stat[1] = ( sumBinHeights() * sumBinHeights() ) / equivalentBinEntries(); if ( abs( equivalentBinEntries() ) > std::numeric_limits<double>::epsilon() )
stat[1] = ( sumBinHeights() * sumBinHeights() ) / equivalentBinEntries();
stat[2] = m_sumwx; stat[2] = m_sumwx;
double meanX = 0;
if ( sumBinHeights() != 0 ) meanX = m_sumwx / sumBinHeights();
stat[3] = ( meanX * meanX + rmsX * rmsX ) * sumBinHeights();
stat[4] = m_sumwy; stat[4] = m_sumwy;
double meanX = 0;
double meanY = 0; double meanY = 0;
if ( sumBinHeights() != 0 ) meanY = m_sumwy / sumBinHeights(); if ( abs( sumBinHeights() ) > std::numeric_limits<double>::epsilon() ) {
meanX = m_sumwx / sumBinHeights();
meanY = m_sumwy / sumBinHeights();
}
stat[3] = ( meanX * meanX + rmsX * rmsX ) * sumBinHeights();
stat[5] = ( meanY * meanY + rmsY * rmsY ) * sumBinHeights(); stat[5] = ( meanY * meanY + rmsY * rmsY ) * sumBinHeights();
stat[6] = 0; stat[6] = 0;
m_rep->PutStats( &stat.front() ); m_rep->PutStats( &stat.front() );
...@@ -231,7 +228,8 @@ void Gaudi::Histogram2D::copyFromAida( const IHistogram2D& h ) { ...@@ -231,7 +228,8 @@ void Gaudi::Histogram2D::copyFromAida( const IHistogram2D& h ) {
// statistics // statistics
double sumw = h.sumBinHeights(); double sumw = h.sumBinHeights();
double sumw2 = 0; double sumw2 = 0;
if ( h.equivalentBinEntries() != 0 ) sumw2 = ( sumw * sumw ) / h.equivalentBinEntries(); if ( abs( h.equivalentBinEntries() ) > std::numeric_limits<double>::epsilon() )
sumw2 = ( sumw * sumw ) / h.equivalentBinEntries();
double sumwx = h.meanX() * h.sumBinHeights(); double sumwx = h.meanX() * h.sumBinHeights();
double sumwx2 = ( h.meanX() * h.meanX() + h.rmsX() * h.rmsX() ) * h.sumBinHeights(); double sumwx2 = ( h.meanX() * h.meanX() + h.rmsX() * h.rmsX() ) * h.sumBinHeights();
double sumwy = h.meanY() * h.sumBinHeights(); double sumwy = h.meanY() * h.sumBinHeights();
...@@ -255,8 +253,3 @@ void Gaudi::Histogram2D::copyFromAida( const IHistogram2D& h ) { ...@@ -255,8 +253,3 @@ void Gaudi::Histogram2D::copyFromAida( const IHistogram2D& h ) {
std::array<double, 11> stat = { { sumw, sumw2, sumwx, sumwx2, sumwy, sumwy2, sumwxy } }; std::array<double, 11> stat = { { sumw, sumw2, sumwx, sumwx2, sumwy, sumwy2, sumwxy } };
m_rep->PutStats( stat.data() ); m_rep->PutStats( stat.data() );
} }
#ifdef __ICC
// re-enable icc remark #1572
# pragma warning( pop )
#endif
...@@ -131,34 +131,30 @@ void* Gaudi::Histogram3D::cast( const std::string& className ) const { ...@@ -131,34 +131,30 @@ void* Gaudi::Histogram3D::cast( const std::string& className ) const {
return nullptr; return nullptr;
} }
#ifdef __ICC
// disable icc remark #1572: floating-point equality and inequality comparisons are unreliable
// The comparison are meant
# pragma warning( push )
# pragma warning( disable : 1572 )
#endif
bool Gaudi::Histogram3D::setRms( double rmsX, double rmsY, double rmsZ ) { bool Gaudi::Histogram3D::setRms( double rmsX, double rmsY, double rmsZ ) {
m_rep->SetEntries( m_sumEntries ); m_rep->SetEntries( m_sumEntries );
std::vector<double> stat( 11 ); std::vector<double> stat( 11 );
// sum weights // sum weights
stat[0] = sumBinHeights(); stat[0] = sumBinHeights();
stat[1] = 0; stat[1] = 0;
if ( equivalentBinEntries() != 0 ) stat[1] = ( sumBinHeights() * sumBinHeights() ) / equivalentBinEntries(); if ( abs( equivalentBinEntries() ) > std::numeric_limits<double>::epsilon() )
stat[1] = ( sumBinHeights() * sumBinHeights() ) / equivalentBinEntries();
stat[2] = m_sumwx; stat[2] = m_sumwx;
double meanX = 0;
if ( sumBinHeights() != 0 ) meanX = m_sumwx / sumBinHeights();
stat[3] = ( meanX * meanX + rmsX * rmsX ) * sumBinHeights();
stat[4] = m_sumwy; stat[4] = m_sumwy;
double meanY = 0;
if ( sumBinHeights() != 0 ) meanY = m_sumwy / sumBinHeights();
stat[5] = ( meanY * meanY + rmsY * rmsY ) * sumBinHeights();
stat[6] = 0; stat[6] = 0;
stat[7] = m_sumwz; stat[7] = m_sumwz;
double meanX = 0;
double meanY = 0;
double meanZ = 0; double meanZ = 0;
if ( sumBinHeights() != 0 ) meanZ = m_sumwz / sumBinHeights(); if ( abs( sumBinHeights() ) > std::numeric_limits<double>::epsilon() ) {
meanX = m_sumwx / sumBinHeights();
meanY = m_sumwy / sumBinHeights();
meanZ = m_sumwz / sumBinHeights();
}
stat[3] = ( meanX * meanX + rmsX * rmsX ) * sumBinHeights();
stat[5] = ( meanY * meanY + rmsY * rmsY ) * sumBinHeights();
stat[8] = ( meanZ * meanZ + rmsZ * rmsZ ) * sumBinHeights(); stat[8] = ( meanZ * meanZ + rmsZ * rmsZ ) * sumBinHeights();
// do not need to use sumwxy sumwxz and sumwyz // do not need to use sumwxy sumwxz and sumwyz
m_rep->PutStats( &stat.front() ); m_rep->PutStats( &stat.front() );
return true; return true;
} }
...@@ -198,7 +194,8 @@ void Gaudi::Histogram3D::copyFromAida( const AIDA::IHistogram3D& h ) { ...@@ -198,7 +194,8 @@ void Gaudi::Histogram3D::copyFromAida( const AIDA::IHistogram3D& h ) {
// statistics // statistics
double sumw = h.sumBinHeights(); double sumw = h.sumBinHeights();
double sumw2 = 0; double sumw2 = 0;
if ( h.equivalentBinEntries() != 0 ) sumw2 = ( sumw * sumw ) / h.equivalentBinEntries(); if ( abs( h.equivalentBinEntries() ) > std::numeric_limits<double>::epsilon() )
sumw2 = ( sumw * sumw ) / h.equivalentBinEntries();
double sumwx = h.meanX() * h.sumBinHeights(); double sumwx = h.meanX() * h.sumBinHeights();
double sumwx2 = ( h.meanX() * h.meanX() + h.rmsX() * h.rmsX() ) * h.sumBinHeights(); double sumwx2 = ( h.meanX() * h.meanX() + h.rmsX() * h.rmsX() ) * h.sumBinHeights();
double sumwy = h.meanY() * h.sumBinHeights(); double sumwy = h.meanY() * h.sumBinHeights();
...@@ -243,8 +240,3 @@ void Gaudi::Histogram3D::copyFromAida( const AIDA::IHistogram3D& h ) { ...@@ -243,8 +240,3 @@ void Gaudi::Histogram3D::copyFromAida( const AIDA::IHistogram3D& h ) {
stat[10] = sumwyz; stat[10] = sumwyz;
m_rep->PutStats( &stat.front() ); m_rep->PutStats( &stat.front() );
} }
#ifdef __ICC
// re-enable icc remark #1572
# pragma warning( pop )
#endif
...@@ -119,6 +119,6 @@ bool Gaudi::Profile1D::setBinContents( int i, int entries, double height, double ...@@ -119,6 +119,6 @@ bool Gaudi::Profile1D::setBinContents( int i, int entries, double height, double
bool Gaudi::Profile1D::fill( double x, double y, double weight ) { bool Gaudi::Profile1D::fill( double x, double y, double weight ) {
// avoid race conditions when filling the profile // avoid race conditions when filling the profile
auto guard = std::scoped_lock{ m_fillSerialization }; auto guard = std::scoped_lock{ m_fillSerialization };
( weight == 1. ) ? m_rep->Fill( x, y ) : m_rep->Fill( x, y, weight ); m_rep->Fill( x, y, weight );
return true; return true;
} }
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment