bugGNU Octave - Bugs: bug #58956, Empty plot for certain xlim / ylim...


bug #58956: Empty plot for certain xlim / ylim settings

Submitter:  Hg200 <hg200>
Submitted:  Fri 14 Aug 2020 06:12:22 PM UTC
Category:  Plotting with OpenGL Severity:  3 - Normal
Priority:  5 - Normal Item Group:  Matlab Compatibility
Status:  Patch Reviewed Assigned to:  None
Originator Name:  hg200 Open/Closed:  * Open
Release:  * other Operating System:  * Any
Fixed Release:  None Planned Release:  None
* Mandatory Fields

Add a New Comment Rich Markup

Jump to the original submission

Sun 10 Apr 2022 06:43:40 PM UTC, comment #32: 

The Internet is not clear about it, but it seems that it depends on the target system whether the data is treated internally as a float or a double. Although the prototype contains double. The code

t = now-1:(1/1440):now;
x = -cos(2*pi*t) + rand(size(t))/10;

looks on my system as in the attached pictures. One image is with the patch file #49866, the other is without.

Hg200 <hg200>
Mon 04 Apr 2022 07:07:01 AM UTC, comment #31: 

At least on the system that I could test, it does make a difference. Afaict, we are using double precision interfaces of the OpenGL function (e.g., `glMultMatrixd` or `glVertex3d`). I'd guess that the transformations are still done in double precision. The result is converted to single precision before it is rendered.
So taking your example from comment #28, it looks like the effective code would look more like the following with your patch:

t = now-1:(1/1440):now;
x = -cos(2*pi*t) + rand(size(t))/10;
t = single(t - fix(now));

Markus Mützel <mmuetzel>
Group administrator
Sun 03 Apr 2022 06:05:04 PM UTC, comment #30: 

