Chasing quantitative biases in neutron imaging with scintillator-camera detectors : a practical method with black body grids

We propose a method for improving the quantification of neutron imaging measurements with scintillator-camera based detectors by correcting for systematic biases introduced by scattered neutrons and other sources such as light reflections in the detector system. This method is fully experimental, using reference measurements with a grid of small black bodies (BB) to measure the bias contributions directly. Using two test samples, one made of lead alloy and having a moderate (20%) neutron transmission and one made of stainless-steel and having a very low (1%) transmission, we evaluated the improvement brought by this method in reducing both the average quantification bias and the uncertainty around this average bias after tomographic reconstruction. The results show that a reduction of the quantification bias of up to one order of magnitude can be obtained. For moderately transparent samples, little sensitivity is observed to the parameters used for the correction. For the more challenging sample with very low transmission, a correct placement of the BB grid is of utmost importance for a successful correction. © 2018 Optical Society of America under the terms of the OSA Open Access Publishing Agreement OCIS codes: (120.0120) Instrumentation, measurement, and metrology; (100.6950) Tomographic image processing; (110.0110) Imaging systems; (120.3940) Metrology; (290.0290) Scattering. References and links 1. M. Strobl, I. Manke, N. Kardjilov, A. Hilger, M. Dawson, and J. Banhart, “Advances in neutron radiography and tomography,” J. Phys. D Appl. Phys. 42(24), 243001 (2009). 2. D. Mannes, F. Schmid, J. Frey, K. Schmidt-Ott, and E. Lehmann, “Combined Neutron and X-ray Imaging for Non-invasive Investigations of Cultural Heritage Objects,” Phys. Procedia 69, 653–660 (2015). 3. C. Gruenzweig, D. Mannes, F. Schmid, and R. Rule “Neutron Imaging: A Non-Destructive Testing Method to Investigate Canned Exhaust After-Treatment System Components for the Three Dimensional Soot, Ash, Urea and Coating Distributions,” 2016. 4. C. Grünzweig, D. Mannes, A. Kaestner, F. Schmid, P. Vontobel, J. Hovind, S. Hartmann, S. Peetermans, and E. Lehmann, “Progress in Industrial Applications using Modern Neutron Imaging Techniques,” Phys. Procedia 43, 231–242 (2013). 5. A. Griesche, E. Dabah, T. Kannengiesser, A. Hilger, N. Kardjilov, I. Manke, and B. Schillinger, “Measuring Hydrogen Distributions in Iron and Steel Using Neutrons,” Phys. Procedia 69, 445–450 (2015). 6. M. Große, M. Steinbrück, J. Stuckert, A. Kastner, and B. Schillinger, “Application of neutron radiography to study material processes during hypothetical severe accidents in nuclear reactors,” J. Mater. Sci. 47(18), 6505– 6512 (2012). 7. P. Boillat, P. Oberholzer, A. Kaestner, R. Siegrist, E. H. Lehmann, G. G. Scherer, and A. Wokaun, “Impact of Water on PEFC Performance Evaluated by Neutron Imaging Combined with Pulsed Helox Operation,” J. Electrochem. Soc. 159(7), F210–F218 (2012). 8. J. P. Owejan, T. A. Trabold, J. J. Gagliardo, D. L. Jacobson, R. N. Carter, D. S. Hussey, and M. Arif, “Voltage instability in a simulated fuel cell stack correlated to cathode water accumulation,” J. Power Sources 171(2), 626–633 (2007). 9. P. Stahl, J. Biesdorf, P. Boillat, J. Kraft, and K. A. Friedrich, “Water Distribution Analysis in the Outer Perimeter Region of Technical PEFC Based on Neutron Radiography,” J. Electrochem. Soc. 162(7), F677–F685 (2015). Vol. 26, No. 12 | 11 Jun 2018 | OPTICS EXPRESS 15769 #327723 https://doi.org/10.1364/OE.26.015769 Journal © 2018 Received 11 Apr 2018; revised 24 May 2018; accepted 25 May 2018; published 7 Jun 2018 10. J. Seweryn, J. Biesdorf, T. J. Schmidt, and P. Boillat, “Communication—Neutron Radiography of the Water/Gas Distribution in the Porous Layers of an Operating Electrolyser,” J. Electrochem. Soc. 163(11), F3009–F3011 (2016). 11. S. Lal, L. D. Poulikakos, I. Jerjen, P. Vontobel, M. N. Partl, D. Derome, and J. Carmeliet, “Investigation of Gravity-Driven Drainage and Forced Convective Drying in a Macroporous Medium Using Neutron Radiography,” Transp. Porous Media 118(1), 119–142 (2017). 12. J. D. Schaap, P. Lehmann, A. Kaestner, P. Vontobel, R. Hassanein, G. Frei, G. H. de Rooij, E. Lehmann, and H. Flühler, “Measuring the effect of structural connectivity on the water dynamics in heterogeneous porous media using speedy neutron tomography,” Adv. Water Resour. 31(9), 1233–1241 (2008). 13. A. Kaestner, R. Hassanein, P. Vontobel, P. Lehmann, J. Schaap, E. Lehmann, and H. Flühler, “Mapping the 3D water dynamics in heterogeneous sands using thermal neutrons,” Chem. Eng. J. 130(2-3), 79–85 (2007). 14. P. Benard, E. Kroener, P. Vontobel, A. Kaestner, and A. Carminati, “Water percolation through the root-soil interface,” Adv. Water Resour. 95, 190–198 (2016). 15. M. Zarebanadkouki, A. Carminati, A. Kaestner, D. Mannes, M. Morgano, S. Peetermans, E. Lehmann, and P. Trtik, “On-the-fly Neutron Tomography of Water Transport into Lupine Roots,” Phys. Procedia 69, 292–298 (2015). 16. H. Justnes, K. Bryhn-Ingebrigtsen, and G. O. Rosvold, “Neutron radiography: an excellent method of measuring water penetration and moisture distribution in cementitious materials,” Adv. Cement Res. 6(22), 67–72 (1994). 17. P. Trtik, B. Münch, W. J. Weiss, A. Kaestner, I. Jerjen, L. Josic, E. Lehmann, and P. Lura, “Release of internal curing water from lightweight aggregates in cement paste investigated by neutron and X-ray tomography,” Nucl. Instrum. Methods Phys. Res. A 651(1), 244–249 (2011). 18. M. Sedighi-Gilani, M. Griffa, D. Mannes, E. Lehmann, J. Carmeliet, and D. Derome, “Visualization and quantification of liquid water transport in softwood by means of neutron radiography,” Int. J. Heat Mass Transfer 55(21-22), 6211–6221 (2012). 19. R. Hassanein, “Correction methods for the quantitative evaluation of thermal neutron tomography,” (ETH Zürich, 2006), p. 149. 20. D. S. Hussey, K. J. Coakley, E. Baltic, and D. L. Jacobson, “Improving quantitative neutron radiography through image restoration,” Nucl. Instrum. Methods Phys. Res. A 729, 316–321 (2013). 21. R. Hassanein, E. Lehmann, and P. Vontobel, “Methods of scattering corrections for quantitative neutron radiography,” Nucl. Instrum. Methods Phys. Res. A 542(1-3), 353–360 (2005). 22. N. Kardjilov, F. de Beer, R. Hassanein, E. Lehmann, and P. Vontobel, “Scattering corrections in neutron radiography using point scattered functions,” Nucl. Instrum. Methods Phys. Res. A 542(1-3), 336–341 (2005). 23. M. H. Hassan, “Point Scattered Function (PScF) for fast neutron radiography,” Nucl. Instrum. Methods Phys. Res. B 267(15), 2545–2549 (2009). 24. S. Liu, T. Bücherl, Y. Zou, S. Wang, Y. Lu, and Z. Guo, “Study on scattering correction in fast neutron tomography at NECTAR facility,” Sci. China Phys. Mech. Astron. 57(2), 244–250 (2014). 25. A. S. Tremsin, N. Kardjilov, M. Dawson, M. Strobl, I. Manke, J. B. McPhate, J. V. Vallerga, O. H. W. Siegmund, and W. B. Feller, “Scatter rejection in quantitative thermal and cold neutron imaging,” Nucl. Instrum. Methods Phys. Res. A 651(1), 145–148 (2011). 26. M. Raventos, R. P. Harti, E. Lehmann, and C. Grünzweig, “A method for neutron scattering quantification and correction applied to neutron imaging,” Phys. Procedia 88, 275–281 (2017). 27. R. Ning, X. Tang, and D. Conover, “X-ray scatter correction algorithm for cone beam CT imaging,” Med. Phys. 31(5), 1195–1202 (2004). 28. E. H. Lehmann, P. Vontobel, and L. Wiezel, “Properties of the radiography facility NEUTRA at SINQ and its potential for use as European reference facility,” Nondestruct. Test. Eval. 16(2-6), 191–202 (2001). 29. R. Franke, “Smooth interpolation of scattered data by local thin plate splines,” Comput. Math. Appl. 8(4), 273–


