Forest damage and its braking effect on the extreme snow avalanche on Mt. Nodanishoji, Japan, in 2021

ABSTRACT

An extreme snow avalanche occurred on Mt. Nodanishoji, Japan, in 2021. This avalanche was the second largest ever documented in Japan and destroyed many trees and structures. Despite the increasing interest in avalanche disaster mitigation effects of forests, there are limited opportunities to obtain actual data sets on avalanches and forests. In this study, we reconstructed the avalanche velocity based on the cedar forest damage and simulated the avalanche flow using the dynamical run-out model Titan2D to reveal the avalanche movement and ascertain the braking effect of the forest on the avalanche. We successfully simulated the whole movement of the extreme avalanche with a horizontal runout distance of 2,800 m from the starting zone at an altitude of 1,700 m to the run-out zone considering the effect of the forest by increasing the frictional parameters in the model. Whereas in the case of no forest, the avalanche spread beyond its actual reach, considering the forest, the run-out distance decreased and was consistent with the observed reach. Comparing the results with and without forests, we ascertained the distinct braking effect of the forest on the extreme avalanche and quantify it in terms of a locally increased friction angle.

KEYWORDS:

• Dry-slab avalanche
• forest damage
• disaster mitigation effect

Introduction

An extreme dry-slab snow avalanche occurred on Mt. Nodanishoji, Gifu prefecture, Japan, on 10 January 2021. Mt. Nodanishoji is located in the northern part of the Hakusan mountain range on the Sea of Japan side of Honshu and receives large quantities of snow due to the influence of the cold winter monsoon caused by high-pressure episodes over Siberia. The annual maximum snow depth was 1.77 m averaged for thirty years (1991–2020) at the nearest site (Shirakawa: 478 m.a.s.l.) of the Automated Meteorological Data Acquisition System (AMeDAS) of the Japanese Meteorological Agency. Our investigation revealed that this avalanche traveled approximately 2,800 m, the second largest on record in Japan, and was categorized as “extremely large” in the European avalanche size classification (European Avalanche Warning Services 2023) and “size 5” in the Canadian classification (McClung and Schaerer 2006). The avalanche destroyed numerous trees and several structures in the run-out zone. We investigated the forest damage and its braking effect on the avalanche in this study.

Forests reduce disasters by preventing the formation of snow avalanches in starting zones and by hindering the progression of flowing avalanches in run-out zones (Brang et al. 2006). The former effect has been explored for a long time (e.g., Ishikawa, Sato, and Kawaguchi 1969; Gubler and Rychetnik 1991; McClung 2003; Gauer 2018), and there has also been an increasing interest in the latter effect, with active research being conducted in recent years (e.g., Teich and Bebi 2009; Feistl et al. 2014). The deceleration mechanism of forests, as well as their protective capacity, differs according to the size of the avalanche, and research on how forests decelerate avalanches was carried out (Feistl et al. 2012; Teich et al. 2012). For small-to-medium avalanches that do not destroy trees, forests work as effective obstacles to decelerate avalanches, and the braking effect was explained by detrainment that extracts the mass of the snow stopped behind trees directly from the avalanche flow (Teich et al. 2014). Here, the release volumes of small-to-medium avalanches were defined to be <10,000 m³ according to the European Avalanche Classification Scale (European Avalanche Warning Services 2023). The detrainment approach was further improved, and a new avalanche–forest interaction model that considered snow properties and forest configurations was developed (Védrine, Li, and Gaume 2022).

However, simulations based on the detrainment approach are limited to small-to-medium avalanches and cannot be applied to large avalanches. As for large avalanches that destroy trees, the avalanche flow is considered to slow down in the forest through more complicated processes than detrainment, interacting with the forest effects on the avalanche flow and the physical processes of the avalanche flow itself, but the processes are not well understood. Trees are easily destroyed by large avalanches, and the protective function of a forest is disabled or at least reduced once a large avalanche has destroyed the forest (Gubler and Rychetnik 1991). Moreover, the energy required to fracture, uproot, and entrain trees was calculated to be small compared with the total potential energy of an avalanche (Bartelt and Stöckli 2001). Therefore, for large avalanches, the braking effect of a forest was considered to be minor (Brang et al. 2006) and was often neglected in the calculations of avalanche dynamics (Feistl et al. 2014).

Nonetheless, a few studies simulated the effect of the forest on large avalanches by a frictional approach that accounts for the resistance of trees by increasing the frictional parameters in the avalanche dynamics program (Christen, Bartelt, and Kowalski 2010; Takeuchi, Nishimura, and Patra 2018). In Norway, the forest has been identified as the problem with one of the greatest potentials for improving avalanche hazard indication maps, and the run-out model considering the braking effect of forests has been developed (Issler, Gleditsch Gisnås, and Domaas 2020; Issler et al. 2023). These frictional approaches have the advantage of simplifying the complicated braking processes of forests for large avalanches. Utilizing this approach, Takeuchi, Nishimura, and Patra (2018) successfully demonstrated a distinct effect of forest decelerating the avalanches and shortening the run-out distance, although many trees were broken in the process. However, this simulation was performed only in the forest, and it was not possible to reproduce the whole process of the avalanche from the occurrence in the starting zone to stopping in the forest, probably because of the complicated flow path. In addition, it is unknown how much the frictional parameter obtained in the previous studies differ depending on avalanches since there are few cases. Although the friction might depend on various elements such as the size of the avalanche, physical characteristics of snow, topography, and vegetation, it is not directly measurable. The only possible way to obtain the friction parameter is to back-calculate well-documented avalanche events. Hence, it is essential to continue efforts to obtain the friction parameter in various avalanche cases.

Field surveys for documenting forest damage caused by avalanches can provide valuable information on the characteristics of avalanches and will improve our understanding of avalanche dynamics and its interaction with forests (Feistl et al. 2015). However, there are few opportunities to investigate large avalanches and their damage to the forest because such events are rare and require extensive fieldwork. In particular, to obtain a data set on snow and forest, it is necessary to visit the survey field several times depending on the snow conditions because the damaged trees are often buried by debris.

This study aims to reveal the movement of the Mt. Nodanishoji avalanche and quantify the forest braking effect on the avalanche. For this purpose, we investigated the extent of the avalanche and tree damage in detail and estimated the thickness of the avalanche flow and its velocity from these data. Then, we numerically simulated the avalanche with the frictional approach by inputting an optimal friction parameter determined based on both the extent of the avalanche and its velocity flowing into the forest. We also evaluated the dependency of the friction parameter on the size of the avalanche and forest condition by comparing it with the previous studies. The Mt. Nodanishoji avalanche is one of the best-suited cases for this purpose because it was the second-largest avalanche recorded in Japan and stopped within the forest after breaking many trees.

