Add early-termination optimization to attenuated mip

[aganders3]

This idea originated from a thread on image.sc about poor framerates for volume rendering napari. The idea is to exit the main raycasting loop early if there is no chance of finding a higher value based on the provided contrast limits. I believe this will not affect the final appearance at all and is only a performance optimization (and my testing supports this).

I don't have a realistic dataset where this makes a truly meaningful difference, but for now you can try this example (adapted from examples/scene/volume.py) with a huge cube of random data with/without this change and the performance difference should be apparent.

import numpy as np

from vispy import app, scene

# Prepare canvas
canvas = scene.SceneCanvas(size=(800, 600), show=True)
canvas.measure_fps()

# Set up a viewbox to display the image with interactive pan/zoom
view = canvas.central_widget.add_view()

# Create the volume visual
volume = scene.visuals.Volume(
    np.random.random((1024, 1024, 1024)).astype(np.float32),
    parent=view.scene,
    attenuation=0.05,
    method="attenuated_mip",
)

view.camera = scene.cameras.TurntableCamera(
    parent=view.scene,
    fov=60,
    name='Turntable',
)


if __name__ == '__main__':
    app.run()

N.B. @brisvag

[jni]

@aganders3 what do you think of my comment on the imagesc thread that the contrast limits should actually apply to the attenuated value, rather than the original data? Basically, my idea is that the attenuation should affect which voxel gets rendered from the ray, and then the mapping of that voxel value to a color gets put through the colormap. Am I thinking about this wrong?

[aganders3]

@aganders3 what do you think of my comment on the imagesc thread that the contrast limits should actually apply to the attenuated value, rather than the original data? Basically, my idea is that the attenuation should affect which voxel gets rendered from the ray, and then the mapping of that voxel value to a color gets put through the colormap. Am I thinking about this wrong?

Yes I want to look into this. I didn't want to hijack the image.sc thread, so was going to start a Zulip chat about it but just haven't made the time yet.

Have you tried it since the changes in #2417? I feel this mode is much more predictable since then but there's surely room for improvement. I still think it makes sense to scale/clamp when calculating the attenuation such that the attenuation is not dependent on the scale of your data (in particular when you have negative values this does not make sense) but I see where you're coming from as well. I'm going to try a few things.

[aganders3]

Maybe I'm getting it wrong, but here is a comparison of what I think you mean (using the value instead of scaled value in the output). This causes some odd artifacts. But you're right that it prevents needing to bump up the lower contrast limit as you increase attenuation. I'll have to think about that a bit to see if there's a way to get the best of both.

scaled (current)

unscaled

[jni]

@aganders3 Interesting, thanks for the link! Totally missed that PR. I'm convinced that I have no idea what the right thing to do is, carry on! 😂

[brisvag]

Love this :)

I believe this will not affect the final appearance at all and is only a performance optimization (and my testing supports this).

You could do a few buffer captures at different levels of attenuation and compare the result between the two implementations to actually make sure this is true. But in principle I agree that it should be the same.

[djhoese]

As far as I can tell the failures are caused by the new ipython/jupyter comm package being released and/or used before it was ready (NotImplementedError). The code has changed quite a bit in the last 24 hours or so: https://github.com/ipython/comm

[aganders3]

You could do a few buffer captures at different levels of attenuation and compare the result between the two implementations to actually make sure this is true. But in principle I agree that it should be the same.

Thanks for the idea - I added this to my demo script and I'm satisfied it's the same 🙂

As far as I can tell the failures are caused by the new ipython/jupyter comm package being released and/or used before it was ready (NotImplementedError). The code has changed quite a bit in the last 24 hours or so: https://github.com/ipython/comm

Anything I can help with or should it just be re-run? It looks like the test used comm=0.1.2 and now 0.1.3 is released and has a relevant fix: ipython/comm#13

[brisvag]

Thanks for the idea - I added this to my demo script and I'm satisfied it's the same slightly_smiling_face

You could even compare the literal pixel values, they're numpy arrays after all :P That's how some of the older tests in vispy are handled in fact. But I don't think it's necessary, to be clear!

[djhoese]

Restarted failing job. Let's see what happens.

[aganders3]

Now that I think about it - the same optimization can be applied to non-attenuated MIP for any fragment that's going to be maxed out (maxval >= clim.y). This is probably less impactful, but may still be worth including. I guess it's already in place here if you set attenuation = 0 to get a standard MIP, so that gives a relative simple way to roughly compare performance (with the overhead that attenuation is still calculated unnecessarily for each step).

[jni]

This is probably less impactful,

I disagree, because MIP is the default rendering that most people use. I expect that the overhead of attenuation is actually significant and you should make this optimisation! (I thought that was already in place...)

[psobolewskiPhD]

Running the example takes considerable time before anything happens (with a 10Gb+ python process spawned), but the difference in performance with this PR is astounding.
Without I get max 4.9 fps, with 60 fps all the time.

[aganders3]

I disagree, because MIP is the default rendering that most people use. I expect that the overhead of attenuation is actually significant and you should make this optimisation! (I thought that was already in place...)

I mean less impactful in that a MIP with reasonable contrast limits will probably not have too many saturated pixels, where the attenuated MIP might have many that reach attenuation. In any case I added this + the same for minip. I also added a flag to skip the fine-step refinement in these cases.

Running the example takes considerable time before anything happens (with a 10Gb+ python process spawned), but the difference in performance with this PR is astounding.
Without I get max 4.9 fps, with 60 fps all the time.

Yes - sorry for the enormous example :p. The random noise likely shows the effect much more than real data, but it is a striking difference!

[aganders3]

Oops, 9e61444 was bad, but I think it's fixed in 3cb973a. With test data the buffer is still identical with these changes, but there could still be some difference with the depth as it's now sometimes returning on the first value above/below the threshold along the ray. I don't think this would affect the output image but might be missing something.

[jni]

there could still be some difference with the depth as it's now sometimes returning on the first value above/below the threshold along the ray.

I think this is actually what users would expect, to be honest?

[brisvag]

As far as I can tell this is ready to go in! @djhoese, should we get this in the next release?

[brisvag]

Sorry for the last minute doubt; I don't think it's warranted to add that refine flag, unless I misunderstand how it works?

[brisvag]

Sorry, forgot to update the review state, LGTM!

[brisvag]

@djhoese feel free to merge if this is good for you!

[djhoese]

LGTM. Thanks @aganders3!