Introduction
Neutron imaging is a powerful non-destructive and non-invasive analysis method complementary to the well-known X-ray imaging [1].Typically, neutron imaging is used in cases where X-rays cannot be applied because of the need to penetrate through an important thickness of high-Z materials, to visualize low-Z materials, or a combination of both.The range of demonstrated uses include, among others, the study of archeological objects [2], the analysis of soot and ash accumulation in particulate filters [3,4], the observation of hydrogen in metals [5,6].The good penetration of neutrons through a variety of materials also allow for fully non-invasive in situ studies, such as the measurement of the water distribution in operating fuel cells [7][8][9] or electrolysers [10], the analyses of water flow in porous media [11][12][13] and of the water uptake by plant roots [14,15], concrete [16,17] or different types of wood [18].
Although quantification in imaging can refer to a broad variety of aspects and analysesin some cases even qualitative results suffice to answer specific questions -, the need for even more precise quantification of the measured attenuation is growing with a growing number of material science applications.Unfortunately, to date, there is still a lack of knowledge about the systematic biases affecting such quantification and an assessment of the measurement accuracy based on statistical deviations only is not satisfactory.Systematic biases can originate from the detection system itself, for example from neutrons impinging either on the detector casing or internal elements such as a mirror which are scattered back towards the neutron sensitive element (e.g.scintillator screen).This has been evidenced by Hassanein [19], who observed that the intensity measured behind a layer fully opaque to neutrons was not zero, and that its intensity depended on the fraction of the field of view which was left free, i.e. the size of the neutron beam.He proposed a correction method where he measured the intensity of this so called background scattering by placing a neutron-opaque object, called 'black body' (BB) in place of the sample, and subtracted this background contribution.One drawback of this approach is that if the sample covers a large fraction of the neutron beam -which is the case in any well optimized setup -the black body used for measurement strongly impacts the background to be measured by blocking most of the neutron beam.More recently, Hussey et al. [20] described a systematic bias introduced by the "long range tail" of the detector transfer function.They were not very specific concerning the physical basis for this event, which might be the same as the one observed by Hassanein.After measuring the transfer function of the detector, they correct for this long range contribution by deconvolution.Although very interesting, this approach suffers from the drawback -as also mentioned by the authors -that the deconvolution correction is based on the assumption of a spatially invariant transfer function, which can be far from the reality for certain detectors, in particular for the popular setup using a mirror at 45° from the detector surface.A further contribution to systematic errors is introduced by the detection of neutrons scattered by the sample itself.To correct for this sample scattering, a method based on a so called point scattered function (PScF) has been proposed [21,22].The PScF is obtained by Monte Carlo simulation on slabs of a given thickness with infinite size in the directions perpendicular to the beam, and the image is restored in an iterative manner by subtracting the calculated scattering component from the measured image and subsequently recalculating the PScF until convergence is found.Although this method demonstrated promising results and was also used at other facilities for correcting tomographies with fast neutrons [23,24], important drawbacks have to be mentioned.First, it is not a fully experimental method, as the need to simulate the scattered components introduced an a priori on the sample behavior.Second, the sample scattering behavior is assumed to be the same as a slab of the same thickness, but with infinite size in the directions perpendicular to the beam.The real scattering behavior might strongly differ from this assumption, in particular for thin samples (e.g.fuel cells) with their longest dimension aligned in the beam direction.Finally, the use of Monte Carlo simulations heavily increases the complexity of the post processing.In this context, a possibility to avoid the bias introduced by scattered neutrons in a fully experimental approach is highly appealing.Tremsin et al. [25] proposed a method where the scattered neutrons are filtered out using neutron absorbing micro-channel plates (MCP) which only transmit neutrons parallel to the beam axis.The thickness of the MCP (1 mm) allows placing them between the detector and the sample without a significant increase of the sample-detector distance.However, the ratio of the channels area to the total area is approximately 50%, and the sample has thus to be exposed twice longer to obtain the same statistics.Another approach consisting in performing two measurements at different distances to evaluate the contribution of the neutrons scattered by the sample was recently proposed by Raventos et al. [26].This approach demonstrated a clear improvement of quantification, though the approximation of the scattered neutron intensity by a Gaussian shape might be suitable only for a limited set of samples.Here, we propose the use of a robust method to measure and correct for both effects (background scattering and sample scattering) in a single step, with a small increase (4%) of the acquisition time.While similar previous research dealt with scattering correction for Xray CT cone beam imaging [27] with flat panel detectors, we propose here a quantitative method suited for neutron imaging and more in general for imaging systems with scintillatorcamera based detectors.The strength of this method is that, unlike corrections based on simulations, it does not require any a priori knowledge of the sample composition and structure.The method uses reference data obtained with a set of small black bodies (BB) arranged in a grid pattern.These BB have a sufficient thickness of a strong neutron absorber to be fully opaque to neutrons and a diameter of only a tiny fraction of the field of view.Thus, the signal measured behind them (see Fig. 1) corresponds to the sum of all contributions other than the direct transmitted beam.Unlike the approach using a large black body covering the whole sample as proposed by Hassanein [19], our measurement setup results in a rather minute attenuation of the beam (approx.3%) during the scattering measurement -avoiding the corresponding measurement errors -and allows measuring the neutrons scattered by the sample itself.The improvements brought by this correction to test samples will be shown here, and the impact of different implementations of the correction strategy will be discussed.The implication of these results for X-ray imaging will also be mentioned.

