The relationship between the physicochemical properties and permeability of the Fluvisols and Eutric Cambisols in the Zagreb aquifer, Croatia

The physicochemical properties of soil and the unsaturated zone can have a large influence on the infiltration of precipitation into an aquifer. Soil permeability presents soil property that can be very variable. The main objectives of this study were to estimate differences in soil permeability in two dominant types of soils—Fluvisols and Eutric Cambisols—in the area of the Zagreb aquifer and determine the relationship between the physicochemical properties and the permeability of the analyzed soils. For this purpose, the particle size distribution, soil water retention curves (SWRCs), hydraulic parameters, and chemical properties of soils (i.e., electrical conductivity (EC), pH, carbonate content, cation exchange capacity (CEC), and total concentration of zinc (Zn)) from six soil profiles were observed. In general, the results show that Fluvisols have a smaller amount of clay and a higher amount of sand. Furthermore, particle size distribution indicates that Eutric Cambisols have smaller permeability and a slightly higher capacity for retention. In Eutric Cambisols, the percentages of clay, silt, and sand generally do not change with depth. On average, Fluvisols and Eutric Cambisols become impermeable when they reach different values of water content and effective saturation. All results suggest that Fluvisols generally have to desaturate more than Eutric Cambisols to become impermeable. The proportions of sand and Ks increases through the depth of all analyzed soils, while CEC, EC, and Zn decrease. The total Zn generally decrease with depth, which can be attributed to the aerodeposition in the surface horizons of the analyzed soils and their higher availability for binding/sorption elements. Generally, it can be seen that most Zn concentrations increase until 80 cm in depth, after which they decrease. This result indicates that, in these specific locations, the groundwater body is not under the influence of a potentially toxic metal, in this case Zn. Statistical analysis shows a strong correlation between Zn concentrations and some soil properties, such as soil texture and CEC. This may point to the prevalence of Zn retention. Furthermore, statistical results show that silt has a higher influence on the permeability of Eutric Cambisols than Fluvisols.


Introduction
The Zagreb aquifer represents the main source of potable water for the inhabitants of the City of Zagreb and part of Zagreb County. This aquifer is designated as a part of the strategic water reserves of the Republic of Croatia. Due to the anthropogenic influence, problems regarding groundwater quality have been observed. Five main groups of contaminants have been identified: nitrates, pesticides, potentially toxic metals, chlorinated aliphatic compounds, and pharmaceuticals [1]. Although different means of exploration regarding groundwater contamination have been undertaken, most efforts were

Study Area
The study area (the extent of the Zagreb aquifer) is situated in the northwest part of the Republic of Croatia (Figure 1). This area includes the City of Zagreb and its surroundings, covering approximately an area of 300 km 2 . In the research area, different types of soil have been observed, and two types dominate: Fluvisols and Eutric Cambisols on Holocene deposits [29][30][31]. Fluvisols cover approximately 43% of the investigated area, while Eutric Cambisols are developed in approximately 32% of the investigated area. Fluvisols are generally permeable, except during dry periods [32], while Eutric Cambisols developed on the Holocene deposits suggest less permeable soil, where percolation into the aquifer is negligible [33]. The investigated soil profiles are determined as Fluvisols and Eutric Cambisols developed on Holocene deposits [31]. Fluvisols are found on the alluvial plains, river fans, and valleys. This type of soil is characterized by weak horizon differentiation, although a distinct Ah-horizon may be present. The texture of these soils is mainly loam (rarely clayey loam) while their structure is crumbly or blocky. Eutric Cambisols combine soils with at least an incipient subsurface soil formation. Transformation of the parent material is evident from the structure's formation and its mostly brownish discoloration, as well as its increasing clay percentage and/or carbonate removal. The texture of these soils is silty clay loam and the soil's structure is mainly granular. The presence of the clayey cambic horizon, with a pedological profile structure of A-(B)w-C, is characteristic for the specified soil [31]. It has been shown that the thickness of the unsaturated zone varies from two meters in the southeastern part of the Zagreb aquifer to eight meters in the northwestern part [34], while in high waters, some areas can be flooded [35]. A great part of the unsaturated zone, especially in the left bank of the Sava River in the urban part of the City of Zagreb, has been disintegrated by anthropogenic influences, such as traffic, aerodeposition, industry, and illegal waste dumps. The Zagreb aquifer consists of the two aquifer layers, which were deposited during the Middle and Upper Pleistocene and Holocene. Holocene deposits are alluvial, while Pleistocene are lacustrine-marshy [36][37][38]. The Middle Pleistocene unit is relatively uniform in composition. While are lacustrine-marshy [36][37][38]. The Middle Pleistocene unit is relatively uniform in composition. While the lower and middle part is predominantly composed of the grey colored sands, the upper part comprises grey colored or red to yellowish-brown mottled silt and clay sized material. During the warmer periods, the tops of these sediments were exposed to pedogenic processes. The Upper Pleistocene unit is characterized by frequent lateral changes of its gravels, sands, silts and clays. The Holocene unit is composed of pale, yellowish-grey, colored gravels and sands, in which limestone cobbles prevail [37,39]. All investigated soil profiles are situated in the Holocene sediments.
The Zagreb aquifer ( Figure 1) is an alluvial unconfined aquifer, which is in direct contact with the Sava River, the main source of recharge. The general groundwater flow direction is from west to east and south-east. Although recharge also occurs from precipitation, the Sava River presents the dominant boundary condition, which is the most prominent in the vicinity of the Sava River [40]. The groundwater levels are declining. The main reasons for the lowering groundwater levels are (1) extensive riverbed erosion due to upstream Sava River regulation and gravel exploitation from the river; (2) embankment of the Sava River, which stopped occasional flooding of the area and potential infiltration to the groundwater; (3) excessive abstraction for municipal and industrial needs; and (4) prolonged periods of drought [1].

