Development of a tomographic reconstruction method for axisymmetric D_α emission profiles in the ITER plasma boundary

Strong stray light produced by reflections from the all-metal first wall (FW) is the main problem for the ITER Main Chamber H-alpha and Visible Spectroscopy diagnostic systems [1, 2]. For hydrogen spectral lines, the source of the stray light is the emission (in the multi-MW range for high power, burning plasma operation) from these spectral lines in the divertor plasma. Predictive plasma boundary simulations of the stationary phases of ITER Q_DT = 10 baseline discharges obtained with a combination of the SOLPS-4.3 [3, 4] and OEDGE [5] codes (the latter being used to extend the SOLPS-4.3 computational grid out to all plasma-facing surfaces), show that the ratio of the D_α power emitted in the main chamber SOL, $P_{{\mathrm{D\alpha }}}^{{\text{SOL}}}$ to that emitted in the divertor $P_{{\mathrm{D\alpha }}}^{{\text{Div}}}$, can vary in the range 0.005 to 0.1. This wide range reflects the uncertainty of the far SOL plasma specification for ITER, itself a reflection of the current far from complete understanding of the cross-field transport in this region. If the broad density shoulders seen on some tokamaks (e.g. JET [6] and ASDEX Upgrade [7]) under the partially detached divertor conditions in which ITER must operate at high power are seen on ITER, then levels of main chamber recycling may lead to higher emission there, increasing the ratio of $P_{{\mathrm{D\alpha }}}^{{\text{SOL}}}$/$P_{{\mathrm{D\alpha }}}^{{\text{Div}}}$. If the far SOL density is much lower, then the ratio will decrease rapidly. This is captured in the simulations used here to generate synthetic diagnostic signals by varying the radial convective velocity, v_perp of the far SOL particle flux in OEDGE (as a proxy for SOL filamentary transport).

A short description of the SOLPS + OEDGE simulation scenarios is given in table 1, which complements table 1 in [2] with three additional columns. The comparison of D_α emission profiles in the two uttermost scenarios 'd' and 'o' is given in figure 1 in [2]. Depending on the simulation scenario and the position of the line of sight (LoS) in the main chamber, the ratio of the useful signal to the stray light varies from 0.001 to 2. To take into account a strong stray light and a strong deviation of neutral atom velocity distribution function from a Maxwellian in the measurement interpretation, a H-alpha Spectroscopy synthetic diagnostic [8] was developed. Now a new version of this synthetic diagnostic which deploys the Raysect ray-tracing library [9] and the Cherab spectroscopy modelling framework [10] for predictive modelling of the stray light, is under development. The Raysect and Cherab codes are complemented with modules for fast volume integration along the ray trajectories over regular spatial grids and fast calculation of ray transfer matrices (geometry matrices): see below and in [11]. With these additions, Raysect and Cherab perform ray tracing with inhomogeneous volumetric light sources defined on a very fine spatial grid, hundreds of times faster than general-purpose illumination design software packages such as LightTools or Zemax OpticStudio previously used in [2, 11, 12].

Zoom In Zoom Out Reset image size

Table 1. Main parameters of SOLPS + OEDGE simulation scenarios from ITER database.

Case	Pedestal type	Far SOL vperp( m s−1)	Far SOL Te( eV)	$P_{{\mathrm{D\alpha }}}^{{\text{Div}}}\,\left( {{\text{kW}}} \right)$	$P_{{\mathrm{D\alpha }}}^{{\text{SOL}}}{\mathrm{\,}}\left( {{\text{kW}}} \right)\,$	$P_{{\mathrm{D\alpha }}}^{{\text{SOL}}}/P_{{\mathrm{D\alpha }}}^{{\text{Div}}}\,$
d	L-mode-like	30	10	129.2	0.693	5.37 × 10–3
e	H-mode-like	30	10	129.6	0.708	5.46 × 10–3
f	L-mode-like	30	20	129.7	0.761	5.88 × 10–3
g	H-mode-like	30	20	130.0	0.780	6.00 × 10–3
h	L-mode-like	65	10	129.7	1.111	8.56 × 10–3
i	H-mode-like	65	10	129.9	1.137	8.75 × 10–3
j	L-mode-like	65	20	130.2	1.230	9.45 × 10–3
k	H-mode-like	65	20	130.4	1.250	9.59 × 10–3
l	L-mode-like	100	10	137.4	13.77	1.00 × 10–1
m	H-mode-like	100	10	138.0	13.76	9.97 × 10–2
n	L-mode-like	100	20	140.1	13.52	9.65 × 10–2
o	H-mode-like	100	20	140.6	13.48	9.59 × 10–2

