Retrieving fields from proton radiography without source profiles

TL;DR

This paper introduces a statistical method for retrieving integrated magnetic and electric fields from proton radiography data without requiring the unmodulated source profile, enhancing robustness and applicability in experimental settings.

Contribution

The authors develop a novel statistical approach that eliminates the need for source profile knowledge in proton radiography, improving field retrieval accuracy.

Findings

Method successfully applied to experimental data

Retrieves both field properties and statistical characteristics

Demonstrates robustness against shot-to-shot variability

Abstract

Proton radiography is a technique in high energy density science to diagnose magnetic and/or electric fields in a plasma by firing a proton beam and detecting its modulated intensity profile on a screen. Current approaches to retrieve the integrated field from the modulated intensity profile require the unmodulated beam intensity profile before the interaction, which is rarely available experimentally due to shot-to-shot variability. In this paper, we present a statistical method to retrieve the integrated field without needing to know the exact source profile. We apply our method to experimental data, showing the robustness of our approach. Our proposed technique allows not only for the retrieval of the path-integrated fields, but also of the statistical properties of the fields.

Equations

$$\mathbf{\alpha }(x_0,y_0)=-\nabla \Phi (x_0,y_0),$$

$$\Phi (x_0,y_0)=-\frac{e}{\sqrt{2mW}}\int \mathbf{A}(x_0,y_0,z_0)\cdot \mathrm{d}{\mathbf{z}}_0.$$

$$\int \mathbf{B}\times \mathrm{d}{\mathbf{z}}_0=-\frac{\sqrt{2mW}}{e}\nabla \Phi .$$

$$x(x_0,y_0)$$

$$y(x_0,y_0)$$

$$I(x,y)=I_0(x,y){\left|\frac{\partial (x,y)}{\partial (x_0,y_0)}\right|}^{-1}.$$

$$\mu =L\alpha /a,$$

$${\nabla }^2\Phi \approx \frac{2W}{e}\left(\frac{I}{I_0}-1\right).$$

$$\mathbb{P}(\Phi |I,I_0)=\delta (\Phi -{\Phi }^{\prime }(I,I_0)),$$

$$\mathbb{P}(\Phi |I)=\int \mathbb{P}(\Phi |I,I_0)\mathbb{P}(I_0)\mathrm{d}{I}_0,$$

$$\mathbb{P}(I_0|\theta )=\mathcal{G}\mathcal{P}(I_0|⟨I_0⟩,\kappa (\cdot ,\cdot |\theta )),$$

$$\mathcal{G}\mathcal{P}(\mathbf{y}|\mu ,\kappa )=\frac{\exp \left[\frac{1}{2}(\mathbf{y}-\mu {)}^T{\mathbf{K}}^{-1}(\mathbf{y}-\mu )\right]}{\sqrt{(2\pi {)}^n\det (\mathbf{K})}},$$

$${\kappa }_{\mathrm{SE}}({\mathbf{x}}_1,{\mathbf{x}}_2|\theta )$$

$${\theta }_{ML}=\arg \underset{{\theta }^{\prime }}{\max }\sum _j\log \mathbb{P}({\hat{I}}_{0j}|{\theta }^{\prime }),$$

$$\mathbb{P}(\Phi |I)=\int \mathbb{P}(\Phi |I,I_0)\mathbb{P}(I_0|\theta )\mathbb{P}(\theta )\mathrm{d}\theta \mathrm{d}{I}_0.$$

Full text

I Introduction

Proton radiography is a popular method to diagnose magnetic and/or electric fields in laser-produced plasmas and inertial confinement fusion experiments Nilson et al. (2006); Tzeferacos et al. (2018); Willingale et al. (2011); Romagnani et al. (2010); Li et al. (2007); Huntington et al. (2015); Borghesi et al. (2002). In proton radiography, a low density proton beam is fired into a plasma, gets deflected by its electric and magnetic fields, and the beam intensity profile is measured on a screen placed at some distance. The modulated intensity of the beam is related to the deflection of the proton beam and thus related to the integrated magnetic and electric field. Therefore, by comparing the unmodulated beam intensity and the modulated intensity one can retrieve the integrated magnetic and electric fields.

