1.

Introduction

Rice is one of the main cereal crops in the world that plays an important role in food security, especially in the context of climate change, environmental pollution, and population growth. It is a daily staple food for a major part of the world’s population, especially in the developing countries across Asia, Latin America, and Africa.¹ This has led to an increasing demand for stable rice supply, especially for major rice-producing countries such as China, Thailand, and Vietnam. This issue requires a system that can provide important information for monitoring and managing rice cultivation such as rice age (the number of days after planting), date of planting (dop), and date of harvest (doh). The information on the rice growth stages is also important for implementing management practices including irrigation systems, pest warning, fertilizer application, procurement and export, and other purposes. Furthermore, climate change, resulting in an increase in extreme weather events such as drought and floods,² can damage rice crops and affect food security and farmer livelihoods at local and regional scales. It is, therefore, necessary to accurately and timely determine how much of the rice area is in any stage at any time or position on the map.

Many published studies have used optical imagery such as Landsat³ to map rice paddy fields at a local scale with either unsupervised or supervised classification methods. High-temporal resolution optical sensors such as AVHRR,⁴ SPOT VEGETATION,^5,6 and MODIS^7–9 were also used to map large-scale rice fields. However, these optical images often have low spatial resolution ($>250\mathrm{m}$) when applied at a regional scale and are frequently affected by cloud cover. Although there are techniques^10–12 to reduce the effects of cloud, but in tropical monsoon regions like Vietnam, quasipermanent cloud cover limits the use of optical data during the rainy season. In contrast to optical sensors, radar sensors with their polarizing characteristics, can capture differences in the vegetation growth stages through their effect on the canopy structure, such as leaf area, plant height, and plant biomass.^13,14 In addition, high-temporal resolution synthetic aperture radar (SAR) systems, such as ALOS-2 (14 days), TerraSAR-X (11 days), and Sentinel-1 (12 or 6 days), can provide more reliable data for monitoring rice growth information.

Radar remote sensing is recognized as an effective tool for mapping rice growing areas.^15–17 The most frequently used SAR data for monitoring common rice parameters are C-band¹⁸ and X-band.¹⁹ The C-band SAR data are particularly suitable for rice monitoring and mapping,^20–22 because of the sensitivity of the C-band backscatter to plant vertical structure and to the inundation status of the underlying soil.²³ Until recently, SAR data used for rice mapping and monitoring were from RADARSAT-1/2,²⁴ ENVISAT ASAR²⁵ (advanced synthetic aperture radar), ERS-1/2²⁶ (European remote sensing), and Sentinel-1.^27,28 However, previous SAR-based methods were not used for large-scale rice mapping due to limited data availability and high cost.^15,16,29 Only, with the launch of Sentinel-1A satellite in 2014,²² users can access free SAR data with spatial resolution of 20 m [interferometric wide swath (WS) mode], 12-day repeat cycle, and WS of 250 km. By 2016, the Sentinel-1B satellite²² was launched to allow images to be recorded with a 6-day cycle when combined with Sentinel-1A. The dense time series of Sentinel-1 data and free access at high spatial resolution provide an opportunity to monitor near-real-time rice growth in the study area. Sentinel-1 image data have been widely used in many applications such as crop classification and area detection,²⁷ height and biomass estimation of rice plants.²⁸ Therefore, rice monitoring over large areas using SAR images with high spatial resolution and near real-time repeat cycles became widely applicable.

In this paper, we develop a method to determine growth stages using Sentinel-1 data over the Mekong Delta. In the region, conventional approaches of determining growth stages are often dominated by expensive and time-consuming surveys and estimation methods, which are based on a limited number of samples at commune and district levels with hidden errors. As stated above, SAR data have the potential to detect rice growth stage.^30,31 This is based on the temporal variation of the radar backscattering of paddy fields, which can provide information on farming practices such as field preparation, planting (sowing and transplanting) and on crop growth from germination to maturity stage.

Using SAR images often pay no attention to the effect of the incidence angle, however, it usually affects the backscatter value.^32,33 Recently, a number of studies have begun to focus on reducing the influence of the incidence angle on SAR images, involving the work of Pathe et al.,³⁴ who used a linear model to normalize the incidence angle suitable for the whole ENVISAT ASAR scene to estimate the soil moisture. Topouzelis and Singha³³ used this method for SAR WS mode data for applications in the oceanography. Nguyen et al.³⁵ also applied this method to standardize the incidence angle for ENVISAT ASAR WS images to monitor rice growing areas and rice crops. In this study, normalization of the incidence angle for Sentinel-1 IW data with three subswaths was performed using the quadratic function. The Mekong Delta is a large area covering almost a Sentinel-1 scene with a swath width of 250 km, which is affected by the incidence angle. Thus it is necessary to consider the influence of incidence angle on the scattering pattern of rice plants to reduce its influence on the area considered in this study.

