Total mass estimate of the January 23, 2018, phreatic eruption of Kusatsu-Shirane Volcano, central Japan

On January 23, 2018, a small phreatic eruption (VEI = 1) occurred at the Motoshirane Pyroclastic Cone Group in the southern part of Kusatsu-Shirane Volcano in central Japan. The eruption ejected ash, lapillus, and volcanic blocks from three newly opened craters: the main crater (MC), west crater (WC), and south crater (SC). Volcanic blocks were deposited up to 0.5 km from each crater. In contrast, the ash released during this eruption fell up to 25 km ENE of the volcano. The total mass of the fall deposit generated by the eruption was estimated using two methods, yielding total masses of 3.4 × 104 t (segment integration method) and 2.4 × 104 t (Weibull fitting method). The calculations indicate that approximately 70% of the fall deposit was located within 0.5 km of the craters, which was mainly attributed to the low height of the eruption plume.


Introduction
At 10:02 on January 23, 2018 (all times in this paper are in JST), a phreatic eruption was initiated and lasted for a few minutes in the northern area of the Motoshirane Pyroclastic Cone Group (MPCG), one of two summit pyroclastic cone groups of Kusatsu-Shirane Volcano in central Japan . The eruption occurred at three newly opened craters, two of which were located only a few hundred meters away from a popular ski resort ( Fig. 1a, b). During the eruption, fall deposits were distributed to the NE and were observed up to 25 km ENE of the source craters (The Joint Research Team for ash fall from Kusatsu-Shirane 2018 eruption 2018a; b). Volcanic blocks from the eruption were ejected up to 0.5 km from the source craters ) and, unfortunately, hit several skiers and ski gondolas. Consequently, one person was killed, and 11 others were injured (Fig. 1a).
Immediately after the 2018 eruption, a research team consisting of researchers from several universities and public institutes surveyed the fall deposits. Preliminary results of this survey and total discharged mass calculations were reported at the 140th Meeting of the Japanese Coordinating Committee for the Prediction of Volcanic Eruptions (3 × 10 4 -5 × 10 4 t; The Joint Research Team for ash fall from Kusatsu-Shirane 2018 eruption 2018a). However, the mass of the fall deposits around the craters could not be measured by the survey team, as entry into the area near the source craters was restricted immediately after the eruption for safety reasons. As a result, the calculation of the discharged mass around the source craters was based on thicknesses estimated from aerial photographs taken before and after the 2018 eruption.
In this paper, we present the entire distribution of fall deposits and total discharged mass estimates for the 2018 eruption, including results from field surveys conducted near the source craters in April, May, and September Open Access *Correspondence: kame@mfri.pref.yamanashi.jp 2 Mount Fuji Research Institute Yamanashi Prefectural Government, 5597-1 Kenmarubi, Kamiyoshida, Fujiyoshida, Yamanashi 403-0005, Japan Full list of author information is available at the end of the article of 2018. Several small craters of the same size as those formed in the 2018 eruption are located in the MPCG, suggesting that eruptions of the same magnitude and type as the 2018 eruption have occurred in the past and are likely to occur in the future. Therefore, the distribution of the fall deposits and the total discharged mass obtained in this study are essential for understanding eruption magnitudes and evaluating hazardous areas of past and future phreatic eruptions in the MPCG.
For discussion purposes, we define the area, where volcanic blocks/lapillus were deposited with ash as the proximal area (< 0.5 km from the source craters, which is roughly equivalent to the area enclosed by the 30,000 g/ m 2 isomass contour described later). The area beyond the 0.5 km limit is referred to as the distal area.

