Thin-film approximate point scattered function and its application to neutron radiography

Thin-film approximate point scattered function and its application to neutron radiography 1 first-author 1 2 corresp *Jiayu_NI@cumt.edu.cn *Jiayu_NI@cumt.edu.cn Jiayu_NI@cumt.edu.cn https://orcid.org/0000-0003-1841-4044 1 2 1 2 1 Department of Nuclear Engineering and New Energy Technology, The College of Engineering and Technology of Chengdu University of Technology, Sichuan Leshan 614000, China Department of Nuclear Engineering and New Energy Technology, The College of Engineering and Technology of Chengdu University of Technology, Sichuan Leshan 614000, China Society for Theoretical and Computational Physics, The College of Engineering and Technology of Chengdu University of Technology, Sichuan Leshan 614000, China Society for Theoretical and Computational Physics, The College of Engineering and Technology of Chengdu University of Technology, Sichuan Leshan 614000, China A simpler and improved utility approximate point scattered function for thin-film converters currently used in neutron photographic devices is proposed as a correction method to produce clearer, more realistic images. The validity of the model was demonstrated through a simulation experiment. Based on the results, an error analysis was carried out, certain corrections were made to the original model, and the final model achieved a very low relative error in the simulation experiment. The model can also be optimized for quantitative neutron photographic analysis using iterative algorithms to obtain realistic neutron photographic images more quickly. At the end of the article, the model is extended to consider the case of energy spectrum hardening by introducing a temperature correction parameter. Neutron radiography Point scattered function Image processing Introduction 1 Neutron radiography is a method for surface or internal diagnosis of a sample to be measured using the properties of neutron attenuation and scattering. Compared to other radiographic techniques, neutrons are more complex because of the manner in which they interact with the nucleus, and the macroscopic cross section does not have a definite relationship with the nuclide species. These properties provide neutron radiography with an unparalleled advantage over other radiographic techniques, such that neutron radiography is currently used in the diagnosis of a wide variety of industrial materials and is of great importance in the metal and aerospace industries[ 1 1 - 3 3 1-3 ]. Simply put, neutron radiography takes advantage of the transmission of neutrons and the absorption of neutrons by the sample being evaluated, among other properties. By comparing the imaging of the detector in the presence or absence of the sample, the physical condition of the sample can be determined. However, there are inherent errors in this technique, mainly due to the scattering of neutrons by the sample, absorption of neutrons by the detector material, and attenuation of neutrons during transmission. In response to this problem, two main correction methods have been proposed: the first is to modify the energy spectrum of the neutron source to produce more realistic final detection results [ 4 4 , 5 5 ], and the second is to introduce various functions to correct for the effect of scattering in neutron radiography [ 6 6 , 7 7 ]. Energy spectrum correction is proposed mainly to counter the neutron energy spectrum hardening problem. This problem becomes more significant as the thickness of the sample to be measured increases to a point at which energy spectrum correction becomes necessary. However, in a broader context, the scattering effect is more direct and important than the energy spectrum hardening effect. Therefore, this paper discusses energy spectrum correction as a secondary correction after scattering correction. For example, the 14 MeV fast neutron radiography (FNR) technique [ 20 20 ], which has been extensively studied in recent years, has promising applications with high contrast sensitivity and other features. However, it has difficulty dealing with the impact of scattering effects, which makes this technique much less effective and restricts its use in practical applications. There are two popular methods for energy spectrum correction: the maximum-likelihood expectation-maximization (MLEM) method, and the quantitative neutron photographic analysis (QNPA) method. The MLEM method is mainly based on least squares and maximum-likelihood estimations to perform statistical correction [ 21 21 ], whereas the QNPA method uses the characteristics of neutron decay to perform iterative correction of the neutron distribution in order to correct the energy spectrum. Therefore, this paper discusses energy spectrum correction as a secondary correction after scattering correction, with the optimization of the scattering effect as the primary objective. For the scattering correction, the "point scattered function" has been proposed [ 8 8 ]. It is essentially the probabilistic density distribution function of the scattered neutrons detected by the detector after the neutron beam passes through the sample to be measured. However, it is not easy to apply directly because of its complex form and the many parameters it contains. Therefore, scattering correction models based on a variety of approximate point-scattered functions have been proposed [ 9 9 , 10 10 ]. In recent years, materials science has advanced rapidly, and an increasing number of neutron detectors can use very thin detection layers as a converter to optimize neutron radiography [ 11 11 - 13 13 11-13 ]. Thus, we propose a new approximation of the point scattered function called the "thin-film approximate point scattered function". This new model has a significant computational advantage and the error is sufficiently low for practical applications. In the following, we derive the thin-film approximate point scattered function from the conventional point scattered function and then perform numerical simulations using software such as Monte Carlo N-particle transport (MCNP) and Mathematica to evaluate the accuracy of the model and to perform error analysis. Thin-film approximate scattering model 1 Thin-film approximate point scattered function 2 We first review the basic theory of the point scattered function, which is based on four basic assumptions: 1. The scattered neutrons are homogeneous. 2. The neutrons are scattered only once in the sample. 3. The source and scattered neutrons have the same attenuation coefficient in the sample. 4. The sample width is large enough to ensure that the scattered neutrons penetrate the bottom of sample (not the sides) and reach the detector. These are reasonable assumptions. Because neutrons are electrically neutral, each should be homogeneous after scattering. The vast majority of neutrons also scatter only once in the sample and therefore the need for multiple scattering corrections is not significant. Figure 1 Figure 1 shows a schematic of neutron radiography and the angular response of the detector. In Fig. 1 Fig. 1 , r is the coordinate parameter of the detector, d is the thickness of the detector converter, θ is the scattering angle, D is the distance from the sample to the detector, T is the thickness of the sample, t is the location where scattering occurs in the sample, and R and l are both geometric parameters. The following relationship exists between the geometric parameters. <math display="block" id="M1"><mrow><mi>cos</mi><mi>θ</mi><mo>=</mo><mfrac><mrow><mi>D</mi><mo>+</mo><mi>t</mi></mrow><mi>R</mi></mfrac><mn>,</mn><mi>R</mi><mo>=</mo><msqrt><mrow><msup><mrow><mrow><mo>(</mo><mrow><mi>t</mi><mo>+</mo><mi>D</mi></mrow><mo>)</mo></mrow></mrow><mn>2</mn></msup><mo>+</mo><msup><mi>r</mi><mn>2</mn></msup></mrow></msqrt><mn>,</mn><mi>l</mi><mo>=</mo><mfrac><mrow><mi>t</mi><mi>R</mi></mrow><mrow><mi>t</mi><mo>+</mo><mi>D</mi></mrow></mfrac></mrow></math> Then, for the detector coordinate r , the incident neutron flux Φ sc ( E ) and the scattered neutron flux of the detector Φ st ( E ) can be expressed in the following form. <math display="block" id="M2"><mtable columnalign="left"><mtr><mtd><mtext>d</mtext><msub><mtext>Φ</mtext><mrow><mtext>sc</mtext></mrow></msub><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><mo>=</mo><mfrac><mrow><mtext>d</mtext><msub><mtext>Φ</mtext><mrow><mtext>st</mtext></mrow></msub><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow></mrow><mrow><mn>4</mn><mi>π</mi><msup><mi>R</mi><mn>2</mn></msup></mrow></mfrac><mo>⋅</mo><msup><mtext>e</mtext><mrow><mo>−</mo><msub><mi>Σ</mi><mtext>t</mtext></msub><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><mi>l</mi></mrow></msup></mtd></mtr><mtr><mtd><mtext>d</mtext><msub><mtext>Φ</mtext><mrow><mtext>st</mtext></mrow></msub><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><mo>=</mo><msub><mtext>Φ</mtext><mn>0</mn></msub><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><mo>⋅</mo><msup><mtext>e</mtext><mrow><mo>−</mo><msub><mi>Σ</mi><mtext>t</mtext></msub><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><mrow><mo>(</mo><mrow><mi>T</mi><mo>−</mo><mi>t</mi></mrow><mo>)</mo></mrow></mrow></msup><mo>⋅</mo><msub><mi>Σ</mi><mtext>s</mtext></msub><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><mtext>d</mtext><mi>t</mi></mtd></mtr></mtable></math> Here, Φ 0 ( E ) is the neutron flux of the source, ∑ s ( E ) and ∑ t ( E ) are the scattering macroscopic cross section and the absorbing macroscopic cross section of the sample, respectively. ∑ det ( E ) is the absorption cross section of the detector material. We have to take into account the neutron attenuation effect, where a single neutron is detected with a probability of <math display="block" id="M3"><mrow><mi>P</mi><mo>=</mo><mn>1</mn><mo>−</mo><msup><mtext>e</mtext><mrow><mo>−</mo><mi>Σ</mi><mi>d</mi><mo>/</mo><mi>cos</mi><mi>θ</mi></mrow></msup><mn>.</mn></mrow></math> The neutron flux at the detector surface is <math display="block" id="M4"><mrow><mtable columnalign="left"><mtr columnalign="left"><mtd columnalign="left"><mrow><mtext>d</mtext><msub><mi>I</mi><mrow><mtext>det</mtext></mrow></msub></mrow></mtd><mtd columnalign="left"><mrow><mo>=</mo><mstyle displaystyle="true"><mrow><msub><mo>∫</mo><mi>E</mi></msub></mrow></mstyle><mi>P</mi><mo>⋅</mo><mtext>d</mtext><msub><mtext>Φ</mtext><mrow><mtext>sc</mtext></mrow></msub><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><mtext>d</mtext><mi>E</mi></mrow></mtd></mtr><mtr columnalign="left"><mtd columnalign="left"><mrow/></mtd><mtd columnalign="left"><mrow><mo>=</mo><mstyle displaystyle="true"><mrow><msub><mo>∫</mo><mi>E</mi></msub></mrow></mstyle><msub><mtext>Φ</mtext><mn>0</mn></msub><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><mfrac><mrow><mi>D</mi><mo>+</mo><mi>t</mi></mrow><mrow><mn>4</mn><mi>π</mi><msup><mi>R</mi><mn>2</mn></msup></mrow></mfrac><mo>⋅</mo><msup><mtext>e</mtext><mrow><mo>−</mo><msub><mi>Σ</mi><mi>t</mi></msub><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><mi>l</mi></mrow></msup><msup><mtext>e</mtext><mrow><mo>−</mo><msub><mi>Σ</mi><mi>t</mi></msub><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><mrow><mo>(</mo><mrow><mi>T</mi><mo>−</mo><mi>t</mi></mrow><mo>)</mo></mrow></mrow></msup></mrow></mtd></mtr><mtr columnalign="left"><mtd columnalign="left"><mrow/></mtd><mtd columnalign="left"><mrow><msub><mi>Σ</mi><mi>s</mi></msub><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><mi>cos</mi><mi>θ</mi><mrow><mo>(</mo><mrow><mn>1</mn><mo>−</mo><msup><mtext>e</mtext><mrow><mo>−</mo><msub><mi>Σ</mi><mrow><mi>d</mi><mi>e</mi><mi>t</mi></mrow></msub><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><mi>d</mi><mo>/</mo><mi>cos</mi><mi>θ</mi></mrow></msup></mrow><mo>)</mo></mrow><mtext>d</mtext><mi>t</mi><mtext>d</mtext><mi>E</mi></mrow></mtd></mtr><mtr columnalign="left"><mtd columnalign="left"><mn>.</mn></mtd><mtd columnalign="left"><mrow/></mtd></mtr></mtable></mrow></math> Then, according to the definition of the point scattered function we have <math display="block" id="M5"><mtable><mtr><mtd><mtext>PScF</mtext><mrow><mo>(</mo><mi>r</mi><mo>)</mo></mrow><mo>=</mo><msub><mi>I</mi><mrow><mtext>det</mtext></mrow></msub><mo>/</mo><mstyle displaystyle="true"><mrow><msub><mo>∫</mo><mi>E</mi></msub></mrow></mstyle><mi>P</mi><mo>⋅</mo><msub><mtext>Φ</mtext><mn>0</mn></msub><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><mtext>d</mtext><mi>E</mi></mtd></mtr><mtr><mtd><mo>=</mo><mfrac><mrow><mstyle displaystyle="true"><mrow><msub><mo>∫</mo><mrow><mi>E</mi><mn>,</mn><mi>t</mi></mrow></msub></mrow></mstyle><msub><mtext>Φ</mtext><mn>0</mn></msub><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><mfrac><mrow><mi>D</mi><mo>+</mo><mi>t</mi></mrow><mrow><mn>4</mn><mi>π</mi><msup><mi>R</mi><mrow><mn>3</mn><mo>/</mo><mn>2</mn></mrow></msup></mrow></mfrac><msup><mtext>e</mtext><mrow><mo>−</mo><msub><mi>Σ</mi><mi>t</mi></msub><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><mrow><mo>(</mo><mrow><mfrac><mi>t</mi><mrow><mi>t</mi><mo>+</mo><mi>D</mi></mrow></mfrac><mi>R</mi><mo>+</mo><mi>T</mi><mo>−</mo><mi>t</mi></mrow><mo>)</mo></mrow></mrow></msup><mrow><mo>(</mo><mrow><mn>1</mn><mo>−</mo><msup><mtext>e</mtext><mrow><mo>−</mo><mfrac><mrow><msub><mi>Σ</mi><mrow><mi>d</mi><mi>e</mi><mi>t</mi></mrow></msub><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><mo>⋅</mo><mi>d</mi></mrow><mrow><mi>cos</mi><mi>θ</mi></mrow></mfrac></mrow></msup></mrow><mo>)</mo></mrow><msub><mi>Σ</mi><mi>s</mi></msub><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><mtext>d</mtext><mi>t</mi><mtext>d</mtext><mi>E</mi></mrow><mrow><mstyle displaystyle="true"><mrow><msub><mo>∫</mo><mi>E</mi></msub></mrow></mstyle><msub><mtext>Φ</mtext><mn>0</mn></msub><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><mo>⋅</mo><mrow><mo>(</mo><mrow><mn>1</mn><mo>−</mo><msup><mtext>e</mtext><mrow><mo>−</mo><msub><mi>Σ</mi><mrow><mi>d</mi><mi>e</mi><mi>t</mi></mrow></msub><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><mo>⋅</mo><mi>d</mi></mrow></msup></mrow><mo>)</mo></mrow><mtext>d</mtext><mi>E</mi></mrow></mfrac><mn>.</mn></mtd></mtr></mtable></math> This function is quite complicated and inconvenient to use. Many studies related to scattering correction have been conducted to approximate this function. For a very thin detector, it can be considered that <math id="M6"><mrow><mi>d</mi><mo>→</mo><mn>0</mn></mrow></math> . Because the series expansion of the exponential function <math id="M7"><mrow><msup><mtext>e</mtext><mrow><mo>−</mo><mstyle displaystyle="true"><mo>∑</mo><mi>d</mi></mstyle></mrow></msup><mo>=</mo><mn>1</mn><mo>−</mo><mstyle displaystyle="true"><mo>∑</mo><mi>d</mi></mstyle><mo>+</mo><mi mathvariant="script">O</mi><mrow><mo>[</mo><mn>2</mn><mo>]</mo></mrow></mrow></math> involves a second-order infinitesimal term <math id="M8"><mrow><mi mathvariant="script">O</mi><mrow><mo>[</mo><mn>2</mn><mo>]</mo></mrow></mrow></math> , we can use this approximation by considering that the main and first-order terms occupy 99% of the exponential function (i.e., d satisfies 1-∑ d =0.99e -∑ d ). This implies that the following approximation can be applied. <math display="block" id="M9"><mrow><mtable columnalign="left"><mtr columnalign="left"><mtd columnalign="left"><mrow/></mtd><mtd columnalign="left"><mrow><mn>1</mn><mo>−</mo><msup><mtext>e</mtext><mrow><mo>−</mo><mfrac><mrow><msub><mi>Σ</mi><mrow><mi>d</mi><mi>e</mi><mi>t</mi></mrow></msub><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><mi>d</mi></mrow><mrow><mi>cos</mi><mi>θ</mi></mrow></mfrac></mrow></msup><mo>→</mo><mfrac><mrow><msub><mi>Σ</mi><mrow><mtext>det</mtext></mrow></msub><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><mi>d</mi></mrow><mrow><mi>cos</mi><mi>θ</mi></mrow></mfrac></mrow></mtd></mtr><mtr columnalign="left"><mtd columnalign="left"><mrow/></mtd><mtd columnalign="left"><mrow><mn>1</mn><mo>−</mo><msup><mtext>e</mtext><mrow><mo>−</mo><msub><mi>Σ</mi><mrow><mi>d</mi><mi>e</mi><mi>t</mi></mrow></msub><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><mi>d</mi></mrow></msup><mo>→</mo><msub><mi>Σ</mi><mrow><mtext>det</mtext></mrow></msub><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><mi>d</mi></mrow></mtd></mtr></mtable></mrow></math> Using the geometric relations and the adjustment order, and writing the new function as TPScF ( r ), we obtain <math display="block" id="M10"><mtable columnalign="left"><mtr><mtd><mtext>TPScF</mtext><mrow><mo>(</mo><mi>r</mi><mo>)</mo></mrow><mo>=</mo><mfrac><mrow><mstyle displaystyle="true"><mrow><msub><mo>∫</mo><mrow><mi>E</mi><mn>,</mn><mi>t</mi></mrow></msub></mrow></mstyle><msub><mtext>Φ</mtext><mn>0</mn></msub><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><mo>⋅</mo><msub><mi>Σ</mi><mrow><mtext>det</mtext></mrow></msub><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><mo>⋅</mo><msub><mi>Σ</mi><mtext>s</mtext></msub><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><mo>⋅</mo><mi>g</mi><mrow><mo>(</mo><mrow><mi>E</mi><mn>,</mn><mi>t</mi></mrow><mo>)</mo></mrow><mtext>d</mtext><mi>t</mi><mtext>d</mtext><mi>E</mi></mrow><mrow><mstyle displaystyle="true"><mrow><msub><mo>∫</mo><mi>E</mi></msub></mrow></mstyle><msub><mtext>Φ</mtext><mn>0</mn></msub><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><mo>⋅</mo><msub><mi>Σ</mi><mrow><mtext>det</mtext></mrow></msub><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><mtext>d</mtext><mi>E</mi></mrow></mfrac></mtd></mtr><mtr><mtd><mi>g</mi><mrow><mo>(</mo><mrow><mi>E</mi><mn>,</mn><mi>t</mi></mrow><mo>)</mo></mrow><mo>=</mo><mfrac><mn>1</mn><mrow><mn>4</mn><mi>π</mi><msqrt><mi>R</mi></msqrt></mrow></mfrac><msup><mtext>e</mtext><mrow><mo>−</mo><msub><mi>Σ</mi><mtext>t</mtext></msub><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><mrow><mo>(</mo><mrow><mfrac><mrow><mi>t</mi><mi>R</mi></mrow><mrow><mi>t</mi><mo>+</mo><mi>D</mi></mrow></mfrac><mo>+</mo><mi>T</mi><mo>−</mo><mi>t</mi></mrow><mo>)</mo></mrow></mrow></msup><mn>.</mn></mtd></mtr></mtable></math> The integral in the numerator of Eq. 7 with respect to t can be written as <math display="block" id="M11"><mrow><mtable columnalign="left"><mtr columnalign="left"><mtd columnalign="left"><mrow/></mtd><mtd columnalign="left"><mrow><mi>F</mi><mrow><mo>(</mo><mrow><mi>E</mi><mn>,</mn><mi>r</mi></mrow><mo>)</mo></mrow><mo>=</mo><mstyle displaystyle="true"><mrow><msubsup><mo>∫</mo><mn>0</mn><mi>T</mi></msubsup></mrow></mstyle><mi>f</mi><mrow><mo>(</mo><mrow><mi>E</mi><mn>,</mn><mi>r</mi><mn>,</mn><mi>t</mi></mrow><mo>)</mo></mrow><mtext>d</mtext><mi>t</mi></mrow></mtd></mtr><mtr columnalign="left"><mtd columnalign="left"><mrow/></mtd><mtd columnalign="left"><mrow><mo>=</mo><mstyle displaystyle="true"><mrow><msubsup><mo>∫</mo><mn>0</mn><mi>T</mi></msubsup></mrow></mstyle><mfrac><mn>1</mn><mrow><mn>4</mn><mi>π</mi><msqrt><mrow><mi>R</mi><mrow><mo>(</mo><mrow><mi>r</mi><mn>,</mn><mi>t</mi></mrow><mo>)</mo></mrow></mrow></msqrt></mrow></mfrac><msup><mtext>e</mtext><mrow><mo>−</mo><msub><mi>Σ</mi><mi>t</mi></msub><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><mrow><mo>(</mo><mrow><mfrac><mrow><mi>t</mi><mi>R</mi><mrow><mo>(</mo><mrow><mi>r</mi><mn>,</mn><mi>t</mi></mrow><mo>)</mo></mrow></mrow><mrow><mi>t</mi><mo>+</mo><mi>D</mi></mrow></mfrac><mo>+</mo><mi>T</mi><mo>−</mo><mi>t</mi></mrow><mo>)</mo></mrow></mrow></msup><mtext>d</mtext><mi>t</mi><mn>.</mn></mrow></mtd></mtr></mtable></mrow></math> This integral can be calculated numerically using the Newton-Cotes method after the associated parameters are given. Finally, the thin-film approximation of the model is <math display="block" id="M12"><mrow><mtext>TPScF</mtext><mrow><mo>(</mo><mi>r</mi><mo>)</mo></mrow><mo>=</mo><mfrac><mrow><mstyle displaystyle="true"><mrow><msub><mo>∫</mo><mi>E</mi></msub></mrow></mstyle><msub><mtext>Φ</mtext><mn>0</mn></msub><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><mo>⋅</mo><msub><mi>Σ</mi><mrow><mtext>det</mtext></mrow></msub><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><mo>⋅</mo><msub><mi>Σ</mi><mi>s</mi></msub><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><mo>⋅</mo><mi>F</mi><mrow><mo>(</mo><mrow><mi>E</mi><mn>,</mn><mi>r</mi></mrow><mo>)</mo></mrow><mtext>d</mtext><mi>E</mi></mrow><mrow><mstyle displaystyle="true"><mrow><msub><mo>∫</mo><mi>E</mi></msub></mrow></mstyle><msub><mtext>Φ</mtext><mn>0</mn></msub><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><mo>⋅</mo><msub><mi>Σ</mi><mrow><mtext>det</mtext></mrow></msub><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><mtext>d</mtext><mi>E</mi></mrow></mfrac><mn>.</mn></mrow></math> Here, the integral of E is related to the specific energy spectrum. Pre-processing of the model 2 To better apply this model, we need to pre-process the model, as described in this section. We focus on Eq. 8 and apply the Newton-Cotes method to obtain a numerical solution [ 14 14 ] using <math display="block" id="M13"><mtable><mtr><mtd><mi>F</mi><mrow><mo>(</mo><mrow><mi>E</mi><mn>,</mn><mi>r</mi></mrow><mo>)</mo></mrow><mo>≈</mo><mfrac><mi>h</mi><mn>3</mn></mfrac><mrow><mo>(</mo><mrow><mi>f</mi><mrow><mo>(</mo><mrow><mi>E</mi><mn>,0</mn></mrow><mo>)</mo></mrow><mo>+</mo><mi>f</mi><mrow><mo>(</mo><mrow><mi>E</mi><mn>,</mn><mi>T</mi></mrow><mo>)</mo></mrow></mrow></mrow></mtd></mtr><mtr><mtd><mo>+</mo><mn>4</mn><mo>⋅</mo><mstyle displaystyle="true"><munderover><mo>∑</mo><mrow><mi>n</mi><mo>=</mo><mn>1</mn></mrow><mi>m</mi></munderover></mstyle><mi>f</mi><mrow><mo>[</mo><mrow><mi>E</mi><mn>,</mn><mrow><mo>(</mo><mrow><mn>2</mn><mi>n</mi><mo>−</mo><mn>1</mn></mrow><mo>)</mo></mrow><mi>h</mi></mrow><mo>]</mo></mrow></mtd></mtr><mtr><mtd><mrow><mrow><mo>+</mo><mn>2</mn><mo>⋅</mo><mstyle displaystyle="true"><munderover><mo>∑</mo><mrow><mi>n</mi><mo>=</mo><mn>1</mn></mrow><mrow><mi>m</mi><mo>−</mo><mn>1</mn></mrow></munderover></mstyle><mi>f</mi><mrow><mo>[</mo><mrow><mi>E</mi><mn>,2</mn><mi>n</mi><mi>h</mi></mrow><mo>]</mo></mrow></mrow><mo>)</mo></mrow><mo>,</mo></mtd></mtr></mtable></math> where h represents the step size and m is obtained by solving the following equation. <math display="block" id="M14"><mrow><mrow><mo>(</mo><mrow><mn>2</mn><mi>n</mi><mo>−</mo><mn>1</mn></mrow><mo>)</mo></mrow><msub><mi>h</mi><mtext>m</mtext></msub><mo>=</mo><mi>T</mi><mn>,</mn><mi>m</mi><mo>=</mo><mrow><mo>[</mo><mrow><msub><mi>h</mi><mtext>m</mtext></msub></mrow><mo>]</mo></mrow></mrow></math> Here, [ x ]represents the Greatest Integer Function, which returns the value of the greatest integer, which is less than or equal to x . The error that results from Eq. 10 is given by the following equation. <math display="block" id="M15"><mrow><mtext>Error</mtext><mrow><mo>(</mo><mi>x</mi><mo>)</mo></mrow><mo>=</mo><mfrac><mrow><mi>T</mi><msup><mi>h</mi><mn>4</mn></msup></mrow><mrow><mn>180</mn></mrow></mfrac><mfrac><mrow><msup><mo>∂</mo><mn>4</mn></msup></mrow><mrow><mo>∂</mo><msup><mi>t</mi><mn>4</mn></msup></mrow></mfrac><mi>f</mi><mrow><mo>(</mo><mrow><mi>E</mi><mn>,</mn><mi>ξ</mi></mrow><mo>)</mo></mrow><mn>,</mn><mi>ξ</mi><mo>∈</mo><mrow><mo>(</mo><mrow><mn>0,</mn><mi>T</mi></mrow><mo>)</mo></mrow></mrow></math> Eq. 10 can be calculated very quickly once the determined thickness T of the sample and the absorption cross section ∑ t ( E ) are known. In addition, the integration of TPScF ( r ) is difficult. Therefore, the approach we adopt is to divide the integration interval into many small intervals as follows. <math display="block" id="M16"><mrow><mtext>TPScF</mtext><mrow><mo>(</mo><mi>r</mi><mo>)</mo></mrow><mo>=</mo><mfrac><mrow><mstyle displaystyle="true"><msubsup><mo>∑</mo><mrow><mi>i</mi><mo>=</mo><mn>0</mn></mrow><mrow><msub><mi>E</mi><mtext>m</mtext></msub></mrow></msubsup><mrow><msub><mtext>Φ</mtext><mn>0</mn></msub><mrow><mo>(</mo><mrow><msub><mi>E</mi><mi>i</mi></msub></mrow><mo>)</mo></mrow><msub><mi>Σ</mi><mrow><mtext>det</mtext></mrow></msub><mrow><mo>(</mo><mrow><msub><mi>E</mi><mi>i</mi></msub></mrow><mo>)</mo></mrow><msub><mi>Σ</mi><mtext>s</mtext></msub><mrow><mo>(</mo><mrow><msub><mi>E</mi><mi>i</mi></msub></mrow><mo>)</mo></mrow><mi>F</mi><mrow><mo>(</mo><mrow><msub><mi>E</mi><mi>i</mi></msub><mn>,</mn><mi>r</mi></mrow><mo>)</mo></mrow></mrow></mstyle></mrow><mrow><mstyle displaystyle="true"><msubsup><mo>∑</mo><mrow><mi>i</mi><mo>=</mo><mn>0</mn></mrow><mrow><msub><mi>E</mi><mtext>m</mtext></msub></mrow></msubsup><mrow><msub><mtext>Φ</mtext><mn>0</mn></msub><mrow><mo>(</mo><mrow><msub><mi>E</mi><mi>i</mi></msub></mrow><mo>)</mo></mrow><msub><mi>Σ</mi><mrow><mtext>det</mtext></mrow></msub><mrow><mo>(</mo><mrow><msub><mi>E</mi><mi>i</mi></msub></mrow><mo>)</mo></mrow></mrow></mstyle></mrow></mfrac></mrow></math> here, <math id="M17"><mrow><mi>E</mi><mo>∈</mo><mrow><mo>(</mo><mrow><msub><mi>E</mi><mn>0</mn></msub><mn>,</mn><msub><mi>E</mi><mtext>m</mtext></msub></mrow><mo>)</mo></mrow></mrow></math> is the range of values of the energy spectrum of the neutron source. Thus, we consider the model from a numerical perspective, making it less abstract and more computable in concrete terms. Next, we use examples to demonstrate the validity of the model and then perform an error analysis. Application example 1 Numerical simulation of scattering correction for water and polyethylene 2 We construct a simple neutron radiography experiment. The neutron converter is made of 100 μ m-thick 6 Li material, and cylindrical metallic aluminum with a radius of 0.6 cm and a height of 1.1 cm is used as the sample, which is located 0.8 cm away from the detector. The experimental setup is shown in Fig. 2 Fig. 2 . As the neutron source in the experiment, we use thermal neutrons ranging in energy from 0.01 eV to 1 eV and whose energy spectrum satisfies Maxwell’s spectrum. The neutron distribution on the converter can be obtained by importing the relevant parameters into MCNP. The results are represented as density maps for PScF ( r ) and TPScF ( r ), in Fig. 1 Fig. 1 . In the figure, r indicates the distance from the detection point to the center of the detector plane. Fig. 3(a) Fig. 3(a) represents the simulation results using PScF ( r ) as a point scattered function, and Fig. 3(b) Fig. 3(b) represents TPScF ( r ) for comparison purposes. The red areas in Fig. 3 Fig. 3 indicate heavy neutron scattering, whereas the blue areas indicate very light neutron scattering. The boundary of the red area lies roughly at 0.6 cm, which is generally in line with the actual condition of the sample. Note that the boundaries of the red areas in each figure are not particularly clear. The causes behind this are the idealization of the neutron energy spectrum and the haphazard nature of the scattering in the sample. However, there are well-established methods to solve this problem, mainly by correcting the energy spectrum and by using computer algorithms to make the edges clearer, thus making neutron imaging more representative of the actual experimental conditions [ 15 15 , 16 16 ]. Comparing Fig. 3(a) Fig. 3(a) with (b), it can be seen that TPScF ( r ) gives an overall smaller result than PScF ( r ), which indicates that an overall error exists between them. The Relative Error Function is therefore defined as <math display="block" id="M18"><mrow><mi>σ</mi><mrow><mo>(</mo><mi>r</mi><mo>)</mo></mrow><mo>=</mo><mfrac><mrow><mrow><mo>|</mo><mrow><mtext>TPScF</mtext><mrow><mo>(</mo><mi>r</mi><mo>)</mo></mrow><mo>−</mo><mtext>PScF</mtext><mrow><mo>(</mo><mi>r</mi><mo>)</mo></mrow></mrow><mo>|</mo></mrow></mrow><mrow><mtext>PScF</mtext><mrow><mo>(</mo><mi>r</mi><mo>)</mo></mrow></mrow></mfrac><mn>.</mn></mrow></math> For the experiment described in Fig. 2 Fig. 2 , the image of the relative error function is mapped in Fig. 4 Fig. 4 This result indicates that the error is very small in the central area of the screen and is larger along the edge areas. The average error over the entire detection screen is approximately 0.306. However, this model is valid because of the similarity between Figs. 3(a) and (b) Figs. 3(a) and (b) , so, we further review the validity of the model through error analysis below. Error analysis 2 The error function presented in Fig. 4 Fig. 4 resembles an exponential function in the radial direction, so we can introduce an exponential function to adjust for the error. First, we define <math display="block" id="M19"><mrow><mtext>NPScF</mtext><mrow><mo>(</mo><mi>r</mi><mo>)</mo></mrow><mo>=</mo><mtext>TPScF</mtext><mrow><mo>(</mo><mi>r</mi><mo>)</mo></mrow><mo>+</mo><mi>C</mi><msup><mtext>e</mtext><mrow><mo>−</mo><mi>a</mi><mi>r</mi></mrow></msup></mrow></math> as the new approximate point scattered function. By fitting the relative error data, the values of the parameters C, a can be obtained, and then we have <math display="block" id="M20"><mrow><mtext>NPScF</mtext><mrow><mo>(</mo><mi>r</mi><mo>)</mo></mrow><mo>=</mo><mtext>TPScF</mtext><mrow><mo>(</mo><mi>r</mi><mo>)</mo></mrow><mo>+</mo><mn>8.7646</mn><mo>×</mo><msup><mrow><mn>10</mn></mrow><mrow><mo>−</mo><mn>5</mn></mrow></msup><msup><mtext>e</mtext><mrow><mn>0.3980</mn><mi>r</mi></mrow></msup><mn>.</mn></mrow></math> Then, to study the error relative to the original scattered function, the distribution is plotted in Fig. 5 Fig. 5 . The figure indicates that the relative errors are at a very low level in most regions, and a maximum value of 17% and a mean value of 2.5% exist across the entire region. This result indicates that a further improvement can be made to the approximation model, which validates that NPScF ( r ) is a very effective point scattered function for the thin-film approximation. Neutron imaging correction 1 Quantitative neutron photographic analysis method 2 Owing to the nature of the interaction between neutrons and matter, the results of neutron imaging are often very inaccurate in actuality. We can reduce the error by using more advanced experimental equipment, and/or we can improve the imaging technique using algorithms to correct an image. One such method, the Quantitative Neutron Photographic Analysis Method uses algorithms to improve imaging technology. Quantitative neutron photographic analysis uses the exponential neutron decay pattern, <math display="block" id="M21"><mrow><msub><mi>I</mi><mrow><mtext>sample</mtext></mrow></msub><mrow><mo>(</mo><mrow><mi>i</mi><mn>,</mn><mi>j</mi></mrow><mo>)</mo></mrow><mo>=</mo><msub><mi>I</mi><mtext>o</mtext></msub><mrow><mo>(</mo><mrow><mi>i</mi><mn>,</mn><mi>j</mi></mrow><mo>)</mo></mrow><msup><mtext>e</mtext><mrow><mo>−</mo><mrow><mo>(</mo><mrow><msub><mi>Σ</mi><mtext>A</mtext></msub><mi>T</mi><mrow><mo>(</mo><mrow><mi>i</mi><mn>,</mn><mi>j</mi></mrow><mo>)</mo></mrow></mrow><mo>)</mo></mrow></mrow></msup><mn>,</mn></mrow></math> where I is the intensity of the neutron, ( i, j ) is the prime point in the image, and ∑ A is the thickness of the sample. The distribution of substances can be determined using Eq. 17. <math display="block" id="M22"><mrow><mi>T</mi><mrow><mo>(</mo><mrow><mi>i</mi><mn>,</mn><mi>j</mi></mrow><mo>)</mo></mrow><mo>=</mo><mo>−</mo><mfrac><mrow><mi>ln</mi><mrow><mo>(</mo><mrow><mfrac><mrow><msub><mi>I</mi><mrow><mtext>sample</mtext></mrow></msub></mrow><mrow><msub><mi>I</mi><mtext>o</mtext></msub></mrow></mfrac></mrow><mo>)</mo></mrow></mrow><mrow><msub><mi>Σ</mi><mtext>A</mtext></msub></mrow></mfrac></mrow></math> This is the principle behind quantitative neutron photographic analysis method. Numerous experimental studies have shown that quantitative neutron photographic analysis is influenced by three main factors. (1). Collimation Ratio The collimation ratio is the ratio of the distance from the neutron source to the sample with respect to the diameter of the output aperture. Because neutrons are electrically neutral, they cannot be made parallel by focusing the neutron beam through the electromagnetic field, and therefore it is not possible to achieve full collimation of the neutron beam using a collimator. (2). Scattered Neutrons The presence of substances with high scattering cross section (e.g., water and polyethylene) in the sample can cause a large number of scattered neutrons to reach the detector, thus distorting and blurring the neutron image. The solution is generally to increase the distance between the sample and the detector. (3). Neutron Energy Spectrum Because incident neutrons are not mono-energetic and neutrons having different energies interact with matter having different cross sections, the images and results obtained using a single cross section parameter can produce large errors. In addition, the energy spectrum of the neutrons scattered by the sample may have changed from that of the neutron source owing to attenuation during propagation. The focus of this study is the correction of neutron images by applying a correction for scattered neutrons. To further investigate the optimization of the proposed model with respect to neutron imaging, the current iterative algorithm for the quantitative analysis of scattered samples is presented in Fig. 7 Fig. 7 . In the algorithm, I su is the normalized neutron intensity, I sc is the scattered neutron distribution function, and k is the number of iterations. The criterion equation <math display="block" id="M23"><mrow><mfrac><mrow><mstyle displaystyle="true"><msub><mo>∑</mo><mi>i</mi></msub><mrow><mstyle displaystyle="true"><msub><mo>∑</mo><mi>j</mi></msub><mrow><mrow><mo>|</mo><mrow><msubsup><mi>I</mi><mrow><mtext>su</mtext></mrow><mrow><mi>k</mi><mo>+</mo><mn>1</mn></mrow></msubsup><mrow><mo>(</mo><mrow><mi>i</mi><mn>,</mn><mi>j</mi></mrow><mo>)</mo></mrow><mo>−</mo><msubsup><mi>I</mi><mrow><mtext>su</mtext></mrow><mi>k</mi></msubsup><mrow><mo>(</mo><mrow><mi>i</mi><mn>,</mn><mi>j</mi></mrow><mo>)</mo></mrow></mrow><mo>|</mo></mrow></mrow></mstyle></mrow></mstyle></mrow><mrow><mstyle displaystyle="true"><msub><mo>∑</mo><mi>i</mi></msub><mrow><mstyle displaystyle="true"><msub><mo>∑</mo><mi>j</mi></msub><mrow><msubsup><mi>I</mi><mrow><mtext>su</mtext></mrow><mi>k</mi></msubsup><mrow><mo>(</mo><mrow><mi>i</mi><mn>,</mn><mi>j</mi></mrow><mo>)</mo></mrow></mrow></mstyle></mrow></mstyle></mrow></mfrac><mo>⩽</mo><mi>ε</mi></mrow></math> is used to select a convergence coefficient ε that terminates the iterative calculations. The above criterion actually represents the relative change in the transmitted neutron intensity calculated for two adjacent iterations. As the algorithm converges, the value of the change should become increasingly smaller, gradually converging toward zero. This algorithm has been tested in a large number of neutron scattering experiments and is highly reliable. This iterative algorithm can be applied to the thin-film point scattering function model developed above. Optimization algorithm based on thin-film approximation models 2 It should be noted that in the iterative algorithm ( Fig. 7 Fig. 7 ), PScF ( r ) is calculated based on the sample thickness, which is later brought into the calculation of the neutron intensity prior to the iterative calculations. However, if we use NPScF ( r ) instead of PScF ( r ), the iterative convergence can be greatly accelerated, as NPScF ( r ) is equivalent to an approximate numerical integration of PScF ( r ). When we applied the iterative algorithm to the experiments in Section 3.1; using this substitution only three iterations were required to obtain satisfactory results, as indicated in Fig. 6 Fig. 6 . The results obtained show the reliability and validity of the iterative algorithm on the one hand, and the faster correction of imaging results using NPScF ( r ) instead of PScF ( r ) on the other. The iterative calculation using PScF ( r ) requires five iterations to converge to something closer to actuality and takes approximately three times as long to run as the calculation using NPScF ( r ). Comparing result (a) using the NPScF ( r ) function with result (b) using the conventional PScF ( r ), it can be observed that (b) produces poor imaging results, is slower to converge, and exhibits greater distortion at the edges, whereas (a) produces good imaging results and slight distortion at the edges. We followed up with similar simulations for trapezoidal objects and cones and found that using NPScF ( r ) instead of PScF ( r ) for iterative calculations not only improved the simulation, but also improved the speed of the algorithm. Moreover, in recent years, a “3D detector” has been introduced for use in neutron photographic techniques [ 18 18 , 19 19 ]. It adds 3D microstructures such as trenches or holes to a planar detector and fills them with a neutron converter, which can significantly improve the detection efficiency. For this class of detectors, the model proposed in this paper is also applicable, and together with the algorithm in Fig. 7 Fig. 7 , yields optimized results. However, the relationship between the geometric parameters of the planar detector ( Fig. 1 Fig. 1 ) must be recalculated. That is, to apply the proposed model to the 3D detector, the relationship between geometric parameters must be solved again. A mathematically feasible approach is to use the symmetry of the 3D detector microstructure to solve by parts. For example, if the 3D microstructure is axisymmetric, then it can be solved by the radius with the center of symmetry being the center of a circle, determining the geometric relationship at a fixed radius, and then solving using the method given in Sect. 2. Effects of energy spectrum hardening 2 For the thermal neutron spectra considered in this study, the cross sections interacting with matter approximately obey the 1/ v law. As neutrons penetrate the sample, low-energy neutrons decay more rapidly than high-energy neutrons, ultimately increasing the proportion of higher-energy neutrons in the neutron energy spectrum, a phenomenon known as energy spectrum hardening. The current solution to this problem is to define an "effective neutron cross section to correct for the energy spectrum. Considering that the macroscopic cross section is no longer a constant, the neutron transmittance is no longer simply exponential in relation to the thickness of the sample, even in the absence of scattered neutrons. Therefore, the energy response of the detector should also be considered. We can define an effective cross section based on the ratio of the incident neutron intensity to the transmitted neutron intensity as <math display="block" id="M24"><mrow><mfrac><mrow><msub><mi>I</mi><mrow><mtext>su</mtext></mrow></msub></mrow><mrow><msub><mi>I</mi><mtext>o</mtext></msub></mrow></mfrac><mo>=</mo><msup><mtext>e</mtext><mrow><mo>−</mo><msub><mi>Σ</mi><mrow><mtext>eff</mtext></mrow></msub><mrow><mo>(</mo><mi>T</mi><mo>)</mo></mrow><mo>⋅</mo><mi>T</mi></mrow></msup><mn>,</mn></mrow></math> where -∑ eff ( T ) represents the effective cross section for a sample of thickness T . The intensity of the neutrons detected by a detector is typically expressed by the following equation: <math display="block" id="M25"><mrow><msub><mi>I</mi><mtext>o</mtext></msub><mo>=</mo><mstyle displaystyle="true"><mrow><msub><mo>∫</mo><mi>E</mi></msub></mrow></mstyle><mtext>Φ</mtext><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><mrow><mo>(</mo><mrow><mn>1</mn><mo>−</mo><msup><mtext>e</mtext><mrow><mo>−</mo><msub><mi>Σ</mi><mrow><mi>det</mi></mrow></msub><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><mo>⋅</mo><mi>d</mi></mrow></msup></mrow><mo>)</mo></mrow><mtext>d</mtext><mi>E</mi><mn>,</mn></mrow></math> and the transmitted neutron intensity is <math display="block" id="M26"><mrow><msub><mi>I</mi><mrow><mtext>su</mtext></mrow></msub><mo>=</mo><mstyle displaystyle="true"><mrow><msub><mo>∫</mo><mi>E</mi></msub></mrow></mstyle><mtext>Φ</mtext><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><msup><mtext>e</mtext><mrow><mo>−</mo><mi>Σ</mi><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><mo>⋅</mo><mi>T</mi></mrow></msup><mrow><mo>(</mo><mrow><mn>1</mn><mo>−</mo><msup><mtext>e</mtext><mrow><mo>−</mo><msub><mi>Σ</mi><mrow><mi>det</mi></mrow></msub><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><mo>⋅</mo><mi>d</mi></mrow></msup></mrow><mo>)</mo></mrow><mtext>d</mtext><mi>E</mi><mn>.</mn></mrow></math> According to the 1/ v law, which is approximately satisfied by thermal neutrons, the conventional neutron cross section can be expressed as <math display="block" id="M27"><mrow><mi>Σ</mi><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><mo>=</mo><msub><mi>Σ</mi><mrow><mtext>th</mtext></mrow></msub><mfrac><mrow><msqrt><mrow><msub><mi>E</mi><mrow><mtext>th</mtext></mrow></msub></mrow></msqrt></mrow><mrow><msqrt><mi>E</mi></msqrt></mrow></mfrac><mn>,</mn></mrow></math> where <math id="M28"><mrow><msub><mi>E</mi><mrow><mtext>th</mtext></mrow></msub><mo>≈</mo><mn>0.025</mn><mtext>eV</mtext></mrow></math> is the energy of a thermal neutron. Considering the detector energy response as a whole, <math id="M29"><mrow><mi>ε</mi><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><mo>=</mo><mn>1</mn><mo>−</mo><msup><mtext>e</mtext><mrow><mstyle displaystyle="true"><msub><mo>∑</mo><mrow><mi>det</mi></mrow></msub><mrow><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><mo>⋅</mo><mi>d</mi></mrow></mstyle></mrow></msup><mo>≈</mo><mstyle displaystyle="true"><msub><mo>∑</mo><mrow><mi>det</mi></mrow></msub><mrow><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><mo>⋅</mo><mi>d</mi></mrow></mstyle></mrow></math> , and combining the above equations yields <math display="block" id="M30"><mrow><msub><mi>Σ</mi><mrow><mtext>eff</mtext></mrow></msub><mrow><mo>(</mo><mi>T</mi><mo>)</mo></mrow><mo>=</mo><mo>−</mo><mfrac><mn>1</mn><mi>T</mi></mfrac><mi>ln</mi><mrow><mo>(</mo><mrow><mfrac><mrow><mstyle displaystyle="true"><mrow><msub><mo>∫</mo><mi>E</mi></msub></mrow></mstyle><mtext>Φ</mtext><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><msup><mtext>e</mtext><mrow><mo>−</mo><msub><mi>Σ</mi><mrow><mtext>th</mtext></mrow></msub><mfrac><mrow><msqrt><mrow><msub><mi>E</mi><mrow><mtext>th</mtext></mrow></msub></mrow></msqrt></mrow><mrow><msqrt><mi>E</mi></msqrt></mrow></mfrac></mrow></msup><mi>ε</mi><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><mtext>d</mtext><mi>E</mi></mrow><mrow><mstyle displaystyle="true"><mrow><msub><mo>∫</mo><mi>E</mi></msub></mrow></mstyle><mtext>Φ</mtext><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><mi>ε</mi><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><mtext>d</mtext><mi>E</mi></mrow></mfrac></mrow><mo>)</mo></mrow><mn>.</mn></mrow></math> However, owing to the difficulties of performing integration operations and in obtaining relevant data, use of a simplified approximation of the above equation is desirable. The approach taken in [ 17 17 ] is to apply a linear fit, but we believe that the error in doing so would be very large. Instead, we can adopt a numerical integration approach, as described in Section 2.2 to approximate the above equation. Again applying the Newton-Cotes method, let <math display="block" id="M31"><mrow><msub><mi>Σ</mi><mrow><mtext>eff</mtext></mrow></msub><mrow><mo>(</mo><mi>T</mi><mo>)</mo></mrow><mo>≈</mo><mfrac><mn>1</mn><mi>T</mi></mfrac><mi>ln</mi><mrow><mo>(</mo><mrow><mfrac><mrow><mn>4</mn><mo>⋅</mo><mstyle displaystyle="true"><msubsup><mo>∑</mo><mrow><mi>n</mi><mo>=</mo><mn>1</mn></mrow><mi>m</mi></msubsup><mrow><mtext>Φ</mtext><mrow><mo>(</mo><mrow><mrow><mo>(</mo><mrow><mn>2</mn><mi>n</mi><mo>−</mo><mn>1</mn></mrow><mo>)</mo></mrow><mi>h</mi></mrow><mo>)</mo></mrow></mrow></mstyle><msup><mtext>e</mtext><mrow><mo>−</mo><msub><mi>Σ</mi><mrow><mtext>th</mtext></mrow></msub><mfrac><mrow><msqrt><mrow><msub><mi>E</mi><mrow><mtext>th</mtext></mrow></msub></mrow></msqrt></mrow><mrow><msqrt><mi>E</mi></msqrt></mrow></mfrac></mrow></msup><mi>ε</mi><mrow><mo>(</mo><mrow><mrow><mo>(</mo><mrow><mn>2</mn><mi>n</mi><mo>−</mo><mn>1</mn></mrow><mo>)</mo></mrow><mi>h</mi></mrow><mo>)</mo></mrow><mo>+</mo><mn>2</mn><mo>⋅</mo><mstyle displaystyle="true"><msubsup><mo>∑</mo><mrow><mi>n</mi><mo>=</mo><mn>1</mn></mrow><mrow><mi>m</mi><mo>−</mo><mn>1</mn></mrow></msubsup><mrow><mtext>Φ</mtext><mrow><mo>(</mo><mrow><mn>2</mn><mi>n</mi><mi>h</mi></mrow><mo>)</mo></mrow></mrow></mstyle><msup><mtext>e</mtext><mrow><mo>−</mo><msub><mi>Σ</mi><mrow><mtext>th</mtext></mrow></msub><mfrac><mrow><msqrt><mrow><msub><mi>E</mi><mrow><mtext>th</mtext></mrow></msub></mrow></msqrt></mrow><mrow><msqrt><mi>E</mi></msqrt></mrow></mfrac></mrow></msup><mi>ε</mi><mrow><mo>(</mo><mrow><mn>2</mn><mi>n</mi><mi>h</mi></mrow><mo>)</mo></mrow></mrow><mrow><mn>4</mn><mo>⋅</mo><mstyle displaystyle="true"><msubsup><mo>∑</mo><mrow><mi>n</mi><mo>=</mo><mn>1</mn></mrow><mi>m</mi></msubsup><mrow><mtext>Φ</mtext><mrow><mo>(</mo><mrow><mrow><mo>(</mo><mrow><mn>2</mn><mi>n</mi><mo>−</mo><mn>1</mn></mrow><mo>)</mo></mrow><mi>h</mi></mrow><mo>)</mo></mrow></mrow></mstyle><mi>ε</mi><mrow><mo>(</mo><mrow><mrow><mo>(</mo><mrow><mn>2</mn><mi>n</mi><mo>−</mo><mn>1</mn></mrow><mo>)</mo></mrow><mi>h</mi></mrow><mo>)</mo></mrow><mo>+</mo><mn>2</mn><mo>⋅</mo><mstyle displaystyle="true"><msubsup><mo>∑</mo><mrow><mi>n</mi><mo>=</mo><mn>1</mn></mrow><mrow><mi>m</mi><mo>−</mo><mn>1</mn></mrow></msubsup><mrow><mtext>Φ</mtext><mrow><mo>(</mo><mrow><mn>2</mn><mi>n</mi><mi>h</mi></mrow><mo>)</mo></mrow></mrow></mstyle><mi>ε</mi><mrow><mo>(</mo><mrow><mn>2</mn><mi>n</mi><mi>h</mi></mrow><mo>)</mo></mrow></mrow></mfrac></mrow><mo>)</mo></mrow><mn>,</mn></mrow></math> where h is the set step size, and note that Φ( E )=0 at both E =0 and <math id="M32"><mrow><mi>E</mi><mo>=</mo><mi>∞</mi></mrow></math> . With the introduction of the effective cross section, the use of an effective cross section to correct the energy spectrum can make the decomposition lines of the spectrum more visible, so that the edges of the sample or large variations in density are more clearly reflected in the image. However, it would be preferable to make certain modifications to NPScF ( r ), i.e., Eq. 15, to achieve an improved energy spectrum correction. Maxwell’s spectrum can be formally defined as Φ ( E ) = aE exp (- E/b ), and the hardening of the energy spectrum is essentially an excess retention of fast neutrons and an excess depletion of slow neutrons, and is therefore reflected in an increase in the parameter b . Considering the specific form of Maxwell’s spectrum, an increase in parameter b is equivalent to an increase in the background temperature. Thus, we can perform the following transformations to account for these factors. <math display="block" id="M33"><mrow><mtext>Φ</mtext><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><mo>→</mo><msub><mtext>Φ</mtext><mi>λ</mi></msub><mrow><mo>(</mo><mi>E</mi><mo>)</mo></mrow><mo>=</mo><mfrac><mrow><mn>2</mn><mi>π</mi><mi>n</mi></mrow><mrow><msup><mrow><mrow><mo>(</mo><mrow><mi>π</mi><mi>k</mi><mi>T</mi></mrow><mo>)</mo></mrow></mrow><mrow><mn>3</mn><mo>/</mo><mn>2</mn></mrow></msup></mrow></mfrac><msqrt><mrow><mfrac><mn>2</mn><mi>m</mi></mfrac></mrow></msqrt><mi>E</mi><msup><mtext>e</mtext><mrow><mo>−</mo><mfrac><mi>E</mi><mrow><mi>k</mi><mrow><mo>(</mo><mrow><mi>T</mi><mo>+</mo><mi>λ</mi></mrow><mo>)</mo></mrow></mrow></mfrac></mrow></msup></mrow></math> The temperature-corrected energy spectrum is then used in the step calculations. Parameter λ can then be determined by fitting the specific energy spectrum to Eq. 26. Substituting the corrected energy spectrum into Eq. 13 and replacing the mean scattering cross section with the effective scattering cross section in Eq. 25 yields <math display="block" id="M34"><mrow><mtext>TPScF</mtext><mrow><mo>(</mo><mi>r</mi><mo>)</mo></mrow><mo>=</mo><mfrac><mrow><mstyle displaystyle="true"><msubsup><mo>∑</mo><mrow><mi>i</mi><mo>=</mo><mn>0</mn></mrow><mrow><msub><mi>E</mi><mtext>m</mtext></msub></mrow></msubsup><mrow><msub><mtext>Φ</mtext><mi>λ</mi></msub></mrow></mstyle><mrow><mo>(</mo><mrow><msub><mi>E</mi><mi>i</mi></msub></mrow><mo>)</mo></mrow><msub><mi>Σ</mi><mrow><mtext>det</mtext></mrow></msub><mrow><mo>(</mo><mrow><msub><mi>E</mi><mi>i</mi></msub></mrow><mo>)</mo></mrow><msub><mi>Σ</mi><mrow><mtext>eff</mtext></mrow></msub><mrow><mo>(</mo><mrow><msub><mi>E</mi><mi>i</mi></msub></mrow><mo>)</mo></mrow><mi>F</mi><mrow><mo>(</mo><mrow><msub><mi>E</mi><mi>i</mi></msub><mn>,</mn><mi>r</mi></mrow><mo>)</mo></mrow></mrow><mrow><mstyle displaystyle="true"><msubsup><mo>∑</mo><mrow><mi>i</mi><mo>=</mo><mn>0</mn></mrow><mrow><msub><mi>E</mi><mtext>m</mtext></msub></mrow></msubsup><mrow><msub><mtext>Φ</mtext><mi>λ</mi></msub></mrow></mstyle><mrow><mo>(</mo><mrow><msub><mi>E</mi><mi>i</mi></msub></mrow><mo>)</mo></mrow><msub><mi>Σ</mi><mrow><mtext>det</mtext></mrow></msub><mrow><mo>(</mo><mrow><msub><mi>E</mi><mi>i</mi></msub></mrow><mo>)</mo></mrow></mrow></mfrac><mn>.</mn></mrow></math> Subsequently, we can correct for the error by applying error analysis, as discussed in Section 3.2, by choosing a suitable exponential function based on the symmetry of the spatial distribution. Conclusion 1 In this study, an approximate point scattered function is investigated for the correction of neutron photographic imaging. The modeling is carried out in Section 2.1, whereas Section 2.2 explains how the algorithm is used to carry out the calculations. In the simulation experiment described in Section 3, the results shown in Fig. 3 Fig. 3 demonstrate the validity of the model. The model was then corrected for the relative error results shown in Fig. 4 Fig. 4 . This correction is reflected in the addition of an exponential function to adjust for the error (see Eq.15), and the parameters are determined using the fitting method to finally obtain Eq.16 to produce satisfactory simulation results. Fig. 5 Fig. 5 shows the validity of the correction, indicating a mean relative error of 2.5% for the entire neutron converter. In Sect. 4, we apply the model to the popular quantitative neutron photographic analysis algorithm and find that replacing the conventional point scattered function with the proposed approximate point scattered function results in better neutron imaging, faster convergence of the algorithm, and a significant reduction in the number of iterations. In addition, we investigated the effect of energy spectrum hardening on neutron photographic imaging, concluding that energy spectrum hardening is equivalent to increasing the background temperature in Maxwell’s spectrum, and therefore a temperature correction parameter can be introduced (Eq. 26). After replacing the mean scattering cross section with the effective scattering cross section and including the temperature correction, a new point scattering function is given by Eq. 27. Finally, the error analysis method applicable to the proposed model is discussed. References C.M. Brenner , S.R. Mirfayzi , D.R. Rusby et al. , Laser-driven x-ray and neutron source development for industrial applications of plasma accelerators . Plasma Phys. Contr. F. 58 ( 1 ), 014039 ( 2015 ). doi: 10.1088/0741-3335/58/1/014039 http://doi.org/10.1088/0741-3335/58/1/014039 Laser-driven x-ray and neutron source development for industrial applications of plasma accelerators A.E. Ersundu , M. Büyükyıldız , M.Ç. Ersundu et al. , The heavy metal oxide glasses within the WO 3 -MoO 3 -TeO 2 system to investigate the shielding properties of radiation applications . Prog. Nucl. Energ. 104 , 280 - 287 ( 2018 ). doi: 10.1016/j.pnucene.2017.10.008 http://doi.org/10.1016/j.pnucene.2017.10.008 The heavy metal oxide glasses within the WO3-MoO3-TeO2 system to investigate the shielding properties of radiation applications B. Wang , S. Zhong , T.L. Lee et al. , Non-destructive testing and evaluation of composite materials/structures: A state-of-the-art review . Adv. Mech. Eng. 12 ( 4 ), 1687 ( 2020 ). doi: 10.1177/1687814020913761 http://doi.org/10.1177/1687814020913761 Non-destructive testing and evaluation of composite materials/structures: A state-of-the-art review S.Q. Liu , T. Bücherl , H. Li et al. , Corrections on energy spectrum and scatterings for fast neutron radiography at NECTAR facility . Chinese Phys. C 37 ( 11 ), 118201 ( 2013 ). doi: 10.1088/1674-1137/37/11/118201 http://doi.org/10.1088/1674-1137/37/11/118201 Corrections on energy spectrum and scatterings for fast neutron radiography at NECTAR facility N. Hachouf , F. Kharfi , A. Boucenna , Characterization and MCNP simulation of neutron energy spectrum shift after transmission through strong absorbing materials and its impact on tomography reconstructed image . Appl. Radiat. Isotopes 70 ( 10 ), 2355 - 2361 ( 2012 ). doi: 10.1016/j.apradiso.2012.06.017 http://doi.org/10.1016/j.apradiso.2012.06.017 Characterization and MCNP simulation of neutron energy spectrum shift after transmission through strong absorbing materials and its impact on tomography reconstructed image R. Hassanein , E. Lehmann , P. Vontobel , Methods of scattering corrections for quantitative neutron radiography . Nucl. Instrum. Meth. A 542 ( 1-3 ), 353 - 360 ( 2005 ). doi: 10.1016/j.nima.2005.01.161 http://doi.org/10.1016/j.nima.2005.01.161 Methods of scattering corrections for quantitative neutron radiography M. Raventós , E.H. Lehmann , M. Boin et al. , A Monte Carlo approach for scattering correction towards quantitative neutron imaging of polycrystals . J. Appl. Crystallogr. 51 ( 2 ), 386 - 394 ( 2018 ). doi: 10.1107/s1600576718001607 http://doi.org/10.1107/s1600576718001607 A Monte Carlo approach for scattering correction towards quantitative neutron imaging of polycrystals N. Kardjilov , F. De Beer , R. Hassanein et al. , Scattering corrections in neutron radiography using point scattered functions . Nucl. Instrum. Meth. A 542 ( 1-3 ), 336 - 341 ( 2005 ). doi: 10.1016/j.nima.2005.01.159 http://doi.org/10.1016/j.nima.2005.01.159 Scattering corrections in neutron radiography using point scattered functions M.H. Hassan , Point Scattered Function (PScF) for fast neutron radiography . Nucl. Instrum. Meth. B 267 ( 15 ), 2545 - 2549 ( 2009 ). doi: 10.1016/j.nimb.2009.05.062 http://doi.org/10.1016/j.nimb.2009.05.062 Point Scattered Function (PScF) for fast neutron radiography X. Xu , X. Xie , H. He et al. , Imaging objects through scattering layers and around corners by retrieval of the scattered point spread function . Opt. Express 25 ( 26 ), 32829 - 32840 ( 2017 ). doi: 10.1364/oe.25.032829 http://doi.org/10.1364/oe.25.032829 Imaging objects through scattering layers and around corners by retrieval of the scattered point spread function L. Fernandez‐Izquierdo , M.G. Reyes‐Banda , X. Mathew et al. , Cesium lead bromide (CsPbBr 3 ) thin‐film‐based solid‐state neutron detector developed by a solution‐free sublimation process . Adv. Mater. Technol-US 5 ( 12 ), 2000534 ( 2020 ). doi: 10.1002/admt.202000534 http://doi.org/10.1002/admt.202000534 Cesium lead bromide (CsPbBr3) thin‐film‐based solid‐state neutron detector developed by a solution‐free sublimation process M.T. Tchouaso , Efficiency limitation of thin film coated planar semiconductor thermal neutron detector . Appl. Radiat. Isotopes 173 , 109716 ( 2021 ). doi: 10.1016/j.apradiso.2021.109716 http://doi.org/10.1016/j.apradiso.2021.109716 Efficiency limitation of thin film coated planar semiconductor thermal neutron detector C. Avila-Avendano , I. Mejia , R. Garcia-Lozano et al. , Integrated thin-film radiation detectors and in-pixel amplification . IEEE T. Electron Dev. 65 ( 9 ), 3809 - 3815 ( 2018 ). doi: 10.1109/ted.2018.2859778 http://doi.org/10.1109/ted.2018.2859778 Integrated thin-film radiation detectors and in-pixel amplification T. Sauer , Numerical Analysis . Pearson 231 ( 2006 ). Numerical Analysis Hassanein, René Karim . Correction methods for the quantitative evaluation of thermal neutron tomography . ETH Zurich ( 2006 ). Correction methods for the quantitative evaluation of thermal neutron tomography G. Akkaya , The calculation of self-shielding correction factors for large samples in 241Am-Be isotopic neutron source . Appl. Radiat. Isotopes 179 , 109990 ( 2022 ). doi: 10.1016/j.apradiso.2021.109990 http://doi.org/10.1016/j.apradiso.2021.109990 The calculation of self-shielding correction factors for large samples in 241Am-Be isotopic neutron source S.C. Capelli , G. Romanelli , An effective hydrogen scattering cross section for time-of-flight neutron experiments with simple organic molecules . J. Appl. Crystallogr. 52 ( 5 ), 1233 - 1237 ( 2019 ). doi: 10.1107/s1600576719011592 http://doi.org/10.1107/s1600576719011592 An effective hydrogen scattering cross section for time-of-flight neutron experiments with simple organic molecules J. Uher , C. Fröjdh , J. Jakůbek et al. , Characterization of 3D thermal neutron semiconductor detectors . Nucl. Instrum. Meth. A 576 ( 1 ), 32 - 37 ( 2007 ). doi: 10.1016/j.nima.2007.01.115 http://doi.org/10.1016/j.nima.2007.01.115 Characterization of 3D thermal neutron semiconductor detectors S. Douglas , L. Steven , J. Bellinger et al. , Present status of microstructured semiconductor neutron detectors . J. Cryst. Growth 99 ( 379 ), ( 2013 ). doi: 10.1016/j.jcrysgro.2012.10.061 http://doi.org/10.1016/j.jcrysgro.2012.10.061 Present status of microstructured semiconductor neutron detectors C.B. Lu , J. Bao , Y. Huang et al. , Contrast sensitivity in 14 MeV fast neutron radiography . Nucl. Sci. Tech. 28 ( 6 ), 78 ( 2017 ). doi: 10.1007/s41365-017-0228-5 http://doi.org/10.1007/s41365-017-0228-5 Contrast sensitivity in 14 MeV fast neutron radiography M.S. Beni , D. Krstic , D. Nikezic et al. , Studies on unfolding energy spectra of neutrons using maximum-likelihood expectation-maximization method . Nucl. Sci. Tech. 30 ( 9 ), 134 ( 2019 ). doi: 10.1007/s41365-019-0662-7 http://doi.org/10.1007/s41365-019-0662-7 Studies on unfolding energy spectra of neutrons using maximum-likelihood expectation-maximization method All authors contributed to the study conception and design. Material preparation, data collection and analysis were performed by Lin-Feng Ye, Wei-Jun Jiang and De-Hong Gao. The first draft of the manuscript was written by Jia-Yu Ni and Jun Qin and all authors commented on previous versions of the manuscript. All authors read and approved the final manuscript. 10.1007/s41365-022-01094-y This work was supported by the College Foundation Project, the College of Engineering and Technology of Chengdu University of Technology (No. C122018029). 2022-06-03 2022-07-23 2022-07-27 2022-09-08 2022-10-24 © China Science Publishing & Media Ltd. (Science Press), Shanghai Institute of Applied Physics, the Chinese Academy of Sciences, ChineseNuclear Society and Springer Nature Singapore Pte Ltd. 2022 This is a post-peer-review, pre-copyedit version of an article published in Nuclear Science and Techniques. The final authenticated version is available online at: http://dx.doi.org/10.1007/s41365-022-01094-y http://dx.doi.org/10.1007/s41365-022-01094-y . Springer Nature or its licensor holds exclusive rights to this article under a publishing agreement with the author(s) or other rightsholder(s); author self-archiving of the accepted manuscript version of this article is solely governed by the terms of such publishing agreement and applicable law. 2022 Nuclear Science and Techniques 09 33