In the following sections, we first describe an outline of the Mt. Nodanishoji avalanche. Then we explain the methods of tree damage investigation, the avalanche velocity estimation based on the investigation, and the numerical simulation based on the estimated velocity. The results of the investigation and the avalanche velocity estimation are described next, followed by a discussion on the friction parameter of the avalanche flow and forest braking effect on the avalanche based on the numerical simulation.

Outline of the Mt. Nodanishoji avalanche in 2021

The dry-slab avalanche occurred around 2:30 a.m. (local time) on 10 January 2021. The avalanche was released during heavy snowfall and had multiple starting zones, according to the report of Japan Avalanche Network (JAN 2021). In this study, we treated the two subavalanches that started from the two largest subzones. Hereafter, we refer to the subavalanches from the northern and southern subzones as N- and S-avalanches, respectively (Figure 1). The run-out zones of the two subavalanches appeared to overlap and could not be separated. It was unclear which avalanche occurred first. Hence, we considered the two subavalanches as the first and second parts of the large avalanche. The starting zones were open without trees on an east–northeast-facing slope located at 1,700 m.a.s.l., just below the ridgeline near the peak of Mt. Nodanishoji. The maximum slab thickness appeared to be approximately 2 m by eye (JAN 2021), and judging from the photo, the entire width of each starting zone appeared to be a failure scarp. The avalanche stopped at 710 to 720 m.a.s.l. with a horizontal run-out distance of approximately 2,800 m (Figure 1). The area above about 1,000 m.a.s.l. is steep and small to medium-sized avalanches seem to occur frequently, and there are almost no trees in the avalanche path. Below 1,000 m.a.s.l., the avalanche flowed along the river, and when the slope became smaller around 800 m.a.s.l., it spread out and proceeded to the left along the river and to the right through the cedar forest. In addition to the cedar forest, there were some broad-leaved deciduous trees in the run-out zone (Figure 1, Figure 2(b)).

Figure 1. Topographic map of the avalanche starting and run-out zones on Mt. Nodanishoji. The starting zones of N- and S-avalanches are shown by the red dashed ellipses. The green hatched areas and the yellow areas indicate the cedar forest and the deciduous broad-leaved forest, respectively. The red dots indicate the observed reach of the avalanche, and the red rectangle and triangles in the run-out zone are the investigation plot of a cedar forest and a row of unbroken cedar trees, respectively. Each small red line represents a fallen tree. This figure was prepared by editing the maps of the Geospatial Information Authority of Japan (GSI).

Figure 2. Aerial view of the avalanche starting and run-out zones on Mt. Nodanishoji. The white arrows in (a) indicate the starting zones of N- and S-avalanches, respectively. The red arrows are the principal direction of avalanche flow.

The avalanche destroyed and carried away a bridge crossing the run-out zone, toppled the poles of a power line, and broke or uprooted many trees (Figure 3). Fortunately, nobody was injured.

Figure 3. Avalanche damage to cedar forest in the run-out zone in 2021. The arrows in the right photo indicate the avalanche starting zones.

Intensive snowfall events had been observed around Mt. Nodanishoji since December 2020. The total precipitation from 30 December to 2 January and from 7 January to the time of the avalanche release on 10 January was 103 and 140.5 mm, respectively, at the nearest AMeDAS site (Shirakawa). We could not observe the snowpack properties and their vertical profiles in the starting zone just after the avalanche, due to bad weather and the risk of avalanches. Therefore, the properties of the snowpack in the starting zone of the avalanche were estimated using the numerical snowpack model, and it was presumed that two heavy snowfall events before the avalanche generated two different weak layers made of precipitation particles and slabs above the weak layers (Katsuyama et al. 2023). We determined the reach of the avalanche through a visual survey in the run-out zone, and the positions were measured by a handheld global navigation satellite system device in the snowmelt season. The mass of debris was estimated by photogrammetry using an unmanned aerial vehicle. The area in the run-out zone and the deposited volume and mass of the debris were more than 32 ha, more than 5.0 × 10⁵ m³, and 2.7 × 10⁵ tons, respectively (Katsushima et al. 2021). This was the second-largest avalanche that was ever recorded in Japan.

Methods

Tree damage investigation

Tree damage is effective information to estimate the avalanche flow and its velocity (Takeuchi et al. 2011; Feistl et al. 2015). Many trees in the run-out zone of the avalanche were damaged, and it was not possible to investigate all of them. Hence, to obtain the extent of tree damage, we selected an investigation plot (10 m wide and 160 m long) parallel to the direction of avalanche flow in a cedar forest, planted in the 1950s, in April, after the snow had almost disappeared. The plot was approximately flat and its upper edge was located at around 760 m.a.s.l. (Figure 1). Cedar (Cryptomeria japonica) is an evergreen conifer that is widely planted in Japan. We evaluated the damage status of every tree in the plot. First, the damage status such as leaning, uprooting, broken, or uninjured was visually determined. Then, we measured the heights of the tree, the uppermost broken branch, and the broken trunk, the diameter of the trunk at the breaking height, the diameter at breast height (DBH), the direction of the leaning or fallen trunk, and the positions of the tree’s root and the broken trunk (Figure 4, Table 1). The height of the uppermost broken branch is an indicator of the avalanche height (thickness). The position of each tree was determined by measuring its distance and direction from the reference tree using a laser distance meter and a compass, respectively.

Figure 4. Measured quantities in the cedar forest.

Table 1. Damage status and measurement items.

Items		Tree height	Diameter at breast height	Height of the broken branch	Breaking height	Breaking diameter	Fallen direction
Upright tree	Uninjured	○	○	—	—	—	—
	Broken branch	○	○	○	—	—	—
Leaning		○	○	○	—	—	○
Uprooted		○	○	—	—	—	○
Fallen trunk with stump		○	—	—	○	○	○
Stump without fallen trunk		—	△	—	○	○	—
Fallen trunk without stump		—	—	—	—	○	○

Note: The circles represent items that were measured, and the triangle represents one that was measured if possible.

There were some cedar trees with the bark partially rubbed away by the avalanche in the lower part of the investigation plot (Figure 5). Because the highest and lowest levels of the peeled bark were considered indicative of the heights of a dense-flow layer of the avalanche and snow surface, respectively, we measured both levels.

