Detection of on-surface objects with an underground radiography detector system using cosmic-ray muons

Hirofumi Fujii, Kazuhiko Hara2,∗, Kohei Hayashi, Hidekazu Kakuno, Hideyo Kodama, Kanetada Nagamine, Kazuyuki Sato, Kotaro Sato, Shin-Hong Kim, Atsuto Suzuki, Kazuki Takahashi, and Fumihiko Takasaki 1High Energy Accelerator Research Organization (KEK), 1-1 Oho, Tsukuba, Ibaraki 305-0801, Japan 2University of Tsukuba, 1-1-1 Tennodai, Tsukuba, Ibaraki 305-8571, Japan 3Tokyo Metropolitan University, 1-1 Minami-Osawa, Hachioji, Tokyo 192-0397, Japan †Present Address: Iwate Prefectural University, 152-52 Suko, Takizawa, Iwate 020-0693, Japan ∗E-mail: hara@hep.px.tsukuba.ac.jp


Introduction
Muon radiography is a promising method for investigating the inner structure of large objects without the need for direct access by using cosmic muons [1][2][3][4][5][6].It is very effective for investigating the status of nuclear debris in the Fukushima Daiichi Reactors.After demonstrating the performance of a remotely located compact (1 × 1 m 2 cross section) muon radiography detector system by imaging the inner structure of a nuclear plant [7] and after modifying the detector systems to make it more robust against environmental radiation, we placed two similar detector systems at the Fukushima Daiichi Nuclear Plant.The Unit-1 Reactor measurement was consistent with a large portion of the fuel material not being in the original fuel-loading zone [8].
The muon radiography detector, which maps the distribution of cosmic muons behind the structure of interest, is limited in detection of the structure when it is at lower elevation angles.To further investigate the status of the fuel debris that has probably fallen to the floor, at the same level as a ground-level detector, a unique solution is to locate the muon radiography detector underground.
At KEK, we performed a series of measurements in which iron blocks on the ground surface were imaged while the detector was located on an underground floor.The iron blocks were typically 2 m PTEP 2017, 053C01 H. Fujii et al. × 2 m× 2 m (2 m cubic), as shown in Fig. 1.The purpose was to investigate whether on-surface objects can be imaged and whether the corresponding iron density can be derived by the same muon radiography detector system used previously [7].

Experimental setup
Figure 2 shows the configuration of the detector system.The detector [7] is composed of two pairs of XY planes 1.5 m apart.Each plane consists of 100 scintillator bars, each of which is 1 cm wide and 1 m long.A 1 mm diameter wavelength-shifting fiber embedded in each bar efficiently extracts the scintillating light to one end of the plane to which a photodetector multipixel photon counter (MPPC) [9] is coupled.A field-programmable gate array (FPGA) data acquisition (DAQ) system checks for coincidence hits in the XY planes within a set time window (32 ns), and the data with coincidence hits in both DAQ systems are recorded.The muon tracks are reconstructed from the two space points defined in the two XY detector units.
The detector was placed on a floor 8 m below the ground surface (see Fig. 1).The front XY unit was 40 cm higher than the rear unit so that the geometrical acceptance of the detector was maximized to the more distant (far) block.The two blocks were placed on the ground surface at z = 16 and 28 m from the detector.The 2 m cubic far iron block was placed outside the building.The closer (near) iron block was either 2 m cubic or 1 m× 2 m× 2 m.The muons that passed the near or far block traversed on average 8.5 or 21 m of soil, respectively, before reaching the detector.The configuration for the far block is similar to that at Fukushima Daiichi with respect to the soil thickness and the muon zenith angle if the detector were buried outside but close to the reactor building.Fuel material at the fuel-loading zone is typically 4 m cubic, so the sizes of the imaged iron blocks were selected accordingly.
Data acquisition began on August 28, 2014, and ended on December 15, 2014.During this period, the configuration of the near block was changed.

2D images of the on-surface objects
Two-dimensional (2D) lateral images obtained during two periods, 5 days with the near block only and 67 days with both blocks are shown in Fig. 3.The images are the lateral muon yield distributions normalized by the muon distribution obtained on the ground surface by taking an image of the sky to correct for the detector acceptance and azimuthal angle dependence of muon flux.The shadow areas indicate the existence of additional material along the muon path.The absorption due to the objects is evident in both images.There is little absorption in the upper part of the images because the soil is detected only below that height.Although the width of the near block in (b) is half that in (a), the near block is clearly seen in (b) and with a more uniformly distributed background.

