Monte Carlo Simulation of γ − γ Correlation Functions

Abstract

$\gamma -\gamma$ correlation functions are mathematical expressions that describe the angular distribution of cascade $\gamma$-rays emitted from an atomic nucleus. Cascade transitions may occur in either a two-step deexcitation or through an excitation-deexcitation process of a particular energy level inside the nucleus. In both cases, the nucleus returns to its ground energy state. Spin and parity of the excited state can be determined experimentally using the asymmetry of the angular distribution of the emitted radiation. $\gamma -\gamma$ correlation functions are only valid for point-like targets and detectors. In the real experiments, however, neither the target nor the detector is point-like. Thus, misassignment of the spin-parity of energy levels may easily take place if only the analytical equations are considered. Here, we develop a new Monte Carlo simulation method of the $\gamma -\gamma$ correlation functions to account for the extended target and detector involved in spin-parity measurements using nuclear resonance fluorescence of nuclei. The proposed simulation tool can handle arbitrary geometries and spin sequences. Additionally, we provide numerical calculations of a parametric study on the influence of the detection geometry on the angular distribution of the emitted $\gamma$-rays. Finally, we benchmark our simulation by comparing the simulation-estimated asymmetry ratios with those measured experimentally. The present simulation can be employed as a kernel of an implementation that simulates the nuclear resonance fluorescence process.

1. Introduction

Resonance fluorescence of an atom or a nucleus is one of the richest sources of our information about the structure of these too-small-to-see systems. Many aspects of quantum theory have been proven by X-ray and $\gamma$-ray spectroscopies. The availability and developments of photon sources with polarization capabilities have provided accurate means to measure the spin and parity of the nuclear energy levels in the laboratory. As an example, the polarized $\gamma$-rays generated by laser Compton scattering (LCS) are being used in nuclear structure studies that focus on the determination of spin and parity [1,2]. The technique of using polarized $\gamma$-ray beams offers an unambiguous estimation of the spin and parity, while the main detection system is simply $\gamma$-ray spectroscopy. Furthermore, many applications rely on the resonance fluorescence, especially in nuclear safeguards and nuclear security [3,4,5].

Owing to the significance of the nuclear resonance fluorescence (NRF), many simulation works were conducted [6,7,8]. For example, Jordan and Warren [7] and Hayakawa et al. [8] developed a Monte Carlo simulation for the NRF process. More recently, Vavrek et al. [9,10] have upgraded the simulation model of Jordan and Warren to enhance its accuracy. They have tested the simulation model, which lowered the discrepancy between the simulation and experiment to approximately $15\%$. However, all of these attempts considered the scattering of unpolarized photons limited by the specific resonances of some nuclei. Other simulation works were performed by Omer and Hajima [11,12] and Turturica et al. [13] on the elastic scattering as a beam background. Nonetheless, the literature is missing a Monte Carlo simulation for $\gamma -\gamma$ correlation functions, which are crucial for the nuclear decays.

Monte Carlo simulation is an indispensable tool to represent the extended volumes of the targets and detectors correctly. Here, we develop a Monte Carlo calculation framework of $\gamma -\gamma$ correlation functions as the basic constituent of the formalism of NRF interactions. The simulation is then exploited to conduct a parametric study to investigate the effect of the detection geometry on the spin-parity assignment of nuclear levels. We demonstrate that at certain geometries, the anisotropic parameters may be changed by a factor of ≥100%, which certainly can lead to misassignment of the spin-parity of the levels. Although our simulation focuses on $\gamma -\gamma$ correlation functions, it can be extended to include other correlation functions, such as $\alpha -\gamma$ or $\beta -\gamma$ correlations.

In the following section, we briefly introduce the formalism of the angular correlation involved in the NRF interaction. Then, we show in detail how we implement a Monte Carlo simulation-based of this formalism. We discuss some examples of the simulation in Section 3. These examples are carefully selected to match real NRF experiments. Finally, we summarize the results of the simulation and provide an outlook for applying the simulation in future research works.

2. Method

2.1. Angular Correlation Formalism

The spins of the ground and excited states and the multipolarity involved in the transition are fundamental quantities constructing the angular correlation function of the emitted photons. The level diagram, shown in Figure 1, is a general case of an excitation involved in an NRF measurement. The initial state $J_i^{\pi }$ represents the ground state of the nucleus, while $J_f^{\pi }$ is the final state stemming from excitation of a multipolarity $L_i$. After populating the final state, it deexcites either to the ground state again or to a low-lying energy level. In the latter case, the nucleus may deexcite to an intermediate level with a spin of $J_m^{\pi }$. In such a situation, the NRF excitation is said to have a probability to deexcite to a low-lying energy level, likely the first excited state. The multipolarity of the deexcitation is $L_f$.

When the nucleus is excited with a well polarization-defined $\gamma$-ray photon, the angular distribution of the re-emitted photons is given by [14]

$$W(\theta ,\varphi )=\sum _{k=0,2}{A}_k\left({\gamma }_i\right)A_k\left({\gamma }_f\right)P_k(cos\theta )+{(\pm )}_{L}cos2\varphi \sum _{k=0,2}{B}_k\left({\gamma }_i\right)A_k\left({\gamma }_f\right)P_k^2(cos\theta ),$$

(1)

where $\theta$ is the angle between the direction of the first radiation ${\gamma }_i$ and the second radiation ${\gamma }_f$, denoted here by the scattering angle. In the NRF process, ${\gamma }_i$ and ${\gamma }_f$ are simply the incident and re-emitted (or scattered) photons, respectively.