Information on rice growth, such as rice age, date of planting, and date of harvest, is important to support irrigation and crop management, monitoring, and planning, including pest warning and fertilizer practices according to each growth stage of rice to ensure rice yield. Therefore, many studies were conducted to quantify these parameters using time-series SAR data. An example of determining the number of days after sowing of short-cycle crops is the work of Yang et al.,³⁶ who analyzed the scattering behavior of oilseed rape using C-band RADARSAT-2 data. Asilo et al.³⁷ used the TerraSAR-X data series to determine flooding/transplanting date of rice plants with good results. Similarly, Hoa et al.³⁸ estimated sowing date of paddy fields using COSMO-SkyMed data with root-mean-square error (RMSE) = 4.3 days. However, most of these studies focused on small areas. For a large region such as the Mekong delta, there are several issues to be addressed before to realize operational applications. The first issue is related to the change in the radar backscatter with radar incidence angle for large swath. The second is related to the diversity of rice cultivation across a large region. To this aim, in this study, an algorithm is developed to estimate growth parameters such as rice age, and planting and harvest dates, using Sentinel-1 IW data for a large region, the Vietnamese Mekong Delta.

In particular, the objectives of this study were (1) to process a time-series image data for reducing noise effect and reducing the effect of incidence angles on a large region; (2) to build algorithms to estimate rice age, date of planting, and date of harvest using time-series data processed for the whole Vietnamese Mekong Delta; and (3) to validate the results of the rice age map in the Mekong Delta region with field data.

The study is conducted in relations with initiatives for establishing communities for support and cooperation in crop monitoring, including rice monitoring, such as Group on Earth Observations Global Agricultural Monitoring Initiative and the Asian Rice Crop Estimation and Monitoring,³⁹ where SAR data plays a major role to facilitate rice monitoring.

2.

Materials and Methods

The study used satellite data acquired from the Sentinel-1 satellite. The collected Sentinel-1 image data covers almost the entire Mekong Delta with swath width of 250 km (Fig. 1). The Sentinel-1 A and B 6-day temporal resolution products were used to create time series of backscattering values with a spatial resolution of 20 m. The products were ground range data (GRD) calculated from digital number (DN), converted to sigma-naught value, and then terrain correction was carried out to remove the effect of the terrain. Ground data were used, including field data for result validation (Fig. 2). Other reference datasets include district and provincial rice area statistics in the Mekong Delta.

Fig. 1

The Mekong Delta overlaid by Sentinel-1 scenes.

Fig. 2

Study design flowchart for extracting rice age information from Sentinel-1 images.

2.1.

Data Description

2.1.1.

Image data

Sentinel-1 is a constellation of two satellites, Sentinel-1A and Sentinel-1B. The data products have a 12-day repeat cycle and spatial resolution of 20 m. Sentinel-1A satellite was launched on April 3, 2014, whereas Sentinel-1B satellite was put into orbit 2 years later, on April 25, 2016.⁴⁰ Since 2016, the image acquisition cycle in the Mekong Delta is 6 days with a spatial resolution of 20 m (interferometric WS mode) and dual polarization data (VV—for vertical transmit and vertical receive, VH—for vertical transmit and horizontal receive). For this research, in order to monitor the rice growth status of the 2017 autumn–winter crop, data from August 03, 2017, to February 17, 2018, were collected with a total of over 30 images captured by both Sentinel-1A and B satellites with the descending mode. All the image data was co-registered based on a reference image so that they have the same coordinate system.⁴¹ SAR data were calibrated to convert DN value to the backscatter coefficient.

2.1.2.

Ground data collection

Depending on the purpose and requirements of the field data, the field measurements were designed at different levels of detail. To build the dataset for this research, the sample fields were designed so that they are located far from each other and spatially distributed throughout the study area. Field data collected include rice variety, maturity stage, harvest date of the previous crop, sowing/transplanting date of the current crop (farmer survey), other types of vegetation, and photographs.

The field data were collected in two groups: (1) the training data comprised measurement data of rice parameters for 2017 autumn–winter crop collected continuously in An Giang province with a 12-day repeat cycle coinciding with image acquisition date of Sentinel-1; the rice parameters comprised of rice variety, sowing/transplanting method, sowing/transplanting and harvest date, growth stages and (2) the validation data by survey for the whole delta from December 26 to 31, 2017, serving the main purpose of evaluating the estimation method of rice age, sowing/transplanting and harvest date.

Field measurements of rice parameters were conducted in 60 sample paddy fields (Fig. 3), with areas ranging from 3000 to $\mathrm{77,000}{\mathrm{m}}^2$. These measurements were conducted from September 9 to December 26, 2017. A total of 10 measurements were performed every 12 days corresponding to image acquisition dates of Sentinel-1A satellite to cover changes in paddy fields during the crop season from tillage to harvest. The selected sample fields have different rice varieties including short- and long-day rice varieties of the Mekong Delta, ranging from 80 to 120 days and three major rice varieties including Nep, OM5451, and DS1, accounting for 90% of sample fields (Table 1). These sample fields were selected to represent the study area and also for their accessibility. Figure 3 shows the locations of sample fields across An Giang province, Mekong Delta.

