

\chapter{X-Ray Dispersive Spectroscopy Methods in HED Science}








%\section{\label{Introduction} Introduction}





The birth of spectroscopy as a scientific discipline can be traced back to Joseph von Fraunhofer's work, where he used diffraction gratings to systematically divide the spectral component of the electromagnetic radiation, dispersing them into distinct regions of space~\cite{born2013principles}. Following this, the invention of the photographic plate in the mid-19th century enabled the first permanent recording of spectral lines, representing the earliest example of a light detector. The following discovery of the photoelectric effect, along with the development of quantum mechanics to support the underlying theory, led to the development of modern detectors, such as photomultiplier tubes (PMTs)~\cite{foord1969use} and later, semiconductor-based detectors like charge-coupled devices (CCDs) and finally hybrid detectors~\cite{alle2016comparison}, which are the most commonly used nowadays and are the detectors considered in this thesis. These detectors consist of individual units, called pixels, arranged in ordered arrays. Each pixel is composed of two distinct components (hence the term "hybrid"): a semiconductor layer that captures photon energy and a readout electronics layer that amplify and quantify the signal produced by the photon capture.

As mentioned in the previous chapter, x-ray spectroscopy is an essential tool to investigate numerous properties (e.g. chemical composition, electronic structure) of HED systems, as extreme thermodynamic conditions can only be maintained experimentally for very short times, necessitating the use of \textit{in situ} photon-based diagnostics. Furthermore, the deep penetration of x-rays allows for the exploration of bulk material properties, while their short wavelengths enable the excitation of inner-shell electrons and detailed probing of material structure, expanding the range of spectroscopic techniques available. Lastly, the typical short pulse durations of x-ray sources (see Introduction) allow for time-resolved studies of dynamic processes, usually driven by other coupled instruments.

X-ray spectroscopy is predominantly performed using dispersive elements such as crystals or gratings, coupled with a detector in various geometries. These elements disperse (or selectively transmit in scanning mode~\cite{muller2015quick}) the spectral components of the radiation emitted by the target, which are then collected by the detector. In the following sections, I will analyse separately two of the most common crystal-based configurations in HED science: the flat crystal and the von H\'amos geometries. Additionally, the focus will be on how to treat camera images characterised by low photon counts, except for small regions around transition lines -- a common occurrence in HED physics. These low photon numbers typically result from the low cross-sections of some spectroscopic processes used in HED physics (e.g., RIXS, which, as discussed in Chapter 3, is a second-order scattering process), short exposure times due to the brief lifespan of HED systems, and limited solid angle coverage driven by the pursuit of high resolution and geometrical constraints within the experimental chambers. If the photon density is low everywhere on the camera and resolution requirements are not stringent, single-photon counting spectroscopy (SPCS)~\cite{basden2003photon} - a method suited to pixel detectors - becomes preferable to dispersive spectroscopy. In SPCS, photon energies are determined based on the signal generated by each photon hit, rather than their position on the detector. Although SPCS is widely used with modern detectors, it has not been employed in this thesis and therefore will not be discussed further.

The objective of this chapter is first to provide a step-by-step guide on how to construct a spectrum from the initial detector images, a resource which, to the best of our knowledge, is not found in existing literature as a single comprehensive document, but rather scattered across various sources. In doing this, I also present methods, which are validated on synthetic and HED experimental data, for enhancing the signal-to-noise ratio (SNR) and the resolution of a camera image, and capable of treating simultaneously regions with high and low density of photons. Although these cases are motivated by the typical features occurring in
HED experiments, the techniques discussed are applicable to other science fields with analogous image characteristics and geometrical setups.

%Despite these features being common and completely general, we will showcase these techniques on experimental  data coming from HED science, after having tested them on synthetic images.





























\section{\label{Constructing a spectrum} Constructing a Spectrum}






\begin{figure*}
    \centering
    \includegraphics[width=\textwidth,keepaspectratio]{Figures/Introduction_figure.pdf}
    \caption[Example of spectrum construction from a camera image]{Example of the construction of a spectrum (right) from the relative camera image $A_{ij}$ (left), which is a 2048x2048 pixel array. The image presents area with different photon density, labelled in the associated spectrum. In this example, the signal is the emission spectrum of warm-dense MgF$_2$ from Ref.~\cite{Gawne_PRE_2023} The bright arcs in the image correspond to the $K\alpha$ (main) and $K\beta$ (secondary) emission lines, coming from different ionisation states of the Mg ions. The different strength of the lines comes from the different ion populations.
    }
    \label{fig:general process}
\end{figure*}

In this section, the techniques used to construct the spectrum from the initial camera image are outlined. For completeness, the section begins with a brief overview of how a semiconductor pixel detector works, which will be beneficial for subsequent discussions. Next, the process of finding the energy map for two of the most common types of spectrometers is detailed. Finally, the process of converting the Analog-to-Digital Units (ADU) map to the distribution of photon hits by means of single photon counting algorithms and image deconvolution techniques is explained. In the following, the underlying physical processes generating the photons detected by the spectrometer is disregarded, focusing only on the problem of constructing a spectrum from the initial matrix $A_{ij}$, whose elements represent the total signal on pixel $(ij)$. Here, the pixels form a two-dimensional array labelled by $i=1,\dots, N_y$ and $j=1,\dots, N_x$, with $N_y$ and $N_x$ representing the number of pixels in the vertical and horizontal directions, respectively. As an example, Fig.~\ref{fig:general process} illustrates the construction, from the corresponding camera image, of the emission spectrum of warm-dense MgF$_2$. This data, taken from Ref.~\cite{Gawne_PRE_2023}, was obtained by irradiating MgF$_2$ with intense x-rays approximately 300 eV above the atomic K-edge of Mg ($\sim$ 1303 eV). As described in~\cite{vinko2012creation}, this radiation generates a population of Mg and F ions in various ionisation states, with their relative abundances determined by the rates of atomic processes. The energy of photons emitted during electron photodecays depends on the ionisation state of the ion involved, yielding a series of spectral lines. The resulting spectrum is divided into different regions based on the signal intensity, which corresponds to the local photon density on the camera. This data is analysed in Section~\ref{Real data} using the image processing methods introduced in the following sections.
 
The construction of the spectrum can be roughly divided into 3 independent and sequential steps:
\begin{enumerate}
    \item \textbf{Computing the energy map of the camera.} This process involves determining the energy associated with each pixel, $E_{ij}$. These energies represent averaged values for each pixel, accounting for effects such as source broadening and finite pixel size. These effects mix the energies of photons hitting a single pixel, ultimately reducing the maximum achievable spectral resolution. Section \ref{Energy calibration} details the techniques for computing this energy map;
    \item \textbf{Estimating the photon hits on the detector.} For a given camera image, $A_{ij}$, this steps calculates the most likely distribution of photon hits on the detector, $N^{\gamma}_{ij}$, providing a good estimate for the number of photons hitting each pixel $(i,j)$. The methods for performing this operation are discussed in Section \ref{Estimate of the photon hits distribution};
    \item \textbf{Constructing the spectrum.} This final step involves placing $N^{\gamma}_{ij}$ into the appropriate energy bin according to the energy map $E_{ij}$, with corrections applied for the varying solid angles associated with each energy level.
\end{enumerate}

\noindent The specifics of how the image is converted to a spectrum depends on the specific geometry of the spectrometer. While all the spectrometers relevant to this work are based on Bragg's law and this is ultimately used to construct the spectrum, in practice specific geometries (such as von H\'amos geometry which is given some attention later) may allow for simpler dispersion equations.









\subsection{\label{CCDs} Pixel Detectors}






\begin{figure}
    \centering
    \includegraphics[width=0.8\textwidth,keepaspectratio]{Figures/CCD Drawing.pdf}
    \caption[Scheme of a pixel detector]{Simplified scheme of a pixel detector section with illustration of charge spreading.
    }
    \label{ccd}
\end{figure}


Pixel detectors are integral components in the detection and measurement of light in various scientific applications, including spectroscopy, astronomy, and medical imaging. These detectors work by converting incoming photons into electronic signals that can be quantified and analysed. A simplified diagram of a pixel detector is shown in Fig.~\ref{ccd}. The general functioning of a pixel detector involves several stages: photon absorption, charge generation and collection, charge transfer, and signal readout. 

\begin{enumerate}
    \item \textbf{Photon absorption.} The process begins when photons enter the detector and strike the photosensitive surface, usually made of a doped or compound semiconductor (e.g. doped silicon, GaAs). At x-ray energies, photons interact with the semiconductor atoms via photoionisation, generating electron-hole pairs. The number of these pairs is approximately proportional to the energy of the photon hitting the surface, with random fluctuations quantified by the Fano factor \cite{lewenkopf2008numerical}. The efficiency in converting the photon energy into electrons energy is instead measured by the quantum efficiency (QE)~\cite{jone1959quantum}, which is essentially the probability that a photon is photoabsorbed, creating a photocurrent.
    \item \textbf{Charge Generation and Collection.} The generated electrons are collected in potential wells created by an array of electrodes on the detector surface. Each pixel in the detector corresponds to one potential well. The number of electrons collected in each well is roughly proportional to the energy contained in the light hitting that pixel, making the detector highly effective for capturing detailed images. However, thermal noise, which is intrinsic to the detector and depends on its temperature, can also generate electron-hole pairs, known as ``dark current'', which disturbs the measurement~\cite{newberry1991signal}. Additionally, the electron cloud generated by a photon absorption can diffuse over multiple pixels, causing image blurring. This diffusion is quantified by the Point Spread Function (PSF)~\cite{du2004effects}, which is a two-dimensional decaying function (typically a Gaussian, but sometimes distorted) with a characteristic width referred to as charge spreading radius ($R_s$)~\cite{abboud2013sub, hagino2019measurement}. The spreading radius can be estimated both theoretically and experimentally ~\cite{gatti1987dynamics, hagino2019measurement}, and depends on various parameters, including detector thickness and incoming photon energy. While reducing the sensor thickness can mitigate this spreading effect, it also decreases the probability of photon absorption, thereby impacting the quantum efficiency. Furthermore, even though $R_s$ is rigorously dependent on the photon energy — affecting the depth at which the photon is absorbed (see Fig. \ref{ccd}) — it is assumed to be constant in this thesis, as our analysis can be restricted to a portion of the detector illuminated by a sufficiently narrow photon energy range, thanks to Bragg spectroscopy.
    \item \textbf{Charge Transfer.} Once the exposure is complete, the collected charge needs to be transferred to the readout electronics. For modern hybrid detectors, this is done \textit{in situ} with the readout electronics connected to each individual pixel. In the case of CCDs, this is achieved through a process known as charge transfer, where the charges are moved sequentially through the CCD's structure to the output register. This transfer is carefully controlled to preserve the spatial information of the image. The loss of electrons occurring in this operation, which is not present in modern detectors, is measured by the so-called Charge Transfer Efficiency (CTE)~\cite{janesick1987charge, epperson1988applications}. For these reasons, hybrid detectors offer faster and more accurate readout than traditional CCDs.
    \item \textbf{Signal Readout.} In the final stage, the collected charges are read out by converting them into a voltage signal. This is done using a charge amplifier that translates the charge into a corresponding voltage. The voltage signal is then sent to an Analog-to-Digital Converter (ADC), which converts it into digital counts known as Analog-to-Digital units (ADU)~\cite{lathi1998modern}. These ADUs represent the total energy of the photons incident on each pixel, combined with noise, and form the basis of the digital image.