Here a tomographic method for reconstructing the axially symmetric D_α emission profile in the ITER SOL and divertor from D_α intensity measurements in the fields of view (FoV) of the ITER divertor impurity monitor (DIM) and Vis/IR TV diagnostics is developed. The method was originally proposed in [11], where it was used for reconstruction of the D_α emission profile in the divertor and predictive modelling of divertor stray light on the LoS in the main chamber observing the central column. Because of low performance of LightTools when applied to this task, the results in [11] were obtained for simplified D_α emission profile.

When the D_α emission profile is reconstructed, it can be used to calculate the useful signal on the diagnostic receivers (LoS arrays, filtered cameras with a wide FoV, etc.), that is needed to assess the hydrogen recycling, namely the influx of neutral hydrogen isotopes in the form of atoms and molecules from the FW into the plasma.

The ray transfer matrix (or geometry matrix) provides transformation of the signal from each individual unit light source to each individual receiver (CCD pixel, etc.), based on a ray-tracing simulation (see [11] for detailed explanation). It acts as a Green's function in similar tasks [13]. Ray transfer matrices can be used for reconstructing the emission profile, using the signals from a finite number of detectors which observe the same area at different angles. The inverse problem for reconstructing the emission profile is the following:

$$\begin{equation}\left\| {{\mathbf{Ax}} - {\mathbf{b}}} \right\|\mathop \to \limits_{\mathbf{x}} {\text{min}},\end{equation} \tag{ 1 }$$

where A is the ray transfer matrix, b is the vector of the signals on individual receivers and x is the sought after vector of emissivities of individual emitters (light sources), which form the spatial profile of emission.

Recently, the ray transfer matrix approach has been used to reconstruct the two-dimensional D_α emission profile in the JET ITER-like Wall (JET-ILW) divertor using the image from the KL11 filtered divertor camera [14, 15]. In this tomographic approach, the bidirectional reflectance distribution function (BRDF) of the FW plasma-facing surfaces is assumed to be known. The problem of obtaining this BRDF has not yet been resolved and a few options are currently under discussion among plasma diagnostic groups. Direct measurement of the BRDF of the FW components before the FW is assembled is helpful but may be insufficient since the light reflection properties are likely to change after the following exposure to plasma. It is likely that some model BRDF with a few free parameters will be applied. These free parameters are to be recovered by fitting the camera images obtained for calibrated light sources (see Figure 6 in [14]) with the ray-traced images. Alternatively, instead of using calibrated light sources, a pulse with a strongly localized emission distribution (e.g. the neon distribution during neon seeding in the divertor) can be performed.

In ITER, reconstruction of D_α emission spatial profile can be performed in three steps.

• (a)  
  Ray transfer matrices are calculated for all the light sources and all the receivers.
• (b)  
  Emissivity of the light sources in divertor is reconstructed using the LoS arrays of the ITER DIM diagnostic and certain segments of images from the filtered cameras of the ITER Vis/IR TV diagnostic.
• (c)  
  The obtained solution for divertor is used as the initial guess to reconstruct the emissivity for all light sources, now using both the DIM receivers and all available filtered camera images of Vis/IR TV.

The divertor is the main source of radiation in the vacuum chamber. However, it is observed with a relatively low number of receivers. Therefore, reconstructing the emission first in divertor, before the rest of the SOL, provides more accurate results than if the emission is reconstructed everywhere in a single step. If the BRDFs of the FW components were known exactly, one could simply subtract the divertor stray light from the measured signals after Step b and then reconstruct the emission in the main chamber SOL using the residual. However, in so far as the BRDFs are known with some errors, subtracting the divertor stray light may lead to a negative residual for some segments of the images. Step c of the algorithm prevents this from happening. Namely if the solution in divertor results in divertor stray light which is higher than the actual signal, this solution is corrected in Step c.

The inverse problem of equation (1) usually tends to overfit, so a regularization is required. In this work, a local smoothness of the emission profile is required. At the same time, a high gradient of emission profile near the divertor targets requires higher spatial resolution in this area. Figure 1 shows a map of 9549 axially symmetric light sources of unknown emissivity, which covers the entire SOL and divertor. The source linear sizes vary from 0.5 cm to 4 cm. The map is defined on a rectangular grid with 0.5 × 0.5 cm² grid cells. The larger light sources cover multiple (up to 64) grid cells. The source cross-sections are square except for those located at the map edges and at the borders of the areas with different sizes of light sources in the divertor. In the second step of the algorithm, the emission distribution is reconstructed only below the dashed line in figure 1, the rest of the emission in the main chamber is neglected.