@mmuetzel: The patch manipulates (double precision) matrices that are passed to OpenGl, which still performs the actual transformation work, in single precision. Therefore, even with the patch, we still have this problem of inaccuracy (see comment #28 the note about "glVertex3d" and "opengl_renderer::draw_line"). The point was to try to see if the numerical accuracy (in OpenGL) could be improved by rearranging the coordinate system. However, this comes at the cost of losing the camera target, which in turn leads to incompatibilities with Matlab. Since we get no big improvement for what we lose, please do not use the code. My intention was to just share some thoughts how we could work (or not work) on the topic.

From this attempt my conclusion is that we cannot easily fix the problem in graphics.cc. And btw, even with a newer Opengl API (for example writing into VBOs) we will still have to deal with this precision problem. IMO regarding this bugreport the best compromise is your patch in bug58956_gl_clipping_plains.patch.

So, although not fully fixed, I suggest to close this thread with a link to the OpenGl collection bug #59418. Maybe we can create a link in bug #59418 to here.

Hg200 <hg200>
Sun 03 Apr 2022 10:59:39 AM UTC, comment #29: 

@hg200: Your patch has been sitting here for a while (more than 1.5 years now). Imho, it would still be very interesting because it would probably avoid the most pressing issues users are seeing when having, e.g., datenums (i.e., pretty large numbers that hardly differ) on one of the axes.
IIUC, the patch basically subtracts the offset of the view target in double precision (before the OpenGL backend converts the coordinates to single precision).

We are currently early in the development cycle for Octave 8. Do you think it would make sense to apply it now (if it still applies...)?

Markus Mützel <mmuetzel>
Group administrator
Mon 02 Nov 2020 12:06:45 AM UTC, comment #28: 

<< With the OpenGL coordinate system centered at the target, do you think that it would be possible to do the rescaling for bug #32980 with these transformation matrices (in a follow up patch)?

I think it is not possible to fix bug #32980 comment 25 by merely changing the view transformation and the projection matrix. If I take the example from comment 25, convert the time to simple precision and try to "normalize" the view by translating the abscissa by a constant, the same artifacts are produced:

t = now-1:(1/1440):now;
x = -cos(2*pi*t) + rand(size(t))/10;
t = single(t) - single(fix(now));

Most likely the data is already unrecoverable when we leave opengl_renderer::draw_line() because eventually glVertex3d() is called.

Hg200 <hg200>
Sun 01 Nov 2020 02:07:48 PM UTC, comment #27: 

Markus, I wouldn't push this patch for now. It is not "clean" to do this transformation, although it is technically allowed for orthographic projection and actually solves some of our problems. But, if someone wants to implement "perspective" (not implemented in Octave yet, but available in Matlab) the shift into the origin will cause problems because of the perspective divide. And: I did a test implementation of perspective (full implementation is tricky because interface between bounding box and view port) in Octave and can observe similar clipping artifacts like in glOrtho, which is not too surprising. But then we have to solve this problem too. By the way, I can reproduce the accuracy problems of "glClipPlane" in an offline Qt demo code that is completely independent from the Octave source code. So we can add glClipPlane as a new and separate problem to the OpenGl Open Point list.

More and more points are coming to light and I think it is a bit difficult to continue the OpenGl discussion here because there are quite some points and questions. I have started to prepare a presentation on the topic (explanation of update_camera, collecting the fall traps with respect to our observations etc.), but could not finish the slide for our last octave developer meeting. If there is interest, we could e.g. talk a bit more generally about the topic in one of the next developer meetings. Or maybe discuss it in a workshop with a reduced number of participants, since not everybody might be interested in this topic.

JWE has the OpenGl back end on his wish list [1]. Currently, if I could make a wish, I would remove all legacy glvertex, write everything to local memory (allocated, double type), scale everything down, convert everything to GLfloat and then pass everything to a VBO and let the vertex and fragment shaders do the rest. But apart from the runtime, I have a doubt that simple normalization will cure us of all the problems we have. It will heal us from axes larger than an order of ^38, but probably not from clipping?

Any comments? Also from other members?

[1]: https://wiki.octave.org/JWE_Project_Ideas

Hg200 <hg200>
Sun 01 Nov 2020 11:27:36 AM UTC, comment #26: 

@hg200: I tested your patch from comment #24. It seems to work very good for me.
I ran the dump plot demos script (in scripts/testfun/private):


All plotted demos look correct to me. In particular, afaict there is no visual difference in the demos for light (in comparison to figures I dumped in June).
I'd like to push your patch to the default branch (Octave 7). Is that ok with you?

With the OpenGL coordinate system centered at the target, do you think that it would be possible to do the rescaling for bug #32980 with these transformation matrices (in a follow up patch)?

Markus Mützel <mmuetzel>
Group administrator
Tue 06 Oct 2020 09:54:18 AM UTC, comment #25: 

@hg200: I'm very happy with the progress so far.
Your analysis sounds reasonable. I didn't try the patch. But from what you describe, it seems to work reasonably well.
I didn't get the point about "cameraposition" though. Why would it make a difference wrt projection if everything (including the camera position and target!) is translated to a different coordinate system?
I also don't understand why lighting would be special?
All of these translations are internal to the rendering. None of them should have any influence on the properties returned by `get(gca)`.
Or was your point the following?: The currently observed deviation between the "cameraposition" in Octave and in Matlab doesn't matter as long as we only support orthographic projection. But it will be important once we'd like to implement perspective projection. That point however applies completely independent of the used coordinate system.
Is the projection done by OpenGL or do we have to manipulate all object coordinates to give the "illusion" of perspective projection? Anyway, this is probably a better topic for a report dedicated to implementing perspective projection.

Do you agree that it would make sense to move from a camera centered coordinate system to a target centered coordinate system?
Would it be better to use an object centered coordinate system? (I guess target centered like you described is the best choice and easier to handle.)

Re comment #23: The camera target is not necessarily the mean value of xlim,ylim,zlim (logarithmically scaled axes!). And it is possible to set a different camera target (see e.g. `demo camlookat`).

Markus Mützel <mmuetzel>
Group administrator
Sun 27 Sep 2020 10:46:34 AM UTC, comment #24: 

Small update:

Concerning comment #22:
It looks like we have two clipping algorithms in the line primitive. One is in "uint8_t clip_code (double x, double y, double z)", which can be found in gl-render.h. This clipping routine is done in "draw_line" just before a line vertex is sent to the OpenGl system via "glBegin(GL_LINE_STRIP)". It prevents invisible points from being sent to the renderer. The other clipping is performed in OpenGl via the OpenGl clipping planes and is used to cut off data when the "pan" symbol in a figure is clicked and the figure is moved around. The latter causes the problems because it unfortunately also cuts payload when the axis are set large. The coordinates are probably given in world coordinates, but since clipping is performed as one of the last operations, the coordinates are internally transformed throughout the entire view matrix most likely.

Concerning comment #20:
Just to see what happens, I moved the camera target to the origin [0,0,0] and reverted the changes in gl-render.cc, i.e. x_zlim is now centered around z=0, which is also the projection plane. Octave still shows the view point like matlab when the command "get(gca)" is executed, but in orthographic projection, the distance between view point and camera target has no meaning since it has no influence on the image. So we can do this manipulation without violating any view point properties. It turns out that we do it better if we change the matrix: I do not see any clipping problems, and also large axes are displayed correctly. "xlim ([-1e42, 3e42])" still does not work, as expected. I don't know how "light" works, but at first sight it doesn't seem to be too bad. I have attached an experimental "proof of concept" patch to see the effects.
However there is a problem: Although not yet implemented in Octave, Matlab also offers a perspective view ("projection" can be either "orthographic" (default) | "perspective"). It should not be too difficult to implement "perspective", although I failed on the first attempt. But if we move the target to the origin because of the normalization problems, we will run into big problems with "perspective", because unlike orthographic projection, the distance between camera target and view point has very big impact in perspective view. I wonder how the Matworks people solve this.

This topic will still require some engineering ...

(file #49866)

Hg200 <hg200>
Sun 20 Sep 2020 06:31:00 PM UTC, comment #23: 

>> Would it be possible to use a different coordinate system for the final rendering? E.g. one that has its origin in the center of the scene that is rendered and the camera position (possibly far out) on the negative z-axis?

IMHO yes but you have take care of some restrictions. First of all it should be emphasized that you do not have to use the functions glOrtho () and gluPerspective () either. Without invoking these functions, the "active region" of the OpenGl renderer becomes the cuboid [-1,+1]^3 and the projection plane becomes z=0. Everything you place in this region will appear on the screen.

For example what we can do is to write a new update_camera() that transforms the camera target (in Octave this is the mean value out of xlim,ylim,zlim) into the origin [0,0,0]. glOrtho () must then be removed from gl-render.cc. The problem with this approach is that we still have to make the transformation from "world coordinates" into "camera coordinates". This is done with the command

m_glfcns.glMultMatrixd (x_mat1.data ());

so the transformation remains in OpenGl and hence will be single precision. My guess is that the limit of approx 3.4e+38 on the axes will remain. And i don't know how "clipping" and how "light" is done in this approach. But it sounds like a nice project for next vacation just to find out what is happening if this change is implemented.

Hg200 <hg200>
Sun 20 Sep 2020 05:52:31 PM UTC, comment #22: 

There might be a similar issue (difference between single precision floating point numbers) in opengl_renderer::set_clipbox:

I don't know which coordinate system is used for the input arguments.
If it is the "camera coordinate system", maybe we need to make sure to extend the clipbox taking the "sparse neighbors" in single precision into account.
If it is the "world coordinate system", maybe we still need to make sure that that problem won't cause issues when the OpenGL routines transform the clipping planes into the camera coordinate system.

At the moment, it looks like it might be increased to each side by 0.1% of the extent of the scene(?). We might need to increase that if the plains are "close in single precision" in the camera coordinate system.

Markus Mützel <mmuetzel>
Group administrator
Sun 20 Sep 2020 04:04:38 PM UTC, comment #21: 

with respect to the line problem: the problem disappears when

void opengl_renderer::set_clipping (bool enable)

is erased. this seems to be the clipbox which is initialized in opengl_renderer::set_clipbox with the axes (plot box) limits. i tend to say this is another opengl precision "feature".

Hg200 <hg200>
Sun 20 Sep 2020 04:01:46 PM UTC, comment #20: 

You are probably right.
It very much looks like we need something similar to that change. And if we should find a better solution, we can still change what we have so far.
So I cleaned up the patch and pushed it here:

You are understanding the applied transformations better than I do. So please tell me if this doesn't make any sense:
Would it be possible to use a different coordinate system for the final rendering? E.g. one that has its origin in the center of the scene that is rendered and the camera position (possibly far out) on the negative z-axis?
That might relax issues with large single precision floating point numbers (whose neighbors are "sparse").

Markus Mützel <mmuetzel>
Group administrator
Sun 20 Sep 2020 11:01:08 AM UTC, comment #19: 

file #49794, while not perfect, is definitely much better than what we currently have implemented.

In the meanwhile I have completed the "reverse engineering" of update_camera (). The object we want to plot is defined in the Octave frame that is also displayed when the "axes ()" command is executed. I would call these coordinates "world coordinates". At the other side OpenGl is an API for rendering vector graphics. Unsurprisingly, the data to be visualized has to be passed to OpenGL in a "normalized format". For example, the “view point” and the "camera target" are on the Z axis. OpenGl provides two additional “built-in” projection matrices with which the projection properties can be augmented. glOrtho () is one of them, but basically the viewing direction is still the z-axis. These are the camera coordinates and they are different from the octave world coordinates.

Two important matrices are calculated in update_camera (). These are "x_gl_mat1" and "x_gl_mat2". The matrix "x_gl_mat1" is what I would call the "view matrix", and it is basically a transformation of the octave world coordinates into camera coordinates. For a better understanding I have put some pictures here:


"x_gl_mat1" is constructed from operation steps 1-8, which are a series of translations, scalings and one rotation. The near / far values ​​(clipping values) are finally calculated using a subset of these steps, i.e. a "normalized" view matrix transformation. The near and far values ​​are obtained by pushing a "unit cube" through this subset to get the "rotated" plotbox in camera coordinates. Min and max values ​​are used to get the max and min Z coordinates of the box.

To fix this bug, my idea was to change the view transformation so that the near and far values ​​are calculated depending on the distance between the view point and the camera target. However, a look at this section of code shows that in the end this is basically the same as manipulating the near and far values ​​in gl-render.cc. We're already doing that here.

Hg200 <hg200>
Mon 14 Sep 2020 06:33:22 AM UTC, comment #18: 

Are lines the only graphics object that disappear for certain axis limits? In that case, it might also be that there is an additional (precision?) bug in the calculation of the z-limits for line graphics objects...

Markus Mützel <mmuetzel>
Group administrator
Sun 13 Sep 2020 05:57:32 PM UTC, comment #17: 

So it looks like the "single precision difference" idea was a red herring.
I wonder what the real limiting factor is here...

The problem with the last command you showed is the long standing bug #32980.
It might be that solving that bug - which in itself is probably quite a large project - would relax the error condition for the bug here...

It is nice to read that you are making some progress on understanding the transformations. Unfortunately, I don't know if there is any kind of documentation apart from what is in the code.

Markus Mützel <mmuetzel>
Group administrator
Sun 13 Sep 2020 02:47:48 PM UTC, comment #16: 

Our posts are crossing, but i am attaching mine anyway to see the different solutions ...

Thanks for the prototype! Just a shoot in the dark to get an idea about the performance: If we keep the normalized difference between near and far around two float epsilons it seems to fix the problem for axes and scatter graphic objects and avoids the opengl errors too. However to have reliable code it is probably advisable to add some more saftey margin.

const double opengl_eps = std::numeric_limits<float>::epsilon ();
const double sum = x_zlim(1) + x_zlim(0);
const double diff = x_zlim(1) - x_zlim(0);

if (2.0 * opengl_eps * sum > diff)
    xZ1 = (0.5 - opengl_eps) * sum;
    xZ2 = (0.5 + opengl_eps) * sum;
    xZ1 = x_zlim(0) - 100.0 * diff;
    xZ2 = x_zlim(1) + 100.0 * diff;

This approach is already much better as what we have. Although i haven't checked all details (e.g. one of x_zlims can become negative but i think it excludes with the true branch), i saw the problem with line plots is not fixed. Sometimes lines are displayed, sometimes they are not. I don't know what the trouble is, but to display lines under all conditions the near / far values must be chosen much wider:

xZ1 = 0.49 * sum;
xZ2 = 0.51 * sum;

To check the problem run following code:

xlim ([-172099058766.5512 797603098357.2810]);
ylim ([-1.535202916065375e+11 2.103152628728371e+12]);
zlim ([-1 1]);
hold on;
line ([-1e+11 5e+11], [1e+12 2e+11]);
plot (0, 0, 'r*');
axis equal;
xlim ([-1.7e+11 7.0e+11]);
ylim ([-1.5e+11 2.1e+12]);
zlim ([-1 1]);
hold on;
line ([-1e+11 5e+11], [1e+12 2e+11]);
plot (0, 0, 'r*');
axis equal;

Also the command

xlim ([-1e42, 3e42]);
warning: opengl_renderer: data values greater than float capacity.  (1) Scale data, or (2) Use gnuplot

still fails. But in Matlab it opens axes. So it looks like there other limitations too. To get rid of all problems we probably need more than one fix. Maybe need to normalize the OpenGl modeltransformation too?

I'm making slow progress with axes::properties::update_camera(). What i understand is that the view point, target and upvector are calculated in the first part. Depending on aspect ratio and displacement from the origin, the OpenGl 4x4 matrices are then calculated by matrix-scalar and -sift operations. The trick is to understand what exactly x_view, x_projection, x_normrender, x_pre is doing. I have translated axes::properties::update_camera() from .c to .m, but when i apply x_gl_mat1 to a random scene i still obtain garbage. Looks like a bigger project for next vacation.

Hg200 <hg200>
Sun 13 Sep 2020 02:21:17 PM UTC, comment #15: 

The attachment got lost when I clicked the "Preview" button.
Here it is.

(file #49794)

Markus Mützel <mmuetzel>
Group administrator
Sun 13 Sep 2020 02:20:17 PM UTC, comment #14: 

The attached patch translates the code outlined in comment #13 to the C++ function.
With it the original example and the examples in bug #54551 render correctly for me. But I haven't tested much more yet.
I'm not sure which value we should chose for ´singlePrecFac´. At the moment it is similarly arbitrary as the ´expansionFac´ for bug #54551.

Markus Mützel <mmuetzel>
Group administrator
Fri 11 Sep 2020 04:40:57 PM UTC, comment #13: 

Maybe we could use some logic like this in the C++ code (written in Octave for simplicity):

Z = [707106781186.04761, 707106781187.04761];

if sum(Z)*eps('single')*5 > diff(Z)*100
  Znew = sum(Z)/2 + sum(Z)*eps('single')*5 * [-1,1];  % 5 is an arbitrary value here. Choose something that works
  Znew = sum(Z)/2 + diff(Z)*100 * [-1,1];  % close to what we currently have for bug #54551 without the ´min´ and ´max´

span = diff(single(Znew))

If I got this right, ´span´ should always be greater than 0 after that for any Z with Z(1)<Z(2).
I'm still struggling to understand why the ´min´ and ´max´ that we currently have are needed.

Markus Mützel <mmuetzel>
Group administrator
Mon 07 Sep 2020 06:30:13 PM UTC, comment #12: 

The single precision thresholds are where the Qt demo code stops displaying the cube and they are also near the critical value for xZ1 and xZ2 in gl-render.cc. Here are some random scenarios i tested:

Octave plot is incorrect (xZ1, xZ2):
>> single(70710676.618654922) - single(70710679.618654922)
ans = 0

Octave plot is correct (xZ1, xZ2):
octave:1> single(56568540.994924024) - single(56568543.994924024)
ans = -4

Simple Qt demo program plot is incorrect (farVal, nearVal):
>> single(4e+7 + 2.0) - single(4e+7 - 2.0)
ans = 0

Simple Qt demo program plot is correct (farVal, nearVal):
>> single(3e+7 + 2.0) - single(3e+7 - 2.0)
ans = 4

But hold on, there is still more! Do you remember fix of bug #54551? Now look at this:

octave:1> daspect ([1, 1, 1]);
octave:2> xlim ([-1e11, 1e11]);
octave:3> ylim ([-1e11, 1e11]);

xZ1 = std::max (-1e6, x_zlim(0)-(x_zlim(1)-x_zlim(0))*100.0);
xZ2 = std::min (1e6, x_zlim(1)+(x_zlim(1)-x_zlim(0))*100.0);
m_glfcns.glOrtho (0, vp(2), vp(3), 0, xZ1, xZ2);

(gdb) print *x_zlim.rep.data@x_zlim.rep.len
$1 = {707106781186.04761, 707106781187.04761}
(gdb) print xZ1
$2 = 707106781086.04761
(gdb) print xZ2
$3 = 1000000

octave:1> single(1000000) - single(707106781086.04761)
ans = -7.0711e+11

Surprise is that the plot in Octave is still invisible ;-) So my conclusion is that the position of the object inside the view frustum plays also a role. It looks like if the object distance to the clipping plane becomes too small over camera object distance that the plot also becomes invisible. This problem is solved for me if i increase xZ1 by a factor of 1.1

So maybe this is a good time to recall the observations:

1.) There is a lmitiation by the floating point expression "farVal-nearVal" that becomes relevant for large xlim/ylim settings (in the order greater 1e+8) and daspect([1,1,1])

2.) The position of the object with respect to the clipping plane over camera object distance also plays a role.