$P_k(cos\theta )$ and $P_k^2(cos\theta )$ are the ordinary and associated Legendre polynomials of the order k respectively. The two summations of Equation (1) over k are over all even numbers for which the transition coefficients, $A_k\left(\gamma \right)$, and linear polarization coefficient, $B_k\left(\gamma \right)$, do not vanish. The factor ${(\pm )}_L$ is the parity coefficient, which is $-1$ for the magnetic dipole transition $M1$, and $+1$ for the electric dipole transition $E1$. $\varphi$ is the angle between the polarization of the incident photon and the scattering plane, which is formed by the momentum direction of the incident and the scattered photons. It should be noted that the first term of Equation (1) is the angular distribution function, $W\left(\theta \right)$, in case that the incident photon is unpolarized. The excitation transition coefficient $A_k\left({\gamma }_i\right)$ is given by

$$A_k\left({\gamma }_i\right)=\frac{F_k\left(L_i{L}_i{j}_i{j}_f\right)+2\delta F_k\left(L_i{L}_i^{\prime }{j}_i{j}_f\right)+{\delta }^2{F}_k\left(L_i^{\prime }{L}^{\prime }{j}_i{j}_f\right)}{1+{\delta }^2}$$

(2)

where $\delta$ is the mixing ratio of the two multipolarities $L_i$ and $L_i^{\prime }=L_i+1$. For a pure dipole or quardupole transition, $\delta =0$. This leads to $A_k\left({\gamma }_i\right)=F_k\left(L_i{L}_i{j}_i{j}_f\right)$ where $F_k\left(L_i{L}_i{j}_i{j}_f\right)$ are the F-coefficients, which may be given as

$$\begin{array}{c}\hfill F_k\left(L_i{L}_i{j}_i{j}_f\right)={\left(-1\right)}^{j_f+j_i-1}\left(2L_i+1\right){\left[(2j_f+1)\left(2k+1\right)\right]}^{1/2}\underset{3\mathrm{j}-\mathrm{symbol}}{\underbrace{\left(\begin{array}{ccc}{L}_i& L_i& k\\ 1& -1& 0\end{array}\right)}}\underset{6\mathrm{j}-\mathrm{symbol}}{\underbrace{\left\{\begin{array}{ccc}{j}_f& j_f& k\\ L_i& L_i& j_f\end{array}\right\}}}\end{array}$$

(3)

General expressions for the 3j-and 6j-symbol can be found in references [15,16], and are given by Equations (4) and (5), respectively.

$$\begin{array}{c}\left(\begin{array}{ccc}{j}_1& j_2& j_3\\ m_1& m_2& m_3\end{array}\right)={\left(-1\right)}^{j_1-j_2-m_3}\times \hfill \\ {\left[\frac{\left(j_1+j_2-j_3\right)!\left(j_1-j_2+j_3\right)!\left(-j_1+j_2+j_3\right)!\left(j_1+m_1\right)!\left(j_1-m_1\right)!\left(j_2+m_2\right)!\left(j_2-m_2\right)!\left(j_3+m_3\right)!\left(j_3-m_3\right)!}{\left(j_1+j_2+j_3\right)!}\right]}^{1/2}\hfill \\ \times \sum _k\left[\frac{{\left(-1\right)}^k}{k!\left(j_1+j_2-j_3-k\right)!\left(j_1-m_1-k\right)!\left(j_2+m_2-k\right)!\left(j_3-j_2+m_1+k\right)!\left(j_3-j_1-m_2+k\right)!}\right].\hfill \end{array}$$

(4)

$$\begin{array}{c}\left\{\begin{array}{ccc}{j}_1& j_2& j_3\\ m_1& m_2& m_3\end{array}\right\}={\left(-1\right)}^{j_1+j_2+m_1+m_2}\Delta \left(j_1{j}_2{j}_3\right)\Delta \left(m_1{m}_2{j}_3\right)\Delta \left(m_1{j}_2{m}_3\right)\Delta \left(j_1{m}_2{m}_3\right)\times \hfill \\ \sum _k\left[\frac{{\left(-1\right)}^k\left(j_1+j_2+m_1+m_2+1-k\right)}{k!\left(j_1+j_2-j_3-k\right)!\left(m_1+m_2-j_3-k\right)!\left(m_1+j_2-m_3-k\right)!\left(-j_1-m_1+j_3+k\right)!\left(-j_2-m_2+j_3+k\right)!}\right],\hfill \end{array}$$

(5)

where

$$\Delta \left(xyz\right)={\left[\frac{\left(x+y-z\right)!\left(x-y+z\right)!\left(-x+y+z\right)!}{\left(x+y+z+1\right)!}\right]}^{1/2}$$

Similarly, one can easily write the deexcitation transition coefficient as $A_k\left({\gamma }_f\right)=F_k\left(L_f{L}_f{j}_f{j}_i\right)$. Finally, to determine the polarization coefficient, it is convenient to use the expression $B_k\left({\gamma }_i\right)=\kappa \left(L_i{L}_i\right)F_k\left(L_i{L}_i{j}_i{j}_f\right)$, where the coefficient $\kappa \left(L_i{L}_i\right)$ can be found in references [14,17].

2.2. Anisotropy

Anisotropy of the angular distribution of the emitted $\gamma$-rays is the key factor to identify the spin of the excited states. The anisotropy is estimated by the experimental observables, which include the asymmetry ratio $\mathcal{A}$ and the angular distribution ratio $\mathcal{R}$. The asymmetry requires two simultaneous measurements at a particular scattering angle. The first measurement is parallel to the polarization plane of the incident photon ($\varphi =0$), while the other one is orthogonal to the polarization plane ($\varphi =\pi /2$). The asymmetry ratio is simply the ratio between the difference and the summation of the two measurements,