The accuracy of this tomographic approach is tested within the framework of the synthetic diagnostic. The synthetic signals are calculated using the emissivity specified on a very fine regular mesh with 2.5 mm step in the R and Z directions, thus keeping the accuracy of the original SOLPS-4.3 + OEDGE mesh everywhere, including the divertor. The D_α emission profiles were obtained for the two-dimensional contour of the FW, and its adjustment to the three-dimensional FW model required the profiles be shrunk by a couple of centimetres to put them inside the three-dimensional vacuum vessel.

First, the inverse problem is solved in an ideal case, namely for the scenario 'o' in table 1 characterised by the highest ratio of $P_{{\mathrm{D\alpha }}}^{{\text{SOL}}}/P_{{\mathrm{D\alpha }}}^{{\text{Div}}}$, and assuming that the light reflective properties of the FW are known exactly. Next, a case is considered in which the synthetic signal is calculated for slightly different BRDFs of the FW components than those that were assumed when calculating the ray transfer matrices. Finally, less favourable cases with a lower value of $P_{{\mathrm{D\alpha }}}^{{\text{SOL}}}/P_{{\mathrm{D\alpha }}}^{{\text{Div}}}$ are considered.

Here, as in [14, 15], the Cook–Torrance micro-facet BRDF model [16] is used to model light reflections from the metal surfaces of the FW. The distribution of microfacets and consequently the deviations in the direction of the normal vector of the surface are characterized by a surface roughness. In the Cook–Torrance model, the surface roughness is a dimensionless parameter responsible for both microfacet distribution and self-shadowing of microfacets (see equation (8) and figure 3 in [14]). The data for the refractive index and extinction coefficient in the Fresnel term are taken from the refractive index database [17]. The dimensionless parameter of surface roughness is taken equal to 0.26 for beryllium (blanket modules in the main chamber), 0.29 for tungsten (divertor plates) and 0.13 for stainless steel (port plugs). The first two values of roughness are currently used for beryllium and tungsten components of JET-ILW in the Cherab module for JET tokamak (see table 1 in [14]). Until real experimental data for the actual ITER components are obtained, these are assumed to be reasonable choices given the similarity of the JET-ILW materials to those which will be used in ITER and the fact that the JET plasma-facing surfaces have been exposed to many hours of real tokamak discharges. Because the definitive final version of the ITER FW three-dimensional model is not yet available, the same FW model of a single 20° diagnostic sector as in [18] is used in this work. The complete FW model is obtained by copying and rotating this sector. The latest models of ITER diagnostic port plugs are used where available, and the general model of a port plug is used in all other sectors.

The inverse problem is first solved in the case when the same model BRDFs are used to compute both the ray transfer matrices and the synthetic signals. In the case of success, this proves that the method works in principle. Ray transfer matrices and synthetic signals are calculated for the following diagnostics:

• 5 LoS arrays of the DIM located in upper port (UP) 1 (1 array), equatorial port (EP) 1 (1 array) and divertor port (DP) 2 (3 arrays).
• 3 filtered cameras of the Vis/IR TV located in EP 12 (2 cameras) and EP 17 (1 camera).
• 3 filtered cameras of H-alpha and Visible Spectroscopy located in EP 11 (2 cameras) and EP 12 (1 camera).

Since both the synthetic and sought after D_α emission profiles are two-dimensional and the wall model is periodic in the toroidal direction (except for the sector-specific port plugs), the toroidal sector number in which the camera is located, is not so important. Each DIM LoS array contains 71 aligned LoS. All cameras are modelled as pinholes, neglecting the image distortions caused by the optics. Note that Raysect is capable of modelling complex optical systems, so when the final specifications for the Vis/IR TV cameras are available, the optics will be taken into account in the model. The image resolution for the filtered cameras is taken equal to 10 pixels per degree of observation angle, e.g. the camera with the FoV of 21.5° × 26.9° has the image resolution of 215 × 269 pixels. For the two filtered cameras of the Vis/IR TV system in EP 12, which observes some parts of divertor, the dynamic range limit is set to 20% of maximum light intensity for better resolution of relatively dark SOL regions. As a result, all pixels with intensity higher than this limit are considered saturated and their respective intensity values are excluded from the vector b in (1). Saturated pixels in the Vis/IR TV images do not prevent reconstruction of the emission in the areas affected by the saturation because the same areas are also observed by DIM detectors. All filtered cameras and LoS arrays are considered cross-calibrated and the calibration error is neglected. Adaptive sampling is used to calculate synthetic signals by adjusting the number of rays cast per pixel so that the standard error does not exceed 10% for 99% of the pixels.