3.) The statement

xZ1 = std::max (-1e6, x_zlim(0)-(x_zlim(1)-x_zlim(0))*100.0);
xZ2 = std::min (1e6, x_zlim(1)+(x_zlim(1)-x_zlim(0))*100.0);

has a problem in the transition region (xlim, ylim is in the order of 1e+5), where "farVal-nearVal" can become zero, small or even reverse. The plot becomes also invisible or causes render errors.

4.) There was bug #54551. Don't forget to fix this

5.) Unanswered question in comment #6:

> Is there a reason the proposed changeset breaks the existing symmetry between lower and upper clipping limits?

Actually, the idea was to increase the distance between farVal and nearVal only for the part of the scene which is in front of the camera, not behind. The attached file #49697 solves problem 3.) but it does not solve problem 1.)+2.).

Hg200 <hg200>
Mon 07 Sep 2020 07:21:01 AM UTC, comment #11: 

That might be far fetched: But most of OpenGL "is" single precision. The interface to glOrtho uses ´double´, but it might be that these values are cast to ´float´ somewhere for its calculations.
Depending on where this conversion happens, it might be necessary that "nearVal" and "farVal" must be different in single precision to avoid a division by 0 (or something similar):

>> single((1e+11 + 2.0) - (1e+11 - 2.0))
ans =
>> single(1e+11 + 2.0) - single(1e+11 - 2.0)
ans =
>> single(1e+6 + 2.0) - single(1e+6 - 2.0)
ans =