Absorption in the soil: determination of the muon mass attenuation length
The first step in evaluating muon absorption by the blocks is to determine the absorption in the soil and concrete.We measured the absorption as a function of the sum of the muon path lengths in the soil, the concrete wall and the concrete floor (1.5 m thick) of the experimental hall.The building wall was not included in this calculation because its effect is important only in very small density length regions (density × thickness <500 kg m −2 ).In addition, the horizontal region of the iron blocks was excluded from the data analysis.
Figure 4 shows the yield N normalized by the yield measured for the sky N 0 , as a function of the linear density × thickness (density length), where the density of the soil and the concrete wall were assumed to be 1.8 and 2.3 g cm −3 , respectively. 1 The figure shows that the muon yield is very dependent on the density length when it is < 4000 kg m −2 , where the contribution of the scattering of low-momentum muons is pronounced.There is a more detailed discussion on this in Sect.3.4.In the density length range given above, the dependence is approximately exponential.Because the  near and far blocks are located near 13000 and 30000 kg m −2 , respectively, and the purpose of the measurement was to evaluate the muon mass attenuation around these points, the logarithm of the normalized yield is fitted to linear functions in the ranges of 4000-16000 and 10000-30000 kg m −2 respectively.The resulting mass attenuation coefficients are μ mass = (1.12 ± 0.05) × 10 −5 m 2 kg −1 for the near block region and (0.694±0.014)×10 −5 m 2 kg −1 for the far block region.The difference in the mass attenuation coefficients between the two blocks was due to the difference in the average muon momentum because the path length of the muons is dependent on the azimuthal angle of the cosmic muons and, thus, on muon momentum.

Absorption by the iron block
We studied the image of the iron block before determining the block density.We compared the horizontal event distributions in two 2 m height ranges projected at z = 16 m from the detector: one (A) about the height of the near block and the other (B) beneath the range of the near block Figure 5 shows the number of events in range A subtracted from that in range B measured over 14 days.The subtraction reference (B) is not the same as that in Fig. 3 because we avoid including the neighboring horizontal region in the evaluation of the horizontal image resolution.The absorption by the 2 m-wide iron block is clearly seen in Fig. 5 between x = 0 and 2 m.The sharpness of the two edges, defined as σ of the error function fitted to the edge distribution, was reconstructed differently, i.e., 10 ± 6 cm at x = 0 m and 16 ± 7 cm at x = 2 m.The values differed because the side wall of the block was projective to the detector at x = 0 m, and there were more tracks traversing a fraction of the non-projective block side wall around x = 2 m.In addition, the number of events at the intermediate region 0 < x < 2 m decreased with x reflecting the reduction of the detector acceptance.Figure 6 shows a similar horizontal image for the far iron block measured over 78 days.Although the statistics were limited by the lower elevation angle, the object is discernible in the plot.
After identifying the locations of the blocks, the absorption due to the iron was evaluated as a function of the path length through them.To minimize the effects of the difference in muon momentum due to the difference in zenith angle and of the detector acceptance, we assumed there was a 2 m cubic "air" block horizontally next to the iron block and at the same distance.For this purpose the block was not centered (see Fig. 3) to allow the "air" region to be set at a horizontally symmetric position.The contributions from the above-mentioned effects can then be assumed to be identical in these two regions.Figure 7 presents separate plots of the numbers of events through the assumed "air" region and through the iron block region as a function of the path length through the assumed air and iron block regions.The modulation in the event distributions is due to the finite width of the scintillator (10 mm) and to the limited number of scintillators (100) per plane: The combination of the finite-width scintillator bars in the two planes is favorable for events at digitized zenith angles, so specific measurable path lengths are not favored.The difference between the two distributions cancels out this digitization effect so only the absorption by the iron block is seen, as in  the inset in Fig. 7.Because the absorption is exponential as expected, the obtained linear attenuation length, μ block (m −1 ) is used with the mass attenuation coefficient μ mass (m 2 kg.−1 ) to calculate the density of the object: The density of the block obtained with this equation as a function of the DAQ time is presented in Fig. 8.The config number indicates the data for the different configurations of the near block: config-1 2 m cubic block offset to the left; config-2, 2 m cubic block offset to the right; and config-3, 1 × 2× 2 m 3 block offset to the right.The far block was in only config-4.The calculated density had an additional uncertainty of 10% because of the uncertainty of the soil density.The density of  the iron block was measured independently as 7.0 g cm −3 .Therefore the calculated density of the objects was consistent with the measured density.The uncertainty of the calculated density improved over time.