\end{enumerate}

\noindent For the following analysis, the relevant detector characteristics include $R_s$, the level of noise of each pixel, and the ADU counts associated to a photon of energy $E$, $ADU_{sp}(E)$, which in general depends upon the pixel. The noise on each pixel is measured simply by running the detector multiple times while the source is blocked (known as ``dark runs''~\cite{mullikin1994methods}), and then looking at the statistics of the empty frames for each pixel across the different shots. The dark runs are taken regularly during an experiment to measure the pixel noise as the detector conditions naturally change over time, and their mean - usually computed excluding tail values to account for anomalies (e.g. stray photons) - is then subtracted from the raw image to remove any bias on the ADU values due to thermal noise. Typically, this subtraction is executed by the facility's software, as part of automatic post-processing. Additionally, it is now common for facilities to provide the detector images, with the background subtracted, in energy units. The ADU to energy unit conversion of detectors, i.e. $ADU_{sp}(E)$, is already regularly measured as part of performance testing, and usually, this is automatically applied to the raw ADU image in a post-processing step. However, in this discussion, I always assume to work with background-subtracted images, but in ADU units, and I show how to compute $ADU_{sp}(E)$ directly from the experimental images. $R_s$ instead, is typically given in the camera documentation, but accounting for the charge spread between pixels remains a task for the user as it is computationally quite intensive, and not always necessary depending on the detector geometry or experiment. Lastly, image post-processing often includes masking or correcting faulty pixels, with special attention to edge pixels, where leakage currents reduce charge collection efficiency~\cite{shor2004edge}.

In the regions with low photon densities, which cover most of the detector in our cases, it is possible to identify a photon hit by searching clusters of neighbouring pixels with ADU counts consistently above the noise. This technique, known as single photon counting, allows one to count the photons on the detector, thereby reducing the uncertainties in estimating $N^{\gamma}_{ij}$ by exploiting the spatial correlation between ADU values of neighbouring pixels. The small regions with medium to high photon densities require instead more elaborated techniques to improve signal quality, which is discussed in the following. Although it is possible to estimate the energy of the photon directly from the ADU count, as done in SPCS~\cite{basden2003photon, morton1968photon}, working in dispersive (or Bragg) spectroscopy enables one to determine the energy of the photons by their positions, facilitating the use of the photon counting techniques.












\section{\label{Energy calibration} Energy Calibration}








For both spectrometers of interest, the energy map depends on the geometrical parameters of our setups (e.g. source-detector distance): 

\begin{eqnarray}
E_{ij} = E(x'_{ij},y'_{ij}; \boldsymbol{\Lambda}),
\label{energy_map}
\end{eqnarray}

\noindent where $x'_{ij}$ and $y'_{ij}$ are the coordinates of the centre of pixel $(ij)$ on the planar detector and $\boldsymbol{\Lambda}$ collects all the geometrical parameters. These parameters are used as fitting variables to find the energy map, ensuring that it is consistent with known experimental lines. The two types of geometries are now examined separately.










\subsection{\label{flat crystal spectrometer}Flat Crystal Spectrometer}









\begin{figure*}
    \centering
    \includegraphics[width=\textwidth,keepaspectratio]{Figures/Flat crystal spec.pdf}
    \caption[Scheme of flat crystal spectrometer geometry]{Scheme of flat crystal (FC) spectrometer geometry with the illustration of source broadening (inset). The laboratory frame is defined by $\hat{z}$ being the normal of the crystal, $\hat{x}$ such that $\boldsymbol{r_0}$ is contained in the $\hat{x}-\hat{z}$ plane, and $\hat{y}$ chosen for the frame to be an orthonormal right-handed triad. The camera frame is defined with $\hat{x}$ and $\hat{y}$ aligned with the pixel array, $\hat{z}$ orthogonal to the camera plane, and the origin $O'$ at the centre of the detector.
    }
    \label{Flat crystal spec}
\end{figure*}


Starting from the assumption that the rocking curve of the reflecting crystal~\cite{zachariasen1967general} is a Dirac delta centred at the Bragg angle, $\theta_B$, and that the sample acts as a point source, thereby neglecting source broadening.
Under these assumptions, using Bragg's law~\cite{zachariasen1967general} and referring to Fig. \ref{Flat crystal spec}, one can derive the following relation between $E_{ij}$, the coordinates of the pixel and the geometrical parameters~\cite{humphries2020isochoric}:

\begin{equation}
\begin{split}
E_{ij} & = E(\theta_{B}(x'_{ij},y'_{ij}; \boldsymbol{\Lambda}) \\
& =  \frac{hc}{2d \sin(\theta_{B}(x'_{ij},y'_{ij}; \boldsymbol{\Lambda}))} \\
& = \frac{hc}{2d}\frac{\|\boldsymbol{r_0}+\boldsymbol{Rr'_{ij}}\|}{\boldsymbol{\hat{z}}\cdot(\boldsymbol{r_0}+\boldsymbol{Rr'_{ij}})} \ \ ,
\label{energy map equation}
\end{split}
\end{equation}