Outline of the geology of the Motoshirane Pyroclastic Cone Group and the 2018 eruption
Three pyroclastic cones (two of which are complexes of multiple pyroclastic cones) are located in the summit area of Kusatsu-Shirane Volcano (Fig. 1a), namely, the MPCG, the Ainomine Pyroclastic Cone (APC; it is undated but is likely Holocene in age based on its well-preserved original morphology), and the Shirane Pyroclastic Cone Group (SPCG; it has a crater lake called "Yugama"), of which the SPCG has been the site of numerous historical phreatic eruptions (Japan Meteorological Agency 2013). In this section, we briefly describe the stratigraphic relationship of the MPCG. A detailed description of the volcanic geology of the MPCG can be found in Ishizaki et al. (2020).
The MPCG, which was the site of the 2018 eruption, consists of five overlapping pyroclastic cones, including Motoshirane-nishi, Older and Younger Motoshirane, Kagamiike (~ 4800 cal BP), and Kagamiike-kita (which is slightly younger than 1500 cal BP), from south to north. Each cone has a large summit crater and several smaller craters, the latter being distributed within the large crater and on the sides of the cone. Ishizaki et al. (2020) named the large craters in each cone: Motoshirane-nishi Crater, Older Motoshirane Crater, Younger Motoshirane Crater, Kagamiike Crater, and Kagamiike-kita Crater (Fig. 1b).
The 2018 eruption began at 10:02 on January 23. Volcanic tremors were recorded for 8 min beginning at 09:59 . According to the record of the change in the tilt , the area around the MPCG uplifted from 10:00 to 10:02 and subsided from 10:02 to 10:10, suggesting that the subsidence was due to the release of eruption products. Skiers near the gondola station saw a black vertically elongated jet that spouted from the western foot of the Kagamiike-kita Pyroclastic Cone (Sankei News 2018). Tourists on a gondola also witnessed volcanic blocks falling and hitting their gondola's roof and windows, followed by ashfall (Kyodo News 2018). The maximum plume height was ~ 5500 m based on radar observations (Sato 2021). A ground-hugging flow displaying hybrid pyroclastic surge and avalanche characteristics also occurred on the northern foot of the Kagamiike-kita Pyroclastic Cone during the early stages of the eruption .
The 2018 eruption formed three craters ( Fig. 1a, b). The largest crater (hereafter referred to as the main crater, or MC), with dimensions of 140 m on the long axis and 20 m on the short axis (Chiba et al. 2018), formed in the northern area of the pre-existing Kagamiike-kita Crater. The other two smaller craters, the west crater (WC) and the south crater (SC), formed at the western foot of the Kagamiike-kita Pyroclastic Cone and on the eastern side of the pre-existing Kagamiike Crater, respectively.
The 2018 eruption products comprised a dark gray clayey ash layer deposited as ash aggregates in the distal area (Fig. 1c) and weakly stratified, matrix-supported mixtures of dark gray clayey ash, lapillus, and volcanic blocks in the proximal area. At location c (4.8 km NE of the MC), the thickness of the fall deposit was < 1 cm, and ash aggregates had diameters reaching 0.5 cm (Fig. 1c). At location d (0.2 km NE of the MC), the thickness of the fall deposit was 10 cm and contained lapillus with diameters reaching 4 cm ( Fig. 1d). At the rim of the MC (Fig. 1e), the thickness of the fall deposit was ~ 180 cm and contained angular to subangular volcanic blocks with diameters reaching 30 cm. Volcanic blocks that were not buried in the ash layer were also deposited, and their distribution limit is shown by the blue dashed line in Fig. 1a. The 2018 eruption products contained abundant clay minerals derived from a high-temperature acidic alteration zone beneath the volcano, and no juvenile materials were observed (Yaguchi et al. 2019).

