Investigation of the limits of high-definition muography for observation of Mt Sakurajima.

A multi-wire proportional chamber-based muo- graphy observatory is under development for the monitoring of the internal structure of Mt Sakurajima in Kyushu, Japan. We investigated the limits of large-scale and high-definition muography. We adjusted the parameters of a modified Gaisser model and found that the spectral index of γ =  - 2.64 and normalization factor of C = 0.66 reproduce more accurately the measured fluxes than the original parameters at large thickness. A thickness and zenith angle-dependent correction is suggested to the measured muon flux due to the energy cut which is introduced to suppress the background particles. The multiple scattering of muons was simulated across the standard rock and sea-level atmosphere up to the distance of 5 km. We found that multiple scattering decreases from 10 mrad to 4 mrad across the rock due to the decrease in the steepness of muon spectra. The multiple scattering falls down to about 2 mrad after the object in the atmosphere due to the increase in observed arrival zenith angles. The 2 m2 sized multi-wire proportional chamber-based Muographic Observation System (MMOS) was operating between February and June 2018. Three tracking systems operated reliably with tracking efficiencies of above 95%. The muon flux has been measured correctly down to 10-3 m-2 sr-1 s-1 The average density map of Mt Sakurajima has been measured with angular resolution of 12 mrad × 12 mrad (spatial resolution of 34 m × 34 m from the distance of 2.8 km). The average density values were found between 1.4 and 2 g cm-3, except at the crater regions where lower densities were observed.This article is part of the Theo Murphy meeting issue 'Cosmic-ray muography'.

Introduction

Muography is a novel technique for the visualization of the internal structure of large objects [1]. The imaging process of muography is similar to radiography: tracking detectors are installed around an object and the absorption rate of penetrated cosmic-ray muons is measured. If the object is modelled, the average density length can be deduced by the comparison between measured and expected fluxes. Measurements from different directions allow tomographic imaging [2]. In addition, the newly developed airborne muography is not restricted by local topology and the imaging resolution of computational axial tomography can be achieved [3]. The applicability of cosmic-ray muography has been demonstrated to explore the internal structure of natural formations [2-11] and for inspection of human-made constructions [12-14].

The Hungary–Japan Joint Muography Observatory is focusing on the investigation and observation of the internal structure of Mt Sakurajima, which is located in Kyushu, Japan (see figure 1). Mt Sakurajima is one of the most active volcanoes with hundreds of explosive eruptions every year. The time-sequential imaging of density variations below the erupting craters could help to understand volcanic activity. Furthermore, the construction of a database with the variation of muon flux across the crater regions would contribute to a complex eruption prediction system which combines different geophysical techniques, such as gas emission [16], seismograph [17], video [18], GPS [19], satellite-based SAR [20] or LIDAR [21].

Figure 1.

The Sakurajima muography campaign. (a) Photograph of Mt Sakurajima with the Craters A, B and Showa [15]. The location and orientation of MMOS are shown by a star and an arrow, respectively. (b) Photograph of the MMOS with three tracking systems in February 2018. (Online version in colour.)

Mt Sakurajima is a difficult target for muography because of the kilometre-thick crater regions, which result in a muon flux below 10−2 m−2 sr−1 s−1. The craters A and B (Minami-dake) have been imaged with the angular resolution of 3 mrad with the first module of the MWPC-based Muographic Observations System (MMOS) [11,22] from the distance of 2.8 km in the southwest direction from the Showa crater in 2017. However, the large-scale extension of MMOS is necessary to provide time-sequential muographic images of Mt Sakurajima. Besides the limited flux of muons, the accuracy of large-scale muography is influenced by the following factors. This paper aims to improve the accuracy of muon flux calculation and measurement for large-scale muography, quantify the effect of multiple scattering of muons on the imaging resolution and give a report on the present status of the Hungary–Japan Joint Muography Observatory at the Sakurajima volcano. The paper is organized as follows. §2 focuses on the energy spectra of muons. §3 investigates the multiple scattering of muons. §4 focuses on the present status of the Sakurajima muography campaign. The results are summarized in §5.

Physical effects: The creation of muons is influenced by the variation of pressure and temperature in the atmosphere [23,24]. The flux of penetrated muons after the investigated object is comparable with the flux of background particles [25], which are also detected from the direction of the object and result in overmeasured flux as well as underestimated density. The background consists of low-energy muons [26,27], soft [28] and hadronic components of cosmic-rays. Besides the backround particles, the multiple scattering of muons across the investigated object causes a blurring effect on the muographic image.

