#1911 invalidating data points with conditional operator does not work when data is binary

[Rainald Koch]

Hello,
this case is similar to https://sourceforge.net/p/gnuplot/bugs/1107/ of Hiroki Motoyoshi in that NaN leads to data points shown at the last value instead of being omitted. The difference to that bug is that the NaN is not read from the binary file but calculated while reading a binary file.

My test data consists of two columns,
$1 is 0 1 2 3 4 5 6 7 8,
$2 is 3 2 1 0 0 0 1 2 3.
I prepared two files containing that data, a binary and a text file.
I want to suppress $1 where $2 == 0.

plot 'binNaNtest.dat' \
binary record=0 format="%2uchar" using 0:1 w l, '' \
binary record=0 format="%2uchar" using 0:2 w l, '' \
binary record=0 format="%2uchar" using 0:(($2>0.5?1.:1/0)*$1) w p;

plot 'txtNaNtest.dat' \
using 0:1 w l, '' \
using 0:2 w l, '' \
using 0:(($2>0.5?1.:1/0)*$1) w p;

It works with the text file, but gnuplot shows three extra points with y=2 when reading the binary file.

[Ethan Merritt]

It works to place the NaN in first column rather than the second:

plot 'binNaNtest' using ($2>0.5? $0 : 1/0) : 1 w lp

[Rainald Koch]

THANKS

[Ethan Merritt]

After some poking about, I conclude this behaviour is an unintended side effect of an intended change. In gnuplot 4 any NaN/undefined field in the input data would result in that point being skipped, although there were inconsistencies between for instance "using 1:2" and "using ($1):($2)" In gnuplot 5 the code was changed to always return a data point but set an UNDEFINED flag if a requested field was NaN or otherwise undefined. The idea was the you might still be able to make use of x and y if, for example, a requested x_error field was not present. Your example shows that the new code has inconsistencies also.

I have patched the binary input code in CVS for both 5.0 and 5.1 so that it acts the same as ascii input for expression evaluations in a using spec. That takes care of this particular problem case but I am pretty sure that other corner cases remain.

[Rainald Koch]

Remotely related, feel free to move: When smoothing (e.g., smooth acsplines) meets undefined points, it seems to independently operate on connected components (which is reasonable for large gaps and can be avoided for small gaps by using dummy values with almost no weight). That splitting, however, leads to an error message "Cannot smooth: no data within fixed xrange!" when a component is zoomed out of sight. Gnuplot should not even try smoothing components out of sight.

set key top left
set print $data
do for [t=-11:211] {
    y = 1000*sin(t*0.015708)
    y = y<0.? 0.: (y + sqrt(2*y)*invnorm(rand(0)))
    print sprintf("%d %5.3g", t, y)
}
unset print
plot "$data" using 1:2 with points, \
"" using 1:2:(0.001) smooth acsplines with lines, \
"" using 1:(  ($2==0.||$1==99)? 1/0: $2):(0.001) smooth acsplines with lines lw 3, \
"" using 1:2:(($2==0.||$1==99)? 1e-8    : 0.001) smooth acsplines with lines lc rgb 'black';

.
Zoom around x = 100 to see the behaviour at the gap and then try to zoom around x = 0.

Edit: +Wikilink

[Ethan Merritt]

Hmm. I'll have to think about that one.

If you are using a gnuplot built from current CVS source (either 5.0 or 5.1) you can avoid this problem by saying

set datafile missing NaN

That tells gnuplot to treat data points with NaN values as if they were not present at all. Maybe when smoothing this should always be the case? As I say, I'll have to think about it.

[Ethan Merritt]

OK, here is the story.

This behaviour is intentional but poorly documented. In fact it is not documented at all that I can find. The "smooth" options treat the input data as a set of N curves separated by one or more points whose value is undefined. Each of these N curves is processed separately and there is no consideration given to continuity from one curve to the next. There are dedicated routines in the code to do exactly this: num_curves() and next_curve() in the file interpol.c. I don't know why the original authors chose to do it that way and I certainly don't know why they failed to document it more clearly.

