Domestic sewage dispersion scenarios as a subsidy to the design of urban sewage systems in the Lower Amazon River, Amapá, Brazil

Study area

The dynamics of the pollutant dispersion was assessed in a representative section of the North Channel of the Amazon River (NCM), between the Vila Nova River and the Igarapé (= creek) do Curiaú. The Santana Channel is located along this stretch of Santana, partially confined by Santana Island and the edge of Macapá, to the Igarapé do Curiaú (Fig. 1). The channel comprises a distance of 50 km in the coastal estuarine region of Amapá (Torres, El-Robrini & Costa, 2018).

In Santana, Companhia das Docas de Santana (CDSA) is equipped with a port area that exports commodities (manganese, wood, eucalyptus, and cellulose, among others) to and from other Amazonian regions. The navigation route is unique for vessels of up to 50 thousand tons through the NCM (Pereira et al., 2014; Cunha et al., 2021; Araújo et al., 2022) and is considered a sanitary risk area close to the urban areas of Macapá and Santana. For example, ships’ ballast water in this zone is considered an environmental and sanitary threat to biological resources because it is an import zone of ballast water (Pereira et al., 2014). This zone also includes an increased risk of oil and derivatives spills (Cunha et al., 2021; Araújo et al., 2022); therefore, it is an area of significant economic and environmental interest in water resource management, biodiversity conservation, and environmental sanitation.

Hydrodynamically, the study stretch is considered fluvial (≈250 km upstream of the North Atlantic Ocean), maintaining freshwater quality, but significantly influenced by variations in semidiurnal estuarine tides and with no detectable saline intrusion (Ward et al., 2013). CDSA tide gauges indicate that the amplitude or water level can vary up to ≈3 m (meso tides), favoring two daily reversals of the Amazon River currents (Less et al., 2021).

This coastal estuarine region can be subdivided into high, medium, and low stretches and is situated in Macapá Bay. The plain of this bay is often interrupted by tertiary formations (Barreiras Group sediments) similar to cliffs, with meanders, residual lakes, and “understories”. The term “understory is related to coastal wetlands, which consist of lagoons and lakes exclusively found in urban areas of Macapá and Santana, influenced or not by the tide, in addition to mangroves that make up the ecosystems of the Macapá bay (Torres, El-Robrini & Costa, 2018).

Climatic data

The climate of the region is classified as Am (monsoon equatorial), with elevated temperatures (never lower than 18 °C) and a pronounced drought period (August to November) (Köppen, 1936). The climatic variation occurs in the Amazon rainy season (December to August), with maximum rainfall between March and May, and less rainy Amazon season (September to November), with minimum rainfall in October (Limberger & Silva, 2016; Marengo & Espinoza, 2016; De Souza et al., 2017).

The pattern of seasonal variation in wind intensity follows a behavior opposite to that of rainfall (Figs. 2A and 2B), meaning it is less intense in the rainy season and more intense in the less rainy season (2010–2018 series). This series was obtained from the meteorological station of Macapá (MACAPA-AP-OMM: 82098)—Meteorological Database for Teaching and Research (Banco de Dados Meteorológicos para Ensino e Pesquisa—BDMEP) of the National Institute of Meteorology of Brazil (Instituto Nacional de Meteorologia do Brasil—INMET).

Experimental determination of liquid discharge (CSA and NCM) and bathymetry between the Igarapé do Curiaú and the Vila Nova River

The stretch that includes the CSA and NCM between the Igarapé do Curiaú (P1) and the Beija-Flor River (P11) was selected to evaluate the dynamics of natural coastal flow under the influence of large bodies of water and similar boundary and initial conditions in open channels.

Experimental campaigns of liquid discharge quantification were performed in the CSA and NCM sections to calibrate the hydrodynamic model. The procedure was performed with an acoustic method using an acoustic doppler current profiler (ADCP, Surveyor M9; SonTek). The ADCP was installed on the side of a boat measuring 22 m long, with the transducers immersed 1.0 m below the water surface. Compass calibration with GPS was performed at the beginning of each field session (M. de Abreu et al., 2020).

