Summary:

Most current land surface assimilation systems are basically single point assimilation. Single point assimilation can easily break the coherent large-scale spatial structures of soil moisture anomaly, which are usually the important land surface factors to influence short-term climate variations over the land. The manuscript entitled “Development and preliminary validation of a land surface image assimilation system based on the common land model” propose an image assimilation method by using the curvelet transform to denoise the observational data with only the primary structural information to be assimilated. Preliminary results showed that this assimilation method can adjust the structures of model soil moisture based on the observed spatial structure characteristics, increasing the spatial similarity of soil moisture between the model and the observation. This image assimilation method shows potential for improving the forecast of short-term climate variability related to soil moisture anomalies. However, benefit of the image assimilation is not well evaluated. That is, the paper should show more results regarding the advantages of the image assimilation over traditional single point assimilation. The paper is generally well built up. However, still this manuscript needs to be improved greatly, especially regarding the issues mentioned above. Efforts should be made to improve the readability. I think this paper can be considered for publication after some issues/questions are resolved/explained.

Response:

Thanks to your valuable comments and suggestions. Following your suggestions, we have revised the manuscript carefully from the beginning to the end. The point-by-point response is listed below according to your specific comments.

Major issues:

1. There are many other image denoising techniques, why use curvelet for land surface images? Curvelets are anisotropic, they have a high directional sensitivity and are very efficient in representing vortex edges. Therefore, the curvelet transform is suitable for geophysical fluids. But what is the argument for choosing it for land surface?

Response:

As the reviewer highlighted, the basis function of curvelet analysis exhibits anisotropic characteristics, thereby demonstrating its exceptional capability in accurately reproducing the rapidly-evolving properties of earth fluids. Although the soil moisture does not exhibit rapid temporal variations, there are many small-scale spatial structures of soil moisture due to the high spatial heterogeneity of soil. Therefore, the curvelet analysis method is selected as a more effective approach to capture the intricate local variations of soil moisture.

Indeed, a variety of mathematical image analysis techniques are available. For instance, the Fourier decomposition method and wavelet analysis method are commonly employed in meteorological research. However, the Fourier analysis method primarily focuses on the average feature of the sequence at different frequencies, and lacks the ability to accurately describe the regional variations. The Wavelet analysis could provide more detailed variation information in the time-frequency domain, but its basis functions with isotropic characteristic limit the ability to accurately represent the characteristics of small-scale spatial variations.

 On the other hand, the curvelet analysis method has been selected to fulfill the requirements of variational data assimilation. The curvelet transform is an observation operator in an image assimilation system. During the process of minimizing the cost function in variational data assimilation, the adjoint function of the observation operator becomes necessary. The adjoint function of curvelet analysis is just its inverse transformation, which proves to be a highly advantageous property for minimizing the cost function in a variational data assimilation system.

1. The argument for choosing the threshold of 0.5 for the curvelet denoising is not convincing enough. Probably different thresholds will lead to different assimilation results. If it is true. How should understand this?

Response:

Thanks for your valuable suggestions. Just as the reviewer pointed out, the image assimilation system determines the spatial structural characteristics of assimilation according to the threshold values, and different threshold values could result in certain variations in the spatial structure of assimilation.

Naturally, a higher threshold can capture more spatial structural features of the observed variables, but the presence of observation errors imposes limitations on its continuous increase. More discussions have been added to the revised manuscript in Line 366-404 to prove that a threshold of 0.5 can not only capture the spatial structure information of observation data, but also mitigate the impact of observational errors. The specific content comprises the following three aspects:

(1) The definition of the threshold σ has been further elaborated in order to provide a clearer rationale for its selection. This detailed description has been incorporated into line 370-371 of the revised manuscript. The specific wording is as follows:

The threshold σ means the modulus of the decomposition coefficient falls within the first 100*σ% percentile. For instance, a value of 0.5 indicates that the mode retaining the top 50% of decomposition coefficient.

(2) By employing the spatial correlation method, we demonstrate that a threshold of 0.5 adequately captures the primary spatial information derived from soil moisture observations, the following discussions have been added to Line 366-404 of the revised manuscript:

The image assimilation system finds the spatial structural characteristics of assimilation according to the threshold values, and different thresholds could result in certain variations in assimilated spatial structure. In order to clarify the spatial structure differences corresponding to different thresholds, the spatial correlation method (Daley, 1991) is employed in this study to elucidate the distinctive characteristics of spatial structure corresponding to varying thresholds.