Maybe we should make sure that "nearVal" and "farVal" differ on the 23-bit mantissa of a single precision floating point number?

Markus Mützel <mmuetzel>
Group administrator
Mon 07 Sep 2020 06:33:06 AM UTC, comment #10: 

Attached is a demo code that is independent of Octave and that demonstrates the same behavior in Qt. The code places a unit cube (length = 1) between "nearVal" and "farVal" and sets "nearVal" and "farVal" in glOrtho() to a large value, while the difference "farVal - nearVal" is chosen to be small in comparison to the object-camera distance. The adjustment can be done in "mainwindow.h":

double nearVal = 1e+6 - 2.0;
double farVal = 1e+6 + 2.0;

If the values are selected as

double nearVal = 1e+11 - 2.0;
double farVal = 1e+11 + 2.0;

the cube is not displayed, which is consistent with the observation made in Octave.

The two relevant commands are:

glOrtho(-aspectratio, aspectratio, -1, 1, nearVal, farVal);

glTranslated(0.0, 0.0, -(farVal - 2.0));

The cube is displayed when the difference "farVal-nearVal" becomes larger. For example, if nearVal is set to zero or set to a negative value (moving the clipping plane behind the camera is explicitly allowed for glOrtho). And / or farVal is selected as a factor greater than 1.0 of the object-camera distance. This observation is also consistent with what i observe in Octave.