In the CSA section, the campaign was carried out between the neap and low spring tides on March 18, 2019, from 6:00 a.m. to 6:30 p.m., comprising a complete semidiurnal tide cycle, with 133 crossings in the channel section, whose average width is ≈650 m. In the NCM section, the campaign was carried out on November 7, 2017 (less rainy season) and March 19, 2019 (rainy season), in 12 crossings lasting approximately one hour each (≈12.5 km) (Fig. 1D).

The refined bathymetric survey of the margin, later incorporated into the hydrodynamic model, was carried out between P₁ and P₁₁ on November 8, 2017 (Fig. 1B). The procedure started at 9:00 a.m. in front of Igarapé do Curiaú and ended at 5:00 p.m. in front of the Vila Nova River, on a “zigzag” trajectory along the bank (Figs. 1C and 1D).

Modeling and hydrodynamic simulation and dispersion of sanitary effluents

SisBaHiA (Rosman, 2018) was used to generate the analysis of predicted tide and hydrodynamic simulation (May, August, and November) and effluent dispersion (May and November) between the Igarapé do Curiaú and Vila Nova River. This software has a 3D or 2DH hydrodynamic circulation model optimized for natural water bodies, to which Eurelian and/or Lagrangian models can also be integrated to describe pollutant transport phenomena (Rosman, 2018). In this study, we opted for the 2DH model (Eqs. (1)–(3)), which represents the hydrodynamic behavior. The Lagrangian model was considered the most suitable for this purpose, represented by the following equations:

(1)${\displaystyle \frac{\partial U}{\partial t}+U\frac{\partial U}{Ux}+V\frac{\partial U}{\partial y}=-g\frac{\partial \zeta }{\partial x}-\frac{gH}{2{\rho }_{°}}\frac{\partial \hat{\rho }}{\partial x}+\frac{1}{H{\rho }_{°}}\left(\frac{\partial \left(H{\hat{\tau }}_{xx}\right)}{\partial x}+\frac{\partial \left(H{\hat{\tau }}_{xy}\right)}{\partial y}\right)+\frac{1}{H{\rho }_{°}}\left({\tau }_x^S-{\tau }_x^B\right)-\frac{1}{H{\rho }_{°}}\left(\frac{\partial S_{xx}}{\partial x}+\frac{\partial S_{xy}}{\partial y}\right)+2\Phi senV-\frac{U}{H}\Sigma q}$ (2)${\displaystyle \frac{\partial V}{\partial t}+U\frac{\partial V}{Ux}+V\frac{\partial V}{\partial y}=-g\frac{\partial \zeta }{\partial y}-\frac{gH}{2{\rho }_{°}}\frac{\partial \hat{\rho }}{\partial x}+\frac{1}{H{\rho }_{°}}\left(\frac{\partial \left(H{\hat{\tau }}_{xy}\right)}{\partial x}+\frac{\partial \left(H{\hat{\tau }}_{xy}\right)}{\partial y}\right)+\frac{1}{H{\rho }_{°}}\left({\tau }_y^S-{\tau }_y^B\right)-\frac{1}{H{\rho }_{°}}\left(\frac{\partial S_{yy}}{\partial y}+\frac{\partial S_{xy}}{\partial x}\right)-2\Phi senU-\frac{V}{H}\Sigma q}$ (3)${\displaystyle \frac{\partial \zeta }{\partial t}+\frac{\partial UH}{Ux}+\frac{\partial VH}{\partial y}=\Sigma q.}$