Fig. 3

Study area and locations of sample fields across An Giang province, Mekong Delta.

Table 1

Information on rice varieties in 60 sample fields for 2017 autumn–winter crop in An Giang province.

No.	Rice variety	Growth cycle (day)	Sample number	Percent (%)
1	Nep	95 to 100	24	40
2	OM5451	90 to 95	22	37
3	DS1	105 to 110	8	13
4	Other varieties	85 to 100	6	10
	Total		60	100

One of the problems in this study pertains to the difference in the backscatter time series between the two cultivation methods, sowing and transplanting. To solve this problem, the backscatter time series of VH polarization of the sample fields were compared the two sample fields with similar sowing/transplanting dates. The sample CT01 with long-day rice variety DS1 having a growth cycle of 105 to 110 days was sown on August 28, 2017, whereas the sample TS24 with short-day rice variety OM5451 having a growth cycle of 90 to 95 days was sown on August 31, 2017. Figure 4 shows a pattern of backscattering values of VH polarization for the sample TS24 from sowing to harvest and the field photographs corresponding to the time series of image acquisition with a 12-day cycle.

Fig. 4

Backscattering pattern of VH polarization for the sample TS24 with sowing method and field photographs collected corresponding to the time of image acquisition. OM5451 sown on August 31, 2017, has a growth cycle of 90 to 95 days.

Figure 5 shows a pattern of backscattering values of VH polarization for the sample CT01 from sowing to harvest and the field photographs corresponding to the time series of image acquisition with a 12-day cycle. Similarly, Fig. 6 shows the fluctuation of backscattering values of the sample CT07 with the same rice variety as the sample CT01 but different in the cultivation method, i.e., the transplanting method. The difference between the sowing and transplanting dates of these two sample fields is no more than two days. The results showed that the scattering patterns of both sowing and transplanting methods recorded an increase in the backscattering values after sowing/transplanting from day 10 onward for the VH polarization. This makes it possible to apply the same method of determining sowing/transplanting dates for both cultivation methods. In this study, the date of sowing/transplanting will be referred to as dop.

Fig. 5

Backscattering pattern of VH polarization for the sample CT01 with sowing method and field photographs collected corresponding to the time of image acquisition for every 12 days. DS1 sown on August 28, 2017, has a growth cycle of 105 to 110 days.

Fig. 6

Backscattering pattern of VH polarization for the sample CT07 with transplanting method and field photographs collected corresponding to the time of image acquisition for every 12 days. DS1 transplanted on August 26, 2017, has a growth cycle of 105 to 110 days.

Field data for the verification of the growth stages—rice age, dop, and doh were collected from December 26 to 31, 2017, across different provinces in the Mekong Delta. In particular, data were collected at 241 sample points of rice fields spatially distributed over the Mekong Delta (Fig. 7). The main objective of this data collection was to evaluate the classification results of the spatial distribution of rice area and rice age in the 2017 autumn–winter crop. In addition, this dataset was also used to investigate the effect of incidence angle on the scattering pattern of the rice field.

Fig. 7

Locations of field sample points over the Mekong Delta surveyed from December 26 to 31, 2017.

2.2.

Time-Series Data Preparation

2.2.1.

Image preprocessing

Sentinel-1 image data were calibrated first and then corrected for terrain with 30-m Shuttle Radar Topography Mission digital elevation model data using Sentinel Application Platform software. Afterward, the radar data were filtered by the multitemporal speckle filter⁴² and the spatial filter to reduce noise. The temporal filter method can improve the Sentinel-1 SAR data time series by reducing noise while preserving as much as possible the delicate structures that existed in the radar image.⁴¹ The pixel size of the image data after preprocessing is 20 m.

2.2.2.

Incidence angle normalization

Since Sentinel-1 image data cover a large area of 250 km, the incidence angle is a factor that affects the radar scattering results. Therefore, incidence angle normalization is an essential requirement for large-scale monitoring applications of rice (i.e., regional level). Sentinel-1A and 1B image data collected in the study area with IW swath mode consist of three subswaths: IW1, IW2, and IW3.⁴³ There is sufficient overlap between the subswaths to ensure continuous ground coverage as provided in GRD products.