Accuracy of muon flux modelling: The expected flux is determined by the integration of altitude corrected muon spectra from minimal energy which is necessary to penetrate the object or by Monte Carlo simulation of air showers across the atmosphere [23]. For the former approach, different muon spectra parametrizations [29-33] have been developed based on measurements performed in different energy ranges, see e.g. [34,35]. These discrepancies result in the differences of 10–30% between the different spectra models [23] and between the expected and measured fluxes [31,36,37]. The refining of muon spectra parametrization using world data could improve the accuracy of the calculated flux.

The detector effect: This depends on the applied technology (e.g. scintillator, nuclear emulsion, gaseous detector) as well as data acquisition (DAQ) and readout systems. Shielding plates are applied between the detector layers to deflect and absorb the background particles. The application of shielding is practically a momentum (energy) cut. The determination of an appropriate energy cut and the optimization of detector arrangement is suggested for each experiment.

Investigation of muon spectra for large-scale imaging

Adjustment of muon spectra parametrization

For large-scale muography, the high-energy region of the muon spectra has to be modelled accurately. These are described with the so-called modified Gaisser parametrization [29,30,33]: where dN/dE dΩ is the differential flux of muons, C = 1 is an adjustable constant factor, γ = − 2.7 is the spectral index inherited from the spectra of primary cosmic-rays. The 115 GeV and 850 GeV are critical energies for pions and Kaons, above these energies, the pions and Kaons interact with the atmosphere before they decay. The correction takes into account the curvature of the atmosphere and extends the zenith angle range of the parametrization [30,38] with the constant parameters of p1 = 0.102573, p2 = 0.068287, p3 = 0.958633, p4 = 0.0407253 and p5 = 0.817285. In addition, low-energy and high-energy muon spectra parametrizations have been introduced in [31,32], respectively.

Various groups calculated and simulated the near-vertical flux of muons based on parametrized spectra for underground experiments and reported 10–30% differences between their results and experimental data [31,36,37]. The original Gaisser parametrization has been refined and C = 0.805, γ =  − 2.643 parameters have been suggested to reproduce accurately the so-called Crouch curve [36]. The comparison of modelled fluxes to near-horizontal experimental data is suggested to improve the accuracy of large-scale muography.

For this study, the muon spectra have been parametrized according to equation (2.1) and integrated from minimum energies based on CSDA range [39] for different thickness of standard rock at different zenith angles. The integrated flux values have been compared to fitted experimental data measured deep underground [40-42]. Figure 2 summarizes the experimental data (black empty dots) and the calculated flux values with the original (blue dashed line) and best fitting (red solid line) parametrizations. The parametrizations presented in [31,32] are also plotted with dashed green line and solid cyan line, respectively. We found that the modified Gaisser parametrization with γ = − 2.64 and C = 0.66 parameters produces a minimal difference between the expected fluxes and experimental data. The average values of the relative flux differences, (Fmeasured − Fexpected)/Fmeasured, were found to be 10.4% and 1.6% for the original (blue dashed line) and for the best fitting (red solid line) parametrization, respectively.

Figure 2.

(a–c) The muon flux is plotted as the function of the zenith angle at different depths for experimental data (black empty dots), for the original (blue dashed line) as well as for optimal (red solid line) parametrization, respectively. Further models are drawn with green and cyan lines. The coloured lines on the calculated points are drawn to guide the eyes.

An energy dependent correction to the flux penetrated muons

The suppression of low-energy background particles is required for large-scale muography. The energy spectra of background particles were found to be different at different measurement sites, e.g. they were below 1 GeV at Showa–Shinzan lava dome [25] and below 5 GeV at La Soufriére volcano [27]. These results suggest target specific investigation of background to determine the required energy cut and optimal detector arrangement. This section aims to investigate the reduction of the flux of penetrated muons due to the introduced energy cut independently from detector arrangement.

A simulation code based on GEANT4.10.01 [43] has been developed for this study. For each simulation, the object was implemented with different thicknesses between 250 and 2000 m standard rock equivalent (m.s.r.e.). The muons were injected near-horizontally (θ0≥70°) across the object and their energy spectra were parametrized according to Reyna and modified Gaisser models. A total of 40 000 muons were generated in each run. The minimum vertex energies of muons were set well below the minimum energies [39], which are required to penetrate the object. The application of lower energy minima was motivated by the energy loss of high-energy (>500 GeV) muons, which is dominated by stochastic processes (pair production, Bremsstrahlung, nuclear interaction) and some of the lower energy muons can also penetrate the object. The standard physics list with all of the relevant electromagnetic processes and muon nuclear processes were applied [44].

