Meteorological and snow-cover data

The meteorological data obtained at both sites were air temperature, relative humidity, wind speed, incoming and outgoing shortwave radiation, net radiation (directly measured) and precipitation. Moreover, variations of snow depth and snow water equivalent (SWE) were measured. The parameters were measured automatically at 1 hour intervals. The datasets used in this study were obtained in 1997/98, 1998/99, 1999/2000, 2000/01 and 2001/02 at NISIS and in 2001/02 at SB.

Snow-pit investigations were carried out at 0900 h at 5 day intervals during 2000/01 and 2001/02 at NISIS, and at 10 day intervals during 2001/02 at SB. The elements examined in the snow-pit investigation were snow depth, stratigraphic factors, snow type, snow density, snow temperature and liquid-water content. Liquid-water content was measured using a dielectric probe (Denoth type). The parameters were measured at 5–10 cm intervals.

Simulation conditions

In order to simulate snow characteristics (e.g. snow type, snow temperature, density, and liquid-water content), SNOWPACK requires the following meteorological data: air temperature, relative humidity, wind speed, shortwave radiation, snow surface temperature and/or incoming long-wave radiation, and snow depth and/or precipitation (Reference Lehning, Bartelt, Brown and FierzLehning and others, 2002b).

In this study, the mass added through snowfall was determined based on the variations in measured snow depth. New-snow density (ρ_new) was treated as the following empirical function obtained from data gathered in the northern part of Honshu island, Japan (Reference KajikawaKajikawa, 1989):

where T_A is air temperature and V_W is wind speed.

Snow surface temperature and incoming and outgoing longwave radiation were not measured at either site, but were estimated based on the following methods: When air temperature was <0°C, snow surface temperature was set to the air temperature; when air temperature was >0°C, snow surface temperature was set to 0°C. Since most air temperatures were >0°C at both sites, we expect only a small influence on the simulation results. Incoming longwave radiation was determined using the measured net radiation and the estimated surface temperature by means of Planck’s law.

Although the conventional SNOWPACK model does not require surface albedo values (Reference Lehning, Bartelt, Brown and FierzLehning and others, 2002b), each daily mean albedo estimated from incoming and outgoing shortwave measurements was introduced into the simulation in order to eliminate errors arising from the albedo estimation in this study.

Comparison algorithms

The simulated profile had much higher resolution and finer structure than did the measured profile, and differences were seen between simulated snow depth and measurements. Thus, for a quantitative comparison of measured and simulated profiles, a mapping procedure was required. We calculated agreement scores for all parameters based on algorithms described in Reference Lehning, Fierz and LundyLehning and others (2001), but the following new equation was introduced to calculate scores for snow temperature, density and water contents:

Here S_Xprofile is the agreement score, N is number of data, and X_{obs i} and X_{mod i} are a set of N predicted and measured data pairs, respectively.

Each score shows fluctuation between 1 and 0, where 1 is the maximum and 0 stands for no agreement.

A study focusing on the statistical validation of the SNOWPACK model was performed in an area having a Montana (U.S.A.) climate (Reference Lundy, Brown, Adams, Birkeland and LehningLundy and others, 2001). Similar methods were applied to the standard evaluation of SNOW-PACK. In this study, the following three statistical parameters were used for discussion:

Mean deviation (D), which indicates the model’s tendency towards overestimation or underestimation,

Root-mean-square error (RMSE), which estimates the expected magnitude of error associated with a model’s prediction,

Pearson’s correlation coefficient (r), which indicates positive linear correlation with positive values and negative linear correlation with negative values.

Here, μ_mod and μ_obs are the means of the predicted and measured datasets, respectively. The correlation coefficient has an upper bound of 1 (perfect positive linear correlation) and lower bound of –1 (negative linear correlation).

Snow depth

Figure 1a–c show the snow-depth measurements and simulations for 2001/02 at NISIS, 2000/01at NISIS and 2001/02 at SB. The simulation in 2001/02 at NISIS was in good agreement with the measurements (Fig. 1a). The snow-cover disappearance dates were 11 March in the measurements and 17 March in the simulations, and the maximum snow-depth error was 16 cm on 7 February. The snow-depth error was defined by subtracting measured snow depth from simulated snow depth. On the other hand, the simulations in 2000/01 at NISIS and 2001/02 at SB showed large differences from the measurements (Fig1b and c). Both simulated snow depths corresponded closely with measurements before mid-February, but the snow-depth errors increased during the course of the melt season due to slower decrease of simulated snow depths. The maximum deviation between measured and simulated snow depth was approximately 60 cm at NISIS and 50 cm at SB.

