In this study, the authors have employed a combination of empirical modelling and measurements to assess the evolution of the cold content of the snowpack in four closely located boreal forests in Québec, Canada. The general methodology of the study is explained sufficiently well and the results are of interest to the forest hydrology and water resource communities, as well as the land surface and snow modelling communities. I feel that the results are worthy of publication, but I would like to see a little more effort towards examining the differences between the snowpack properties and CC estimates between the sites.

Specific Comments:

In Figure 5, can the authors infer from the detailed snow pit data whether the underestimated SWE is caused by excessive sublimation or melt in the model or underestimated precipitation? Are the latent heat flux measurements of any help here at two sites?

In Figure 6, the density of the top layer may be underestimated if the density of fresh falling snow is underestimated, or if densified snow falling from the canopy is not accounted for. If density is underestimated then CC will be too small even if temperature is correct. Site 1 which has a larger fraction of the canopy buried shows a smaller underestimation of density, suggesting that canopy effects are stronger at Sites 2, 3 and 4.

Would Figures 7 and 10 be more informative if small x-y plots were employed, with the r² or bias value shown? This would also provide the ability to show nonlinear or curvilinear relationships.

Section 4.2: Site A1 is more exposed to turbulence, solar radiation and longwave losses, and so probably has more rapid temperature changes as weather systems change. I believe the authors are capable of adding more interpretation here. Can they infer snowfall, sublimation, decoupling, and canopy effects? I like the discussion about Site A3 on lines 334-341. I see times in Figure 2 when site A3 has the coldest temperatures. Can the authors relate this to periods with stable conditions and/or low wind speeds to bolster this discussion?

Figure 12 and related discussion: Are the authors comparing the formulae for the density of fresh falling snow against the density found in the top 10 cm of the snowpack? I would only consider this a valid comparison if the snow survey were conducted immediately after a snowfall event and before unloading or density changes had taken place.

Line 364-6: If the snow density estimates were disastrous when employing the Brun model, why did the authors not change to a different model or to fixed values of 100 kg/m³ for fresh snow?

Minor points and Corrections:

Line 50: I would phrase this as "plays a central role in the timing of snowmelt". "Delaying" makes it seem like the snow is not melting at the correct time.

Line 59: I would replace “resorted to” with “employed”. “Resorted to” implies that snow pits are not a good method. They are labour intensive, but can provide much information.

Line 62-64: Is “Of note, a slight contrast was observed by Seligman et al. (2014), who reported that the contribution of spring snow storms to CC had a smaller impact on delaying snowmelt than the porous space from dry fresh snow” intended to mean, that in spring storms the snow is near 0°C and so adds the minimum possible heat content based on its mass, whereas the pore space in cold dry low density snow results in a low thermal conductivity which delays snowmelt more than the cold content of the warm spring snow? If so, a few extra words would make that more clear.

Line 64-65: “However, Jennings et al. (2018) reported shifts in the onset of snowmelt by 5.7 h and 6.7 h at alpine and subalpine sites, respectively, when CC at 6:00AM was less than 0 MJ m⁻².” Shifts from when to when, caused by what? By definition, CC must be ≤ 0 MJ m⁻² but surely the amount would affect the timing of snowmelt.

Line 81: I would delete “our” unless the field observations were provided to a different group.

Line 83: I would change “resort to” to “employ” or “use”.

Line 91: I would change “forested” to “forest”.

Line 93: I would change “was” to “were”.

Line 136: I would change “enabling the stability of the prognostic modelled variables” to “ensuring an uninterrupted time series of the prognostic variables”.

Line 152: Change “weighting” to “weighing”.

Line 152-3: I would acknowledge that there can be significant differences in snowfall for some snowfall events over a 4 km distance. It is probably not a consistent spatial bias, however, so the methodology is acceptable as long as there is such an acknowledgement.

Figure 2: Could site A1 be a brighter blue to aid in seeing it as distinct from site A4? The dots are otherwise fine as long as none are completely hidden. If any are completely hidden, perhaps some can have no fill or have the sizes made smaller.

Figure 11: I am not sure what the light blue shading in the caption is referring to.

The authors detail the evolution of snowpack cold content at 4 sites with differing vegetation cover properties in Quebec, Canada. They do this using a combination of field observations and model output over one field season at their research site. They detail differences in cold content across the sites and present extensive model validation plots and figures. There is a sufficient amount of in-depth analysis and it’s clear the authors have devoted a great deal of time and care to the plots, tables, and manuscript. However, there are two major shortcomings/revisions that I feel need improvement before acceptance should be considered:

1. Given the large number of multi-layer models available, why use a single-layer model and then attempt to split the layers? It’s particularly unclear to me how the Andreadis method is applied. The results make it seem like the authors only track 3 layers, but then subsequent sections break out layers in 10 cm increments. Do the layers stay consistent or are they combined and divided as would be done in a multi-layer model? The methods text says simulated SWE from CLASS is used, but then the equations are all for density. Do the final simulated CC series use the thermocouple data or snowpack temperature output from CLASS? If the latter, how is temperature reallocated? I think a schematic showing the complete workflow for the reconstructed CC time series would benefit the readers, in addition to providing much more detail in the methods text.
2. Scientific contribution. I’m left wondering what the major novel contribution of the paper is, which means it needs further work to be accepted in The Cryosphere. Overall, my feeling is the authors can devote less results text to model validation and output (particularly the analysis of the 3-layer scheme as I note in my line-by-line comments below) and spend more time on process-based analysis. For example, what happens during the rain-on-snow events from an energy balance perspective? Why are there differences in cold content development across the sites? Is it the effect of canopy cover on snow accumulation and/or the snowpack energy balance? Given the field observations and modeling, I feel there is a lot more the authors could unpack that give us, the readers, deeper information on the processes governing cold content development and removal at these particular sites. I provide further recommendations in this regard in the line-by-line comments below.

Line-by-line comments:

Line 9: Surface melt can begin before snowpack cold content = 0.

Line 37: The snowpack energy balance predates USACE (1956). For example, Angström’s work from the 1910s on the radiation budget and Western Snow Conference papers from the 1940s published in Transactions of the American Geophysical Union. I’m sure there’s plenty more going back in time.

Line 37: DHSVM was developed with a single-layer snowpack, but now uses a two-layer formulation (Wigmosta et al., 2002): Wigmosta, M. S., Nijssen, B., Storck, P., & Lettenmaier, D. P. (2002). The distributed hydrology soil vegetation model. Mathematical models of small watershed hydrology and applications, 7-42.

Line 46: ρ_w is not in the equation

Line 55: Change snow surveys to snow pits.

Line 55: Reconsider “tedious and demanding.” I'd probably say time-consuming (I find digging snow pits to be quite enjoyable, not tedious).

Line 56: You can note example datasets here. E.g. the Williams and Morse Niwot LTER data: Williams, M., J. Morse, and Niwot Ridge LTER. 2020. Snow cover profile data for Niwot Ridge and Green Lakes Valley, 1993 - ongoing. ver 15. Environmental Data Initiative. https://doi.org/10.6073/pasta/a5fca9d02a4a6a0744cc0d0ffccacd09 (Accessed 2021-05-06).

Lines 56–59: This reads very similarly to the opening sentence of paragraph 2 in Jennings et al. (2018). Please note as such.

Line 59: Please replace all uses of resort/resorted with a more appropriate word (e.g., used, leveraged, utilized, etc.). Resort typically has a negative connotation (i.e., there was nothing else we could use so we had to use this.).

Lines 62–64: This was noted in Jennings et al. (2018) (please add text as such).

Lines 65–66: This is again very similar to the language in the Jennings et al. (2018) manuscript. Please be careful with paraphrasing and providing citations for paraphrased work.

Line 68: Jennings et al. (2018) includes a forested site (the subalpine) and compares it to a higher, open site (the alpine).

Line 69: Remove “obviously”

Lines 78–79: “model articulated around the energy balance” > rephrase for clarity.

Line 84: “even models such as this are not free of biases” > rephrase because all models have some form of bias. You could highlight the previously cited snow model comparison literature.

Line 87: What is “acceptable” in this context? Most readers, myself included, will not have a good intuition of what large and small SSA values are.

Line 89: Remove “obvious” and rephrase sentence. I think you're referring to cold content modeling specifically, but there are many forest snow model studies.

Line 105 (figure 5 caption and throughout): Consider changing the names A1 through A4 to more meaningful terms or abbreviations. There are three forest cover classes (sapling, juvenile and mature). You could use some variation thereof and specify which mature forest site is denser.

Lines 117–124: This could use greater explanation. What instruments were used, what was the vertical frequency of density and temperature sampling, etc? How were the pit values used to impute the missing data?

Lines 128–130: “A simple approach would be to interpolate the density values extracted from snowpits, but this would be incomplete and error-prone given their limited number and absence early in the season.” This seems overly subjective, especially considering the low r² values for the hybrid approach.

Lines 136–137: Further explanation needed. For example, what are the prognostic variables?

Line 147: Change Table 2 to Table 1

Table 2: It doesn’t seem like all sites should be marked as having PPT measurements when the text says it was measured at one location 4km away?

Line 165 (Section 2.2.3): Please see my major revision comment #1. This section needs significant improvement in terms of clarity and specificity.

Lines 169–171: Text says minimum, but equation uses maximum.