$$\mathcal{A}\left(\theta \right)=\frac{W\left(\theta ,0\right)-W\left(\theta ,\pi /2\right)}{W\left(\theta ,0\right)+W\left(\theta ,\pi /2\right)}$$

(6)

Unlike the asymmetry ratio, the angular distribution ratio requires four simultaneous measurements with two of them at a scattering angle of ${\theta }_1$ parallel and orthogonal to the polarization plane, while the other two follow the same manner at a different scattering angle, ${\theta }_2$. $\mathcal{R}\left({\theta }_1,{\theta }_2\right)$ is then defined by ratio of the summations over the parallel direction to the summation of the orthogonal direction,

$$\mathcal{R}\left({\theta }_1,{\theta }_2\right)=\frac{W\left({\theta }_1,0\right)+W\left({\theta }_1,\pi /2\right)}{W\left({\theta }_2,0\right)+W\left({\theta }_2,\pi /2\right)}$$

(7)

Interestingly, the symmetry properties of Legendre polynomials dictate that $\mathcal{R}\left({\theta }_1,{\theta }_2\right)$ is independent on the azimuth angle $\varphi$. Consequently, it can be used either in the case that the incident beam is linearly polarized [18] or in the case that the incident beam in unpolarized [19]. In the latter case, estimating $\mathcal{R}\left({\theta }_1,{\theta }_2\right)$ requires only two simultaneous measurements at two different angles within the scattering plane.

We calculated the angular correlation functions (see Appendix A for the final formulas) using the method illustrated in the previous subsection. For the even-even nuclei, the spin of the ground state is $j_i=0^{+}$. The allowed transitions involved in the NRF interaction are limited to dipole and quadrupole. The probability of magnetic quadrupole transition, $M2$, is too small to be considered. Therefore, the spin of the excited state is $j_f=1^{\pm }$ for dipole or $j_f=2^{+}$ for quadrupole transition. Ideal asymmetries and radiation patterns of $M1$ and $E2$ transitions are shown in Figure 2. The $E1$ transition is identical to the $M1$ transition but rotated by 90. Therefore, $M1$ and $E1$ transitions can be easily distinguished by $\mathcal{A}\left(\theta =\pi /2\right)$ which is $+1$ for $M1$ and $-1$ for $E1$. However, $\mathcal{A}\left(\theta =\pi /2\right)$ cannot be used to distinguish $M1$ from $E2$ because $\mathcal{A}\left(\theta =\pi /2\right)=+1$ for both. In this case, $\mathcal{A}\left(\theta =3\pi /4\right)=-1\left(+\frac{1}{3}\right)$ for quadrupole (dipole) can be exploited. Similarly, $\mathcal{R}\left(\theta \right)$ can distinguish dipoles from quadrupole transitions, but this is not the case for $M1$ and $E1$.

The situation is more complicated with odd nuclei, which are characterized by half-integer ground state spin. Figure 3 demonstrates $\mathcal{A}\left(\theta \right)$ and patterns of four possible transitions of a nucleus with a ground state spin of $\frac{5}{2}$. The allowed spins of the excited state include $\frac{3}{2}$, $\frac{5}{2}$, or $\frac{7}{2}$ for dipole and $\frac{9}{2}$ for quadrupole. Except for the spin sequence $\frac{5}{2}\to \frac{3}{2}\to \frac{5}{2}$ with almost isotropic angular distribution, it is easy to distinguish the parity of the dipole transition. In contrast, the sensitivity of $\mathcal{A}\left(\theta \right)$ and $\mathcal{R}\left(\theta \right)$ to differentiate the dipole from the quadrupole transitions of the excited state is lower than the even-even nuclei. For example, $\mathcal{A}\left(\theta =\pi /2\right)$ for the spin sequence of $\frac{5}{2}\to \frac{5}{2}\to \frac{5}{2}$ is $0.302$, while it is $0.349$ for the spin sequence of $\frac{5}{2}\to \frac{9}{2}\to \frac{5}{2}$. Furthermore, $\mathcal{A}\left(\theta =3\pi /4\right)$ for the two spin sequences $\frac{5}{2}\to \frac{7}{2}\to \frac{5}{2}$ and $\frac{5}{2}\to \frac{9}{2}\to \frac{5}{2}$ is equal. Nonetheless, the closer the angular distributions of the spin sequences, the more probability to misassign the spin and parity of a certain level and the greater the need to accurately account for the geometrical conditions by Monte Carlo simulation.

2.3. Monte Carlo Simulation

Geant4 [20] is a powerful Monte Carlo simulation toolkit dedicated to simulating particle transport in matter. In Geant4, all particle interactions with matter are tracked. Parameters of the particle, such as position, momentum, and energy, can be retrieved at almost every step from the generation of the particle to the point where the particle is made to dump. Our procedure for the simulation model of $\gamma -\gamma$ correlation functions consisted of three major steps. The first step was to prepare the prerequisites for the simulation. Then, we modeled the event generator. Finally, we modeled the detection geometry and computed the asymmetries as a function of the scattering geometry and spin sequence.

2.3.1. Prerequisites

In the design of the simulation model, we considered that the geometrical parameters of the NRF process were known prior to the simulation. The parameters include the following:

• Spin sequence of the expected NRF transition, which means the spin of the ground and excited states.
• Dimensions of the target and the detectors.
• Detection geometry, including target and detector sizes and their positions with respect to each other.