In the two equations, U is defined as the velocity on the x-axis (m/s); V is the velocity on the y-axis (m s⁻¹);—ζ represents the free surface elevation (m); H is the depth of the water column (m); ρ_o is the average density of water (kg m⁻³); ρ_r is the density of the reference water (kg m⁻³); g is the acceleration of gravity (m s⁻²); τ_ij is the turbulent stress tensor; (i,j) represents the indices on the horizontal (x; y) plane; τ_i$\stackrel{ˆ}{}$S and τ_i$\stackrel{ˆ}{}$B are the wind stresses on the water surface, and the lower friction stress, respectively (kg m⁻¹s⁻²); 2 Φsen θU and 2 Φ θsenV represent the Coriolis accelerations; Φ is the angular rotational velocity of the Earth (rad s-1); θ is the latitude angle; ∑q represents the water balance that considers atmospheric rainfall, infiltration, and evaporation; S_ij represents the effect of the radiation stresses (Cunha et al., 2006; Rosman, 2018).

In the Lagrangian model of SisBaHiA, contaminant sources are represented by number of particles released in the source region at regular time intervals. The randomly arranged particles in the source region are then carried by the currents computed by the hydrodynamic (calibrated) model. The position of any particle at the next instant is defined by Pⁿ⁺¹, approximated by a second-order Taylor Series expansion. Pⁿ⁺¹ iscalculated from the previously known previous position (Pⁿ) (Eq. (4)): (4)${\displaystyle P^{n+1}=P^n+\Delta t\frac{d{P}^n}{dt}+\frac{\Delta t^2}{2!}\frac{d^2{P}^n}{d{t}^2}}$

For effluents discharged from any source, the amount of mass (Ma) of a given species (e.g., biological oxygen demand—BOD or CFU) is represented by a particle entering the modeled domain, which is given by Eq. (5) (Rosman, 2018): (5)${\displaystyle M_a=\frac{Q{C}_{a}x\Delta \tau }{N_P}.}$

Q represents the effluent discharge from the source, Ca is the concentration of the substance a present in the discharge from the source, and NP is the number of particles entering the domain per time step Δt. The dimensions of the source region are defined in such a way that the concentration in the source region is equal to that observed at the end of the initial dilution process, that is, within the near field of the mixture of the contaminant plume (Rosman, 2018).

First-order kinetics of the Lagrangian model

SisBaHiA considers that the total amount of effluent load, Q_T, discharged in a given source region per time interval Δτ, is defined by: (6)${\displaystyle QT=QQ/s\times \Delta \tau }$

If NP particles are launched in each Δτ, the initial amount of each particle, m_o, will be: (7)${\displaystyle m_o=\frac{QT}{N_P}=\frac{QQ/s\times \Delta \tau }{N_P}}$

Thus, it is assumed that over time, the remaining amount in each particle, m_(tv), is a function of its lifetime, “tv”. That is, it is possible to specify kinetic reactions, R_(tv), that alter the initial amount of each particle as follows: (8)${\displaystyle m\left(t_v\right)=m_{o}R\left(t_v\right)}$

First-order kinetic reactions (exponential decay) may be conveniently specified with the prescription of parameter T ₉₀. That is, the time required for the decay of 90% of the value of the initial concentration (m_o) or reduction of an order of magnitude. In this way, the kinetic reaction is given as: (9)${\displaystyle R\left(t_v\right)=exp\left(-k_d{t}_v\right)}$ Where the reaction constant K_d(decay) is calculated as a function of T ₉₀ as: (10)${\displaystyle K_d=-ln\left(\frac{0.1}{T_{90}}\right)}$

Initial contour conditions and the Lagrangian model

We can use two types of boundary conditions in SisBaHiA: 1) closed boundary nodes for the imposition of water discharge values and 2) open boundary nodes for free surface elevation values (Cunha et al., 2006) (Figs. 1B and 1D).