The backscatter values are not only affected by the land cover but also the incidence angle. In order to detect changes in the backscattering due to changes in surface state, it is necessary to eliminate the effect of the incidence angles by normalizing the incidence angle effect. This approach has been presented in this study of estimating soil moisture from ERS⁴⁴ and ENVISAT ASAR data,³⁴ and mapping crops from ASAR WS data.³⁵ The effect of the incidence angle on the backscatter values of the VV and VH polarizations need to be evaluated for different land cover types, and in particular for rice, the angular effect can change with the plant growth stages. This is because of the different angular variations of the backscatter of areas with surface or volume scattering. However, for Sentinel-1 over the Mekong delta with the incidence angle range of 31 deg to 46 deg, the angular behavior of the volume scattering, which can be modeled as a cosine function of the incidence angle, will have a small change, whereas the surface scattering, depending on the surface roughness, can have a large variation across the range. For rice fields which correspond to surface scattering at the beginning of the season and volume scattering at the peak growth season (and in between, surface-volume interaction), a simple approach is to neglect the variation due to volume scattering, in this case, surface scattering variation can be normalized. The water surface was selected for surface scattering because it is the most easily extracted homogeneous object in the study area where the two main branches of the Mekong River, namely Tien and Hau Rivers, flow across the study area. Figure 8 shows the spatial distribution of the incidence angle across the Mekong Delta.

Fig. 8

Spatial distribution of incidence angles in Sentinel-1 image over the Mekong Delta.

The results of the relationship between the scattering values of the water surface of the VH polarization are shown in Fig. 9(a). It indicates that the change of the incidence angle has an effect on the scattering values in 3 subswaths. This will affect the result of the developed algorithm to extract information of the planting date or the rice age in large areas with a large change in the incidence angle. Water surface value changes continuously from $-21$ to $-34\mathrm{dB}$ depending on the change of the incidence angle. Therefore, it is necessary to normalize the incidence angle to eliminate its effect on the Sentinel-1 images.

Fig. 9

The backscatter values of the (a) VH and (c) VV polarizations of the water surface affected by the incidence angle, and the backscatter values of the (b) VH and (d) VV polarizations of the water surface after incidence angle normalization.

To normalize the influence of the incidence angle, this study proposed an empirical formula presented as Eq. (1). The formula normalizes the backscatter coefficient at any incidence angle to the backscatter coefficient at an incidence angle of 41 deg. The backscatter value was normalized to an incidence angle of 41 deg because this is the incidence angle at the 60 sample fields in An Giang province, and to facilitate comparison with the measured field data. The empirical equation is expressed as

Eq. (1)

$${\sigma }_{\mathrm{nor}}={\sigma }_0-{\sigma }_{\phi }+{\sigma }_{41},$$

where ${\sigma }_{\mathrm{nor}}$ is the normalized backscatter coefficient, ${\sigma }_0$ is the unnormalized backscatter coefficient, ${\sigma }_{\phi }=y(\phi )$ is the backscatter coefficient for modeling the incidence angle, and ${\sigma }_{41}$ is the backscatter coefficient at an incidence of 41 deg.

The equation for modeling the incidence angle of the VH polarization can be divided into three parts corresponding to different incidence angles of different subswaths. The function $y^{\mathrm{VH}(\phi )}$ of the polarization is expressed in detail in the following equation:

Eq. (2)