We introduced a new class within the frame of the Geant4 toolkit; namely, FCoefficients class. The class returns F-coefficients based on the spin sequence, according to Equation (3). It should be emphasized that we tested the FCoefficients class for the F-coefficients with all constituents, e.g., 3j- and 6j-symbols, and compared them with the tabulated data of references [15,16]. We found no differences between our coded calculations and the reference values.

2.3.2. Event Generator

We developed a new class, which is inherited from the Geant4 event generation class; namely, G4PrimaryGeneratorAction. In this class, we introduced two functions aimed at generating the position and momentum of photons emitted from the target. One function was dedicated to generating the dipole transition patterns depending on the desired spin sequence, while the other function was responsible for the generation of quadrupole patterns. In both functions, we employed the rejection sampling method to sample the momentum of the emitted photons randomly.

The shape and size of the target were controlled by the G4PrimaryGeneratorAction class. The positions from which the emitted photons started were homogeneously sampled over the whole size of the target. In the present study, we considered a cylindrical target of radius $R_t$ and length ${\lambda }_t$ whose axis coincides with the z-axis. However, the target shape may be changed to match any specific experimental setup.

We adapted the G4PrimaryGeneratorAction class such that the emitted particles were of the geantino virtual photons. Hence, there were no interactions considered through the virtual photons’ transport. The virtual photons, geantinos, are identical to the mathematical vectors. The particles were generated at the beginning of the simulation and continued to move in straight lines until they were counted by a detector located at the desired position.

2.3.3. Detection Geometry

In our simulation, we employed the G4DetectorConstruction class to build a detection geometry similar to those usually used in typical NRF measurements [1,18]. In such detection geometry, we positioned four detectors, each of radius $R_d$ and length ${\lambda }_d$, such that their surfaces were at a distance d from the center of the target. Furthermore, the center of the target was located at the center of the detection system. The four detectors were positioned at $\left(\theta ,\varphi \right)$ of $\left(\pi /2,0\right)$, $\left(\pi /2,\pi /2\right)$, $\left(3\pi /4,\pi \right)$, and $\left(3\pi /4,3\pi /2\right)$. The numbers of geantinos arriving at each detector, $C\left(\theta ,\varphi \right)$, were counted at the end of the simulation. Then, the simulation-expected values of the asymmetry ratios ${\mathcal{A}}_{\mathbf{Sim}}\left(\pi /2\right)$ and ${\mathcal{A}}_{\mathbf{Sim}}\left(3\pi /4\right)$, and the angular distribution ratio ${\mathcal{R}}_{\mathbf{Sim}}\left(\pi /2,3\pi /4\right)$, were calculated from

$${\mathcal{A}}_{\mathbf{Sim}}\left(\pi /2\right)=\frac{C\left(\pi /2,0\right)-C\left(\pi /2,\pi /2\right)}{C\left(\pi /2,0\right)+C\left(\pi /2,\pi /2\right)}$$

(8)

$${\mathcal{A}}_{\mathbf{Sim}}\left(3\pi /4\right)=\frac{C\left(3\pi /4,\pi \right)-C\left(3\pi /4,3\pi /2\right)}{C\left(3\pi /4,\pi \right)+C\left(3\pi /4,3\pi /2\right)}$$

(9)

$${\mathcal{R}}_{\mathbf{Sim}}\left(\pi /2,3\pi /4\right)=\frac{C\left(\pi /2,0\right)+C\left(\pi /2,\pi /2\right)}{C\left(3\pi /4,\pi \right)+C\left(3\pi /4,3\pi /2\right)}$$

(10)

The relative statistical uncertainty at a $95\%$ confidence level of each count was calculated from the following relation [21],

$$\frac{\Delta C(\theta ,\varphi )}{C(\theta ,\varphi )}=1.96{\left[\frac{N-C(\theta ,\varphi )}{N\times C(\theta ,\varphi )}\right]}^{1/2}$$

(11)

where N is the number of the simulated events. In fact, depending on the geometric conditions and the spin sequence of the transition, the number of the simulated events ranged from $1\times {10}^7$ to $5\times {10}^9$ such that the relative uncertainty was ≤1%.

3. Results and Discussion

Our method described in the previous section can generally be used for arbitrary spin sequences and detection geometries. Here we show the results of Monte Carlo simulation of four examples. Two of them are for the dipole and quadrupole transitions of even-even nuclei as examples of ground states with integer spins. The spin sequences $0\to 1\to 0$ and $0\to 2\to 0$ are discussed. The other two examples illustrate how the detection geometry affects the observed quantities for nuclei of ground state with a spin of $\frac{5}{2}$, as example for ground states with half-integer spins. The spin sequences $\frac{5}{2}\to \frac{5}{2}\to \frac{5}{2}$ and $\frac{5}{2}\to \frac{9}{2}\to \frac{5}{2}$ are described.

3.1. Integer-Spin Ground States