So yes, if you need to treat the full set of points as describing a single curve despite some of them being undefined then the recent option set datafile missing NaN does exactly what you want. I can't think of a good work-around for previously released versions of the program that don't provide this option.

This should all be added to the documentation. Or if it is mentioned somewhere that I didn't find, that section should be referenced by the entries for "smooth", "data", and so on.

[Rainald Koch]

I'm using version 5.0 rc2 and know the effect of the option set datafile missing "NaN" on drawing lines (connecting points over gaps). My current data has large gaps in it, so that connecting points over gaps would not be appropriate. I nevertheless checked the behaviour of smoothing patchy data with that option and found the gap bridged with a wide-swing 3rd-order polynomial -- nice but inappropriate, and only after changing NaN to 1/0.

set datafile missing "NaN"
set print $data
do for [t=-100:100] {
    y = t*t < 100. ? NaN : (rand(0)+rand(0))
    print sprintf("%d %5.3g", t, y)
}
unset print
plot "$data" using 1:2 with lines, \
"" using 1:2:(0.001) smooth acsplines with lines, \
"" using 1:2:(0.1  ) smooth acsplines with lines;

The code above yields a plot as if set datafile missing "NaN" had no effect at all.
I then changed NaN to 1/0 also in this demo code and got a plot over only [-100:-10] and a weird error message. Boiled-down code reproducing that error is

do for [t=0:2] {
     y = t > 1.5 ? 1/0 : t
     print t
}

This is what happens:

gnuplot> do for [t=0:2] {
more>      y = t > 1.5 ? 1/0 : t
more>      print t
more> }
0
1

more> ;     y = t > 1.5 ? 1/0 : t;     print t;
                ^
      undefined value

Any idea?

[Ethan Merritt]

The "undefined value" refers to the entire expression being evaluated. The caret points to the start of the expression. When the iteration parameter t reaches 2, the expression evaluates to 1/0 and you get the same error message you would have gotten from

gnuplot> y = 1/0
             ^
         undefined value

Taking a step back...

I am always interested in correcting bugs or other problems with gnuplot, so I thank you for this report and test case. Having said that, if you are really dealing with data that looks like the points generated by your original test script then I think a spline fit of any sort is the wrong tool for the job. It seems that you want to fit a smooth curve through a fairly dense scattering of points. If you have some idea of the form of the curve desired then "fit" is probably the way to go rather than "smooth".

For the data generated by your latest test script where there is a big gap between the two clusters of points, the current default behaviour of "smooth" actually looks more reasonable: split into two data sets and fit each separately. Perhaps the program would do better to split the data based on looking for gaps (large change in x coordinate) rather than looking for undefined points?

[Rainald Koch]

So the latest bug is that NaN written to $data (as string?) is not recovered as a NaN on reading?

The result of "smooth acsplines" without "set datafile missing 'NaN'" is already fine for my data, except the bug occurring when zooming a component out of view. Again: Gnuplot should not even try smoothing components out of sight.

Another bug showed up while searching for the one above: acsplines recalculates on the zoomed range with results depending on the width in x -- try a width of ~10:

set print $data
do for [t=0:999] {
    y = invnorm(rand(0))
    print sprintf("%5.3g", y)
}
unset print
plot "$data" using 0:1 with points, \
          "" using 0:1:(100.) smooth acsplines with lines;

It seems that the provided weights, which I think should be taken as absolute, are somehow 'relativized'. Compare the manual, "the absolute magnitude of the weights determines..."

[Ethan Merritt]

So the latest bug is that NaN written to $data (as string?) is not recovered as a NaN on reading?

So far as I know that works correctly. Do you have an example where it fails? If you are talking about string matching, are you sure that it was written as "NaN" and not "nan"? String input is case-sensitive; numerical input is not.

