cherab / core

The core source repository for the Cherab project.
https://www.cherab.info
Other
45 stars 24 forks source link

How to measure emission of plasma using voxels and the sensitivity matrix with bolometers? #470

Closed leferi99 closed 4 weeks ago

leferi99 commented 1 month ago

I am following the demo "Performing Inversions of Bolometer Measurements Using the Voxel Framework". I have the voxels and sensitivity matrix generated already, and it is stated in the demo that:

If many emission profiles are being inverted for the same bolometer and voxel grid geometries, it would be substantially faster to sample the emission function on the voxel grid and then multiply this by the sensitivity matrix (this is done in this demo script to produce the plot of the phantom).

And the part in the parenthesis refers to this line I believe:

phantom_samples = voxel_grid.emissivities_from_function(emission_function_3d)

Now what if my emission function is not predefined this way, but instead I have a radiating plasma? How to do the voxel observation in that case? I didn't find a suitable method for the Plasma to return a 3D emission function and neither did I find a suitable method for ToroidalVoxelGrid or VoxelCollection to get the emissivity from a plasma.

vsnever commented 1 month ago

Hi @leferi99, You can create a Python 3D function, which returns isotropic, spectrally-unresolved emission from plasma. Suppose you already have a Plasma object and the emission models connected to it.

from raysect.core.math import Point3D, Vector3D
from raysect.optical import Spectrum

...

spectrum = Spectrum(min_wavelength, max_wavelength, bins)
direction = Vector3D(0, 0, 1)

def emission_function_3d(x, y, z):
    emission = 0
    for model in plasma.models:
        point = Point3D(x, y, z)
        emission += model.emission(point, direction, spectrum.new_spectrum()).total()

    return emission

phantom_samples = voxel_grid.emissivities_from_function(emission_function_3d)

Note that this function returns emission in Plasma coordinate space, so perform transformations if necessary. The units are W m-3 sr-1. Also, the emission returned by the model must be axisymmetric, because voxel_grid.emissivities_from_function() will use only x and z coordinates.

leferi99 commented 1 month ago

Thank you @vsnever! I didn't specify but I would need this spectrally resolved. If I leave out the .total() at the end of model.emission(...) would that work?

(I am using BolometerFoil objects for the observation but I specify the pipeline to be SpectralPowerPipeline to be able to model AXUV diodes and the degradation of their spectral response.)

vsnever commented 1 month ago

voxel_grid.emissivities_from_function() wraps the emission_function_3d into Raysect's Function3D which returns a scalar, but you can sample spectral emission in a loop:

from raysect.core.math import Point3D, Vector3D
from raysect.optical import Spectrum

...

spectrum = Spectrum(min_wavelength, max_wavelength, bins)
direction = Vector3D(0, 0, 1)

phantom_samples = np.zeros(bins)

for i in range(bins):

    def emission_function_3d(x, y, z):
        emission = 0
        for model in plasma.models:
            point = Point3D(x, y, z)
            emission += model.emission(point, direction, spectrum.new_spectrum()).samples[i]

        return emission

    phantom_samples[i] = voxel_grid.emissivities_from_function(emission_function_3d)

Note that the units of the emission_function_3d will be W m-3 sr-1 nm-1 in this case.

leferi99 commented 1 month ago

Okay, thank you, I will test this, and if I don't have more questions, I will close this issue in a couple of days.