Earlier attempts to reconstruct fields from proton radiography utilized Poisson’s equation solver Kugland et al. (2012) and a diffusion model Graziani et al. (2017). However, those approaches use linear approximations of the proton radiography forward model, while most of the interesting interactions happen in the non-linear regime. Later, Graziani et al. Graziani et al. (2017) and Kasim et al. Kasim et al. (2017) independently realized that integrated field retrieval with proton radiography is a subclass of the optimal transport problem first posed in 1781 Monge (1781). Establishing the connection between proton radiography and optimal transport opens up numerous algorithms for proton radiography field reconstruction even in the non-linear regime. This has been shown by Kasim et al. Kasim et al. (2017) and Bott et al. Bott et al. (2017) by applying off-the-shelf optimal transport algorithms Aurenhammer et al. (1998); Sulman et al. (2011) to proton radiography, achieving remarkable results in some parameter regions.

Despite their successful implementations, current algorithms to solve the proton radiography inverse problem need the exact knowledge of the unmodulated beam intensity profile (source profile) as well as the modulated beam intensity profile. Unfortunately, there is no viable non-destructive method to simultaneously capture the proton beam intensity profiles before and after interacting with the plasma. Moreover, the source profiles are different from shot-to-shot by some considerable amount Manuel et al. (2012) making it harder to determine the source profile when the measurement is performed. Only the statistical properties of the source profile variability can be determined reliably, by measuring the beam intensity multiple times without any field interaction Manuel et al. (2012); Rygg et al. (2015).

In this paper, we introduce a statistical method in retrieving the integrated field only by using the statistics of the source profiles, without the exact knowledge of the source profile. Besides retrieving the integrated field, the method also provides the probability distribution of the retrieved field allowing one to perform a statistical analysis of the results, such as checking if the measurement suffers from the inverse problem instability Kasim et al. (2018).

II Theory

A schematic of a typical proton radiography set-up can be seen in Figure 1. In a typical case, a proton beam is fired from a source through an object with magnetic and/or electric fields. The fields deflect the trajectory of the beam, forming an intensity modulation on the screen behind the object.

Only cases with a single type of dominant field, either magnetic or electric field, is considered in this paper. It is also assumed that the deflection is small enough so that most of the deflected protons reach the screen. It is noted that the deflection can still be large enough to reach the non-linear regime.

If a proton beam is fired through a plasma with magnetic field, the beam at coordinate $(x_{0},y_{0})$ on the object plane will be deflected to an angle

[TABLE]

for $\alpha\ll 1$, where $\Phi$ is the integrated field given by Kugland et al. (2012)

[TABLE]

The integrated field, $\Phi$, depends on the charge $e$ and mass $m$ of each proton, the kinetic energy of the beam $W$ and the magnetic vector potential in the object $\mathbf{A}$. The equation for cases where the electric field is dominant was presented in Kugland et al. (2012).

The integrated field, $\Phi$, is an important parameter in the retrieval because by knowing $\Phi$ it is straight-forward to retrieve the integrated magnetic field,

[TABLE]

In the cases where the effect from electric field is dominant, the integrated transversal electric field can be retrieved by following the similar steps as above. For both cases, it is the integrated field profile, $\Phi$, that needs to be retrieved before obtaining the integrated magnetic or electric fields.

If the distance between the object and the screen is $L$, the beam at $(x_{0},y_{0})$ on the object plane is mapped to the coordinate $(x,y)$ on the screen where they are given by

[TABLE]

assuming the plasma size is much smaller than $L$ and the beam source is collimated. If the beam is diverging from a point source at distance $l$ from the object, then $x_{0}$ and $y_{0}$ are replaced by $x_{0}\rightarrow x_{0}(1+L/l)$ and $y_{0}\rightarrow y_{0}(1+L/l)$ in the equations above.