Gnuplot should not even try smoothing components out of sight.

I guess that most people would expect the opposite: that the curve will not change shape just because you zoom in on an interior segment. If the program re-fits to ignore points outside the current zoom window then the curve will [probably] not retain its previous y value at the left/right borders of the zoom window.

[Rainald Koch]

I am sorry for the confusion. Too many bugs for a single thread. Let's make a tree out of it, first leaf (closed, case sensitive):

Do you have an example where it fails?

Given above, search for:

as if set datafile missing "NaN" had no effect

Edit: closed.

[Ethan Merritt]

I see the expected effect with all of version 5.0.0, version 5.0.5 and current CVS. So I guess I am not reproducing what you get or else I misunderstand the test code you want me to use.

[Rainald Koch]

I copy-paste the test code without its first line and get the expected result. I then copy-paste the whole test code and get the same result, which is unexpected. Both the data and the approximation should be connected, by a straight line and a 3rd-order poly, respectively.

[Ethan Merritt]

I can't reproduce that in any gnuplot version I have here.
2 attached plots produced by gnuplot version 5.0.0

set print $data
do for [t=-100:100] {
    y = t*t < 100. ? NaN : (rand(0)+rand(0))
    print sprintf("%d %5.3g", t, y)
}
unset print
set output 'default.png'
plot "$data" using 1:2 with lines, \
"" using 1:2:(0.001) smooth acsplines with lines, \
"" using 1:2:(0.1  ) smooth acsplines with lines
set datafile missing "nan"
set output 'set_missing_nan.png'
replot

[Rainald Koch]

The difference between our codes is "NaN" versus "nan".

set datafile missing "nan" is working, because that is the string written for the special variable NaN (as I just learned by writing to the terminal window instead of to $data).

I had set "NaN" following your suggestion 2 days ago, which is in line with the section 'Set datafile missing' in the manual,

invalid data (e.g. "NaN", 1/0.)

which should be corrected.

Edit: What about using the string set via set datafile missing "<string>" for both reading and writing?</string>

[Ethan Merritt]

Whether it is written in formatted output as "nan" or "NaN" depends on the OS compiler or libc implementation of sprintf(). Similarly it is OS-dependentwhether that string is accepted as valid floating point formatted input. That is a different matter than string matching.

Gnuplot does have an OS-independent formatted output routine gprintf() that is roughly the same as standard sprintf() for a single output value, with the addition of gnuplot-specific formats. See "help gprintf". If you replace the sprintf in your test script with
print gprintf("%g ",t), gprintf("%g ",y)
the program will consistently output NaN (capital Ns) even if your local libc sprintf would have produced lower-case "nan".

Note that my original suggestion did not place quotes around NaN:
set datafile missing NaN
with no quotes is a recent variant that your gnuplot executable doesn't support. This recent variant does catch the strings "nan" or "NaN" in an input file because both are valid representations of floating point not-a-number. This is different from the mechanism available in earlier gnuplot versions that only compared string values. The new variant was added partly to reduce this kind of ambiguity.

[Rainald Koch]

2nd leaf (closed):

I guess that most people would expect the opposite: that the curve will not change shape just because you zoom in on an interior segment. If the program re-fits to ignore points outside the current zoom window then the curve will [probably] not retain its previous y value at the left/right borders of the zoom window.

So do I, but your "if" is almost happening: The program refits. It does so for unknown reason. It does not ignore points (changing only y values at the border), but changes all y values as if the weights were scaled up. Please zoom deeply into my latest example above.

Edit: I suspect now that the change in y values -- in extreme cases, a minimum occurs in a zoomed view where a maximum showed up in the full view -- is merely due to undersampling in the less magnified views, where the approximation of the spline approximation by linear segments is missing some structure. This is a feature, not a bug.

[Rainald Koch]