\noindent where $\theta_{B}$ is the Bragg angle associated with pixel $(i,j)$, $h$ is the Planck constant, $c$ is the speed of light and $d$ is the spacing between the reflecting planes of the crystal. In this thesis, the first-order reflection ($n=1$) is considered, but the generalisation of equations to any arbitrary order $n$ is straightforward. The geometrical quantities $\boldsymbol{r_0}$, $\boldsymbol{r'_{ij}}$, $\boldsymbol{R}$, are given by the following expressions (see Fig. \ref{Flat crystal spec}):

\begin{align}
& \boldsymbol{r_0} = D \ (\cos\bar{\theta}_B \ , \ 0 \ , \ \sin\bar{\theta}_B) \ , \\
& \boldsymbol{r'_{ij}} =  (x'_{ij} \ , \ y'_{ij} \ , \ 0) \\
& \boldsymbol{R} = \boldsymbol{R}_x(\alpha_x)\boldsymbol{R}_y(\alpha_y)\boldsymbol{R}_z(\alpha_z) \ .
\label{Geometrical quantities}
\end{align}

\noindent Here, $\boldsymbol{R}$ is the rotation matrix that specifies the orientation of the camera relative to the laboratory frame, parametrised with the Euler's angles $(\alpha_x, \alpha_y, \alpha_z)$~\cite{goldstein:mechanics}. Such angles, along with the distance $D$ from the source to the camera origin $O'$ and the Bragg angle $\bar{\theta}_B$ associated to $O'$, constitute the geometrical parameters $\boldsymbol{\Lambda}$ for fitting and determining the energy map. In order to do this, one first identifies, in terms of pixels coordinates, a bright curve $S_{exp}$ on the experimental image $A_{ij}$ associated with an emission line of known energy $\bar{E}$:

\begin{equation}
\begin{split}
& S_{exp} = \biggl\{(i, f(i)) \ , \ i = 1, \dots, N_y \biggr\} \\ 
& \mathrm{with}  \ f(i) = \max_{N_1 \leq j \leq N_2} A_{ij} \ , 
\label{experimental line1}  
\end{split}
\end{equation}

\noindent where $N_1$ and $N_2$ are properly chosen to isolate the emission line of energy $\bar{E}$. A more sophisticated and accurate method is to fit, for each row, the emission line to an analytic function and to take the peak of the fit as the line position. This approach is less susceptible to noise and allows for sub-pixel resolution, by constructing the energy map on a finer sub-pixel grid. However, for the purposes of this thesis, the method described by Eq. (\ref{experimental line1}) is adopted. As a second step, the relative theoretical contour $\{(ij)| \ E_{ij} = \bar{E}\}$ is computed using Eq. (\ref{energy map equation}) with a guess for $\boldsymbol{\Lambda}$:

\begin{equation}
\begin{split}
& S_{th}(\boldsymbol{\Lambda}) = \biggl\{(i, g(i;\boldsymbol{\Lambda})) \ , \ i = 1, \dots, N_y \biggr\} \\ 
& \mathrm{with}  \ g(i;\boldsymbol{\Lambda}) = \min_{j} |E_{ij}(\boldsymbol{\Lambda})-\bar{E}| \ . 
\label{experimental line}  
\end{split}
\end{equation}

\noindent $S_{exp}$ and $S_{th}$ are parametrised with $i$ as independent variable because of the expected shape for these curves. Finally, one calculates a loss function $\mathcal{L}$, dependent on the geometrical parameters $\boldsymbol{\Lambda}$, measuring a distance between the experimental line $S_{exp}$ and the theoretically-computed contour $S_{th}$:

\begin{equation}
\mathcal{L}(\boldsymbol{\Lambda}) = \frac{1}{N_y}\sum_{i=1}^{N_y} (f(i)-g(i;\boldsymbol{\Lambda}))^2
    \label{loss function}
\end{equation}

\noindent At this point, $\mathcal{L}$ is minimised over the space of the geometrical parameters with standard optimisation routines~\cite{gao2012implementing}. The loss function can include multiple emission lines by summing the loss functions of each individual line, as previously defined. Each line in the summation can also be weighted proportionally to its intensity. Figure \ref{FC_calibration} illustrates an example of this calibration for the Mg spectrum discussed earlier. Here, the geometrical parameters of the setup are determined by fitting the $K\alpha$ lines of the Mg ion in its fundamental charge state, and they are as follows:

\begin{equation}
\begin{split}
& (\alpha_x, \alpha_y, \alpha_z) = (-1.9^\circ,55.8^\circ,0^\circ) \\
& D =  79.5 \ \mathrm{mm} \\
& \bar\theta_B = 34.7^\circ\ \ .
\label{Mg parameters}
\end{split}
\end{equation}

\noindent As it is observed from these values, it is common practise to work with the camera orthogonal to $\boldsymbol{r_0}$ in order to maximise the solid angle detected by the detector, as well as to minimise any grazing incidence broadening in the sensor layer (see Section \ref{sec:Resolution}). In this work, the uncertainties in the $\boldsymbol{\Lambda}$ parameters arising from the fitting process are disregarded, as their impact on the final energy map is minimal. This is justified by the fact that these errors are typically small, especially when multiple emission lines are used for fitting. Consequently, their contribution to the energy map uncertainties is negligible compared to the broadening effects discussed in Section \ref{sec:Resolution}.

\begin{figure}[t]
    \centering
    \includegraphics[width=0.8\textwidth,keepaspectratio]{Figures/FC calibration_v2.pdf}
    \caption[Example of flat crystal spectrometer calibration]{Detector image of MgF$_2$ emission spectrum with the theoretically-computed energy contour (red) fitted to a Mg $K\alpha$ line.
    }
    \label{FC_calibration}
\end{figure}





















%  Note that we took a first-order Bragg reflection, which is the typical case and therefore the only one considered in this work. 








\subsection{\label{Von Hamos spectrometer}Von H\'amos Spectrometer}









\begin{figure*}
    \centering
    \includegraphics[width=\textwidth,keepaspectratio]{Figures/Von Hamos spec.pdf}
    \caption[Scheme of the von H\'amos spectrometer]{Scheme of the von H\'amos spectrometer with an illustration of mosaic focusing (inset). The mosaicity of the crystal increases the area that reflects a certain photon energy. For the correct functioning of the spectrometer, the crystal must be placed halfway between the source and the centre of the detector. As in the flat crystal case, $\bar{\theta}_B$ is defined as the Bragg angle associated with the centre of the detector.
    }
    \label{VH_scheme}
\end{figure*}

In Von H\'amos (VH) geometry~\cite{vonHamos_1932}, an axis (the dispersion axis) is drawn between the source and the detector, and the detector surface lies along this line (see Fig.~\ref{VH_scheme}). A cylindrically-bent crystal is then placed parallel to this axis, with the crystal positioned at its radius of curvature below this axis. The crystal therefore acts to focus the reflected rays back on to the dispersion axis at a distance twice as a far away in the dispersive direction. This results in a simple dispersion relationship between the photon energy $E$, the position from the source to the front of the camera along the dispersion axis $L$, and the position on camera along the dispersion axis $X_{\rm cam}$:
\begin{equation}
    E = \frac{hc}{2d}\sqrt{1 + \frac{(L+X_{\rm cam})^2}{4R^2} } \, ,
\label{VH dispersion}
\end{equation}

\noindent where $R$ is the radius of curvature of the crystal, and the remaining symbols have the same meaning as before.
The calibration is also quite simple and involves the optimisation of $R$ and $L$ using analogous techniques to those described for the flat crystal geometry. Bright emission lines with known energies are well-suited to this task. Two sets of lines that cover the detector range are ideal for the calibration and results in a very well-constrained dispersion relationship. An example of an experimental calibration is shown in Fig.~\ref{VH_calibration}.
% [...] but considering only a narrow stripe around the dispersion axis rather than the entire detector (see Fig. \ref{VH_calibration_total.pdf}).
\begin{figure}[t]
    \centering
    \includegraphics[width=0.8\textwidth]{Figures/VH_calibration_total.pdf}
    \caption[Example of von H\'amos spectrometer calibration] {Detector image of Fe emission using a von H\'amos spectrometer from the experiment reported in Ref.~\cite{forte2024resonant} with a zoom (inset) around the central dispersion axis, indicated by the red dotted line. The two vertical bright lines in the inset correspond to the isoenergy contours of the Fe K-alpha lines used for the calibrations. In reality, the greater the energy separation between these reference lines, the more accurate the calibration becomes. Generally, the central dispersion axis may be slightly tilted relative to the detector's edge. This tilt is accounted for by adjusting the region of interest while moving along the dispersive direction. However, Eq. (\ref{VH dispersion}) remains valid, as corrections due to the tilt are typically negligible.}
    \label{VH_calibration}
\end{figure}

The primary benefit of a von H\'amos spectrometer is its collection efficiency. By focusing all the rays onto a single axis instead of arcs as with the flat crystal spectrometer, it is easier to detect events above the noise.
As the bending process deforms the lattice of a perfect crystal, its rocking curve becomes extended, which degrades the intrinsic resolution of the crystal~\cite{taupin1964theorie,uschmann1993x}. It is therefore common to use mosaic crystals as the dispersive crystal to further increase the collection efficiency of the spectrometer.
Common examples of such crystals are highly annealed pyrolytic graphite (HAPG), highly oriented pyrolytic graphite (HOPG) and Lithium Fluoride (LiF). A mosaic crystal is one that consists of small perfect crystallites with their normal vectors randomly distributed around the crystal surface normal. As a result, in contrast to perfect crystals where reflections only take place around the Bragg angle within the rocking curve, the Bragg condition can be satisfied anywhere on a mosaic crystal by finding a crystallite appropriately angled to the normal (see Fig.~\ref{VH_scheme}). This effect increases the area of the crystal reflecting a certain photon energy on a specific spot of the detector (mosaic focusing). The end result is that von H\'amos spectrometers have very high collection efficiencies, and so are often employed in situations where photon numbers are very low, such as in x-ray Thomson scattering and resonant inelastic x-ray scattering experiments~\cite{gawne2024effectsmosaiccrystalinstrument,humphries2020probing,forte2024resonant}.

As the x-rays are focused onto the dispersion axis, the spectrum is produced by integrating over a few pixels along the non-dispersive direction. The number of pixels used is chosen to maximise the signal-to-resolution ratio (see~\cite{preston}). Outside this region of interest (RoI), unfocused x-rays are disregarded, as different photon energies overlap on the same pixels.

The drawback to mosaic crystal VH spectrometers is that their resolution tends to be worse than for flat crystal spectrometers. The bulk of this is due to the mosaicity of the crystal, which both enhances depth broadening and itself contributes directly to the broadening. As a result, photons of a specific energy are sent to different locations on the detector rather than to a precise, nominal point, producing a quite broad, X-shaped instrument function, as demonstrated in ray-tracing simulations~\cite{gawne2024effectsmosaiccrystalinstrument}. This broadening leads to the aforementioned energy mixing outside the central RoI. In the experimental image presented here, however, these characteristic shapes are not visible, as broadening is dominated by defocusing, which deforms the energy contours towards those of a FC spectrometer, creating small and thick bright arcs (see Fig. \ref{VH_calibration}). Focal aberrations are in fact an additional factor that degrades the resolution~\cite{2012_Zastrau_Focal}, but this effect can be reduced by using a relatively narrow crystal~\cite{2012_Zastrau_Focal,preston}. Typical resolutions of von H\'amos spectrometers range from 1 to 10s eV~\cite{preston,glenzer2016matter,MacDonald_PoP_2021}.



%$forming for each energy characteristic banana-shaped regions, as a result of the range of transverse angles at which photons hit the mosaic cylindrical crystal~\cite{szlachetko2012hamos} (see the bright areas produced by the Fe K-alpha lines in Fig.~\ref{VH_calibration}).

% which creates banana-shaped regions resembling the energy contours of the flat crystal spectrometers

% LCLS run:340













\subsection{\label{sec:Resolution}Broadening Effects}









So far it has been assumed that there is a one-to-one relationship between the photon energy and the position it lands on the detector. However, in practice there are numerous effects that result in photons of the same energy being detected at different points on the detector away from their nominal position.
As the spectrometer is calibrated so that spatial points on the detector correspond to specific energies, there is no way to distinguish between a photon arriving at its nominal position, and a photon with a different energy landing at the same point due to the broadening effects.
This point spread function, also known as the instrument function (IF), means that the features of the measured spectrum appear broader than they actually are.

While high resolution spectra are generally desirable, there is a trade-off between resolution and signal intensity as the very effects that result in broadening also tend to result in higher spectrometer efficiencies. For example, a broader rocking curve reduces resolution but increases crystal reflectivity, while a larger pixel size leads to greater quantisation errors yet allows more photons to be detected by a single pixel, yielding a better SNR.

Broadening effects are generally complex to account for, as the instrument function (IF) can be strongly geometry-dependent and varies between experiments and setups. Furthermore, the IF rigorously depends also on the photon energies incident on the crystal~\cite{gawne2024effectsmosaiccrystalinstrument}. In other words, the broadening of the spectrum due to the IF cannot be expressed as a simple convolution, but must be formulated as a Fredholm integral equation (FIE):

\begin{equation}
I(E) = \int K(E, E_0) I_{orig}(E_0) dE_0,
    \label{FIE}
\end{equation}

\noindent where $I$ and $I_{orig}$ represent the measured and original spectra, respectively, and the kernel $K(E,E_0)$ is the energy-dependent instrument function. Methods for determining $K$ are outlined in Ref.~\cite{gawne2024effectsmosaiccrystalinstrument} and rely on a combination of quasi-elastic signal measurements and theoretical modelling of the IF. Once $K$ is known, $I_{orig}$ can be retrieved using standard techniques for solving FIEs. Even in cases where the dependence of the IF on the incident energy could be neglected -- making the calculation of $I_{orig}$ a simpler deconvolution problem -- finite spectral windows and experimental noise introduce instability in the deconvolution process. As a result, a small amount of noise in $I$ can cause significant distortions in $I_{orig}$.

Here, the main effects that deteriorate the resolution of the measured spectra are summarised.

\begin{itemize}
    \item \textbf{Finite pixel size and quantisation error:} Pixels are rectangles with sides of finite length $L_p$, causing each pixel to collect photons within an energy range of order $O(\frac{2hc}{d}\frac{\cos\bar\theta_B}{\sin^2\bar\theta_B}\frac{L_p}{D})$.
    In other words, the energy calibration of a pixel $E_{ij}$ actually represents an average energy.
    For the cameras used in the experiments, typical $L_p$ values are in the range $\sim$10-100~$\mu$m~\cite{mozzanica2018jungfrau}, which is typically small enough to measure a smooth spectrum with energy bins $\lesssim 2$~eV. The energy bin width $\Delta E$ must in fact be chosen as small as possible to reduce the quantisation error ($\sim\Delta E/\sqrt{12}$), but large enough to contain several $E_{ij}$ in order to obtain smooth spectra. As an aside, there is a trade-off to be considered in pixel size -- while smaller pixels have less pixel broadening, they also suffer from increased charge spreading between pixels. If the charge of a photon is spread between too many pixels, it may become difficult to detect above the noise. As a result, there is a balance between the possible energy bin size and the ability to reliably detect photons above the noise. These two broadening effects, unlike the ones discussed next, degrade resolution without affecting the correspondence between photon energy and its position on the detector.
    \item \textbf{Source broadening (SB):} As the incident beam has a finite size, photons leave the target from a region on the sample surface with an effective size $D_s$, typically ranging from several nm to dozens of $\mu$m (this region is often approximated as a circle, although imprints of the beam show that the FEL beam spatial profile is decidedly non-trivial~\cite{Chalupsky_Opt_2010}). Additionally, as samples have finite thickness, photons also emerge from different points from within the sample.
    Consequently, photons of the same energy emitted from different points in the sample reach the detector at different points, representative of the shape of the source, which results in photons of different energies landing on the same point, and therefore, in the measured spectrum being broadened (see Fig.~\ref{Flat crystal spec}). The energy uncertainty on each pixel, due to this effect, is of order $O(\frac{2hc}{d}\frac{D_s \cos\bar\theta_B}{D})$.
    \item \textbf{Depth Broadening (DB):} A photon can travel a finite depth into a crystal before being reflected back out. As a result, it reflects at a point displaced from the nominal position it would reflect on the surface, and therefore hits the detector at a different point. The DB resembles an exponential function as it represents the probability a photon reflects from a certain depth in the crystal. For perfect crystal, this is typically an extremely narrow function and does not contribute significantly to the broadening. However, for mosaic crystals the broadening is quite substantial (several eV) due to the distribution of the crystallites~\cite{gawne2024effectsmosaiccrystalinstrument}, and is made even broader by non-negligible multi-reflection paths within the crystal~\cite{Schlesiger_JAC_2017}. When the path travelled inside the crystal is much smaller than its attenuation length - a necessary condition to maintain the crystal's reflectivity - the SB is proportional to $t/sin(\theta_B)$, with $t$ representing the crystal thickness. Since $t$ cannot be reduced excessively without compromising reflectivity, it is important to avoid working at overly shallow incidence angles. The same consideration applies to the detector, where small incidence angles influence the depth and, depending on $L_p$, the pixel position of the photoabsorption, thereby impacting pixel detection and increasing charge spreading~\cite{dusi2003study}.
    \item \textbf{Finite width of the rocking curve:} Another source of energy broadening is the finite width of the crystal rocking curve (RC), $D_W$\footnote{This notation stems from its common name in dynamical diffraction theory, the \textit{Darwin width}.}. The RC allows for a photon that does not satisfy the Bragg condition to nevertheless be reflected from the crystal, albeit with a lower reflectivity~\cite{zachariasen1967general}. For perfect crystals, their intrinsic RC is calculated using dynamic diffraction theory, and it is an extremely narrow function~\cite{zachariasen1967general}. Typical values for $D_W$ are dozens of $\mu$rad. For mosaic crystals, the mosaic distribution of the crystallites also contributes to the RC, and is about an order of magnitude larger than the intrinsic RC of the crystallites~\cite{Gerlach_JAC_2015}.
    The energy uncertainty per pixel from this effect is then simply of order $O(\frac{2hc}{d}D_W)$. $D_W$ is generally weakly dependent on the photon energy~\cite{Gerlach_JAC_2015}, and is often assumed to be constant for typical spectrometer energy windows.
\end{itemize}

\noindent For the geometrical setup of the Mg experiment, summarised by the parameters provided in~(\ref{Mg parameters}), with typical values for $L_p$, $D_s$ and $D_W$, the energy resolution losses produced by these effects are all around $0.5-1$ eV. Normally in HED science, these resolution losses are of secondary importance, as other factors such as source bandwidth dominate~\cite{feldhaus1997possible}. Nevertheless, they must be considered if one wants to reach high resolutions (sub eV). These resolutions can be achieved by means of specific experimental setups, characterised by large source-detector distances and high resolution crystals with very narrow RC~\cite{wollenweber2021high, mcbride2018setup}. However, these features also tend to decrease the SNR of the single-shot images. It is therefore useful to contrive alternative methods that mitigate these resolution losses, without compromising the signal intensity. One potential solution to overcome the resolution limitations imposed by the finite pixel size is to analyse the charge distribution between neighbouring pixels. Techniques used to manage charge spreading, based primarily on centroid calculations, have been successfully employed in the past to handle such correlations and achieve sub-pixel resolutions~\cite{abboud2013sub, hagino2019measurement}. Nevertheless, it is important to note that these methods are applicable only when the pixel dimensions are sufficiently small compared to $R_s$ and in low photon-density regimes. Section~\ref{Estimate of the photon hits distribution} introduces an alternative approach to treat these correlations, based on deconvolution algorithms and capable of handling also cases with large photon densities on the detector.

The other three effects can instead be mitigated by studying the instrument function of the experimental setup and then solving the integral equation mentioned earlier. Assuming the pixel size effect is negligible or that it has already been accounted for, the broadening due to the spectrometer instrument function can be determined by irradiating the target with a monochromated source and comparing the elastic scattering measured on the spectrometer to a high-resolution beam spectrometer. For this method to work, a material with very weak inelastic scattering in the vicinity of the elastic scattering must be chosen -- good choice include plastics such as polymethyl methacrylate and polypropylene, or metallic glasses~\cite{Voigt_PoP_20201,gawne2024effectsmosaiccrystalinstrument}.
%The monochromatisation must be sufficient to neglect the broadening due to the source bandwidth.
Another common way to determine the resolution of a setup is from emission lines with known lineshapes; e.g. H\"olzer \textit{et al.}'s $3d$ metal transition lines~\cite{Hoelzer_PRA_1997}. In any case, the FWHM of the instrument function used provides an estimate of the spectrometer resolution.
























\section{\label{Experimental estimation of camera parameters}Experimental Estimation of Camera Parameters}








A fundamental quantity for converting the ADU map to the photon distribution hits is the RMS of the camera read noise, $\sigma_{N}$, which measures the level of thermal noise in the detector. This can be measured pixel-wise by simply running the camera without the source and calculating the variance of the pixel values during these dark runs. However, for simplicity, this quantity is assumed to be constant across the camera, which is a reasonable approximation for the typical detectors employed~\cite{mozzanica2018jungfrau, carini2016epix100}. In these cases, $\sigma_{N}$ is directly inferable from the histogram of the camera images with signal present, whose shape is now discussed.

This work deals with detector images that have undergone dark image subtraction, so that the expected value of the ADU count for a pixel not hit by any photon is zero. Additionally, for our cases, the analysis of the histogram is simplified because, in Bragg spectroscopy, photons energies are determined by their position on the detector. Consequently, the camera can be subdivided into regions (e.g. rectangles that tile the entire camera) sufficiently small to approximate $ADU_{sp}(E)$, the ADU counts produced by a photon of energy $E$, as constant ($ADU_{sp}$) in each region. The shape of the histogram for a single region is in general a series of peaks centred approximately in the following sequence of values:
\begin{eqnarray}
P_{n} = \zeta n ADU_{sp}, \ \ \ \ \ n = 0,1,2,\dots
\label{ADU sequence}
\end{eqnarray}

\noindent where $\zeta$ is a constant smaller than 1 (typically around 0.8-0.9), whose origin is discussed below. The n-th peak contains the pixels which have been likely been hit by $n$ photons. In reality, only the initial peaks of this sequence appear in the histogram, depending on the density of photon hits, which is measured by the so-called fill fraction $\lambda$. This important parameter is defined as:
\begin{eqnarray}
 \lambda = \frac{\# \ \textrm{photons on the camera}}{\# \ \textrm{pixels in the camera}} \ \ .
\label{fill fraction}
\end{eqnarray}

\noindent Each peak is then broadened by the intrinsic noise of the detector and the charge spreading, quantified by $\sigma_{N}$ and $R_s$, respectively. The magnitude of these broadening effects determines whether the peaks are merged together or remain distinct. Additionally, charge spreading systematically shifts the peaks, initially centred at $n ADU_{sp}$, towards lower ADU values, as some of the photon energy is deposited across multiple pixels. The constant $\zeta$ in Eq. (\ref{ADU sequence}) accounts for this shift, decreasing as $R_s$ increases. To illustrate these effects, Figs.~\ref{histogram}~(a-d) present histograms generated from synthetic camera images with various levels of noise and fill fractions. More details on the fabrication of these synthetic images are provided in Section~\ref{Synthetic data}. Alongside these synthetic data, an experimental histogram for MgF$_2$ is included in Fig.~\ref{histogram}e for data from Ref.~\cite{Gawne_PRE_2023}.

\begin{figure}[t]
    \centering
    \includegraphics[width=0.6\columnwidth,keepaspectratio]{Figures/histograms.pdf}
    \caption[Histograms of synthetic and experimental images]{(a-d) Synthetic histograms of 700x700 images for various levels of noise and fill fractions. $ADU_{sp}$ is taken equal to 120. (e) Histogram of the experimental image reported in Fig.~\ref{fig:general process}, acquired using a PIXIS-XF camera. This camera employs 16-bit analog-to-digital converters, enabling quantisation into 65,536 disrete signal levels.}
    \label{histogram}
\end{figure}



If $ADU_{sp}(E)$ is not provided by the facility, it can be directly inferred from the position of the first peak in the series given by Eq. (\ref{ADU sequence}), which corresponds to single photons hitting pixels on the camera.
Similarly, from the width of the peak centred around 0, one can estimate $\sigma_{N}$. Intuitively, the more the peaks merge together, the more challenging it becomes to accurately determine these quantities. In such cases, one approach to extract these two parameters is to fit the experimental histograms to a simulated one. The latter can be constructed (see Section \ref{Synthetic data}) given certain values for $ADU_{sp}$, $\sigma_N$, and $\lambda$, which can be estimated from the expected number of photons, and $R_s$, usually provided in the camera documentation. 


















 







\section{\label{Estimate of the photon hits distribution}Estimate of the Photon Hits Distribution}








The objective of this section is to describe how to compute the distribution of photon hits, $N^{\gamma}_{ij}$, from the camera image $A_{ij}$, thereby eliminating the distortions due to thermal noise and charge spreading. Depending on the spectrum of photons hitting the camera, the image has regions with different local fill fractions (see Fig.~\ref{synthetic images}), which must be treated with distinct methods. Therefore, to process the image, the following two techniques are applied:

\begin{enumerate}
    \item \textbf{Photon counting algorithms:} To address regions with small fill fractions ($\lambda \lesssim 0.1$), where single photon clusters are present, photon counting algorithms are used. These algorithms can identify the single photon clusters, exploiting the correlation between neighbouring pixels to reduce the uncertainty on the photon detection.
    \item \textbf{Image deconvolution techniques:} To treat the regions with large fill fractions ($\lambda \gtrsim 1$), where agglomerated clusters are present (see next section), the image is processed with the Richardson-Lucy deconvolution method~\cite{richardson1972bayesian}.
\end{enumerate}

\noindent As discussed below, these two operations can be carried out sequentially, improving the signal quality across the entire range of $\lambda$, including also intermediate fill fractions. Following the discussion in Section~\ref{Experimental estimation of camera parameters}, the camera parameters $ADU_{sp}$, $\sigma_N$ and $R_s$ are considered known quantities in the following analysis.










\subsection{\label{Photon counting algorithms} Clustering Algorithms}





\begin{figure*}
\begin{center}
    \includegraphics[width=0.5\textwidth,keepaspectratio]{Figures/layout.pdf}
    \caption[Layouts for a single photon cluster]{The 4 layouts for a single photon cluster. The brighter pixel represents the pixel hit by the photon. Each of the last 3 can rotate by $90^{\circ}$, $180^{\circ}$ and $270^{\circ}$ to produce 13 layouts in total.
    }
\end{center}
\label{clusters pattern} 
\end{figure*}


If one supposes, as it is often the case, that the charge produced by a photon does not spread further than a pixel side length, then the shapes of single-photon clusters (i.e. the pixels in which the charged produced by the photon has leaked) are few, and given in Fig. 2..
As there are only 13 possible shapes, a simple clustering algorithm can be employed to identify these clusters and re-accumulate the intensity back on to the central pixel. The benefit of doing so is that this makes the pixel brighter relative to the background, allowing for more of these events to be included once the final thresholding is performed.

The first step to identify the photon hits from the ADU map is thresholding the camera image. Although different choices are possible for the thresholding value, here, 1.5 $\sigma_{N}$ is chosen, as a good compromise to cut the noise and retain the signal:
\begin{equation}
A^{th}_{ij} = 
\begin{cases}
A_{ij} \ \ \ \ \ A_{ij}\ge 1.5 \ \sigma_{N} \\
0 \ \ \ \ \ \ \ \ A_{ij}< 1.5 \ \sigma_{N}
\end{cases}
\label{thresholding}
\end{equation}

\noindent The thresholding value can be chosen in a systematic way by analysing their performances over multiple synthetic images for a given set of parameters ($ADU_{sp}$, $\sigma_N$ and $R_s$). This method, adopted for this thesis, is briefly described in the following section. After the thresholding, the pixels with an ADU value different from zero are the ones that most likely received a contribution to their ADU by a photon hit. 
One may also set pixels exceeding some upper threshold to zero to mask out possible faulty pixels or high intensity events that do not correspond to real signal. 

Next, all the clusters are located and it is determined whether they have fall into one of two cases:
\begin{enumerate}
    \item \textbf{Single-photon clusters.}
    If the shape of the cluster $C$ matches that in Fig. 2.8, then the charge spread is considered to originate from the brightest pixel, and all the ADU values of the cluster are accumulated onto its brightest pixel ($\bar{i},\bar{j}$):
    \begin{equation}
        \begin{split}
            & ADU_{\bar{i},\bar{j}}' = \sum_{ij\in \textrm{C}}A_{ij}^{th} \\
            & A_{ij}'=0 \ \ \textrm{for } (i,j)\in \textrm{C and } (i,j) \ne (\bar{i},\bar{j}) \ .
        \end{split}
        \label{cluster sum}
    \end{equation}
    \item\textbf{Agglomerated clusters.}
    In this case the ADU values are not modified:
    \begin{equation}
        A_{ij}' = A_{ij}^{th}
        \label{cluster sum 2}
    \end{equation}
\end{enumerate}

\noindent After the clustering algorithm is applied, the image is thresholded once more with a higher threshold, assuming that any events that fall below the threshold did not originate from a single photon event. Using the method described in the following section, $(ADU_{sp} - \ \sigma_N)$ is chosen as the second threshold, which results in discarding roughly 20\% of the single photons ($A'_{ij}$ has a larger variance than $A_{ij}$), but eliminates almost entirely noise fluctuations. This second threshold can also be systematically determined using the method outlined in the next section. Charge accumulation over the brightest pixel can also be carried out by considering the centroid of the single-photon cluster, enabling subpixel resolution~\cite{abboud2013sub}. However, this more advanced process is beyond the scope of this thesis.

Fig.~\ref{algorithms} presents an example of this operation using a synthetic camera image (see Section~\ref{Synthetic data} for more details). As it is observed, agglomerated clusters remain unchanged by this algorithm and are addressed in the following section.

\begin{figure*}
    \centering
    \includegraphics[width=0.85\textwidth,keepaspectratio]{Figures/Algorithms.pdf}
    \caption[Illustration of the various steps in the camera image processing]{Illustration of the various steps in the camera image processing. (a) Map of photon hits on the camera for a small $\lambda$, to which thermal noise and charge spreading have been applied to obtain the camera image (b). The same image, processed with the clustering counting algorithm (c), and then with the deconvolution algorithm (d). In the figure, an example of a single-photon cluster and of an agglomerated cluster are indicated. All camera images, except the first initial photon map, are in ADU units.}
    \label{algorithms}
\end{figure*}

Exploiting the correlation of ADU values across neighbouring pixels enables one to boost the signal. As an example, consider three adjacent pixels arranged in an L-shape, with ADU values between the first and second thresholds. If the correlation is ignored, and consequently, the ADU values are not accumulated on the brightest pixel, all three pixels are zeroed after the second thresholding, resulting in a loss of signal that likely corresponds to a photon. In fact, a simple analysis shows that the probability of three L-shaped pixels with values above 1.5 $\sigma_N$ being caused by a photon is significantly higher than by noise. 

In addition to this, thresholding mitigates the noise-related fluctuations in the reconstructed spectrum. The noise reduction is approximately proportional to $\mathcal{P}(N^{\gamma}_{ij} \geq 1) < 1$, where $\mathcal{P}(N^{\gamma}_{ij} \geq 1)$ is the probability that at least one photon lands on pixel $(i,j)$. This result is derived by considering the variance of the random variables $A_{ij}\approx N^{\gamma}_{ij}ADU_{sp} + \eta_{ij}$ (see Section \ref{Probability distribution of the pixel values}) and $A^{th}_{ij}$ (as defined earlier), where $\eta_{ij}$ represents the random noise on pixel ($i,j$), which is assumed to be Gaussian-distributed with mean 0 and variance $\sigma_N$: $\eta_{ij}\sim \mathcal{N}(0, \sigma_N)$. Firstly, fixing the number of photons $N^{\gamma}_{ij}$, one finds:





\begin{equation}
    \begin{split}
            & \mathbb{V}(A_{ij}|N^{\gamma}_{ij}) = \mathbb{V}(\eta_{ij}) = \sigma_N^2 \ \ \ \ \ \ \ \forall N^{\gamma}_{ij}, \\
            & \mathbb{V}(A^{th}_{ij}|N^{\gamma}_{ij}) =  
            \begin{cases}
            \sigma_N^2 \ \ \ \ \ A_{ij}\ge \kappa \sigma_{N} \\
            0 \ \ \ \ \ \ \ \ A_{ij}< \kappa \sigma_{N}
             \end{cases} ,
    \end{split}
    \label{bonus formula}
\end{equation}
where $\mathbb{V}(\cdot)$ represents the variance and $\kappa \sigma_N$ is a generic threshold proportional to $\sigma_N$. Now, taking the mean values $\mathbb{E}$ over the random variable $N_{ij}^{\gamma}$, one obtains:

\begin{equation}
    \begin{split}
            & \mathbb{E}(\mathbb{V}(A_{ij})) =  \sigma_N^2 \mathcal{P}(N^{\gamma}_{ij}=0) + \sigma_N^2 \mathcal{P}(N^{\gamma}_{ij}=1) + \cdots = \sigma_N^2 \ , \\
            & \mathbb{E}(\mathbb{V}(A^{th}_{ij})) = 0 \cdot\mathcal{P}(N^{\gamma}_{ij}=0) + \sigma_N^2 \mathcal{P}(N^{\gamma}_{ij}\ge1) < \sigma_N^2,
    \end{split}
    \label{bonus formula 1}
\end{equation}
having taken as an approximation that $\mathcal{P}(A_{ij}<\kappa\sigma_N|N^{\gamma}_{ij}>1) = 0$ and $\mathcal{P}(A_{ij}>\kappa \sigma_N|N^{\gamma}_{ij}=0) = 0$. Summing the ADU values over the energy contours preserves the ratio between the variance of the raw and thresholded values.








\subsubsection{\label{Determination of thresholding values}Determination of Thresholding Values}





The two thresholding values can be determined using the following procedure. First, one generates $N_s$ synthetic images (see also Section \ref{Synthetic data}) for specified parameters ($ADU_{sp}$, $\sigma_N$, $R_s$), and for a small fill fraction (e.g. $\lambda=0.1$), so that single-photon clusters are predominant. It is important to note that generating these synthetic images requires knowledge of the underlying spectrum. For experimental applications, a reasonably accurate approximation of the underlying spectrum can be obtained from the experimental images prior to optimizing the thresholding values. Each of these synthetic images is then processed with the photon clustering algorithm, applying a trial pair of thresholding values, and divided by $ADU_{sp}$ to reconstruct the photon hits map $\{N^{rec}\}$. The reconstructed photon maps are then compared to the original photon map $\{N^{orig}\}$, and a global loss function is computed:

\begin{equation}
\mathcal{L}oss = \frac{1}{N_s} \sum_{l = 1}^{N_s} ||\{N^{rec}\}_l - \{N^{orig}\}||_F \ \ ,
\label{appendix 0}
\end{equation}

\noindent where $\{N^{rec}\}_l$ and $\{N^{orig}\}$ are two matrices and $||\cdot||_F$ is the Frobenius norm. The optimal threshold pair is obtained as the one that minimises this loss function.













\subsection{\label{Images devonvolution}Images Deconvolution}











To address the charge spreading and the noise contributions in regions with large $\lambda$, the Richardson-Lucy (RL) deconvolution algorithm~\cite{richardson1972bayesian} is applied to the matrix $A'_{ij}$. This method is an iterative procedure, based on Bayesian inference, that computes the most likely expected values $V_{ij}$ of the Poissonian and independent variables $X_{ij}$, given a realisation $R_{ij}$, which has been blurred out by a known kernel function $H$ and corrupted by random noise. In formulas:

\begin{equation}
\begin{split}
    & R_{ij} = (X \ast H)_{ij} + \eta_{ij}, \\
    & X_{ij} \sim \textrm{Poisson}(V_{ij}), \ \ \eta_{ij} \sim \mathcal{N}(0,\sigma^2_N) \\
    & R_{ij} \xrightarrow{\text{RL algorithm}} V_{ij}
    \label{RL formulas}    
\end{split}
\end{equation}

\noindent where `$\ast$' denotes the convolution symbol, and $X_{ij}$ and $\eta_{ij}$ are uncorrelated random variables for each $i$ and $j$. It is intuitive to see that this problem is analogous to our case, where $X_{ij}$ is the number of photons landing on pixel $(i,j)$, $V_{ij}$ is the intensity of the spectrum at that pixel energy $I(E_{ij})$, $H$ is an effective charge spreading function (defined in the next section), $\eta_{ij}$ is the thermal noise on pixel $(i,j)$ and $R_{ij}$ is the ADU map. Therefore, defining $A_{ij}^{proc}$ as:

\begin{equation}
A'_{ij} \xrightarrow{\text{RL algorithm}} A_{ij}^{proc} \ ,
\label{RL deconvolution}
\end{equation}

\noindent it constitutes an estimate for $I(E_{ij})$, the quantity of interest. A more rigorous justification for the application of this method is presented in the next section.

Fig.~\ref{algorithms}d illustrates the result of applying the RL algorithm to the camera image already treated with the photon counting algorithms.
As observed, the algorithm decomposes the agglomerated clusters into the most likely distribution of photon hits (still in ADU units), leaving the single-photon clusters unchanged. The comparison between the reconstructed map $A_{ij}^{proc}$ and the original photon hits is excellent, with only very few mismatches.

Finally, the camera units are converted to photon numbers, which provides the intensity in the measured spectrum:
\begin{equation}
    N^{\gamma}_{ij} = \frac{A_{ij}^{proc}}{ADU_{sp}} \ , %\frac{1}{E_{ij}} \, 
\label{photon number}
\end{equation}



%where the first step converts the ADU values on each pixel to energy units, and the second step converts the energy units on camera to the equivalent number of photons measured on that pixel.









\subsubsection{\label{Probability distribution of the pixel values}Probability Distribution of the ADU Counts}








Let us focus on regions with high local $\lambda$, where the photon counting algorithm leaves pixel values unchanged. While this discussion can be extended to include thresholding and charge accumulation operations, doing so adds complexity without altering the fundamental concepts. First, consider the signal on a sufficiently fine subgrid $k,l = 1,2, ...$ (see also Section \ref{Synthetic data}):

\begin{equation}
    A_{kl} = \sum_{k'l'} ADU_{sp} (N^{\gamma}_{k'l'} + \beta_{k' l'}) H_{k-k' \ l-l'} \ \ \ \ \ \ ,
    \label{extra1}
\end{equation}

\noindent where $H$ is the charge spreading function, and $N^{\gamma}_{k'l'} \sim {\rm Poisson}(I_{k'l'})$, with $I_{k'l'}$ denoting the spectral intensity at subpixel $(k'l')$. The term $\beta_{k' l'} \sim \mathcal{N}(0,f N^{\gamma}_{k'l'})$ is a compound probability distribution accounting for random fluctuations in e-h pairs production, quantified by the Fano factor $f$. This factor is usually small (e.g. $f=0.12$ in Silicon~\cite{humphries2020isochoric}), and therefore this noise source is often neglected. The value $A_{ij}$ is obtained by summing $A_{kl}$ over the subpixels within each pixel, and then adding the thermal noise specific to the pixel:

\begin{equation}
    A_{ij} = \sum_{k,l \in (i,j)} A_{kl} + \eta_{ij} \ \ ,
    \label{extra2}
\end{equation}

\noindent where $k,l\in(i,j)$ indicates that subpixel $(k,l)$ lies in pixel $(i,j)$. Using the properties of the Poisson distribution and neglecting the e-h pairs production noise, one now writes:

\begin{equation}
\begin{split}
      & A_{ij} = \sum_{k,l \in (i,j)} \sum_{k'l'} ADU_{sp}{\rm Poisson}(I_{k'l'}) H_{k-k' \ l-l'} +  \eta_{ij} \\
      & = ADU_{sp} \ {\rm Poisson}\biggl(\sum_{k,l \in (i,j)} \sum_{k'l'} I_{k'l'} H_{k-k' \ l-l'}\biggr) + \eta_{ij} \\
      & = ADU_{sp} \ {\rm Poisson}\biggl(\sum_{i' j'} \sum_{k'l' \in (i',j')} I_{k'l'} \sum_{k,l \in (i,j)} H_{k-k' \ l-l'}\biggr) + \eta_{ij} \\
      & \approx ADU_{sp} \ {\rm Poisson}\biggl(\sum_{i' j'} \bar{I}_{i'j'} \sum_{k'l' \in (i',j')} \sum_{k,l \in (i,j)} H_{k-k' \ l-l'}\biggr) +  \eta_{ij} \\
      & = ADU_{sp} \ {\rm Poisson}\biggl(\sum_{i' j'} I_{i'j'} \tilde{H}_{i-i' \ j-j'}\biggr) + \eta_{ij} \\
      & = ADU_{sp} \ \sum_{i' j'} {\rm Poisson} (I_{i'j'}) \tilde{H}_{i-i' \ j-j'} + \eta_{ij} \ \ ,
\end{split}
    \label{extra3}
\end{equation}

\noindent where $I_{k'l'}$ was reasonably assumed to be constant within each pixel and equal to $\bar{I}_{i'j'}$. Here, I defined $I_{i'j'} = \sum_{k'l' \in (i',j')} I_{k'l'}\approx N_p \bar{I}_{i' j'}$, with $N_p$ representing the number of subpixels per pixel, and an effective charge spreading function as:

\begin{equation}
    \tilde{H}_{i-i' \ j-j'} = \frac{1}{N_p}\sum_{i' k' \in (i' j')} \sum_{i k \in (i j)} H_{k-k' \ l-l'} \ \ . 
    \label{extra4}
\end{equation}

\noindent The translation invariance of $\tilde{H}$ follows directly from the inherent translation invariance of $H$. Equation (\ref{extra3}) proves that the operation defined by Eq. (\ref{RL deconvolution}), combined with the division in Eq. (\ref{photon number}), yields an estimate for $I_{ij}$.

To account for the e-h pairs production noise as well, one can perform Bayesian inference on the camera image before applying the Richardson-Lucy algorithm. Defining $\Chi_{ij} = ADU_{sp} \sum_{i' j'} {\rm Poisson} (I_{i'j'}) \tilde{H}_{i-i' \ j-j'}$, from Eq. (\ref{extra1}) one has:

\begin{equation}
\begin{split}
    & A_{ij} = \Chi_{ij} + \eta_{ij} + ADU_{sp}\sum_{kl \in (ij)} \sum_{k'l'} \beta_{k'l'} H_{k-k' \ l-l'} \\
    & = \Chi_{ij} + \mathcal{N}(0,\sigma_N^2) + \mathcal{N}(0, f ADU^2_{sp}\sum_{kl \in (ij)} \sum_{k'l'} N^{\gamma}_{k'l'} H^2_{k-k' \ l-l'}) \\
    & = \Chi_{ij} + \mathcal{N}(0, \sigma_N^2+f ADU^2_{sp} \sum_{kl \in (ij)} \sum_{k'l'} N^{\gamma}_{k'l'} H^2_{k-k' \ l-l'}) \\
    & \approx \Chi_{ij} + \mathcal{N}(0, \sigma_N^2+f \Chi_{ij}ADU_{sp}) \ \ , 
\end{split}
\label{extra5}   
\end{equation}

\noindent where the Gaussian distribution properties are used, and the charge spreading is ignored as a zero-order approximation ($H_{0,0} = 1$ and 0 elsewhere). This approximation can be refined by expressing $N_{kl}^{\gamma}$ in terms of $\Chi_{ij}, \Chi_{i+1,j}, \Chi_{i,j+1},...$, using the relations $\Chi_{ij} = ADU_{sp}\sum_{kl \in(i,j)}\sum_{k'l'}N^{\gamma}_{k'l'} H_{k-k' \ l-l'}$ ($H$ is known). While Eq. (\ref{extra5}) could incorporate thresholding and charge accumulation, doing so would significantly complicate the expressions without bringing large benefits. From Eq. (\ref{extra5}), one obtains:

\begin{equation}
    \mathcal{P}(A_{ij}|\Chi_{ij}) = \frac{\exp\biggl(-\frac{1}{2}\frac{(A_{ij}-\Chi_{ij})^2}{\sigma_N^2+ f\Chi_{ij}ADU_{sp}}\biggr)}{\sqrt{(2\pi(\sigma_N^2+ f\Chi_{ij}ADU_{sp})}} \ \ .
\end{equation}

\noindent Next, by Bayesian inference and assuming a flat prior for $\Chi_{ij}$, one estimates $\Chi_{ij}$ from a given ADU map as:

\begin{equation}
    \langle \Chi_{ij} \rangle = \int_{\mathbf{R}} d\Chi_{ij} \Chi_{ij} \mathcal{P}(\Chi_{ij}|A_{ij})  \ \  {\rm with} \ \mathcal{P}(\Chi_{ij}|A_{ij}) = \frac{\mathcal{P}(A_{ij}|\Chi_{ij})\mathcal{P}(\Chi_{ij})}{\int_{\mathbf{R}} d\Chi_{ij} \mathcal{P}(A_{ij}|\Chi_{ij})}  \ \ ,
\end{equation}

\noindent where the integrals are convergent and can be evaluated numerically. This approach allows one also to estimate the uncertainty on $\Chi_{ij}$ as the standard deviation of $\mathcal{P}(\Chi_{ij}|A_{ij})$. The RL method can then be applied to $\langle \Chi_{ij} \rangle$, providing an estimate for $I_{ij}$, along with an error that comes from the uncertainty on $\langle \Chi_{ij} \rangle$ propagated through the RL algorithm.


















\section{\label{Construction of the spectrum and solid angle corrections}Construction of the Spectrum and Solid Angle Corrections}








The final step is to construct an estimate of the photon emission spectrum, $S(E)$, in the direction of the detector, based on the photon hit distribution $N^{\gamma}_{ij}$ obtained in the previous section. In the following, the emission spectrum is assumed to be independent of the unit direction $\boldsymbol{\hat{k}}$ in the solid angle detected by the camera $\Omega_{det}$. If this assumption does not hold, the following discussion remains largely unaffected, except that the calculated spectrum would represent values averaged over all detected directions for a given energy. $S(E)$ describes the probability that, given a source spectrum, a photon of energy $E$ is scattered in the direction observed by the camera. Its units are sr$^{-1}$eV$^{-1}$. The measured spectrum $I(E_l)$, which represents the probability of detecting a photon within a given energy interval, can then be expressed as follows:

\begin{equation}
\begin{split}
& I(E_l) = \int_{\Omega_{det}}d\Omega(\boldsymbol{\hat{k}}) \int_{E_l-\Delta E/2}^{E_l+\Delta E/2} dE \ S'(E,\boldsymbol{\hat{k}}) \\
& l=1,\dots,N_b \ ,  
\end{split}
\label{measured spectrum} 
\end{equation}

\noindent where $E$ and $\boldsymbol{\hat{k}}$ are the energy and unit wave vector of the photons, $\Delta E$ is the energy resolution, and $N_b$ is the number of energy bins. Here, $S'(E,\boldsymbol{\hat{k}})$ is the emission spectrum $S(E)$ modified by the action of the Bragg crystal, summarised in the shape function $F(E,\boldsymbol{\hat{k}})$:

\begin{equation}
S'(E, \boldsymbol{\hat{k}}) = S(E)F(E,\boldsymbol{\hat{k}}) \ .
\label{emissionspectrum}
\end{equation}

\noindent Integrating Eq. (\ref{measured spectrum}) over the energies, one obtains:

\begin{equation}
I(E_l) = \int_{\Omega_{det}}d\Omega(\boldsymbol{\hat{k}}) \ \bar{S}'(\boldsymbol{\hat{k}}; E_l) \ .
\label{passage 1}
\end{equation}

\noindent Given the energy-dispersive effect of the crystal on the radiation, and assuming that crystal broadening (including depth and mosaicity effects) is sufficiently small, $\bar{S}'(\boldsymbol{\hat{k}}; E_l)$  (l=1,...,$N_b$) can be approximated in the following way:

\begin{equation}
\begin{split}
& \bar{S}'(\boldsymbol{\hat{k}}; E_l) = \bar{S}(E_l)=\int_{E_l-\Delta E/2}^{E_l+\Delta E/2} dE \ S(E) \\ 
& \forall k \ | \ E_l = \argmin_{E_{ij}}|E(\boldsymbol{\hat{k}})-E_{ij}| \ ,
\end{split}
\label{passage 2}
\end{equation}

\noindent and equal to zero for all the other values of $\boldsymbol{\hat{k}}$. Here, $E(\boldsymbol{\hat{k}})$ is given by Eq. (\ref{energy map equation}), with:

\begin{equation}
\boldsymbol{\hat{k}} = \frac{\boldsymbol{r_0}+\boldsymbol{Rr'}}{\|\boldsymbol{r_0}+\boldsymbol{Rr'}\|} \ , \ \boldsymbol{r'} = (x',y',0) \ .
\label{passage 3}
\end{equation}

\noindent Now, the integration over the solid angle in Eq. (\ref{passage 1}) becomes:

\begin{equation}
I(E_l) = S(E_l)\cdot \sum_{ij|E_{ij}\in \Delta (E_l)} \Omega_{ij} \ \ . 
\label{passage 4}
\end{equation}

\noindent Having denoted the interval $(E_l - \Delta E/2, E_l - \Delta E/2)$ as $\Delta (E_l)$, and $\Omega_{ij}$ the solid angle associated with pixel $ij$. This solid angle  for the flat crystal spectrometer is then computed through the following formulas:




\begin{align}
& (\hat{z}')_{\rm lab} = \boldsymbol{R}(\hat{z}')_{\rm cam} \ \ \ \ \ \ {\rm with} \ \ (\hat{z}')_{\rm cam} = (0,0,1)  \\
& \boldsymbol{x}_0 = \hat{z}'(\hat{z}'\cdot \boldsymbol{r}_0) \ \ \ \ \ \ \ \ \ \ {\rm with} \ \ (\boldsymbol{r}_0)_{\rm lab} = D(\cos\bar{\theta}_B,0,\sin\bar{\theta}_B) \\
& |\boldsymbol{x}_0| = |\hat{z}'\cdot \boldsymbol{r}_0| = D' \ \ \ \longrightarrow  \ \ (\boldsymbol{x}_0)_{\rm cam} = (0,0,D') \\
& \boldsymbol{x}'' = \boldsymbol{r}_0 - \boldsymbol{x}_0 + \boldsymbol{x}' \ \ \ \ \ \longrightarrow \ \ (\boldsymbol{x}'')_{\rm cam} = (x'',y'', 0) \\
& {\rm with} \ \ (\boldsymbol{r}_0)_{\rm cam} = \boldsymbol{R}^{-1} (\boldsymbol{r}_0)_{\rm lab}  \ \ {\rm and} \ \ (\boldsymbol{x}')_{\rm cam} = (x',y', 0) \\ 
& \boldsymbol{x} = \boldsymbol{x}_0 + \boldsymbol{x}''  \ \ \ \ \ \ \ \ \ \ \ \ \longrightarrow \ \ (\boldsymbol{x})_{\rm cam}=(x'',y'',D') \\
& \Omega(x,y) = \frac{dx dy \cos(\beta)}{|\boldsymbol{x}|^2} = \frac{dx dy |\boldsymbol{x}_0|}{|\boldsymbol{x}|^3} = \frac{D'}{(x''^{2}+y''^{2}+ D'^{2})^{3/2}}
\label{passage 5}
\end{align}

\noindent where the relative quantities are defined in Fig. \ref{solid angles}. Similar reasoning can be carried out for computing the solid angles in the von H\'amos geometry.

\begin{figure}[t]
    \centering
    \includegraphics[width=0.8\textwidth,keepaspectratio]{Figures/solid angles.pdf}
    \caption[Scheme for the solid angles calculations in the flat crystal spectrometer case]{Scheme for the solid angles calculations in the flat crystal spectrometer case.
    }
    \label{solid angles}
\end{figure}
On the other hand, $I(E_l)$ can be estimated (see also Section \ref{error calculations}) from the distribution of photon hits $N^{\gamma}_{ij}$ as:

\begin{equation}
\tilde{I}(E_l) = \sum_{ij|E_{ij}\in \Delta E_l} N^{\gamma}_{ij} \ \ .
\label{passage 6}
\end{equation}

\noindent Putting Eq.~(\ref{passage 4}) and~(\ref{passage 6}) together, one finally obtains an estimate $\tilde{S}(E_l)$ for the emission spectrum:

\begin{equation}
\tilde{S}(E_l) = \frac{\sum_{ij|E_{ij}\in \Delta E_l} N^{\gamma}_{ij}}{\sum_{ij|E_{ij}\in \Delta (E_l)} \Omega_{ij}} \ \ .
\label{passage 7}
\end{equation}

\noindent When summing the number of photons over the energy contours associated to these bins, the fluctuations due to the thermal noise and Poissonian statistics are reduced, thanks to the large numbers of pixels.














\section{\label{sec3} Application to Synthetic and Experimental Data}






In this section, the discussed methods are applied to both synthetic and experimental data. The synthetic data are first utilised to evaluate the algorithms' performance in enhancing resolution and SNR for different conditions. The experimental data then demonstrate how these methods can amplify specific spectral features that would otherwise be weak, showcasing the efficacy of these techniques in real situations. 
 




\subsection{\label{Synthetic data}Synthetic Data}




In order to test these techniques, some artificial camera images were generated through the following steps (illustrated in Fig. \ref{synthetic images}):

\begin{enumerate}
    \item A 40x40 array is created to represent the camera pixels and a finer 4000x4000 subgrid to depict its surface. For simplicity, the energy contours on the camera are vertical, corresponding to the columns of the pixel matrix.
    \item A number of photons, specified by the selected $\lambda$, are then sent onto the subgrid. Their energy, which determines the pixel column they land on, is sampled from the spectrum $S(E_l)$ shown in Fig. \ref{synthetic images}a and then converted to an ADU value, taking ADU$_{sp}=120$.
    \item Charge spreading is applied by convolving the subgrid image with a gaussian kernel of width $R_s$. The ADU values of the subgrid belonging to the same pixel are then summed together to obtain the actual camera image (see Fig \ref{synthetic images}c)
    \item Finally, gaussian random noise with standard deviation $\sigma_N$ is added to the image (Fig. \ref{synthetic images}d).
\end{enumerate}

\begin{figure*}
    \centering
    \includegraphics[width=0.8\textwidth,keepaspectratio]{Figures/synthetic images.pdf}
    \caption[Illustration of the steps for the generation of synthetic images]{(a) The normalised spectrum used for producing the synthetic camera images. (b) The charge deposited on the camera surface, with the pixel grid depicted. (c) The ADU map on the pixel grid. (d) The final camera image with the noise applied.
    }
    \label{synthetic images}
\end{figure*}

\noindent The performances of these techniques were evaluated for different values of $\lambda$, $R_s$ and $\sigma_N$ by computing the $\mathcal{L}_2$ distance between $S(E_l)$ and the reconstructed spectrum $\tilde{S}(E_l)$. The latter is obtained by integrating the processed camera images along its columns and then normalising it for comparison with $S(E_l)$. Fig. \ref{performances} compares, for $\lambda = 0.1$ and $\lambda = 1$, the original spectrum with the spectra reconstructed using different methods to treat the camera image: only simple thresholding, only the photon counting algorithm, only the deconvolution algorithm, and the sequential application of photon counting and deconvolution techniques (`Hybrid' mode), as discussed in the previous sections. Panel (a) shows that, for small $\lambda$ values, the photon counting algorithm alone outperforms the deconvolution algorithm alone and the simple thresholding. This is because, in regions with very few photons—where ADU counts are mostly produced by noise—the photon counting algorithm achieves higher accuracy, outweighing the contributions from regions with a stronger signal. In contrast, for large $\lambda$ values (panel (b)), the photon counting algorithm and the simple thresholding offer little improvement over the raw image, whereas the deconvolution algorithm is able to treat better regions with high $\lambda$, enhancing and better resolving the closely spaced spectral features (as seen in the two spectral lines at the lower end of the energy range). The sequential application of these two techniques (`Hybrid' mode) combines their benefits, yielding the best performances for both values of $\lambda$. It is noted that the applied method can have a substantial effect on the extracted spectral features, such as spectral linewidths or line ratios, which are often used as key diagnostics in experiments. For the implementation of the RL deconvolution, the Python’s Scikit-image library was used~\cite{scikit}. The number of iterations used in the RL algorithm was selected from a range of values by determining which provided the best match between the processed and original spectrum. The optimal value, found to be 100, was also applied to the analysis of experimental data, considering comparable levels of noise and charge spreading to the synthetic images.

\begin{figure*}[t]
    \centering
    \includegraphics[width=\textwidth,keepaspectratio]{Figures/performances.pdf}
    \caption[Comparison of performances for different image processing methods]{The original spectrum compared, for $\lambda = 0.1$ (a) and $\lambda = 1$ (b), to the spectra reconstructed with different methods: only the photon counting algorithm (`Photon counting'), only the deconvolution algorithm (`Richardson-Lucy'), both the techniques combined (`Hybrid'), none of them (`Simple Thresholding'). The camera parameters are $R_s=0.4$, $\sigma_N=30$.
    }
    \label{performances}
\end{figure*}

To evaluate performance, the quality of the spectrum reconstructed from the `Hybrid'-processed image was studied by comparing it to that obtained through `Simple thresholding' of the raw image, as a function of the filling fraction. The effect of changing the noise level is included in this study, as the ratio $\lambda/\sigma_N$ is the dominant parameter determining the quality of the reconstruction. This analysis was conducted for two values of the charge spreading, $R_s=0.3$ and $R_s=0.6$ in pixel units, and its results are reported in Fig. \ref{performances3}.

As anticipated, the `Hybrid'-processed image consistently yields a higher quality spectrum than the thresholded image across all $\lambda$ values and for each $R_s$. The smallest difference occurs at intermediate $\lambda$ values, where neither the photon counting algorithm nor the deconvolution algorithm are particularly effective. Additionally, the quality of the reconstruction improves for smaller $R_s$ and larger $\lambda$, reaching an asymptote as $\lambda\to\infty$. This behaviour is expected, as a larger number of photons on the camera reduces the Poissonian error and the relative strength of noise, while a smaller $R_s$ facilitates the deconvolution process and the identification of the thresholding value, better separating the 0th and 1st peak in the histogram of photon counts. 

\begin{figure}[t]
    \centering
    \includegraphics[width=0.8\textwidth,keepaspectratio]{Figures/performances3_v2.pdf}
    \caption[Performances of the proposed image processing method for various fill fractions]{$\mathcal{L}_2$ distance between the original spectrum and the one reconstructed from the processed and thresholded image as a function of the filling fraction, for two values of $R_s$. Note that the asymptote value is not zero, as discrepancies from charge spreading and finite pixel size persist in the limit $\lambda \to \infty$.
    }
    \label{performances3}
\end{figure}











\subsection{\label{Real data}Experimental Data}







This section shows the application of this analysis approach to a real camera image and the resulting improvements in the spectrum obtained. Figure \ref{real_data} reports the MgF$_2$ emission spectrum computed either with simple thresholding or with the image processing described in Section \ref{Estimate of the photon hits distribution}. The parameters used in these operations ($ADU_{sp}$, $\sigma_N$, $\lambda$) were determined by fitting a synthetic histogram to the experimental one (see Fig. \ref{histogram}e). To generate the synthetic histogram, I used an estimated spectrum computed by integrating the camera signal over the energy contours defined by the energy map, and the parameters of the PIXIS-XF camera (the detector employed for this experiment), which were provided in the relevant documentation~\cite{yaffe1997x}. Using this documentation, the fitting process obtained $R_s = 1.7$ in pixel unit, $ADU_{sp} = 43$ ADU and $\sigma_N=4$ ADU, which were applied in the processing of the experimental images. For the thresholding, the values specified in Section \ref{Photon counting algorithms} were used. The error bars for the two spectra were computed using Poissonian inference, and neglecting noise-related errors, drastically reduced by thresholding and integration over the energy contours. Details of these error calculations are provided in the next section.

\begin{figure*}[t]
\centering
    \includegraphics[width=0.8\textwidth,keepaspectratio]{Figures/Performances_real.pdf}
    \caption[Comparison of performances for different image processing methods on experimental data]{Emission spectrum of MgF$_2$ obtained with simple thresholding (red), and processing the images with the photon counting and deconvolution methods (blue). Error bounds are indicated by shaded regions. The peaks in the spectrum correspond to the $K\alpha$ (from 1250 to 1370 eV) and $K\beta$ (from 1430 to 1500 eV) emission lines of the different ionisation states of the Mg ions.}
    \label{real_data}
\end{figure*}
The enhancement of features is particularly noticeable for the narrowest and most intense ones, such as the low-energy peaks around 1260 eV. I believe that such improvements are mainly due to the reduction of the background noise, and they could therefore be essential when the physical information of the probed system is contained in weaker features (e.g. Raman scattering peaks). Moreover, the enhanced resolution achieved through the deconvolution algorithm for the narrowest and most intense peaks has significant implications for HED science, where line ratios and widths play a critical role. In this example, line ratios are essential for determining relative ion populations and, consequently, plasma temperature~\cite{ciricosta2016detailed}. Temperature is also commonly inferred using x-ray Thomson scattering via detailed balance, where precise line ratio measurements are equally crucial~\cite{doppner2009temperature}. On the other hand, accurate estimation of line widths can be important for determining internal electromagnetic field intensities and collision dynamic through Zeeman, Stark and collisional broadening~\cite{ferri2022stark, doron2014determination, gornushkin1999line,PhysRevE.109.045204}. For the photon counting algorithm to be effective, it must be applied to single-shot images rather than averaged ones, where single photon clusters are reduced by summation. This requirement, combined with the typical size of a camera image (here 2048x2048 pixels), makes the process time-consuming for real-time analysis. However, this procedure can be performed during post-analysis, where time constraints are less demanding. The current implementation takes approximately 8 minutes to process a single 2048×2048 image, which remains quite slow even for post-analysis. Improving the efficiency of these algorithms would greatly enhance their practical usability. Additionally, the computation can be parallelised by distributing different regions of the image across multiple processors, significantly reducing processing time.








\subsubsection{\label{error calculations} Error calculations}





To compute the uncertainties on the experimental spectra, one begins by calculating the best estimate of the parameter $I$ (also the mean) of a Poisson distribution, given an observation $n$:

\begin{equation}
N \sim Poisson(I) \  \Rightarrow \ \mathcal{P}(N = n|I) = \frac{I^n e^{-I}}{n!}.
\label{appendix 1}
\end{equation}

\noindent Assuming a flat prior distribution for $I$ over a large interval $[0,M]$ (since $I$ cannot be infinitely large), and noting that $\mathcal{P}(n)=\int_0^M \mathcal{P}(n|I)  \mathcal{P}(I) = \int_0^M dI \frac{I^n e^{-I}}{n!} \frac{1}{M} \approx \frac{1}{M}$ if $M$ is sufficiently large, from Bayes' theorem one finds:

\begin{equation}
\mathcal{P}(I|n) = \frac{\mathcal{P}(n|I)\mathcal{P}(I)}{\mathcal{P}(n)} = \mathcal{P}(n|I) \ \ .
\label{appendix 2}
\end{equation}

\noindent From this equation, one can calculate the best estimate for $I$ as:

\begin{equation}
\langle I \rangle (n) = \int_0^M dI I \mathcal{P}(I|n) \int_0^M dI I \frac{I^n e^{-I}}{n!} \approx n+1 \ \ .
\label{appendix 3}
\end{equation}

\noindent Similarly, the error on this estimate is determined as the standard deviation of $\mathcal{P}(I|n)$:

\begin{equation}
\sqrt{\langle I^2 - \langle I \rangle^2 \rangle} (n) \approx \sqrt{n+1} \ \ .
\label{appendix 4}
\end{equation}

\noindent To analyse our case, $\tilde{I}(E_l)= \sum_{ij|E_{ij}\in \Delta E_l} N^{\gamma}_{ij}$ is considered as a realisation of the random variable $I(E_l) \sim {\rm Poisson}(\Omega(E_l) S(E_l))$. Here, the contribution to $I(E_l)$ from noise, which is generally small due to the thresholding operations and the summation over many pixels, is disregarded. Additionally, the error propagation through processing steps is not considered, taking instead directly $N^{\gamma}_{ij}$ as the estimate of the photon hits map. Under these assumptions, the spectrum and its error, for each energy bin, are given by:
\begin{equation}
\langle S(E_l) \rangle = \frac{\tilde{I}(E_l)+1}{\Omega(E_l)} \ \  \sigma(S(E_l)) = \frac{\sqrt{\tilde{I}(E_l)+1}}{\Omega(E_l)} \ \ .
\label{appendix 5}
\end{equation}

\noindent As a final consideration, it is observed that the error normalised to the signal, $\ensuremath{\frac{\sigma(S(E_l))}{\langle S(E_l) \rangle} = \frac{1}{\sqrt{\tilde{I}(E_l) + 1}}}$, correctly decreases as the number of photons per bin increases.









\section{\label{sec4} Summary}






This chapter provided a detailed guide for constructing spectra from the relative camera images, for two of the most common geometries employed in HED science spectroscopy. Specifically, I examined the computation of the energy map and the solid angle correction for the pixel array, the estimation of important camera parameters such as the thermal noise intensity, and presented techniques to estimate the distribution of photon hits from the raw camera data, with the aim of increasing the SNR. These image processing techniques are based on the sequential combination of clustering algorithms and deconvolution methods for simultaneously treating regions with different local fill fractions. Application of these methods to both synthetic and experimental data demonstrated notable improvements in SNR and feature clarity, proving their potential in contexts where the measured physical features are weak and the precise determination of line widths and intensities is crucial. Future work could explore introducing in this framework the capacity of multiple photons-cluster identification, or expanding the variety of allowed cluster shapes. Additionally, it would be useful to integrate this work with studies on the instrument function to combine the improvement of SNR with resolution enhancement. 
The experimental data presented in this thesis were collected using detectors and setups where pixel broadening is negligible compared to IF broadening. As a result, deconvolution techniques were not applied to process the experimental images. Photon counting algorithms were instead applied to the camera images from the flat crystal spectrometer. However, for the von H\'amos data, photon counting was not necessary due to the high fill fraction in the region of interest. The analysis of the von H\'amos data, which constitutes the primary focus of this work, relies solely on simple thresholding.




% Maybe the last sentence is for the conclusion chapter