The hourly soil moisture data from ERA5-Land from May 1 to 30, 2016 is selected for analysis. The threshold σ means the modulus of the decomposition coefficient falls within the first 100*σ% percentile. For instance, a value of 0.5 indicates that the mode retaining the top 50% of decomposition coefficient. The original image can be reconstructed by selecting different threshold ranges, namely (0,0.01], (0.01,0.03], (0.03,0.05], (0.05,0.1], (0.1,0.2], (0.2,0.3], (0.3,0.4], (0.4,0.5], (0.5,0.6], (0.6,0.7], (0.7,0.8], (0.8,0.9] and (0.9,1.0]. The correlation coefficient between each grid point and its neighboring grid points can be obtained based on the reconstructed time series of each grid point. The spatial structural characteristics of different scales in the reconstructed images could be quantitatively expressed by the average correlation coefficients corresponding to different grid point distances.

The mean correlation coefficient corresponding to grid point distance is illustrated in Figure 4. As can be seen, the variation characteristics of the inter-grid correlation coefficient of the original soil moisture represented by the black line with respect to the grid distance. The average correlation coefficient can exceed 0.5 within a radius of 200 km, while maintaining above 0.4 within a radius of 300 km. The distance corresponding to high correlation coefficients represents the characteristics of consistent changes in soil moisture within a similar range, that is, soil moisture has the characteristics of spatial structure at the corresponding scale.

When the threshold value is 0.01, the average correlation curve exhibits a similar change in correlation coefficient of the original variable, thereby indicating that the curvelet coefficient corresponding to this threshold value effectively reproduces the large-scale spatial structure. The spatial structure scale represented by the corresponding curvelet transformation reconstruction results decrease as the threshold value increases, leading to a rapid decrease in the correlation coefficient with increasing distance. The curvelet reconstruction results with different threshold intervals represent the structural characteristics of different horizontal scales, while the cumulative threshold can well represent the spatial structural characteristics of soil moisture variables represented by the selected threshold in the assimilation. The average correlation coefficient of the cumulative threshold is depicted in Figure 4b. As can be seen, the top 10% of curvelet coefficients can effectively replicate the spatial correlation characteristics of soil moisture variables. The results also indicate that the variations in threshold values have minimal impact on the assimilated spatial structure when the threshold value exceeds 0.1.

Figure 4: Variation curves of the average correlation coefficient between grid points with the distance in the reconstructed ERA5-Land hourly soil moisture image of the study area from May 1 to 30, 2016, which is reconstructed based on the curvelet coefficients of (a) different threshold intervals and (b) cumulative thresholds.

(3) According to the stochastic characteristics of observation errors, we conducted an analysis on the probability distribution properties of the reconstructed residuals and found that a threshold value of 0.5 effectively mitigates the impact of observation errors.

Naturally, a higher threshold can effectively capture more spatial structural features of the observed variables, but the presence of observation errors imposes limitations on its continuous increase. The observational error is typically characterized by stochastic fluctuations. When the discrepancy between the reconstructed results and the original variables exhibits random variation characteristics, it can be inferred that the observation information eliminated by the threshold method primarily consists of observation errors.

To better clarify the statistical characteristics of the reconstruction errors under different thresholds, Figure 5 shows the probability density distribution curves of the reconstruction errors for 100 reconstructed fields at different thresholds. For the error at the threshold of 0.5, the skewness coefficient of the probability density distribution curve is 0.00 and the kurtosis coefficient is 0.38, indicating the curve is close to the standard normal distribution curve (the skewness and kurtosis coefficients are all 0). With the gradual increase of threshold value, although the reconstruction error decreases, the residual error is mainly concentrated in the range of smaller values, and the curve shows a "sharp peak" distribution. Considering that the observation errors are mostly random errors, it is reasonable to believe that the reconstruction errors at the threshold of 0.5 are mainly observation errors, which also implies this threshold is good for the purpose of denoising the observation images.

Figure 5: Probability density distributions of 100 reconstructed errors under different thresholds. The magenta dashed line represents standard normal distribution. The red, blue, black, green and orange solid lines represent threshold values of 0.4, 0.5, 0.6, 0.7 and 0.8, respectively.

The revised manuscript now includes the newly added Figure 4, which has been incorporated into the section on observational error analysis along with its corresponding discussion.

And the relevant reference has been added to the revised manuscript:

Daley, R.: Atmospheric Data Analysis, Cambridge University Press, Cambridge, 1991.

1. Don’t understand why there is no error covariance matrixes involved in the term J_1 in Equation (5).

