Skip to content

Introduce indexed ray transfer APIs and tests#503

Open
munechika-koyo wants to merge 9 commits into
cherab:developmentfrom
munechika-koyo:feature/add-new-raytransfer
Open

Introduce indexed ray transfer APIs and tests#503
munechika-koyo wants to merge 9 commits into
cherab:developmentfrom
munechika-koyo:feature/add-new-raytransfer

Conversation

@munechika-koyo

@munechika-koyo munechika-koyo commented Jun 18, 2026

Copy link
Copy Markdown
Member

Summary

This PR introduces a general index-function-based ray transfer API by adding indexed emitter and integrator classes, replacing mesh-specific naming with functionality-based naming.

Key Changes

Unit Test

Test Purpose Setup Verification Why it matters
test_evaluate_function Validate that IndexedRayTransferEmitter works correctly with NumericalIntegrator and maps contributions to the correct bins via index_function. A 3x3x3 Box domain is used with bins=27. A diagonal ray crosses the volume. The index function maps in-domain points to 0..26 and returns -1 outside. Only bins 0, 13, and 26 are non-zero, each with path-length contribution sqrt(3). The output spectrum matches the expected vector with atol=0.001. Confirms correct geometric integration and bin assignment for the index-function-based workflow.
test_default_integrator Confirm default integrator behavior when no integrator is explicitly passed. IndexedRayTransferEmitter is created without an integrator argument under the same ray/volume setup as above. The emitter uses IndexedRayTransferIntegrator by default, and the resulting spectrum matches the same expected vector (atol=0.001). Guarantees safe default behavior and avoids mandatory integrator wiring for users.
test_discrete3dmesh_as_index_function Validate that Discrete3DMesh can be used as an equivalent index function source. A Discrete3DMesh is built from a 4x4x4 vertex grid over 3x3x3 cells; each cube is split into 6 tetrahedra. Cell values follow the same indexing rule as the reference index function. Representative points across all 27 cells match the reference mapping; outside-domain points return -1; ray-transfer spectra from mesh-based and function-based indexing are equal within atol=0.001. Demonstrates implementation-agnostic design and compatibility with tetrahedral mesh indexing in practical ray-transfer use.

Executed:

python -m unittest cherab.tools.tests.test_raytransfer.TestIndexedRayTransferEmitter -v

Result:

  • test_default_integrator: ok
  • test_discrete3dmesh_as_index_function: ok
  • test_evaluate_function: ok
  • Ran 3 tests, all passed.

Example Usage

from raysect.optical import World
from cherab.tools.raytransfer import IndexedRayTransferEmitter

def index_func(x, y, z):
    if x < 0:
        return 0
    return 1

world = World()
material = IndexedRayTransferEmitter(index_func, bins=2)

Compatibility and Risk

  • Scope is limited to ray transfer emitter/integrator API naming and related tests.
  • Runtime behavior is validated by focused unit tests for indexed evaluation and Discrete3DMesh integration.

Reviewer Notes

Checklist

  • API implementation updated
  • Cython declaration file updated
  • Unit tests added/updated
  • Targeted tests executed in pixi test environment

Benchmark