Test samples
Metal blocks were used as test samples, the first one made of lead alloy (96% Pb, 4% Sb) and having x-, y-, z-dimensions of 50 x 15 x 100 mm 3 , the second one made of stainless-steel (EN 1.4404 / AISI 316L) and having dimensions of 40 x 10 x 150 mm 3 .

Measurement setup
Each sample was placed on a rotation table with the rotation axis parallel to the longest dimension of the sample.In all results, the 0° orientation refers to the position where the second longest dimension (50 mm or 40 mm) is parallel to the detector surface.The distance between the rotation axis and the detector surface was set to 32 mm.The samples were imaged at 625 different orientations homogeneously distributed over a full rotation (360°).The experiments were conducted at the NEUTRA beam line [28] of the Paul Scherrer Institut (PSI).The neutron image formed behind the sample was converted to light using a 100 µm thick 6 LiF/ZnS scintillator screen, and the image was captured using a digital camera (Andor Neo sCMOS, 2160 x 2560 pixels) with an exposure time of 8 s.The field of view was 145 x 172 mm 2 , resulting in a pixel size of 67 µm.

Reference attenuation coefficient measurements
The energy averaged attenuation coefficients can in principle be computed provided that the exact composition and crystalline structure of the sample is known.However, we deemed more practical to measure them directly for the NEUTRA spectrum using the reference measurements described as follows.In these measurements, a pencil beam (7 x 10 mm 2 ) was used so that the measured sample would completely cover it, and the transmission was measured as the ratio between the beam intensity with and without the sample, from the average of 20 exposures with a duration of 10 s each.The neutron imaging detector was placed 5 m away from the measured sample, so that the contribution from scattered neutrons was clearly negligible.An evacuated flight tube between the sample and detector allowed to minimize the losses over this distance.All scattering contributions occurring within the detector being proportional to the beam intensity, they were considered not to impact the transmission measurement.Based on the measured transmission and samples thicknesses, the obtained reference attenuation coefficients were 0.335 ± 0.01 cm −1 for the lead sample and 1.08 ± 0.025 cm −1 for the stainless-steel sample.The uncertainties on the attenuation coefficients were obtained from an evaluated ± 1% uncertainty in the intensity measurements (derived from the standard deviation of the single exposures) and an assumed ± 0.1 mm uncertainty in the sample thickness measurement.

