User:Bill Flanagan/Notebook/OWW Projects/2009/12/04/Latex-Doc

%Edited by Elaine, 29/09/09 \NeedsTeXFormat{LaTeX2e} \documentclass[annals,prodtf,twocolumn]{igs} \usepackage{igsnatbib} \usepackage[dvipdf]{graphicx} \usepackage{color} \usepackage{subfigure} \usepackage{amsmath, amsthm, amssymb} \usepackage{igsupmath} \usepackage{stfloats} \usepackage{longtable} \usepackage{wasysym} \usepackage{array} \usepackage{psfrag} \usepackage{lineno} %\usepackage[english]{babel} %\usepackage{latexsym} %\usepackage{amsfonts} %\usepackage{tabularx} %\usepackage{epsfig} %\usepackage{floatflt} %\usepackage{rotating} %\usepackage{float} %\usepackage{lineno} %\usepackage{cancel} %\usepackage{psfrag}

\setcounter{secnumdepth}{2} \begin{document} \volume{54} \issueno{} \pubyear{2009} \title[Virutal network for estimating daily new snow water equivalent]{A virtual network for estimating daily new snow water equivalent and snow depth in the Swiss Alps} \author[Egli and others]{Luca EGLI,$^1$\footnote{Present address: WSL Institute for Snow and Avalanche Research SLF, Fl\"{u}elastrasse 11, 7260 Davos-Dorf, Switzerland}~ Tobias JONAS$^1$ and Jean-Marie BETTEMS$^2$}  \affiliation{%  $^1$WSL Institute for Snow and Avalanche Research SLF, Davos, Switzerland\\  E-mail: egli@slf.ch\\    $^2$Federal Office of Meteorology and Climatology MeteoSwiss, Z\"urich, Switzerland} \abstract{Daily new snow water equivalent (HNW) and snow depth (HS) are of significant practical importance in cryospheric sciences such as snow hydrology and avalanche formation. In this study we present a virtual network (VN) for estimating HNW and HS on a regular mesh over Switzerland with a grid size of 7\,km. The method is based on the HNW output data of the numerical weather prediction model (COSMO-7), \textbf{[AUTHOR: IS COSMO an acronym? Please define]} driving an external accumulation/melting routine. The verification of the VN shows that, on average, HNW can be estimated with a mean systematic bias close to 0 and an averaged absolute accuracy of 4.01\,mm. The results are equivalent to the performance observed when comparing different automatic HNW point estimations with manual reference measurements. However, at the local scale, HS derived by the VN may significantly deviate from corresponding point measurements. We argue that the VN presented here may introduce promising cost-effective options as input for spatially distributed snow hydrological and avalanche risk management applications in the Swiss Alps. } \maketitle \section{Introduction} In cryospheric sciences such as snow hydrology, avalanche formation/dynamics and snow climatology, daily new snow water equivalent (HNW, mm) and snow depth (HS, m) are important measurement categories. Together with wind speed, HNW is of significant practical importance for avalanche hazard estimation \citep{mcclung93,egli08} and represents a key input parameter for spatially distributed snow models \citep[e.g.][]{lehning06,liston06}. The monitoring and modelling of water resources, represented by HS (respectively, SWE), \textbf{[AUTHOR: What is SWE? (snow water equivalent?) Is HNW meant? Please clarify!]} is the basic aim of snow hydrology \citep[e.g.][]{grayson00,anderton02}. In particular, estimates of the spatial and temporal distribution of HS (or SWE) \textbf{[AUTHOR: See previous comment]} is essential \citep{luce98,skaugen07} for water resources management \citep[e.g.][]{schaefli07} and flood prevention.

Since the early 1990s, \textbf{[AQ: Have reworded; please confirm correct]} various expensive networks of manual and automatic point HS (and occasionaly HNW) measurements have been set up in the Swiss Alps. Automatic methods have the advantage of providing data at high temporal resolution in terrain that is difficult to access. \cite{durand93} presented a first analysis method for relevant meteorological parameters for snow models, while \cite{eglietal09} investigated different methods to reveal a general feasibility benchmark in automatic estimation of HNW.

In this study, we present a virtual network (VN) to estimate HNW and HS on a regular grid (grid size 7\,km) covering the entire Swiss Alps. The method is based on the HNW output data of the numerical weather prediction model (COSMO-7). \textbf{[AUTHOR: Is COSMO an acronym? Please define]} The raw snow precipitation output of COSMO-7 drives a simple external accumulation/melting routine that allows HNW and HS to be computed for each grid cell. This routine does not account for individual physical processes (e.g. settling of snow cover and fresh snow). Instead, it evaluates bulk accumulation and melting rates based on COSMO-7 HNW data calibrated using observed HS data from existing snow monitoring networks.

We believe that the VN presented here may introduce promising cost-effective options as input for spatially distributed snow hydrological models and avalanche risk management applications. If a numerical weather model and HS point measurement stations are available, a VN network can be developed for other mountainous regions that are not as densely instrumented as the Swiss Alps.

\section{Data}\label{Data} \subsection{COSMO-7 model output data} A numerical weather prediction system is operated by MeteoSwiss (www.meteoswiss.ch) for a wide range of applications. This system is based on the COSMO-7 model (www.cosmo-model.org), which is a primitive equation model with {non-hydrostatic, fully compressible dynamics}. The prognostic variables include precipitation (separated into snow and rain) and many other meteorological parameters. COSMO-7 is used in two modes: (1) a free forecast mode, where the temporal evolution of the atmosphere is computed without any constraints other than the lateral boundary conditions from the driving model (European Centre for Medium-range Weather Forecasts, www.ecmwf.int/research), and (2) an assimilation mode to produce the best gridded representation of the atmosphere by blending available current observations with the model dynamic (the so-called data assimilation process). Note that no observed precipitation is assimilated in COSMO-7. The data assimilation process is not a simple interpolation of available observations, but instead produces a consistent \textbf{[AUTHOR: Is continuous (or constant) meant?]} state of the full atmosphere.

In this study we have chosen to use data from the assimilation mode of COSMO-7 between September 2005 and February 2008. For each day, snow precipitation (i.e. solid precipitation only) cumulated from 0000\,UTC to 2400\,UTC for the automatic reference stations (and from 0800\,UTC to 0800\,UTC on the following day for the manual reference stations) \textbf{[AUTHOR: Is this correct?]} have been extracted from COSMO-7 and used to derive $\mathrm{HNW_{VN.RAW}}$ at the surface. Introduced abbreviations are summarized in Table~\ref{tab1}.

%table 1 here

\begin{table}%table1 \centering \caption{Contingency table for the validation of the HNW estimations of the virtual network. The point measurements at the control sites $\mathrm{HNW_{MEAS}}$ are observed and the estimations of the virtual network $\mathrm{HNW_{VN}}$ are forecast}\label{tab1} \begin{tabular}{@{}lcc}\hline & \multicolumn{2}{c}{$\mathrm{HNW_{MEAS}}$ (observed)}  \\ $\mathrm{HNW_{VN}}$ (forecast) & NO & YES\\ \hline NO  & a & b\\ YES & c & d\\ \hline \end{tabular} \end{table}

\subsection{Point measurements} %figure 1 here

\begin{figure*}%fig1 \centering \includegraphics[width=90mm, angle=-90]{fig_1.ps} \caption{Measurement stations: 141 point measurements of snow depths in the Alpine region over Switzerland. The stations are located at altitudes 800--3130\,m\,a.s.l., where the three different elevation zones are indicated with squares, points and triangles. Reference stations are marked with a cross.}\label{fig1} \end{figure*}

For this study, 141 stations spread over the entire range of the Swiss Alps (Fig.~\ref{fig1}) were available to provide daily measurements of HS during the period September 2005 to February 2008. The sites have been carefully selected in flat open terrain, with as little wind influence as possible, to ensure a regional representativeness of the HS/HNW estimation \citep{egli08}. In order to derive HNW from HS measurements, we used a simple parameterization proposed by \cite{eglietal09}: \begin{eqnarray} && \mathrm{HNW_{MEAS}}(\delta \mathrm{HS_{24h}})\,(\mathrm{mm}) \nonumber \\ &=& \begin{cases} 1+1.09\,\delta \mathrm{HS_{MEAS.24h}}\,(\mathrm{cm}) & \text{for~} \delta \mathrm{HS_{MEAS.24h}} >0 \\ 0 & \text{for~} \delta \mathrm{HS_{MEAS.24h}} \leqslant 0, \end{cases}\label{eqn1} \end{eqnarray} where $\delta \mathrm{HS_{MEAS.24h}}$ denotes the 24\,h difference between an HS reading (measured here in centimetres) and the respective observations from the previous day. Possible limitations of using Equation~(\ref{eqn1}) will be discussed later. \textbf{[AUTHOR: Please provide section number]} Note that for automatic measurements (for ENET and IMIS \textbf{[AUTHOR: Please define acronyms]} network, see \citealp{rhyner02}), HS estimations were provided at 0000\,am while manual HS readings (for observations, see \citealp{marty08}) at around 0800\,am. The values of HS were later checked for plausibility and single missing/faulty values were corrected manually by interpolation. Finally, HNW was only calculated for the period 1 November--30 April for each season (analogous to \citealp{eglietal09}), while HS is considered during the entire seasons of 2005--08 if HS$>$0.

\section{Methodology} \subsection{Virtual network (VN)} \subsubsection{HNW estimated by VN} The cumulated snow precipitation output of COSMO-7 ($\mathrm{HNW_{VN.RAW}}$) constituted the basis of the VN. Thirty-five snow stations (of 141) were set aside for adjusting $\mathrm{HNW_{VN.RAW}}$ to measured HS ($\mathrm{HNW_{MEAS}}$). These stations are referred to as reference stations (Fig.~\ref{fig1}, crosses), while the remaining 106 stations are referred to as control stations. All stations were separated into three different elevation bands (800--1500\,m, 1501--2200\,m and 2201--3200\,m; see legend of Fig.~\ref{fig1}).

Reference stations were selected to cover all regions and elevation bands approximately evenly, in the horizontal as well as in the vertical space. COSMO-7 mesh points were assigned to corresponding snow stations by searching the grid point with centre nearest to the station. The horizontal station coordinates deviated about $\pm$2.7\,km in average (maximum $7\,\mathrm{km}\times \sqrt{2}/2=4.9$\,km) from the centre of $\mathrm{HNW_{VN.RAW}}$, while the vertical difference between stations and respective $\mathrm{HNW_{VN.RAW}}$ grid cells deviated by about $\pm$250\,m \citep{schaub07}. The VN therefore covers an area of approximately 7\,km\,$\times$\,7\,km on the surface (horizontal). In comparison with the point measurement stations, a vertical distention \textbf{[AUTHOR: Is resolution meant? accuracy?]} of about 500\,m is provided.

The correction routine for $\mathrm{HNW_{VN.RAW}}$ first identified possible outliers of $\mathrm{HNW_{VN.RAW}}$ if $\mathrm{HNW_{VN.RAW}}$ exceeded the maximum of $\mathrm{HNW_{MEAS}^{r}}$ at the reference stations (index r). Outliers were replaced by the maximum of $\mathrm{HNW_{MEAS}^{r}}$ for each elevation band. Secondly, a correction matrix $\mathrm{CorrMat^{r}}$ was calculated from the difference between $\mathrm{HNW_{MEAS}}$ and $\mathrm{HNW_{VN.RAW}}$ for each reference station and for each elevation band: \begin{equation} \mathrm{CorrMat^{r}} = \mathrm{HNW_{MEAS}^{r}} - \mathrm{HNW_{VN.RAW}^{r}}. \end{equation} For all remaining locations (index c) of the respective elevation band, the values of the correction matrix were then interpolated by a inverse distance interpolation: \begin{equation} \mathrm{CorrMat^{c}} = \frac{\sum_{i=1}^{N}w_{i}\times \mathrm{CorrMat^{r}}_{i}}{\sum_{i=1}^{N} w_{i}}, \end{equation} where \begin{equation} w_{i} = \frac{1}{d(|r_{i}-c|)} \end{equation} is the inverse distance weighting function. $d$ represents the Euclidean metric distance operator and $N$ is the total number of reference locations ($N=35$ here). The corrected COSMO-7 HNW estimation for the control stations ($\mathrm{HNW_{VN.CORR}^{c}}$) were then defined by: \begin{equation} \mathrm{HNW_{VN.CORR}^{c}} = \mathrm{CorrMat^{c}} + \mathrm{HNW_{VN.RAW}^{c}} \end{equation} where negative values of $\mathrm{HNW_{VN.CORR}^{c}}$ were set to zero.

\subsubsection{HS estimation by the VN} HS is basically calculated from the summation of HNW during the accumulation period, while a simple melting procedure is used to calculated HS during the ablation period. To calculate a temporary term $\mathrm{TEMP.\mathrm{HS_{VN}^{cr}}}(t)$ for the COSMO-7 grid points, the values of $\mathrm{HNW_{VN.CORR}^{cr}}$ were added iteratively day by day: \begin{equation} \mathrm{TEMP.\mathrm{HS_{VN}^{cr}}}(t) = \mathrm{HNW_{VN.CORR}^{cr}}(t) + \mathrm{HS_{VN}^{cr}}(t-1). \end{equation}

Note that instead of a simple cumulation of $\mathrm{HNW_{VN.CORR}^{cr}}(t)$ for each day, this iterative procedure also allows HS to be determined for snowfall events during the melting period. However, at this point $\mathrm{TEMP.\mathrm{HS_{VN}^{cr}}}(t)$ does not indicate a real snow depth measure. In order to convert $\mathrm{TEMP.\mathrm{HS_{VN}^{cr}}}(t)$ to corresponding HS data, we employed data from the reference stations to calibrate $\mathrm{TEMP.\mathrm{HS_{VN}^{cr}}}(t)$ using a simple linear parameterization: \begin{equation} \mathrm{HS_{VN}^{r}}(t) = \mathrm{TEMP}\times \mathrm{HS_{MEAS}^{r}}(t) \times \mathrm{slope} + \mathrm{intercept}.\label{eqn7} \end{equation} Again, such a parameterization was obtained for each elevation band separately. The resulting parameters of Equation~(\ref{eqn7}) (slope and intercept) were subsequently used to convert TEMP.HS to HS at the control stations: \begin{equation} \mathrm{HS_{VN}^{c}}(t) = \mathrm{TEMP}\times \mathrm{HS_{VN}^{c}}(t) \times \mathrm{slope} + \mathrm{intercept}. \end{equation}

The above procedure was applied for all days (considered here as the period of accumulation) except during ablation, which is defined as all days if $\mathrm{mean(HS_{MEAS}^{r}}(t)) - \mathrm{mean(HS_{MEAS}^{r}}(t-4)) < 0$ and $\mathrm{mean(HS_{MEAS}^{r}}(t)) < 0.75 \times \max(\mathrm{mean(HS_{MEAS}^{r}}(t)))$. Note that this criterion does not take into account any information about the physical processes of melting (e.g. liquid water content). It approximately defines the point where the mean snow depth at the reference stations for each elevation band decreases due to melting conditions \citep[c.f.][]{egli09}. In case of melting conditions, we calculated a melting matrix from $\mathrm{HS^{r}_{MEAS}}(t)- \mathrm{HS^{r}_{MEAS}}(t-1)$ for every elevation band at the reference stations: \begin{equation} \mathrm{MeltMat^{r}}(t) = \mathrm{HS^{r}_{MEAS}}(t)- \mathrm{HS^{r}_{MEAS}}(t-1). \end{equation}

We calculated melt rates for all remaining locations using the inverse distance weighting method, in the same way as for Equations~(3) and (4): \begin{equation} \mathrm{MeltMat^{c}}(t) = \frac{\sum_{i=1}^{N}w_{i}\times \mathrm{MeltMat^{r}}_{i}(t)}{\sum_{i=1}^{N} w_{i}}, \end{equation} where \begin{equation} w_{i} = \frac{1}{d(|r_{i}-c|)}. \end{equation}

Finally, HS (for melting conditions) was derived for the control stations as \begin{equation} \mathrm{HS_{VN}^{c}}(t) = \mathrm{MeltMat^{c}}(t) + \mathrm{HS_{VN}^{c}}(t-1). \end{equation}

Note that at the very end of the winter season the reference stations would eventually meltout, preventing any melting rates from being derived. In this case, the final melting rates before meltout have been used for the interpolation. In addition, negative values of $\mathrm{HS_{VN}}$ have been set to 0.

\subsection{Parameters of verification}\label{sect3.2} \subsubsection{HNW verification} In order to validate the performance of $\mathrm{HNW_{VN}}$ at the control stations (in the following the index $c$ is not included) by comparing to the point measurements $\mathrm{HNW_{MEAS}}$, and to compare the results to the performance of other automatic methods, parameters of comparison were used \citep{egli09}. The parameters are briefly summarized below.

%table 2 here

\begin{table*}%table2 \centering \caption{Comparative statistics to assess the performance of VN relative to competing methods tested in \cite{eglietal09}; notation is described in section~\ref{sect3.2} \textbf{[AUTHOR: Will the reader know what all these methods are/how acronyms are defined?]}}\label{tab2} \begin{tabular}{lcccccc} \hline Method & $n_\mathrm{stations}$ & $n_\mathrm{days}$ & $\overline{\delta \mathrm{HNW}}$& $\sigma({\delta \mathrm{HNW}})$ & $R^{2}_{\log}$ & Ranking points \\ \vspace{4.5pt} \\ & & & mm & mm & \\ \hline {\bf VN.CORR} & {\bf 106} & {\bf 484} & {\bf --0.03 ($\pm$ 0.35)} & {\bf 4.01 ($\pm$ 0.62)} & {\bf 0.78 ($\pm$ 0.037)} & {\bf 51} \\ GAUGE & 1 & 725 & --1.15 & 4.12 & 0.89 & 50 \\ SNOWPILLOW & 1 & 719 & 0.6 & 4.41 & 0.79 & 49 \\ SNOWPACK & 1 & 726 & 0.6 & 5.07 & 0.82 & 48 \\ SIMPLE-HNW & 1 & 723 & --0.68 & 5.04 & 0.8 & 45 \\ {\bf VN.RAW} & {\bf 106} & {\bf 484} & {\bf 0.22 ($\pm$ 0.47)} & {\bf 4.89 ($\pm$ 0.65)} & {\bf 0.72 ($\pm$ 0.035)} & {\bf 30} \\ COSMO forecast & 1 & 722 & 0.09 & 6.59 & 0.74 & 30 \\ SNOWPOWER & 1 & 348 & 2.59 & 15.37 & 0.41 & 5 \\ \hline \end{tabular} \end{table*}

\begin{enumerate} \item Systematic bias $\overline{\delta \mathrm{HNW}}$: \begin{equation} \overline{\delta \mathrm{HNW}} = \mathrm{mean(HNW_{VN}-HNW_{MEAS})}. \end{equation}

\item The absolute accuracy $\sigma(\delta \mathrm{HNW})$: \begin{equation} \sigma(\delta \mathrm{HNW}) = \mathrm{standard~deviation(HNW_{VN}-HNW_{MEAS})}. \end{equation}

\item $R^{2}_{\log}$: The coefficient of correlation ($R^{2}_{\log}$) between log-transformed $\mathrm{HNW_{VN}}$ and $\mathrm{HNW_{MEAS}}$.

\item POD and FAR: Parameters adopted from severe weather forecast theory \citep{murphy86} have been used to detect certain classes of snowfall events, namely the Probability of Detection (POD) and the False Alarm Rate (FAR). To calculate POD and FAR, we used contingency tables (Table~\ref{tab2}) for four classes of HNW intensity (Snow-NoSnow: $\mathrm{HNW_{MEAS}}\geqslant 0$\,mm  and $<1$\,mm; Low: $\mathrm{HNW_{MEAS}}\geqslant 1$\,mm and $< 15$\,mm; Medium: $\mathrm{HNW_{MEAS}}\geqslant 15$\,mm and $<30$\,mm; High: $\mathrm{HNW_{MEAS}}\geqslant 30$\,mm).

The POD and FAR were calculated for every class using \begin{equation} \mathrm{POD} = \frac{d}{b+d} \end{equation} and \begin{equation} \mathrm{FAR} = \frac{c}{c+d}, \end{equation} where $b, c$ and $d$ were derived from the contingency table (Table~\ref{tab2}).

\item Ranking points: A ranking point scale for the overall assessment relative to other HNW estimation methods was tested in \cite{eglietal09}. For each comparative measure as described above, the best performance was rated with 7 points and the lowest performance with 0 points. The different ranking points for each parameter and class were added. \end{enumerate}

\subsubsection{HS verification} HS data from the VN ($\mathrm{HS_{VN}}$) was validated against $\mathrm{HS_{MEAS}}$ by applying the following parameters of comparison \citep{egli09}, restricted to the control stations.

\begin{enumerate} \item Systematic bias ($\alpha$): The percent systematic bias $\alpha$ was determined by a least squares fit to \begin{equation} \mathrm{HS_{VN}} = \alpha\times \mathrm{HS_{MEAS}}, \end{equation} where the standard error of $\alpha$ was also determined ($\mathrm{error}_{\alpha}$). Here, $\mathrm{error}_{\alpha}$ is a measure of the spread of the data around the fitting line.

\item Absolute error (RMSE): The averaged absolute difference between $\mathrm{HS_{VN}}$ and $\mathrm{HS_{MEAS}}$, the root mean squared error (RMSE), was included in the comparison statistics: \begin{equation} \mathrm{RMSE} = \sqrt{\frac{1}{N}\sum_{i=1}^{N}\left(\mathrm{HS_{VN}}^i- \mathrm{HS_{MEAS}}^i\right)^2} \end{equation} where $i$ indexes one day of the analysed periods if either $\mathrm{HS_{VN}}^i>0$ or $\mathrm{HS_{MEAS}}^i>0$. Since HS is investigated instead of SWE \citep{eglietal09}, a direct comparison of these studies by a ranking point score is not appropriate. \end{enumerate}

\section{Results and discussion} \subsection{HNW}\label{sect4.1} Results of the validation of $\mathrm{HNW_{VN}}$ and $\mathrm{HNW_{MEAS}}$ are listed in Table~\ref{tab3} and presented in Figure~\ref{fig2}. The parameters of comparison were calculated for both $\mathrm{HNW_{VN.CORR}}$ and $\mathrm{HNW_{VN.RAW}}$ compared to $\mathrm{HNW_{MEAS}}$. For reference, they are listed together with corresponding results of other HNW estimation methods which were tested in \cite{eglietal09}. The parameters for VN.CORR and VN.RAW are averaged over all 106 control stations (Table~\ref{tab3}: $n_\mathrm{stations}$), where the lines below and above the end of the vertical bars in Figure~\ref{fig2} (as well as the numbers in brackets in Table~\ref{tab3}) represent $\pm$half of the standard deviation.

%table 3 here

\begin{table}[h]%table3 \centering \caption{Table of nomenclature}\label{tab3} \begin{tabular}{ll} \hline {Abbreviation} & {Description} \\ \hline HNW & Daily new snow water equivalent \\ HS & Snow depth \\ $\mathrm{HNW_{MEAS}}$ & Observed daily new snow water equivalent \\ VN.RAW & Uncorrected COSMO-7 output \\ VN.CORR & Corrected COSMO-7 (virtual network)\\ $\overline{\delta \mathrm{HNW}}$ & Systematic bias \\ $\sigma(\delta \mathrm{HNW})$ & Absolute accuracy \\ POD & Probability of detection \\ FAR & False alarm rate \\ \hline \end{tabular} \end{table}

%figure 2 here

\begin{figure*}%fig2 \centering \includegraphics[width=90mm, angle=-90]{fig_2.ps} \caption{POD and FAR values are presented for the four different classes of HNW intensity (diamond, point, triangle, star). The values are listed either for single point measurements or the average of 106 control stations, where the error bars denote the deviation from the mean.}\label{fig2} \end{figure*}

Table~\ref{tab3} summarizes the parameters for HNW comparison ($\overline{\delta \mathrm{HNW}}$, $\sigma(\delta \mathrm{HNW})$ and $R^2_{\log}$) \textbf{[AUTHOR: $R^2_{\log}$ not in Table 3]} in order of the highest ranking points. While $\overline{\delta \mathrm{HNW}}<0$ indicates a general underestimation of $\mathrm{HNW_{VN}}$ to $\mathrm{HNW_{MEAS}}$, $\overline{\delta \mathrm{HNW}} >0$ denote a statistical overestimation. The results for VN.CORR and VN.RAW reveal a systematic bias close to 0. However, the considerably large deviation around the mean values of VN.CORR and VN.RAW indicates that some stations exhibit a large under-/overestimation of $\mathrm{HNW_{VN}}$.

A more detailed analysis of the frequency distribution showed that the individual $\overline{\delta \mathrm{HNW}}$s at each control station are approximately equally distributed around the mean value of --0.03\,mm. This indicates an arbitrary character which stems from the comparison of the point measurements $\mathrm{HNW_{MEAS}}$ with data representing a much larger domain of estimation for $\mathrm{HNW_{VN.CORR}}$ or $\mathrm{HNW_{VN.RAW}}$. This is supported by the fact that the stations are randomly distributed in a cuboid of 7\,km\,$\times$\,7\,km by 500\,m \citep{schaub07}, which also explains why the deviations are generally similar for VN.CORR and VN.RAW.

For VN.CORR, $\sigma(\delta \mathrm{HNW})$ demonstrates considerably small values for the mean ($\sigma({\delta \mathrm{HNW}})=4.01$\,mm) with a deviation from the mean of about $\pm$0.62. This low absolute error $\sigma({\delta \mathrm{HNW}})$ implies that, in principle, $\mathrm{HNW_{VN.CORR}}$ can also be calibrated for each single station of the measuring network. The systematic bias of each station can therefore be reduced, resulting in a lower deviation ($\pm$) in $\overline{\delta \mathrm{HNW}}$. However, due to the comparison of a point measurement with a cuboid area, a calibration of $\mathrm{HNW_{VN.CORR}}$ is not meaningful. In fact, the values of $\mathrm{HNW_{VN.CORR}}$ may be considered as an average measurement of HNW over a large area, therefore representing a complementary network to the point measurements.

Taking the performance of SNOWPACK for $\sigma({\delta \mathrm{HNW}})$ as a benchmark, about 82\% of the 106 control stations of the VN exhibit a $\sigma(\delta \mathrm{HNW})< 5.07$\,mm. The absolute accuracy between $\mathrm{HNW_{MEAS}}$ and $\mathrm{HNW_{VN.CORR}}$ is therefore comparable to the operational SNOWPACK-HNW estimations for manual reference measurements investigated at one single location. Note that SNOWPACK calculations are in operational use for Swiss avalanche warning \citep{rhyner02} and provide different types of snowpack properties \citep{lehning02}. However, they are restricted to the locations of the IMIS network \citep{rhyner02}. In contrast, VN.CORR constitutes a complementary network over the entire Swiss Alps (although restricted to the areas of HNW and HS estimation).

The POD/FAR statistics show that VN.CORR performs for all classes of intensities with better POD/FAR values than VN.RAW. Note that POD denotes the percentage of $\mathrm{HNW_{MEAS}}$, also measured by $\mathrm{HNW_{VN}}$, while FAR represents the percentage of $\mathrm{HNW_{VN}}$ which are not observed by $\mathrm{HNW_{MEAS}}$. VN.CORR has similar features to the POD/FAR statistics of SNOWPACK, in particular for the two highest classes of intensities (Medium and High) and for the Snow-NoSnow class, most important for avalanche risk management decisions. \textbf{[AUTHOR: Have reworded; please confirm correct]} The good performance of VN.CORR in POD/FAR statistics with respect to the Snow-NoSnow class also implies that COSMO-7 is capable of distinguishing between snow and rain. Note that results do not deteriorate if evaluated for the three elevation bands separately (data not shown).

Furthermore, the calibration routine to derive VN.CORR implies that at the reference station $\mathrm{HNW_{VN.CORR}} = \mathrm{HNW_{MEAS}}$, independent of the processes leading to $\mathrm{HNW_{MEAS}}$ such as combinations of solid, liquid precipitation and evaporation.

The low performance of VN.CORR in POD ($\sim$0.45) and FAR ($\sim$0.35) for the Low class diverges from the SNOWPACK results. In this class, VN.CORR is comparable to the performance of the SIMPLE-HNW method. This is unsurprising given the fact that HNW are determined using SIMPLE-HNW at the control stations. This limits the potential of our approach to deal with situations in which the change in snow depth is strongly influenced by both snowfall and settling. While Equation~(\ref{eqn1}) implicitly includes the effect of settling of new snow and the entire snow cover, it may become inaccurate after the first day of a multi-day snowfall event. In such situations VN.CORR therefore may not perform as well as a physical snowpack model. Note also that the difference in the time of measurements between the manual and automatic measurements (8\,h; see section~\ref{Data}) may lead to erroneous HNW estimates in a few cases (e.g. if a significant precipitation event occurs between 0000 and 0800). However, an HNW analysis of manual and automatic stations separately (data not shown) revealed very similar results with respect to POD/FAR statistics, which implies that this effect is minor.

%figure 3 here

\begin{figure*}%fig3 \centering \includegraphics[width=140mm]{fig_3.eps} \psfrag{Date} \caption{The temporal development of snow depths derived by the virtual network (black lines) and the point measurements (grey lines) for three selected stations. The systematic bias ($\alpha$) and the abolute error (RMSE) are indicated in the legend. \textbf{[AUTHOR: Can a new version of better quality be provided?]}}\label{fig3} \end{figure*}

Overall, VN.CORR demonstrates an advantage over VN.RAW. In particular, the calculation described by Equations~(1)--(5) demonstrated the largest improvement. The method of removing the outliers has a minor impact on the results of VN.CORR. Moreover, in comparison to the other methods, VN.CORR yields the highest number of ranking points and therefore constitutes an adequate automatic estimation for HNW over the entire Swiss Alps.

\subsection{HS}

Figure~\ref{fig3} shows the temporal development of the modelled HS ($\mathrm{HS_{VN}}$, black lines) and the measured HS ($\mathrm{HS_{MEAS}}$, grey lines) during the season of winter 2006/07 for three locations of automatic measurement stations (ELM~2, DAV~2 and DAV~3). A temporal development of HS is first displayed where $\mathrm{HS_{VN}}$ describes a curve consistently below $\mathrm{HS_{MEAS}}$ (ELM~2), resulting in a systematic underestimation of $\alpha= 0.6$ and a mean absolute error (RMSE) of 0.42\,m. Secondly, a trajectory is shown (DAV~3) where the systematic bias ($\alpha$) is close to 1 and RMSE\,=\,0.19\,m, indicating a strong congruence between the curves. Finally, a station (DAV~2) is shown which describes a path constantly above $\mathrm{HS_{MEAS}}$, resulting in an $\alpha=1.54$ and RMSE\,=\,0.45\,m (i.e. similar to that of ELM~2). Despite the fact that $\mathrm{HS_{VN}}$ shows considerable deviation from $\mathrm{HS_{MEAS}}$ for some stations, all stations generally reproduce the temporal progression qualitatively well. This is particularly true for the peaks during accumulation and the following settling (produced by Equations~(7) and (8)). The same result is obtained for the ablation derived by Equations~(9)--(12).

Basically, the ablation period is generated by the VN only by means of differences of HS, where the melt rate matrix (Equation~(8)) actually refers to a combination of some fresh snow accumulation, settling and melting. As the peaks are attributed to snow accumulation by daily new snowfall ($\mathrm{HNW_{VN.CORR}}$), the deviation from $\mathrm{HS_{MEAS}}$ stem from possible under-/overestimations of $\mathrm{HNW_{VN.CORR}}$, as discussed in section~\ref{sect4.1}. Analogous to the results and discussion of the systematic bias of $\mathrm{HNW_{VN.CORR}}$, the specific values of $\alpha$ for each station are approximately equally distributed around the mean value of $\overline{\alpha}=1.06 \pm 0.286$. Accordingly, while $\alpha$ is reasonably close to 1, $\overline{\mathrm{error}_{\alpha}}=0.022 \pm 0.026$ and $\overline{\mathrm{RMSE}}=0.40 \pm 0.12$\,m are considerably large (see Table~\ref{tab4}). %table 4 here

\begin{table*}%table 4 \centering \caption{The parameters of $\mathrm{HS_{VN}}$ evaluation statistics ($\alpha$, $\mathrm{error}_{\alpha}$ and RMSE) as described in section~\ref{sect3.2} calculated either for the analysis per station averaged over 106 stations ($\mathrm{HS_{VN}}$ per station), or for the analysis of $\mathrm{HS_{VN}}$ averaged over all stations (average of $\mathrm{HS_{VN}}$)}\label{tab4} \begin{tabular}{lccccc} \hline & $n_\mathrm{stations}$ & $n_\mathrm{days}$ & $\overline{\alpha}$ & $\overline{error_{\alpha}}$ & $\overline{RMSE}$ \\ $\mathrm{HS_{VN}}$ per station & 106 & 782 ($\pm$ 73)  & 1.06 ($\pm$ 0.286) & 0.022 ($\pm$ 0.026) & 0.40 ($\pm$ 0.12)   \\ \hline & $n_\mathrm{stations}$ & $n_\mathrm{days}$ & $\alpha~ of~ \overline{HS}$ & $error_{\alpha}~ of~ \overline{HS}$  & $RMSE~ of~ \overline{HS}$ \\ Average of $\mathrm{HS_{VN}}$ & 106 & 809 & 0.96  & 0.004  & 0.08 \\ \hline \end{tabular} \end{table*}

\subsection{Spatial analysis} Since the analysis of the single point measurements of $\mathrm{HS_{VN}}$ is limited due to its random nature, an investigation of $\mathrm{HS_{VN}}$ and $\mathrm{HNW_{VN.CORR}}$ considering the Swiss Alps as an entire area is also considered. First, the mean of all snow depths of point and VN estimations ($\overline{\mathrm{HS_{MEAS}}}$ and $\overline{\mathrm}$) and their standard deviation ($\sigma(\mathrm{HS_{MEAS}})$ and $\sigma(\mathrm{HS_{VN}})$) are calculated. Figure~\ref{fig4} displays the evolution of $\sigma(\mathrm{HS_{VN}})$ with $\overline{\mathrm}$ (black symbols) and $\sigma(\mathrm)$ with $\overline{\mathrm{HS_{MEAS}}}$ (grey symbols), where each point of the trajectory represents one day of the year 2006/07. The curves show a characteristic hysteretic dynamics as discussed in \cite{egli09}. The trajectories of the period of accumulation (points) and ablation (triangles) are clearly separated. The quasi-linear increase of $\sigma(\mathrm{HS_{VN}})$ with increasing $\overline{\mathrm}$ highlights that accumulation of snow leads to an increase in the differences between sites. The path during ablation, on the other hand, is mainly attributed to the spatial distribution of melting rates in the Alpine region.

%figure 4 here

\begin{figure*}%fig4 \centering \includegraphics[width=100mm, angle=-90]{fig_4.ps} \caption{The evolution of the mean snow depths ($\overline{\mathrm{HS}}$) and the standard deviation of snow depths ($\sigma{(\mathrm{HS})}$) derived by the virtual network ($\mathrm{HS_{VN}}$, black symbols) and the point measurements ($\mathrm{HS_{MEAS}}$, grey symbols). The trajectory during the period of accumulation is indicated by points; the trajectory during ablation is indicated by triangles.}\label{fig4} \end{figure*}

It has been shown \citep{egli09} that the curves of the hysteretic dynamics are similar between years and therefore characterize the seasonal development of HS in the Swiss Alps. Figure~\ref{fig4} shows that both the periods of accumulation and ablation are well reproduced by $\mathrm{HS_{VN}}$ when compared to $\mathrm{HS_{MEAS}}$. The seasonal characteristic of the HS development is therefore also displayed by the VN. We therefore speculate that $\mathrm{HS_{VN}}$ correctly reproduces the total amount of snow depths during a season at Swiss Alpine scale. Moreover, when evaluating $\alpha$, $\mathrm{error}_{\alpha}$ and RMSE for $\overline{\mathrm{HS_{MEAS}}}$ and $\overline{\mathrm}$ (HS is averaged day by day over all stations), the percentage bias is as small as --4\% (Table~\ref{tab4}).

Furthermore, the variogram analysis of HNW \citep{egli08} showed that the statistical correlation length of HNW over the entire Swiss Alps is about 50--60\,km. The same correlation length was found in \citep{schaub07}, analysing the raw HNW output of COSMO ($\mathrm{HNW_{VN.RAW}}$). This additionally supports the assumption that $\mathrm{HNW_{VN}}$ reproduces the regional precipitation patterns of the Swiss Alps appropriately and can be used as an input parameter for spatially distributed snow models.

Finally, the POD statistics of point HNW measurements as a function of the distance between two measuring sites \citep{egli08} showed that, for the smallest scale (5--10\,km), the POD is about 60\% for HNW intensities $>30$\,mm (High). The same value has been derived from a comparison of $\mathrm{HNW_{VN.CORR}}$ and $\mathrm{HNW_{MEAS}}$. Again, a better performance of the VN HNW estimation in the POD statistics of two point measurement stations about 7 km distant would have been astonishing, since the large area of $\mathrm{HNW_{VN.CORR}}$ covers a grid cell of 7\,km\,$\times$\,7\,km. As the POD decreases rapidly with the distance from a point measurement station \citep{egli08}, the VN may also be applied for avalanche risk management warning for the locations between the point measurements where no HNW estimation is available.

\section{Conclusion and outlook} In this study, the HNW output of the numerical weather prediction model COSMO-7 was coupled with a simple snow accumulation/melting model in order to provide HNW/HS grids of 7\,km resolution for the entire Swiss Alps. The results for the $\mathrm{HNW_{VN}}$ estimation showed that, on average, its performance is comparable to the performance of different automatic point measurements. $\mathrm{HNW_{VN}}$ therefore represents a complementary network for avalanche risk management applications and can also be used as input data for spatially distributed snow hydrological models where HNW is required.

This is also the case for the HS estimation by the VN, where the statistical dynamics of the mean of HS ($\overline{\mathrm}$) and its standard deviation ($\sigma(HS_{VN})$) are congruent to the measured dynamics. As a consequence, the estimations of HS by the VN over the entire Swiss Alps may be used to estimate the total amount of snow and snow water equivalent for larger catchments. At a local scale, however, a direct comparison of the point measurement to the corresponding grid cell of the VN may result in considerable deviation. This may stem from the fact that a measurement over a large area (7\,km\,$\times$\,7\,km for the VN) is compared to a single point measurement.

Future effort is necessary to investigate the spatial scales over which the VN is capable of providing a good performance for estimation of HNW/HS. For this purpose we will extend the application to COSMO-2, a new high-resolution version of COSMO available from February 2008 for a mesh size of 2.2\,km. Additionally, the differences in the VN between the assimilation mode and the forecast mode of COSMO-2 will be investigated. Finally, the principle of the VN presented here can be applied to other regions which are less densely equipped with automatic or manual point measurement stations than the Swiss Alps. If point measurement stations and a numerical weather prediction model are available, HNW and HS can be estimated with a refined accuracy for spatially larger extended regions for avalanche risk and water resources management.

\section{Acknowledgements} The authors thank D.~Schaub and F.~Schubiger for their support with the data acquisition and processing and the external reviewer E.~Brun for his constructive suggestions which helped to improve the manuscript. % \begin{verbatim} \begin{thebibliography}{9} \expandafter\ifx\csname natexlab\endcsname\relax \def\natexlab#1{#1}\fi \expandafter\ifx\csname selectlanguage\endcsname\relax \def\selectlanguage#1{\relax}\fi

\bibitem[Anderton and others, 2002]{anderton02} Anderton, S.P., S.M. White and B. Alvera. 2002. Micro-scale spatial variability and the timing of snow melt runoff in a high mountain catchment. \textit{J. Hydrol.}, \textbf{268}(1--4), 158--176.

\bibitem[Durand and others, 1993]{durand93} Durand, Y., E. Brun, L. M\'{e}rindol, G. Guyomarc'h, B. Lesaffre and E. Martin. 1993. A meteorological estimation of relevant parameters for snow models. \textit{Ann. Glaciol.}, \textbf{18}, 65--71.

\bibitem[Egli, 2008]{egli08} Egli, L. 2008. Spatial variability of new snow amounts derived from a dense network of Alpine automatic stations. \textit{Ann. Glaciol.}, \textbf{49}, 51--55.

\bibitem[Egli and Jonas, 2009]{egli09} Egli, L. and T. Jonas. 2009. Hysteretic dynamics of seasonal snow depth distribution in the Swiss Alps. \textit{Geophys. Res. Lett.}, \textbf{36}(2), L02501. (10.1029/2008GL035545.)

\bibitem[Egli and others, 2009]{eglietal09} Egli, L., T. Jonas and R. Meister. 2009. Comparison of different automatic methods for estimating snow water equivalent. \textit{Cold Reg. Sci. Technol.}, \textbf{57}(2--3), 107--115.

\bibitem[Grayson and Bl\"{o}schl, 2001]{grayson00} Grayson, R. and G. Bl\"{o}schl, \textit{eds}. 2001. \textit{Spatial patterns in catchment hydrology: observations and modelling.} Cambridge, etc., Cambridge University Press.

\bibitem[Lehning and others, 2002]{lehning02} Lehning, M., P. Bartelt, B. Brown, C. Fierz and P. Satyawali. 2002. A physical SNOWPACK model for the Swiss avalanche warning. Part II: snow microstructure. \textit{Cold Reg. Sci. Technol.}, \textbf{35}(3), 147--167.

\bibitem[Lehning and others, 2006]{lehning06} Lehning, M., I. V\"{o}lksch, D. Gustafsson, T.A. Nguyen, M. St\"{a}hli and M. Zappa. 2006. ALPINE3D: a detailed model of mountain surface processes and its application to snow hydrology. \textit{Hydrol. Process.}, \textbf{20}(10), 2111--2128.

\bibitem[Liston and Elder, 2006]{liston06} Liston, G.E. and K. Elder. 2006. A distributed snow-evolution modeling system (SnowModel). \textit{J. Hydromet.}, \textbf{7}(6), 1259--1276.

\bibitem[Luce and others, 1998]{luce98} Luce, C.H., D.G. Tarboton and K.R. Cooley. 1998. The influence of the spatial distribution of snow on basin-averaged snowmelt. \textit{Hydrol. Process.}, \textbf{12}(10--11), 1671--1683.

\bibitem[Marty, 2008]{marty08} Marty, C. 2008. Regime shift of snow days in Switzerland. \textit{Geophys. Res. Lett.}, \textbf{35}(12), 12501. (10.1029/2008GL033998.)

\bibitem[McClung and Schaerer, 1993]{mcclung93} McClung, D.M. and P.A. Schaerer. 1993. \textit{The avalanche handbook.} Seattle, WA, The Mountaineers.

\bibitem[Murphy and Winkler, 1987]{murphy86} Murphy, A.H. and R.L. Winkler. 1987. A general framework for forecast verification. \textit{Mon. Weather Rev.}, \textbf{115}(7), 1330--1338.

\bibitem[Rhyner and others, 2002]{rhyner02} Rhyner, J. \textit{and 7 others}. 2002. Avalanche warning Switzerland -- consequences of the avalanche winter 1999. \textit{In} Stevens, J.R., \textit{ed.} \textit{Proceedings of the International Snow Science Workshop 2002, 30 September--3 October 2002, Penticton, British Columbia}. Victoria, B.C., B.C. Ministry of Transportation. Snow Avalanche Programs. \textbf{[AUTHOR: please provide page nos.]}.

\bibitem[Schaefli and others, 2007]{schaefli07} Schaefli, B., B. Hingray and A. Musy. 2007. Climate change and hydropower production in the Swiss Alps: quantification of potential impacts and related modelling uncertainties. \textit{Hydrol. Earth Syst. Sci.}, \textbf{11}(3), 1191--1205.

\bibitem[Schaub, 2007]{schaub07} Schaub, D. 2007. Vergleich realer IMIS Wetterdaten gegen COSMO-7 Modelldaten. (Diplomarbeit thesis, Z\"{u}rcher Hochschule f\"{u}r Angewandte Wissenschaft.) \textbf{[AUTHOR: is this a completed thesis, with an awarded degree? Please check institution]}

\bibitem[Skaugen, 2007]{skaugen07} Skaugen, T. 2007. Modelling the spatial variability of snow water equivalent at the catchment scale. \textit{Hydrol. Earth Syst. Sci.}, \textbf{11}(5), 1543--1550.

\end{thebibliography}

\end{document}