Adding / subtracting a small value to a large value in floating point arithmetic is dangerous, but 1e+11+2.0 should be fine. In any case, we might also want to draw axes with 1e+22 rather than 1e+11, and then the problem becomes more obvious.

What is really happening in Octave depends now on the modelview- and projection matrix calculated in graphics.cc (update_camera). But according to GDB it looks like the "farVal-nearVal" difference is not proportional to the distance between the object and the camera. If no one can provide information about it, I will dig more into graphics.cc now.

(file #49750)

Hg200 <hg200>
Sat 29 Aug 2020 08:00:22 AM UTC, comment #9: 

Thank you for your feedback. Speaking of alternative bugfixes, I don't think that just correcting the camera position is enough to solve the whole problem.

If I undo the patch from bug #54551 and add

set(gca, 'CameraPosition', [0.014392034052167, 0.134014656653948, 2.411976487699734]*1e13);

after "axis equals" in three_body_planar.m (see comment #8), the axis and scatter objects are displayed correctly, but the trajectories (= line plots) in the first figure, before the loop begins are still invisible. To see the demo correctly, you may take a fresh copy of three_body_planar.m, without modification and remove "axis equal" and uncomment the following lines:

dz = 1/sqrt((xrange(2)-xrange(1))^2+(yrange(2)-yrange(1))^2);
daspect ([1, 1, dz]);

Not showing the lines isn't a smoothing artifact, since

set (gcf, 'graphicssmoothing', 'off');

does not make the lines appear (references: bug #57336, bug #56979, bug #54554).

I still have a doubt it's a problem when (far+near)/(far-near) in glOrtho gets too big. Question: Is there a documentation on "axes::properties::update_camera", or do we have to analyze the source code? I want to understand the transformation matrices but it costs much time to go through the code.

Hg200 <hg200>
Fri 28 Aug 2020 04:11:03 PM UTC, comment #8: 

I made the change you described in comment #7. Here again as a patch to make it easier to find those lines:

diff -r 6dbd32dd2a5f libinterp/corefcn/gl-render.cc
--- a/libinterp/corefcn/gl-render.cc        Fri Aug 28 14:50:41 2020 +0200
+++ b/libinterp/corefcn/gl-render.cc        Fri Aug 28 17:47:08 2020 +0200
@@ -1279,8 +1279,10 @@

     Matrix x_zlim = props.get_transform_zlim ();

-    xZ1 = std::max (-1e6, x_zlim(0)-(x_zlim(1)-x_zlim(0))*100.0);
-    xZ2 = std::min (1e6, x_zlim(1)+(x_zlim(1)-x_zlim(0))*100.0);
+    // xZ1 = std::max (-1e6, x_zlim(0)-(x_zlim(1)-x_zlim(0))*100.0);
+    // xZ2 = std::min (1e6, x_zlim(1)+(x_zlim(1)-x_zlim(0))*100.0);
+    xZ1 = x_zlim(0)-(x_zlim(1)-x_zlim(0))/2;
+    xZ2 = x_zlim(1)+(x_zlim(1)-x_zlim(0))/2;

     Matrix x_mat1 = props.get_opengl_matrix_1 ();
     Matrix x_mat2 = props.get_opengl_matrix_2 ();

Without that change, I couldn't make the axes in your example appear in any test.
CC'ing Pantxo. Maybe he still remembers why that change was done in the first place for bug #54551.

With that change, I can the OpenGL errors (while the axes have disappeared). Octave crashes eventually for me, too.

If I add the following line after "axes equal" in your script, I see no OpenGL errors, the axes don't disappear and Octave doesn't crash:

  set(gca, 'CameraPosition', [0.014392034052167, 0.134014656653948, 2.411976487699734]*1e13);

That is the camera position Matlab calculates at that point in the script.

For reference, the camera position Octave calculates at that point in the script (if it isn't set manually):

>> format long
>> get(gca, 'CameraPosition')
ans =

   1.439510455551360e+11   1.340147276008146e+12   1.685352061831762e+13

It is again the z-component that differs. So there seems something off in how we calculate the camera position.

Markus Mützel <mmuetzel>
Group administrator
Mon 24 Aug 2020 05:09:53 PM UTC, comment #7: 

Thank you for your response. Regarding the segfault: I have made following changes to hg id f0414ee0fefe and made a fresh new build.

    //xZ1 = std::max (-1e6, x_zlim(0)-(x_zlim(1)-x_zlim(0))*100.0);
    //xZ2 = std::min (1e6, x_zlim(1)+(x_zlim(1)-x_zlim(0))*100.0);
    xZ1 = x_zlim(0)-(x_zlim(1)-x_zlim(0))/2;
    xZ2 = x_zlim(1)+(x_zlim(1)-x_zlim(0))/2;

The demo code that forces the segfault is attached to this report (run "three_body_planar.m" and wait about a minute). Forgive me, it's not a minimum working example, but it should be not to hard to read. Please notice that there are OpenGL error messages in the backtrace. I admit that the conclusion that the segfault is precisely caused by a poorly shaped projection matrix is ​​debatable. Could be also a problem how the OpenGl error is managed. But those OpenGl erros occur especially for large x,y settings at the same time we set the axis to equal.

For the other points mentioned in comment #6 i will respond later ;-)

Thread 1 "lt-octave-gui" received signal SIGSEGV, Segmentation fault.
0x00007fe790c82850 in vtable for charNDArray () from /home/.../libbgui/.libs/liboctgui.so.6
(gdb) bt
#0  0x00007fe790c82850 in vtable for charNDArray () at /home/.../libgui/.libs/liboctgui.so.6
#1  0x00007fe79094b1e8 in std::_Sp_counted_base<(__gnu_cxx::_Lock_policy)2>::_M_release() (this=0x7fe76493d690) at /usr/include/c++/9/bits/shared_ptr_base.h:155
#2  0x00007fe79094acd7 in std::__shared_count<(__gnu_cxx::_Lock_policy)2>::~__shared_count() (this=0x7ffe7cf5a5e8, __in_chrg=<optimized out>) at /usr/include/c++/9/bits/shared_ptr_base.h:730
#3  0x00007fe78fc75e8e in std::__shared_ptr<octave::stack_frame, (__gnu_cxx::_Lock_policy)2>::~__shared_ptr() (this=0x7ffe7cf5a5e0, __in_chrg=<optimized out>) at /usr/include/c++/9/bits/shared_ptr_base.h:1169
#4  0x00007fe78feab8e7 in std::__shared_ptr<octave::stack_frame, (__gnu_cxx::_Lock_policy)2>::operator=(std::__shared_ptr<octave::stack_frame, (__gnu_cxx::_Lock_policy)2>&&) (this=0x7ffe7cf5a640, __r=...) at /usr/include/c++/9/bits/shared_ptr_base.h:1265
#5  0x00007fe78feaae92 in std::shared_ptr<octave::stack_frame>::operator=(std::shared_ptr<octave::stack_frame>&&) (this=0x7ffe7cf5a640, __r=...) at /usr/include/c++/9/bits/shared_ptr.h:335
#6  0x00007fe78fea6cf5 in octave::call_stack::find_current_user_frame() const (this=0x7fe764004148) at libinterp/corefcn/call-stack.cc:488
#7  0x00007fe78fea71ad in octave::call_stack::backtrace_frames[abi:cxx11](long&) const (this=0x7fe764004148, curr_user_frame=@0x7ffe7cf5a818: -1) at libinterp/corefcn/call-stack.cc:636
#8  0x00007fe78fea7399 in octave::call_stack::backtrace_info[abi:cxx11](long&, bool) const (this=0x7fe764004148, curr_user_frame=@0x7ffe7cf5a818: -1, print_subfn=true) at libinterp/corefcn/call-stack.cc:673
#9  0x00007fe78fea7617 in octave::call_stack::backtrace_info[abi:cxx11]() const (this=0x7fe764004148) at libinterp/corefcn/call-stack.cc:695
#10 0x00007fe78fdf9cf8 in octave::tree_evaluator::backtrace_info[abi:cxx11]() const (this=0x7fe764004030) at libinterp/parse-tree/pt-eval.cc:2061
#11 0x00007fe78fdf9dcc in octave::tree_evaluator::backtrace_message[abi:cxx11]() const (this=0x7fe764004030) at libinterp/parse-tree/pt-eval.cc:2083
#12 0x00007fe78ff49c87 in octave::error_system::vwarning(char const*, char const*, char const*, __va_list_tag*) (this=0x7fe7640033a8, name=0x7fe7903b8a6a "warning", id=0x7fe7903b8a54 "", fmt=0x7fe7903bc588 "opengl_renderer: Error '%s' (%d) occurred drawing '%s' object", args=0x7ffe7cf5ade8)
    at libinterp/corefcn/error.cc:558
#13 0x00007fe78ff4a31c in octave::error_system::vwarning(char const*, char const*, __va_list_tag*) (this=0x7fe7640033a8, id=0x7fe7903b8a54 "", fmt=0x7fe7903bc588 "opengl_renderer: Error '%s' (%d) occurred drawing '%s' object", args=0x7ffe7cf5ade8) at libinterp/corefcn/error.cc:613
#14 0x00007fe78ff47a7e in vwarning(const char *, const char *, typedef __va_list_tag __va_list_tag *) (id=0x7fe7903b8a54 "", fmt=0x7fe7903bc588 "opengl_renderer: Error '%s' (%d) occurred drawing '%s' object", args=0x7ffe7cf5ade8) at libinterp/corefcn/error.cc:115
#15 0x00007fe78ff4d81e in vwarning(char const*, __va_list_tag*) (fmt=0x7fe7903bc588 "opengl_renderer: Error '%s' (%d) occurred drawing '%s' object", args=0x7ffe7cf5ade8) at libinterp/corefcn/error.cc:1069
#16 0x00007fe78ff4d8c0 in warning(char const*, ...) (fmt=0x7fe7903bc588 "opengl_renderer: Error '%s' (%d) occurred drawing '%s' object") at libinterp/corefcn/error.cc:1077
#17 0x00007fe78ff9ad2d in octave::opengl_renderer::draw(graphics_object const&, bool) (this=0x16d1210, go=..., toplevel=true) at libinterp/corefcn/gl-render.cc:772
#18 0x00007fe78ffa5483 in octave::opengl_renderer::draw_axes_children(axes::properties const&) (this=0x16d1210, props=...) at libinterp/corefcn/gl-render.cc:2258
#19 0x00007fe78ffa5c97 in octave::opengl_renderer::draw_axes(axes::properties const&) (this=0x16d1210, props=...) at libinterp/corefcn/gl-render.cc:2346
#20 0x00007fe78ff9a0d0 in octave::opengl_renderer::draw(graphics_object const&, bool) (this=0x16d1210, go=..., toplevel=false) at libinterp/corefcn/gl-render.cc:730
#21 0x00007fe78ffb5010 in octave::opengl_renderer::draw(Matrix const&, bool) (this=0x16d1210, hlist=..., toplevel=false) at libinterp/corefcn/gl-render.cc:4082
#22 0x00007fe78ff9b3a3 in octave::opengl_renderer::draw_figure(figure::properties const&) (this=0x16d1210, props=...) at libinterp/corefcn/gl-render.cc:797
#23 0x00007fe78ff9a019 in octave::opengl_renderer::draw(graphics_object const&, bool) (this=0x16d1210, go=..., toplevel=true) at libinterp/corefcn/gl-render.cc:728
#24 0x00007fe79096af59 in QtHandles::GLCanvas::draw(octave_handle const&) (this=0x16d1180, gh=...) at libgui/graphics/GLCanvas.cc:87
#25 0x00007fe790954779 in QtHandles::Canvas::canvasPaintEvent() (this=0x16d11b0) at libgui/graphics/Canvas.cc:284
#26 0x00007fe79096c2a2 in QtHandles::GLCanvas::paintGL() (this=0x16d1180) at libgui/graphics/GLCanvas.cc:241
#27 0x00007fe78e3ce835 in QOpenGLWidgetPrivate::invokeUserPaint() (this=0x13ad020) at ../../include/QtGui/../../src/gui/kernel/qopenglcontext.h:152

(file #49702)

Hg200 <hg200>
Mon 24 Aug 2020 08:46:26 AM UTC, comment #6: 

Very good detective work.

I haven't tried your patch yet.

Tbh, I don't understand the reason for the segfault you see without the min and max. Why do we need these (arbitrary) limits? Where does the segfault occur?
The min and max were introduced to solve bug #54551. Maybe that can give a clue:

What is the reason for the segfault? Just a wild guess: Must the difference between the z-clipping limits be larger that a certain fraction of its average/minimum/maximum value?
The current approach limits the distance between the clipping planes. The new approach extends it. Why do both not cause a segfault, but removing min and max does?

Is there a reason the proposed changeset breaks the existing symmetry between lower and upper clipping limits?
Would it make sense to extent the z-clipping limits symmetrically instead?

"xZ2" is the same in both branches. Maybe it could be moved outside the conditional blocks?

The z-component of the camera position after running the first example in comment #0 differs by a factor of approx. 20 in Matlab and Octave.
Matlab R2020a:

>> plot (0, 0, 'ro');
daspect ([1, 1, 1]);
ylim ([-1e5, 3e5]);
xlim ([-1e5, 1e5]);
get(gca, 'CameraPosition')

ans =

   1.0e+05 *

         0    1.0000    1.1180

Octave 5.2.0:

>> plot (0, 0, 'ro');
>> daspect ([1, 1, 1]);
>> ylim ([-1e5, 3e5]);
>> xlim ([-1e5, 1e5]);
>> get(gca, 'CameraPosition')
ans =

         0.00000    100000.00000   2236068.47751

That might be the reason why the original example exposes this bug (like you already wrote in comment #3).

But I also think that we need to make sure that we don't clip important parts from the scene (in particular not the plain with the axes in 2d plots).
IIUC, xZ2 (farVal) should not be smaller than xZ1 (nearVal). That happens in the original example (see comment #3). So there most probably is an error with the current logic.

Markus Mützel <mmuetzel>
Group administrator
Sat 22 Aug 2020 04:59:33 PM UTC, comment #5: 

The answer is probably yes and no. Using min, max for the calculation of xZ1, xZ2 is most likely motivated by the fact that the projection matrix becomes ​​undefined for a large object distance over a small clipping plane gap. For example if we remove the min, max statements and enter large xlim, ylim values, Octave segfaults:


However with this min, max implementation we run into multiple problems: For example the scene can become invisible. But it is also possible to jam the renderer if xlim, ylim is chosen unfortunate. The reason for the latter case is that the difference between the truncated far and near value can become zero. This is not allowed. Example:

daspect ([1, 1, 1]);
xlim ([-1e5, 1e5]);
ylim ([-1e5, 2.46456585e5]);
warning: opengl_renderer: Error 'invalid value' (1281) occurred drawing 'text' object

The attached test code does circumvent the above problems, but as mentioned in comment #3 it does not solve the displaying issue for larger axis settings. I guess at this point OpenGl limitations come in.

Comments are welcome.

(file #49697)

Hg200 <hg200>
Sun 16 Aug 2020 08:28:10 PM UTC, comment #4: 

Could this be another manifestation of the notorious single-precision-in-OpenGL bug?
Esp. because the gnuplot toolkit does it right.

Philip Nienhuis <philipnienhuis>
Group Member
Sun 16 Aug 2020 05:03:37 PM UTC, comment #3: 

Okay, then i think it's a bug related to plotting, because in a numerical tool like Octave one can't expect axis limits to always be less than the order of +-3e5.

1.) If we switch to

graphics_toolkit ('gnuplot')

everything is displayed correctly. If we switch to fltk, similar problem occurs:

graphics_toolkit ('fltk')

2.) When we compare the output of

get (gca)

between Octave and Matlab R2018, the camera position and some other properties differ. It looks like Octave rather correlates with Matlab R2013 than R2018. Anyway, I believe there's basically nothing wrong with it, although the camera position has an influence on the bug.

3.) I've searched via trial and error. For example in opengl_renderer::setup_opengl_transformation, the near and far values ​​for the example from comment 0 are at the clipping edge:

xZ1 = std::max (-1e6, x_zlim(0)-(x_zlim(1)-x_zlim(0))*100.0);
xZ2 = std::min (1e6, x_zlim(1)+(x_zlim(1)-x_zlim(0))*100.0);
m_glfcns.glOrtho (0, vp(2), vp(3), 0, xZ1, xZ2);
(gdb) print *vp.rep.data@vp.rep.len
$1 = {0, 0, 560, 420}
(gdb) print *x_zlim.rep.data@x_zlim.rep.len
$2 = {1118033.4887610751, 1118034.4887610751}
(gdb) print xZ1
$4 = 1117933.4887610751
(gdb) print xZ2
$3 = 1000000

glOrtho produces a parallel projection. If we allow a wider near/far interval it seems to solve the problem, although it does not solve the problem for all variations of larger axis (i.e. 1e12). This and also axes::properties::update_camera look like a starting point to me.

Hg200 <hg200>
Fri 14 Aug 2020 06:30:28 PM UTC, comment #2: 

can reproduce with af302efce502 (stable)


Dmitri A. Sergatskov <dasergatskov>
Fri 14 Aug 2020 06:27:26 PM UTC, comment #1: 

Yes I can reproduce the problem.
Ubuntu 18.4 and octave 5.2.1

Doug Stewart <dastew>
Fri 14 Aug 2020 06:12:22 PM UTC, original submission:  

The following code creates an empty plot. I.e. the axes and additional graphic objects are not displayed. Affects Linux development branch, Windows MXE Octave 5.2.0 and an intermediate version of the stable MXE Octave branch.

In Matlab 2018, the "dot" and the axis are visible:

plot (0, 0, "ro");
daspect ([1, 1, 1]);
ylim ([-1e5, 3e5]);
xlim ([-1e5, 1e5]);

We can "heal" the plot by doing the following adjustment:

plot (0, 0, "ro");
daspect ([1, 1, 1e6]);
ylim ([-1e5, 3e5]);
xlim ([-1e5, 1e5]);

Sometimes but not always this also "heals" the plot:

plot (0, 0, "ro");
daspect ([1, 1, 1]);
ylim ([-1e5, 1e5]);
xlim ([-1e5, 1e5]);

And in this variant, the plot is empty too:

hold on;
axis equal;
ylim ([-1e5, 3e5]);
xlim ([-1e5, 1e5]);
plot (0, 0, "ro");

I don't get any gl-renderer warnings.

Please: Can anybody can reproduce?

Hg200 <hg200>


(Note: upload size limit is set to 16384 kB, after insertion of the required escape characters.)

Attach Files:

Attached Files
file #53065:  comparison.png added by hg200 (25KiB - image/png)
file #49794:  bug58956_gl_clipping_plains.patch added by mmuetzel (2KiB - application/octet-stream)
file #49750:  opengl_test.tar added by hg200 (10KiB - application/x-tar)
file #49702:  three_body_planar.m added by hg200 (3KiB - text/x-objcsrc)
file #49697:  gl-render_tests.diff added by hg200 (948B - text/x-patch)


Depends on the following items: None found

Items that depend on this one: None found


Carbon-Copy List
  • -email is unavailable- added by mmuetzel
  • -email is unavailable- added by mmuetzel (Posted a comment)
  • -email is unavailable- added by dasergatskov (Posted a comment)
  • -email is unavailable- added by dastew (Posted a comment)
  • -email is unavailable- added by hg200 (Submitted the item)

    There are 0 votes so far. Votes easily highlight which items people would like to see resolved in priority, independently of the priority of the item set by tracker managers.

    Only group members can vote.


    Follow 9 latest changes.

    Date Changed by Updated Field Previous Value => Replaced by
    2022-04-10 hg200 Attached File- Added comparison.png, #53065
    2020-11-01 mmuetzel StatusPatch Submitted Patch Reviewed
    2020-09-27 hg200 Attached File- Added patch_camera_target_experimental.diff, #49866
    2020-09-13 mmuetzel Attached File- Added bug58956_gl_clipping_plains.patch, #49794
    2020-09-13 mmuetzel StatusNone Patch Submitted
    2020-09-07 hg200 Attached File- Added opengl_test.tar, #49750
    2020-08-28 mmuetzel Carbon-Copy- Added pantxo
    2020-08-24 hg200 Attached File- Added three_body_planar.m, #49702
    2020-08-22 hg200 Attached File- Added gl-render_tests.diff, #49697

    Back to the top

    Powered by Savane 3.13-6ae7.
    Corresponding source code