BB grid measurements
For the measurement of scattered neutrons and other background contributions, we used an arrangement of 5x5 cylindrical black bodies (BB) made of 10 B 4 C (2.5 mm diameter, 3 mm length) with a center-to-center spacing of 25 mm -a period which was chosen as a compromise providing a reasonable spatial resolution for the background measurement while having only a small portion of the beam covered by absorbers and thus minimize the impact of the grid on the intensity of the scattered neutrons.A grid shaped frame cut out of a 10 mm thick aluminum plate was used to hold the neutron absorbers.The neutron absorbers cover 0.8% of the area and the holding frame has an average neutron attenuation of 2%, resulting in an average transmission of 97.2%.The grid was placed at a distance of 160 mm from the detector to avoid any disturbance from the neutrons scattered by the grid itself.Taking into account the beam divergence (L/D ~ = 350) in this experiment, the geometric blurring of the BB absorbers was approximately 0.5 mm, and smaller absorbers (e.g.1.5 mm) could in principle also be used.However, the BB grid was designed as a generic grid which can also be used with beams of higher divergence.For the correction process defined as the base case (a suggestion of the standard correction parameters), the grid was positioned with a vertical line of BB coinciding with the sample rotation axis, and measurements with the sample in 25 different orientations distributed over 360° were acquired.These 25 measurements were conducted as a single block of measurements rather than "interleaved" in the main tomography sequence.In this case, measurements of the background were obtained at the BB positions only (with a spacing of 25 mm) and for selected sample orientations only.The values for missing positions and orientations were obtained by interpolation (see "image processing").The following variations around this base case were realized: 1) Displacing the BB grid in two directions perpendicular to the beam axis, with 4 horizontal and 4 vertical steps resulting in 16 different positions.This allowed a measurement of the background with a finer spacing of 6.25 mm and was done with the sample at 25 different orientations.This variation is referred to as "fine grid" measurement later on, opposed to the "coarse grid" of the base case.
2) Measuring with one single grid position (as in the base case), but with all 625 different sample orientations.
3) Measuring with the grid shifted horizontally by 12.5 mm compared to the base case, so that the sample rotation axis was between two vertical lines of BB.This was done with the sample in 25 different orientations.This variation is referred to as "shifted grid" measurement later on, opposed to the "centered grid" of the base case.For all cases, measurements of the BB grid without the sample (called "open beam with BB") were performed.The same BB grid positions were used for these measurements as for the respective measurements with the sample.
To avoid too much influence of scintillator afterglow effects, a delay of 3 minutes was observed after each movement of the BB grid.