Methods
The 2018 eruption occurred during the winter. Consequently, the fall deposits were deposited on top of the snow and were quickly covered and packed in by the falling snow (Fig. 1c, d). Due to legal restrictions, access within 2 km of the craters was restricted for three months following the eruption, and the distal fall deposits were primarily sampled on January 23-25 and 29, February 28 and March 8, 2018. Most of these deposits, although well preserved in the snow, were thin and mixed with the snow, making it difficult to measure their thicknesses with reliable accuracy. Therefore, in this study, samples of the fall deposits were collected with reference to the method described by Scott and McGimsey (1994), including the under-and overlying snow layers.
Depending on the amount of the fall deposit at a given site, samples were collected from surface areas varying from 0.04 m 2 to 0.25 m 2 (to make them easier to carry, a small area was sampled if the amount of the samples was large, while a large area was sampled if the amount was small). In the laboratory, the mixed snow and ash samples were thawed in a beaker at room temperature (~ 20 °C), and the water was either filtered (if the amount of material was small) or drained from the beaker (if the amount of material was large). The residues left on the filter paper or in the beaker were then dried in an oven (approximately 40-60 °C) and weighed. The weight of the residue was divided by the sample area to calculate the mass per unit area (MPUA; g/m 2 ). Using the same process, samples were also collected from the roofs of the gondolas, which were stopped due to a power outage immediately after the start of the eruption and later stored at the gondola station. These fall deposits were collected from areas covering 0.01 m 2 due to the limited area available for sampling. The location of each gondola when the power went out was determined from the operation schedule provided by the Kusatsu Tourism Corporation (green circles in Fig. 2b). Field surveys within 2 km of the source craters were conducted on April 21, May 5 and 10-11, and September 11, 2018, once the safety restrictions were lifted. Proximal deposits around the MC were too thick (reaching ~ 2 m) to calculate the MPUA using the Scott and McGimsey (1994) method used for the distal deposits. Therefore, we converted the thickness data to MPUA using a density of 700 kg/m 3 for the fall deposits. This density value was based on average density measurements of fall deposits at two sites located 0.2 km NE and W of the MC (609 kg/m 3 and 778 kg/m 3 , respectively). The thicknesses of the deposits around the MC were measured only at the crater rim by our survey. Therefore, we used the thickness data estimated by Chiba et al. (2018) as the thicknesses of the eruption products around the crater. These data were obtained by subtracting the elevation of the surface before the eruption from that of the surface after the eruption. Isomass contours of 5, 20, 50, 100, 250, 500, 1000, 5000, 10,000, 30,000, 45,000, 60,000, 700,000, and 1,260,000 g/ m 2 were drawn (Fig. 2) based on the MPUA values obtained at each sampling point, and the areas enclosed within each contour (including the contour lines of the SC) were calculated using QGIS software. We calculated the total discharged mass using three methods. The first method was the segment integration method (Takarada et al. 2001). Because this method uses multiple massarea segments, more available isomass data produce more accurate mass values. For this reason, the segment integration method has been used for recent small-scale eruptions in Japan, such as the 2008 Shinmoedake eruption (Geshi et al. 2010) and the 2014 Ontake eruption (Takarada et al. 2016). The isomass values were plotted against the isomass areas (Fig. 3a), and the contours were subdivided into eight segment lines. The equations of each regression line were used to integrate the total mass. The second method was the Weibull fitting method (Bonadonna and Costa 2012), which can reproduce the gradual thinning of fall products even if less data are obtained, and it does not depend on choosing arbitrary segments or arbitrary extremes for integration. Here, the MPUA (g/m 2 ) and the square root of the area (km 1/2 ) of each mass contour are plotted (Fig. 3b). In this method, the total mass was obtained by integrating the Weibull fit of these point data. The third method was the empirical formula V = 12.2 TA (Hayakawa 1985), where V is the bulk volume, T is the thickness of an isopach, and A is the area enclosed by the isopach. We calculated the total discharged mass using 14 isomass contour lines and a density of 700 kg/m 3 to convert the MPUA to thickness data.