Fig. 1. Comparisons between snow-depth measurements and simulations. Solid lines indicate measurement, and dashed lines indicate simulation. (a) NISIS in 2001/02; (b) NISIS in 2000/01; (c) SB in 2001/02.

Figure 2 shows the relationship between maximum snow depth and maximum snow-depth error in 1997/98, 1998/99, 1999/2000, 2000/01, 2001/02 at NISIS, and in 2001/02 at SB. It is apparent from this figure that maximum error increases with maximum snow depth when SNOWPACK is applied to wet- and heavy-snowfall areas inJapan. It is likely that snow-depth errors significantly affect snow characteristics. Thus, only simulations of snow characteristics in 2001/02 at NISIS were used for a detailed analysis because the simulated snow depths agreed closely with measurements.

Fig. 2. Comparisons between maximum snow depth and maximum snow-depth error.

Quantitative profile comparison

Figure 3 summarizes the individual and overall agreement scores for NISIS in 2001/02. The average total score was 0.75, which is smaller than the scores reported in the European Alps (Reference Lehning, Bartelt, Brown and FierzLehning and others, 2002b). The scores for snow temperature and snow type showed high values and small fluctuations, while those for snow density and liquid-water content showed large fluctuations (Fig. 3; Table 1). These results imply that the methods for simulating density and liquid-water content need to be improved.

Fig. 3. Time series of individual and overall agreement scores between the measurement profile and simulations at NISIS in 2001/02.

In the following, detailed analyses are carried out to account for the differences between simulation and measurements.

Snow type and snow temperature

Figure 4 shows comparisons of snow-type simulations and measurements for winter 2001/02 at NISIS. Measured snow was usually wet even during mid-winter because of melt-water, and the simulation showed similar trends; that is, melt forms appeared in all seasons.

Fig. 4. Comparisons between snow-type measurement and simulations in 2001/02 at NISIS.

One reason for the decreasing agreement score is that the simulated speed of the change of (lightly) compacted snow into melt forms was slower than in the measurements. Another reason is that the model occasionally constructed unrealistic snow types. Faceted particles appeared on 28 December, 15 January and 25 February in the simulations, but such a snow type was never observed at NISIS during winter 2001/02.

Most of the snow-temperature measurements and simulations were 0°C. Negative temperatures were measured on only three days during the snow-pit survey, whereas seven days had negative temperatures in the simulation. Thus SNOWPACK sometimes simulated lower temperatures than those measured.

In the model, the factor of snow-grain roundness (sphericity) is treated as a function of temperature gradient (Reference Lehning, Bartelt, Brown, Fierz and SatyawaliLehning and others, 2002a). When the temperature gradient is <5°Cm^–1, the individual snow grains’ sphericity increases and their shapes become rounder. Conversely, sphericity decreases and snow grains become more angular when the temperature gradient is >5°Cm^–1. Thus the appearance of unrealistically faceted particles is caused by the larger simulated temperature gradients associated with the lower temperatures.

Density and liquid-water content

Table 2 summarizes the results of the descriptive statistical analysis of density and liquid-water content in 2001/02 at NISIS. All measured density values and all measured liquid-water values which were not zero were used to compare with the simulations at the corresponding height.

It is apparent from the RMSE of 93.5 kgm^–3 that SNOW-PACK had difficulty predicting snow-cover density. A mean deviation (D) of –55.7 kgm^–3 indicates that most of the prediction error was due to underestimation of the density.

Figure 5a shows comparisons between the density measurements and simulations. It indicates that simulated densities were overestimated when they were 5200 kgm^–3, and underestimated when >200 kgm^–3. These trends correspond well with the results obtained in the Montana climate (Reference Lundy, Brown, Adams, Birkeland and LehningLundy and others, 2001). One of the reasons why a larger density appeared at a small density may be the overestimation of snow depth. In the model, the amount of new snow is estimated based on the difference between simulated and measured snow depth (Reference Lehning, Bartelt, Brown and FierzLehning and others, 2002b), so overestimation of snow depth causes underestimation of the amount of new snow.