PSF measurements
To measure the point spread function (PSF) of the detection system, measurements with a very small "pencil" beam (7 x 5 mm 2 ) and with a medium sized (25 x 25 mm 2 ) beam were recorded.Shielding blocks fully opaque to neutrons (5 cm thick polyethylene with 5% B 4 C) were used to limit the beam and an additional measurement with the shielding blocks fully blocking the beam was used as reference background.These measurements do not focus on the short range (< 1 mm) part of the PSF corresponding to the detector intrinsic resolution, but on the middle range (< 10 mm) and the long range tails of the PSF.

Image processing
The images were filtered to remove outliers corresponding to high energy photons directly hitting the detector, and were corrected for the signal offset in the camera (due to the dark current and readout shift) as well as for beam intensity fluctuations.The middle range tail of the PSF (see below) was corrected for by deconvolution using direct inversion in the Fourier domain.This correction is named "middle range deconvolution" and is applied in the base case.In specific cases identified throughout the manuscript, this deconvolution was omitted for comparison purposes.The correction steps described above were applied to each single image including the sample projections, all images with the BB grid and the open beam images (measurement of the beam and detector sensitivity spatial distribution, also called "flat field" measurement).The values measured in the position of each grid dot are considered the local background caused by the outlined biasing scattering effects.These values were obtained by taking the median value (in order to reject any remaining outlier) of the pixels contained in a circle in the center of the dot and having half of its diameter.The background values for the whole field of view, i.e. including the intermediate positions between the BB grid dots were obtained by interpolation using a thin plate splines (TPS) algorithm [29].The processing data flow is illustrated in Fig. 2.

Tomographic reconstruction
The reconstruction was performed using the PSI in-house developed MuhRec software [30], with an additional module allowing for the subtraction of externally computed scattered background images from each sample projection image, as well as from the open beam image.The voxel size of the reconstructed image is the same as the pixel size in the projections, i.e. 67 µm.For comparison purposes, nine different combinations of reconstruction parameters were used as described in Table 1.The "base case" was defined as being the case no.6 because it corresponds to the configuration which will likely be used routinely, owing to its low overhead on measurement time.