Lines 206–207: Are these differences computed within the snowpack at each site per time or computed across all sites?

Line 213 (and throughout text): Please change amplitude to magnitude.

Figure 3: I like the amount of info conveyed by this figure, but find the color scale to be confusing. Because it's representing a single, non-divergent variable, you would be better suited by using the shade of single color. Also, please add "Layer cold content" to the scale bar.

Line 222: Is this the maximum difference or difference in maximum CC? Table 3 shows a larger difference between A1 and A4.

Table 3: The minimum appears to occur in late March (early April?) for site A2 as seen in Fig 4.

Line 232: Please use less subjective terms (for example, "CLASS produces low mean biases in...").

Line 233: Somewhat interestingly, the negative SWE bias would indicate the model has a cold bias in snowpack temperature in order to get CC so accurately. Please include snowpack temperature in the validation.

Line 239: It appears you're producing a 3-layer scheme, but this is not specified in the methods section. Please add that along with information on how the layers are defined/combined/separated (see major revision comment #1).

Lines 240–242: You need to clarify when you're using snow temperature from the thermocouple arrays versus from the snow pit or CLASS. It's not clear here or in the methods section. I think it would be worth adding some material to the methods and again to the results so the readers are certain when observed versus simulated values are being used.

Line 249: This statement seems incorrect. If you're using the snow temperature from the thermocouple, then most of the error is coming from the density estimates and layering scheme.

Figure 6 and preceding text: I’m leaning towards getting rid of the arbitrary 3-layer scheme and only validating/describing the 10 cm layer results shown below. It’s unclear what information/utility the 3-layer scheme provides given that most of the findings are based on the observations and the 10 cm layer discretization.

Figure 7: Please change to a table. The color scale provides the same information as the numbers, but the numbers on their own are easier to interpret.

Lines 264–265: Were the rain-on-snow events missing because no snow pits were dug at those times or the snow pits suggested different changes to cold content than the simulations? Either way, this needs to be clarified.

Figure 8: Same comment as figure 3. Consider changing the color ramp and adding “Layer cold content” to the scale bar.

Lines 273–274: This assertion is presented without supporting evidence. In line with my major revision comment #2, this would be an ideal place to test some process-based hypotheses.

Figure 9: Why is simulated cold content plotted with observed snow depth? It seems like observed snow depth should be plotted with earlier observational figures. Also, keep the rain-on-snow shading consistent with previous figure.

Lines 281–285: These results need more unpacking and their associated methods need to be moved to the methods section. You should also be careful with “positive” and “negative” correlations here. Most cold content values have been discussed in terms of their magnitude, which may lead readers to think cold content declines (i.e. approaches 0) when air temperature decreases. Additionally, simulated and observed values need to be noted explicitly in the text along with the period of comparison (is CC correlated with 30-min air temperature, daily air temperature, or average air temperature to date?).

Lines 286–291: Similar to my above comment, this needs much further explanation. This section could become an important part of the paper (and its novel contribution to the field) if you can further evaluate forest cover differences and their quantitative effect on cold content evolution. For example, you could include an assessment of snowpack energy balance differences and/or changes in snow accumulation as caused by forest cover in both observations and the model.

Figure 10: Please change to table (same comment as figure 7).

Line 298–299: Figure 11 does not show cold content.

Figure 11: Need to clarify if these are simulated or observed values. There is no light-blue shading in the plots. Also, the color ramp should not be divergent as the values they represent are not. Consider using gradation of single color.

Lines 304–309: Figure 4 does not show mass, only cold content. It might be worth further evaluating SWE, depth, and cold content over time in the results sections. Also note that the average winter temperatures in Jennings et al. (2018) were ~4°C cooler at the alpine site (colder frozen mass and more of it led to greater CC in the alpine).

Lines 318–350: I like this section, but I feel like the paper would have a greater impact if there was a greater reliance on results versus discussion when comparing the sites (please see my comment on lines 286–291). For example, you discuss cold air pooling here, but don’t provide data. Why not add this to the results section with data from the air temperature sensors at the different sites? If these data don’t support the hypothesis, then it can be removed. Data from Jennings et al. indicated that the energy balance was typically positive when snow was not actively accumulating. However, there were exceptions at night as a result of radiative cooling from the snowpack. You could provide a comparison from your sites here by providing energy balance output from CLASS in the results.

Figure 12: This is an unfair comparison. The density of freshly fallen snow is not comparable to the density of the top 10 cm of a snowpack. There’s no fresh snow I know of that falls at 400 kg m^-3.

Lines 362–366: Please see comment above. Consider removing plot and text unless the analysis is changed to provide a more important validation of new snow.

Line 397: Please add link in revised manuscript.