Figure 5. A cedar stand with the bark rubbed away by the avalanche.

Although not included in the investigation plot, there was a row of relatively thick and tall cedar trees whose trunks remained unbroken, so we also measured the heights of the tree and the uppermost broken branch and the DBH of these twenty-three cedar trees. The row of cedar trees can be seen in Figure 2, and their positions are shown in Figure 1.

Avalanche velocity estimation

We estimated the avalanche velocity assuming the following two layers: a dense-flow layer at the bottom with a high density and a powder snow layer with a low density, based on the tree damage. The avalanche velocity was estimated from the bending stress of the broken cedar trunks based on on-site measurements. Here, we assumed the avalanche broke a trunk when the bending stress exceeded the tree’s modulus of rupture (MOR), which depends on the tree species and varies due to genetic factors and environmental conditions even within the same species. Hence, to obtain MOR excluding such dependencies, we measured the Young’s modulus of the cedar stands in the investigation plot by a nondestructive method using strain gauges (Koizumi 1987, 1988). We obtained an MOR of 49.7 × 10⁶ Pa (SD = 10.6 ×10⁶ Pa) for the cedar stands in this region using the following relationship because the MOR (${\sigma }_b$) correlates closely with the Young’s modulus ($\mathit{E}$, in Pa; Sawada 1985):

(1) ${\sigma }_b=6.96\times {10}^{-3}\times \mathrm{E}.$(1)

Assuming that the trunk is a column with diameter $\mathit{D}$ (m), the maximum bending stress, $\mathit{\sigma }$ (Pa), is expressed as

(2) $\mathit{\sigma }=\frac{M}{Z},$(2)

where $\mathit{M}$ (N m) and $\mathit{Z}$ (m³) are the bending moment and the section modulus, respectively. The section modulus of a circle is

(3) $\mathrm{Z}=\frac{\pi }{32}{D}^3.$(3)

The force of the avalanche on the trunk, $\mathit{F}$ (N), is expressed as

(4) $\mathrm{F}=\frac{1}{2}{C}_d\mathit{\rho A}{v}^2,$(4)

where $C_d$ is the drag coefficient, $\mathit{\rho }$ (kg m⁻³) is the avalanche density, $\mathit{v}$ (m s⁻¹) is the avalanche velocity, and $\mathit{A}$ (m²) is the area of the vertical section. The drag coefficient was set to 1.2, following Kashiyama (1967), Johnson, Ramey, and O’Hagan (1982), and Shi-igai (1993). Assuming a two-layer avalanche, the sum of the bending moments from the two layers, $M_1$ (N m) and $M_2$ (N m), caused by forces $F_1$ (N) and $F_2$ (N), operates on the breaking height, $h_b$ (m), of the trunk (Figure 6(a)). For the broken trunks, the force applied to the tree crown (F₃) was not included because we could not obtain any data relevant to the size of the crown (e.g., the original tree height and DBH), and $\mathit{A}$ was regarded as $\mathit{Dh}$ $(\mathit{h}$ [m] is flow thickness). The avalanche was assumed to flow on the snow surface, whose height, $h_s$, was set to 1.7 m based on the lowest level of peeled bark; details are described in the section on estimation of avalanche layer thickness. The force on trees by the snowpack support was negligible. Then, from Equation (4)(4) $\mathrm{F}=\frac{1}{2}{C}_d\mathit{\rho A}{v}^2,$(4) , the bending moments operating at the breaking height are expressed as

(5) $\left[\begin{array}{c}{M}_1\\ M_2\end{array}\right]=\left[\begin{array}{c}\frac{1}{4}{C}_d{\rho }_1\mathit{D}{h}_1{v}^2\left(h_1+2h_s-2h_b\right)\\ \mathit{\hspace{0.17em}}\\ \frac{1}{4}{C}_d{\rho }_2\mathit{D}{h}_2{v}^2\left(h_2+2h_1+2h_s-2h_b\right)\end{array}\right],$(5)

Figure 6. (a) Schematic of load on trunks from an avalanche. (b) The relationship between tree height and height of the avalanche.

where $h_1\left(=2.1\mathit{m}\right)$ and $h_2$ are the thickness of the dense-flow and the powder snow layer, respectively, which is obtained from the difference between the highest and lowest levels of the peeled bark and the height of the broken branch, $\mathit{H}$ (m). The densities of the dense-flow layer ${\rho }_1$ and the powder snow layer ${\rho }_2$ were unknown, and they were assumed as 300 and 3 (kg m⁻³), respectively (Issler 2003); we also tried calculations with ${\rho }_1$ set to 200 and 400 (kg m⁻³) and ${\rho }_2$ set to 10 (kg m⁻³).

For the lower part of the investigation plot where no trunk rupture was observed, DBH values of leaning and upright trees were used instead, and the bending moments operating at the ground level ($h_b$ = 0) were calculated to estimate the avalanche velocity. In these calculations, the bending moment $M_3$ (N m) caused by the force $F_3$ (N) applied to the tree crown was also added depending on the tree’s height (Figure 6(b)). The area of the vertical section of the crown was estimated assuming the isosceles triangle crown section and the crown spread $\mathit{W}$ (m) and crown length $\mathit{L}$ (m) were roughly estimated from DBH (m) using the empirical formulas (6) and (7) (Ishida 2016):

(6) $\frac{W}{2}=6.06\cdot \mathit{DBH}+0.216$(6)

(7) $\frac{L}{\mathit{HT}}=0.726\cdot \mathit{DBH}-0.0126$(7)

where $\mathit{HT}$ (m) is the tree height. The flow density upper than the powder snow layer was regarded as the same as air density (1.2 kg m⁻³). In the case of (i) in Figure 6(b), $\mathit{HT}-\mathit{L}\ge \mathit{H}$, assuming the position of the center of gravity is half of the layer thickness $h_3$, the bending moment $M_3$ is expressed as

(8) $\begin{array}{rl}& M_3=\frac{1}{4}{C}_d{\rho }_3{v}^2\left\{\frac{\mathit{WL}}{2}+\left(\mathit{HT}-\mathit{L}-\mathit{H}\right)\mathit{D}\right\}\\ & \left(\mathit{H}+\mathit{HT}-2h_b\right).\end{array}$(8)

In the case of (ii) and (iii), the force applied to the crown was taken into consideration in the same way, but there was no case of (iv) in the investigation plot.

Simulation methods