Estimate of total discharged mass
The isomass contour map produced in this study (Fig. 2a) indicates that the dispersal axis of the fall deposits was oriented NE in the proximal area and E in the distal area. The areas enclosed by each contour are shown in Table 1. The MPUA values obtained along the dispersal axis decreased exponentially with increasing distance from the MC: 60,897 g/m 2 at 0.2 km, 11,221 g/m 2 at 1 km, 823 g/m 2 at 2.8 km, 511 g/m 2 at 4.8 km, 89 g/m 2 at 6.3 km, 48 g/m 2 at 9.0 km, 25 g/m 2 at 11.9 km, and 8 g/m 2 at 23.3 km (dark blue text in Fig. 2).
The total discharged mass of the 2018 eruption was calculated as 3.4 × 10 4 t (3.4 × 10 7 kg; Fig. 3a) using the segment integration method (Takarada et al. 2001); 2.4 × 10 4 t (2.4 × 10 7 kg; Fig. 3b) using the Weibull fitting method with the calculation parameters θ = 100.15313, λ = 0.29515, and n = 0.73597 (see Bonadonna and Costa 2012 for details); and 2.7 × 10 3 t (using the 5 g/m 2 isomass contour) to 9.4 × 10 4 t (using the 45,000 g/m 2 isomass contour) (Table 1) using the empirical formula developed by Hayakawa (1985). The mass obtained using the Weibull fitting method was ~ 70% of the mass obtained using the segment integration method. The fitting line calculated using the Weibull method agreed with the values measured in the vicinity of the source craters but was much lower than the values measured in the distal areas. This discrepancy in the distal areas may have caused an underestimation of the total discharged mass. In contrast, for the near-source region, the mass may have been overestimated, because the thickness of the mixture of volcanic blocks/lapillus and ash deposits was converted to MPUA and used for the calculations. The mass discharged from the SC was calculated as 175 t using the segment integration method, which accounts for only 0.5% of the mass discharged by the 2018 eruption. The reason for the 30-fold difference in the calculation results using the empirical formula developed by Hayakawa (1985) is that the value of TA is not constant. For the same reason, the method of Legros (2000) (V = 3.69 TA) should be applied to small eruptions with caution, because the results can vary widely depending on which contour is used. As Takarada et al. (2016) noted, the single isopach method (Hayakawa 1985;Legros 2000) may not be suitable for small-scale phreatic eruptions. As such, the results obtained using the single isopach method have been omitted from the remainder of this paper. The magnitude (M = log 10 (erupted mass, kg) − 7; Hayakawa 1993) of the 2018 eruption calculated using the values obtained from the segment integration and Weibull methods was 0.4-0.5. The 2018 eruption began at 10:02 and continued until 10:10 . The duration of the eruption was 8 min, and the total mass discharged was approximately 2.4 × 10 4 -3.4 × 10 4 t. Therefore, the flux rate of the eruption was approximately 5.0 × 10 4 -7.0 × 10 4 kg/s. This value is similar to or slightly larger than the value obtained for the 2014 Ontake eruption (4.2 × 10 4 kg/s; VEI = 2, M = 2.0; Takarada et al. 2016).
One of the most important conclusions of our calculations is that most of the fall deposits generated by the 2018 eruption were deposited near the source craters. The proximal deposits accounted for approximately 50-70% of the total discharged mass (~ 2.4 × 10 4 t using the segment integration method and ~ 1.2 × 10 4 t using the Weibull fitting method, Fig. 3) from the 2018 eruption. Similar results have been obtained for other phreatic eruption deposits in Japan, including the 2014 Ontake eruption (Takarada et al. 2016) and the 2008 Shinmoedake eruption (Geshi et al. 2010). Geshi et al. (2010) found that more than 70% of the fall deposits generated by the phreatic eruption at Mt. Shinmoedake (Kirishima Volcano, Kyushu) on August 22, 2008, was deposited within 1 km of the crater, which they attributed to the low eruption plume height. As in the 2008 Shinmoedake eruption, ~ 70% (Fig. 3a) of the total mass ejected during the 2018 Kusatsu-Shirane eruption was deposited as relatively coarse pyroclasts in the proximal region. The distribution of fall deposits near the source craters is consistent with the low altitude of the 2018 eruption plume, which was just above the crater (Sato 2021).

Interpretation of the distribution of fall deposits
Parameters that can affect the distribution of fall deposits include plume height and intensity, wind direction and speed, and vent inclination (e.g., Bonadonna et al. 2015;Carey and Bursik 2015;Houghton et al. 2017); however, the distribution of the fall deposits produced by the 2018 eruption was likely influenced primarily by plume height and wind direction. The 2018 eruption occurred in three newly opened craters, and the isomass contour map (Fig. 2a) indicates that most fall deposits were ejected from the MC. In contrast to the relatively widespread distribution of fall deposits that originated from the MC, the fall deposits that originated from the SC and WC were limited to areas near the source craters (Fig. 2b). Based on the extent of the isomass contours, we infer that the eruption from the MC was the largest, while the eruption from the SC was the smallest in magnitude. Although the eruptions from each crater varied in magnitude and the distribution of the fall deposits from the MC and WC partially overlapped (arrowhead in Fig. 2b), the main axes of the proximal deposits from each crater were likely identical and oriented northeast (Fig. 2b). The direction of the main axes of the proximal deposits around (See figure on next page.) Fig. 3 Total discharged and cumulative mass estimates for the 2018 eruption of Kusatsu-Shirane Volcano. a Log-log plot of mass/area (g/m 2 ) versus area (m 2 ) of the isomass map and cumulative mass (%) using the segment integration method (Takarada et al. 2001). The regression lines are subdivided into eight segments (shown in different colors), and each segment was used in the calculation. The estimated masses of each segment are shown in the lower part of the figure. Black solid circles indicate field data. b Semilog plot of mass/area (g/m 2 ) versus square root area (km) of the isomass map and the cumulative mass (%) using the Weibull fitting method (Bonadonna and Costa 2012). Black solid circles and dashed lines indicate field data and the fitting curve, respectively Kametani et al. Earth, Planets and Space (2021)  ). The distribution limit of the volcanic blocks also extended from the craters in a NE direction and overlapped with the distribution of the fall deposits near the source (Fig. 1a), which indicates that the trajectories of the volcanic blocks (which were mostly < 20 cm in diameter) were also affected by near-surface wind advection. In contrast, the direction of the main axis of the distal deposits shifted eastward (Fig. 2a), which coincides with the wind direction at high altitudes. The Japan Meteorological Agency radiosonde data indicate a 20-30 m/s westerly wind at high altitudes (2000-5000 m a.s.l at WJM and TTN; Fig. 1a) prior to the 2018 eruption (9:00 am) (Japan Meteorological Agency 2021).