To demonstrate the effect of the detection geometry on the angular distribution in an NRF measurement, we performed a large number of simulations at different geometrical conditions. Figure 4 shows how the asymmetry and angular distribution ratios of the dipole transition corresponding to the spin sequence $0\to 1\to 0$ are affected by the detector radius and target-to-detector distance. Generally, the asymmetry ratios decreases as the target-to-detector distance decreases. This is because when the distance between the target and the detector decreases, the solid angle that the target subtends the detector increases. Thus, the number of photons counted within the detector are not only those photons emitted at an angle of $\theta$ but also photons emitted at $\theta \pm \Delta \theta$, where $\Delta \theta$ is the angular acceptance of the detector as seen by the target. However, the decrease of ${\mathcal{A}}_{\mathbf{Sim}}(\pi /2)$ is faster than ${\mathcal{A}}_{\mathbf{Sim}}(3\pi /4)$ when decreasing the target-to-detector distance. For example, when the detector radius is 30 mm, ${\mathcal{A}}_{\mathbf{Sim}}(\pi /2)$ decreases to $80\%$ of the ideal value when the d changes from 5 to 30 mm, while ${\mathcal{A}}_{\mathbf{Sim}}(3\pi /4)$ decreases only to $87\%$ of the ideal value by the change over the same distance. This difference arises from the fact that at $\pi /2$, there is a maximum value of ${\mathcal{A}}_{\mathbf{Sim}}\left(\theta \right)$, which means ${\mathcal{A}}_{\mathbf{Sim}}\left(\theta \right)$ has lower values at either side; that is, ${\mathcal{A}}_{\mathbf{Sim}}(\theta <\pi /2)<{\mathcal{A}}_{\mathbf{Sim}}(\theta =\pi /2)>{\mathcal{A}}_{\mathbf{Sim}}(\theta >\pi /2)$, as shown in Figure 2. On the contrary, at $3\pi /4$ there is no maxima which leads to ${\mathcal{A}}_{\mathbf{Sim}}(\theta <3\pi /4)<{\mathcal{A}}_{\mathbf{Sim}}(\theta =3\pi /4)<{\mathcal{A}}_{\mathbf{Sim}}(\theta <\pi /2)$. In practice, an NRF experiment may be conducted using a close geometry with small target-to-detector distance and big detector aiming at enhancing the counting rates. For such close geometry ($d=6$ mm and $R_d=60$ mm) ${\mathcal{A}}_{\mathbf{Sim}}(\theta =\pi /2)$ drops to $53\%$ of its ideal value, while ${\mathcal{A}}_{\mathbf{Sim}}(\theta =3\pi /4)$ drops to $63\%$ of its ideal value. This reduction reflects the significance of taking the extended detector into consideration.

In contrast to the asymmetry ratios, the angular distribution ratio ${\mathcal{R}}_{\mathbf{Sim}}(\pi /2,3\pi /4)$ increases when the distance d gets shorter. From the definition of ${\mathcal{R}}_{\mathbf{Sim}}(\pi /2,3\pi /4)$, Equations (7) and (10), the second term of the numerator, $W(\pi /2,\pi /2)$, vanishes for the dipole transition with spin sequence $0\to 1\to 0$. But for the extended detector, $W(\pi /2,\pi /2)\ne 0$ and gets larger for larger solid angles. Therefore, ${\mathcal{R}}_{\mathbf{Sim}}(\pi /2,3\pi /4)$ increases at larger solid angle (smaller target-to-detector distance).

To demonstrate the effect of the dimensions of the target, we carried out some simulations at different target lengths and radii while fixing the detector radius at 45 mm. Figure 5 shows the variation of the asymmetry and angular distribution ratios when the target length is changed from 5 to 70 mm. It should be emphasized that, in real experiments, the target length is selected based on the density of the target and the expected counting rates of the detectors. For dense targets, e.g., Pb or Ta, several-millimeter-long targets are usually sufficient to produce reasonable counting rates in NRF experiments. However, if the target is of low density, e.g., Al or Sn, several-centimeter-long targets are usually used. Therefore, the range of target lengths shown in Figure 5 represents a wide range of values that are used by the experiment. The target length decreases ${\mathcal{A}}_{\mathbf{Sim}}(\pi /2)$ by a factor of $15\%$ over the whole range of target lengths considered in the present study. The effect of target length is smaller at a greater target-to-detector distance (decrease of $15\%$ and $11\%$ at $d=7$ and $d=12$ cm, respectively). ${\mathcal{R}}_{\mathbf{Sim}}(\pi /2,3\pi /4)$ decreases as well with the increase of the target length by a factor of $17\%$ at $d=7$ and $8\%$ at $d=12$ cm. Opposing ${\mathcal{A}}_{\mathbf{Sim}}(\pi /2)$, ${\mathcal{A}}_{\mathbf{Sim}}(3\pi /4)$ increases by a factor of approximately $20\%$, almost independent of d, with the increase of the target length. This may be explained in terms of the inclined position of the detector located at $3\pi /4$ with respect to the target.

Figure 6 demonstrates that the target radius has almost no effect on the asymmetry and the angular distribution ratios. Only a change of <3.5% occurs when the target radius changes from $2.5$ to 15 mm. Actually, this change is too small to be observed experimentally because the experimental uncertainties accompanying NRF experiments are in the order of 5–8%. The negligible effect of the target radius may be due to the symmetry of the target shape with respect to the radiation pattern of the transition. Further study is required to investigate the effect of the target shape. The present Monte Carlo simulation can be used to extract information about the effect of the target shape.

Similarly, the effects of the geometrical conditions on the angular distribution of the quadrupole transition with a spin sequence of $0\to 2\to 0$ are shown in Figure 7, Figure 8 and Figure 9. As indicated in Figure 7, both asymmetry ratios depart from the ideal values with shorter target-to-detector distances, following the same trend of the dipole transition. However, ${\mathcal{A}}_{\mathbf{Sim}}(3\pi /4)$ decreases at much faster pace than ${\mathcal{A}}_{\mathbf{Sim}}(\pi /2)$, which is just the opposite to the case of spin sequence $0\to 1\to 0$. Repeatedly, this can be accounted for by the aid of Figure 2 in which the asymmetry of the quadrupole transition exhibits almost a flat-top near $\pi /2$ while it changes rapidly near the minimum at $3\pi /4$. For the close geometry ${\mathcal{A}}_{\mathbf{Sim}}(\pi /2)$, it drops to $41\%$ of its ideal value. The reduction is greater than that which occurs for the dipole radiation. This can be understood in terms of the versatility of directions that the quadrupole transition is emitted in. Note that the ideal dipole transition is emitted within the scattering plane while the quadrupole transition is emitted in two orthogonal planes.

