1. Introduction
Water has been studied for a long time for a number of reasons including its importance in terms of water resources that directly relate to human life activities. Therefore, in hydrology research, studies have been continuously conducted to understand the water circulation mechanism and to efficiently manage limited water resources. Along with development of technology, there has been a breakthrough in hydrology. However, the amounts of snow-melt and evapotranspiration remain critical missing points in the water cycle. Because the amount of snow-melt can be utilized not only as a study of the circulation mechanism of water in hydrology but also as an indirect indicator for monitoring global warming, studies regarding snow-melt have been conducted in the environmental field during recent years. To estimate the amount of snow-melt, the amount of snowfall must be known, which can be calculated indirectly from snow cover area and snow depth. Because it is not possible to measure at all points in a large area, a method of estimation using a limited number of Automated Weather Source (AWS) data is mainly used. However, there are many limitations because observational stations are sparse and irregularly distributed. To overcome these limitations, a remote-sensing method has been proposed as a solution. However, it is difficult to estimate the depth of snow compared to the easily estimated snow cover area.
Several researchers have proposed a linear sow depth estimation algorithm that is a spectral polarization difference (SPD) algorithm based on the relation between the microwave brightness temperature difference and snow depth (Chang et al., 1987; Foster et al., 1997; Aschbacher, 1989). The microwave brightness temperature difference is calculated from the difference between 18 GHz and 36 GHz. However, such a linear algorithm has a limitation in that the accuracy greatly deteriorates beyond a certain depth. It is also noted that the accuracy of the estimated snow depth varies greatly according to the characteristics of the region. To overcome these limitations, Davis et al. (1993) conducted a study to estimate the parameters affecting the correlation between microwave brightness temperature difference and snow depth, and identified that grain size, density, etc. should be considered in snow depth estimation. Goita et al. (2003) noted that the two variables (microwave brightness temperature difference and snow depth) have different relationships depending on the cover type and type of vegetation in the area where the snow exists. Tedesco et al. (2004) estimated the snow depth of by applying the artificial neural network (ANN) method to SSM/I (Special Sensor Microwave/Imager) data, and validated the result compared to an existing algorithm using R-squared and root mean square error (RMSE). The R-squared and RMSE were 0.73 and 14.28 cm of the depth, respectively. Zahir and Mahdi (2015) estimated snow depth by applying genetic algorithm and Support Vector Regression (SVR) to SSM/I radar images. As a result, an RMSE of 6.97 cm and an R-squared of 0.9 were derived. Liang et al. (2015) additionally used the visible band and near infrared (NIR) band of a MODIS image with SSM/I data. The result was an RMSE of 6.21 cm and an R-squared of 0.87. As previously mentioned, research using microwaves has been extensively conducted; however, has been noted that that the coverage area is narrow and the time resolution is low.
Optical satellites such as MODIS, Geostationary Operational Environmental Satellite - R Series, and Himawari-8 can cover a large area and have a short-period time resolution, which can solve the problem noted as a limitation in microwave images. However, unlike microwaves, the wavelength of visible range irradiated from the sun can only penetrate a few centimeters of the snow. Thus, the depth must be indirectly estimated. It is known that snow up to a certain thickness gradually increases the reflectivity of the surface Baker et al. (1991), because the snow on the surface reflects sunlight. Romanov and Tarpley (2004) estimated the snow depth within 30 cm in an unforested plain area by calculating the relationship between snow depth and Snow Cover Fraction (SCF) which represents the portion of snow-covered area in the pixel. Romanov and Tarpley (2007) estimated the depth of snow in a forested area as well as a plain area considering the forest type which is a factor affecting the SCF. The validation result of the estimated snow depth using the derived empirical equation of the study was approximately 10 cm of RMSE.
As a result, the most important process for estimating snow depth using optical images is to obtain the SCF for each pixel unit. Appel and Salomonson (2002) suggested estimating the SCF using a quadratic equation with NDSI or NDVI as independent variables. Romanov et al. (2003) proposed a method of estimating the SCF through a simple equation with the reflectance of land and snow. Salomonson and Appel (2004) proposed an NDSI-based first-order equation and the results of an empirically derived equation were validated using Landsat 7 Enhanced Thematic Mapper Plus (ETM+) images. Lin et al. (2012) proposed to calculate the SCF from each variable such as NDSI, Ratio Snow Index (RSI), and Difference Snow Index (DSI) using an estimation equation based on an exponential function. The result found that the RSI estimation shows the highest R-square value (0.83) through validation using Landsat 7 ETM+.
It is quite significant that the proposed equation obtained information that is difficult to estimate by indirectly calculating the SCF using an index correlated with the SCF. However, existing equations have a limitation in that if the value of the independent variable is less than a certain value, the SCF tends to be underestimated, and if it is more than a certain value, the SCF tends to be overestimated. In addition, although there are any factors affecting the SCF, existing equations only consider one factor to estimate the SCF. Therefore, in this paper, we 1) suggest an equation for estimating a more precise SCF by solving these problems, 2) develop a more accurate estimation equation for snow depth which is optimized for the Republic of Korea using the relationship between the calculated SCF and measured snow depth, and 3) suggest a method which conducts post-processing using the calculated snow depth based on the relative error of the image for improving the accuracy of the calculated snow depth map.
2. Study Area and Data
The Republic of Korea is a peninsular country between latitude 33-39 north degrees and longitude 124-130 east degrees. Because it is in the mid-latitude region, the four seasons of spring, summer, autumn, and winter are clearly distinct. In the Republic of Korea, there are different meteorological characteristics in each region by season because the three sides of the country are surrounded by the sea and there are many mountainous areas. Winter typically starts at the end of November and continues until the end of February of the following year. During this period, snowfall can be observed in all parts of the country and it is as little as 1 cm and as much as approximately 1 m. Because the country stretches from north to south, the temperature difference between the north and south is quite significant and snowfall is more frequently observed in the north. I fact, Gangwon-do in the high-latitude region, a large amount of snow is annually observed. Whereas, the observational days of snow are not only relatively rare but also limited in amount in Jeolla-do and Gyeongsang-do, which are at further south compared to Gangwon-do. Therefore, the pattern of snowfall is various. In addition, 60% of the Republic of Korea’s land comprises forests and there is considerable rough mountainous terrain. Therefore, it is necessary to derive the SCF estimation and snow depth estimation equations that can be applied not only to the Republic of Korea but also other regions that have complex regional characteristics. In this study, snow depth was estimated from MODIS images taken of the Republic of Korea. Fig. 1 shows the geographical location of the study area. The red marks in Fig. 1 indicate the locations of manned or unmanned weather stations that record the in-situ snow depth.
Fig. 1. Location of the study area; weather stations are marked with red dots.
The purpose of this study was to develop an algorithm to estimate snow depth in the Republic of Korea using MODIS images. MODIS is a sensor that can acquire 36 bands from visible to infrared, with resolutions of 250 m, 500 m, or 1 km, and is onboard the Aqua and Terra. Its re-visit period is very short, approximately one day, and it has a swath width of approximately 2,330 km, which is suitable for observing snow cover in a wide area. In this study, a snow depth algorithm was designed based on the relationship between the snow depth and SCF. For this purpose, the red and NIR bands were used to calculate the NDVI, and green, and short wave infrared (SWIR) bands were used to calculate the NDSI. In addition, Landsat images with a resolution of 30 m were used together to estimate the pixel unit truth SCF of the MODIS images. The Landsat images included band information from the visible to the infrared region and had a 16-day re-visit period. Therefore, we constructed the data for deriving the estimation equation by acquiring the Landsat-8 image where the snow cover was observed in the study area and the MODIS image of the same date as the Landsat-8 image was also acquired. The constructed data included 38 Landsat images and 27 MODIS images for a total of 38 pairs. It is noted that the Landsat images can be replaced into the Sentinel-2 images that have a 5-day revisit period. To derive the coefficient of the equation of snow depth from the SCF and to validate the accuracy of the estimated snow depth map using the derived equation, we used in-situ snow depth data provided from meteorological stations of the Republic of Korea. Depth is recorded every hourly in 0.1-cm units. To minimize the error in the time difference between the measured value and the estimated value of the snow depth, the snowfall data recorded at the temporally nearest the acquisition time of the MODIS image were used.
3. Method
In this study, there were two main steps(Fig. 2). The first step was the part to derive the coefficients of the equation for estimating the SCF and snow depth, and the second step was correction of the snow depth map calculated from the MODIS image using the acquired equation in step 1. For the first step, the coefficients of the equation for estimating the SCF map of the MODIS image from the NDSI and NDVI values were derived using the Landsat top-of-atmosphere (TOA) images and the MODIS images as input data. In this study, we derived the equation based on the Gaussian function. Using the derived equation, an SCF map of MODIS was developed from the NDSI and NDVI maps. Then, an estimation equation of snow depth was derived based on the relationship between the calculated SCF map and the in-situ snow depth data. The second step was to correct the calculated snow depth value to improve the accuracy of the snow depth map. Assuming that the conditions in a single image are the same, the relative difference of each pixel in the calculated snow depth map had a high correlation with the actual difference. Using this, it was possible to derive a linear equation that can post-process the calculated snow depth map. The final snow depth map was developed using the obtained equation, and a validation was performed using this map compared to in-situ values.
Fig. 2. Detailed work flow in developing the snow depth map from MODIS images.
In this study, reflectance values of Landsat and MODIS images were calculated and used through preprocessing. In addition, clouds were masked out prior to the processing because they is only cause error during the calculating step of the true SCF by distorting the results classifying the Landsat image into snow and non-snow pixels, but they also affect the calculation of the NDSI value in the MODIS image. For convenience of calculation, the MODIS images used in this study were resampled to have a resolution of 990 m. This was done to ensure that the same area covered by the two images would fit in each resolution when calculating the true SCF value of a MODIS image from a Landsat image with a resolution of 30 m. The ratio of the Landsat pixel (spatial resolution: 30 m) to the MODIS pixel (spatial resolution: 1000 m, resampled resolution: 990 m) covering the same area was expressed as an integer ratio.
1) Generation of truth SCF map
The SCF is the proportion of snow-covered area in a pixel. Assuming that there is one pixel in a MODIS image with 1 km of spatial resolution, it covers an area of 1 km2. The SCF is the ratio of the area occupied by snow-covered areas of 1 km2. Assuming that Fig. 3 covers a pixel of MODIS with a resolution of 990 m, it can cover 33 Landsat pixels in the row and column directions(1089 total pixels in total), respectively. The SCF of the MODIS pixel can be calculated using the number of pixels classified as snow among the 1089 Landsat pixels. The green-colored pixels represent a snow-free area and the blue-colored pixels represent a snow-covered area in Fig. 3. The number of snow-covered pixels is 312. Therefore, the SCF of Fig. 3 is approximately 0.29.
Fig. 3. Landsat image clip covering one pixel of a MODIS image. Snow-free and snow-covered pixels are represented by the green and blue colors, respectively.
To generate the true SCF map of MODIS, a snow cover map should first be developed. To do so, an NDSI map was created using equation 1 in the cloud-masked Landsat image. The NDSI was calculated using the green band and the SWIR band. In the case of the green band, the reflectance value was high in the snowy region but the SWIR had a low reflectance value. Therefore, a snow-covered pixel had a high value of NDSI and the opposite was the case for a snow-free pixel as follows:
\(N D S I=\frac{\text { Reflectance Green }-\text { Reflectance } S W I R}{\text { Reflectance_Green }+\text { Reflectance } S W I R},\) (1)
where the green and SWIR reflectance values of the Landsat-8 image correspond to bands 3 and 6, respectively. Based on the NDSI map and spectral characteristics, we developed a snow cover map, which was divided into snow-covered pixels and snow-free pixels. Then, a true SCF map was produced using a 33 by 33 matrix to have the same resolution as the MODIS image of 990-m resolution.
2) Deriving the SCF estimation equation
The SCF map of the image can be calculated using the difference in spatial resolution. However, it is impossible to obtain the high-resolution image of the same time and regions of the MODIS image because of its vast swath. To solve this problem, studies have been conducted to indirectly calculate the SCF using indices acquired from the image. The NDSI is known to have the most dominant relationship with SCF. Therefore, the equation for estimating the SCF using the NDSI has been suggested in many previous studies. However, although the NDSI is the most critical variable for obtaining the SCF, other variables also affect the SCF. Therefore, in this study, we proposed an estimation equation that used the NDVI together with the NDSI considering that there are many forests in the Republic of Korea. The NDVI is a vegetation sensitive index that can be calculated using the red band and the NIR band. The following is the equation for calculating the NDVI:
\(N D V I=\frac{\text { Reflectance }_{-} \text {NIR }-\text { Reflectance } \text { Red }}{\text { Reflectance_NIR }+\text { Reflectance_Red }} \text { , }\) (2)
where the red and NIR reflectance values of MODIS correspond to bands 4 and 5, respectively.
In addition, when deriving the estimation equation, previous studies have suggested an equation based on linear, quadratic, and exponential functions. However, in this study we derived the equation based on the Gaussian distribution. The Gaussian distribution enable us to estimate the snow depth more precisely.
In the case of the MODIS image used during this step, only the same area as the Landsat image was clipped. To derive the coefficients of the equation, the equation should be derived empirically based on the relationship between the variables (NDSI and NDVI) and the true SCF. However, because the vast swath of the MODIS image is too wide, distortion becomes severe toward the edge of the image. Because of this distortion, the Landsat image and the MODIS subset image do not exactly coincide, which is an error factor in deriving the correct estimation equation. Therefore, in this study, co-registration of the MODIS image to the Landsat image was performed using a cross-correlation method. The normal coefficient between the two images was calculated by applying moving window to the NDSI map of the Landsat image and the MODIS image. A polynomial equation was derived using the correlation coefficient calculated for each pixel and reconstructing the image using it. The same equation was applied to the NDVI image to produce an NDVI subset image that exactly spatially matches the Landsat image. Then, the coefficients of the equation were derived to estimate the SCF using the relationship between the two indices and the true SCF. In this case, the images without the pixel of the water body area were used because the NDSI value was high and the NDVI value was low in areas similar to the snow area, which may affect the result.
3) Deriving the snow depth estimation equation and the correction of value
The obtained SCF estimation equation was applied to the NDSI and NDVI map to calculate the SCF map. Then, the coefficient of the equation for estimating the snow depth was derived using the in-situ data and the value of the calculated SCF map. Because the SCF and the snow depth are known to be related to the exponential function, the coefficient was derived based on this equation.
Then, we performed a post-correction to increase the accuracy of the calculated snow depth map. Generally, in the Republic of Korea, it is rare that more than 10 cm of snowfall is observed. Therefore, the maximum value of the derived estimation of the snow depth equation was relatively low. This caused a large error if snow depths of more than 10 cm were actually observed, thus, post-correction must be performed. Assume that the condition of the image is nearly constant and there are two points A and B. When the points A and B are independently considered, the difference between the measured value and the calculated value does not have a correlation, but the result of subtracting the measured value of B from that of point A and the result of subtracting the calculated value of B from that of A have a correlation. Using this, we derived a linear equation that corrected the calculated value and generated the final snow depth map using this equation. Then, the accuracy of the final result was validated by comparing it to the measured value.
4. Results
Fig. 4 shows a coastline area which make it easy to compare before and after processing of the co-registration step. Fig. 4(a) is the original position of the MODIS and Landsat images. The green-colored area is represented by a pseudo color image of MODIS and beneath it is the Landsat image covering the same area. As can be seen from Fig. 4(a), both satellite images do not exactly correspond although they cover the same region. An NDSI map of both images has a correlation because it was acquired on the same day. Thus, we calculated the correlation coefficient from the two NDSI images and derived a co-registration polynomial equation from it. Using this polynomial equation, the MODIS image of Fig. 4(a) was transformed to Fig. 4(b). Prior to the processing, the MODIS pixel does not cover the same area as the Landsat pixel. The coastline of MODIS was shifted in an inland direction as shown in Fig. 4(a) but it corresponds after processing as shown in Fig. 4(b). After processing of the co-registration step, deriving an accurate fitting equation was possible because the position of both images covered exactly the same area which means the value of SCF from the NDSI or NDVI was accurately plotted.
Fig. 4. Before and after image of co-registration processing. (a) Before the co-registration the coastlines of the two images are not matched well. (b) After the co-registration processing the coastlines match each other.
The co-registered NDSI and NDVI maps of the MODIS image were compared to the same position of the true SCF map. Fig. 5 shows the relationship of the NDSI with the SCF and the NDVI with the SCF accumulated from 38 pairs of images via scattergram. In both scattergrams, the SCF has a value of from 0 to 1, the NDSI value is between -0.4 and 1, and the value of NDVI is between -0.2 and 0.6. As the NDSI increases, the SCF tends to increase, and as the NDVI increases, the SCF tends to decrease. This means that the SCF increases when the influence of the snow is dominant in the pixel value and that the SCF value decreases even though there is the same snowfall in a region where the vegetation ratio is high.
Fig. 5. Scatter grams of variables and SCF. (a) Relationship of NDSI and SCF. (b) Relationship of NDVI and SCF.
To accurately grasp the relationship between the two variables and the SCF, the scattergram was reconstructed using a representative value for each interval, as shown in Fig. 6. The gray dot in Fig. 6 shows the distribution of the whole data and the black dot is the representative value in the corresponding interval. The whole data are divided into 1000 sections according to values. Then, the median value belonging to each section was used as the representative value. Reviewing the tendency of the graph expressed by the black dot, it can be seen that the SCF gradually increases while the NDSI increases from 0.2 to 0.4, and in the section where the NDSI is 0.4 or greater, the slope of the graph sharply increases and the rate of change in the SCF increases. This tendency continues to an NDSI value of approximately 0.65 Thereafter, the slope gradually decreased, such that the increase in the SCF according to the increase in the NDSI decreased. In the case of the NDVI, there is slight decrease in the SCF while increasing to 0.1 of NDVI, but the value rapidly decreases between 0.1 and 0.15. Then, the value gradually decreases with an increase in the NDVI until 0.2. It is more appropriate to explain this increase/ decrease pattern based on the Gaussian function rather than the previously suggested equation based on a linear, secondary, or exponential function. Therefore, in this study, we proposed an equation to estimate the SCF using the Gaussian function based an equation considering both the NDSI and NDVI.
Fig. 6. Scatter diagrams: (a) NDSI vs. SCF and (b) NDVI vs. SCF.
Fig. 7 shows the relationship between the NDSI and SCF based on each equation. Fig. 7(a) shows the relationship using the previously proposed linear equation, and Fig. 7(b) shows the relationship using the exponential function. Finally, Fig. 7(c) shows the relationship between the NDSI and SCF in the form of a Gaussian function proposed n this study. Equations derived from a linear equation, an exponential function, and a Gaussian function correspond to Eq. 3, Eq. 4, and Eq. 5, respectively, as follows:
Fig. 7. Derived fitting graph of the NDSI to SCF based on each equation. (a) Based on a linear function. (b) Based on an exponential function. (c) Based on a Gaussian function. The red lines denote the fitted functions.
\(S C F=\left\{\begin{array}{cl} e^{-18.16(N D S I-0.73)^{2}} & N D S I \leq 0.73, \\ 1 & N D S I>0.73 \end{array}\right.\) (5)
When the correspondence between the representative value of the interval and the defined formula (red line) was compared, the linear equation was well matched in the interval where the SCF sharply increases, but tended to be under-estimated and over-estimated at small and large NDSI values, respectively. I the case of estimating with an exponential function, it can be seen that the tendency to underestimate increases with increasing NDSI value. However, the Gaussian-based equation had the highest matching with the interval representative value in the whole range. Therefore, in this study, we used the Gaussian-based equation to estimate the SCF using the NDSI and NDVI as variables.
In this manner, the equations for calculating the SCF from the NDSI and the equation for calculating the SCF from the NDVI were derived (Fig. 8). The equation using the NDVI is defined as Eq. 6 and defined to be 1 below a certain value (0.07). Likewise, the NDSI(0.73) greater than a predetermined value was defined to have a value of 1. The horizontal bar of Fig. 8 means the standard deviation range of the NDSI and NDVI from which the SCF value can be calculated. Both the DSI and NDVI have a large standard deviation range in a small value range of the SCF, but the range gradually decreases as the value of the SCF increases. That is, the greater the snow accumulation, the more accurate the SCF estimation.
Fig. 8. SCF fitting graph using the NDSI and NDVI; the horizontal bar represents the standard deviation of the independent variable: (a) Fitting graph using the NDSI. (b) Fitting graph using the NDVI.
\(S C F=\left\{\begin{array}{cl} 1 & N D V I \leq 0.07 \\ e^{32037(N D V I-0.07)^{2}} & N D V I \leq 007 \end{array}\right.\) (6)
The equation for estimating the SCF based on the two equations (Eq. 5 and Eq. 6) can be differently defined according to the region of the NDSI and NDVI as Eq. 7. As the equation shows, the ratio of the NDSI and NDVI affecting the determination of the SCF value was approximately 6:4. If the NDSI is greater than 0.68, the value of SCF is determined by the NDVI value and shows a value of at least 0.72. In addition, if the value of the NDVI is less than 0.06, the value of the SCF will be at least 0.28 according to the value of NDSI as follows:
\(\begin{array}{c} S C F= \\ \left(\begin{array}{c} 0.58 e^{23.1(N D S I-0.68)^{2}}+0.42 e^{\left(-286.68(\mathrm{NDV} \cdot 0.06)^{2}\right.} \\ N D S I \leq 0.68, N D V \geq 0.06 \end{array}\right. \\ 0.58 e^{23.1(N D S I-0.68)^{2}+0.42} \\ \quad N D S I \leq 0.68, N D V I<0.06 \\ 0.58+0.28 e^{286.68(N D V I-0.06)^{2}} \\ N D S I>0.68, N D V \geq 0.06 \\ 1 \\ N D S I>0.68, N D V I<0.06) \end{array}\) (7)
The results of estimating the SCF using the derived equations are summarized in Table 1. As can be seen, the equations based on a linear function and an exponential function show a similar graph shape, and the RMSE obtained by validating the estimated SCF using the two equations has a similar value. The results of the SCF computation using the Gaussian function-based equation show better results than those of the other equation. An RMSE of 0.26 was recorded using the two equations in the previous study, while an RMSE of 0.24 was recorded using the proposed equation in this study. This was interpreted that the overall trend is best explained by the proposed formula, but the standard deviation of the NDSI and NDVI for each SCF interval is large, such that there is no significant difference in the calculated RMSE.
Table 1. RMSE calculated according to the fitting function
The relationship between the SCF value calculated using the proposed SCF estimation equation and the measured snow depth is shown in Fig. 9. In general, the snow depth is known to increase in the form of an exponential function as the SCF increases; thus, we derived the exponential function based on the equation of snow calculation. Therefore, the equation for calculating the snow depth from the SCF is as follows. In the Republic of Korea, snowfall is generally less than 10 cm and it is difficult to observe snow because of fog, clouds, etc. because the area where snowfall of 10 cm or more occurs is mostly mountainous. Therefore, the depth of the snow in the constructed data was not as deep: thus, the maximum depth of snow that can be calculated in the derived equation was limited as follows:
Fig. 9. Relationship of the calculated SCF and snow depth. The red line is an equation that calculates the snow depth from the SCF.
\(\mathrm{SD}=6.95 *\left(e^{0.67 * \mathrm{SCF}}-1\right)\) (7)
When the value of the snow depth was calculated using Eq. 8 and compared to the measured value, the accuracy was an RMSE of approximately 3.57 cm. However, because the maximum value to be calculated was limited, the difference between the in-situ and calculated snow depths was increasingly greater as the measured value increased (Fig. 10(a)). When considering the value of the whole image, it was difficult to correct this difference. However, because the condition can be assumed to be same in a single image, the error between the measured value and the calculated value has a similar bias for each pixel. The correction equation can be derived by using the relationship of the pixels. Therefore, the correction equation was derived from the relationship between all of the measured values of the image and the corresponding calculated values of the pixel. Post-processing of the calculated snow depth was performed using this relationship. In the case of the corrected image, the RMSE (2.98 cm) was further reduced compared to that of the pre-correction map (3.55 cm). It was possible to provide a more accurate snow depth map because it was possible to correct not only the value at the in-situ location but also all the pixel values of the image.
Fig. 10. Comparison of the in-situ and calculated values before and after post-processing. () Before postprocessing. (b) After post-processing.
Fig. 11 shows the intermediate output and the final output according to the processing flow of the algorithm developed in this study. Fig. 11(a) shows the MODIS RGB image of the Republic of Korea acquired on January 8, 2013. The snow depth map was developed by applying the algorithm to the NDSI and NDVI map of this image. Fig. 11(b), Fig. 11(c), and Fig. 11(d) correspond to the calculated SCF, snow depth, and corrected snow depth maps, respectively. The RMSE of the corresponding image was approximately 5.19 cm before correction but 3.32 cm after correction.
Fig. 11. Intermediate output and final output using the proposed algorithm. (a) RGB image of the study area. (b) Estimated SCF map. (c) Calculated snow depth map. (d) Corrected snow depth map
5. Conclusions
This study has presented an advanced methodology using optical images to estimate snow depth. In the case of optical images, the snow depth is estimated using indirect methods because the reflectivity does not directly change according to the amount of snowfall in the visible region with a short wavelength. Because the SCF had a positive correlation with snow depth, it could be used to estimate the snow depth. Therefore, accurate production of the SCF map of the image is the most important step in estimating the depth of snow from optical images. Because the SCF has characteristics that are difficult to directly calculate in a vast area and the value varies depending on the characteristics of the surface, many studies have been conducted to estimate the value using the index obtained from the image. In this paper, we proposed a Gaussian-function-based estimation equation using the NDSI and NDVI as independent variables. When comparing the results of the proposed equation to the results of the existing equation, the underestimation and overestimation regions were significantly reduced. Therefore, it is possible to estimate an SCF value near the true SCF value. The proposed equation also showed the best results in the calculation of the RMSE of the estimated SCF. However, although the proposed equation most closely matched the increasing and decreasing trends of the SCF value, it can be confirmed that the difference was not large. This was interpreted as a problem caused by the large fluctuation range of the NDSI and the NDVI. In other words, because the NDSI and NDVI were not correctly calculated, they are the primary cause of the errors in estimating the SCF value, and this problem seems to affect the RMSE calculation as well. Therefore, to estimate a more accurate SCF, it is important to first obtain an accurate reflectance value for each pixel. In the case of optical images, it is difficult to correct the effect of the atmosphere. Therefore, a more advanced atmospheric correction technique is needed because the remaining atmospheric effect effects the reflectance value. In addiion, because the NDSI value is high in a water region, an error may occur because of a water body existing in a pixel or adjacent pixels. There is also a factor that may cause error in the estimation of the true SCF. Snow-covered pixels and snow-free pixels may be misclassified when producing a snow cover map, which is an error factor when deriving the SCF estimation equation because it distorts the true SCF value. In addition, because pixels such as haze, thin clouds, and fog are difficult to mask, it is not possible to exclude the possibility that a pixel having a distorted value exists in the image.
The snow depth estimation algorithm proposed in this study resulted in an RMSE of approximately 3.55 cm. Considering that the Republic of Korea does not have abundant snowfall, the result is not perfect. This is a limitation in the algorithm that indirectly estimates the snow depth. Also, the in-situ data used for the validation is point data, while the depth of the snow obtained from the image is calculated at a resolution of 1 km. In addition, if snow removal is performed or there is heavy snowfall in small areas, the error may occur because the pixel value does not properly reflect the snow depth. In this study, a final accuracy of the RMSE of approximately 2.98 cm was obtained through post-correction. It is expected that a more accurate correction value can be obtained if the number of measurement points available for the image increases.
As a result, this study proposed an optimized algorithm to estimate the depth of snow in the Republic of Korea and its accuracy is satisfactory. Following this study, we can summarize a direction for follow-up study to improve the accuracy as follows: 1) calculate the reflectance value of the pixel using an advanced atmospheric correction technique. 2) study additional variables that can be considered in estimating the SCF, and 3) design a more accurate algorithm through fusion with a microwave image. It is expected that it is possible to provide snow depth data required in various areas of research by simultaneously using the results of this study and in-situ data, and it is expected that the quality of the data will be further improved through subsequent studies.
Acknowledgement
This work was supported by the “Development of Scene Analysis & Surface Algorithms” project, funded by the Electronics and Telecommunications Research Institute, which is a subproject of the “Development of Geostationary Meteorological Satellite Ground Segment (NMSC-2017-01)” program funded by the National Meteorological Satellite Center (NMSC) of the Korea Meteorological Administration.
References
- Appel, I.L. and V.V. Salomonson, 2002. Estimate of fractional snow cover using MODIS data, Proc. of 2002 International Geoscience and Remote Sensing Symposium, Toronto, Canada, Jun. 24- 28, vol. 5, pp. 3044-3046.
- Aschbacher, J., 1989. Land surface studies and atmospheric effects by satellite microwave radiometry, University of Innsbruck, Innsbruck, Austria.
- Baker, D.G., R.H. Skaggs, and D.L. Ruschy, 1991. Snow depth required to mask the underlying surface, Journal of Applied Meteorology, 30(3): 387-392. https://doi.org/10.1175/1520-0450(1991)030<0387:SDRTMT>2.0.CO;2
- Chang, A.T.C., J.L. Foster, and D.K. Hall, 1987. Nimbus-7 SMMR derived global snow cover parameters, Annals of Glacialogy, 9: 39-44. https://doi.org/10.1017/S0260305500000355
- Davis, D.T., Z. Chen, L. Tsang, J.N. Hwang, and A.T. Chang, 1993. Retrieval of snow parameters by iterative inversion of a neural network, IEEE Transactions on Geoscience and Remote Sensing, 31(4): 842-852. https://doi.org/10.1109/36.239907
- Foster, J.L., A.T.C. Chang, and D.K. Hall, 1997. Comparison of snow mass estimates from a prototype passive micro wave snow algorithm, a revised algorithm and a snow depth climatology, Remote Sensing of Environment, 62(2): 132-142. https://doi.org/10.1016/S0034-4257(97)00085-0
- Goita, K., A.E. Walker, B.E. Goodison, 2003. Algorithm development for the estimation of snow water equivalent in the boreal forest using passive microwave data, International Journal of Remote Sensing, 24(5): 1097-1102. https://doi.org/10.1080/0143116021000044805
- Liang, J., X. Liu, K. Huang, X. Li, X. Shi,Y. Chen, and J. Li, 2015. Improved snow depth retrieval by integrating microwave brightness temperature and visible/infrared reflectance, Remote Sensing of Environment, 156: 500-509. https://doi.org/10.1016/j.rse.2014.10.016
- Lin, J., X. Feng, P. Xiao, H. Li, J. Wang, and Y. Li, 2012. Comparison of snow indexes in estimating snow cover fraction in a mountainous area in northwestern China, IEEE Geoscience and Remote Sensing Letters, 9(4): 725-729. https://doi.org/10.1109/LGRS.2011.2179634
- Romanov, P. and D.Tarpley, 2007. Enhanced algorithm for estimating snow depth from geostationary satellites, Remote Sensing of Environment, 108(1): 97-110. https://doi.org/10.1016/j.rse.2006.11.013
- Romanov, P. and D. Tarpley, 2004. Estimation of snow depth over open prairie environments using GOES imager observations, Hydrological Processes, 18(6): 1073-1087. https://doi.org/10.1002/hyp.5508
- Romanov, P., D. Tarpley, G. Gutman, and T. Carroll, 2003. Mapping and monitoring of the snow cover fraction over North America, Journal of Geophysical Research, 108: D16.
- Salomonson, V.V. and I. Appel, 2004. Estimating fractional snow cover from MODIS using the normalized difference snow index, Remote Sensing of Environment, 89(3): 351-360. https://doi.org/10.1016/j.rse.2003.10.016
- Tedesco, M., J. Pulliainen, M. Takala, M. Hallikainen, and P. Pampaloni, 2004. Artificial neural network-based techniques for the retrieval of SWE and snow depth from SSM/I data, Remote Sensing of Environments, 90(1): 76-85. https://doi.org/10.1016/j.rse.2003.12.002
- Zahir, N. and H. Mahdi, 2015. Snow Depth Estimation Using Time Series Passive Microwave Imagery via Genetically Support Vector Regression (case Study Urmia Lake Basin), Proc. of International Conference on Sensors & Models in Remote Sensing & Photogrammetry, Kish Island, Iran, Nov. 23-25, vol. 40, pp. 555-558.