No
Grid positions Projections count with BB (1)   Centered/shifted (2)  Use of deconvolution Overhead (3)     (1) Projections count for the background measurements with the sample and the BB grid.The projection count for the measurement with the sample only was 625 in all cases. (2)"Centered" means that one vertical line of dots coincides with the rotation axis, "Shifted" means that the rotation axis is in the middle position between two vertical lines of dots. (3)Additional measurement time for all BB measurements (with and without the sample), as a fraction of the time used for the set of 625 sample projections and 10 open beam measurements.The additional time of 3 minutes used after removing the BB to avoid the afterglow is not included here, as it does not scale with the exposure time.
The obtained 3D data set was further processed to obtain the average value of the measured attenuation coefficient within the sample volume, as well as the spread of values around this average due to systematic errors.This processing was done as follows: A threshold of half of the theoretical attenuation coefficient was applied to identify the voxels within the sample, and the holes due to outliers were removed using a closing operation in 3 dimensions.All voxels within a distance of 5 pixels of the border were discarded to eliminate the impact of blurring.The average of the remaining voxel values was computed as the average measured attenuation coefficient.Then cubes of 8x8x8 voxels were averaged together, keeping only the cubes fully inside the sample.The spread due to systematic error was obtained after discarding the 0.5% lowest and 0.5% highest values of these 8x8x8 cubes, corresponding to a 99% confidence interval.The reason for using average values of 8x8x8 voxel cubes instead of the single voxel values is the need to reduce sufficiently the impact of the statistical spread (due to noise) to identify the effect of systematic errors.The analysis of the pencil beam measurements indicates that different physical phenomena contribute to the scattered background.On the long range (Fig. 3(a)), a strong vertical distribution is observed with a higher intensity in the upper part of the image.We attribute this specific feature as follows: part of the light produced by the scintillator screen in one point is reflected by the mirror to illuminate another position of the scintillator (Fig. 3(e)).The relatively sharp transition between the upper and lower part of the image corresponds to the limit between the regions where this is or isn't geometrically possible.The intensity measured in the lower part of the image is attributed to neutrons which are back-scattered by the mirror (Fig. 3(d)).The small pencil beam measurement (Fig. 3(b)) allows emphasizing the middle range (0-20 mm) contribution to the scattered background.This contribution -unlike the long-range interactions described above -can be well represented by a convolution operation, with the best results obtained using a Lorentzian PSF (Fig. 3(c)).We attribute this local effect to the scattering of neutrons by the aluminum substrate of the scintillator screen (Fig. 3(f)), and by the polymeric binder in the scintillator screen.As mentioned in the experimental section, we decided to correct this specific effect by applying deconvolution with the PSF calibrated on the small pencil beam measurement.The measurement of the scattered background using the BB grid with the open beam (Fig. 4) gives results which are consistent with the pencil beam measurements.A background of 2% of the nominal beam intensity is measured in the lower part of the image not subject to the light reflection effect described above.A value of 4% is measured in the upper half of the image subject to this effect.Finally, the absence of deconvolution correction leads to a measured background further increased by 2%.We conclude that the three physical processes described above (neutron scattering by the mirror, neutron scattering by the scintillator, and light reflection off the mirror) have approximately equal contributions to the measured background.For the sake of simplicity, we will still refer to these contributions as the scattered background, although a part of it is caused by light reflection.As shown in Fig. 4(c), the measurements with one single position of the BB grid -which is much easier to implement -nearly perfectly catches the spatial profile of the background, with only a small discrepancy in the lowest 10 mm of the image.We will further on refer to the use of a "coarse grid" when only one position of the grid is used, and to a "fine grid" when the information from all 16 positions of the grid are used.When placing a sample in the beam, the total scattering component is affected by two effects.First, it is increased by the addition of neutrons scattered by the sample and impinging on the detector.Second, it is decreased because the attenuation of the beam proportionally reduces the effects described in the previous section.The first clearly dominates over the background reduction effect for the lead block sample (Fig. 5) with a measured total scattering component as high as 11% of the nominal beam intensity for the 0° orientation (to be compared to the 4% measured with the open beam only).The total scattering images (Fig. 5(b)) do not have sharp features and are logically well approximated by the coarse grid interpolations.The lower amount of sample scattered neutrons measured for the 90° orientation can be explained by a self-shielding effect: neutrons scattered in the forward direction with this orientation have a relatively high probability of being scattered a second time, and therefore there is a higher chance that neutrons exit the sample in a direction perpendicular to the beam axis and therefore outside of the field of view.This also explains the slight depression of the total background visible directly behind the sample for the 90° orientation.This depression is not perfectly fitted by the coarse grid interpolation, but the sample transmission (> 20%) is high enough for this to have a rather low impact on the corrected thickness profiles (Fig. 5(d)).The situation is different for the more challenging stainless-steel sample (Fig. 6), which has a very low transmission in the 90° orientation.Not only the depression of the total background behind the sample is stronger as for the lead sample, but the measured sample intensity is very close to the background and even a minimal interpolation error has a major impact.This is particularly true when the coarse grid is placed in an unfavorable position with no BB position overlapping with the sample.In this case, the interpolation results in an overestimation of the background and correspondingly to a significant bias in the thickness measurement.If the coarse grid is well placed (i.e.BB positions overlapping with the sample in the 90° orientation), the effect is much less pronounced as the positions with the highest interpolation errors are outside the sample.Finally, it is worth noting that the measured thicknesses significantly differ if the corrections described here (middle range PSF deconvolution and BB grid correction) are not applied.In particular, when the sample is in diagonal position, we do not obtain the expected triangular profile if we do not apply the black body correction.This first qualitative indication of the importance of applying the correction will be completed by a quantitative analysis of the tomography reconstruction.Our correction method is based on measuring the total scattering for a limited set of 25 projections evenly distributed over a full rotation.This represents only 4% of the images used for a full tomographic data set with 625 projections and therefore induces only an insignificant overhead in experimental time.However, for comparison, we also recorded the total scattering for the full 625 projections.The comparison shown in Fig. 7 indicates that there is little difference between the full measurement and the linear interpolation of the limited measurement.The most important divergence is seen at the 90° position, where the sharp minimum falls exactly between sampling positions.This induced an evaluation error of approximately 0.5% (see insert in Fig. 7(b)), which is expected to have a negligible impact in the case of the lead sample but might represent a significant bias for the stainless-steel sample having a very low transmission at the 90° orientation.