In the closed limit, a net discharge of the Amazon River was imposed (collected at the Óbidos station and experimental data in the NCM), and two of its main tributaries (the Xingu and Tapajós rivers) based on a measurement campaign carried out at the CSA and NCM, previously described or on experimental data available in the literature (Silva, 2009; Barros & Rosman, 2018), such as those of the National Water Agency (Agência Nacional de Águas (ANA), 2017) of the Óbidos station (code 00155001). In addition, rainfall levels, direction, and wind intensity were obtained from the Macapá weather station.

In the open limit, a series of water surface elevations was imposed to represent the tidal waves considering the six main tidal components (M2, S2, N2, K1, O1, and M4) imposed on nodes of the open limit mesh (Fig. 1B).

In the Lagrangian model, effluent discharge points (Fig. 1C), estimated raw sewage flow, and their respective concentrations were defined (Table 1). In addition, due to its environmental and sanitary importance to the coast, control points or imaginary reference lines were defined along the coast (RF at 300 m away from the riverbanks—RF-300; RF at 600 m away from the river banks RF-600), in addition to the point at the water collection station (CAESA) of Macapá. For example, RF-300 has been commonly used as a precautionary and safety line in outfall projects (RF-300 m) (Metcalf & Eddy, 1991; Rodrigues, 2012). Justification for this choice is the prior knowledge that displacement of pollutant plumes is concentrated in this zone and that there is a tendency to remain close to the riverbanks and for long distances greater than 50 km (Cunha et al., 2021; Araújo et al., 2022).

The choice of effluent launch points was based on an experimental campaign (November 2017), and Sewage Atlas reports available on the website of the National Water Resources Information System (Sistema Nacional de Informações sobre Recursos Hídricos—SNIRH) and the National Water Agency (Agência Nacional de Águas (ANA), 2017). These reports define specific raw sewage launch points (Figs. 1C and 1D) for what is known as the Metropolitan Zone of Macapá—MZM: Macapá (Igarapé do Curiaú), Santana (Vila Nova River), and Mazagão (Beija-Flor tidal channel) according to their population numbers.

A maximum concentration value of 10⁸ (CFU/100 mL of the sample) was used for the pollutant (thermotolerant CF or CFU) defined at the launch points imposed on the model. For CFU, the lifetime (T) or decay time of 90% of the particles (T₉₀) equal to T₁ of 21,600 s (6 h) was defined (Barros et al., 2015; Von Sperling, 2017). In addition, simulations were performed for T₉₀ equal to T₂ = 16,200s (4.5 h) and T₃ = 27,000s (7.5 h) to analyze sensitivity, self-depuration capacity, dynamics of pollutant dispersion in the analysis sections in NCM and CSA, and variations of the results in relation to decay.

Computational mesh

The computational domain covers an area of approximately 3 × 10⁷ km², containing 2,238 homogeneously distributed elements and totaling 10,627 nodes (Figs. 1B and 1D). The bathymetric data inserted in the mesh were obtained from the nautical charts of the Brazilian Navy [Directorate of Hydrography of the Navy (Diretoria de Hidrografia da Marinha—DHN); Hydrography Center of the Navy (Centro de Hidrografia da Marinha”—CHM)] (Ministério da Defesa do Brasil/Marinha do Brasil, 2019). They were then integrated with local experimental data and refined along the coast near Macapá and Santana, between the Igarapé do Curiaú and Vila Nova River (November 2017), generating approximately 19,000 additional points in the mesh.

Model calibration and verification in estimating pollutant dispersion

Data of the tidal rise (tide gauge station—CDSA) and the hydrodynamic data measured with those simulated by the model were compared to calibrate the hydrodynamic model. Flow and velocity measurements performed in the field (ebb and flood) were statistically compared with the model results, observing the acceptable error ranges of the parameters. When the errors exceeded acceptable statistical limits, adjustments to the parameters were made until their calibration (Hsu et al., 1999; Liu et al., 2008; Fossati & Piedra-Cuevai, 2013; M. de Abreu et al., 2020; Cunha et al., 2021; Araújo et al., 2022).