$$y^{\mathrm{VH}}(\phi )=\{\begin{array}{l}0.1527{\phi }^2-10.464\phi +154.53;41.5\le \phi <46.0\\ 0.4274{\phi }^2-33.654\phi +635.11;36.5\le \phi <41.5\\ 1.1847{\phi }^2-103.58\phi +2233.9;31.0\le \phi <36.5\end{array},$$

Eq. (3)

$$y^{\mathrm{VV}}(\phi )=-0.559\phi +1.9368.$$

The results of incidence angle normalization of the VH polarization are presented in Fig. 9(b). The findings indicate that the scattering values of the water surface after normalization are almost stable with the change of the incidence angle. Figure 9(c) shows that the relationship between the scattering coefficient of the VV polarization and the incidence angle is almost linear. Thus the equation for modeling the incidence angle of the VV polarization will be simpler with different incidence angles. The function $y^{\mathrm{VV}(\phi )}$ of the VV polarization is presented in detail in Eq. (3). The results of incidence angle normalization of the VV polarization are illustrated in Fig. 9(d). It shows that the backscatter values of the water surface after normalization have a stable range of values with respect to the incidence angle.

2.2.3.

Fast Fourier transform filter

Although the use of multitemporal data effectively reduces the amount of speckle noise, the time series of scattering values still show significant differences related to morphology or land cover changes. Interpolation of data gaps and smoothing to reduce noise are often encountered when processing time-series data.

The phenomenon of noise due to sudden anomalous changes to each pixel is detected by performing an initial harmonic analysis based on a period of $\sim 5$ to 6 months, which is the minimum time required to achieve satisfactory results for a rice crop. The harmonic analysis is based on fast Fourier transformations which have a strong impact on noise phenomena in a time series.⁴⁵ The Fourier transform allows decomposing a time-dependent periodic signal into the frequency domain, through which frequency information is expressed as the constituent sine and cosine functions.⁴⁶ Thereby, the Fourier transform decomposes different components of the frequency of the time series into different harmonic waves, and each harmonic wave is characterized by a unique amplitude and phase for the data series.⁴⁷ According to Moody and Johnson,⁴⁸ for plant morphology, the main signal concentrated on low-order harmonic waves and high-frequency signals are usually noise. In the case of rice fields, for the VH backscatter time series, the main signals from orders of 1 to 4 for a 5- to 6-month period were investigated, while the anomalous noise was the components in higher order. The analysis results for the VH data series for the main subjects on the image over a 5-month period showed that this filtering method strongly affects water surface, perennial trees, urban areas, and aquaculture lands, which have relatively stable scattering values. This helps to eliminate sudden anomalous variations (Fig. 10). Figure 10 also shows that the scattering pattern of rice fields is almost not affected by this filter method.

Fig. 10

VH backscatter time series of land covers before and after FFT filtering, in which the data series has been processed with multitemporal speckle filter and incidence angle normalization.

2.2.4.

Smoothing filter

Fast Fourier transform (FFT) filtering has no significant effect on the backscatter curve of the paddy fields. Therefore, a more efficient filtering method for multitemporal scattering of rice fields should be considered. This can be performed by applying a Savitzky–Golay (SG) smoothing filter⁴⁹ programmed in interactive data language (IDL). The SG algorithm can generally be described as a moving window filter that uses linear least squares regression to fit a data series.⁵⁰ The implementation of SG filter has the ability to smooth the data series and increase the signal-to-noise ratio without distorting the signal. However, the decisive factor lies in the choice of the correct window size, which is especially important for large areas with very irregular soil surface morphology over time, as rice fields in this case. Large window sizes lead to smoother curves, though its ability for monitoring the short-term changes is more limited. It is obvious that the rapid change in morphology, as the much changes in scattering values at growth stages, is not fully preserved with SG filter if a large window size is applied. In contrast, a small window size can overcome this problem, often resulting in local noise reduction of rice fields and creating a smooth time series that mimics the typical backscattering pattern of a rice crop. As a result, a SG filter with a small window size (three observations) was applied to the pixels of rice fields in the VH polarization data series of 5 to 6 months (Fig. 11).

Fig. 11

VH backscatter time series for rice fields in the following cases: multitemporal speckle filter, incidence angle normalization, FFT, and smoothing filter.

2.3.

Extraction of Phenological Parameters

The expected rice parameters extracted from the Sentinel-1 image data include planting date (sowing/transplanting), harvest date, and the number of days after planting (rice age). Biological data are extracted for each growing season based on a time series that has been interpolated and smoothed. This information will be extracted from the data series through conditional functions programmed in IDL. Thresholds for sowing/transplanting (Thr_dop) and harvest (Thr_doh) dates were determined based on statistical values from field-surveyed samples. The beginning and end of the growth cycle have been determined by special points in the multitemporal data series.

Accurate determination of rice age in agriculture is a result of changes in cultivation practices and the growth stages. This study investigated the backscattering model of three different rice varieties with different growth cycles including: OM5451 (90 to 95 days), Nep (95 to 100 days), and DS1 (Taiwan rice variety, 105 to 110 days). The changes in the scattering coefficients of VH, VV polarizations, and VH/VV ratios of OM5451, Nep and DS1 are described in Figs. 12Fig. 13–14, respectively.

Fig. 12

Multitemporal radar scattering coefficient of OM5451 for (a) VH, (b) VV polarizations, and (c) VH/VV ratio of Sentinel-1 images for 2017 autumn–winter crop. This rice variety has a growth cycle of 90 to 95 days.

Fig. 13

Multitemporal radar scattering coefficient of Nep for (a) VH, (b) VV polarizations, and (c) VH/VV ratio of Sentinel-1 images for 2017 autumn–winter crop. This rice variety has a growth cycle of 95 to 100 days.

Fig. 14

Multitemporal radar scattering coefficient of DS1 for (a) VH, (b) VV polarizations, and (c) VH/VV ratio of Sentinel-1 images of 2017 autumn–winter crop. This rice variety has a growth cycle of 105 to 110 days.

Figures 12–14 reveal that sowing/transplanting and harvest dates can be determined using VH polarization. The scattering value of VH polarization increases with the increase in the number of days after sowing/transplanting, especially after about 20 days for all long-cycle and short-cycle rice varieties and this was the basis for the method of determining sowing/transplanting dates. To determine the sowing/transplanting date, the study used VH-polarized data series, where ${\sigma }_{(i)}^{\mathrm{VH}}$ is the backscatter value with $i=1,\dots ,n$, and $n$ is the number of observations in the data series. Let $\text{day}(i)$ be the image acquisition date of the $i$’th image in the data series, as coded in Algorithm 1.

Algorithm 1

Determining the date of planting (dop).

for $i=1$, $n-4$ do

if {$[{\sigma }_{(i+4)}^{\mathrm{VH}}-{\sigma }_{(i)}^{\mathrm{VH}}]$ >= Thr _dop} and [${\sigma }_{(i)}^{\mathrm{VH}}$ < −18 dB] and [${\sigma }_{(i+4)}^{\mathrm{VH}}$ > −21 dB] then begin

$\text{dop}=\text{day}(i)-3$

i_dop = i

endif

endfor

All short-to-long rice varieties show that the common harvest date is always equal to or longer than the growth cycle. This is explained by the fact that rice fields are usually planted in each subregion with the same irrigation or dyke system and with nearly the same sowing/transplanting date. This facilitates soil preparation, irrigation, and harvest in each area. Therefore, in each subregion with the same irrigation system, the harvest dates will be almost the same and so some early ripening fields are ready to be harvested but not yet harvested. To determine the harvest date, the study also used the VH-polarized data series as shown in Algorithm 2.

Algorithm 2

Determining the date of harvest (dop) and days after planting (dap).

if {[day(n) – dop >= 70] and [day(n − 2) − day(i_dop + 10) > 0]} then begin

for i = i_dop +10, n − 2 do

if {[${\sigma }_{(i)}^{\mathrm{VH}}$ − ${\sigma }_{(i+2)}^{\mathrm{VH}}$] >= Thr_doh} and [${\sigma }_{(i)}^{\mathrm{VH}}$ > −18 dB)] then begin