Tomographic reconstruction
The results of tomographic reconstruction are shown in Fig. 8 (lead alloy sample) and Fig. 9 (Stainless-steel sample).Here, only the base case for the black body correction (use of the middle range deconvolution, 25 projections for the black body measurements, coarse grid centered on the object) is shown.They are displayed as variations of the measured local attenuation coefficient relatively to the measured sample average.Thus, the results presented in these two figures do not show the average quantification bias, which will be discussed later.From the observation of the figures, it is obvious that the corrections proposed here leads to a significantly improved result better reflecting the actual homogeneous material.When the corrections are not applied, a higher attenuation coefficient is measured in the vicinity of the borders of the samples than in their center.This can logically be understood by the fact that the negative measurement bias introduced by an unaccounted background is stronger when large thicknesses are involved.Thus, the points in the center of the sample, for which the neutron transmission path through the object is -on average-the largest for the entire range of orientations, suffer from a larger negative bias than the points in the border.This explanation holds only for the distribution along the x and y axes.The z-axis distribution is rather to be attributed to the spatial distribution of the background components themselves.It is worth noting that this type of measurement artifacts (higher density measured at the edges than in the center of the samples) is similar to the artifact induced by beam hardening (reduction of the effective attenuation coefficient for large thickness samples due to the wavelength dependency of the cross section).However, we can easily exclude any significant impact of beam hardening for these materials which have a fairly low dependency of the cross section on the neutron wavelength, in particular for the thermal range used here.The qualitative improvements brought by the BB correction seem evident from these results, but have been further analyzed in order to quantify the impact of this correction.The average quantification bias observed in the reconstructed object in comparison to the real attenuation coefficient of the respective samples is shown in the top sections of Figs. 10 and 11.Additionally, the relative variations around this average value are shown in the bottom sections of these figures.The interval of relative variations was constructed as a 99% confidence interval around the average value (meaning the 0.5% smallest and largest values were discarded).To mitigate the impact of noise on this measurement, we did not use single voxel values, but the average value of 8x8x8 voxels cubes.The different correction strategies as outlined in Table 1 are applied, including the base case, which is suggested as standard correction.For both samples, when the black body correction is not applied, a significant quantification bias of up to 25% is observed.Applying the middle range deconvolution alone only brings marginal improvement.For the lead sample, the application of the black body correction reduces the measurement bias to approximately 5%, nearly regardless the correction parameters.  1.The baseline for the attenuation coefficient corresponds to the real attenuation coefficient of 0.335 cm −1 for this sample.The error bars and relative variations are based on a 99% confidence interval on the average value of all 8x8x8 voxel cubes (equivalent to volumes of 0.54 x 0.54 x 0.54 mm 3 ) contained in the sample.
On the contrary, the use of an adequate correction including a correct middle range deconvolution is of the utmost importance for the stainless-steel sample.The base case gives an excellent quantitative result, with an average bias of less than 3% and relative variations of less than 5%.Omitting the middle range deconvolution, however, has a significant negative impact on the results, leading to a significantly enhanced bias and relative variation of values.The application of the fine grid correction (grid in 16 different positions) and of the black body correction with 625 projections for the black body measurements both bring a slight reduction of the relative variations, though the much higher experimental effort is not justified in the light of the benefits for most applications.Finally, the use of a misplaced black body grid results in a massive overcorrection rendering the results even worse than the noncorrected ones -at least for the case where the deconvolution is omitted.This is in line with the results shown for single projections (Fig. 6), where the misplacement of the grid has a significant impact on the measurement regions of very low transmission.This is not an issue for plain samples roughly featuring a central symmetry respective to the rotation axis such as the blocks used here, because the lowest transmission region will obviously coincide with the rotation axis and the correct placement of the grid is evident.The situation will be very different for samples with more complex shapes/compositions, such as hollow or dissymmetrical, where the very low transmission regions can be decentered and changing as a function of the object orientation.For such samples, different measurement strategies can be considered.One of them would be moving the grid as a function of the orientation so as to provide the best possible coverage of the low transmission regions.Another approach to be considered is the use of iterative reconstruction algorithms, where, unlike in the filtered back projection used here, regions of too low transmission could be considered as unreliable and excluded from the calculus.It is worth noting that the method shown here is not only applicable to tomography, but can also be applied to single radiographies or to time series.However, in the latter case, the possibility of an additional "accelerated" experiment (such as the measurement with only 25 projections instead of 625 shown here for tomography) is excluded.There are a few options to deal with this issue.One of them, applicable for example to fuel cell research, is to measure the scattered background with the sample in a reference state and to assume that the sample changes over time do not induce a significant change to this background.For cases where the sample changes do affect the scattered background, a practical option is to perform the measurement with the BB grid always present, at the cost of a few sacrificial locations which are hidden behind the black bodies.1.The baseline for the attenuation coefficient corresponds to the real attenuation coefficient of 1.08 cm −1 for this sample.The error bars and relative variations are based on a 99% confidence interval on the average value of all 8x8x8 voxel cubes (equivalent to volumes of 0.54 x 0.54 x 0.54 mm 3 ) contained in the sample.