At the boundary of the open limit in the computational domain region, the harmonic constituents (M2, S2, N2, K1, O1, and M4) were considered and observed in the study region after the analysis of tidal elevations (Figs. 1B and 1D). The value of the constituents was imposed on the nodes of the open limit using adjustment of coefficients (in amplitude and phase) for each harmonic constituent (Liu et al., 2008; Nzualo, Gallo & Vinzon, 2018).

For model calibration and due to quality and consistency, the tidal rise data of 2016 from the CDSA were used as a reference. The objective was to reproduce the same harmonic components of the local tide with the smallest possible error and use the tide projection of SisBaHia. Thus, the tidal behavior predictions for 2019 were made. To estimate the reliability of the projections, comparisons were made between data projected by the model and those observed in 2017 and 2018.

Statistical analysis

The significance level evaluated between the experimental results and those of the predicted tides generated in SisBaHiA was tested, comparing the tidal levels recorded in March and November 2017 and May and August 2018. The modeled hydrodynamic results for May, August, and November 2019 were compared with the same corresponding predicted tides for 2019. Thus, flows between the tidal cycles (experimental and simulated) were statistically tested by successive comparisons with results simulated by the model for March 2019.

Quantification of the efficiency of the mathematical models was performed using the Nash-Sutcliffe (NSE) methods (Eq. (11)), Pearson’s correlation, R², and concordance index (d) (Eq. (12)). These statistical indicators were selected according to their different practical applications in hydrodynamic modeling, such as testing the accuracy and ability of models to represent physical reality (Devkota & Fang, 2015; Bakken, King & Alfredsen, 2016; Skhakhfa & Ouerdachi, 2016; Haddout & Maslouhi, 2018).

(11)${\displaystyle {\mathrm{E}}_{\mathrm{NS}}=1-\frac{\sum _{\mathrm{i}=1}^{\mathrm{n}}{\left({\mathrm{O}}_{\mathrm{i}}-{\mathrm{P}}_{\mathrm{i}}\right)}^2}{\sum _{\mathrm{i}=1}^{\mathrm{n}}{\left({\mathrm{O}}_{\mathrm{i}}-\overline{\mathrm{O}}\right)}^2}}$ (12)${\displaystyle d=1-\frac{\sum {\left|\mathrm{P}-\mathrm{O}\right|}^2}{\sum {\left(\left|\mathrm{P}-\overline{\mathrm{O}}\right|+\left|\mathrm{O}-\overline{\mathrm{O}}\right|\right)}^2}}$

where Oi and O represent the observed tidal rises and Pi and P represent the values predicted by the model. The results of the model efficiency can vary between 1 and −∞ in Eq. (11) (Nash-Sutcliffe efficiency), where 1 indicates the model that perfectly represents the observed data. Values greater than 0.75 indicate that the model is capable of representing what was experimentally observed (Al-Asadi & Duan, 2015).

Similar to the NSE efficiency test, the concordance index test (d) compares data generated by the model (P) with those observed (O): if d = 1, the model is considered perfectly coherent; if d = 0, the model shows a lack of ability to represent the reality of the observed data set. The minimum acceptable value for d would be 0.75.

Experimental results in the Santana Channel(CSA) and North Channel of the Amazon River—Macapá (NCM)

The results of experimental flow measurements performed in CSA and NCM (March 2019) and NCM (November 2017) were represented by the liquid discharge curves (Supplement-S1). In March 2019, the maximum ebb value in the CSA section was 22,729 m³s⁻¹, and the maximum flood value was −13,381 m³s⁻¹. In the rainier period (March 2019) in the NCM, values higher than those measured in November were obtained, with a maximum of 254,944 m³s⁻¹ in the ebb and −149,653 m³s⁻¹ in flood. In the less rainy season (November 2017), the values of 208,637 m³s⁻¹ (ebb) and −195,307 m³s⁻¹ (flood) were obtained, respectively. For both channels, the reversal period (in flood) was significantly shorter (between 3 h 20 min and 3 h 40 min) than the ebb period (between 8 h 20 min and 8 h 40 min).