Interestingly, at $d=6$ mm and $R_d=60$ mm, ${\mathcal{A}}_{\mathbf{Sim}}(3\pi /4)$ is changed from the ideal value of $-1$ (point detector) to $+0.1$ (extended detector) which represents a change of a factor of $110\%$. More importantly, this value is very close to that of the dipole transition under the same conditions, ${\mathcal{A}}_{\mathbf{Sim}}{(3\pi /4)|}_{0\to 1\to 0}=0.2$. This indicates that unless the counting statistics are sufficiently high, it is possible to misassign the spin-parity. Note that ${\mathcal{A}}_{\mathbf{Sim}}(\pi /2)$ cannot be used to differentiate dipoles and quadrupoles, as mentioned earlier. Furthermore, it would be more difficult to use the angular distribution ratio in this case too. This is because, at the same conditions mentioned above, ${\mathcal{R}}_{\mathbf{Sim}}{(\pi /2,3\pi /4)|}_{0\to 1\to 0}=0.79$, which is very close to ${\mathcal{R}}_{\mathbf{Sim}}{(\pi /2,3\pi /4)|}_{0\to 2\to 0}=0.9$. It is clear that it would be too difficult to assign spin-parity using the angular distribution functions only correctly. Indeed, in such close geometries, Monte Carlo simulation is necessary.

Unlike the dipole case, ${\mathcal{R}}_{\mathbf{Sim}}(\pi /2,3\pi /4)$ decreases when the target-to-detector distance decreases. This can be interpreted in terms of the difference between the radiation patterns of the dipole and quadrupole transitions. The radiation pattern of the dipole is limited to the scattering plane. But the radiation pattern of the quadrupole includes the orthogonal plane in addition to the scattering plane. Consequently, $W(\pi /2,\pi /2)$ does not vanish and is affected by the target-to-detector distance similar to $W(\pi /2,0)$. The effect of the target length becomes more apparent for the quadrupole case. As shown in Figure 8, changing the target length from 5 mm to 70 mm causes an approximate change of $30\%$ in ${\mathcal{A}}_{\mathbf{Sim}}(\pi /2)$ when $d=7$ cm and $18\%$ at $d=12$ cm, which is close to twice the change occurred for the dipole transition at fixed detector radius. Furthermore, ${\mathcal{A}}_{\mathbf{Sim}}{(3\pi /4)|}_{R_d=60\mathrm{mm}}$ changes from $-0.14$ to $0.07$ over the considered range of target lengths. Additionally, the change of the target length decreases ${\mathcal{R}}_{\mathbf{Sim}}(\pi /2,3\pi /4)$ by a factor of approximately $30\%$. Finally, the quadrupole transition exhibits identical behavior to that of the dipole transition regarding the effect of the target radius. The change over the target radius only alters the asymmetry ratios by ≤2.7%.

3.2. Half-Integer-Spin Ground States

The effects of the geometrical parameters for the dipole transition with spin sequence of $\frac{5}{2}\to \frac{5}{2}\to \frac{5}{2}$ follow the same trend of the even-even nuclei as is indicated in Figure 9, Figure 10, Figure 11 and Figure 12. The asymmetry ratios decrease with decreasing the target-to-detector distance, and also with increasing the detector radius. The opposite occurs for the angular distribution ratio. Generally, it is clear that ${\mathcal{A}}_{\mathbf{Sim}}\left(\theta \right)\le {\mathcal{A}}_{\mathbf{Ideal}}\left(\theta \right)$, while ${\mathcal{R}}_{\mathbf{Sim}}({\theta }_1,{\theta }_2)\ge {\mathcal{R}}_{\mathbf{Ideal}}({\theta }_1,{\theta }_2)$ for the to dipole transitions regardless the spin sequence. However, in odd nuclei, the ideal values of the asymmetry and angular distribution ratios are often lower than those of even-even nuclei. Moreover, the selection rules of dipole transitions in odd nuclei increase the possibilities of excitation to various levels with different spins. The ground state spin of $\frac{5}{2}$ may be excited to a state of spin $\frac{3}{2}$, which is almost isotropic.

The quadrupole transition with a spin sequence of $\frac{5}{2}\to \frac{9}{2}\to \frac{5}{2}$ manifests a rather complicated dependence on the geometrical parameters, as indicated in Figure 13, Figure 14 and Figure 15. While ${\mathcal{A}}_{\mathbf{Sim}}(\pi /2)$ is similar to previously discussed quadrupole transition of even-even nuclei, ${\mathcal{A}}_{\mathbf{Sim}}(3\pi /4)$ shows a dramatic difference in comparison with the case of spin sequence $0\to 2\to 0$. At small distances, ${\mathcal{A}}_{\mathbf{Sim}}(3\pi /4)$ increases with increasing target-to-detector distance. Then, at a distance that depends on the detector radius, it decreases with increasing the target-to-detector distance. This behavior is more distinctive than the case of $0\to 2\to 0$, which showed a monotonic decrease with the target-to-detector distance. The irregular trend of ${\mathcal{A}}_{\mathbf{Sim}}(3\pi /4)$ may be explained by nature of the quadrupole transition associated with spin sequence $\frac{5}{2}\to \frac{9}{2}\to \frac{5}{2}$. There is a minimum value of the angular correlation function at $W(3\pi /4,0)$, (see for example Equation (A6)). The effect of this minimum value may be propagated at large solid angles when the detector is close to the target.