doh = day(i)

i_doh = i

endif

endfor

endif

Days after planting (dap) is calculated by the formula as follows:

dap = day(i) – dop

3.

Results

3.1.

Effect of Incidence Angle on Backscatter Pattern of Rice Fields

The application of the normalized method for backscatter measurements ensures the equivalence of the scattering signal of similar object types located at different incidence angles in the images. Figures 15(a)–15(f) illustrate the scatter patterns of different incidence angles for the two cases of original and normalized values of the incidence angle of backscattering for rice, perennial plants, and water surface. The comparison between the normalized scattering values and the incidence angle of 41 deg in the sample points surveyed in An Giang province [Fig. 15(d)] shows that there is not much difference between the scattering data series before and after normalization. In contrast, Figs. 15(a), 15(b), and 15(e) show a significant difference between the scattering data series before and after incidence angle normalization.

Fig. 15

Multitemporal backscattering coefficients of rice, water surface, and perennial trees for VH polarization of Sentinel-1 images in 2017 autumn–winter crop before and after incidence angle normalization. Backscattering coefficients at different incidence angles include (a) 32 deg, (b) 35 deg, (c) 37 deg, (d) 40 deg, (e) 43 deg, and (f) lines showing the backscatter coefficients of the water surface at different incidence angles.

3.2.

Spatio-Temporal Growth Stages of Rice Paddy Cultivation in the Mekong Delta

Most rice areas are planted in the Vietnamese Mekong Delta with two or three crops per year, namely winter–spring, summer–autumn, and autumn–winter crops. The rice age map of 2017 autumn–winter crop estimated from the Sentinel-1 image data series is shown in Fig. 16. Throughout the study area, the map shows that the rice age of different growing regions in the Mekong Delta is different, corresponding to the different growth stages or sowing/transplanting dates. Rice age was classified into 12 information layers with a 10-day interval after sowing/transplanting, corresponding to the longest rice variety planted in this area (nearly 120 days), although common rice varieties in the Mekong Delta are 85 to 105 days. There were some areas where the late autumn–winter crop has not been harvested yet on December 31, 2017, such as those in Bac Lieu and Soc Trang provinces. At the same time, there were some areas where the 2017 autumn–winter crop has been harvested but have not been sown/transplanted for new crops as in some parts of An Giang, Long An, and Tra Vinh. On the southern coast of the Mekong Delta, rice fields were sparse, especially visible in the southern parts of Ca Mau, Bac Lieu, and Soc Trang provinces.

Fig. 16

Rice age map on December 31, 2017, in the Mekong Delta derived from Sentinel-1 multitemporal data from August 03 to December 31, 2017.

3.3.

Validation of Rice Age Map

A part of the estimation results for rice age was the determination of sowing/transplanting and harvest dates for rice paddy fields. Figure 17 presents the results of the comparison between sowing/transplanting dates estimated from VH-polarized series and field-surveyed data. The time to determine the rice age from the image data is on December 31, 2017, and the time to collect field information on sowing/transplanting and harvesting dates is from December 26 to 31, 2017. The comparison results show the ability of Sentinel 1 data to accurately estimate sowing/transplanting and harvest dates of rice fields as indicated by high coefficients of determination, respectively, $R^2=0.917$ and $R^2=0.697$ (Fig. 17).

Fig. 17

