Evaluating the Relationship between Field Aerodynamic Roughness and the Modis Brdf, Ndvi, and Wind Speed over Grassland

Aerodynamic roughness (AR) is an important parameter that influences the momentum and energy exchange between the earth's surface and the atmosphere. In this study, profile wind data observed during the vegetation growing period (April–September) in 2013 and 2014 at the A'rou grassland station, which is in the upstream of the Heihe River Basin (HRB), were used to determine the relationship between the field AR and the Moderate-resolution Imaging Spectroradiometer (MODIS) near-infrared (NIR) bi-directional reflectance distribution function (BRDF) R index, the normalized difference vegetation index (NDVI), and a combination of these indices. In addition, the relationship between the average wind speed at a height of 1 m and the field AR is also presented. The results indicate that the MODIS NIR BRDF_R index and the NDVI are both sensitive indicators of the AR over grassland (R 2 : 0.5228 for NIR BRDF_R; R 2 : 0.579 for NDVI). Moreover, the combined index shows a significantly increased R 2 value of 0.721, which is close to the result inferred from the wind speed (R 2 : 0.7411). The proposed remote sensing-based combination index (CI) has the potential for use in evaluations of the AR over grasslands during growing season and its sensitivity can reach levels that are comparable to considering the effects of wind speed, which usually requires ground-based observations.


Introduction
Aerodynamic roughness (AR) plays a significant role in characterizing the interface between the earth's surface and the atmosphere.Developing an effective remote sensing-based AR estimation method will promote research on the temporal and spatial process of evapotranspiration, hydrology, land surface processes, and sand transport [1,2].Most of the current research on remote sensing-based AR estimations is focused on farmland, forest, and desert regions.Researchers have used the normalized difference vegetation index (NDVI), the airborne Light Detection and Ranging (LIDAR), multi-angle observations, radar backscattering coefficients, and a morphology method to estimate the AR over these land surfaces [3][4][5][6][7][8][9][10][11][12].
For grassland, recent work has primarily concentrated on ground-based experimental data analyses.Stull (1999) provided the referenced AR value for different land covers, including prairie grasslands [13].Ma et al. (2002) analyzed in-situ turbulent data over the grassy marshland surface of the Tibetan Plateau [14].Yang et al. (2003) investigated the AR at three sites-including Anduo PBL, Naqu, and MS3478 PAM, which are covered with sparse and short grasses-in the monsoon season based on wind and temperature profiles [15].However, thematic research on remote sensing-based AR estimations over grassland is seriously lacking, which restricts associated research on the interaction between grasslands and the atmosphere at a regional level.Although Brown et al. (2012) used airborne LIDAR data to estimate the AR in a mixed grassland prairie by the empirical relationship between wind profile-derived AR and airborne LIDAR-derived vegetation height [16], acquiring airborne LIDAR terrain coverage is an apparent limitation for large-scale applications.
This paper analyzed the relationship between the field AR over grasslands which was determined based on the profile wind data measured at the A'rou station, which is in the upper stream of the Heihe River Basin (HRB), and the temporal Moderate-resolution Imaging Spectroradiometer (MODIS) near-infrared (NIR) bi-directional reflectance distribution function (BRDF) R index, MODIS NDVI, a combination of the BRDF and NDVI indices, and the average wind speed during the vegetation growing period of 2013-2014.The proposed combination index (CI) in this paper provides a new scenario for remote sensing-based AR estimations over grasslands during growing season.

Research Site and Experiments
The research site lies in the Babao River sub-basin, which is part of the upper reaches of the HRB and has an area of approximately 2452 km 2 and an elevation ranging from 2640 to 5000 m.The annual average air temperature and precipitation are 1 • C and 270 to 600 mm, respectively.Grassland is the main form of vegetation cover and is mainly in the form of alpine meadow [17,18]     AR estimations over grassland is seriously lacking, which restricts associated research on the interaction between grasslands and the atmosphere at a regional level.Although Brown et al. (2012) used airborne LIDAR data to estimate the AR in a mixed grassland prairie by the empirical relationship between wind profile-derived AR and airborne LIDAR-derived vegetation height [16], acquiring airborne LIDAR terrain coverage is an apparent limitation for large-scale applications.This paper analyzed the relationship between the field AR over grasslands which was determined based on the profile wind data measured at the A'rou station, which is in the upper stream of the Heihe River Basin (HRB), and the temporal Moderate-resolution Imaging Spectroradiometer (MODIS) near-infrared (NIR) bi-directional reflectance distribution function (BRDF) R index, MODIS NDVI, a combination of the BRDF and NDVI indices, and the average wind speed during the vegetation growing period of 2013-2014.The proposed combination index (CI) in this paper provides a new scenario for remote sensing-based AR estimations over grasslands during growing season.