The detailed velocity profiles of the currents measured with the ADCP, CSA (Supplement-S2a), and NCM (Supplement-S2b) are naturally interconnected (Fig. 1D). Therefore, they present interdependent hydrodynamic and bathymetric characteristics. The different colors of each vertical profile indicate vertical and transverse variations of the instantaneous velocities along the sections of each channel. When these profiles are integrated, it is possible to quantify maximum speed values in channels up to ≈2 m s⁻¹ (light and warmer colors).

The characteristics of the channels are very different, with average widths of ≈ 600–700 and depths of 20–35 m to CSA. The average width of the NCM is 12.3 km, and depth varies between 10–38 m. In the NCM section, three “sub-channels” normally influence the currents’ lateral, transverse, and longitudinal distribution, resulting in a high complexity of the 3D hydrodynamic flow profile.

Predicted observed tide and measured and simulated liquid discharge

The SisBaHiA predicted model provided the estimated tidal rise data (May to November 2019). To confirm the accuracy of the data predicted by the model, we performed statistical efficiency tests to compare the data observed at the CDSA station in 2017 and 2018 (M. de Abreu et al., 2020) (Table 2).

The statistical tests indicated reliable comparative results at the CDSA station and were used to reproduce numerical tidal elevation outputs without needing experimental data for March, May, August, and November 2019. Subsequently, these results were statistically compared with the increase in the tide in CDSA, which was used as a local tide gauge reference, and the levels generated by the hydrodynamic model of SisBaHiA.

The analysis indicated that the hydrodynamic model satisfactorily corresponds to the tidal elevations of the CDSA tide gauge station. Therefore, we correlated the efficiency of the simulated hydrodynamic model (variation of the experimental liquid discharge) of March 2019, both in the CSA and the NCM calibration section. The adjustment values of R² for the simulated and observed hydrodynamic behavior were greater than 0.95. The other statistical tests (Pearson, “d”, and NSE) for 03/18/2019 and 03/19/2019 indicated values higher than 0.9. The same occurred for values of experimental and simulated flood flow and velocity intervals with regard to CSA (03/18/2019) and NCM (03/19/2019) (Table 3).

Hydrodynamic simulation model

In simulations of the hydrodynamic model for the CSA and NCM sections in 2019, the moments of maximum flow (May) and maximum flood (November) were selected as references for simulating pollutant dispersion using the Lagrangian Model during the spring tide (Supplement-S3a–S3b).

The model-derived current values in NCM varied between 0.134 m/s and 3.34 m/s, with peak values concentrated in the center of the channel. In the CSA, speeds ranged from 0.133 m/s to 1.47 m/s, with higher speeds observed throughout the reference section (Supplement-S4).

Similar analyses were conducted for November 2019 during both spring and flood tide. The current velocity values in the reference section during the flood in November varied between 0.138 m/s and 3.44 m/s in NCM and 0.25 m/s and 1.36 m/s in CSA, respectively (Supplement-S5).

Dispersion and concentration of pollutants on the coast of Macapá and Santana

Spatiotemporal variations in the behavior of pollutant plumes showed a tendency of concentration and distribution along the Macapá and Santana shores. According to the Lagrangian model, after 1, 10, 20, and 30 days from the start of simulations, and considering a hypothetical decay rate of 21,600 s (May/maximum and November/minimum flows), the results of the dispersive processes of pollutants (total thermotolerant coliforms = CFU) are detailed in scalar fields (concentrations) (Figs. 3 and 4). Figures 3A–3B shows the fields of variation of the CFU concentration and geometric variation of the plumes for May (Fig. 3A) and November (Fig. 3B). The pollutant plume projected more toward the center of the channel in up to ≈2 km (NCM) in November compared to May.