Response:

As the reviewer emphasized, the data assimilation typically involves the covariance of observation error and background error. But in this study, the image term in the cost function serves only as a weak constraint to adjust the spatial structure of the analysis field within the image assimilation system, thereby the error covariance is not necessary. Additionally, as explained in response to question 2, we have elucidated in detail that the significance of setting a threshold value for effectively filtering out erroneous information from the observed image.

Minor comments:

1. No information of the used atmospheric forcing data.

Response:

The information of used atmospheric forcing data is given in Line 166:

Atmospheric forcing conditions provide constraints on land-surface models. The quality of atmospheric forcing data greatly affects the ability of land surface models to realistically simulate land surface conditions. The atmospheric forcing dataset used to drive the CoLM in this study includes the downward short-wave solar radiation at surface, downward long-wave radiation, near-surface air temperature, specific humidity, precipitation rate, surface atmospheric pressure, U-component wind speed, and V-component wind speed. It has a temporal resolution of three hours (at 0000 UTC, 0300 UTC, 0600 UTC, etc.) and the spatial resolution is T62 (about 1.875°) (Qian et al., 2006). The forcing dataset was derived through combining observation-based analyses of monthly precipitation and surface air temperature with intramonthly variations from the National Centers for Environmental Prediction-National Center for Atmospheric Research (NCEP-NCAR) reanalysis. To correct the spurious long-term changes and biases in the NCEP-NCAR reanalysis precipitation, surface air temperature, and solar radiation fields, Qian et al. (2006) combined the intramonthly variations from the NCEP-NCAR 6 hourly reanalysis with monthly time series derived from station records of temperature and precipitation. It is shown that the CLM3 reproduces many aspects of the long-term mean, annual cycle, interannual and decadal variations when it was forced by this dataset.

We have added the relevant reference to the revised manuscript：

Qian, T., Dai, A., Trenberth, K. E., and Oleson, K. W.: Simulation of global land surface conditions from 1948 to 2004. Part I: Forcing data and evaluations, J. Hydrometeorol., 7, 953– 975, doi:10.1175/JHM540.1., 2006.

1. Figure 9: the correlations generally decline until the middle of July and then increase, how to understand this?

Response:

The occurrence of this phenomenon is attributed to the amount of precipitation in the driving data. To elucidate this matter, we superimpose the temporal variation of precipitation in the forced data within the study period (indicated by gray shading) on Figure 10 of the original manuscript in Line 398.

The following discussions have been added to Line 593-600 of the revised manuscript:

It is important to note that the SCC exhibits a clear temporal variation, which does not necessarily imply a time-varying assimilation effect. This can be attributed to the dominant influence of precipitation on the changes in the SCC. Hence, Figure 13 also includes the hourly total precipitation (represented by grey bars) in the model domain. The changes in precipitation exhibit a strong correlation with the SCC. From May 16 to June 15, there is minimal precipitation, corresponding to sustained high SCCs of soil moisture (red line) after assimilation. Subsequently, as precipitation increases, the SCC gradually diminishes. From August 15 to September 1, the SCC exhibits an inverse variation with decreasing precipitation.

Figure 13: Hourly variations of the spatial correlation coefficient of the surface (red and blue solid lines) and subsurface (black and gray solid lines) soil moisture between the observations and the experiments with (black and red solid lines) and without (black and gray solid lines) image assimilation, and the precipitation in the forced data was indicated by gray shading. After the vertical dashed line, it is the prediction period.

1. line 62-67: the sentences are ambiguous and hard to follow, please clarify to be concise and accurate.

Response:

Thanks for your valuable suggestions. We have revised these sentences to make them clear and concise.

The sentence “However, in reality, the observation quality varies sharply across regions, and the strong spatial heterogeneity of soil variables also tends to cause large spatial variations in the accuracy of surface variables simulated by the land surface model (Li, 2013; Li et al., 2020b). This leads to the regional differences in the accuracy of the estimations of observation error and background error in the single-column assimilation, and ultimately causes discontinuities in the spatial structure of the anomalies in the analyzed soil moisture fields.” has been revised as follows.

Due to the non-uniform spatial distribution of precipitation, as well as the heterogeneous spatial distribution of soil properties, land cover types and topographic elevations, there are significant variations in the spatial distribution of soil moisture (Tian et al., 2021). The estimation of soil moisture by the land surface model is adversely impacted by the uncertainties in atmospheric forcing, model dynamics and parameterization, leading to significant spatial variations in the accuracy of simulated surface variables (Li, 2013; Li et al., 2020b). Furthermore, there are regional differences in the accuracy of the estimation of the observation error and the background error resulting from the single column assimilation, which ultimately contribute to the discontinuity of the abnormal spatial structure in the analyzed soil moisture field.