Research Site and Experiments
The research site lies in the Babao River sub-basin, which is part of the upper reaches of the HRB and has an area of approximately 2452 km 2 and an elevation ranging from 2640 to 5000 m.The annual average air temperature and precipitation are 1 °C and 270 to 600 mm, respectively.Grassland is the main form of vegetation cover and is mainly in the form of alpine meadow [17,18].The superstation is in A'rou Village of Qilian County, Qinghai Province (100.4643°E, 38.0473° N, 3033 m), as shown in Figure 1.The experimental site features comprehensive, precise, and multi-scale observations of land surface and hydrological processes in a cold region.Intensive and long-term profile meteorological observations of the wind, air temperature, and humidity at the study site were conducted at six layers: 1 m, 2 m, 5 m, 10 m, 15 m, and 25 m.The specific manufacturer and model of each observation instrument are shown in Table 1.

Methodology
Different remote sensing indices are introduced in this section, including the BRDF_R index, NDVI, and the CI.Furthermore, an in-situ AR calculation method based on profile wind data is also presented.

MODIS BRDF_R Index
The Ross-Li-Maignan model predicts the reflectance as the sum of three terms [9,10]: where θ s and θ v are the sun and viewing zenith angles, respectively; φ is the difference between the azimuth angles of the sun direction and the viewing direction; F 1 is the volume-scattering kernel based on the Ross-Thick function corrected for the hot-spot process, F 2 is the geometric kernel based on the Li-sparse-reciprocal function [9,10]; k 0 is the isotropic scattering component, which is equal to the bi-directional reflectance for a view zenith angle of θ v = 0 and a solar zenith angle of θ s = 0; k 2 is the coefficient of the Li-sparse-reciprocal geometric scattering kernel F 2 , which is derived from surface scattering and the geometric shadow casting theory; and k 1 is the coefficient for the Ross-Thick volume-scattering kernel F 1 , which is derived from radiative transfer models.The parameters k 0 , k 1 , and k 2 are included in the MODIS MCD43 B1 data and were provided by the MODIS working team.We used MODIS remote sensing data because they represent standard data products for the BRDF model parameters and the NDVI at a spatial resolution of 1 km with the frequent simultaneous observations of the Terra and Aqua satellite constellation.Furthermore, because of the homogeneity of the underlying surface at A'rou station, as shown in Figure 1, the 1 km pixel size did not introduce the apparent uncertainty of mixed pixels to the analysis and the results presented in this study.In the future, we will attempt to develop the BRDF model parameters using Belgian Proba-V satellite data, which have a spatial resolution of 100 m, in order to reduce the influence of mixed pixel.
The BRDF_R index, which reflects the geometric roughness of the land surface, is calculated as follows [12]: The eight-day composited MODIS MCD43 B1 data covering the experimental site in 2013 and 2014 were downloaded from the National Aeronautics and Space Administration (NASA) website (http://ladsweb.nascom.nasa.gov).The visible, near-infrared, and shortwave broad bands cover 0.3 to 0.7 µm, 0.7 to 5.0 µm, and 0.3 to 5.0 µm in sequence, respectively, and have a spatial resolution of 1 km.In this paper, the MODIS BRDF_R index in the near-infrared band was used to conduct the analysis according to the conclusion [12] that near-infrared and shortwave-infrared BRDF parameters show similar sensitivities for the AR to those of the visible band.The data were pre-processed, which included format transformation, band extraction, projection transformation, subset formation, Julian day transformation, Digital Number (DN) value type transformation, R factor computation, and singular value elimination.

NDVI
The NDVI was calculated based on the reflectance of the near-infrared and red bands.The formula is as follows: where ρ nir and ρ red are the reflectance from the near-infrared and red bands, respectively.MODIS Level 1B calibrated radiances data with a spatial resolution of 250 m (MYD02QKM) for 2013 and 2014 were also downloaded from the NASA website at the above address and then pre-processed using a data processing chain, including radiometric, geometric, atmospheric corrections, and cloud detection and NDVI calculation [21].

Combination Index
The CI was calculated based on the MODIS NIR BRDF_R index and NDVI using the following formula: where BRDF_R refers to the MODIS NIR BRDF_R index, which was calculated using Formula (2).The NDVI was calculated using Formula (3).

Field Aerodynamic Roughness Calculation
The profile data were pre-processed and quality controlled by removing singular and repeated values exceeding the range of physical possibility via date format unification and data arrangement into a standard Excel format for ease of use.The original wind profile data were collected at 10-min intervals, which resulted in 144 data records each day.The effects of the underlying surface were minimized by selecting the wind and atmospheric profiles above the landscape for the field AR calculations.The 10-min profile data were first used to calculate the instantaneous 10-min AR and then averaged into daily and eight-day average values for consistency with the remote sensing data.
Meteorological profile time series data from the wind profile stations were processed to acquire the field AR based on Monin-Obukhov similarity theory (MOST).Following Zhou et al., the wind velocities and air temperatures measured at different heights can be expressed as follows [12,22]: where u is the wind velocity; u * is the friction velocity, k is the von Kármán constant (equal to 0.4); z is the height at which the wind velocity and air temperature are observed; θ is the potential air temperature of each layer, which is defined with local surface pressure used as standard pressure; θ * is the friction temperature; z om is the aerodynamic roughness; d is the zero-plane displacement; z oh is the thermal roughness length; and θ 0 is the potential temperature near the surface.In addition, L is the Obukhov length, which is calculated as follows [12]: where g is the acceleration caused by gravity (equal to 9.8 m/s 2 ), T is the mean air temperature, p is the ambient air pressure, and P 0 is 101 kPa.The Richardson Number (R i ) is used to determine the atmospheric stability near the surface.When R i ranges from −0.03 to 0.03, it exhibits neutral stratification.When R i is less than −0.03, it exhibits unstable stratification.When R i is greater than 0.03, it exhibits stable stratification.
where ∂z is the height difference between the maximal observation height and the minimal height above the vegetation; ∂u is the wind velocity difference of the two heights; θ v is the virtual potential temperature at the observation height; and ψ m and ψ h are the stability correction functions for momentum and sensible heat transfer, respectively, which are determined based on atmospheric stratification.These functions for unstable conditions are expressed as follows [22]: x = (1 ) These functions for stable conditions are expressed as follows [23]: The profile data at the six levels of 1 m, 2 m, 5 m, 10 m, 15 m, and 25 m were all used in the calculation.The aerodynamic roughness z om and displacement height d are concurrently optimized and determined through iteration.The d was allowed to increase from 0.1 to 1 m in intervals of 0.01 m considering the maximum grass height around 1.1 m there.In the iteration, u * and θ * were given initial values, and then z om , u * , z oh , and θ * were determined for each given d.The L was then repeatedly calculated using new values of u * and θ * from the second step until the difference in L between two consecutive iterations was less than 0.01.For every time period, a series of AR and d estimates and corresponding fitting coefficients and correlation coefficients between the fitted and measured wind velocity and temperature data were determined.The AR and d with the highest correlation coefficient was selected as the optimal result for the temporal [12].

Results and Analysis
In this section, an example of the original observation data at various layers is firstly shown to demonstrate the consistency of the data.Afterwards, the change trend of the indicators, including the MODIS BRDF_R and NDVI, the average wind speed at 1 m height, and the in-situ AR, were presented on the eight-day average scale to ensure the synchronization of the ground measurements and remote sensing observations, followed by the relationship analysis between the in-situ AR and these remote sensing indicators.

Consistency Analysis of the Observed Data
An example of the original observed wind speed and air temperature data recorded at 10-min intervals on 1 April 2013 is shown in Figure 2a,b respectively, for heights of 1 m, 2 m, 5 m, 10 m, 15 m, and 25 m.Both the wind speed and air temperature display consistent variation trends at different layers over time, thus confirming the reliability of the original experimental instrument and data.At most times within a given day, a higher altitude is correlated with greater wind speeds.The difference in wind speeds among the various layers is relatively apparent.The wind speed itself shows relatively high fluctuations and randomness.At most times within a given day, the difference in air temperature among different layers is relatively apparent, and the variation trends differ as the observation height increases.During the night, the air temperature increases as the height increases, whereas during the day, the air temperature decreases as the height increases.The air temperature itself shows a clear and stable curve with a single peak.

Change Trend Analysis of the Indicators
The change profiles of the MODIS NIR BRDF_R index, NDVI, average wind speed, and in-situ AR from April to September in 2013-2014 are shown in Figure 3.
The four parameters of the MODIS NIR BRDF_R index, MODIS NDVI, average wind speed, and in-situ AR all show apparent variations in vegetation growth for the monitored period from 2013 to 2014.Variations in the NDVI (average value: 0.442; variation range: 0.113 to 0.671) directly reflect the

Relationship Analysis
The relationships between the in-situ AR and the MODIS NIR BRDF_R index, MODIS NDVI, and the CI are shown in this section.
The MODIS NIR BRDF_R index is negatively correlated with the in-situ AR and has an R 2 value of 0.5228, and the relationship is in the form of an exponential power function that reflects the geometric structure information for the vegetation [9][10][11].The MODIS NDVI is positively correlated with the in-situ AR and has an R 2 value of 0.579, which is comparable to that of the BRDF_R index, and the relationship is in the form of a power exponent that primarily reflects variations in the vegetation leaves [6,7].The results in Figure 4 show that the NIR BRDF_R index and NDVI are both sensitive indicators for characterizing the AR variation process over grassland.The four parameters of the MODIS NIR BRDF_R index, MODIS NDVI, average wind speed, and in-situ AR all show apparent variations in vegetation growth for the monitored period from 2013 to 2014.Variations in the NDVI (average value: 0.442; variation range: 0.113 to 0.671) directly reflect the growing process of the grass, and the variation trends in the NIR BRDF_R index (average value: 0.067; variation range: 0.009 to 0.159) are opposite to that of the NDVI. Figure 2 indicates that as the grass grows, the NDVI and AR will increase and the NIR BRDF_R index and average wind speed will decrease.The AR (average value: 0.0101 m; variation range: 0.0028 to 0.0234 m) shows similar variation trend as the NDVI, and the average wind speed (average value: 2.27 m/s; variation range: 1.785 to 2.967 m/s) shows consistent variation trend with the NIR BRDF_R index.

Relationship Analysis
The relationships between the in-situ AR and the MODIS NIR BRDF_R index, MODIS NDVI, and the CI are shown in this section.
The MODIS NIR BRDF_R index is negatively correlated with the in-situ AR and has an R 2 value of 0.5228, and the relationship is in the form of an exponential power function that reflects the geometric structure information for the vegetation [9][10][11].The MODIS NDVI is positively correlated with the in-situ AR and has an R 2 value of 0.579, which is comparable to that of the BRDF_R index, and the relationship is in the form of a power exponent that primarily reflects variations in the vegetation leaves [6,7].The results in Figure 4 show that the NIR BRDF_R index and NDVI are both sensitive indicators for characterizing the AR variation process over grassland.The CI based on the NIR BRDF_R and the NDVI shows a stronger relationship with the in-situ AR, which is also in the form of an exponential power function, as shown in Figure 5.As the CI increases below 2.0, the AR decreases rapidly.With CI values between 2.0 and 6.0, the AR decreases more slowly before converging.The R 2 value of 0.721 and the much higher aggregation degree indicate an obvious improvement of the fit compared with the result in Figure 4. Thus, the combination of the NIR BRDF_R index and the NDVI can provide a more accurate characterization of the AR dynamics related to grass growth.The proposed CI in this paper will provide significant guidance towards remote sensing-based AR estimations over grassland during growing season.The CI based on the NIR BRDF_R and the NDVI shows a stronger relationship with the in-situ AR, which is also in the form of an exponential power function, as shown in Figure 5.As the CI increases below 2.0, the AR decreases rapidly.With CI values between 2.0 and 6.0, the AR decreases more slowly before converging.The R 2 value of 0.721 and the much higher aggregation degree indicate an obvious improvement of the fit compared with the result in Figure 4. Thus, the combination of the NIR BRDF_R index and the NDVI can provide a more accurate characterization of the AR dynamics related to grass growth.The proposed CI in this paper will provide significant guidance towards remote sensing-based AR estimations over grassland during growing season.The CI based on the NIR BRDF_R and the NDVI shows a stronger relationship with the in-situ AR, which is also in the form of an exponential power function, as shown in Figure 5.As the CI increases below 2.0, the AR decreases rapidly.With CI values between 2.0 and 6.0, the AR decreases more slowly before converging.The R 2 value of 0.721 and the much higher aggregation degree indicate an obvious improvement of the fit compared with the result in Figure 4. Thus, the combination of the NIR BRDF_R index and the NDVI can provide a more accurate characterization of the AR dynamics related to grass growth.The proposed CI in this paper will provide significant guidance towards remote sensing-based AR estimations over grassland during growing season.

Discussion
The AR results in this paper were compared with the results over grassland from the literature.In the present study, the average AR value at the A'rou grassland station was 0.0101 m, and the AR values ranged of 0.0028 to 0.0234 m during the growing period.Stull (1999) provided a general reference AR value for prairie grassland of approximately 0.03 m [13].Ma et al. (2002) studied the AR over the grassy marshland surface of the Tibetan Plateau using in-situ turbulent data and obtained an AR value of 0.00466 m at Anduo for grass heights of 5 cm and a value of 0.0139 m at the HPAM for grass heights of 15 cm [14].Yang et al. (2003) determined the AR at Anduo PBL and Naqu, which are covered with sparse and short grasses, in the monsoon season based on wind and temperature profiles, and the order of magnitude of the AR was from 0.001 to 0.01 m [15].Thus, the AR results in the present paper are comparable to previously published in-situ results.
The relationship between the average wind speed at the height of 1 m and the in-situ AR on the eight-day scale (Figure 6) indicates that the average wind speed has a significant influence on the AR over grassland [7,22], with R 2 values reaching 0.7411 in the form of an exponential power function, resulting from Equation ( 5).However, the wind speed over land surfaces is difficult to obtain from remote sensing; therefore, field meteorological observations are conducted at the single-site level and spatial interpolations are performed for the limited number of ground stations at the spatial-distribution level.Therefore, wind speed is likely a limiting factor for the use of remote sensing applications in large regions.The findings in the present study indicate that a reasonable combination of the MODIS NIR BRDF_R index and NDVI can obtain a good fit and present an R 2 value of 0.721, which is close to the value obtained using wind speeds.The CI provides a significant new method of estimating temporal AR values over grassland via remote sensing.However, there is also some limitation for the CI, as it was developed under grassland during its growing season with NDVI greater than 0. As a result, it is probably not suitable for poorly vegetated areas.In the next step, we will attempt to develop another index to estimate AR for sparse or non-vegetated areas.

Discussion
The AR results in this paper were compared with the results over grassland from the literature.In the present study, the average AR value at the A'rou grassland station was 0.0101 m, and the AR values ranged of 0.0028 to 0.0234 m during the growing period.Stull (1999) provided a general reference AR value for prairie grassland of approximately 0.03 m [13].Ma et al. (2002) studied the AR over the grassy marshland surface of the Tibetan Plateau using in-situ turbulent data and obtained an AR value of 0.00466 m at Anduo for grass heights of 5 cm and a value of 0.0139 m at the HPAM for grass heights of 15 cm [14].Yang et al. (2003) determined the AR at Anduo PBL and Naqu, which are covered with sparse and short grasses, in the monsoon season based on wind and temperature profiles, and the order of magnitude of the AR was from 0.001 to 0.01 m [15].Thus, the AR results in the present paper are comparable to previously published in-situ results.
The relationship between the average wind speed at the height of 1 m and the in-situ AR on the eight-day scale (Figure 6) indicates that the average wind speed has a significant influence on the AR over grassland [7,22], with R 2 values reaching 0.7411 in the form of an exponential power function, resulting from Equation ( 5).However, the wind speed over land surfaces is difficult to obtain from remote sensing; therefore, field meteorological observations are conducted at the single-site level and spatial interpolations are performed for the limited number of ground stations at the spatialdistribution level.Therefore, wind speed is likely a limiting factor for the use of remote sensing applications in large regions.The findings in the present study indicate that a reasonable combination of the MODIS NIR BRDF_R index and NDVI can obtain a good fit and present an R 2 value of 0.721, which is close to the value obtained using wind speeds.The CI provides a significant new method of estimating temporal AR values over grassland via remote sensing.However, there is also some limitation for the CI, as it was developed under grassland during its growing season with NDVI greater than 0. As a result, it is probably not suitable for poorly vegetated areas.In the next step, we will attempt to develop another index to estimate AR for sparse or non-vegetated areas.

Conclusions
Two remote sensing-based indices, the MODIS NIR BRDF_R index and the NDVI, are sensitive to the AR over grassland surfaces (R 2 : 0.5228 for NIR BRDF_R; R 2 : 0.579 for NDVI), and the combination of these indices can significantly increase the R 2 to 0.721, which is consistent with the relationship between the wind speed and the field AR (R 2 : 0.7411).The proposed CI from remote sensing technology can accurately calculate AR values at a level comparable to considering the effects of wind speed, which is still difficult to acquire from remote sensing.The findings in this paper have powerful application potential for estimating temporal AR values over grassland surfaces during the growing season based on remote sensing data.

Conclusions
Two remote sensing-based indices, the MODIS NIR BRDF_R index and the NDVI, are sensitive to the AR over grassland surfaces (R 2 : 0.5228 for NIR BRDF_R; R 2 : 0.579 for NDVI), and the combination of these indices can significantly increase the R 2 to 0.721, which is consistent with the relationship between the wind speed and the field AR (R 2 : 0.7411).The proposed CI from remote sensing technology can accurately calculate AR values at a level comparable to considering the effects of wind speed, which is still difficult to acquire from remote sensing.The findings in this paper have powerful application potential for estimating temporal AR values over grassland surfaces during the growing season based on remote sensing data.
. The superstation is in A'rou Village of Qilian County, Qinghai Province (100.4643• E, 38.0473 • N, 3033 m), as shown in Figure 1.The experimental site features comprehensive, precise, and multi-scale observations of land surface and hydrological processes in a cold region.Intensive and long-term profile meteorological observations of the wind, air temperature, and humidity at the study site were conducted at six layers: 1 m, 2 m, 5 m, 10 m, 15 m, and 25 m.The specific manufacturer and model of each observation instrument are shown in

Figure 1 .
Figure 1.Landscape of the study site and observation equipment studied at A'rou station.Figure 1. Landscape of the study site and observation equipment studied at A'rou station.

Figure 1 .
Figure 1.Landscape of the study site and observation equipment studied at A'rou station.Figure 1. Landscape of the study site and observation equipment studied at A'rou station.
Atmosphere 2017, 8,16 6 of 11 temperature among different layers is relatively apparent, and the variation trends differ as the observation height increases.During the night, the air temperature increases as the height increases, whereas during the day, the air temperature decreases as the height increases.The air temperature itself shows a clear and stable curve with a single peak.

Figure 4 .
Figure 4. Relationship between the MODIS NIR BRDF_R index and the field AR (a) and between the MODIS NDVI and the field AR (b).

Figure 5 .
Figure 5. Relationship between the combination index (CI) and the field AR.

Figure 4 .
Figure 4. Relationship between the MODIS NIR BRDF_R index and the field AR (a) and between the MODIS NDVI and the field AR (b).

Figure 4 .
Figure 4. Relationship between the MODIS NIR BRDF_R index and the field AR (a) and between the MODIS NDVI and the field AR (b).

Figure 5 .
Figure 5. Relationship between the combination index (CI) and the field AR.

Figure 5 .
Figure 5. Relationship between the combination index (CI) and the field AR.

Figure 6 .
Figure 6.Relationship between the average wind speed and the AR at A'rou station.

Figure 6 .
Figure 6.Relationship between the average wind speed and the AR at A'rou station.