The avalanche flow was numerically simulated using the dynamics program Titan2D v3.0.0 (Pitman et al. 2003; Patra et al. 2005) to assess the braking effect of the forest on the avalanche. Titan2D is an open-source code that was developed to simulate geological mass flows over natural terrain, and it has also been used in snow avalanche simulations (Mori et al. 2018; Takeuchi, Nishimura, and Patra 2018). This model solves the depth-averaged mass and momentum conservation equations for a shallow granular flow representing the dense flow. It assumes the bed friction to follow the Mohr-Coulomb failure criterion. The internal Coulomb force due to velocity gradients in the direction transverse to the flow is expressed in terms of the internal friction angle φ, but the value of φ was found to exert a minor effect on the simulation result and was set to φ = 20° in our simulations. Details of solution techniques are described in Pitman et al. (2003) and Patra et al. (2005).

Dry snow avalanches exhibit three flow regimes, including the fluidized regime, which has intermediate density and very high mobility (Sovilla, McElwaine, and Köhler 2018; Issler, Gauer et al. 2020), but this flow regime cannot be modeled directly in Titan2D. Also, large dry snow avalanches form a head and tail, and the velocity is maximum near the avalanche front generally (Gubler 1987; Nishimura and Ito 1997), but such velocity and thickness distributions cannot be reproduced by Titan2D: in the simulation by Titan2D, the thickness tends to be larger around the center of the avalanche, and the velocities are also larger there. Therefore, unavoidably, we regarded the velocities in the avalanche to be uniform at the average velocity. The velocity of the avalanche front was also assumed to be the average velocity.

The starting zone of the N-avalanche was approximated as a circular area with a diameter of 300 m. Similarly, the starting zone of the S-avalanche was approximated as an elliptical area with major and minor axes of 360 and 240 m, respectively, based on the topography and the photographs of the starting zones taken immediately after the avalanche (JAN 2021). We confirmed that the position of the starting zone did not affect the simulation result substantially. The maximum initial thickness of each avalanche was regarded as 2 m based on the photographs and the snow thickness above the weak layer estimated by the snowpack model (Katsuyama et al. 2023). Then, the volume of each flux source is approximately similar, 1.4 × 10⁵ m³, and the sum was about 2.8 × 10⁵m³. The fact that the observed deposited volume of more than 5.0 × 10⁵ m³ was nearly twice as large was presumed to be due to entrainment of snow by which the deposited volume was increased compared to the starting one. Because our simulations did not account for entrainment, we increased the starting volumes of both by inputting the maximum initial thickness as 4 m instead. The total volume was increased to 5.6 × 10⁵ m³, exceeding the observed deposited volume. This is plausible because the deposits of large avalanches tend to be much denser than the released slab. We performed both types of simulations with initial thicknesses of 2 and 4 m and the results were compared.

Although Titan2D dynamically obtains the required topographic data by incorporating a direct connection to geographic information system databases, vegetation such as forests is not considered. Therefore, in this study, forested areas were assigned a larger bed friction angle than nonforested areas to mimic the flow resistance due to the trees. Tree resistance is caused by hydrodynamic drag while the tree is intact and thereafter by fracturing the trunk, overturning, and entrainment of tree debris (Bartelt and Stöckli 2001). Although trees decelerate the avalanche through physical processes different from bed friction, we substituted bed friction for tree resistance and assumed that the combined resistance of bed friction and trees was equivalent to the bed friction angle in the cedar forest (δ_f).

The values of bed friction angle in nonforested and forested areas were determined by trial and error as described below. The optimal bed friction angle in a nonforested area (δ) was selected considering the velocity at which the avalanche flows into the investigation plot in the cedar forest and the row of unbroken cedar trees as well as the observed run-out zone.

The bed friction angle in the cedar forest (δ_f) was examined separately using the local simulations limited to the immediate surroundings of the investigation plot for the following two reasons: (1) the distribution of velocity of the snow avalanche flow cannot be reproduced well by Titan2D and (2) the calculation grid of Titan2D becomes coarser as the avalanche spreads (due to dynamic mesh adaptation) and the velocity variation in the investigation plot of the cedar forest was difficult to extract in the whole-scale simulation from the starting zone. The local simulation was started from the upper edge of the forest, giving an initial velocity moving in the direction of the fallen trees (78° from the north). The initial velocity was input as the same velocity at which the avalanche released from the starting zone flowed into the investigation plot by the simulation. The initial pile was approximated as an ellipse with major and minor axes of 30 and 20 m, respectively, to cover the width of the investigation plot. The maximum initial thickness of the pile was input as 2 m based on the thickness of the dense-flow layer. The bed friction angle (δ_f) was selected so that the velocity was consistent with velocities estimated from the bending stress of the cedar trunks and that the run-out distance coincided with the observed run-out.

We mainly investigated damage to the evergreen cedar forest only in the investigation plot in this study. We could not obtain any data on the deciduous broadleaved forest out of the investigation plot and, hence, the numerical simulation could not consider the deciduous broadleaved forest around the cedar forest.

Results

Extent of tree damage

A total of 247 trees within the plot were investigated, of which 198 (80 percent) were cedar trees and the remaining were deciduous broadleaved trees. We primarily used the data from cedar trees to estimate the avalanche flow in this study. As a result of aggregation for cedar trees, we found that 23 trees were uninjured, 43 were upright with only the branches broken, 10 were leaning, 79 were uprooted, 27 were stumps without fallen trunks, 15 were fallen trunks without stumps, and 1 tree was a fallen trunk with a stump. The broken trunks were swept away and piled with the uprooted trees, so it was too crowded and not possible to find their original stumps. This is the reason for the difference in the number of stumps without fallen trunks and fallen trunks without stumps. The stand density of the cedar forest was 1,144 ha⁻¹ from the number of cedar trees within the plot (183) divided by the area of the investigation plot (10 m × 160 m). The number of fallen trunks without stumps was not included because their original positions were unknown. The height of cedar trees ranged from 4.7 to 35.2 m, with an average height of 18.7 m. The DBH ranged from 0.11 to 0.83 m, with an average of 0.31 m. The trees fell toward 78° from the north on average, indicating the avalanche direction.

The position and damage status of each cedar tree in the investigation plot are shown in Figure 7. The status of the trees showed obvious changes around 80 m from the upper edge of the investigation plot; that is, most trees on the upstream side were uprooted or had trunk ruptures with no upright trees, whereas on the downstream side, almost all but a few leaning trees were upright. This indicates that the avalanche slowed down as it flowed through the forest.

Figure 7. Cedar trees and their damage status in the investigation plot.