We performed further simulations with an asymmetry ratio at a different angle, which was ${\mathcal{A}}_{\mathbf{Sim}}(2\pi /3)$, where the angular correlation function shows no minima. The results of our simulations are depicted in

Table 1. A comparison among ${\mathcal{A}}_{\mathbf{Sim}}(\pi /2)$, ${\mathcal{A}}_{\mathbf{Sim}}(2\pi /3)$, and ${\mathcal{A}}_{\mathbf{Sim}}(3\pi /4)$ for the spin sequence $\frac{5}{2}\to \frac{9}{2}\to \frac{5}{2}$ as affected by the target-to-detector distance d. The other parameters are fixed such that $R_d=60$ mm, ${\lambda }_t=10$ mm, and $R_t=10$ mm. Numbers in parentheses are the uncertainties referred in the least significant digit(s).

d cm	${\mathcal{A}}_{\mathbf{Sim}}(\mathit{\pi }/2)$	${\mathcal{A}}_{\mathbf{Sim}}(2\mathit{\pi }/3)$	${\mathcal{A}}_{\mathbf{Sim}}(3\mathit{\pi }/4)$
6	$0.18216$$\left(6\right)$	$0.12650$$\left(8\right)$	$0.07810$$\left(3\right)$
7	$0.20651$$\left(7\right)$	$0.13845$$\left(7\right)$	$0.08184$$\left(4\right)$
8	$0.22689$$\left(3\right)$	$0.14784$$\left(3\right)$	$0.08409$$\left(4\right)$
9	$0.24395$$\left(7\right)$	$0.15535$$\left(9\right)$	$0.08532$$\left(6\right)$
10	$0.25809$$\left(8\right)$	$0.16121$$\left(06\right)$	$0.08599$$\left(5\right)$
12	$0.27981$$\left(11\right)$	$0.16971$$\left(14\right)$	$0.08635$$\left(8\right)$
14	$0.29502$$\left(14\right)$	$0.17549$$\left(10\right)$	$0.08614$$\left(9\right)$
16	$0.30585$$\left(16\right)$	$0.17946$$\left(13\right)$	$0.08564$$\left(18\right)$
18	$0.31378$$\left(22\right)$	$0.18224$$\left(3\right)$	$0.08529$$\left(10\right)$
20	$0.31993$$\left(18\right)$	$0.18380$$\left(87\right)$	$0.08499$$\left(10\right)$
25	$0.32970$$\left(19\right)$	$0.18765$$\left(25\right)$	$0.08441$$\left(82\right)$
30	$0.33528$$\left(15\right)$	$0.18925$$\left(13\right)$	$0.08388$$\left(14\right)$
Ideal value	$0.3489$	$0.1937$	$0.0826$

, along with those of ${\mathcal{A}}_{\mathbf{Sim}}(\pi /2)$ and ${\mathcal{A}}_{\mathbf{Sim}}(3\pi /4)$. The data of Table 1 demonstrates that both ${\mathcal{A}}_{\mathbf{Sim}}(\pi /2)$ and ${\mathcal{A}}_{\mathbf{Sim}}(2\pi /3)$ show the same trend, while only ${\mathcal{A}}_{\mathbf{Sim}}(3\pi /4)$ is different. As a result, it is obvious that ${\mathcal{A}}_{\mathbf{Sim}}(3\pi /4)$ is a special case, and therefore, its irregular trend may be reasonable. The minimum value of the angular correlation function at $3\pi /4$ also affects the angular distribution ratio for the present quadrupole transition, as inferred from Figure 13. The effects of the target length and target radius are similar to even-even nuclei, either for the case of the dipole or the quadrupole transitions.

The irregular behavior of ${\mathcal{A}}_{\mathbf{Sim}}(3\pi /4)$ affirms the significance of the Monte Carlo simulation as a beneficial tool that strongly affects the success of experiments involving spin-parity assignments. The spin-parity assignment of odd nuclei is a difficult task using NRF interaction [19]. Only a few experimental works have investigated odd nuclei, such as Rusev et al. [22] and Shizuma et al. [1]. However, the present simulation tool opens the door to study the spin and parity of odd nuclei by providing a precise estimation of the required geometry prior to the experiment. Furthermore, the present Monte Carlo method of the $\gamma -\gamma$ correlation functions has the advantages of simplicity, generality, and, at the same time, accuracy over previously implemented methods based on mapping [23] or parametrization [24,25]. It worth noting that the present study is limited to pure transitions. An upgrade would be required to implement a simulation tool that can handle the phenomenon of multipole mixing.

3.3. Benchmarking against Experimental Data

To validate our simulation, we compared the simulation-estimated asymmetry ratios with those determined experimentally for ${}^{52}$Cr and ${}^{40}$Ar as examples of even nuclei, and for ${}^{27}$Al and ${}^{11}$B as examples of odd nuclei. In each case, we could extract the dimensions of the target and the detector from the corresponding reports, and then we precisely simulated the experimental setups. In a recent study, Wilhelmy et al. [26] assigned 35 parity quantum numbers of ${}^{52}$Cr to either $M1$ or $E1$ transition. The asymmetry ratios ${\mathcal{A}}_{\mathbf{Exp}}(\pi /2)$ averaged $0.878\left(313\right)$ overall assignments. The simulation-estimated ${\mathcal{A}}_{\mathbf{Sim}}(\pi /2)$ for the corresponding setup is $0.888$, which agrees with measured asymmetry ratios within the experimental uncertainty.

