www.mdpi.com/journal/remotesensing Article Estimating Flow Resistance of Wetlands Using SAR Images and Interaction Models

The inability to monitor wetland drag coefficients at a regional scale is rooted in the difficulty to determine vegetation structure from remote sensing data. Based on the fact that the backscattering coefficient is sensitive to marsh vegetation structure, this paper presents a methodology to estimate the drag coefficient from a combination of SAR images, interaction models and ancillary data. We use as test case a severe fire event occurred in the Parana River Delta (Argentina) at the beginning of 2008, when 10% of the herbaceous vegetation was burned up. A map of the reduction of the wetland drag coefficient is presented.


Introduction
Wetlands, particularly those in floodplains and in coastal areas, contribute to flood control by storing and decreasing the velocity of excess water during heavy rainfall.Wetland vegetation also provides a natural barrier to fast moving water and therefore aids in flood speed reduction.The result of wetland activity during floods is often to decrease damage to surrounding areas (i.e., nearby cities).
Significant advances have been made to gain a better understanding of flow phenomena in floodplains and natural wetlands.However, the majority of these studies are based on laboratory experiments with simple artificial vegetation roughness, whereas in reality natural vegetation exhibits a wide variety of forms and spatial distributions and also a considerable variability in time and space [1,2].
The estimation of the vegetation drag coefficient constitutes a measure of flow resistance and is a critical variable for the development of reliable hydrodynamic models of wetlands [1].However, to date there is no operational technique to estimate the spatial distribution of the vegetation drag coefficient in natural wetlands.The fact that this variable can be considered a complex function of both vegetation geometry and density [1], provides arguments to explain why there are difficulties in obtaining spatially distributed information about this variable.These difficulties are related to the unavailability of spatial information on vegetation structural conditions.
In previous works [3,4] it was shown that significant information about marsh condition can be achieved using active microwave remote sensing.In particular it was found that for some marsh species, the backscattering coefficient σ 0 measured by a Synthetic Aperture Radar (SAR) is sensitive to plant structure and condition.Moreover, it was shown that the σ 0 can be considered a function of both marsh condition and hydrological regime [2,[5][6][7].The purpose of this paper is to investigate changes in flow resistance in the floodplain of the Paraná River Delta Region in Argentina, caused by changes in vegetation density due to severe burnings occurred in herbaceous vegetation wetlands.The methodology proposed here estimates the vegetation drag coefficient using as input: (1) Envisat ASAR images; (2) SACC/MMRS optical images; (3) fieldwork and (4) a microwave interaction model.The result obtained is a map of the vegetation drag coefficient of the wetland marshes after the fire event.

Study Area
The Paraná River Delta Region in Argentina is a large (17,000 km 2 ) wetland located in the last portion of the Paraná floodplain and near major cities such as Buenos Aires.One of the most common disturbances of this ecosystem is the damage or destruction of herbaceous vegetation by intentional or natural fires.This occurs because marshes are both the most extended autochthonous vegetation and the least valued by society [8].The most common species is the junco (California bulrush, Shoenoplectus Californicus), which covers up to 15% of the wetland area (~2,900 km 2 ).Eighteen field sites inside the Junco marsh area have been regularly monitored since 2005 [9], in order to obtain multitemporal biogeophysical data of the marsh status.Standing biomass, water level, junco marsh plant density, junco height, junco radius and junco dry matter density are systematically measured [10].Two junco marsh soil status are considered, normal condition (nonflooded) when the soil has a permanent layer of water of about 5 cm, and flooded condition, when water level rises and juncos can be partially or completely flooded.

Fire Event
At the beginning of 2008, a series of fires burned large areas of the herbaceous vegetation of the Paraná River Delta Region.The consequences of these disturbances in the future behaviour of the wetland are still under discussion [11], but the magnitude and the duration of this event indicates that it was a major disturbance for this ecosystem.
From the hydrodynamic point of view, the most important expected consequence is the modification of the wetland flow regime in periods of high water levels of the Paraná River.This is due to the fact that this wetland works as a huge reservoir capable of accepting the excess of water from this river and reducing its average momentum [12].As a whole, it has the overall effect of increasing the time the water remains as fresh water, before being released into the sea.This has two very important and positive effects: since effective discharge and water resilient time increase, flood water level peak decreases and total sedimentation increases [13].
In order to characterize the future hydrodynamic behaviour of the wetland after the fire event, there are two important magnitudes to be estimated: marsh water storage capacity and marsh hydraulic conductivity.The determination of marsh water storage capacity consists on the estimation of the volume of water inside the wetland as a function of time, and can be obtained using microwave remote sensing data [3,14].This magnitude was previously estimated in the Paraná River wetland for flooded and nonflooded conditions using field data and ENVISAT ASAR images [4,15].Characterization of the hydraulic conductivity is a much more difficult task, which requires the estimation of the drag coefficient of vegetated areas of different places of the watershed [2].

Satellite Data Sets and Fieldwork
The methodology described in the next sections is based on the exploitation of optical and radar data.Data from the optical Argentine system SACC/MMRS was used to derive a land cover map of the area and a map of the burnt area [9,11].Data from ENVISAT ASAR was used for the estimation of marsh plant density in burnt areas.Main characteristics of these two systems are summarized in Table 1.

Methodology
The main objective of this work is to provide an estimate of the spatial distribution of the vegetation drag coefficient of the Paraná River wetland marshes after the fire event of the summer of 2008.To this end, we need to: (1) Determine marsh burned area; (2) Determine vegetation condition prior and after fire event; (3) Derive information about the drag coefficient, using the vegetation condition as input.
To determine marsh burned area (1), our methodology was based on the supervised classification of SAC-C MMRS images.This system presents a good trade off between resolution and swath coverage and was successfully used in previous works to produce landcover maps of the Paraná River Delta and to produce the most recent map of the area [9].This landcover map gives reliable information about Parana's River Delta status prior to the fire events and will be used as a base map to determine the junco marsh basal area (Figure 1).To determine the burned area, we used a more recent map derived from a new SACC MMRS image collected after the fire events ended.The classification scheme was based on a progressive segmentation of Band 4 of MMRS (795-835 nm) [11].This new map provides the necessary information about marsh condition after the fire events.
As stated before, marsh σ 0 is sensitive to vegetation structure and hydrological condition; therefore, to monitor marsh condition after the fire event (2), we used ENVISAT ASAR images acquired in Wide Swath (WS) mode at HH polarization [16].This operation mode was chosen because it covers the entire area with one image and adequate resolution (75 m).In addition, the HH polarization maximizes σ 0 dynamic range and sensitivity to vegetation condition [3,6].ASAR images corresponding to a time prior to the burning events and after them are shown in Figure 2.
At this point, two different problems arise: the estimation of vegetation structure from σ 0 and the estimation of the drag coefficient using this information.

Estimation of the Reduction of Junco Plant Density
All our test sites were totally burned during the fire events.From photographs obtained during the field campaigns (Figure 3), it can be seen that all the standing plants were removed by the fire, only leaving a residual value of 15 plants/m 2 .
Figure 4 shows the HH σ 0 of junco marsh patches extracted from ASAR images collected prior and after the burning events.For the second date, σ 0 's of both burned and unburned areas are reported.For all the cases, the samples were extracted from eighteen patches of 2 m × 2 m randomly selected from six large homogeneous junco marsh areas (~0.25 Km 2 , see Figure 1).More details of the fieldwork procedure and the data statistical analysis are given in [10] and [15].From the figure, it can be seen that the values of the backscattering coefficient change ~6 dB from normal condition (green) to burned condition (red).This increase in the observed σ 0 values with increasing vegetation density is not unexpected, and it was also observed for natural vegetation in [17] and for crops in [18].Since the images ENL is ~21, this change cannot be attributed to speckle fluctuations, and should be related to changes in the phenological/hydrological condition of the marsh.Since the hydrological condition was similar for the two dates, these changes should be linked with the decrease of junco biomass.This fact is supported by hydrological data (Figure 5) and by the comparison of HH σ 0 box plots of not-burned junco marshes for the two dates (Figure 4).Furthermore, the observed σ 0 values of both nonburned junco marsh conditions (prior and after the fire events) corresponds to the observed σ 0 values of normal (nonflooded) Junco marsh condition previously reported in [6].This result is also expected, since both image acquisitions correspond to normal junco marsh condition (not flooded and not burned).
Figure 5. Water level at Zárate city, used to infer water level inside field sites (see [6]).The dates corresponding to ASAR acquisitions are indicated.This observed decreased in HH σ 0 , though not related to a change in hydrological condition, is expected, since the marsh burned condition corresponds to a mean value of JPD ~15 plants/m 2 and the nonburned condition to a mean value of ~85 plants/m 2 (Figure 4).Therefore, there is evidence of a correlation between junco marsh HH σ 0 and JPD.This correlation between σ 0 and JPD is expected from theoretical studies [19] and was also found in [10], where the regrowth of junco marshes using data from ERS-2 SAR was monitored.
Since the entire set of field sites was destroyed by fire, only extreme values of JPD junco marshes patches were observed.But there is strong evidence derived from optical images and field work that there are other areas in the wetland where the burning was not total, i.e., 15 plants/m 2 < JPD < 85 plants/m 2 .To build a continuous map, it is necessary to estimate JPD from HH σ 0 for all the intermediate values of JPD, for which there are not available ground truth.Instead of relying on a simplified fit (linear, polynomial) between the mean values of both extreme observation conditions, an interaction model was used to simulate the variation of junco marsh σ 0 values as a function of JPD.This model is capable of simulating the interaction of electromagnetic waves with vegetation and to estimate the junco marsh σ 0 values, provided their geometric and dielectric properties are given as inputs [6].This model has the advantage of being already tested for junco marshes in different (and extreme) environmental conditions [4,6,15].Furthermore, this approach has the advantage of being based on physical hypothesis about the vegetation structure and the interaction between vegetation and electromagnetic radiation.

The interaction model
The scattering processes have been simulated using the electromagnetic model developed at Tor Vergata University [16] and adapted to junco marshes in [6].The model the soil as a homogeneous half space.The vegetation is described as an ensemble of discrete dielectric elements.Important details are given in [6], while the main points are summarized below.
For the soil, the bistatic scattering coefficient is computed by means of the Integral Equation Model [20].In our case, the backscattering of the soil is expected to be very low, since flooded wetland soil roughness is very low [6].Junco plants are described as dielectric cylinders.The "infinite length" approximation [21] is used to compute the bistatic scattering cross section and the extinction cross section.Radius, height and orientation distributions are given on the basis of available ground truth.Junco plant permittivity is related to junco plant moisture by the semiempirical formula given in [22].The scattering contributions due to the soil and the single vegetation elements are then combined to include multiple scattering effects of any order [19].In combining soil and vegetation effects, both specular soil-vegetation reflection and diffuse soil-vegetation multiple scattering processes are considered.In our case, the first contribution is dominant.The main input values to the model are summarized in Table 2.
This microwave interaction model was successfully validated in previous works and it was used as a key component of a water level retrieval scheme in junco marshes [6,10,15].Using this model, and the inputs of Table 2, it is possible to simulate σ 0 values of the junco marsh for intermediate JPD values.Figure 6 shows the result of this simulation.2.
Several remarks about the simulation results in Figure 6 are listed below: 1.The interaction model generally reproduces the observed σ 0 values of junco marsh corresponding to extreme values of JPD (box plots red and green).2. The interaction model predicts a rapid increase of the total σ 0 as function of JPD, related to an increase in junco marsh radar cross section reaching a maximum at about 40 plants/m 2 .Then, the HH σ 0 trend becomes flat and shows almost no sensitivity to JPD changes.This can be explained considering that the junco marsh total σ 0 can be considered a balance between soil-junco double bounce interaction and Junco extinction [6].An increase in JPD means additional junco plants available for backscattering, but also more junco plants for extinction.The complex balance between these two magnitudes depends strongly on vegetation structure and dielectric properties (see [6] for a complete discussion about junco marsh scattering behaviour).3. Finally, this simulation can be used to estimate remaining JPD as a function of HH σ 0 , at least for low values of JPD.

Estimation of the Junco Marsh Hydraulic Conductivity
Now that we have developed a procedure to estimate JPD as a function of σ 0 , the following step is to relate JPD to the vegetation drag coefficient.Reference [23], states that there are two components of the terrain that are potentially able to reduce the momentum of the water: soil and vegetation.However, there is experimental evidence that in systems dominated by dense vegetation, the drag force exerted by the vegetation is dominant [23].Consequently, the reduction of vegetation density should have a very strong impact on wetland hydrodynamic behaviour.The decrease of flow velocity can be related to the loss of fluid linear momentum caused by the force exerted by the vegetation using the horizontally averaged momentum equation [24]: where ∂H/∂x is the mean slope of the channel, g is the acceleration of gravity, (uw) is the averaged Reynolds stress, z the vertical coordinate, C d is the vegetation drag coefficient, u 0 is the magnitude of the velocity field of the fluid and A the reference area of the body projected on a plane normal to u, the velocity field.The reference area A depends explicitly on vegetation geometry.Structurally, junco shoots are long (mean height ~180 cm), cylindrical plants (mean diameter ~0.7 cm).They are located inside islands surrounded by levees, where water velocity is typically small (~0.1 cm/s).Therefore, junco shoots will be modelled as vertical cylinders, and the main task will be to estimate the drag coefficient of a circular cylinder immersed in a flow of water with an average velocity u 0 relatively small.
The problem of estimating the vegetation hydraulic resistance given the vegetation characteristics is not new, and it is still open.This is mainly related to the failure of simplistic theoretical formulations and the intrinsic difficulty related to experimental measurements [2].The most important parameter dictating the hydrodynamic behaviour of a system is the Reynolds number, a dimensionless number defined as the ratio of inertial forces to viscous forces used to quantify the relative importance of these two for a given problem.Mathematically: where a is the cylinder radius and υ the water dynamical viscosity.In this case, R ≈ 70; therefore, we are working in the low Reynolds number regime, where flux is usually laminar.Many experiments have been carried out for different types of vegetation in different Reynolds number regimes in order to measure the drag coefficient of leafless bushes (i.e., junco marshes) [25,26].Reference [25] presented a semiempirical formula for the computation of the drag of a plant inside a vegetation patch, based on experimental studies with cylindrical elements.He showed that the drag of a cylinder inside a cylinder array for the low Reynolds number regime can be estimated as: where C d∞ is the drag coefficient of a single cylinder in an ideal 2D flow, d is the diameter of the cylinder and a x and a y are the longitudinal and lateral distances between the cylinders.Reference [23] notes that the two terms of (3) represent the effect of cylinder blockage and the effect of free surface, respectively, the two most important contributions to vegetation drag in low Reynolds number condition for vertical vegetation.It is important to note that: 1. a x and a y depend implicitly on vegetation spatial distribution and 2. the only dynamical-dependent input of Equation ( 3) is C d∞ , the drag coefficient of a single cylinder.
At this point, two separate problems arise: the estimation of C d∞ , and the estimation of a x and a y .Since junco shoot height is large compared with its radius, they are usually considered as infinite cylinders and the key parameter to estimate is the force per unit length of cylinder.The following simplifying assumptions are taken: 1. Junco shoots can be considered as rigid circular cylinders.Actually, junco shoots are flexible, but its flexibility is only important at large water velocities [23].2. The flux is stationary (∂u/∂t = 0).This means that the input flux to the channel does not vary with time.3. The Reynolds number of the problem is low.These are all reasonable assumptions, which were easily checked using historical data from our field work and are usually adopted by other authors for similar environmental conditions [2].If these assumptions apply, it can be shown that the drag coefficient of a circular cylinder is on the range 0.8-1.5 [21,24].Following [2], we will assume C d∞ = 1.
The other parameters required to estimate the overall junco patch drag coefficient are a x and a y , the longitudinal and lateral distances between junco plants.These values were measured for all the field sites and their dependence as a function of JPD is shown in Figure 7.This plot summarizes the measurements of mean a x , a y for 18 different junco patches representing a range of JPD from 5-110 plants/m -2 .

Figure 7.
Plot of the mean values of a x (black) and a y (red) as a function of JPD.A linear fit of mean a x , a y as a function of JPD is also presented (R 2 (a x ) = 0.94, R 2 (a y ) = 0.92).
Since there is an evident (and expected) linear dependency between a x , a y and JPD, a linear fit (Equation ( 4)) is proposed: It can be seen that there is almost no difference between a x and a y values, and this is consistent with a random junco plant spatial distribution, a reasonable hypothesis in natural vegetation.Therefore, both variables will be considered as identical and will be estimated as a function of JPD using (4).

Estimation of Junco Marsh Patch Drag Coefficient
Model results analysis (Figure 6) shows that due to the structural/dielectric nature of junco marsh ecosystem, the sensitivity of C band HH σ 0 with JPD is only important at low JPD values, where increasing biomass has as a net effect an increase in the overall HH σ 0 .After a value of near 40 plants/m 2 (corresponding to a σ 0 value ~-6 dB), the increase of JPD has a negligible effect on junco marsh σ 0 .This implies that using this frequency and polarization, we can only expect to discriminate between different low JPD conditions.But this information is important, since these are the most disturbed patches and the most important ones for the hydrodynamic modelling of the watershed.
Given that the model simulations of Figure 6 predict that the HH σ 0 trend as function of JPD is roughly linear from ~15 plants/m 2 (~-11 dB) to ~40 plants/m 2 (~-6 dB), a simple segmentation of the image was carried out in order to establish 7 classes of JPD values.The segmentation step was chosen to be of the order of the expected image radiometric accuracy.Therefore, the segmentation step was set to 1 dB.Using this segmentation, we can relate σ 0 values to JPD values.Also using Equation (3) and Equation ( 4), we can also relate JPD values to junco plant drag coefficient values, and ultimately to junco marsh drag coefficient values, since for low Reynolds number the drag of the whole patch can be well approximated by the sum of the corrected drags of the single plants [2].The proposed segmentation of the backscattering coefficient is shown in Table 3.

Results
In order to derive a map of the vegetation drag coefficient from SAR data, the methodology proposed in (4.3) was applied to the ASAR image of 2008.09.01 (corresponding to the burned condition) in the junco marsh burned area.Figure 8 shows a map that represents the junco marsh drag coefficient of the Paraná River Delta by 2008.09.01 derived using the methodology developed in the previous sections.
It is important to address the overall consistency of the map.There are three types of errors to be considered: (i) errors due to speckle fluctuations; (ii) errors in the interaction model that estimates JPD from σ 0 and (iii) errors in the hydrodynamic model that estimates C d for a junco patch given its JPD and ancillary data.The speckle fluctuation (i) should be ~1 dB with a confidence interval of 85%, due to the ASAR image acquisition mode.The interaction model was intensively validated for junco marshes in different and extreme environmental conditions, and an overall error (ii) of 1 dB was estimated from the comparison with the SAR image and field data (see [6] for details).Finally, the errors in the hydrodynamic modelling depend strongly on the validity of the hypothesis presented on section 4.2.As stated before, these are reasonable assumptions usually adopted by other authors in wetland marshes, and most of them were checked experimentally in junco marshes [10].

Discussion
Several conclusions can be extracted from the results presented in Figure 8: • The junco drag coefficient map presents areas with low and high spatial homogeneity.This is typical of burning events in wetlands; where the access to fuel and the presence of water usually leads to a heterogeneous burning pattern.This heterogeneity is not related to the speckle phenomena, since (i) ENL is high and (ii) the speckle phenomenon is expected to affect the whole image.
• Areas where islands are larger (and cattle raising is more common), were the most affected both in intensity and homogeneity (C d ≈ 25 or less, red in the map).
• Areas in the middle part of the Paraná River Delta were less homogeneously affected, due to the intrinsic heterogeneity of the vegetation and the water channels present in the area.
• Values presented in blue (class coastal) correspond to coastal junco marsh patches located at the island line coast that present a higher backscattering coefficient due to a different flood state.Therefore, they will be excluded from our analysis.
• Table 3 shows that most of the burned area corresponds to completely destroyed junco marsh patches (<15 plants/m 2 ), as expected from the information obtained using optical images and field data.This implies that fire events reduced drastically the wetland overall drag coefficient in this area.This should have important effects in the future behaviour of Paraná's River seasonal flooding, reducing water resilient time and increasing water level peak.

Conclusions
In this work we presented a novel methodology to estimate changes in flow resistance in wetlands marshes using a remote sensing approach.The developed methodology was applied for estimation of the decrease of the vegetation drag coefficient in the Paraná River wetland, analyzing a local and recent burning event.We obtained estimations of burned junco marshes area in the wetland, using optical images and fieldwork.Then, we estimated junco marsh plant density using SAR images and a microwave interaction model to make a physically-based fit of the observed data.Finally, we estimated junco marsh hydraulic conductivity using the semiempirical formula introduced in [25] and in situ experimental data obtained from our field sites in the area.The estimated values of drag coefficient are comparable to those obtained in the literature for similar vegetation characteristics in comparable conditions [2,23,27].
Due to the fact that the backscattering coefficient measured at C Band and HH polarization is only sensitive to changes in JPD for low values of JPD (see Figure 6), the methodology described in this paper is best suited to discriminate among areas of low JPD conditions.Nevertheless, the information related to these heavily disturbed areas (low JPD is related to low drag coefficient) is usually the most important one for the differential hydrodynamic analysis of the ecosystem.
It is important to point out that the map was generated using a physically-based fit of the available observations, that limits this map to only a few different classes of drag coefficient values.This limitation is related to the observing system technology (frequency, polarization, incidence angle) and the interaction mechanism between the radar waves and the vegetation.An alternative approach based on a linear segmentation of the ASAR image would imply to assign the high σ 0 values (~-6 dB) to medium and high values of JPD (~80 JPD).This would imply an overestimation of JPD values (and of the vegetation drag coefficient).
In this work, only σ 0 at HH polarization was modelled, because to cover with a single image the whole Paraná River Delta, a Wide Swath image with only one polarization was used.Nevertheless, it is important to mention that, with more advanced systems that are able to measure the full polarimetric signature of vegetation, a more complex and robust estimation can be performed [28].
There are other less extended ecosystems in the Paraná River Delta that suffer frequent burning events and are important for the hydrodynamic modelling of the watershed.Unfortunately, there are no interaction models available for the vegetation present in these ecosystems; therefore no reliable relation between σ 0 and vegetation density can be assessed.This problem can be solved by the development and validation of interaction models for other types of herbaceous vegetation.This requires expertise in electromagnetic modelling and access to ad hoc fieldwork [19].
Summarizing, a quantitative methodology to derive the reduction on marsh drag coefficient in natural watersheds was introduced.This methodology is based on remote sensing data synergy (optical, radar) and the use of interaction models.This constitutes a necessary step for interpretation and prediction of the complex behaviour of the backscattering coefficient of marshes.

Figure 1 .
Figure 1.Landcover map of the Paraná River Delta at the beginning of 2008 [9].Field sites areas are marked in red.

Figure 2 .
Figure 2. Envisat ASAR Wide Swath image of the Paraná River Delta prior to the fire events (2008.02.04, left) and after them (2008.09.01, right).Burned junco marsh areas are shown in red, while selected not burned areas are shown in green.

Figure 3 .
Figure 3. Photographs of one of the test sites, prior to the burning events (2008.02.04, left) and after them (2008.09.01, right).A strong decrease in the emerged biomass was observed.Also, it can be seen that the soil condition is similar in both dates.

Figure 4 .
Figure 4. Box plots of C HH σ 0 values of junco marsh sites in normal condition (green), and for the burned condition (red).The associated mean junco plant density (JPD, plants/m 2 ) measured in the field sites is also informed.

Figure 6 .
Figure 6.Box plots of C HH σ 0 values of junco marsh sites in normal condition (green) and burned condition (red) superimposed with a model simulation of the HH σ 0 values of the junco marsh as a function of JPD.The input model parameters used in this simulation summarized in Table2.

Figure 8 .
Figure 8. Map of junco marsh drag coefficient derived from the ASAR image of 2008.09.01 using the methodology described in (4.3).Only burned junco marsh area is shown.The segmentation criterion is summarized in Table 3. Classes' colours and names are informed in the figure legend.C d range values are informed in the figure legend (Coastal, <25, 25-34, 34-43, 43-54, 54-65, 65-76, and >76).

Table 1 .
Characteristics of satellite data sets.

Table 3 .
Segmentation proposed to derive JPD from σ 0 values and the corresponding single plant and marsh drag coefficient.