Estimation of avalanche layer thickness based on tree damage

The factors required for the calculation of avalanche velocity, such as the snow depth in the plot, the height and thickness of each avalanche layer, the breaking height of trunks, and the trunk diameters at the breaking height, were determined based on on-site measurements. Figure 8 shows the levels of the peeled bark of cedar trees. The lowest levels that indicate snow depth range from 1.3 to 2.3 m, with an average of 1.7 m. The highest levels that indicate the height of a dense-snow layer range from 3.1 to 4.9 m, with an average of 3.8 m. Consequently, the thickness of the dense-flow layer was estimated to be 2.1 m. Considering the possibility that the original snow depth had decreased due to erosion of the avalanche in the run-out zone when estimating the avalanche velocity, we also tried calculating the case where the snow depth was large.

Figure 8. The levels of the peeled bark of cedar trees. The dotted lines indicate the mean value of each level.

The height of the broken branch, which implies the height of the avalanche, is shown in Figure 9, with the tree height also shown for comparison. The first and second highest broken branches were 13.7 and 13.2 m, respectively, whereas most of them coincided with the level of the peeled bark (1.7–3.8 m). In particular, in the lower part of the plot, branches at a high level were not broken. These findings suggest that the avalanche consisted of two layers, viz. a dense-flow layer at the bottom with a high density and a powder snow layer with a low density, and the branches were broken by the dense flow in the lower part of the plot in which the avalanche velocity decreased. Hence, the height of the powder snow layer was regarded as 13.5 m based on the average of the top two heights.

Figure 9. The height of the broken branches. Tree height is also shown for comparison. H = height of powder snow layer; h_s + h₁ = height of dense-snow layer; h_s = height of snow.

Estimation of avalanche velocity in the cedar forest

The avalanche velocities required to break the cedar trunks were calculated using the breaking diameter and breaking height of each trunk. For the lower part of the investigation plot where no trunk rupture was observed, DBH values were used instead, and the avalanche velocity was estimated from the bending moments operating on the ground level. These diameters for every tree and the breaking height of trunks are shown in Figure 10. The uprooted trees were excepted because they were not used for the avalanche velocity estimation. As mentioned earlier, the status of trees changed around 80 m from the upper edge, with no upright trees observed on the upstream side, whereas almost all trees were unbroken on the downstream side (Figure 10(a)). One of the reasons is that the trunks were thick in the 80 to 120 m area from the upper edge, but the reason why even the thin trunks did not break could be because the avalanche slowed down and its destructive force decreased.

Figure 10. (a) The breaking diameter and DBH of leaning and upright trees. (b) The breaking height of the trunks. Height of snow estimated from the level of the peeled bark is also shown for comparison. Note: We counted twenty-seven stumps without fallen trunks but only twenty-six are shown because we could not access one stump and obtained no data on its tree.

The breaking height of the trunks was lower than the height of snow; that is, all trunks were broken in the snowpack (Figure 10(b)). If the trunk diameter were constant, it should break at the ground where the bending moment is maximum. However, in reality, because the lower the tree trunk is, the thicker it is, the height at which the bending stress was maximum was determined due to the individual tree’s shape. Hence, the breaking heights differed among the trees.

The avalanche velocities estimated from the bending stress of broken trunks are shown in Figure 11. The velocities required to break the unbroken trees at the ground level are also depicted in the figure. For the broken trunks, their DBH was difficult to estimate and the force applied to the tree crown could not be taken into account. Considering the tree crown, the velocities required to break trunks were reduced. The rate of decrease varied depending on the tree size (trunk diameter and height of each tree). In calculations for unbroken trees, the velocities required to break trunks with a diameter of about 0.3 m were reduced by about 10 percent on average when the tree crown was taken into account compared to the calculations without a tree crown. Accordingly, in Figure 11, the range of 10 percent decrease is expressed by each error bar. The avalanche velocity required to break the trunks was estimated to be at least 20 to 22 m s⁻¹. It can be inferred that the avalanche flowed into the cedar forest at a speed of >20 to 22 m s⁻¹ and flowed approximately 80 m through the forest while breaking the trunks at a velocity higher than the solid line shown in the figure and subsequently moved >60 m at a velocity lower than the dashed line. The end of the investigation plot was 160 m from the upper edge, and even slight traces of the avalanche disappeared beyond there; hence, the estimation of the avalanche velocities is considered to be reasonable. Based on the critical turning moment of cedar trees from Todo et al. (2015), the required velocities for uprooting were estimated to be lower than those for breaking trunks, and the velocity of avalanche flowing into the cedar forest inferred above was not affected even if the uprooted trees were considered.

Figure 11. Avalanche velocities estimated from the bending stress of broken trunks (×). The range of 10 percent decrease is expressed by each error bar. The velocities required to break the upright trees at the ground level are also shown (○).

In the estimation of the avalanche velocities, the densities of the dense-flow layer ${\rho }_1$ and the powder snow layer ${\rho }_2$ were assumed as 300 and 3 (kg m⁻³), respectively, and we tried calculating by changing the densities. When the density ${\rho }_1$ was changed to 200 (kg m⁻³), the required avalanche velocities increased by about 10 to 20 percent, and when ${\rho }_1$ was changed to 400 (kg m⁻³), they decreased by about 10 percent. When the density ${\rho }_2$ was changed to 10 (kg m⁻³), they decreased by about 15 percent on average for trunk-broken trees and 20 percent on average for upright trees.

Assuming that the snow surface had lowered due to erosion, we also tried calculating the avalanche velocities when the snow depth was large; that is, the thickness of the dense-snow layer was small. As a result of increasing the snow depth by 0.5 or 1.0 m, the velocities were increased by about 10 or 20 percent on average for trunk-broken trees, respectively, and 0 or 10 percent for the upright trees.

The avalanche velocities required to break the upright trees in the cedar row a short distance outside the investigation plot were also calculated. There were twenty-three cedar trees, and the direction of the row was almost the same as the avalanche flow. The diameters of the tree trunks were larger upstream of the row, and the required velocities to break the trees were estimated to be 29 to 34 m s⁻¹ for the five trees from the upper edge of the row with an average of 31 m s⁻¹.

Discussion

The starting zones and courses of avalanches