Conclusions
We demonstrated the use of a robust experimental method for the correction of quantification biases in neutron imaging introduced not only by scattered neutrons, but also by other background contributions appearing in scintillator-camera detectors.This method is based on an additional reference measurement using a grid of black bodies fully opaque to neutrons in order to measure the biases introduced by scattered neutrons and other effects such as light reflections.It implies a reasonable experimental effort with approximately 4% of additional exposure time and does not require any a priori information on the sample scattering behavior.Our results are in agreement with the previous findings reported by Ning et al. [27], in which however scattering correction is applied for X-ray CT cone beam reconstruction with flat panel detectors.Here, we also specifically discuss the implication of scintillatorcamera based detectors.Due to the similarities between the issues induced by neutron and Xray scattering in CT reconstruction, these results are expected to be relevant as well for X-ray imaging with scintillator-camera detectors, such as used for example for synchrotron imaging.The results on tomographies show a clear quantification improvement for an "easy" test sample made of lead and having a transmission of more than 20% in all orientations.The necessity of performing such corrections is even more evident on a "more challenging" stainless-steel sample having a transmission as low as 1% in given orientations.While, for the lead sample, the correction brought similar improvements regardless of the parameters applied, we observed that much more care has to be taken for low transmission samples such as the stainless-steel test object.In particular, misplacing the black body grid can, in this case, even result in a worse bias than when neglecting the correction.When correctly applied, the method was shown not only to reduce the average quantification bias by up to one order of magnitude, but also to significantly reduce the uncertainty around this average bias.

Fig. 1 .
Fig. 1.(a) Illustration of the measurement setup showing the position of the black body grid (BB grid).(b) Detail showing two examples of bias contributions, namely the scattering contribution from the sample and the mirror.

Fig. 2 .
Fig. 2. Simplified data flow diagram of the background correction process.The initial correction steps (filtering, deconvolution, camera offset and intensity correction) are not shown for the sake of clarity.

Fig. 3 .
Fig. 3. (a-b) images of the pencil beam measurements (two different sizes with the blue line showing the beam boundaries).(c) Horizontal intensity profile of the small pencil beam with Lorentzian and Gaussian fit.(d-f) Illustration of the background contributions occurring in the detector: Neutron scattering by the mirror (d), light reflection on the mirror (e) and neutron scattering by the scintillator substrate.

Fig. 4 .
Fig. 4. (a) Example measurement of the open beam with BB grid.(e) Resulting interpolated background image using the BB grid at 16 different positions.(f) Vertical intensity profile of the measured background.Fine grid corresponds to 16 different BB grid positions and coarse grid to one single position.

Fig. 5 .
Fig. 5. (a) Radiographs of the lead sample at different orientations.(b) Total scattering images interpolated from the 16 positions BB grid scanning ("fine grid" measurements).(c) Horizontal profiles of the total scattering at the middle of the image (only the lowest values of the sample transmission are shown).(d) Measured thickness after correcting for the scattering.Abbreviations: CC = coarse grid, centered on sample.CS = coarse grid, shifted.F = Fine grid.

Fig. 6 .
Fig. 6.(a) Radiographs of the stainless-steel sample at different orientations.(b) Total scattering interpolated from the 16 positions BB grid scanning ("fine grid" measurements).(c) Horizontal profiles of the total scattering at the middle of the image (only the lowest values of the sample transmission are shown).(d) Measured thickness after correcting for the scattering.Abbreviations: CC = coarse grid, centered on sample.CS = coarse grid, shifted.F = Fine grid.

Fig. 7 .
Fig. 7. Measured total scattering as a function of the sample orientation (superimposed measurements for the full 625 projections measurements and the limited 25 projections measurements).(a) Lead alloy sample (b) Stainless-steel sample (insert showing the detail of the interpolation error at 90°).

Fig. 8 .
Fig. 8. Result of tomographic reconstruction for the lead sample without (a) and with (b) the corrections (base case).The grayscale intensity shows the deviations from the average attenuation coefficient measured for the entire sample.(c) Horizontal profiles over the region marked with dashed lines.

Fig. 9 .
Fig. 9. Result of tomographic reconstruction for the stainless-steel sample without (a) and with (b) the corrections (base case).The grayscale intensity shows the deviations to the average attenuation coefficient measured for the entire sample.(c) Horizontal profiles over the region marked with dashed lines.

Fig. 10 .
Fig. 10.Average measured attenuation coefficient (top) and relative variations (bottom) of lead sample for all correction methods presented in Table1.The baseline for the attenuation coefficient corresponds to the real attenuation coefficient of 0.335 cm −1 for this sample.The error bars and relative variations are based on a 99% confidence interval on the average value of all 8x8x8 voxel cubes (equivalent to volumes of 0.54 x 0.54 x 0.54 mm 3 ) contained in the sample.

Fig. 11 .
Fig. 11.Average measured attenuation coefficient (top) and relative variations (bottom) of stainless-steel sample for all correction methods presented in Table1.The baseline for the attenuation coefficient corresponds to the real attenuation coefficient of 1.08 cm −1 for this sample.The error bars and relative variations are based on a 99% confidence interval on the average value of all 8x8x8 voxel cubes (equivalent to volumes of 0.54 x 0.54 x 0.54 mm 3 ) contained in the sample.