Materials and Methods
Each soil profile was classified according to the IUSS Working Group WRB [31]. Three soil profiles were determined to be the Fluvisols, while three were determined to be the Eutric Cambisols developed on Holocene deposits. Two profiles were located in the left bank of the Sava River, while four of them were located on the right bank of the Sava River ( Figure 1). A detailed description of the soil profiles can be found in previous research [32,33,[41][42][43]. Profile 1 is located in the eastern part of the Zagreb aquifer near the future regional well field Kosnica in the Fluvisols. The particle size analysis of the soil profile was investigated up to 210 cm deep in six different horizons. Soil particle size analysis was made by the pipette-method with wet sieving and sedimentation after dispersion with sodium-pyrophosphate (Na4P 2 O 7 , c = 0.4M). Profile 2 is located in the area of the Stara Loza well field in the western part of the aquifer (also determined to be Fluvisols). The particle size analysis was investigated up to a 170 cm depth in the six different soil horizons. Profile 3 is located in the area of the well field Petruševec, which is located in the Fluvisols till the depth of 190 cm in the eight different soil horizons. Profiles 4 to 6 are located in the area of the well field Velika Gorica, while the investigated depth varies from 90 to 120 cm in the 3 to 6 soil horizons. The particle size analysis of all soil profiles (by depth) is shown in Table 1. In sum, 34 soil samples were collected from the soil profiles. The soil samples were air-dried and passed through a 2 mm sieve for laboratory analysis (i.e., determination of electrical conductivity, pH, carbonate content, cation exchange capacity (CEC) and total concentration of zinc). Electrical conductivity of the soil was measured in the water with a 1:5 soil to water ratio using a Multi 340i WTW EC meter. Soil pH was measured in the water with a 1:5 soil to water ratio [44] using a Multi 340i WTW pH meter. Carbonate content was determined gas-volumetrically using a Scheibler apparatus [45]. CEC was determined using a barium chloride solution [46]. The total concentration of Zn was analyzed by flame atomic absorption spectrometry (AAnalyst 700, Perkin Elmer) after totally dissolving the samples using an aquaregia solution. Soil hydraulic parameters were estimated using Rosetta Lite software [47] and are shown in Table 2. Theoretical SWRCs and the estimated values of the unsaturated hydraulic conductivities [26] were constructed and calculated using Equations (1)-(4): where θ r is the residual water content (cm 3 /cm 3 ), θ s is the saturated water content (cm 3 /cm 3 ), h b is the bubbling pressure (cm), K(θ) is the unsaturated hydraulic conductivity (cm/day), K s is the saturated hydraulic conductivity (cm/day), S e is the effective saturation (%), ψ is the matrix potential (cm), while α (1/m), n, and m are parameters that depend on the slope of the curve [26,48]. For statistical analysis, the basic statistical parameters (i.e., average value, standard deviation, coefficient of variation, and correlation coefficient) were used. Further, correlation analyses was done between the clay, silt, sand content, estimated saturated hydraulic conductivity, and chemical properties of the investigated soils. Firstly, a correlation analysis was done on all available data, and after it was done separately for data from the Fluvisols and Eutric Cambisols. All statistically significant values are marked in red. The chemical properties of the investigated soils are shown in Table 3.

Relationships Between the Soil Hydraulic Data
Statistical insight into the observed data has shown that the investigated soils in the area of the Zagreb aquifer contain, on average, 17.91% clay, 64.06% silt, and 18.03% sand, with the highest standard deviation and coefficient of variation for silt (Table 4). However, if the soils are investigated separately, it can be seen that the Fluvisols have a smaller amount of clay and a higher amount of sand, while the coefficients of variations are very similar. The first insight into the particle size distribution suggests that Eutric Cambisols should have a smaller permeability and a slightly higher capacity for retention due to the calculated hydraulic parameters (Table 4). This is consistent with the average values of the saturated hydraulic conductivity which are more than 50% higher in the Fluvisols. If the data are evaluated graphically (Figure 2), it can be seen that in the Fluvisol percentage of sand increases with depth, while the silt percentage continues to decrease. Soil texture (percentages of clay, silt, and sand) generally does not change with depth in the Eutric Cambisols.  In the unsaturated zone, hydraulic conductivity decreases with a decrease in water content. Increased Ks values in the surface horizons of the Fluvisols could be due to the presence of root channels and macroporosities. Similar conclusions were reported in some of the previous research [49,50]. The authors found that the root systems increased water flow, as indicated by the higher Ks, macroporosities, and water recharge rates of the Ap horizon. The lowest Ks values were determined in the Eutric Cambisols on the Holocene deposits. This is due to the soil texture and close contact with gravels and sands. Similar results were reported by Wegehenkel et al. [51].
It is very important to monitor the water content in time because it can show the variability of unsaturated hydraulic conductivity and related permeability [32]. However, this kind of data can In the unsaturated zone, hydraulic conductivity decreases with a decrease in water content. Increased K s values in the surface horizons of the Fluvisols could be due to the presence of root channels and macroporosities. Similar conclusions were reported in some of the previous research [49,50]. The authors found that the root systems increased water flow, as indicated by the higher K s , macroporosities, and water recharge rates of the Ap horizon. The lowest K s values were determined in the Eutric Cambisols on the Holocene deposits. This is due to the soil texture and close contact with gravels and sands. Similar results were reported by Wegehenkel et al. [51].
It is very important to monitor the water content in time because it can show the variability of unsaturated hydraulic conductivity and related permeability [32]. However, this kind of data can usually be used only when the pedological pits are excavated and different instruments are installed [52][53][54]. Additionally, water content can also be defined when sampling is carried out [33], but, in that case, unsaturated hydraulic conductivity can be only estimated during the sampling time. However, differences in the permeability between various materials can also be estimated by the creating theoretical SWRCs and by inspecting the range of values of the unsaturated hydraulic conductivities. According to Ružičić et al. [32], water content varies from 39.74% to 46.41% in the investigated Fluvisols developed in the area of the Zagreb aquifer. Borek and Bogdał [55] reported a similar distribution of water content in the Fluvisols. Figures 3-8 show the SWRCs and the full range (for the water content from θ r to θ s ) of the possible unsaturated hydraulic conductivities of all the analyzed soil profiles, together with the full range of the possible effective saturation. It can be seen that the soil depths where the highest percentage of sand is defined show faster desaturation (i.e., their slopes are steeper). Furthermore, they also have smaller values for the possible residual water content. On the other hand, when high percentages of silt and clay are recorded, higher values are obtained for the residual water content and porosity. In addition, it is obvious that in soils with a higher sand percentage, the start of desaturation occurs under a smaller pressure. Unsaturated hydraulic conductivities can also vary widely-from about 50 cm/day up to less than 1 × 10 -32 cm/day, depending on the level of saturation. Kercheva et al. [56] reported similar hydraulic soil properties in the Fluvisols. In Profile 1 (Figure 3, Table 5), it can be seen that the impermeable boundary is reached when the water content is 0.16 and 0.18, and the effective saturation is between 29% and 41%. In Profile 2 ( Figure 4, Table 5), the impermeable boundary is reached when the water content is between 0.19 and 0.3, and the effective saturation is between 31% and 50%. Similar numbers can be seen in Profile 3 ( Figure 5, Table 5), where the impermeable boundary is reached when the water content is between 0.17 and 0.21, and the effective saturation is between 31% and 36%. In general, in the Eutric Cambisols, higher water content (i.e., smaller level of desaturation) is necessary for to establish the impermeable boundary. In Profile 4 ( Figure 6, Table 5), the water content between 0.20 and 0.22 and effective saturation between 32% and 36% are necessary to reach the impermeable boundary. In Profile 5 ( Figure 7, Table 5), the water content between 0.22 and 0.30 generates impermeability, and the effective saturation is between 40% and 51%. Similar values are obtained for Profile 6 ( Figure 8, Table 5), where the impermeable boundary is reached when the water content is between 0.23 and 0.30, and effective saturation is between 42% and 49%. Generally, it can be seen that the soil profiles in the two most dominant soil types in the area of the Zagreb aquifer have different particle size distributions that generate different levels of retention. On average, Fluvisols become impermeable when they reach a water content of 0.19 with an effective saturation of 36%, while Eutric Cambisols become impermeable when their water content reaches 0.25 with an effective saturation of 42% ( Figure 9). All results suggest that effective infiltration of the precipitation would be smaller if the Eutric Cambisols were developed, which is consistent with the study that found Fluvisols to generally have good permeability with low water retention capacity and Eutric Cambisols to have good water retention [57]. become impermeable when their water content reaches 0.25 with an effective saturation of 42% ( Figure 9). All results suggest that effective infiltration of the precipitation would be smaller if the Eutric Cambisols were developed, which is consistent with the study that found Fluvisols to generally have good permeability with low water retention capacity and Eutric Cambisols to have good water retention [57].              Figure 8. SWRCs, effective saturations, and unsaturated hydraulic conductivities in soil profile 6. Figure 9. Difference in the SWRCs, effective saturations, and unsaturated hydraulic conductivities between the Fluvisols and Eutric Cambisols. Table 5. Values of the water content and effective saturation that generate impermeability in the investigated soil profiles.

Relationship Between Water Movement and Physicochemical Properties of Soils
The chemical properties of the analyzed soils are shown in Table 3. The highest CEC and EC values are determined at almost all soil subsurface depths (Table 3). Due to their higher organic matter, clay minerals, and iron oxyhydroxides, subsurface soil horizons-in almost all cases-show higher CEC and EC values. Correlation matrices are shown in Table 6. All statistically significant results are marked red (p <0.05). Results were evaluated both together and separately for the Fluvisols and Eutric Cambisols. When all the data are correlated, much more statistically significant correlation coefficients are obtained. Soil electrical conductivity may influence the CEC values [58]. However, in this research, correlation analyses did not produce statistically significant results between the CEC and EC. A strong relationship between Zn concentrations, CEC, and water content in soils was observed [59]. When a statistically significant correlation analysis of Zn with all other variables is evaluated, it can be seen that, in the Fluvisols, Zn is only negatively correlated with pH and with EC in the Eutric Cambisols. When all the data are correlated together, Zn concentrations have a statistically significant positive correlation with CEC and clay, and a negative correlation with pH, carbonate content, sand, and K s . The correlation coefficients suggest that Zn concentrations are more closely related to the pH, carbonate content, and sand components, and higher Zn concentrations can be expected in soils with lower pH, less carbonate content, and fewer sand components. According to Mertens and Smolders [60], Zn concentration is increased as the pH decreases in soils. Table 6. Correlation matrices between the chemical properties, particle size distribution, and saturated hydraulic conductivity (correlation coefficients marked red are statistically significant).  Although a statistically significant moderate positive correlation between Zn and CEC is present, this correlation is not recorded when the correlation analysis is made separately. A very low negative correlation between Zn and CEC is present in the Fluvisols, and a very low positive correlation is present in the Eutric Cambisols. Due to increase of sand content in the Fluvisols, it seems that Zn more commonly binds to fine particles in the Eutric Cambisols than in the Fluvisols. In some research, positive correlations between Zn and siltous materials have been reported [61]. This correlation is not very pronounced in the present study, although in all three cases, there is a very low-moderate positive correlation between the silt components and Zn, which indicates that Zn is partially bounded in the silt component of the analyzed soils. The total Zn generally decreases with depth, which can be attributed to the aerodeposition in the surface horizons of the analyzed soils and their higher availability for binding/sorption elements. A similar distribution of Zn through the profiles of the Fluvisols and Eutric Cambisols have been found [62][63][64]. In some studies, a highly positive correlation between EC and siltous material has been pointed out [17]. Electrical conductivity generally decreases with depth, but no correlation between silt and EC can be seen in the Eutric Cambisols. On the contrary, in the Fluvisols, a moderately positive correlation between silt and EC is present, which is consistent with the moderately negative correlation between the EC and sand components. In general, the proportion of sand and K s values increase through the depth of all analyzed soils, as CEC, EC, and Zn are decrease. Results also suggest that particle size distribution in the Fluvisols has a greater influence on CEC values, carbonate content, and EC, which is mostly related to the silt and sand components. Changes in pH values along the soil profile can be associated with a change of carbonates, while the alkalinity of the soil can be attributed to the presence of magnesium and calcium carbonates [41], which generally coincides with the correlation results in this research. Furthermore, when all data are correlated, the results suggest that silt is an unessential parameter that does not have any influence on hydraulic conductivities. This is evident within the correlation results for the Fluvisols. However, in the Eutric Cambisols, the silt components have a very positive correlation, which could be the consequence of the lower average content of the sand and the higher average content of the clay than in the Fluvisols. This impact is confirmed by the moderately negative correlation of the saturated hydraulic conductivity with the sand. It seems that in soils with high clay and low sand content, silt can have a very high impact on the soil's permeability. All results indicate that a statistical analysis should be used to evaluate the different levels, which could point to different results and suggest different conclusions.

Conclusions
The Fluvisols and Eutric Cambisols that developed on the Holocene deposits represent the two most dominant types of soils in the area of the Zagreb aquifer. Generally, it can be seen that the soil profiles in the two most dominant soil types in the area of the Zagreb aquifer have different particle size distributions, which generate different retention and possible infiltration parameters. The sand percentage in the Fluvisols increases with depth, while the silt percentage decreases. This can be explained with the movement of the Sava River in the past and formation of various meanders with bigger particles. In average, the Fluvisols become impermeable when they reach a water content of 0.19 and an effective saturation of 36%, while the Eutric Cambisols become impermeable when their water content reaches 0.25, with an effective saturation of 42%. The increased K s values in the surface horizons of the Fluvisols could be due to the presence of root channels and macroporosities. In addition, this study revealed a strong relationship between soil permeability and some chemical properties of the analyzed soils. Most of the highest Zn concentrations were determined in the first 30 cm of the analyzed soils, while the lowest concentrations were observed in the deepest soil depths, which can be attributed to the aerodeposition in the surface horizons of the analyzed soils and their higher availability for binding/sorption elements. Generally, it can be seen that most Zn concentrations increase to 80 cm depth and decrease afterward. In general, the statistical analyses showed strong correlations between the Zn concentrations and some physicochemical soil properties. Zn is negatively correlated with the pH, carbonate content, the sand component and K s . The Zn concentration increases as pH decreases in the analyzed soils. In addition, a negative correlation was observed between the total Zn and sand content, which can be explained by Zn's availability to bind to smaller particle size materials, in this specific case, siltous material. The results suggest that sand represents a very important variable in the estimation of the transport of Zn through the soils in the area of the Zagreb aquifer. However, this result should be tested in more detail in future research. Furthermore, the statistical results indicate that silt is an unessential parameter that does not have a big influence on the values of hydraulic conductivities in the Fluvisols but can be very important for estimating permeability in the Eutric Cambisols. All results suggest that percolation of the surface and soil water from precipitation, as well as Zn transport, is possible in larger quantities where the Fluvisols are developed. This is confirmed with the physicochemical properties and Zn availability of the Fluvisols. All results suggest that water movement or percolation is related with the physicochemical properties of soils, especially the particle size distribution of the analyzed soil. A detailed chemical analysis of the potentially toxic metals, which are also derived from the agricultural activities (e.g., Cd, Cu) and the creation of water flow and transport models using field measurements in the analyzed soils should be of primary interest in future research. Future research should alsofocus on monitoring the physical and chemical properties of the soil and the unsaturated zone in time and in the different parts of the aquifer, but also their interpretations at different scales. In addition, future research needs to consider the distribution of soil cracks and voids, which can also be important for water percolation with contaminants in soils.
Author Contributions: Z.K. have participated in the different previous studies, the data of which were used in this article. S.R. and Z.K. developed the concept of the research and wrote most of the article. Z.K., D.P., and L.B. produced the statistical analysis. D.P. made Figure 1 in ArcMap. Z.K., L.B., and D.P. participated in the technical editing of the figures and article. L.M. participated in the laboratory analysis. All authors participated in the writing of the article.
Funding: This research received no external funding.