Undetectable NaN when using quiet_NaN() and -Ofast

user_1418351 0
Viewed 185

I am writing a classic nanmean function with OpenCV. I try to emulate MatLab's nanmean by default behaviour (i.e. nanmean reduce on the first dimension). I generate a matrix of random size which can be CV_32F or CV_64F with up to 4 channels. I fill it with random values following a uniform law. Then I assign some values to Nan using std::numerical_limits<float> (if CV_32F, double otherwise) :: quiet_NaN();

During the debugging step I was looking for an issue and I print the following:

T v = *it_src;
std::cout<<"v before: "<<v<<" "<<std::isnan(v)<<" "<<cvIsNaN(v)<<" "<<(v==v)<<" "<<" "<<std::isinf(v)<<" "<<((v+v)==v)<<std::endl;

T is template type, it can be either float or double, nothing else.

The output is:

v before: nan 0 0 1 0 1

So the value is "nan" but neither "std::isnan" nor "CvIsNan" can detect it. The comparison feature from IEEE 754 (if v is a Nan then v == v should be false) failed (v == v return true). The only thing that works is the last check ((v+v) == v).

I have few questions:

  • Simply why?
  • Is where this issue come from and how to fix it?
  • Can SIMD instructions also be concerned by this?

#include <opencv2/core.hpp>
#include <iostream>

using namespace cv;

int main()
{
    // Initialization
    int rows(theRNG().uniform(10,21)), cols(theRNG().uniform(10,21));
    int cn(theRNG().uniform(1, 5)); // [1, 5[ -> [1,4]
    
    Mat src(rows, cols, CV_32FC(cn));
    
    theRNG().fill(src, RNG::UNIFORM, 0, 10000);
    
    float ratio = theRNG().uniform(0.1f, 0.8f);
    int nb_points = saturate_cast<int>(src.total() * ratio);
    
    for(int i=0;i<nb_points;i++)
    {
     int x = theRNG().uniform(0, cols);
     int y = theRNG().uniform(0, rows);
     int z = theRNG().uniform(0, cn);
    
     src.ptr<float>(y,x)[z] = std::numeric_limits<float>::quiet_NaN();
    
    }
    
    
    Mat dst = Mat::zeros(1, cols, src.type());
    
    // Computation (reduce mean with omition of the Nan value over the first axis).
    const size_t src_step1 = src.step1();
    
    for(int c=0; c<cols; c++)
    {
        // The default constructor of Scalar_ initialize the elements of the attribute "val" to 0.
        Scalar sum;
        Scalar_<int> cnt;
    
        const float* it_src = src.ptr<float>(0,c);
    
        for(int r=0; r<rows; r++, it_src+=src_step1)
            for(int i=0;i<cn;i++)
            {
                float v = it_src[i];
    
                std::cout<<"v before: "<<v<<" "<<std::isnan(v)<<" "<<cvIsNaN(v)<<" "<<(v==v)<<" "<<" "<<std::isinf(v)<<" "<<((v+v)==v)<<std::endl;
    
                if(!std::isnan(v)) // Failing
                {
                    sum[i]+= saturate_cast<double>(v);
                    cnt[i]++;
                }
            }
    
        for(int i=0; i<cn;i++)
        {
            float den = saturate_cast<float>(cnt[i]);
            if(den==0.f)
                den = 1.f;
            dst.ptr<float>(0, c)[i] = saturate_cast<float>(sum[i]) / den;
        }
    }

 return 0;

 }
nan c++
1 Answers

You mention in the comments that you were using -Ofast, and this was causing the issue. To understand why this is, we start by looking at the GCC documentation for options that control optimizations. Here it lists the following options that are turned on by -Ofast:

It turns on -ffast-math, -fallow-store-data-races and the Fortran-specific -fstack-arrays, unless -fmax-stack-var-size is specified, and -fno-protect-parens.

For your case, looking at the description of -ffast-math seems relevant, as it sets the following set of options:

Sets the options -fno-math-errno, -funsafe-math-optimizations, -ffinite-math-only, -fno-rounding-math, -fno-signaling-nans, -fcx-limited-range and -fexcess-precision=fast.

I didn't look at what each of all these options do, but -ffinite-math-only specifically seems like it would be incompatible with any program that cares about NaN or Inf, as it is documented as:

Allow optimizations for floating-point arithmetic that assume that arguments and results are not NaNs or +-Infs.

This option is not turned on by any -O option since it can result in incorrect output for programs that depend on an exact implementation of IEEE or ISO rules/specifications for math functions. It may, however, yield faster code for programs that do not require the guarantees of these specifications. [emphasis added]

Since your use case clearly includes detecting NaN, it would seem ill-advised to use this flag.

Simple example program that demonstrates the issue, when compile with GCC 11.3 using the -Ofast option (Example on Godbolt because I don't have access to GCC right now). With -O3 nothing is printed as expected, and with -Ofast it prints Brought to you by -Ofast:

#include <limits>
#include <cstdio>
#include <cmath>

int main() {
    auto value = std::numeric_limits<float>::quiet_NaN();

    if(!std::isnan(value)) {
        std::puts("Brought to you by -Ofast");
    }
}
user_3270710
Related