Unlike [11] where the general optimisation method of simulated annealing was used to solve the problem in equation (1), here the simultaneous algebraic reconstruction technique (SART) method [19] is used, which was specially created to solve tomographic problems with a limited number of measurements. A high-performance implementation of the SART method for graphical processing units has been developed for this work and added to the Cherab framework.

The divertor D_α emission distribution is reconstructed using five DIM receivers and the segments of images from the Vis/IR TV filtered cameras in EP 12. The results of useful signal recovery for 4 out of 5 DIM receivers in simulation scenario 'o' (see table 1) are shown in figure 2. The LoS of the 5th DIM receiver pass between the divertor cassettes on both the low-field side (LFS) and the high-field side (HFS). As a result, the signal on this receiver does not contain the reflected light. Similar results, but now for the segments of the filtered camera images of Vis/IR TV cameras in EP 12 are given in figure 3. Evidently, the recovered useful signal fits the true one with a high accuracy for both the DIM LoS and the Vis/IR TV cameras. This confirms the result of [11] that with the help of the tomographic approach, the stray light can be filtered out for the receivers which observe the divertor.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Figure 4 shows the comparison of the reconstructed and true divertor D_α emission profiles in scenario 'o'. Here the emission in the rest of the SOL is neglected, which corresponds to the 2nd step of the algorithm described in section 2. Despite successful recovery of the useful signal, the reconstructed spatial distribution of emission is different from the true emission profile. Note, however, that the profiles are plotted with a logarithmic scale and thus that the emission in the regions of greatest discrepancy is two orders of magnitude lower than that at the strike points. Despite this discrepancy, both the location of strike points and the X-point position can be identified correctly from the reconstructed emission. The recovered value of the total D_α power emitted in divertor is 8% higher than the true value.

Zoom In Zoom Out Reset image size

When recovering the axisymmetric emission profile, the filtered cameras of H-alpha and Visible Spectroscopy do not provide any additional information compared to the Vis/IR TV system. They can be useful for reconstructing the three-dimensional emission profiles, or as a backup to the Vis/IR TV cameras. There is, in addition, no need to use the entire Vis/IR TV image: only one half of symmetric images is required. The results of useful signal recovery for segments of images from the three filtered Vis/IR TV cameras in EP 12 and 17 are given in figure 5. Even for the cameras observing the main chamber, where the divertor stray light is high, the recovered useful signal is a good fit to the true value.

Zoom In Zoom Out Reset image size

Figures 6(d) and (h) compare the reconstructed and true D_α emission profiles in scenario 'o'. Note the decrease of the recovered $P_{{\mathrm{D\alpha }}}^{{\text{Div}}}$ when the emission in main chamber is taken into account. In the LFS SOL, the D_α emission profile can be recovered only along the poloidal direction since, unlike the HFS SOL, this region is not observed at oblique angles. Considering that the scenario 'o' is characterised by a very strong D_α emission in main chamber and that the BRDF uncertainty is neglected, figure 6(d) shows the best possible reconstruction of the D_α emission profile that can be achieved using the capabilities of Vis/IR TV and DIM diagnostics without applying any model functions or complex (e.g. machine-learning based) regularization rules to the emission profile. In the other scenarios in table 1 with lower ratios of $P_{{\mathrm{D\alpha }}}^{{\text{SOL}}}/P_{{\mathrm{D\alpha }}}^{{\text{Div}}}$, or when the model BRDFs do not perfectly reproduce the real ones, or if the cameras are not well calibrated, the results will be worse.

Zoom In Zoom Out Reset image size

Table 2. Values of dimensionless roughness, r, albedo at normal incident angle, ρ₀, and albedo averaged over the incident angle, $\langle \rho \rangle$, used to calculate synthetic signals.

		Beryllium			Tungsten			Stainless steel