The appendix in this paper (https://doi.org/10.1063/5.0225703) compared the raytransfer of Discrete3D meshes with that of regular grids, showing that the geometry matrix calculation for rectangular grids' raytransfer was much faster than the Discrete3D one.

@jacklovell jacklovell left a comment

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Very cool. If I understand this correctly, it's a generalisation of the existing ray transfer objects where a callable function replaces the voxel map, and this function returns the voxel index at a given point in space. As well as the 3D application in https://doi.org/10.1063/5.0225703 (which you should cite somewhere in the documentation by the way), I can see the possibility of application to axisymmetric voxels of arbitrary poloidal cross section without having to approximate them with a rectangular grid.

I think a demo would be highly beneficial, as it's hard to see a concrete use case from the docstrings alone. Perhaps you could adapt Vlad's Space Invaders demos to showcase the new tools.

Comment thread cherab/tools/raytransfer/emitters.pyx
Comment thread cherab/tools/raytransfer/emitters.pyx
@munechika-koyo

Copy link
Copy Markdown
Member Author

Thank you for your review!
I've been thinking about what makes an effective demo script.
Currently, my ideas are using index functions like:

  • analytical geometry with sin, cos, circle, etc.
  • using mesh geometry (triangulation for 2D, tetrahedralization for 3D, e.g., Stanford Bunny data used in raysect).

I will try to work out these ideas with AI assistance for now.
Any idea of yours is really helpful.

@jacklovell

Copy link
Copy Markdown
Member

Analytically-defined voxels will make a nice demo, agreed. For example, flux-aligned axisymmetric voxels for the generomak equilibrium defined in ($\Delta \psi$, $\Delta \beta$) space for normalised flux $\Psi$ and poloidal angle $\beta$. The ToroidalVoxelGrid could handle this but I expect the ray transfer framework to be more performant.

And yes, voxels defined by triangular or tetrahedral meshes would also be a good illustration, showing the ability to extend beyond regular rectilinear grids.

@Mateasek Mateasek left a comment

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Hi @munechika-koyo , this is a very nice generalisation of the RT framework! I have one important question to raise.

Changing the indexing function from Function3D to Function6D would allow to pass also direction information which would make the new IRT framwrok applicable to also anisotropic radiation.

integrator = integrator or IndexedRayTransferIntegrator(step=integration_step)
super().__init__(integrator=integrator)

self.index_function = autowrap_function3d(index_function)

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Making index_function a Function3D makes the IndexedRay framework specific to isotropic radiation. If you made it Function6D which now we have in Cherab, the framework would be generalised to anisotropic applications because you could also pass direction vector components. That is something I'm be very interested in. From the point of view of the code it shouldn't be a large change. What do you think @munechika-koyo?

@munechika-koyo munechika-koyo Jul 22, 2026

Copy link
Copy Markdown
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

That sounds really interesting!
It might also be useful for distributional tomography.
I think we need to add a new API or refactor the existing one to accept Function6D as an index function for the IndexedRayTransferEmitter.emission_function() and IndexedRayTransferIntegrator.integrate() methods.
These changes don't seem straightforward to me, so should we handle them in a separate PR after finalizing this one?

In addition, any good ideas for your demo using Function6D would be really helpful.

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Maybe I'm missing something, but what else do you need to change except the Function3D calls, i.e. in IndexedRayTransferIntegrator.integrate you change <int>index_function(x, y, z) to <int>index_function(x, y, z, dx, dy, dz) and then you need to change the object type and etc.. Maybe I'm missing something.

@Mateasek Mateasek mentioned this pull request Jul 22, 2026
12 tasks

@Mateasek Mateasek left a comment

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thanks @munechika-koyo for your reaction. It made me look at the Raytransfer framework previously added by @vsnever and the IndexedRaytransfer framework you did in this PR from a broader perspecitve.

Let me now express how I see the framework and what it should do. I think that RayTransfer in general is there to discretise continuous space into one dimensional integer indexes. Then it uses the integration to add a sensitivity to the integer bins.

What does the current RayTransfer framework does? Its specified by a single aspect:

  • The mapping between a cartesian space and a 1D array of integers. There are now two distinct mappings for which two sets of classes were made. That is cylindrical and cartesian mappings.

Now what your contribution does is, and correct me if I'm wrong, that it allows user to define the mapping function. This from my point of view generalises the previous approach. Now, we can have a single set of integrator and emitter, which can be used for the old raytransfer functionality if you provide the right discretising mapping functions.

Now I'm going to propose what I think should be done, but please tell me if I'm missing something and I would also like @jacklovell and @skuba31 to think about this, since it would be a major change:

  1. The IndexedRaytransfer approach I described above will replace the older Raytransfer because essentially it can do exactly the same thing and we should not be doubling our functionality.
  2. The most general indexing function (and replacements for the old framework) together with inverse indexing funcitons will be added to tools to give users the possibility to do some basic and most applications out of the box.
  3. Extend the indexing function to accept Function6D to make this framework even more general.

I'm sorry to be proposing such a major overhaul of this contribution but there wasn't any previous discussion in an issue, so there was no place to do it. I actually had the approach I described above in my head for a long time and I was very happy to see that you actually did it @munechika-koyo, but I think it could be done in a more general and systematic way which would be actually very powerful. If you think about it, it will allow users to apply discretisation in most of the "information dimensions" Raysect passes to the integrator with the ray which position and direction. This would allow to do anisotropic contribution matrices for example. Or you could decide that you implement a discretisation function which uses x, y, z and wavelength (and you just don't use 2 parameters in the 6D function evaluation).What do you think?

integrator = integrator or IndexedRayTransferIntegrator(step=integration_step)
super().__init__(integrator=integrator)

self.index_function = autowrap_function3d(index_function)

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Maybe I'm missing something, but what else do you need to change except the Function3D calls, i.e. in IndexedRayTransferIntegrator.integrate you change <int>index_function(x, y, z) to <int>index_function(x, y, z, dx, dy, dz) and then you need to change the object type and etc.. Maybe I'm missing something.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

3 participants