The avalanche flow was simulated using a numerical model. As described before, the thickness of the avalanche in its starting zone was estimated to be 2 m, but we also set a thickness of 4 m in the simulation to make the volume close to the observed deposited volume. Both the gravitational and frictional forces were proportional to the volume in Titan2D, whereas the hydrostatic pressure force grew with the square of the flow depth. Therefore, we expected a larger initial volume to have a minor effect on the velocity and the trajectory of the flow, but it would increase the run-out zone due to pressure-driven spreading. This was borne out by the simulations (see Supplementary Materials, Figure S1). We discuss the braking effect of the forest using the simulation with an initial thickness of 4 m in this study.

The simulation results with bed friction angles of 13° for the N- and S-avalanches are shown in Figure 12. There were almost no trees along the path of the avalanche in areas higher than about 1,000 m.a.s.l. due to the steepness and the frequent occurrence of small-sized avalanches, but we found several fallen trees near the ridge and on the slope in the photographs from the unmanned aerial vehicle or Google Maps. Each red cross in the figure represents the position of a fallen tree. We assumed that the trees were felled by the avalanches, though we could not investigate them on-site. However, according to the simulations, the N-avalanche did not pass through the ridge marked with ellipse A, and the S-avalanche did not pass through the slope marked with ellipse B. The results for δ of 12° and 14° were similar. Accordingly, we concluded that both avalanches must have felled trees in stands A and B.

Figure 12. Simulation results with bed friction angle (δ) of 13° for the (a) N-avalanche and (b) S-avalanche, respectively. The extent of an avalanche with a thickness of >0.05 m is shown by the blue at 20 to 40 seconds after starting. Each red cross represents the position of a fallen tree. The results with δ of 12° and 14° are similar. The area shown in this figure is indicated in the map in Figure 1

Bed friction angle of the avalanches

When two avalanches occur in succession, which one has smaller bed friction depends on the situation. However, regarding the N- and S-avalanches, we assumed that the bed friction was smaller in the subsequent avalanche. We could not obtain data on which avalanche occurred first. The volumes of the avalanches were probably similar because both starting zone areas were similar. The slopes where the avalanches occurred lack vegetation and are prone to erosion, so a large amount of sediment is generated during the season without snow. Therefore, three large erosion control dams have been built in the river where the avalanches passed. Also, because the avalanches occurred relatively early in winter, the shrubs in the avalanche path may not yet have been buried by snow. Thus, the avalanche that occurred first met large resistance by these dams and shrubs. The subsequent avalanche was considered to flow while meeting smaller resistance after the dams and shrubs had been buried by the previous avalanche. Accordingly, it is plausible that the bed friction was smaller in the subsequent avalanche because of these events.

The velocity variations of the N- and S-avalanches with bed friction angles (δ) of 12°, 13°, and 14° from start to deposit are shown in Figure 13. The avalanche reached the upper edge of the investigation plot at a time between the two plus marks and reached the row of unbroken cedar trees at a time between the two circles. The black dashed lines indicate the velocity required to break the cedar trunks in the investigation plot (20 to 22 m s⁻¹; see Figure 11) and the red dashed line indicates the velocity required to break the trunks in the row of cedar trees (29 m s⁻¹). It was found that the N-avalanche reached the row of cedar trees at the same time or immediately after reaching the investigation plot, whereas the S-avalanche reached there before reaching the investigation plot. The N-avalanche with δ of 12° was found to flow into the investigation plot and the row of cedar trees at the same time at the velocity (28 m s⁻¹) between the red and black dashed lines: it could destroy the cedar trees in the investigation plot but not the cedars in the row. When δ was 13°, its velocity flowing into the investigation plot (22 m s⁻¹) was the same as the black dashed line, and the avalanche could barely fell the stand. Consequently, assuming that the N-avalanche destroyed the cedar forest, the optimal δ was 12° to 13° in our simulation. On the other hand, it was found that the δ of the S-avalanche should be larger than 12° because, with a δ of 12°, its velocity upon reaching the row of cedars was larger than the red dashed line: it was possible to fell the stands in the cedar row, contrary to the actual state. When δ was 13°, it could barely fell the stands in the investigation plot but not in the cedar row, and when δ is larger than 13°, the cedar forest was left intact.

Figure 13. Variations in the average velocity of the avalanches. We unavoidably regarded these averages as equal to the velocities in the front of avalanches, because Titan2D cannot reproduce the velocity distribution of avalanches plausibly. The blue, black, and green lines represent the bed friction angles (δ) of 12°, 13°, and 14°, respectively. The plus marks and circles indicate the velocity at which the avalanche flowed into the investigation plot and reached the unbroken row of the cedar stands, respectively. The black dashed lines indicate 22 and 20 m s⁻¹. The red dashed line indicates 29 m s⁻¹. (a) N-avalanche and (b) S-avalanche.

The simulation results at the time just before and after each avalanche reached the investigation plot are shown in Figure 14. As described above, we found that the N-avalanche reached the row of cedars immediately after the investigation plot, whereas the S-avalanche reached the row of cedars before reaching the investigation plot. The flow direction of both avalanches coincided with that of the fallen direction of trees in the plot. The figure shows the case for δ = 13°, though the result for δ = 12° was similar.

Figure 14. The simulation results at the time just before and after the avalanche reached the investigation plot (the red rectangle). The red triangles in the center are the row of unbroken cedar stands. The arrows indicate flow direction. (a) The extent of the N-avalanche at 54 and 57 seconds after starting with a bed friction angle of 13°. (b) The extent of the S-avalanche at 54 and 58 seconds after starting with a bed friction angle of 13°.

The run-out zone of the avalanche simulation is shown in Figure 15. Comparing the reach of the left arm of the N-avalanche, which was not affected by the cedar forest, with the observed run-out, it was found that δ = 12° to 13° is appropriate (Figure 15(a-1), Figure 15(a-2)). For the S avalanche, δ = 13° was found to be too small because the left arm could exceed the observed run-out, and δ = 13° to 14° was assumed to be reasonable (Figure 15(a-3), Figure 15(a-4)). Consequently, the optimal δ of the N-avalanche was smaller than that of the S-avalanche, and it is unlikely that the S-avalanche destroyed the cedar forest.

Figure 15. (a) The run-out zone of the avalanche simulations without forest. The blue area shows the extent of an avalanche with a thickness of >0.05 m. The red dots represent the observed avalanche reach. (a-1) and (a-2) are the results of N-avalanche simulations with bed friction angle (δ) of 12° and 13°, respectively, for the entire area assuming no forest. (a-3) and (a-4) are the results of S-avalanche simulations with δ of 13° and 14°, respectively, for the entire area. (b) The run-out zone of the avalanche simulations with cedar forest (the hatched areas). (b-1) and (b-2) are the results of N-avalanche simulations with δ of 12° and 13°, respectively, in the nonforested area and with δ_f of 27° and 19°, respectively, in the cedar forest. (b-3) and (b-4) are the results of S-avalanche simulations with δ of 13° in the nonforested area and with δ_f of 27° and 19°, respectively, in the cedar forest.