#	Difference with assumed roughness	r	ρ0	$\langle \rho \rangle$	r	ρ0	$\langle \rho \rangle$	r	ρ0	$\langle \rho \rangle$
1	50% lower	0.130	0.530	0.508	0.145	0.491	0.475	0.065	0.529	0.526
2	25% lower	0.195	0.514	0.489	0.218	0.472	0.455	0.098	0.526	0.528
3	15% lower	0.221	0.506	0.481	0.247	0.463	0.446	0.111	0.524	0.515
4	exact match	0.260	0.491	0.469	0.290	0.446	0.433	0.13	0.521	0.510
5	15% higher	0.299	0.475	0.456	0.334	0.428	0.419	0.15	0.517	0.504
6	25% higher	0.325	0.463	0.447	0.363	0.414	0.409	0.163	0.514	0.500
7	50% higher	0.390	0.431	0.425	0.435	0.379	0.386	0.195	0.506	0.491

The results of previous section are obtained under the assumption that the FW BRDFs are known exactly or, put another way, the same BRDFs were used to calculate both the ray transfer matrices and the synthetic signals. This section analyses the case when the assumed BRDFs, used to calculate the ray transfer matrices, are different from the true values used to calculate the synthetic signals. Because the Cook–Torrance BRDF is characterised by a single parameter (the surface roughness) it is important to vary this value in order to test the sensitivity of the reconstructions to this parameter. To estimate the error caused by using the ray transfer matrices calculated for the incorrect value of surface roughness, the synthetic signals are calculated for different sets of roughness, r, considered here as possible variants of the unknown true value of surface roughness, while the matrices are kept unchanged. Table 2 shows the values of r used in these sets, along with the resulting values of albedo at normal incident angle, ρ₀, and albedo averaged over the incident angle, $\langle \rho \rangle$, both at 656.1 nm (the wavelength of the D_α emission). Figures showing the comparison of Cook–Torrance BRDFs for beryllium with different values of roughness and the respective albedos at 656.1 nm as functions of incident angle, as well as those giving the synthetic signals on DIM LoS arrays and Vis/IR TV filtered cameras for different sets of roughness from table 2, are available online as a supplementary material.

The comparison of reconstructed emission profiles in cases when the true surface roughness is different from the assumed value is shown in figure 6. One can see that if the true surface roughness is higher than assumed, the reconstruction is more accurate than if it is lower. Similar figures but for the cases with lower $P_{{\mathrm{D\alpha }}}^{{\text{SOL}}}/P_{{\mathrm{D\alpha }}}^{{\text{Div}}}$ ratio can be found in the supplementary material.

Because the main goal of this tomographic approach is to recover the useful signal in the FoV of spectroscopic diagnostics, the error is analysed for that recovery. The H-alpha and Visible Spectroscopy diagnostic is the one most heavily affected by the stray light problem. Figure 7 compares the true and recovered useful signals as a function of vertical coordinate in the H-alpha Spectroscopy diagnostic FoV for simulation scenario 'o'. For the tangential FoV from EP 12, the results for the HFS and LFS SOL are given separately. The relative error of signal recovery is given below each plot.

Zoom In Zoom Out Reset image size

Evidently, the error varies greatly not only between different FoV but also within the same FoV. Relatively low errors are achieved near the top of the main chamber, and near the lower region of the central column. This can be explained as follows: radiation from divertor impacts the FW panels at the top of the main chamber at near normal incidence so that the resulting distribution of reflected radiation depends weakly on the surface roughness. Tangential views (such as the left view of the Vis/IR TV filtered camera in EP 12) provide relatively high signal-to-stray-light ratio in the HFS SOL near the lower inboard panels, so the D_α emission profile can be reconstructed in that area with a high accuracy, even if the assumed BRDF differs much from the true value.

The results presented thus far have been focused on simulation scenario 'o', with the largest value of the ratio of main chamber to divertor emission. However, it is important for the evaluation of diagnostic performance, to see how the error of the useful signal recovery depends on the calculation scenario, and especially on the ratio $P_{{\mathrm{D\alpha }}}^{{\text{SOL}}}/P_{{\mathrm{D\alpha }}}^{{\text{Div}}}$. Unfortunately, the set of simulation scenarios does not uniformly cover the whole range between the minimum and maximum values of the $P_{{\mathrm{D\alpha }}}^{{\text{SOL}}}/P_{{\mathrm{D\alpha }}}^{{\text{Div}}}$ ratio (e.g. cases are missing between $P_{{\mathrm{D\alpha }}}^{{\text{SOL}}}/P_{{\mathrm{D\alpha }}}^{{\text{Div}}}\,$= 0.0096 and 0.096, corresponding to the two groups of scenarios 'h'–'k' and 'l–o' in table 1). To fill the gap nine additional D_α emission profiles have been constructed by combining the emission profiles of the scenarios 'i' and 'o' linearly with variable weights. This yields additional cases with $P_{{\mathrm{D\alpha }}}^{{\text{SOL}}}/P_{{\mathrm{D\alpha }}}^{{\text{Div}}}$ ratio varying gradually from 0.018 to 0.088. The inverse problem of D_α emission reconstruction has been solved for all of the 21 cases which results when combining those in table 1 with the newly generated values and for all sets of surface roughness from table 2.