Using the results of previous experiments on ${}^{40}$Ar [27,28], the effect of the target length could be appropriately benchmarked. In these experiments, a target of pressurized ${}^{40}$Ar of length 12 cm and radius 6 cm was exploited in NRF measurements while the detector radius was 32 mm. The first experiment on ${}^{40}$Ar reported that 17 dipole transitions with average ${\mathcal{A}}_{\mathbf{Exp}}(\pi /2)=0.544\left(126\right)$, while the second measurement resulted in assignments of 18 dipole transitions with ${\mathcal{A}}_{\mathbf{Exp}}(\pi /2)=0.638\left(171\right)$. Our simulation showed that ${\mathcal{A}}_{\mathbf{Sim}}(\pi /2)=0.616$ and $0.619$ for the first and second experiments, respectively. The small difference between the simulation of the two experiments is due to the difference between the diameters of collimators used to define the $\gamma$-ray beam size ($2.54$ cm for the first experiment and 1.91 cm for the second experiment). Nonetheless, the asymmetry ratios predicted by the simulation agree within the experimental uncertainty with the measured counterparts.

Shizuma et al. [1] have recently assigned the spin and parity of the 2.982 and 3.004 MeV levels of ${}^{27}$Al. The first level was assigned to be $M1$ transition via a spin sequence of $\frac{5}{2}\to \frac{3}{2}\to \frac{5}{2}$ with ${\mathcal{A}}_{\mathbf{Exp}}(\pi /2)=0.02\left(5\right)$. The simulation-calculated value of this level is ${\mathcal{A}}_{\mathbf{Sim}}(\pi /2)=0.015$. For the 3.004 MeV level, an $E2$ transition was identified, having ${\mathcal{A}}_{\mathbf{Exp}}(\pi /2)=0.34\left(5\right)$ with a spin sequence of $\frac{5}{2}\to \frac{9}{2}\to \frac{5}{2}$. This demonstrates a reasonable consistency with the value predicted by the proposed simulation (${\mathcal{A}}_{\mathbf{Sim}}(\pi /2)=0.331$). Furthermore, the angular distribution ratio ${\mathcal{R}}_{\mathbf{Exp}}(\pi /2,3\pi /4)$ was determined from the experiment amounted as ${\mathcal{R}}_{\mathbf{Exp}}(\pi /2,3\pi /4)=0.98\left(7\right)$ for the spin sequence $\frac{5}{2}\to \frac{3}{2}\to \frac{5}{2}$ and ${\mathcal{R}}_{\mathbf{Exp}}(\pi /2,3\pi /4)=0.95\left(7\right)$ for the spin sequence $\frac{5}{2}\to \frac{9}{2}\to \frac{5}{2}$. Again, the simulation-estimated angular distribution ratio is ${\mathcal{R}}_{\mathbf{Sim}}(\pi /2,3\pi /4)=0.97$ and ${\mathcal{R}}_{\mathbf{Sim}}(\pi /2,3\pi /4)=0.91$ for dipole and quadruple transitions, respectively.

Another example of the odd nuclei is the experiment performed by Rusev et al. [22] on ${}^{11}$B. In this experiment the level at $8.916$ MeV was assigned as $M1$ transition with ${\mathcal{A}}_{\mathbf{Exp}}(\pi /2)=0.215\left(20\right)$. The spin sequence of this transition was found to be $\frac{3}{2}\to \frac{5}{2}\to \frac{3}{2}$. Our simulation resulted in a consistent asymmetry ratio (${\mathcal{A}}_{\mathbf{Sim}}(\pi /2)=0.223$) for the given spin sequence.

A general feature of spin-parity measurements is the symmetry of the detection geometry. Usually, two identical detectors are positioned at the same distance from the center of the target, with one of them being parallel to the polarization plane while the other is perpendicular to the polarization plane. In such geometry, NRF photons emitted from the target travel the same distance within the target and through the target-to-detector distance before arriving at the detector. Therefore, the effects of self-attenuation of photons within the target cancel each other out when asymmetry and angular distribution ratios are calculated. This scenario enabled the proposed simulation to calculate the asymmetry and angular distribution ratios without the need to track the propagation of photons through the target. The consistency between the experimental observation of the asymmetry and angular distribution ratios with the simulated-estimated values strengthens the fact that the proposed simulation can be used as a planning utility prior to the experiment. Another utilization of the proposed simulation is to estimate the corrections of angular distribution functions after the experiment as an indispensable part of the data analysis procedure. More importantly, the present simulation can be employed to sample the final states through a complete implementation that simulates the whole NRF interaction. In this case, the angular distribution functions are invoked by the simulation code only for the interacting nuclei. The simulation run within the Geant4 toolkit will then track the propagation of NRF photons in detail.

4. Conclusions

In summary, we developed a new Monte Carlo simulation of the $\gamma -\gamma$ correlation functions. Our simulation is capable of accounting for the effects of detection geometry encountered in the spin-parity assignments using NRF. The simulation can handle arbitrary shapes and sizes of the target and detector, and spin sequences. We showed that the spin-parity assignment is very sensitive to detection geometry. Consequently, incorrect spin-parity, particularity for odd-nuclei, may be assigned unless the detection geometry is carefully planned using the proposed method. The simulation provides a tool to plan and predict the sufficient statistics for correct assignment of the spins and parities of the nuclear energy states. In addition, we provided simulation examples for different geometrical conditions and spin sequences. The effects of the target and detector sizes and target-to-detector distance may lead to a deviation of the order of ≥100% of the ideal values of the $\gamma -\gamma$ correlation functions. Finally, we validated our simulation model by comparing the results of our simulations with the observable quantities involved in spin-parity measurements using NRF; namely, asymmetry and angular distribution ratios. Our validation against experimental data included a variety of detection geometries and spin sequences.