Figure 3 shows the ratio of the number of muons above the applied energy cut to the total number of penetrated muons, N( > Ecut)/N, as the function of energy cut for various thicknesses and zenith angles. The application of an energy cut of a few GeV decreases the number of penetrated muons by only a few per cent, which is consistent with the results of [25,27]. Therefore, large-scale muography can be performed with a muon flux correction of a few per cent. The N( > Ecut)/N depends on thickness because of the decrease in the steepness of muon spectra at larger energies (relatively more energetic muons penetrate across thicker objects). A slight zenith angle dependence is also observed.

Figure 3.

The ratio of the number of detected muons to the number of penetrated muons as the function of energy cut.

Investigation of muon scattering for high-definition muography

The spatial resolution of a muographic experiment, Δximaging, is given by the angular resolution of MOS, ΔαMOS and the object-detector distance, Lair: Δximaging ≈ ΔαMOS × Lair. Therefore, the desired imaging resolution and the accessibility of the object (e.g. it is prohibited to approach Mt Sakurajima closer than 2.5 km) determine the expected angular resolution of MOS. Newly developed, versatile and robust detector technologies allow imaging resolution of below 10 mrad [11,14]. This section aims to investigate the multiple scattering of muons and quantify the upper limit of the angular resolution of muography.

The muons are deflected across the transversed object due to a series of Coulomb scatterings on the nuclei of the medium. The distribution of scattering angles after a small part of medium (0.001 < L/X0 < 100) follows a zero-centred Gaussian distribution [45] with a long tail and its root mean squared (RMS) value, σMS, is approximated as follows: where E is the energy of muons and L/X0 is their path-length across the medium in radiation length units. The total scatterings of muons across the medium are given by the sum of small scatters along their path. The muons are losing continuously their energy across the medium, thus their scattering angles are increasing continuously across the medium.

The multiple scattering cause the observed directions of penetrated muons (red arrow in figure 4) to differ from the effective paths (blue arrow in figure 4) along which they carry information about the amount of transversed material. The multiple scattering of muons in the air, σair and across the tracking system, σMOS, have to be taken into account for long-range muography. The angular resolution of muography is suggested to be chosen above the total multiple scattering to avoid a smearing effect on the images: We note that the background particles also smear the muographic images. In the above criterion, we assume that the background particles are suppressed in the tracking system.

Figure 4.

The schematic drawing of land-based muography. The multiple scattering of a muon across the object causes that its observed trajectory (red arrow) differs from the effective one (blue arrow) which carries information about the amount of material along the path of muon across the object. (Online version in colour.)