Fig. 5. Comparisons between measurement and simulation in 2001/02 at NISIS: (a) density; (b) liquid-water content.

The comparison of predicted and observed liquid-water content was problematic. The simple water-transport routine employed in SNOWPACK works with a fixed volumetric residual water content of 0.025 (Reference Bartelt and LehningBartelt and Lehning, 2002), which results in a poor description of quantitative values of layer liquid-water percentage per weight. The water-transport routine is, however, computationally effective and produces a reasonable wetting front and runoff prediction (Reference EtcheversEtchevers and others, 2004). The poor results of Table 2 with a RMSE of 5.0% have therefore to be interpreted with caution. A mean deviation (D) of –0.21% indicates that the model slightly underestimated liquid-water content. Figure 5b shows comparisons between the measured and calculated liquid-water content. The graph and the low correlation coefficient (r = 0.21) suggest that there was only a very weak linear correlation. Liquid-water contents were sometimes estimated as 0% in the simulation although liquid water was detected in the measurements. These errors may be caused by the underestimation of snow temperature.

Main source of snow-depth error

Decrease in snow depth is caused mainly by two factors: snowpack settlement and snowmelt. Figure 6 shows the simulated ratio between snow-depth decrease due to snow-melt and that due to snow settlement. Snowmelt was estimated by subtracting simulated settlement from simulated snow-depth change.

Fig. 6. Comparison between snow settlement and snowmelt in 2001/02 at NISIS.

The figure reveals that snowmelt was the dominant factor in the snow-depth decrease throughout the winter. It is thus inferred that snow-depth errors, which appeared at the end of winter, were caused by the underestimation of snow-melt. Therefore, in the following, we examine the surface energy balance for a melting period (20 February–10 March).

The energy balance at the snow surface can be estimated by means of the bulk method:

Here, Q_e is the energy balance at the surface, Q_s is sensible-heat flux, Q_l is latent-heat flux and Q_n is net radiation. Q_n is measured directly, while sensible- and latent-heat fluxes are estimated using the following equations:

Here ρ_a is the density of air, c_pa is the heat capacity, T_s is snow surface temperature and P_a is air pressure. L_w/i are the latent heat for vaporization and sublimation, respectively, is the saturation vapor pressure over water and ice, and rH is relative humidity. Assuming that the snow temperature was 0°C during this period and all of Q_e was consumed in snowmelting, the value of bulk coefficient α can be determined to fit the realistic change of SWE (Fig. 7). Its value was 2.3610^–3 in this study. It is emphasized that this method does not say anything about the true value of the surface heat fluxes but we have adopted a constant value of the exchange coefficient such that the decrease in snow depth and SWE becomes as large as in the observation. This becomes important in the analysis presented below.

Fig. 7. Comparison of changes of snow weight between measurement and bulk method. Solid line indicates measurement, and dashed line indicates bulk method.

Sensible- and latent-heat fluxes in SNOWPACK are estimated based on Monin–Obukhov similarity theory with the assumption of a neutral atmospheric surface layer (Reference Lehning, Bartelt, Brown and FierzLehning and others, 2002b). For our simulations, the value of roughness length adopted was 7.0610^–4 m.

Figure 8a–d show the relations between our “best-fit” bulk method and the simulation of energy balance at the snow surface, sensible-heat flux, latent-heat flux and net radiation at each hour during this period.

Fig. 8. Comparisons between bulk method and simulation in 2001/02 at NISIS. Net radiation was measured directly. (a) Energy balance at surface; (b) sensible-heat flux; (c) latent-heat flux; (d) net radiation.

Table 3 summarizes the results of the descriptive statistical analysis for energy balance at the surface and each heat flux during this period. The comparison shows that on average an additional amount of –11.9Wm^–2 would be required to reproduce the snow-mass decrease as observed. This is consistent with results obtained in the Snow Models Intercomparison Project(SnowMIP; Reference EtcheversEtchevers and others, 2004).