An analytical equation to determine the intensity on the screen is given by

[TABLE]

The term $I_{0}(x,y)$ denotes the intensity profile without any deflections from the object and the term $|\partial(x,y)/\partial(x_{0},y_{0})|$ is the absolute determinant of the Jacobian matrix of $(x,y)$ with respect to $(x_{0},y_{0})$. If the deflection is large enough, the determinant of the Jacobian can be close to zero at some points and the intensity at the corresponding positions can reach a very high value. These are known as caustics.

One way to describe the deflection strength in proton radiography is using a dimensionless variable, as introduced by Kugland et al. Kugland et al. (2012), which is defined as

[TABLE]

with $L$ the distance from the system to the screen, $\alpha$ the magnitude of the deflection angle, and $a$ the size of the field perturbation.

If $\mu\ll 1$, equation 6 can be linearized and the integrated field can be obtained by solving Poisson’s equation Kugland et al. (2012),

[TABLE]

The two integrals in Poisson’s equation solver amplify the lower frequency components of $I/I_{0}$ by a factor proportional to $1/k^{2}$ with the wavenumber $k$.

When $\mu$ exceeds a certain value, $\mu\geq\mu_{c}$, caustics form and parts of the beam cross each other. In this region, the relation between $\Phi$ and $I$ is no longer injective Kugland et al. (2012); Bott et al. (2017). This means that there are multiple profiles of $\Phi$ that correspond to the same intensity profile $I$.

If $\mu$ is on the order of unity and less than $\mu_{c}$, it gets into the so-called non-linear injective regime. In this regime, Poisson’s equation is no longer accurate, but the relation between $\Phi$ and $I$ is still injective. The integrated field, $\Phi$, can be obtained by solving the optimal transport problem: find a way in “transporting” the protons from the source profile to resemble the modulated profile with total squared displacement as small as possible. Numerous algorithms are available to solve this type of problem, for example Aurenhammer et al. (1998); Sulman et al. (2011).

For the rest of the paper, it is assumed that the proton deflection is always in the injective domain where there is no beam crossing. Field reconstruction in the non-injective domain is beyond the scope of this paper.

III Method

Optimal transport based algorithms to retrieve the integrated field from a source profile and a modulated intensity profile are deterministic. Therefore, the probability of getting $\Phi$ from $I$ and $I_{0}$ can be written as:

[TABLE]

where $\delta(\cdot)$ is the Dirac delta and $\Phi^{\prime}(I,I_{0})$ is the integrated field profile retrieved using the retrieval algorithms given the modulated intensity profile $I$, and the source profile $I_{0}$.

When we have no exact knowledge of the source profile $I_{0}$, one thing that we can do is to marginalize it by integrating the term $I_{0}$ over its probability distribution,

[TABLE]

with $\mathbb{P}(\Phi|I,I_{0})$ given by equation 9. From the equation above, we can get the probability distribution of the integrated field $\Phi$ by generating multiple source profile samples according to $\mathbb{P}(I_{0})$ and apply the retrieval algorithm for every generated profile. In this paper we use the algorithm from Sulman et al. (2011) to retrieve the integrated field from a modulated intensity profile and a source profile.

Now the problem is shifted to determining the probability distribution of the source profile, $\mathbb{P}(I_{0})$. Based on Manuel et al. (2012), the proton profile from a DD and D3He fusion source deviates by several percent from the mean profile and the deviations of two nearby points are correlated. The deviation is larger for proton sources generated by target normal sheath acceleration (TNSA) Warwick et al. (2017). Correlated deviations cannot be captured by a standard normal distribution because it assumes uncorrelated deviations. One way to capture correlated deviations is by representing them using a Gaussian Process (GP) Williams and Rasmussen (1996); Rasmussen (2003):

[TABLE]