Although we could not obtain any data indicating which avalanche destroyed the cedar forest, we conjecture based on the simulations that it was the N-avalanche. In that case, there are two possible scenarios. In scenario 1, the S-avalanche occurred first and left the cedar stand intact, which was subsequently felled by the N-avalanche. In the second scenario, the N-avalanche occurred first and felled the cedar stand. The S-avalanche released thereafter and its right-hand side flowed over the already destroyed and probably buried cedar trees. We do not have pertinent observations indicating which of the two scenarios is more likely. However, assuming that the bed friction was larger in the previous avalanche as described above, the S-avalanche was supposed to precede; that is, the first scenario is more likely.

These δ values of 12° to 14° were almost similar to but a little smaller than those of the 2008 avalanche in the Makunosawa valley, for which the bed friction angle in the nonforested area was 13° to 14° (Takeuchi, Nishimura, and Patra 2018). This is plausible because both avalanches were large dry-slab avalanches and occurred with similar air temperatures, and the elevations, slope directions, and slope angles of the starting zones were also similar. Because the Mt. Nodanishoji avalanche was larger and two subavalanches occurred in succession, δ became a little smaller and it is also consistent with the trend of δ decreasing with increasing avalanche size (Izumi 1985).

Bed friction angle in the forested area

As shown in Figure 13(a), the N-avalanche’s velocity flowing into the investigation plot of the cedar forest was simulated to be 28 and 22 m s⁻¹ for the bed friction angles (δ) of 12° and 13°, respectively. Accordingly, a local simulation in which avalanches initiated from the position of the upper edge of the investigation plot at initial velocities of 28 and 22 m s⁻¹ was carried out. When the bed friction angle in the forested area (δ_f) was set to 26° to 27° (18°–19°), the avalanche, which initiated at an initial velocity of 28 (22) m s⁻¹, was found to slow down as estimated from the bending stress of trees (Figure 16). The avalanche was simulated to flow in the first half at a velocity approximately higher than the cross marks shown in the figure and slow down to a velocity lower than the circles in the second half. It terminated about 100 to 114 m from the upper edge of the forest, within the investigation plot. Accordingly, the bed friction angle in the cedar forest (δ_f) was selected to be 26/27° to 18/19° in our simulation. This δ_f value was similar to that of the 2008 avalanche in the Makunosawa valley, for which the bed friction angle in the cedar forest was 25° (Takeuchi, Nishimura, and Patra 2018). The difference is probably due to the differences in the forest conditions such as stem diameter, root strength, and stand density. To determine the bed friction angle in forests based on the forest conditions, we still need to acquire more data covering a wider range of avalanche parameters and forest conditions.

Figure 16. Comparison between the simulated velocity with bed friction angles of 27°/26° and 19°/18° (black and blue solid lines) and 12° and 13° (dashed black and blue lines). Velocities estimated from the bending stress of broken trunks (×) and the velocities required to break the upright trees at the ground level (○) are also shown for comparison. Each error bar expresses the range of a 10 percent decrease. The end of the investigation plot is shown by the arrow (↑).

Braking effect of the forest on the avalanche

To quantify the braking effect of the cedar forest on the avalanche, we initiated the avalanche at an initial velocity of 28 (22) m s⁻¹ from the position of the upper edge of the forest with a constant bed friction angle of 12° (13°) assuming no forest. Our result showed that the avalanche decelerated slowly and flowed farther than the observed stopping point (Figure 16). The avalanche would have traveled 30 to 270 m farther than the observed run-out if the cedar forest had not existed.