Figure 8a further reveals that the energy balance at the surface in the model corresponded to the measurements when values were 5200 Wm^–2, while the underestimations appeared with larger values. Although negative mean deviations (D) in both net radiation and latent-heat flux indicate that the simulations tended toward underestimation, each high correlation coefficient r (0.98 in net radiation; 0.96 in latent-heat flux) indicates accurate reconstruction. Simulated net radiation calculated from the balance of incoming and outgoing short- and longwave radiation agreed closely with the measurements because in this study the values of daily mean measured albedo were introduced into the model (Fig. 8d). Simulated latent-heat fluxes were very close to the measurements, although there was some underestimation at negative values and overestimation at large values (Fig. 8c).

On the other hand, it is apparent from the large RMSE of 27.4 Wm^–2 and low correlation coefficient (r =0.64) that SNOWPACK simulations of the sensible-heat flux did not correspond to the “best-fit” bulk method. Negative mean deviation (D) indicates underestimation of sensible-heat flux. Figure 8b reveals that most simulated sensible-heat fluxes were smaller than those independently determined by the “best-fit” bulk method, and the differences between them increased as the sensible-heat fluxes increased. This discrepancy is explained in the following way: SNOW-PACK uses an implicit numerical scheme to calculate the temperature distribution (Reference Bartelt and LehningBartelt and Lehning, 2002). This scheme adapts the surface temperature during numerical iterations to the temperature of the new time-step. Because the scheme is implicit, the surface fluxes which depend on the surface temperature, i.e. particularly the sensible-heat flux and to a smaller degree net longwave radiation and latent-heat flux, are adapted too. This scheme is correct for a non-phase-changing snow cover. However, it becomes inaccurate for a phase-changing snow cover because melt processes are treated explicitly after the temperature distribution is calculated. Therefore, when the temperature distribution of a new time-step shows temperatures >0°C, the temperatures are set back to 0°C and the additional energy is used to melt snow. In reality, however, the snow surface temperature stays at 0°C, and higher fluxes of heat to the snow surface are occurring than calculated by SNOWPACK because its numerical surface temperature is >0°C. This scheme is responsible for surface heat fluxes being underestimated for melt situations. Currently, first tests are being performed with a numerical scheme that avoids this inaccuracy. Since the error depends on the numerical calculation time-step, a sensitivity calculation has been performed with a time-step of 5min instead of the default 15 min. The total simulated sensible-heat flux during the period became almost the same as the “best-fit” bulk method.

The values of sensible- and latent-heat fluxes in SNOW-PACK depend also on the roughness length z0; so simulations were carried out using various values of z0 (Table 4). Although maximum snow-depth errors decreased with increases in z₀, >15 cm errors still remain even in the case of z₀ = 2 cm.

Improvement of compressive viscosity

Underestimations of large density seen in Figure 5a imply the possibility that snow densification was smaller than observed.

The compressive viscosity of wet snow is generally less well known than that of dry snow, and the snow cover at NISIS was usually wet throughout the winter. Although the compressive viscosity of wet snow is treated as a function of liquid-water content, which was obtained in the Swiss Alps (cf. equation (31) by Reference Lehning, Bartelt, Brown, Fierz and SatyawaliLehning and others, 2002a), it is possible that SNOWPACK is not suitable for the wet-snow zone in Japan.

Reference KojimaKojima (1967) found that the relation between the compressive viscosity (η) and dry densities (ρ_dry) can be applied to wet snow by subtracting the contribution of free water from wet-snow density having a free-water content of 55%. Endo and others (2002) then described the depth and the density of hourly new snow in the wet-snow areas of Japan using the following equation:

We attempted to introduce Equation (9) in SNOWPACK, assuming that it holds when the liquid-water content is 45%.

Figure 9 shows the density scores obtained through simulation using the conventional formula and Equation (9). The latter became slightly higher than the former. Therefore, the introduction of Equation (9) into SNOW-PACK may be a necessary improvement when the model is applied to warm heavy-snowfall areas. However, the improvement is small compared to the absolute error. Therefore, it appears that problems in the simulation of liquid-water content (Fig. 5) cause most of the error in the wet-snow settlement simulation.

Fig. 9. Comparison of agreement score of density between conventional model and improvement model. Solid line indicates the simulation result by the conventional model. Dashed line indicates the simulation result by the improvement model with Equation (9).