where $\mathbb{P}(I_{0}|\theta)$ is the parameterized prior of the source profile with hyperparameters $\theta$, $\langle I_{0}\rangle$ is the mean intensity for the source profile, and $\kappa(\cdot,\cdot|\theta)$ is a kernel function that represents the correlation of two points with hyperparameters $\theta$. The source profile’s mean intensity should be equal to the modulated profile’s mean intensity.

The explicit expression for GP for an $n$-elements multivariate random variable, $\mathbf{y}\in\mathbb{R}^{n}$, is given by

[TABLE]

where each element in the covariance matrix $\mathbf{K}$ is the value of the kernel $\kappa(\cdot,\cdot)$. In this case, the variable $\mathbf{y}$ is the vectorized source profile $I_{0}$, and $n$ is the number of pixels. The inputs to the kernel are the position of each pixel in $I_{0}$.

One of the most commonly used kernels is the squared exponential kernel,

[TABLE]

where $\theta=\{\sigma,d\}$ are the hyperparameters of the kernel, $\sigma$ is the expected deviation, and $d$ is the correlated distance. Increasing $\sigma$ will increase the standard deviation from the mean of the profiles. Increasing $d$ will make the deviation less oscillatory. To illustrate the effect of $\sigma$ and $d$, figure 2 shows some samples taken from the Gaussian Process distribution with various values of $\sigma$ and $d$.

If source profiles are available experimentally, then one can fit the hyperparameters $\theta$ to the experimental source profiles by finding the maximum log-likelihood,

[TABLE]

where $\hat{I}_{0j}$ is the $j$-th source profile obtained experimentally. However, if there is no experimental data available, one should make a reasonable assumption over the prior probability of the hyperparameters, $\mathbb{P}(\theta)$. This changes Eq. (10) to

[TABLE]

One way to choose the prior distribution of the hyperparameters when there is no experimental data available is to set it based on known references, for example Manuel et al. (2012); Rygg et al. (2015). Another way is to choose a weak prior of the hyperparameters such as Jeffreys prior Paulo et al. (2005), a uniform, or a log uniform prior for certain range.

We present a summary of our method in algorithm 1.

IV Results

IV.1 Simulations

The first test on the proposed method is to retrieve integrated magnetic field statistics from modulated intensity profiles generated from various source profiles and a single integrated magnetic field profile. The source profiles in the test cases have variation about 10% from the mean value to follow the findings in Rygg et al. (2015) and have various correlation lengths. The proton radiography set up follows the experimental set up in Tzeferacos et al. (2018). The distance between the source and the object is 1 cm, object-to-screen distance is 27 cm, and the modulated intensity profile of 15 MeV proton is captured on a screen with size $(10\times 10)\ \mathrm{cm}^{2}$. The magnetic field profile is chosen to give non-linear deflections on the proton beam.

In retrieving the magnetic field in the simulations, we only use the modulated intensity profiles without using the source profiles. Because the actual source profiles are not known in the retrieval, we need to make a prior assumption about the correlation length. One good prior assumption of the correlation length is log-uniform distribution from 10% to 100% of the screen size. The range was chosen to capture the variation of source profile that could give considerably different integrated field profiles. It is to obtain an upper estimate of the variance and to avoid over-confidence. The log-uniform distribution is to make it more focused on the shorter correlation length where the largest variation of the retrieved integrated field profiles can be found.

The retrieval results for various source profiles are shown in Figure 3. From the bottom row in the figure, we can see that although in some cases the average retrieved profiles are quite different from the actual profile, the difference can be captured well inside the 95% confidence interval, except near the boundaries. This shows that the method can capture the uncertainty of the retrieved profiles relatively well. Worse results on the boundaries are because the inversion algorithm Sulman et al. (2011) forces the boundaries to have constant values which reduces the retrieval accuracy near the boundaries.

The mean retrieved profiles are very close to the actual integrated magnetic field profiles when the correlation length is either very short or very long. From the figure, we can say that the average retrieved profiles is relatively far from the actual profile when the correlation length is about 10% of the screen size (i.e. Figure 3(c)). When the correlation length is increased from 10%, the average retrieved profiles are getting closer to the actual profile and getting really close when it is about 80% of the screen size or more. This is because when the correlation length is very short, the source profile is similar to white noise which can be handled very well by the inversion algorithm. When the correlation length is very large, it becomes similar to a uniform source profile.