3rd leaf:
This bug is not about points but about connected components. These are the items counted in num_curves(). If (at least) one of them is completely zoomed out of view, the zoomed view is not generated and an error message hints at the bug:

Cannot smooth: no data within fixed xrange!

[Ethan Merritt]

That results from this compilation choice in the source to interpol.c

    /* HBB 20010727: Sample only across the actual x range, not the

     * full range of input data */
#if SAMPLE_CSPLINES_TO_FULL_RANGE 
    xstart = this_points[0].x;
    xend = this_points[num_points - 1].x;
#else
    xstart = GPMAX(this_points[0].x, sxmin);
    xend = GPMIN(this_points[num_points - 1].x, sxmax);

    if (xstart >= xend)
        int_error(plot->token,
                  "Cannot smooth: no data within fixed xrange!");
#endif

I am not familiar with the history of this code. Since Hans-Bernhard added the comment, maybe he can reconstruct the reason the default has ever since then been to limit the sample range.

Would you favor going back to the alternate code path, or changing the error to a non-fatal warning, or some third action?

[Rainald Koch]

There are two aspects, the proper density of samples and the guard against sampling invalid regions.

On the first aspect: The two options for the sample density in the code above are both reasonable at times. A third approach would be distributing samples_1 samples equally over all visible spline intervals (3rd-order polynomials), which would make a large difference in case of unevenly spaced data and/or unevenly distributed weights, see attached plot.
Having more than the appropriate number of samples, however, is no problem on current hardware. Therefore, instead of complicating the code, I suggest to go back to the alternate code path, remove the switch, and increase the default value of samples_1 to 1000. That would also fix the problem.

An alternative patch (visible to less many users) copes with the second aspect. Remove the switch leaving the current code for xstart and xend, and move the guard:

    xstart = GPMAX(this_points[0].x, sxmin);
    xend = GPMIN(this_points[num_points - 1].x, sxmax);

    xdiff = (xend - xstart) / (samples_1 - 1);

    if (xdiff > 0.) for (i = 0; i < samples_1; i++) {

I have not checked the calling code whether zero samples for a connected component is ok.

[Hans-Bernhard Broeker]

I wrote that code, added the comment, and even documented the behaviour it achieves, some 15 years ago. That sample is the implementation of this paragraph of "help smooth":

If autoscale is not in effect, and the smooth option is either acspline or cspline, the sampling of the generated curve is done across the intersection of the x range covered by the input data and the fixed abscissa range as defined by set xrange.

Before, gnuplot would try to either sample splines outside the range supported by their input data, or produce samples outside, sometimes way outside, the current viewport, causing all kinds of weird, sometimes broken output. But yes, I'll agree that it may interact poorly with the even older behaviour of smoothing and sampling interrupted polylines one fragment at a time.

[Ethan Merritt]

After staring at the code and running a few tests, my inclination is to leave the code the way it is but remove the early exit via int_error(). The data points making up the segment that triggers this code path are necessarily out of the current plot range so they will not appear in the plot no matter what values are assigned. So far as I can see it is sufficient to mark them OUTRANGE and return directly rather than via int_error().

Do you see any problem with this:

#if SAMPLE_CSPLINES_TO_FULL_RANGE
    xstart = this_points[0].x;
    xend = this_points[num_points - 1].x;
#else
    xstart = GPMAX(this_points[0].x, sxmin);
    xend = GPMIN(this_points[num_points - 1].x, sxmax);

    if (xstart >= xend) {
        /* This whole segment is outside current xrange. */
        /* Not a fatal error, we won't draw it at all. */
        for (i = 0; i < samples_1; i++)
            dest[i].type = OUTRANGE;
        return;
    }
#endif

[Rainald Koch]

I see no problem.
I would omit the comments. The first one just states what is obvious from the code (type = OUTRANGE for all samples) and the second one comments on the previous code and the change, not on a difficulty of the problem to solve.

[Hans-Bernhard Broeker]

Do you see any problem with this:

No. Looking good.