Comparison of sowing/transplanting and harvest dates estimated from Sentinel-1 data on December 31, 2017, and field-surveyed data from December 26 to 31, 2017, in the Mekong Delta.

To assess the ability to estimate the rice age close to real time, the dataset from August 03 to December 31, 2017 was used to produce the rice age map on December 31, 2017. Figure 18 describes the differences in the stages of rice age, depicting the difference between the estimated days after sowing/transplanting and field data collected on December 26 to 30, 2017. It is easy to see the potential application of the multitemporal Sentinel-1 satellite data for near-real-time rice age monitoring. The differences between the actual collected dates and estimated dates correspond to $\mathrm{RMSE}=7.4$ days, $\mathrm{B}\mathrm{i}\mathrm{a}\mathrm{s}=-0.7$, $\mathrm{MAE}=6.2$ days, and the coefficient of determination between these two quantities is $R^2=0.918$ ($n=241$). However, it is also noted that the ability to accurately identify rice age ranging from 0 to 20 days is low. This is due to the limitation of the method using backscatter when determining low biomass rice fields that are $<20$ days old after planting. The result of the rice age map compared with field photos is presented in Fig. 19.

Fig. 18

Comparison of rice age estimated from Sentinel-1 data on December 31, 2017, and field-surveyed data from December 26 to 31, 2017, in the Mekong Delta.

Fig. 19

The result of the rice age map and field photos taken from December 26 to 31, 2017, in the Mekong Delta.

4.

Discussion

Remote sensing data at a single time often does not meet the needs of agricultural monitoring. Therefore, the use of SAR data series is becoming more and more popular in multitemporal analysis methods. Sentinel-1 SAR data, which are freely available, becomes a reliable source of building a multitemporal dataset with the frequency required for agricultural applications. The Sentinel-1 with IW swath mode allows images to be captured over a large swath width of 250 km and high spatial resolution of 20 m. However, to derive applications for a large region such as the Mekong delta requires data preprocessing steps such as incidence angle normalization in addition to reduction of speckle noise.

In this study focused on the rice backscatter time series, the incidence angle was normalized by applying the nonlinear quadratic function to the VH polarization for the water surface to match the backscattering of the surface scattering in rice fields, which has a much more important angular variation than that of the volume scattering (Fig. 9). The proposed method resulted in the similar backscatter temporal variation of the rice fields across incidence angles after normalization as shown in Fig. 15. After normalization, the effect of rice growth stages is a major contributor to backscatter changes in a time series.

With regards to SAR data processing, it is essential to implement multitemporal data processing methods to improve data quality and significantly reduce speckle noise. The time series data acquired from Sentinel-1A and B will play a vital role in providing high-quality data in terms of temporal and spatial resolution. This will enable the use of pixel-based approaches, which are simpler and less time-consuming than other approaches.^52,53 Our processing algorithm has created a time series of backscatter data where the noise has been filtered by FFT (Fig. 10) and smoothing (Fig. 11) techniques. Thus this allows evaluating the changes in scattering values of the predominant rice growing areas, determining the parameters of rice age, as well as the planting and harvest dates. If an area with minor changes in actual sowing/transplanting and harvest dates, time-smoothing methods can allow better identification, classification, and detection of rice areas without losing too much time information.

Analysis of multitemporal high-resolution SAR data, especially using the backscatter coefficient (${\sigma }^0$) from C-band, is highly suited for detecting flooded rice paddies.¹³ This is particularly advantageous in the tropical monsoon regions in Asia where most of the world’s rice is produced in large areas and in cloudy conditions.^17,25 Today, we know that the applicability of SAR data is of great potential including the use of the Sentinel-1 data archive, which is free, with extensive coverage, and high temporal resolution. These advantages allow for the creation of consistent and high-resolution rice age maps with nationwide coverage. High temporal resolution (6 days) enables near real-time monitoring of agricultural information such as rice age with an appropriate interval. For example, the plan in the Sentinel-1A and B mission will allow creating monthly or weekly datasets of rice age nationwide. In addition, the proposed rice age monitoring algorithm allows the use of high-resolution SAR time series to provide information on the growth stages of rice, especially such as the time of sowing/transplanting and harvest in the growing areas. This will facilitate and provide important information to the agricultural management organizations.

The results of this research are very encouraging in comparison to other studies^37,38 in terms of quality and consistency despite differences in data used and study area. This study proposed an algorithm for estimating planting date (sowing/transplanting) for Sentinel-1 C-band IW images for a large area with better results, i.e., with $R^2=0.917$ and $\mathrm{RMSE}=6.2$ days ($n=153$). This proves that the determination algorithm of sowing/transplanting dates using multitemporal data series can be applied in practice with RMSE nearly equivalent to the temporal resolution of Sentinel-1 images.