And the new reference has been added to the revised manuscript:

Tian, S., Renzullo, L. J., Pipunic, R. C., Lerat, J., Sharples, W., Donnelly, C.: Satellite soil moisture data assimilation for improved operational continental water balance prediction, Hydrol. Earth Syst. Sci., 25(8): 4567-4584, https://doi.org/10.5194/hess-25-4567-2021, 2021.

1. Section 3.1: It is indicated in line 230 that resolution of the soil moisture reanalysis data is 31 km, while line 235 states “increased to 9 km”.

Response:

Thanks for your valuable suggestion. The description of line 230 of the original manuscript has now been revised to “9 km”. In line 235, “increased to 9 km” means the resolution of ERA5-land is increased from 31 km (original resolution of ERA5) to 9 km.

1. line 252: To demonstrate the benefit of the image assimilation and evaluate its advantages over traditional single point assimilation, if set J_O=0 in equation (5), authors should do one more set of experiments performing single point assimilation. Another option is to first do J(x) = J_B + J_O as conventional assimilation, and then do J(x) = J_B + J_O + J_I to see the benefit of the image assimilation.

Response:

Thanks to your valuable suggestion. According to your opinions, we have added a specialized section in the revised manuscript to facilitate a comparative analysis of the disparity in the effectiveness between the prevailing single point assimilation method and image assimilation method. The following analysis results have been added into the revised manuscript. Please refer to Line 419-458 for detailed information.

3.4 The influence of image assimilation constraints

Two sets of ideal experiments are designed to validate the impact of image assimilation and evaluate its superiority over traditional single column assimilation in adjusting the spatial distribution structure of soil moisture. The ideal observational data for assimilation is the ERA5-Land reanalysis soil moisture. The first set corresponds to the conventional assimilation experiment, where , as described by Equation (5). Another set is image assimilation experiment, where  in Equation (5), indicating that .

The process of data assimilation entails leveraging the discrepancy between observed data and background field, in conjunction with a priori knowledge of observation error and background error, to derive an analysis field that closely approximates the true value. The primary challenge in single column assimilation lies in acquiring precise prior information regarding observation error. The spatial distribution of observation error for a specific single column assimilation experiment is illustrated in Fig. 6. In consideration of the necessity for an ideal experiment, it is assumed that the observation error outside the China region is negligible, while a significant error is presumed within the China region, so as to emphasize the impact of observation error on assimilation results.

Figure 6: Spatial distribution of observation errors.

The spatial distributions of soil moisture for the ideal observation data and different experiments at 0000 UTC on May 1, 2016 are given in Figure 7. The spatial distribution of surface soil moisture in ERA5-Land is illustrated in Fig. 7a. The northern Siberian region of the selected area exhibits a relatively high soil moisture content overall, with a ring-shaped distinct wet zone in the northwest. The central region stretching from Xinjiang to western Mongolia is a significant arid area. However, the soil moisture in the Tianshan Mountains is wet. The soil moisture of the Qinghai-Tibet Plateau region gradually decreases from west to east. The soil moisture in southern Qinghai, Hunan and Jiangxi is characterized by high level of saturation, while Gansu, Ningxia and Hebei experience relatively arid soil conditions. Figure 7b is the distribution of soil moisture in the control experiment (background field). It is evident that there are significant disparities in the spatial distribution of soil moisture when compared with the reanalysis data. In the control experiment, a dry region extends from west to east in the northern area of Lake Baikal, while eastern Kazakhstan and central Inner Mongolia also exhibit arid conditions.

Figure 7c shows the results of the single point assimilation experiment. The observation error outside the China region is relatively minimal, indicating a strong correspondence between the analysis field and the observation data, and the overall distribution also exhibits a high degree of conformity with the observation. The analysis field in China region, however, closely resembles the background field. Nevertheless, there is a significant disparity between the observed soil moisture and that of the background field, indicating a lack of adjustment based on observed information.

Figure 7d is the assimilation results of the image assimilation ideal experiment. It is evident that image assimilation effectively adjusts the distribution pattern of soil moisture. The above-mentioned characteristics of moist soil moisture in the northwest region of the observation field, the arid region of Xinjiang and Mongolia, and the humid region of the Tianshan Mountains are all well reflected in the analysis field.