Comparison with historical phreatic eruptions of Kusatsu-Shirane Volcano and other volcanoes in Japan
The total discharged masses, VEIs, and magnitudes of historical eruptions at Kusatsu-Shirane Volcano are presented in Table 2. For the 1982and 1937-1939eruptions, masses were described (Aramaki and Hayakawa 1983Minakami et al. 1942);however, for the 1976, 1932, 1902, and 1882 eruptions, only the volumes of the eruption products were described (Shimozuru et al. 1978;Hayakawa 1999; Japan Meteorological Agency 2013); thus, we converted these volumes to masses using a density of 700 kg/m 3 , which is the density used elsewhere in this study. The minimum and maximum masses of these eruption products were 2.8 × 10 3 -3.4 × 10 3 t and 5.0 × 10 6 t, respectively. These values were not calculated using the same method, because the methods used for the 1982, 1976, and 1932 eruption calculations were not described in their respective manuscripts. The 1902 and 1882 eruption products were calculated using the method of Hayakawa (1985), while the 1937-1939 eruption products were calculated using the source crater topography. Therefore, the results of the calculations using different methods should be used cautiously for comparisons, but the discharged mass, VEI, and magnitude of the 2018 eruption were relatively small compared to other eruptions observed over the past 140 years at Kusatsu-Shirane Volcano (Table 2). Compared to other phreatic eruptions that have occurred in Japan during the last 20 years (Table 2), the 2018 eruption was slightly larger than the 2008 Meakandake eruption (1.2 × 10 4 t; Ishimaru et al. 2009) but smaller than the 2008 Shinmoedake eruption (2.0 × 10 5 t; Geshi et al. 2010), the 2014 Ontake eruption (8.9 × 10 5 -1.2 × 10 6 t; Takarada et al. 2016), and the   (2013) 1902 Kusatsu-Shirane 5.0 × 10 5 1 1.7 Hayakawa (1999) 1882 Kusatsu-Shirane 5.0 × 10 6 2 2.7 Hayakawa (1999) 2015 Kuchinoerabujima eruption (1.0 × 10 5 t; Tajima et al. 2015).

Conclusions
Constructing the detailed isomass contour map of the 2018 phreatic eruption at the MPCG was possible, because we conducted field surveys during the snowy season immediately after the eruption and could obtain accurate data, even in distal areas. The detailed distributions of the fall deposits produced by small-scale phreatic eruptions, such as this one, are useful for estimating the distributions and total discharged masses of similar past eruptions. The main dispersal axis of the 2018 eruption was oriented in the NE-E direction. The total mass discharged by the eruption was estimated as 3.4 × 10 4 t using the segment integration method and 2.4 × 10 4 t using the Weibull fitting method. The calculated masses were smaller than those of phreatic eruptions from Kusatsu-Shirane Volcano that occurred during the past 140 years. The values estimated in this study were similar to the preliminary values obtained (3.0 × 10 4 -5.0 × 10 4 t; The Joint Research Team for ash fall from Kusatsu-Shirane 2018 eruption 2018a). This may be because the preliminary data included a value from the location close to the source (collected on the gondola roofs), namely, ~ 0.5 km from the source, and the assumed value near the source was not significantly different from the actual value. In small phreatic eruptions, most of the eruption products are often deposited near the source area, which greatly affects the calculations. Therefore, it is essential to obtain field survey data near the crater for accurate discharged mass estimates. In the case of the MPCG, its thermal activity and eruption history have been poorly clarified, so ski resorts have been built and operated near the crater without adequate safety measures. The 2018 eruption was a prominent example of this lack of safety, as even a small-scale phreatic eruption can cause a disastrous situation in the near-source area due to the ejection of large amounts of pyroclastic material, including volcanic blocks. This feature of phreatic eruptions is important when potential hazards are evaluated for resort areas near active volcanoes. Such resorts may need to conduct a geological survey of the volcanic eruption history, including very small eruptions of less than 1.0 × 10 4 t, if land is used for ski resorts or in some other way near the crater of an active volcano.