The average error of useful signal recovery weighted by the true useful signal is calculated for all FoV of the H-alpha spectroscopy diagnostic:

$$\begin{equation}{E_{\text{w}}} = \mathop \sum \limits_p \frac{{S_p^{\text{T}}}}{{\mathop \sum \nolimits_{\mathop p\limits^` } S_{\mathop p\limits^` }^{\text{T}}}}\frac{{\left| {S_p^{\text{R}} - S_p^{\text{T}}} \right|}}{{S_p^{\text{T}}}} = \frac{{\mathop \sum \nolimits_p \left| {S_p^{\text{R}} - S_p^{\text{T}}} \right|}}{{\mathop \sum \nolimits_{\mathop p\limits^` } S_{\mathop p\limits^` }^{\text{T}}}}.\end{equation} \tag{ 2 }$$

Here, $S_p^{\text{T}}$ is the true useful signal, $S_p^{\text{R}}$ is the recovered useful signal, $p$ and $\mathop p^`$ are the pixel indices. The sum is taken over all pixels of the image, as in the case of the EP 11 top and bottom FoV, or over the segments of the image, as for the EP 12 FoV. Because the relative error is quite different in the top and bottom parts of the EP 12 FoV, the values of E_w are calculated separately not only for the LFS and HFS regions, but also for the top and bottom parts of the image. Thus, six values of E_w are calculated for each simulation scenario for each set of surface roughness, which yields a total of 882 E_w values. All are shown in figure 8 in six separate plots.

Zoom In Zoom Out Reset image size

The ITER Measurement Requirements stipulate an accuracy of 30% for the recovery of the influx of neutral hydrogen isotopes. Obviously, the influx cannot be recovered more accurately than the useful signal from which it will then be derived. Excluding the most unfavourable case when the true surface roughness is 50% lower than the assumed value, the useful signal can be recovered with 30% accuracy in the bottom parts of the EP 11 and EP 12 FoV already at $P_{{\mathrm{D\alpha }}}^{{\text{SOL}}}/P_{{\mathrm{D\alpha }}}^{{\text{Div}}}$ = 0.06. For the bottom HFS part of the EP 12 FoV, the same accuracy can be achieved even at $P_{{\mathrm{D\alpha }}}^{{\text{SOL}}}/P_{{\mathrm{D\alpha }}}^{{\text{Div}}}$ = 0.025. Unfortunately, recovering the useful signal with 30% accuracy in the simulation scenarios 'd', 'e', 'f' and 'g' is impossible anywhere unless the BRDFs are known exactly, which is unrealistic. It may also be noted that 30% accuracy for the useful signal is unlikely to be sufficient for recovery of the influx with the same accuracy since obtaining the influx of neutral hydrogen isotopes from the measured Balmer-alpha signal is a complex problem, requiring additional measurements either of molecular spectra (so called S/XB + D/XB method [20]), or high-resolution Balmer-alpha spectra [21–23]. Moreover, it is not yet clear if the highest achievable accuracy of 7%–10% possible only if the FW BRDF is known exactly, will be sufficient.

Even with the roughness parameter equal to 1, the Cook–Torrance BRDF cannot describe an ideal diffuse reflector, such as a Lambertian surface. In previous analyses of light reflections in the ITER main chamber [2, 11], the BRDF was composed of diffusive and Gaussian-broadened specular components, and the stray light level was obtained for various ratios of diffusive to specular reflections. Similar analysis can be performed here with a generalised two-component model BRDF, which combines the Cook–Torrance and the Lambertian parts:

$$\begin{equation}f = \left( {1 - {w_{\text{L}}}} \right){f_{{\text{C}} - {\text{T}}}}\left( {\lambda ,r,{n_{\text{i}}},{n_{\text{o}}}} \right) + {w_{\text{L}}}\rho \left( {r,\lambda } \right)/\pi .\end{equation} \tag{ 3 }$$

Here, f_C–T is the Cook–Torrance BRDF, λ is the wavelength, n _i and n _o are the directions of incident and outgoing radiation respectively and w_L is the weight of Lambertian component. The surface albedo averaged over the incident angle, $\rho \left( {r,\lambda } \right)$, is defined by equation (4):

$$\begin{equation}\rho \left( {r,\lambda } \right) = \frac{1}{{2\pi }}\mathop \int \limits_{{{{\Omega }}_{\text{i}}}}^{} \mathop \int \limits_{{{{\Omega }}_{\text{o}}}}^{} {f_{{\text{C}} - {\text{T}}}}\left( {\lambda ,r,{n_{\text{i}}},{n_{\text{o}}}} \right)cos\left( {{\theta _{\text{o}}}} \right)d{{{\Omega }}_{\text{o}}}d{{{\Omega }}_{\text{i}}},\end{equation} \tag{ 4 }$$

where the integrals are taken over the solid angles Ω_o and Ω_i in their respective hemispheres and θ_o is the angle between the outgoing direction and the surface normal.

Assume that the Lambertian weight, w_L, in equation (3) is not known exactly. To obtain the impact of the error in determining the Lambertian weight on the accuracy of recovery of the useful signal, r must be fixed and only w_L varied. Unlike the case with variable r, changing the value of w_L in equation (3) does not affect $\rho \left( {r,\lambda } \right)$. Synthetic signals are calculated for seven different values of w_L: 0.25, 0.375, 0.425, 0.5, 0.575, 0.625, 0.75, while the ray transfer matrices are calculated for w_L = 0.5. As before (based on values from the JET-ILW), the surface roughnesses are fixed in all calculations at 0.26 for beryllium, 0.29 for tungsten and 0.13 for stainless steel.

Figures showing the synthetic signals on the DIM LoS arrays and Vis/IR TV filtered cameras for varying w_L for the simulation scenario 'o' are available online as a supplementary material. Likewise for respective comparisons of reconstructed D_α emission profiles similar to figure 6 and comparison of recovered useful signals in the H-alpha Spectroscopy diagnostic FoV similar to figure 7. Both the DIM signals and the Vis/IR TV images demonstrate much weaker dependence on w_L than on r. The reconstructed emission profiles also look much closer to one another. Even when the true Lambertian weight is a factor 2 lower than the assumed value, the emission profile is reconstructed surprisingly well. Even better agreement is found in the recovered useful signals. The relative errors exceed 50% only in the areas where the useful signal is low. Finally, the average weighted errors (equation (2)) of the useful signal recovery in the H-alpha spectroscopy FoV shown in figure 9 are much smaller than those shown in figure 8. Even if the assumed Lambertian weight is 50% lower or higher than the true value, the accuracy of 30% can be achieved when recovering the signals in the HFS SOL. However, these results are obtained for moderate surface roughness. In the case of lower r( and sharper BRDF), the dependence of the synthetic signals on the Lambertian weight may be stronger.

Zoom In Zoom Out Reset image size

The use of the Raysect and Cherab numerical codes in the ITER H-alpha synthetic diagnostic has significantly expanded its capabilities in terms of both the speed and the accuracy of simulations. It would take a great deal longer (and may not even be feasible) to complete the analysis presented here using the general-purpose illumination design software that has been previously used for stray light analysis in ITER. For example, calculating a ray transfer matrix for a set of light from figure 1 and a filtered camera with image resolution of 200 × 300 pixels will take about a month with LightTools on a single PC, because it requires turning on the light sources one at a time. The Raysect and Cherab allow to calculate such ray transfer matrix for all light sources simultaneously, and this calculation takes only a few hours on a PC.

The tomographic approach for filtering out the stray light in the measured signal is promising but requires some techniques (not yet developed) for measuring or reconstructing the BRDFs of the FW during the ITER operation. The light reflection properties are expected to change under plasma exposure, so it will be necessary to recalibrate the ray transfer matrices from time to time. In this work, the same BRDF has been used for all plasma-facing components made of the same material. It is, however, possible that after some time of operation, those plasma wetted areas subject to net erosion will become more reflective compared to shadowed areas (in which net deposition may also occur). As a result, the parameters of the model BRDF may have to adopt different values for different segments of the FW (e.g. not only on different FW panels, but also on the convex and recessed segments of the same panel). While it is not a problem to apply an inhomogeneous BRDF over the plasma-facing surfaces, its reconstruction may be a difficult task. There is the additional question (also discussed in [14]) as to whether the Cook–Torrance BRDF (or a generalised BRDF (3)) originally introduced for computer graphics is an adequate choice for modelling the reflection of the light from plasma-exposed materials. Measurements of BRDF [24] on tungsten-coated graphite samples from the WEST tokamak, and for the beryllium samples similar to those used in the JET-ILW, show little in common with the Cook–Torrance model.

Another issue which may restrict the applicability of the method described here is the necessity of cross-calibration of the signals measured with different receivers by different diagnostics. The Vis/IR TV filtered cameras must be cross-calibrated not only between themselves, but also with DIM LoS arrays and with the filtered cameras of the H-alpha and Visible Spectroscopy system. This problem has not been analysed in this work since all receivers were considered ideally calibrated. While the cross-calibration of different receivers is a feasible task it would require a joint effort of the three teams responsible for three different diagnostics.

Here, the proposed method is shown to work in principle. For discharges with a ratio $P_{{\mathrm{D\alpha }}}^{{\text{SOL}}}/P_{{\mathrm{D\alpha }}}^{{\text{Div}}}$ of D_α power emitted in the main chamber SOL to that emitted in the divertor, $>$ 0.06, the approach allows useful signal in the HFS SOL and in the lower region of the LFS SOL to be recovered in the FoV of the H-alpha and Visible Spectroscopy diagnostic with 30% accuracy, even if the true BRDFs are slightly different from the assumed values. Namely, if the true surface roughness, described by the dimensionless roughness parameter in the Cook–Torrance BRDF model is up to 15% lower or up to 25% higher than the assumed value. In the lower HFS region of the EP 12 FoV, the same accuracy can be achieved for plasmas in which $P_{{\mathrm{D\alpha }}}^{{\text{SOL}}}/P_{{\mathrm{D\alpha }}}^{{\text{Div}}}$ > 0.025.

Note that due to the complexity of the respective inverse problem, the 30% accuracy for the useful signal recovery is insufficient for recovery of the neutral hydrogen influx with the 30% accuracy stipulated by ITER Measurement Requirements. The tomographic method developed here should be incorporated into an extended synthetic H-alpha diagnostic which returns the achievable accuracy of the neutral hydrogenic influx measurement in the ITER main chamber using all available data. The latter should include the high-resolution spectroscopy, the wavelength-integrated tomography of Balmer-alpha spectral lines, the ongoing modification of the Ballistic model [21, 22] and the modified SXB method [23].

In [25] the divertor stray light in the signal on the LoS viewing an inner beryllium limiter in the JET-ILW was modelled using only high-resolution spectroscopy data. Taking into account the recent results in [15], the tomographic approach similar to that proposed for ITER should also be used on JET to filter out the stray light in the signals measured by spectroscopic systems and filtered cameras observing the main chamber.

It is important to note that besides the recovery of the useful signal, reconstruction of the Balmer-alpha emission profile itself is key for benchmarking plasma boundary codes such as the SOLPS.

As a final remark, it is important to recall that the ITER main chamber spectroscopy system is not just a monitor of D_α emission. It also has the key, perhaps even more important function, of providing measurements of the main chamber emission of Be atoms and ions. Such measurements provide an estimate of the FW Be gross erosion and thus, through comparison with modelling, of the Be migration pathways. This is essential regulatory information for estimates of tritium fuel inventory build-up and of a potentially important dust source which may be produced from the growth of Be co-deposits. A similar tomographic approach to that described here, aimed at reconstruction of three-dimensional emission profiles of Be atoms and ions in the ITER SOL is under development. The accuracy of this reconstruction will also be tested in the synthetic diagnostic framework with the results of very recent modelling of global Be migration in ITER with the ERO2.0 code [18] used as synthetic experimental data. Once developed, a tomographic reconstruction of three-dimensional emission profiles may be applied to any plasma species if needed, including the hydrogen isotopes.

The authors are grateful to H-alpha and Visible Spectroscopy team, Divertor Impurity Monitor team, Vis/IR TV team, all developers and contributors of Raysect and Cherab codes and personally to S Kajita, H Natsume, M Kocan, S Lisgo and A Loarte for helpful discussions. R I K thanks the ITER Organization for the hospitality and support during his internship in the second half of 2019, when some parts of this work were completed.

Some calculations for this work were carried out on the ITER HPC cluster.

The views and opinions expressed herein do not necessarily reflect those of the ITER Organization.