The sea-level atmosphere and a tracking system were also implemented in the above presented simulations (see in §2b) to study the multiple scattering of muons. The muons were tracked from the vertex to the end of the tracking system and their position vector, momentum direction vectors, energies were recorded in every 100 m across large-scale (≥500 m.s.r.e.) objects and in every 10 m across medium-scale ( ≤ 250 m.s.r.e.) objects, respectively. The geometrical arrangement of the simulations is presented in figure 4. An object made up from standard rock (i.e. mixture of MgCO3 and CaCO3 with density of 2.65 g cm−3) and a 2 m length MOS with eight MWPCs and five 2 cm thick lead plates were installed. The atmosphere was composed of N2 (78.1%), O2 (21%) and Ar (0.9%) and the pressure and temperature were set according to sea-level data (https://ccmc.gsfc.nasa.gov/modelweb/atmos/us_standard.html (last view 27/07/2018).). The average terrestrial magnetic field of B = 45 μT was also applied.

A track-by-track analysis has been performed on the simulated data to quantify the multiple scattering of muons. The analysis was initiated with the calculation of the projected scattering angles of muons. These were defined by the differences between the expected and observed arrival zenith angles. The projected scattering angles before the MOS (Δθobject+air), across the atmosphere (Δθair) and across the MOS (ΔθMOS) are described by the following equations: where ph and pv denote the horizontal and vertical components of the momentum direction vector; xh and xv denote the horizontal and vertical components of position coordinates (see in figure 4).

Finally, the RMS values were calculated for the central 98% of the scattering angle distributions to determine the multiple scattering values for different locations before the tracking system. Figures 5 and 6 show the multiple scattering of muons after different objects up to the distance of 5 km. The results are summarized as follows. The systematics of the simulations have also been investigated. The systematic errors caused by the Earth's magnetic field and physical properties of the atmosphere were found to be negligible. The chemical composition (Z, A) of rock was found to be a significant source of systematics errors. Figure 7 shows the multiple scattering of 75° muons after standard rock and SiO2 with the same density length of 1000 m.s.r.e. and 3000 m.s.r.e., respectively. The differences between the different compositions were found to be 5–10%. The observed differences are caused by the pair production and Bremsstrahlung energy loss processes which are scaled with Z2/A.

Figure 5.

The multiple scattering of muons across medium-scale ( ≤ 250 m.s.r.e.) objects. The ratio of horizontal component to vertical component of momentum direction vector (a–c), the multiple scattering before the MOS (d–f) and across the air (g–i) are plotted as the function of object-MOS distance.

Figure 6.

The multiple scattering of muons across large-scale objects. The ratio of horizontal component to vertical component of momentum direction vector (pv/ph), the multiple scattering before the MOS (σobject+air) and across the air (σair) are plotted as the function of object-MOS distance for 500 m.s.r.e. (a), 1000 m.s.r.e. (b), 1500 m.s.r.e. (c), 2000 m.s.r.e. (d), 2500 m.s.r.e. (e) and 3000 m.s.r.e. (f), respectively.

Figure 7.

Comparison of the multiple scattering of 75° muons after standard rock and SiO2 with the same density length. (Online version in colour.)

Both expected scattering angle of muons (pv/ph) and multiple scattering (σobject+air) decreases with the rock thickness from 15 mrad to 5 mrad and from 10 mrad to 4 mrad, respectively. The observed decrease in multiple scattering is caused by the decrease in the steepness of muon spectra. We note that tracking systems are placed in limited space in case of underground muography, thus multiple scattering can be approximated with σobject+air(Lair = 0). Based on this observation, the angular resolution of 10 mrad is suggested for portable devices applied for underground muography.

The multiple scattering (σobject+air) falls after an object in the atmosphere due to the increase of the observed arrival zenith angle (Δxv/Δxh). Consequently, it is optimal to install the MOS at a larger distance from the object to decrease the effect of multiple scattering. In the case of medium-sized objects, there is an optimal distance between 0.1 and 1 km where the effect of multiple scattering is minimal, and beyond this optimal distance the contribution of air (σair) dominates.

As it is expected, σair is increasing with the distance from the object. This effect and its contribution to the total multiple scattering (σobject+air) decreases with the thickness of the object due to the decrease in the steepness of muon spectra.

The difference between the same type of quantities under different zenith angles is less than 10%, thus multiple scattering across large-scale objects has a slight dependence on the zenith angle of muons.

The multiple scattering of detected muons across the MOS (σMOS) is found to be ∼3.5 mrad, which is consistent with our previous work [11]. The σMOS depends on the positional resolution of tracking layers, the length of the tracking system and the amount of shielding plates, and thus it should be determined for each experiment.

The MWPC-based Muography Observatory at Mt Sakurajima

An MWPC-based Muographic Observation System has been under development at the Sakurajima volcano since January 2017. The first high-definition images of Minami-dake are presented in [11]. The MMOS has been extended with two tracking systems to the sensitive surface of 2 m2 in February 2018. Each MMOS has at least five MWPCs and five 2 cm thick lead plates. There are 20 MWPCs in the three tracking systems. The length of each tracking system is 2 m. The tracking systems are housed in steel containers and isolated from the environment with Plexiglass walls. A non-flammable and non-toxic Ar-CO2 gas mixture (Ar 80: CO2 20) is supplied parallel to the tracking systems with the flow of 2 l h−1. The electrical supply is provided by +110 V AC which is converted to +12 V DC. Each tracking system is protected by a +12 V (50 Ah) lead gel battery against power cuts. The total power consumption of MMOS is found below 30 watts (below 10 Watts per MMOS).

The tracking systems are operated independently from each other by microcomputer-based DAQ and control systems. The data acquisition process and the system are briefly presented in [11,22]. The MWPCs are supplied with a high-voltage of +1750 V to achieve a reasonable gas gain of above 103. The trigger is provided by the triple coincidence of MWPCs. The front-end cards amplify the analogue signals and transfer the bits to the DAQ with shift registers. The data are compressed and transfered via a virtual private network to a remote server in every hour. Track reconstruction and track-by-track analysis are performed based on the analysis presented in [11]. The track rates, muon fluxes, are monitored on the remote server. Figure 8 shows the variations of tracking efficiencies (εtracking) and environmental parameters during the period from 2 February to 1 June 2018. The two missing periods correspond to two technical stops between 19 and 24 February and between 24 and 26 of April. The detectors operated reliably with ε > 95% during the measurement period of 4 months, except MMOS01 which suffered from electrical noise due to the daily variation of humidity.

Figure 8.

The variation of tracking efficiencies (εtracking) of the three tracking systems (a), the temperature (b) and the humidity (c) during the data taking period between 2 February and 1 June. (Online version in colour.)

All tracking systems have to be oriented to the same direction within the imaging resolution before merging of their fluxes. An off-line correction method was developed for the orientation of tracking systems. First of all, the fluxes were calculated for each MMOS with the angular bin size of 6 mrad × 6 mrad. Thereafter, each flux was compared to the flux of the reference tracking system (MMOS03) and the flux residuals were calculated in different detector orientations. The residuals were found minimal when the tracking systems were oriented to the direction of MMOS03. An azimuthal rotation of 12 mrad and 30 mrad was applied to MMOS01 and MMOS02, respectively.

The expected fluxes have been calculated with different densities using the modified Gaisser model (see in equation (2.1)) and the thickness map of Mt Sakurajima. Figure 9a shows the elevation angle dependence of the muon flux at the azimuth angle of 0 mrad with the angular bin size of 12 mrad. The flux of muons has been measured correctly down to 10−3 m−2 sr−1 s−1, where the contribution of background particles became dominant. The density map of Mt Sakurajima was also measured with the angular resolution of 12 mrad × 12 mrad (see in figure 9b). This angular resolution corresponds to the spatial resolution of 34 m × 34 m from the distance of 2.8 km. The average density values of Mt Sakurajima were observed to be between 1.4 g cm−3 and 2 g cm−3. The lower density regions (blue patches) show the structure of Minami-dake in the top region of Mt Sakurajima. The high-density region observed in Crater A is under investigation. The high-density regions in the downward-going ridge of the volcano are caused by the overestimation of muon flux by the high-energy spectra model. The relative errors of density values were found to be 10–20% between the dashed black lines. The Showa crater erupted only three times during the presented measurement period, thus these data can be used as a reference data set for time-sequential measurements of muon fluxes and average densities in the crater region.

Figure 9.

(a) The elevation angle dependence of muon flux at the horizontal angle of 0 mrad. (b) The density map of Mt Sakurajima with the angular bins of 12 mrad × 12 mrad (34 m × 34 m from the distance of 2.8 km).

Summary

An MWPC-based muography observatory is under development for the monitoring of the internal structure of Mt Sakurajima. This paper aimed to investigate the limits of large-scale and high-definition muography and present the status of the Sakurajima measurement campaign.

The energy spectra of muons have been investigated. We found that the modified Gaisser model reproduces more accurately the earlier experimental data at large rock thicknesses with γ = − 2.64 and C = 0.66 parameters than the original parametrization. The penetration of muons was simulated across standard rock to investigate the reduction of muon flux due to the energy cut applied to suppress the background particles. We found that the number of penetrated muons decreases a few per cent with the application of energy cut of few GeV. A thickness and angular-dependent correction of muon flux is recommended to experiments which apply energy cuts above few GeV. The multiple scattering of muons has been investigated to determine the upper limit of angular resolution of muography. We found that multiple scattering decreases from 10 mrad to 4 mrad with rock thickness due the decrease in the steepness of energy spectra and it falls after an object in the the atmosphere to about 2 mrad due to the increase of the observed arrival zenith angles.

The 2 m2 sized MMOS operated reliably with tracking efficiencies of above 95% between February and June 2018. The flux of penetrated muons was measured correctly down to 10−3 m−2 sr−1 s−1. The craters of Mt Sakurajima were visualized with spatial resolution of 34 m × 34 m from a distance of 2.8 km. The average density values of Mt Sakurajima were found between 1.4 g cm−3 and 2 g cm−3. The Showa crater erupted only three times during the measurement period, thus these data can be used as a reference dataset for time-sequential imaging of flux and density variations during active periods. The improvement and industrialization of MWPCs are ongoing to construct a large MMOS and allow time-sequential muography of Mt Sakurajima.