The accuracy evaluation of harvest date estimated by the algorithm highly correlates with the field surveyed date with $R^2=0.697$ and $\mathrm{RMSE}=5.7$ days ($n=88$). Several studies^13,54 show that the C-band SAR backscattering reaches its maximum at the ripening stage, which is attributed by angular scattering from upright structures of plants and water surface. Torbick et al.⁵⁵ used Sentinel-1 time-series analysis to determine the area of rice harvested in Myanmar and compared it with the area survey statistics, obtaining an $R^2=0.78$. This study shows the potential for determining harvest dates with time-series analysis of rice growth stages. This study considers a fuller evaluation, in which our proposed method seems to be more promising for practical applications. At the same time, this study also showed that the VH polarized scattering values decreased after the rice fields were harvested as the basis for determining the harvest date.

To evaluate the rice age map on December 31, 2017, using Sentinel-1 image dataset for the last 5 months from August 3 to December 31, 2017, this study has used the field dataset corresponding to the time of creating the rice age map from December 26 to 31, 2019. The results obtained have a strong correlation between the estimated rice age and the actual rice age with an $R^2=0.918$ and $\mathrm{RMSE}=7.4$ ($n=241$). These outcomes show that the image dataset used over the last 5 months was able to create a rice age map with appropriate accuracy. In addition, the study also found that if the rice fields had days $<20$ days after planting, it was unlikely to be identified on radar images. This can be easily seen in Figs. 12–14, where the backscattering values of rice $<20$ days old are very low corresponding to the scattering of wet soil or water surface, because rice plants at this stage are still young with low biomass. At the same time, the research results show that the sowing and transplanting methods have an effect on the scattering pattern of rice plants (Figs. 4–6), but it does not greatly affect the algorithm for determining sowing/transplanting date and on the evaluation results with field data.

However, because the algorithm for determining the rice age is based on an empirical model, it is developed and verified with specific field data of rice cultivation technique in the Mekong Delta. Hence, to apply in other areas, the algorithm needs to be verified before adoption. Even though, the general research methods developed here can be applied to all regions of the world. The study also showed similar results to previous studies which suggested that C-band SAR data could track the changes in the growth stages of agricultural crops.^18,56 The potential for practical application of SAR data to determine sowing/transplanting and harvest dates, rice age, and growth was demonstrated in this study. The differences in sowing/transplanting dates across the study area led to differences in crop growth among rice fields in each region. This is due to the influence of irrigation conditions, dike systems, and local environmental conditions and hence having different cropping systems. Therefore, the proposed rice age monitoring method has taken the advantages of SAR data to detect and map rice growth.

5.

Conclusions

The feasibility of using C-band SAR data to monitor rice age, planting and harvest dates of the rice fields was presented. C-band Sentinel-1A and B image-series data with IW swath mode acquisition was collected in the Mekong Delta of Vietnam from August 03, 2017, to February 17, 2018. The data were preprocessed and filtered carefully with FFT and SG filters. The image-series data were then normalized to the same incidence angle for the entire image to reduce the effect of the incidence angle on the radar backscatter value. The study developed an algorithm to determine rice age, planting and harvest dates suited to the processed data series. This study has significantly contributed to the understanding of the response of polarized SAR data to the development of rice plants according to the growth cycle. The results show the high sensitivity of VH polarization when determining rice growth status parameters.

The results show the applicability of the rice age determination algorithm, which can be used to monitor sowing/transplanting and harvesting dates during the growth stages with data series that have been preprocessed. Furthermore, the results of the analysis show that changes in backscatter values are more sensitive to changes in soil conditions and growth stages, which is important when applying algorithms to data time series. Therefore, the information about rice age corresponding to the growth stages of the rice field is important for agricultural management such as irrigation management, fertilizer application or pest prevention according to the growth stages of rice fields. In addition, the results are also able to provide important information for the procurement of agricultural products in each region, and to minimize the negative effects of the market on food security.

The results also emphasize the importance of preprocessing, noise filtering of time-series data, and incidence angle normalization based on SAR data for crop monitoring on large areas where the effect of incidence angles is significant. The Mekong Delta has very different sowing/transplanting dates, depending on the local conditions of the region and the rice crop, which leads to difficulties in statistical and agricultural management according to traditional methods. The results not only help determine the rice age in a timely manner with near real-time data, but it can also assist in improving management plans in agricultural farming, as well as allow farmers/managers to develop optimal strategies for managing agricultural land to ensure productivity.

The results of this study are based on the analysis of C-band SAR data-series and field-surveyed data to propose appropriate rice age monitoring algorithms. As the approach used in this study is based on pixels, the mixed pixels of rice fields may cause confusion or estimation errors—a potential constraint of the method. However, the Sentinel-1 image data have a pixel size of 10 m, which greatly reduces this effect. From the application perspective, it is necessary to evaluate the method more comprehensively for other rice seasons in a year or in another areas. In future studies, it is necessary to investigate the application of this method to an area with different farming conditions and the adaptation of this method to other agricultural crops.