On May 1st (Fig. 3C), the longitudinal extension of the plume was ≈63 km longer than in November. Despite being shorter longitudinally (≈52.2 km) (Fig. 3D), it tended to expand with more intensity laterally.

We found a greater expansion of the pollutant plume toward the center of the NCM ≈3 km wide between P1 (Igarapé do Curiaú) and P4 (downstream of the Igarapé das Pedrinhas) when compared to the 1st and 10th day after pollutant dispersion began (Figs. 4A and 4B). Spatially, the pollutant plume had higher concentrations near P5 (downstream of the Igarapé da Fortaleza) and between P2 (Channel of Jandiá) and P4 (upstream of the Igarapé das Pedrinhas), in addition to increasing its width by up to approximately 3.5 km in the transverse direction of the NCM.

Thirty days after pollutant dispersion simulations started, the increase in the CFU concentration level was also observed in May (Fig. 4C), and the maximum peak concentration value reached between all simulations was 1.9 × 10⁶CFU/100 mL. In the selected period, NCM and CSA were in the flood tide phase, and the pollutant plume had a greater extension area at P4 (Igarapé das Pedrinhas). Thirty days after the release of pollutants in November (Fig. 4D), the maximum concentration peak reached was 1.5 × 10⁶. In November, this same order of magnitude was observed in all simulations.

Regarding P10 (Vila Nova River), the concentration values varied with the tide period, between 200 and 1,000 CFU/100 mL. However, at P11, concentrations remained above 2,000 CFU/100 mL between the Vila Nova River and the launch point without directly contributing to the pollutant plume released into the Amazon River.

Table 4 describes the variations in effluent concentration for the chosen points or along the coastal reference line (P-300, P-600, and CAESA), considering 1, 10, 20, and 30 days after the release of the pollutants. The criteria T₉₀ = 21,600s, 27,000s, and 16,200s were adopted, respectively.

At T₁, there were significant seasonal differences between the average concentration permanence values at the reference geographical imaginary points or lines for the two periods (May and November). The average decay time of CFU concentration in May presented 40% of the average value of November at the CAESA point (hotspot), 54% in relation to the RF-300, and 53% in relation to the RF-600. Spatially for T₁, the concentration variations between the RF-300 and RF-600 references decreased significantly concerning the distance from the margins. The average of the RF-600 values represented only 28% of the average of the RF-300 values (Table 4).

For T₂, the differences between average concentrations in RF-300 and RF-600 behaved differently than in the T₁ simulations. For example, the CFU concentration averages obtained for RF-300 and RF-600 were higher in May (40%) than in November (18%) (Table 4).

At T₂ (27,00s), the concentration of CFU/100 mL in May had the highest averages compared to the T₃ and T₁ decay. In November, only the average concentration at the CAESA point was higher in relation to T₃ and T₁. Only the average concentration in the RF-300 reference in May, along with the T₂ decay time, was higher than in November when compared with T₁ and T₃ (Table 4).

Possibly due to a shorter decay time (in theory, with higher clearance rates), the average concentrations of T₃ (16,200 s) showed lower values than those observed for T₁ and T₂. As expected, they had similar seasonal and spatial behavior. In May, the CAESA reference presented only 40% of the concentration in November, and the RF-300 reference had 55% in relation to November. The behavior was different regarding the RF-600 reference, which displayed a concentration in May that was 22% higher than in November (Table 4).

High values of overall average concentrations in the respective references were found for T₂, with 6,319 CFU/100 mL. Still, similarly, the overall average concentration value was 4,720 CFU/100 mL for T₁ and 3,813 CFU/100 mL for T₃. As expected, the RF-600 reference had the lowest overall average pollutant concentration (Table 4). The RF-300 reference had the highest concentration for T₁ and T₃. Therefore, only for T₂ the overall average pollutant concentration was higher, at the CAESA point (Table 4).