Understanding the mass attenuation distribution
The mass attenuation coefficient derived from Fig. 4 is a crucial result in this analysis obtained without relying on simulations.For calculation or simulation, the momentum distribution of cosmic muons is indispensable, but for muons at a large zenith angle, only the measurement [10] was performed at θ = 75°and momentum above 1 GeV/c.
The available muon flux data [10] as a function of momentum was convoluted with the muon range [11] followed by calculation of the attenuation in the soil.Figure 9 compares the measured attenuation with the calculated attenuation.Because the range of the zenith angle of the data (range >6000 kg m −2 ) in Fig. 9 (56°to 68°) is not covered in Ref. [10], the zenith-angle dependence mentioned in Ref. [10] was not incorporated in the present calculation.The agreement between the calculated and the measured attenuations is reasonable, which suggests that muon attenuation above 6000 kg m −2 can be explained consistently by the cosmic muon flux in Ref. [10].The data points below 6000 kg m −2 are due to muons with momentum <1 GeV/c.The change in slope in this region could be attributed to the enhanced muon flux relative to the simple extrapolated muon flux from Ref. [10].

Fig. 1 .
Fig. 1.Schematic of the experimental setup.An object (a 2 m cubic, or a 1 × 2 × 2 m 3 iron block) was located 16 m from the detector (near configuration) and another 2 m cubic iron block was located 28 m from the detector (far configuration).

Fig. 2 .
Fig. 2. Schematic of the detector system consisting of two XY plane units that are 1.5 m apart, each plane containing 100 scintillator bars.The basic detector unit, a scintillator bar embedded with a wavelength-shifting fiber is shown on the right.The scales are as indicated.

Fig. 3 .
Fig. 3. Images of muon yield distributions obtained (a) with a 2 m cubic near iron block only over 5 days, and (b) with a 1 × 2× 2 m 3 near block and a 2 m cubic far block over 67 days.The absorption due to the object is indicated by the shadow and the locations of the blocks are indicated by the lines.In (a), the location of "air" is indicated by the dash-lined box.This is subtracted from the absorption observed in the iron block region to obtain the absorption due to only the iron block.

Fig. 4 .
Fig. 4. Determination of the mass attenuation coefficient of the soil and concrete.The muon yield normalized by that for the sky is plotted as a function of the sum of the density length of the soil and concrete along the path of the muon.Data were collected over 15 days.The arrows indicate the locations of the blocks, and the data points were obtained for the region excluding the blocks.The dependence of the yield on the density length around the near and far iron blocks was approximately linear as shown by the solid and dashed lines.

Fig. 5 .
Fig. 5. Number of events observed through the 2 m height range below the near object after subtracting the number of events observed through the 2 m height range containing the near object as a function of the horizontal position projected at z = 16 m.The data were collected over 14 days.The 2 m cubic iron block spanned from 0 to 2 m.

Fig. 6 .
Fig. 6.Number of events observed through the 2 m height range below the far object after subtracting the number of events observed through the 2 m height range containing the far object as a function of the horizontal position projected at z = 28 m.The data were collected over 78 days.The 2 m cubic iron block spanned from 0 to 2 m.

Fig. 7 .
Fig.7.Normalized muon yield as a function of path length through the iron block (blue circles) and through the assumed air block (red points).Inset: the difference between the two yields is fitted to a linear function to extract the attenuation length in the iron block.

Fig. 8 .
Fig. 8. Calculated density of the iron blocks as a function of the data acquisition time.Config-4 included detection of the far block, whereas the other three configurations included only the near block object.The error bars are statistical only and the density has a systematic uncertainty of 10%.