IV.2 Real experimental data

We apply our method in retrieving the integrated field as well as the uncertainty to an experiment of dynamo amplification of magnetic fields Tzeferacos et al. (2018). In the experiment, the authors used proton radiography to retrieve the integrated magnetic field using a proton beam generated from a D3He fusion capsule. One of the measured beam intensity profiles is shown in figure 4.

Because there is no information on the source profile statistics from this experiment, we take $\sigma$ to be $\sim$10%, close to the $\sim$13% standard deviation for proton beams generated from a D3He fusion capsule given in Rygg et al. (2015). For the correlated length $d$, we use a prior distribution which is log-uniform from 0.7 cm to 7 cm as in the previous subsection.

In each iteration, we choose a value of $d$ from the prior distribution above, generate a source profile from the Gaussian Process distribution using the chosen hyperparameters, retrieve the integrated field profile using the generated source, and put the integrated field profile into the pool of samples. Figure 5(a,b) shows the mean and standard deviation of the integrated field from 15,000 samples. From the figure, we can see that with $\sigma=10\%$, the integrated magnetic field has standard deviation of about $\sim$20%. The error magnification factor of only 2 is similar to the error propagated by a quadratic equation.

Besides calculating the mean and standard deviation of the integrated field, we also calculated the power spectra from the integrated field samples. This is shown in figure 5(c), along with the Kolmogorov’s power law ($S(k)\propto k^{-5/3}$) for reference. From the figure, we can see that the power spectrum uncertainty is high at small and large wavenumbers. High uncertainty at small wavenumbers is due to the nature of the retrieval algorithm that amplifies low wavenumber elements (i.e. amplification $\sim 1/k^{2}$). High uncertainty at large wavenumbers is due to limited precision in retrieving features with small amplitudes which happen to be at large wavenumber.

Figure 5(d) shows the power of the wavenumber in relation to the spectrum. Here we can see that the sample distribution peaks at $\sim k^{-1.5}$ which is shallower than Kolmogorov’s power law. We note that the “true” spectrum of the turbulent magnetic field in the experiment is expected to be even shallower than $\sim k^{-1.5}$ due to lack of injectivity at small scales that manifests diffusively to steepen the measured power spectrum, as shown in Tzeferacos et al. (2018). One needs to include more uncertainty factors to reach a meaningful conclusion, such as non-injectivity of crossing beams and uncertainty in proton energy, which are not the focus of this paper.

V Conclusions

We have presented a statistical method to retrieve integrated fields from proton radiography without knowing the exact source profile. The probability distribution of the integrated field is obtained by marginalizing out the probability distribution of the source profile. The distribution of source profiles can be obtained by collecting a number of proton radiography samples without any interaction with electric and magnetic fields, or by making a weak prior assumption. The method has been applied to an experiment to retrieve the integrated magnetic field and its statistics, showing the robustness of the proposed approach.

Acknowledgement

S.M.V. is grateful for support from the Royal Society. M.F.K. and S.M.V. acknowledge support from the UK EPSRC grant EP/P015794/1. G.G. and A.F.A.B. acknowledges support from AWE plc., and the UK EPSRC (EP/M022331/1 and EP/N014472/1). This work was supported in part from the U.S. DOE under Cooperative Agreement DE-NA0001944 to the University of Rochester, Field Work Proposal 57789 to Argonne National Laboratory, and grants DE-NA0002724, DE-NA0003605, and DE-SC0016566 to the University of Chicago; and from the National Science Foundation under grant PHY-1619573. We acknowledge support from DOE NNSA under subcontracts B632670 and 536203 with the University of Chicago.

Code availability

The code associated with this paper can be found on: https://github.com/OxfordHED/proton-radiography-no-source/