Furthermore, we compared the extent of the N-avalanche that flowed down from the actual starting zone at an altitude of 1,700 m with and without forest in the runout zone. Titan2D v3.0.0 can be run with a spatially variable friction angle. Therefore, to confirm the braking effect of the cedar forest, the extent of the N-avalanche was simulated by setting different bed friction angles in the nonforested area and the cedar forest. As described previously, the bed friction angle δ was 12° (13°) for the former, and δ_f was 27° (19°) for the latter. The hatched area in Figure 15(b) indicates the cedar forest with the larger bed friction angle. The results showed that the reach of the N-avalanche coincided well with the actual reach (Figure 15(b-1), Figure 15(b-2)), and the braking effect of the cedar forest was ascertained by comparing the simulations assuming no forest with a constant δ = 12° or 13° for the entire area (Figure 15(a-1), Figure 15(a-2_). Whereas in the case of no forest, the avalanche spread beyond its actual reach, considering the cedar forest, the run-out distance decreased and came to be consistent with the observed run-out. Figure 15(b-1) and Figure 15(b-2) show the case with a δ_f value of 27° or 19°; a very similar result was obtained using δ_f of 26° or 18°. In the case of the S-avalanche, the cedar forest seems a little less effective and it was found that δ_f = 19° is too small and δ_f = 27° was necessary at least for the one with δ = 13° (Figure 15(b-3), Figure 15(b-4)). When δ is close to 14° or δ_f is close to 27°, the S-avalanche may not reach or be close to the observed run-out, and there is no contradiction to the first scenario described above. Comparisons of the forest effect in the temporal changes in the N- and S-avalanche extents in the run-out zone are shown in Supplementary Figure S2.

A first attempt at quantifying the braking effect of forests

To create realistic avalanche hazard indication maps, the braking effect of forests should be taken into account in a run-out model in any avalanche-prone country with forests. In this process, it is essential to quantify the braking effect of forests based on the forest characteristics (tree density, trunk diameter, tree species, etc.) in the relevant area. We tried to quantify the braking force exerted by the trees and compared them with the extra contribution to the bed friction coefficient obtained by Titan2D.

As a working hypothesis, one may assume that the retarding force per unit forested area is proportional to the mean drag force on a tree ($F_t$) times the number of trees per unit area ($\mathit{n}$) while the trunks are upright. This neglects possible interactions between the flow disturbances from neighboring trees. With this, the equivalent averaged bed shear stress (${\tau }_f$) follows from Equation (4)(4) $\mathrm{F}=\frac{1}{2}{C}_d\mathit{\rho A}{v}^2,$(4) as

(9) ${\tau }_f=F_t\mathit{n}=\frac{1}{2}{C}_d\mathit{\rho Dh}{v}^2{n}_{\mathrm{hor}}cos\mathit{\theta }.$(9)

The number of trees per unit horizontal area ($n_{\mathrm{hor}}$) was converted by the factor $cos\mathit{\theta }$ to the tree density per unit oblique area (n), but the cedar forest in the run-out zone of Mt. Nodanishoji avalanche is approximately flat ($\mathit{\theta }\approx 0^{\circ }$), so we assumed that $\mathit{n}\approx n_{\mathrm{hor}}$. This force gives an enhancement of the effective friction coefficient $\mathit{\mu }\left(=tan\mathit{\delta }\right)$ of

(10) $\mathrm{\Delta \mu }=\frac{{\tau }_f}{\mathit{\rho gh}cos\mathit{\theta }}=\frac{1}{2\mathit{g}}{C}_{\mathit{d}}\mathit{Dv}{ }^2{n}_{\mathrm{hor}}.$(10)

We tried to obtain $\mathrm{\Delta }\mathit{\mu }$ based on the data on stands observed in the investigation plot in the run-out zone. Calculations were performed separately for the upstream area (0–80 m from the upper edge), where almost all trees were broken or uprooted, and the downstream area (80–160 m), where almost all but a few leaning trees remained upright. $\mathrm{\Delta \mu }$ largely depends on velocity. Here, we calculated $\mathrm{\Delta \mu }$ by using the velocities of the solid and dashed lines in Figure 11, depending on the tree’s position (distance from the forest edge). The number of uprooted trees and other trees than cedar was included in these calculations. As a result, we obtained $\overline{\mathrm{\Delta \mu }}\approx 1.1$ for the upper half of the plot and $\overline{\mathrm{\Delta \mu }}\approx 0.1$ for the lower half (Supplementary Table S1), whereas the difference between the bed friction coefficients for the forested and nonforested areas was 0.3 ($tan{27}^{\circ }-tan{12}^{\circ }$) or 0.1 ($tan{19}^{\circ }-tan{13}^{\circ }$) from the simulation and was close to the $\overline{\mathrm{\Delta \mu }}$ in the downstream area. In the upstream area, $\overline{\mathrm{\Delta }\mathit{\mu }}$ may have been overestimated considerably because many trees probably were broken very quickly.

Conclusions

Based on the field survey and the numerical simulations, we concluded that the extreme avalanche occurred in the following scenario. The avalanche had two major starting zones, and the S-avalanche occurred first and left the cedar forest unbroken, which was subsequently destroyed by the N-avalanche. The optimal bed friction angles in the nonforested area for the N- and S-avalanches were 12° to 13° and 13° to 14°, respectively. The N-avalanche flowed into the cedar forest at approximately 28 to 22 m s⁻¹ and flowed in the forest while breaking many cedar trunks in the first half. Regarding the combined resistance of bed friction and trees as a bed friction angle in the cedar forest, δ_f, the optimal δ_f was found to be 26/27° to 18/19°.

Setting the different optimal bed friction angles to 12° to 13° in the nonforested area and to 26/27° to 18/19° in the cedar forest, we successfully simulated the whole movement of the N-avalanche with the horizontal runout distance of 2,800 m. In contrast, with a constant bed friction angle of 12° to 13° in the entire area assuming no forest, the N-avalanche spread beyond the actual reach. Consequently, the braking effect of the forest effectively mitigated the damage by the extreme avalanche that destroyed trees. This result supports previous research (Takeuchi, Nishimura, and Patra 2018).

In this study, we used Titan2D, a model for dense granular flows. A simplified model like Titan2D has inherent limitations. The most important issue is that Titan2D cannot model flow regimes of snow avalanches, including the fluidized regime, which has intermediate density and very high mobility occurring in dry snow avalanches. Instead, we had to choose very low bed friction like 12° to 13° to capture the run-out distance of the head. Also, distributions of thickness and velocity of the avalanche flow are very different from those of real snow avalanches. Therefore, unavoidably, thickness distribution could not be used, and the velocity was assumed to be uniform at the average velocity in this study. The starting and deposited volumes were almost equal in the simulation because entrainment was not considered. To approximate the observed deposited volume, we simulated the flow with a larger fracture depth than inferred from the photographs and meteorological conditions. In addition to the velocity distribution by Titan2D being different from that of a real snow avalanche, the calculation grid became coarser as the avalanche expanded, so we had to simulate separately to obtain detailed velocity variations in the investigation plot of the cedar forest in the run-out zone. Although the deciduous forest around the cedar forest could also have influenced the movement of the avalanche, the simulation did not account for the forest except the cedar forest.

When the avalanche velocities were estimated from the bending stress of broken trunks, the force applied to the tree crown could not fully be taken into account due to a lack of information. In the estimation of the avalanche velocities, we tried calculating by changing the densities of snow layers and the original snow depth, which might include uncertainties, and confirmed the sensitivities to the calculated velocities; the effects were found to be 10 to 20 percent.

As mentioned above, although there are limitations and uncertainties associated with the model, we could show that the forest had a distinct braking effect even against the extreme avalanche, based on detailed observation data. It is necessary to improve avalanche run-out models so that they can take into account the braking effect of forests, but acquisition of forest data such as stand density, trunk diameter, and tree species through on-site surveys is limited because it requires a great deal of effort and time. Evaluating the forest structure and species composition using remote sensing (Breidenbach et al. 2021; Stritih 2021) is expected to be useful. Some parameters characterizing forest stands based on satellite images and high-resolution aerial LiDaR data are available in Norway and were effective in quantifying the effect of forests over large areas (Issler, Gleditsch Gisnås, and Domaas 2020; Issler et al. 2023). In other countries and regions, the availability of forest data over large areas using remote sensing is necessary and should be effective in reducing avalanche disasters. It is also important to improve avalanche run-out models so that the forest data can be fully utilized.

Acknowledgments

The field investigation was conducted with the cooperation of the Gifu Prefecture Forestry Bureau Takayama Office and the TOYOTA Shirakawa-Go Eco-Institute. The first author thanks Dr. Kae Tsunematsu of Yamagata University for the kind advice on Titan2D. The authors are sincerely grateful to the reviewers for their many useful suggestions and very kind encouragement.

Disclosure statement

No potential conflict of interest was reported by the author(s).

Supplementary material

Supplemental data for this article can be accessed online at https://doi.org/10.1080/15230430.2024.2327652.

Additional information

Funding

This work was supported by the Forestry and Forest Products Research Institute.