Open AccessArticle

The Dunes of Belvedere–San Marco of Aquileia: Integrating High-Resolution Digital Terrain Models and Multispectral Images with Ground-Penetrating Radar Survey to Map the Largest System of Continental Dunes of Northern Italy

Federica Vanzani

¹,

Alessandro Fontana

^1,*

Livio Ronchi

¹,

Jacopo Boaga

Veronica Chiarini

and

Irka Hajdas

Department of Geosciences, University of Padova, Via G. Gradenigo 6, IT-35131 Padova, Italy

Ion Beam Lab, ETH Zurich, Otto-Stern-Weg 5, 8093 Zurich, Switzerland

Author to whom correspondence should be addressed.

Remote Sens. 2024, 16(5), 765; https://doi.org/10.3390/rs16050765

Submission received: 12 January 2024 / Revised: 15 February 2024 / Accepted: 20 February 2024 / Published: 22 February 2024

(This article belongs to the Special Issue Remote Sensing and Geophysics Methods for Geomorphology Research)

Download

Browse Figures

Versions Notes

Abstract

The interpretation of high-resolution remote-sensed data (i.e., LiDAR-derived DTMs, aerial photos and satellite images), compared with ground-penetrating radar surveys, historical cartography, geomorphological surveys and stratigraphic data, allowed us to map a large system of dunes near the Grado-Marano Lagoon (NE Italy) and reconstruct its evolution. Remote sensing investigations allowed us to recognize, map and interpret the sandy reliefs as a field of continental aeolian landforms extending for over 15 km² and consisting of parabolic dunes elongated in the WSW direction. Radar soundings, together with the description of stratigraphic sections and cores, documented the internal clinostratification of the dunes, supporting their aeolian origin. Radiocarbon dating documents that the dunes formed 22 ka ago, at the end of the Last Glacial Maximum, and probably evolved until the first part of the Late Glacial, when vegetation was scarce. The landforms were fed by the sands blown from a paleochannel of Isonzo River flowing eastward of the dune’s field and blown by Bora. This is a very strong katabatic wind, still characterizing the area, but that was likely much stronger during last glaciation, when it was probably sustained by a stronger wind pattern in Central Europe.

Keywords:

light detection and ranging (LiDAR); aeolian geomorphology; parabolic dunes; Friulian plain; Bora wind; Last Glacial Maximum

1. Introduction

Trieste and its Gulf represent the northernmost sector of the Adriatic Sea and are characterized by the Bora, one of the strongest winds in the Mediterranean Basin. Bora is a katabatic wind blowing from ENE, with gusts reaching peak velocities between 30 and 50 m/s [1]. Despite the persistence of this wind, related aeolian landforms are rare in this area and have been discovered only recently. An example is represented by the remnants of fossil dunes connected to the action of Bora, which were documented under the seafloor in front of the Trieste harbor [2]. Another group of relict landforms attributed to aeolian processes is located near Aquileia, along the margin of the Grado-Marano Lagoon. It corresponds to several small sandy reliefs called the dunes of Belvedere–San Marco (Figure 1). Their origin has been intensely debated (ref. [3] and reference therein) and their study is the object of this paper.

Aeolian landforms are generally present along the Italian Peninsula as coastal dunes, along marine beaches and deltas, formed during the Holocene or during previous interglacials periods in a few rare locations. Less frequently, the sedimentary bodies accumulated during marine lowstand conditions, when part of the continental shelf was exposed and prone to deflation. This is the case of the aeolianites of Isola d’Elba [4] and some other deposits located in Sicily, e.g., ref. [5,6] and Sardinia [7].

Figure 1. (a) Square b indicates the location of the study area, south of the city of Aquileia, while c indicates the San Giuliano Isle. The basemap is a hillshade elaboration of the Shuttle Radar Topography Mission (SRTM) data [8]. (b,c) show the study area with the localities’ names used in this work. Red rectangles indicate the locations of Figure 6 (Rectangle 1), Figure 7 (Rectangle 2) and Figure 8 (Rectangle 3). Yellow rectangles indicate the location of Figure 4a (Rectangle 4), Figure 4b (Rectangle 5), Figure 4c (Rectangle 6) and Figure 4d (Rectangle 7). The white star indicates the location of the core BLG1. The green star indicates the location of the outcrop in Figure 9. The purple stars indicate the location of cemented sand samples in Figure 10.

In northern Italy’s inland territory, aeolian landforms are almost completely absent in the alluvial plains. In contrast, aeolian deposits are rare, and they usually consist of thin bodies of loess that sedimented during the Pleistocene glacial periods on raised terraces or over the mountain slopes facing the plains [9,10].

The situation in northern Italy strongly contrasts with Central and Northern Europe, where continental aeolian landforms and deposits are widely diffused, leaving an important overprint on the present landscape. The so-called “Northern European Sand Belt” characterizes large portions between Great Britain and Russia, where Pleistocene aeolian dunes are well documented in the lowlands, as in the Rhine–Meuse Delta [11,12,13]. In Central Europe, abundant aeolian landforms and deposits are found in Hungary, Serbia, the Czech Republic and Romania [14,15,16,17].

In northern Italy, the only system of continental dunes that has been mapped in detail is located in the Venetian Plain, west of the Euganean Hills and near the city of Este, but its chronology and genesis are unclear [18]. In this work, we analyze in detail the much larger dunes field of Belvedere–San Marco, which corresponds to the place where, according to the local legend, the evangelist Saint Mark made landfall during his voyage from Alexandria of Egypt to Christianize the regions of the Northern Adriatic (hence the toponym). These landforms have represented one of the main unresolved geomorphological questions of the northeastern Italian alluvial plains for over 150 years. They have been studied since the 19th century, and several contrasting ideas were put forward about their origin and age, mainly spanning from the hypothesis of Holocene or Pleistocene coastal dunes [19,20,21] to the ones proposed by the geologist Ruggero Marocco. He discussed their possible relations with the remnants of Pleistocene fluvial ridges [22], with the occurrence of a tectonic terrace in the area [23], and, for the first time, he proposed the possibility that they were continental dunes, whose evolution might have been related to the reworking of ancient fluvial channel deposits [24,25]. The dunes of Belvedere–San Marco were recently reported in the map of the geological units of the plain of the Friuli–Venezia Giulia Region at the 1:250,000 scale [3] as continental dunes formed during the Last Glacial Maximum (LGM). In this document, they were mapped as a single unit due to their small cartographic scale, and no details about their planar morphology, internal geometry and stratigraphy were included. Thus, our work aims to analyze the Belvedere–San Marco dunes, adopting a multidisciplinary approach involving the investigation of the geomorphological, stratigraphic and sedimentological aspects to understand the origin and chronology of these peculiar reliefs. In particular, we used elaborations of DTMs, different remote-sensed images and historical maps to reconstruct their planform. Moreover, we compared new stratigraphic data with geophysical information to infer their internal characteristics, age and evolution processes. These landforms have also been compared with other fields of dunes in Europe, to highlight similarities and differences.

2. Setting

The study area is located along the inner margin of the Grado-Marano Lagoon, in the coastal plain of Isonzo River, where large sectors of the former lagoon and palustrine environments were reclaimed and transformed into farmland between the 18th and 20th centuries [26]. This area is part of the distal portion of the Friulian Plain, corresponding to the foreland basin of the southeastern Alps, where the bedrock is found at a depth between 250 and 350 m [27]. Alluvial, deltaic and marine processes interacted during the late Quaternary, and a very significant sedimentary aggradation occurred during the Last Glacial Maximum (LGM, 28,000–19,000 years ago; ref. [28]). In that period, large glaciers occupied the Alpine catchments. They fed the major fluvial systems with extraordinary liquid and solid discharge, leading to the construction of alluvial megafans like those of the Isonzo and Torre rivers, which formed the study area [29]. These large fluvial systems were mainly characterized by the deposition of fine sediments (silt and clay) in their distal portion, where the groundwater table is now found within 1–3 m of the surface [29].

Due to the glacial conditions, the sea level was 120 m lower than today and the marine lowstand exposed the continental shelf of the Adriatic for over 300 km [10]. According to paleobotanic investigations and recent climate modelling, the LGM alluvial plain of NE Italy was characterized by significantly cooler temperatures and the presence of steppe-like environments (refs. [10,30] and references therein).

During the final phase of the LGM, the Alpine glaciers started retreating upstream along the Alpine valleys, causing a massive drop in the sedimentary fluxes and triggering the formation of incised fluvial valleys in the alluvial plain [29]. In the Isonzo catchment, this condition was further enhanced by landslides in the main valley near Bovec (Slovenia), which dammed and trapped the downstream transport of coarse material almost completely [31]. Between the end of the LGM and the early Holocene, the above-described conditions led to the confinement of fluvial deposition inside the incised valleys also in the distal plain, where they were up to 1 km wide and 20 m deep [29]. At the bottom of these depressed landforms, the confined flux allowed the transport to the current coastal plain of gravel and sand eroded in the apical portion of the alluvial plain (ref. [32] and references therein), which is documented in the area west of Aquileia, cf. [25].

A significant change occurred about 8000 years ago, when post-LGM marine transgression led the sea to reach a position slightly comparable to the present one and triggered the formation of deltas and lagoons [33,34,35,36]. The marine highstand promoted sedimentation in the distal plain: the final portions of the fluvial incised valleys were progressively filled with estuarine and fluvial deposits and further expanded over the external alluvial plain [29]. This phase led to the formation of the coastal plain of Isonzo River, which consists of a coastal wedge of lagoon and alluvial deposits whose maximum thickness is about 10 m near Grado and tends to become thinner landwards [33,35]. Thus, south of the Aquileia–Monfalcone latitude, the surface of the LGM alluvial plain has been largely buried by the Late Holocene sedimentation, and the related channel belts of Isonzo River and Torre Torrent are characterized by wide fluvial ridges elevating 1 to 3 m over the surroundings [23,25]. Presently, Torre is a tributary of the Isonzo River, but in the past this stream reached the lagoon directly, as witnessed by traces of its paleochannels documented west of Aquileia [23,37].

The dunes of Belvedere–San Marco rise from the landscape for their topographic evidence south of Aquileia, reaching a maximum elevation of about 10 m above sea level (asl) in the highest sites (near the church of San Marco and the Belvedere Camp Ground Isola Centenera). The extent of the studied dunes is about 15 m², with a maximum north–south length of 4 km, and a maximum east–west length of 6 km. However, if the most western dune, located in the Isle of San Giuliano, is considered, the east–west length reaches 10 km. These dune ridges generally display an ENE–WSW alignment and are largely surrounded by recent lagoon deposits, except along the northern side, where Holocene alluvial deposits are present (refs. [3,24,34]; Figure 1a). Remnants of former partly cemented sandy reliefs have also been found inside the Grado Lagoon, as documented on the islands of Mottaron, Volpera, Gorgo and San Giuliano, which have been largely eroded in the last centuries because of sand mining [33]. Until the first part of the 20th century, the pedological profiles above the dunes in the area near Beligna and Belvedere were still characterized by the occurrence of some rubified horizons with an accumulation of illuviation clays and oxides, especially in the depressions between the ridges [21]. In the deeper soil horizons up to some meters of depth, decimetric carbonate concretions and cemented layers are still frequent and are locally called “maltoni”. This term derives from the accretive version of the Italian word “malta”, which means mortar, and was given because of the cementation of the sand and their large dimensions. Anyhow, the modern reclamation strongly altered the original topography, and the continuous ploughing and levelling practices over the sandy ridges led to an intense reworking of the soil, often involving the unweathered sediments. Today, the plough layer is enriched with abundant fragments of carbonate concretions that partly hinder its cultivation.

The sands of the dunes have been tentatively dated through thermoluminescence dating, with results spanning between 87,000 and 176,000 years before present [23], but these ages turned out to be too old to fit with the occurrence in the subsoil of coastal deposits of the Last Interglacial at a depth of 35–60 m [28,38].

3. Materials and Methods

This research was designed with a multidisciplinary approach, combining the results of remote sensing investigations with the direct information collected in the field and lab analyses of sediment samples. The elaboration of LiDAR (Light Detection and Ranging), derived digital terrain models (DTM) and the analyses and interpretation of aerial photos and satellite images have been compared with historical cartography and coupled with the results provided by the study of new stratigraphic cores, field geological surveys and the acquisition of geophysical profiles with GPR (ground-penetrating radar). Moreover, some selected samples of large pedogenetic concretions (“maltoni”) were investigated in detail for their sedimentary properties, while some organic layers have been radiocarbon dated.

The first phase of this work consisted of the use of remote sensing techniques to produce a detailed map of the study area for the reconstruction of the sandy ridges’ planform patterns in multiscale views, and eventually to relate them to specific geomorphic processes. Moreover, these analyses were essential to guide the following field activities through the selection of specific sites to be investigated through field surveys and corings. In fact, the research considered a lowland area, in which the presence of the dunes often corresponds to subtle topographic variations and small soil and sediment color differences that are generally not visible in the field.

Photo interpretation was carried out using QGIS software (version 3.22.5) to combine the information obtained from conventional aerial photos and multispectral satellite images. We considered the aerial photos taken during flights performed in 1938, 1945, 1954, 1984 and 1993, which are available from the IGM (Istituto Geografico Militare) portal [39]. Their specific properties are listed in Table 1. They are all characterized by a high geometric resolution and a low ground distortion, thanks to the flat topography of the alluvial plain. Most of the photos were already available in the digital version, as they were previously acquired with a scanner at 800 dpi resolution, but for some of them, especially for the flights of 1954 and 1993, we acquired the printed version of the photos of interest with a scanner at 600 dpi, and we georeferenced them in ArcGIS Pro. It is worth noting that the oldest datasets (1938 and 1945) allowed us to study the traces of some dunes which nowadays have been completely eroded or almost vanished because of the intensive agricultural practices affecting the area.

Satellite images were selected from the Sentinel-2, Aster and Planet images [40] (Table 2) and were compared with available images in the Google Earth databases (from Ikonos, Quickbird and GeoEye) and Esri basemaps. Apart from the latter, the satellite geometric resolution is quite low (see Table 2) but, by stretching the RGB color composite images, we were able to emphasize the spectral signature of geomorphological visible traces. To quantify the specific spectral signature of the dunes, we performed a boxplot analysis on Aster and Sentinel-2 imagery. To differentiate the Belvedere–San Marco dunes from the surrounding alluvial and lagoon deposits forming the coastal plain, we focused on the signatures of the visible, near-infrared and short-wave infrared bands. A box plot, or box and whiskers plot, is commonly used to show the distribution of a quantitative variable in different categories. The boxes show the lower quartile (25% of observations below it), upper quartile (75% of the observations below it) and the median, thus representing the interquartile range, corresponding to half of the observations. The whiskers in the plot show the rest of the dataset was contained in 1.5 interquartile ranges, except for the outliers, represented by separate points [41,42]. This analysis was carried out in Python (version 3.10) environment, using the Seaborn library.

One of the most striking characteristics of the Belvedere–San Marco dunes is their topographic relief, which strongly contrasts with the surrounding distal plain, generally below the elevation of 0 m asl. We used data obtained from LiDAR datasets accessible through the geoportal of the Friuli Venezia Giulia Region [43], where two different datasets are available. The first was obtained in 2006–2010, through a ALTM Optech 3033 sensor (manufactured by Teledyne Optech, Vaughan, ON, Canada), mounted on a AS350 helicopter. The point density is 4–5 points/m² and the final result is a DTM with 1 m resolution and a horizontal accuracy of 10–80 cm. The second dataset was acquired in 2017–2020 using a RIEGL LMS-Q780 sensor (manufactured by Riegl, Horn, Austria), mounted on a helicopter. The point density is about 16 points/m² on plain regions and the obtained DTM has a resolution of 0.5 m. In both datasets, the DTM has already been filtered for vegetation and anthropogenic constructions, and the nominal vertical accuracy over the whole region is ±15 cm, although in the study area it is generally much more accurate. The high resolution was of critical importance in this case study, as the topography of the dunes has been severely impacted by levelling and ploughing practices, and in some sectors no more than a few decimeters of relief are left. Several stretching techniques, such as slope and aspect analysis, were used to emphasize the height difference. In the QGIS software (version 3.22.5), through the plugin “Terrain Profile”, it was possible to obtain the topographic profiles, both transverse and parallel to the dune ridges (Figure 2).

To further emphasize the morphologies of these ridges, we employed the Relief Visualization Toolbox, developed by the Research Centre of the Slovenian Academy of Sciences and Arts. This toolbox allows us to implement different visualization techniques, and it has been specifically designed to highlight small-scale topographic features [44,45]. Through this toolbox we calculated the raster images obtained for hillshading, slope, local dominance, simple local relief model (SLRM) and multi-scale relief model (MSRM). Local dominance depicts how dominant an observed standing on a given pixel would be in the surrounding area [46].

This value is calculated as the average steepness of the angle at which the observer looks down on the land surface [47]. Local dominance should be higher for pixels on ridges and slopes and lower for depressions. It can be useful in the detection of very subtle reliefs, such as the dunes’ ridges. Local relief models are generally used to emphasize small and local morphological features. They are obtained by calculating the trend with a low-pass filter (smooth filter) and subtracting it from the original relief model [47]. In the case of the multi-scale relief model, several low-pass filters with diverse kernel sizes were applied, in order to emphasize small morphological features with different scales [48]. It is important to keep in mind that, despite preserving the height proportions, SLRM and MSRM do not provide a direct representation of the topographic differences, but they provide an artificial method to visualize changes in elevation.

After obtaining different visualizations through the Relief Visualization Toolbox, we took a further step by performing a principal component analysis (PCA) in GRASS GIS (version 7.8). This is a multivariate statistical method which is often used to reduce the correlation between quantitative variables by rotating the original space and detecting new uncorrelated variables; these explain most variations in the dataset, called principal components [49]. In this case study, we performed PCA on the original DTM, hillshading, slope, local dominance, SLRM and MSRM. We then combined the three principal components explaining the majority of the variance into a RGB composite image (Figure 3).

Historical maps were analyzed to supply data about the changes in morphology and land use in the research area over the last few centuries, aiming to reconstruct the specific areas where the dunes have been reworked or even levelled by anthropogenic activities before or during the 20th century (Table 3). Indeed, large portions of the sandy ridges were not cultivated in the past because of their very poor soils, and they were covered by woods until the beginning of the last century. The maps produced from 1780 to 1887 by the Austrian State, and later by the Austro-Hungarian Empire during the so-called military surveys, have been the main source of information. Most of these documents have a 1:28,800 scale and correspond to the First Military Survey (1785), the Second Military Survey (1822) and the Third Military Survey (1869), which are available on the geoportal Arcanum Maps [50]. The same scale characterizes the “Topographisch-geometrische Kriegskarte von dem Herzogthum Venedig” (Military topographic–geometric map of the Duchy of Venice), surveyed between 1798 and 1805 under the supervision of Anton von Zach, whose hard copy reproduction was published in Rossi [51]. Other historical geometric maps used in this work are the so-called “tavolette”, produced at 1:25,000 scales in 1917, 1927 and 1949 by the IGM [39].

To verify the observations provided by remote sensing analysis and to investigate the sedimentological and stratigraphical characteristics of the soils and sediments of the studied dunes, field surveys campaigns were carried out between 2020 and 2022. Fieldwork was organized between November and May, since during these months the fields are ploughed, providing the best visibility of bare soil. Moreover, in that season, parts of ditches and canals are cleaned or renewed and, along their slopes, long temporary stratigraphic sections are exposed. The activity consisted of the validation of the geomorphological features mapped through remote sensing and of the documentation of outcrops, to check the sedimentological and stratigraphic characteristics of soils and sediments forming the dunes. In particular, the collected data consisted of the description of grain size, Munsel Charts color, stratigraphic layering, presence of sedimentary structures, presence, dimensions, abundance and typology of the pedogenetic features (i.e., carbonatic concretions, Fe-Mn nodules) and fossils (i.e., shells and paleobotanic remains).

The stratigraphic information has been improved with the use of a hand auger, which allowed us to sample the first meters of the subsoil. In particular, we used an Eijkelkap probe with bayonet extension rods and an Edelmann combination-type head above the groundwater table and a 1 m gouge below it. Between the dune ridges, the water table is generally found at about 1–2 m of depth.

Several stratigraphic data were obtained by analyzing both the archive of the Consorzio di Bonifica Pianura Friulana and the databases of some professional geologists, in which they collected the logs of stratigraphic cores and penetrometric tests, originally carried out for building constructions and territorial planning. To correlate the surface deposits to a wider regional stratigraphic setting, we considered the stratigraphic core BLG1, which was drilled in the area of Beligna in 2008, along the slope of one of the highest preserved sand ridges (45°44’54.37”N–13°21’59.06”E; surface elevation 2.2 m asl; Figure 1b). The core was carried out in the framework of the project VECTOR and reached a depth of 60 m from the surface [38,52]. It was possible to directly sample some of the organic-rich layers of this core that were used for AMS (Accelerator Mass Spectrometry) radiocarbon dating at the Ion Beam Laboratory at ETH (Eidgenössische Technische Hochschule) in Zurich (Switzerland). The results were calibrated to 2 sigma using Calib Rev 8.0.1 (IntCal20.14c; ref. [53]).

During the field activities, several samples of loose sediment and 12 samples of cemented sands were collected in the area of Muson and Beligna and they were analyzed in the laboratory for their petrographic and sedimentologic properties. In particular, we obtained some petrographic thin sections of 30 μm from cemented sands and colored them with alizarina, to highlight and distinguish the dolostone clasts from the limestone ones. Moreover, two samples of cemented sands were analyzed through photogrammetry, in order to produce 2D zenithal photos, and 3D images, using the Agisoft software (version 1.8.5). To visualize them, the CloudCompare software (version 2.11.3) was used (see Supplementary Materials).

In May 2021, a campaign of ground-penetrating radar (GPR) was carried out in the areas of the Belvedere Camp Ground (Isola Centenera) and near the church of San Marco along selected profiles, both perpendicular and longitudinal to the direction of the dune ridges. For this geophysics technique, which is based on the propagation of electro-magnetic waves [54], we used a monostatic GPR, where the sender and receiver antennas were mounted in a manually pushed trolley. As generally used for geological purposes cf. [55], we acquired the profiles with an IDS-MF Hexagon GPR system with a double frequency of 200 and 600 MHz, aiming to investigate up to some meters deep with a decimetric resolution. The reflection of the radar pulses happens at the interface between two materials with different electro-magnetic impedance (ε) and, thus, if the velocity of the reflected pulses is known, it is possible to obtain a pseudo-depth based on the arrival time. This allowed us to see possible internal structures and stratifications. Data acquired through GPR were plotted in radargrams, which have the length of the line on the x-axis, and the pseudo-depth on the y-axis. The visibility of the obtained profiles was corrected and improved, using filters to eliminate noise, and was adjusted to follow the true topography of the surface.

4. Results

4.1. Remote Sensing, Light Detection and Ranging (LiDAR) and Comparison with Historical Cartography

The results obtained from the photointerpretation of remote-sensed images and the topographic analyses of DTMs were integrated into the GIS environment to generate a map of the entire study area. The mapped landforms were compared with the maps produced between the 18th and 20th centuries. Among the aerial photos, the ones captured in the 1993 flight were the most significant source of data and, in particular, the frames numbered 1087, 1088 and 1089.

In the photos and the visible wavelength images, the dunes traces are mainly recognizable as soilmarks, characterized by yellowish and whitish areas with very high tones, which appear in strong contrast with the darker color of the surroundings, where alluvial and lagoon sediments are present (Figure 4). The high visibility of the dune ridges is mainly related to the occurrence of unweathered carbonatic sands that have been exposed to the surface in recent decades because of the recent ploughing and/or levelling, which caused the erosion of the topsoil horizons. Moreover, sand terrains have a higher permeability in contrast with the finer sediments of the distal plain and they are, thus, rapidly drained, while the surrounding deposits have a significant moisture content or can be even waterlogged for rather long periods in winter.

Photointerpretation was mostly performed on soilmarks for detecting the boundaries of the sandy areas, whereas the cropmarks were less effective and gave some interesting information only in very limited zones (Figure 4d). In particular, in some high-resolution images, they evidenced the occurrence of a dense pattern of parallel traces, with decimetric width, forming an acute angle between 15 and 35° with the direction of the dune ridge they belong to. As indicated by the field data described below, these tracks likely correspond to the intersection of the clinostratified layer of the dunes with the ground surface, where the ridges have been strongly levelled, leading the plough layer to also affect the unweathered sediments.

Similarly, ridge traces were mapped from satellite images, elaborated through brightness and contrast stretching to emphasize the marks of sandy sediments. Various RGB composite images were produced, mainly using visible and NIR bands. Satellite images were further analyzed to study the spectral signature, using box plots within the Seaborn library in Python.

Box plots allow us to observe the variation in a specific band of the reflectance of an area. The reflectance of the dunes was compared with that of the surrounding plain using SWIR, NIR and visible bands on Sentinel-2 and Aster images.

The box plots show that the spectral signature of dunes is lower than the alluvial plain in the visible bands (band 3 of Sentinel-2), but it progressively increases and exceeds the reflectance of the surrounding plain in red edge and near-infrared bands (band 2 and 3 of Aster; band 6 of Sentinel-2). In SWIR intervals (band 9 of Aster; band 11 of Sentinel-2), the reflectance of the dunes is much higher than that of the surrounding plain (Figure 5). This characteristic spectral signature of the dunes in SWIR bands was important in identifying the majority of these landforms; however, due to the low resolution of the Aster SWIR band (30 m) and the Sentinel-2 SWIR bands (20 m), it is not sufficient to reconstruct the planar morphology of the ridges in detail. Moreover, satellite images do not allow us to recognize dunes where dense vegetation cover is present, as is the case in the San Marco and Belvedere areas. Therefore, the careful integration of historical maps, LiDAR-derived DTMs, satellite images and aerial photos was crucial to reconstruct the planar morphology of the dune field (see Figure 6, Figure 7 and Figure 8).

In general, aerial photos and satellite images were particularly useful in detecting the sandy reliefs where they had been almost completely levelled by human activity, such as in Due Leoni, Beligna and the eastern zones near Viola and Domine. DTM analysis, in contrast, was crucial in mapping the most preserved sectors of the dunes, which are those with the highest relief. In San Marco and in Belvedere Pineta, the topography reaches 10 m asl, versus −1 m asl of the surrounding plain. Furthermore, these areas are also heavily vegetated, making it difficult to reconstruct their morphology through aerial and satellite images. From the DTM, we calculated the hillshade and multi-directional hillshade, and various profiles, both transverse and parallel to the dune ridges, were extracted. The information obtained from aerial photos, satellite and DTM data was combined to produce a single map illustrating the specific planar morphology of the dunes (Figure 6, Figure 7 and Figure 8).

4.2. Field Survey, Stratigraphic Cores and Sediment Samples

One of the most interesting stratigraphic sections was found during the geomorphological and geological surveys in the area of Muson, in the Due Leoni farm (45°44′17.9′’N 13°21′29.7′’E; Figure 1b), where the renovation of a ditch exposed a section of 80 cm of depth from the surface, but over 250 m in length (Figure 9). The first 40 cm corresponded to the plough layer that reworked silty sandy sediments, related to brackish palustrine deposits characterizing the area before the modern reclamation. At 45 cm of depth, an abrupt unconformity was observed, separating lagoon facies from a sandy clinostratified bed. On this unconformity, a lot of sandy cemented concretions were found, with a downfacing lumpy side (Figure 10). Some of them were in their original position, not reworked by anthropogenic activity. The clinostratification was clearly evidenced by the presence of lighter and darker laminae, with a thickness ranging from 0.1 to 3 cm. The laminae were dipping in a northwesterly direction with an angle of 25–28° for a length of at least 100 m along the section. A close-up of the clinostratified layers is represented in Figure 9c, where the different laminae are labelled from 1a to 1f. The different colors reflected a difference in the grain size; that is, coarser for darker layers (>300 μm) and finer for the lighter ones (<180 μm). Sub-millimetric layering was also recognized, but we also frequently found laminae consisting of very coarse sands and fine gravels, where the coarser clasts reached up to 6 mm in diameter. The sediments forming the dunes were characterized by a high grain size variability, ranging from silt to fine gravels (up to 6 mm), and dramatic changes were frequently documented within a few centimeters. As clearly evident in Figure 11, even in the same clinoform layer, the clasts of sands were generally not well sorted and the clasts often displayed medium to low values of roundness. Thus, for the studied dunes, the grain size distribution curves were not a very useful tool for evidencing the aeolian origin of the deposits. This peculiar situation can be easily explained by the immediate vicinity of the sediment source, likely deflated from river channels. The granules did not experience a long transport to allow a selection of grain sizes and morphometries.

The samples of cemented sands collected in the study area were made of the same sediment forming the dunes, and they were commonly found on the topographic surface, over the ridges. They have already been described and photographed in literature by Comel and Marocco [21,24] and they correspond to pedogenetic concretions formed as a result of the dissolution of carbonates in the upper horizons of the soil and their percolation and reprecipitation in the lower horizons (Bk/Ck horizons, according to the pedologic classifications cf. [56]. The cementation of the layers was not continuous and the formation of concretions was strongly influenced by the grain size differences, generally occurring in the less permeable layers. A complex morphology often characterized by lumps was a typical feature of the cemented sand and, in some cases, it seemed to be related to bioturbation. In fact, the different grain sizes and permeability of the burrows generated by roots and animals might have favored the precipitation of carbonates along pre-existing conduits. This might also explain the lumpy and/or elongated shape of the collected samples, which can be found at over 2 m of depth from the surface, as also photographed by Comel in Belvedere Pineta, and Marocco on the Gorgo and Villanova isles (Figure 12). Today, these cemented samples are commonly disseminated in field surfaces, broken in pieces by ploughs, making the inference of their original orientation quite challenging. The samples we collected were a few decimeters in size and presented with a lumpy and elongated shape. Based on the sampled area, they can be differentiated into two main types. The ones collected near the Due Leoni farm (45°44′41.9′’N 13°22′45.1′’E, Figure 1b) were about 10 cm in size and presented with marks of intense bioturbation. Differently, the samples collected near Muson (45°44′22.9′’N 13°22′23.3′’E, Figure 1b), were 20–30 cm in size and, in many of them, it was still possible to identify a clinostratification, which was not completely erased by bioturbation and pedogenetic processes (Figure 10). This suggests that many samples originally formed deeper than 1 m in the soil profile, where bioturbation was less present, and they have been brought to the surface by agricultural activities only in recent decades. The clinostratification was very similar to the one observed in the cross-section, but the scattered distribution of the samples made it generally impossible to infer the original direction and inclination.

As already documented by Comel [21] and Marocco [23,24], the composition of the original sediments forming the dunes is dominated by clasts derived from limestones and dolostones (60–80%), reflecting rocks outcropping in the Isonzo River catchment basin. In the thin sections, it is also clearly visible that the clasts are not well selected, and they are very angular. Among the larger granules, the clasts of chert were quite common (Figure 11).

The BLG1 borehole was drilled in the lower portion of the slope of one of the longest sandy ridges in the northern part of the study area and allowed to reach the base of that dune, which locally has a maximum elevation of 9.2 m asl (Figure 8d).

Holocene colluvial deposits related to the erosion of part of the dune were present in the first 2.5 m, covering an inceptisol with a centimetric fragment of prehistoric pottery. The sedimentary unit corresponding to the dunes was found between 2.8 m and 7.20 m and consisted of medium to coarse sands (Figure 13). They were marked at the base by the contact with sandy silts interpreted as alluvial deposits. The new radiocarbon dates allowed us to demonstrate that the deposition of the dune likely occurred around 22.1–21.4 ka cal BP, because the ages of the organic layers below turned out to have the same geochronological age of the entisol found at 5.50 m of depth, and probably evolved during a temporary stasis of sedimentation. It is worth noting that, in the core BLG1, the base of the LGM alluvial deposits was found at approximately 29 m of depth, whereas, only a few meters deeper, the sediments already had an age close to the limit of radiocarbon methods.

These new geochronological data (Table 4) suggest that the thermo-luminescence (TL) ages obtained by Marocco [23] and ranging between 76 and 127 ka ago are not reliable, especially in the light of the chronostratigraphic framework available for this sector of the Friuli Plain [28]. Thus, the hypothesis of important structural activity in the area, causing the uplift of tectonic terraces where the dune ridges are found, can be discarded.

East of the studied dune’s field (Figure 1), the surface geological units consisted of recent lagoon deposits and of thick Medieval floods related to an ancient channel belt of Isonzo, corresponding to the present Isonzato River [3,23]. In this sector, the stratigraphic data are few, but east of the Tiel Canal, near Domine zone and east of the site of Viola, some new hand boreholes and few stratigraphic logs stored in the archives of the Consorzio di Bonifica Pianura Friulana indicate that the LGM alluvial deposits were eroded; moreover, in the subsoil, a sequence of Holocene lagoon and fluvial sediments was present up to a depth between 10 and 15 m, where some meters of gravelly sandy deposits were found [57]. These data support the idea that a fluvial incised valley had been formed east of the dune’s field since the end of LGM by the Isonzo River and was later progressively filled during the Holocene.

4.3. Ground-Penetrating Radar

GPR surveys were carried out on the dunes which appeared to be the most preserved ones, owing the highest potential to preserve internal structures and able to provide an efficient geophysical response, because of the occurrence of well-drained sands.

The investigated dunes are located in the area of Belvedere Pineta and San Marco, which are still partly covered by woods. The most noteworthy obstacle we faced was the presence of the tree roots, which appear as parabolic signals in the radargrams (Figure 6 and Figure 7). Fourteen profiles were acquired in Belvedere Pineta and six in San Marco, with a total length of 2 km in both sites combined. In each area, we planned lines transverse and parallel to the dunes’ ridges. In all radargrams, we observed a clear clinostratification in the first 5 m below the surface, displaying various apparent inclinations, with a maximum value of 30° (Figure 14). The dip direction of the layers also varied, with profiles BEL1 and MARC5, which are lines transverse to the north ridge of the respective dunes with a NW dip, while BEL9, parallel to the frontal lobe, has a south-westerly dip direction. Anyhow, as the pseudo-sections were not migrated for the velocity values of the layers, it was not possible to detect, in detail, the actual dip value of the layers. Similar clinostratifications have already been observed in aeolian dunes radargrams documented in the literature, e.g., [58,59].

5. Discussion

Through the integration of remote sensing, field surveys, stratigraphic data and geophysical investigations, we obtained a detailed map depicting the planar morphology of the Belvedere–San Marco reliefs and the reconstruction of their internal structure (Figure 15). Considering their grain size composition, the presence of oriented clinostratification and the general pattern consisting of a complex system of linear ridges joining westwards in frontal lobes with a width ranging between 50 and 250 m, it is reasonable to infer an aeolian origin for the studied dunes. In the literature, this pattern is typical of the parabolic dunes, which are usually made of almost parallel linear and narrow ridges ending in a frontal lobe in the downwind side [60,61,62,63]; in between the ridges, a depressed zone of deflation or palustrine environment can often be found. More specifically, on the basis of location and morphological characteristics, we can distinguish five sectors in the dunes system (Figure 16): Beligna, Muson, San Marco, Viola and Belvedere.

The isles of Volperassa, Volpera, Villanova, Gorgo, Montaron and San Giuliano, which are currently part of the Grado Lagoon, are also part of the dune field; however, because of their patchy occurrence, they have not been associated with a specific morphological unit.

The Beligna unit consists of the northernmost dunes, which are those with the highest length to width ratio between the two stoss-side ridges, ranging from five to seven. In the literature, parabolic dunes with a ratio higher than three are classified as elongated parabolic dunes, also called hairpin dunes [60]. Generally, their main characteristic is to have migrated downwind during their formation for a rather significant distance. The Muson unit, which is almost completely included inside the Due Leoni estate, is the one characterized by the largest extension. The samples of cemented sands analyzed in this research were collected from this unit, and, also, the stratigraphic section presented in the results paragraph was described in this area. The length to width ratio is three, but the wide and composite shape of the frontal lobes gives them a more complex morphology. Similar characteristics are also displayed by the San Marco unit, located along the margin of the lagoon and clearly separated by the previous sector; the best-preserved ridge of this unit is the one reaching the church of San Marco, which is characterized by a linear development of over 2.5 km.

The Viola unit is the easternmost sector of the dune system and one of the less preserved, because of the recent anthropogenic reworking, but also due to the natural erosive processes related to Holocene landscape evolution.

Currently, the area is crossed by the Tiel Canal, and we observe various traces of past lagoon channels and palustrine creeks, which have partially eroded the dune field, reworking and levelling the original morphology. The Belvedere unit is the southernmost one and consists of two sectors corresponding to Centenera Isle and Domine Isle. As suggested by their names and historical maps, these areas were islands surrounded by shallow lagoons until modern reclamation. It is likely that their isolated location might be one of the reasons for their exceptional preservation. In this sector the dunes have a very specific en-echelon morphology, which is made by side-by-side frontal lobes with a typically low length to width ratio. This pattern has been associated with vegetation–sand interactions [64], as plants can modify wind direction, but also hinder the migration of the dunes. In some studies, en-echelon morphologies have also been correlated to modern and past steppe-like environments [14,65].

Considering the specific planforms presented by the dunes, it is possible to infer that the Belvedere–San Marco system is made of different sectors, which have been evidenced in Figure 16. The differences are highlighted by changes in length, width and shape of ridges and they have been separated in five areas, from Beligna (in the north) to Muson, San Marco and Belvedere (in the south). So, the spatial variation seems to be defined as a north–south variation, ranging from elongated to complex and to en-echelon dunes. This trend can be explained by changes in different variables like sediment supply, migration, wind intensity, degree of vegetation cover, soil moisture and chronology. The eastern sectors, showing the lowest amount of migration, might have been the last ones to be deposited, given the presence of other older dunes, on the west side, that would have inhibited their further migration. The morphology of a specific dune can be directly or indirectly influenced by climatic-related factors and can thus provide important information about climate conditions at the time of its formation. As an example, parabolic dunes have been repeatedly associated with the presence of vegetation [60,61].

A study on the White Sands dunes in New Mexico [64] observed that, when a barchan dune (i.e., a type of dune typical of unidirectional winds) migrates in a vegetated area, it becomes a parabolic dune. In other cases [65], dunes formed in a semi-vegetated environment develop a parabolic morphology. Thus, the reconstruction of the planar morphology of a system of paleodunes can reveal important details about the type of vegetation at the time of deposition. Specifically, the parabolic typology and the en-echelon shape of the dunes near Belvedere reveal that the study area was likely characterized by steppe-like vegetation during the formation and migration of the dunes. This observation highlights the key role that remote sensing can have, not only in mapping a specific morphology, but also in indirectly revealing some paleoclimatic aspects.

Nevertheless, it is important to highlight that the severe impact that anthropogenic activities had on the study area, and the consequent degradation of the elevation and the planform of the sand ridges, makes the reconstruction of the dune field prone to some bias, when using only remote sensing.

The photointerpretation of aerial images was mainly based on soilmarks for the detection of the dune ridge boundaries. However, the cropmarks supplied some interesting information only in very limited zones, where the differential growth of the plants highlighted the traces related to the decimetric clinoform layers forming the dunes. This characteristic was also well documented by the stratigraphic sections and new radargrams, showing that the internal structure of the reliefs is characterized by pervasive clinostratification, from a small (<1 cm) to a large scale (>10 cm). The apparent inclination is 20–30°, and the grain size of the laminae varies from silty fine sands to fine gravels, with clasts up to 4–6 mm. We infer that coarser grain-sized laminae might be associated with moments during which the wind was stronger, while finer grain-sized laminae, rich in silt, would be associated with low energy phases, when only dust could be accumulated, in the form of loess. The clinostratification observed in the radargrams also confirms that a similar internal structure is present until a significant depth, at least 5 m below the surface. As documented by Figure 12, the sedimentary structures detected by GPR are very sounding with the deeper stratigraphic sections described by previous authors [21,22] and with the stratigraphic cores carried out over the dune deposits. Moreover, the depth of the aeolian sand bodies has also been confirmed by the analyses of surface waves velocities, Electrical Resistivity Tomography profiles (ERT) and geotechnical tests carried out over the dunes near San Marco and the mouth of the Natissa River, where these geophysical surveys indicated that the sands forming the dunes are documented until 7–8 m of depth and they are partly cemented for the first 2–5 m [66].

Besides the well-evoluted pedological profiles described by Comel [21], which already suggested the old age of the surface of the dunes compared with the Holocene coastal dunes of the northern Adriatic, the new geochronological data presented in this work allowed us to clearly assess the age of formation of some of the sand ridges. In particular, the radiocarbon dates in the core BLG1 document that the ages at the base and at the top of the dune gave almost identical results, testifying that the northernmost dunes of the Beligna unit accumulated between 22.1 and 21.4 ka cal BP; thus, in the final part of the LGM. In fact, the aeolian sands sedimented over the alluvial plain, which experienced its last aggradation after a value contained in the interval 21.8.–21.2 ka cal BP, but most parts of the dunes already existed before an age contained in the time span between 21.8 and 21.4 ka cal BP. The chronological resolution of the radiocarbon is too coarse for solving the real interval in which aeolian deposition occurred but, according to modern analogues [59,61,63], it seems a consistent hypothesis that the dune sampled by the core BLG1 formed in the span of few years/decades or maximum of few centuries.

This information directly constrains the age of a limited sector of the investigated field of dunes, but it is likely that, also, the other aeolian ridges mapped in the area have a rather comparable age. In particular, in the final part of the LGM, the climate was still fully glacial in terms of temperature, atmospheric circulation (ref. [30] and references therein, ref. [67]) and for the related scarce vegetation cover [68]. These conditions did not last after the first part of the Late Glacial, as conifer woods are documented by paleobotanic data in the distal Friulian Plain at 15 ka cal BP (ref. [32] and references therein, refs. [69,70]), and their presence would have hampered the aeolian deflation. It is worth noting that river-bordering dunes have also been recognized near Trieste, formed between 22 and 17.5 ka cal BP [2]. This interval fits nicely with the final part of the LGM, when, in NE Italy, the glaciers started to withdraw and also triggered an important change in the sedimentary dynamic of the related alluvial megafans [29].

The plan pattern of the dunes, combined with their sedimentological and stratigraphic characteristics, clearly points to their aeolian origin and allows us to exclude the hypothesis that they were formed by wind reworking former fluvial ridges with parallel direction, as hypothesized originally by Lenardon and Marocco [22] and also recently recalled by Fontana et al. [3]. In particular, in the study area, the LGM and early Holocene channel belts of Isonzo and Torre streams display a completely different pattern, characterized by divergent downstream directions originated from the apex of their megafans and ranging between ESE and WSW [3,23]. The direction of the reconstructed dunes is thus almost perpendicular to those channel belts and, therefore, not parallel to any of them. Moreover, the parabolic typology of the dunes testify that they migrated from ENE to WSW, implying that the source of the deflated sediment was located eastward from them and the wind was blowing from a direction of about 60°E, parallel to the dunes of Belvedere–San Marco. This coincides with the current Bora wind, which today blows over the Trieste Gulf from the ENE sector and can reach peak velocities of 30–40 m/s [1], suggesting that a similar wind, but probably stronger, persisted during the LGM.

Considering the volume of sands forming the dune’s field, which extends in continuity for over 15 km² and forms ridges with a thickness of up to 15 m, an important source of deflation was necessary, with high sediment availability, at a close distance and mainly in the form of loose sand. In continental aeolian systems, these characteristics are usually satisfied by near fluvial or lacustrine environments during low-water stages, e.g., [13,61,71]. Our interpretation is that the dunes were formed by the Bora wind carrying the sands from the braided channel of Isonzo, flowing just east of the area, along an incised fluvial valley. According to preliminary stratigraphic data [57], east of the zone of Comugna, this feature was at least 1000 m wide and with gravelly sandy braided channels at a depth between 10 and 15 m. It is likely that this channel belt entrenchment formed during the LGM peak and led that surface to be free of floods and, thus, suitable for the accumulation of aeolian sands blown from the incised valley. This geomorphological configuration would supply enough loose sand to feed the aeolian processes, and even the fine sandy gravels that form some of the clinoform layers characterizing the dunes of Belvedere–San Marco. It is also possible that the incised valley was rather shallower at the end of the LGM, when it started to form, while it further deepened during the Late Glacial, as demonstrated for other fluvial incisions of NE Italy cf. [32].

The link between the sand ridges and a local source of sediments corresponding to an ancient paleochannel of the Isonzo River has already been proposed by Lenardon and Marocco [22], on the basis of petrographic investigations of heavy minerals. In particular, they analyzed sand samples collected from three sectors of the dunes (Beligna, San Marco and Belvedere) and the present channels of the Torre Torrent and Isonzo rivers. The most noteworthy difference in composition between the sediments of the two streams is the percentage of picotite, which is a chrome-rich spinel. At present, this mineral is abundant in the sediments of Torre (24.6%), whereas it is very low in Isonzo sediments upstream of the confluence of the rivers (5.4%). Sands from the Belvedere–San Marco dunes are also characterized by a low percentage or absence of picotite (3.7% in Beligna, but not present in San Marco and Belvedere). These mineralogic analyses document that the dunes have been formed by sands deflated by a paleochannel of Isonzo in a period in which the Torre Torrent was not its tributary. This paleohydrographic setting differs from the modern one and is consistent with the LGM situation, when Torre and Isonzo were still forming two separated alluvial megafans, cf. [72].

Considering that the presence of a large continental aeolian dune field is a novelty in alluvial plain of NE Italy, we decided to look at other known case studies (Figure 17). In northern Italy, the only dune field known before the studied one is located in Este, in the Veneto Plain, as described by Rizzetto [18]. That dune field seems to be made by barchanoid–transverse dunes; however, those landforms have been almost completely levelled by anthropogenic activity and, thus, the reconstruction of their planar morphology is limited. In that case study, both the age of deposition and the direction of the formative wind are not clear, but a NNE–SSW direction is hypothesized as dominant. Anyhow, it is interesting to highlight the proximity of LGM channel belts of the Adige River, with a WNW–SES direction, which could have provided the source of blown sands.

From a European perspective, wide fields of parabolic dunes are documented in the European Sand Belt, a region spanning from the Netherlands to Germany, Denmark, and Poland, where important systems of continental dunes and loess deposits exist, which also mainly formed during the LGM or the Late Glacial, e.g., refs. [67,73,74]. In particular, detailed geomorphological reconstructions of riverine dunes exist in the Netherlands, where the use of LiDAR, stratigraphic and pollen investigations allowed researchers to characterize the morphological pattern, the age and the climate existing during the formation of dune systems along braided river channels. Several examples of these landforms have been documented along former channel belts of the Rhine River [75,76] and the Scheldt River [13]. In most of these cases, the aeolian deposits accumulated over abandoned fluvial terraces, a few meters (i.e., 2–20 m) above the river channel active during their formation and from where the sands have been deflated.

Figure 17. Comparison with wind directions (orange arrows) of the aeolian landforms between Central Europe and the northern Adriatic in Moravia [17], the Transdanubian area [14], Deliblato [16], Walachia [77] and Este [18]. The basemap is hillshade, derived from SRTM data [8].

Compared to our case study, a rather close area with widespread occurrences of continental aeolian landforms is the Great Hungarian Plain and, in general, the Carpathian Basin. This region is characterized by a dry and windy climate which promotes aeolian activity, even in the present day, e.g., refs. [78,79].

Aeolian deposits cover vast sectors of the plain and less intense agriculture has led to a better preservation of the reliefs in some areas. The largest aeolian landforms were shaped before the first part of the Holocene, with the most spectacular case being the Mega-Yardangs, south of the Transdanubian Range in Hungary. These are wind erosional landforms, formed by winds carving over 100 m deep corridors in the Upper Miocene bedrock [14]. Interesting examples of large dune fields are present in the area of Nyírség in Eastern Hungary [15], and the Deliblato area, in the Vojvodina Plain of Serbia [16]. Both these fields have a Late Pleistocene–Early Holocene age, and the ridges display a typical parabolic morphology. The Nyírség deposits were formed by winds blowing from NW–N, while the Deliblato field was formed by winds from the SE. The large quantity of loose sands necessary for supplying the formation of the aeolian landforms were fed by the channels of important streams, as the Tisza River for the Nyírség and the Danube River for the Deliblato field.

Downstream of the Danube River, other large systems of relict elongated continental dunes are also present in the Walachia Plain, in SW Romania, as documented in the geomorphological map [77]. In particular, an evident parabolic pattern is displayed in the aerial images of the dune ridges above the fluvial terraces bounding the river channel between Calafat and Corabia. Another system is present in Walachia along the eastern edge of the terrace bordering the Jiu River, downstream of Craiova. In these cases, the landforms were active between the LGM and early Holocene and the winds were blowing from the WNW.

A recent study has characterized a complex system of transverse dunes called the Moravian Sahara dune field, in the Czech Republic, that was formed between the LGM and the Late Glacial [17]. In this case, the sediment was supplied by a braided paleochannel system of the Morava River and the sand ridges accumulated over a fluvial terrace 20 m above the river channel; the direction of the winds was N–S for the older dunes and W–E for the youngest ones. According to the authors, the Moravian Sahara dune field is part of the same wind circulation system that affected, during the last glaciation, the Great Hungarian Plain and was mainly sustained by the baric gradient existing from the Eurasian ice sheets towards the medium latitudes.

Considering the strong winds that shaped the aeolian landforms in Central Europe and their connection with the environmental conditions occurring during the LGM and Late Glacial, it is likely that the formation of the Belvedere–San Marco dunes, together with the ones detected near Trieste by Ronchi et al. [2], represents relict traces of the same wind pattern. Today, the Bora has very strong gusts that can reach 30–40 m/s, due to the orographic setting and the katabatic gradient existing from the Karst Plateau to the Gulf of Trieste. Anyhow, the topography of the region has not changed since the LGM and, even if the study area is rather far from the Carpathian Basin, the local orographic situation could have fostered the strength of the winds reaching the Gulf of Trieste during the last glaciation. Thus, the parabolic dunes existing between Aquileia and Grado could be interpreted as the product of a “glacial Bora”, which was probably blowing more frequently and stronger than today, with the capability to transport even fine gravels. Glacial conditions were also a key factor in their formation because of the necessity of a very limited vegetation cover in favoring the aeolian transportation of the sediment and, moreover, because of the availability, in the distal alluvial plain, of wide braided river channels to supply the deflated sands.

The characterization of the dunes of Belvedere–San Marco testifies that the formation of aeolian landforms was not negligible during the LGM and probably part of the Late Glacial in the alluvial plain of NE Italy and, more generally, in northern Italy. Our research demonstrates that the combined use of remote sensing methodologies and geophysical surveys, accompanied by some fieldwork and stratigraphic data, could support a detailed reconstruction of the pattern of the dune ridges, helping in the assessment of the overall conditions that led to the aeolian processes and the deposition of the related landforms. We are confident that the results of this work can be applied also in other areas of northern Italy and lead to the recognition of other dune fields, maybe of smaller dimensions. Similar geomorphological features were also present in the alluvial plain occupying the floor of the Adriatic Sea before the post-LGM marine transgression but, probably, they have been largely dismantled and partly recycled during the sea-level rise that led to the present coastline.

6. Conclusions

This work is focused on the Belvedere–San Marco sandy reliefs. The use of complementary and various techniques, ranging from high-resolution aerial photos, satellite images, DTM, historical cartography and geophysical and field surveys allowed us to obtain a comprehensive understanding of the study area, to decipher the dunes origin and their chronology.

The Belvedere–San Marco dunes extend for over 15 km² and they represent the largest continental aeolian system ever identified in the Italian Peninsula. Their planar morphology is parabolic, with ridges up to 3 km long, joining in frontal lobes in the west. The internal structure is characterized by a clinostratification, and the total thickness of the dunes can reach 15 m in some ridges;
The deposition of the dunes started at the end of the LGM (22.1–21.4 ka cal BP) and possibly continued during the first part of the Late Glacial. The radiocarbon samples indicate that the depositional events were short-lived;
Chronology and parabolic shapes suggest that, during the deposition, the climate was characterized by low temperatures, less precipitation, and sparse, steppe-like vegetation. The orientation of the dunes permits us to link them with the Bora wind, an ENE katabatic wind, and one of the strongest winds in the Italian Peninsula;
The sandy material of the reliefs has a local origin, probably supplied by a braided paleochannel, located east of the study area;
We observed many similarities with other continental aeolian dunes located in Central and Northern Europe, which were partially deposited in the same time period or in periglacial conditions;
This work opens new perspectives on the role and presence of other large aeolian landforms in Northern Italy.

Supplementary Materials

The following supporting information can be downloaded at: https://www.mdpi.com/article/10.3390/rs16050765/s1. Suppl1_DEM_Dunes-Belvedere-SanMarco.kmz; Suppl2_3D-Concretion_sample_DueLeoni-2.pdf.

Author Contributions

Conceptualization, A.F. and F.V.; data analysis, F.V. and A.F.; field surveys, F.V., A.F., L.R. and J.B.; writing—original draft preparation, F.V.; writing—reviewing and editing, F.V., A.F., V.C., L.R., J.B. and I.H.; funding acquisition, A.F. All authors have read and agreed to the published version of the manuscript.

Funding

This research was funded by the University of Padova, project DOR2023139/20.

Data Availability Statement

All relevant data are available in the main text or in the Supplementary Materials.

Acknowledgments

The authors thank the owners of the Due Leoni farmland, the owners of the San Marco croplands and the owners and manager of the Belvedere Pineta Camp Ground, for their collaboration in letting us carry out field surveys and GPR surveys. We also wish to thank Consorzio di Bonifica Pianura Friulana for making available their database of stratigraphic logs. We thank Planet “Education and Research Program” for providing free access to the imagery data used in this work. Finally, we thank Cristina Stefani and Giovanni Monegato for the analysis of the thin sections and the discussion about the petrography of the samples, as well as Fabrizio Antonioli and Giorgio Fontolan for the availability of the core BLG. This paper is dedicated to the memory of Ruggero Marocco, who devoted significant efforts to improve the knowledge of the study area.

Conflicts of Interest

The authors declare no conflicts of interest. The funders had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript; or in the decision to publish the results.

References

OSMER. Il Clima del Friuli Venezia Giulia; Regione Friuli Venezia Giulia: Trieste, Italy, 2016; 27p. [Google Scholar]
Ronchi, L.; Fontana, A.; Novak, A.; Correggiari, A.; Poglajen, S. Late-Quaternary Evolution of the Semi-Confined Alluvial Megafan of Isonzo River (Northern Adriatic): Where the Fluvial System of the Southern Alps Meets the Karst. Geosciences 2023, 13, 135. [Google Scholar] [CrossRef]
Fontana, A.; Monegato, G.; Rossato, S.; Poli, M.E.; Furlani, S.; Stefani, C. Carta Delle Unità Geologiche della Pianura del Friuli Venezia Giulia alla Scale 1:150,000 e Note Illustrative; Regione Autonoma Friuli Venezia Giulia—Servizio Geologico: Trieste, Italy, 2019; 80p. [Google Scholar]
D’Orefice, M.; Graciotti, R. Geomorphological and sedimentological features on the aeolian deposits in the western sector of the Mt. Calamita promontory (Elba Island, Italy). Mem. Descr. Carta Geol. d’It. 2008, 78, 113–126. [Google Scholar]
Di Maggio, C.; Agate, M.; Contino, A.; Basilone, L.; Catalano, R. Unità a limiti inconformi utilizzate per la cartografia dei depositi quaternari nei fogli CARG della Sicilia nord-occidentale. Il Quat. Ital. J. Quat. Sci. 2009, 22, 345–364. [Google Scholar]
Ferranti, L.; Burrato, P.; Sechi, D.; Andreucci, S.; Pepe, F.; Pascucci, V. Late Quaternary coastal uplift of southwestern Sicily, central Mediterranean sea. Qua. Sci. Rev. 2021, 255, 106812. [Google Scholar] [CrossRef]
Melis, R.T.; Di Gregorio, F.; Panizza, V. The Coastal Dunes of Sardinia; Landscape Response to Climate and Sea Level Changes. In Landscapes and Landforms of Italy; Soldati, M., Marchetti, M., Eds.; Springer: Cham, Switzerland, 2017; pp. 365–376. [Google Scholar]
SRTM 90m DEM Digital Elevation Database. Available online: https://srtm.csi.cgiar.org (accessed on 21 December 2023).
Cremaschi, M. The LOESS in Northern and Central Italy: A Loess Basin between the Alps and the Mediterranean Region; Excursion Guidebook, September–October 1988; Editrice Gutemberg: Milano, Italy, 1990; 187p. [Google Scholar]
Peresani, M.; Monegato, G.; Ravazzi, C.; Bertola, S.; Margaritora, D.; Breda, M.; Fontana, A.; Fontana, F.; Janković, I.; Karavanić, I.; et al. Hunter-gatherers across the great Adriatic-Po region during the Last Glacial Maximum: Environmental and cultural dynamics. Quat. Int. 2021, 581–582, 28–163. [Google Scholar] [CrossRef]
Koster, E.A. The “European Aeolian Sand Belt”: Geoconservation of Drift Sand Landscapes. Geoheritage 2009, 1, 93–110. [Google Scholar] [CrossRef]
Pierik, H.J.; van Lanen, R.J.; Gouw-Bouman, M.T.; Groenewoudt, B.J.; Wallinga, J.; Hoek, W.Z. Controls on late-Holocene drift-sand dynamics: The dominant role of human pressure in the Netherlands. Holocene 2018, 28, 1361–1381. [Google Scholar] [CrossRef] [PubMed]
Kasse, C.; Woolderink, H.A.G.; Kloos, M.E.; Hoek, W.Z. Source-bordering aeolian dune formation along the Scheldt River (southern Netherlands—Northern Belgium) was caused by Younger Dryas cooling, high river gradient and southwesterly summer winds. Neth. J. Geosci. 2020, 99, e13. [Google Scholar] [CrossRef]
Sebe, K.; Csillag, G.; Ruszkiczay-Rüdiger, Z.; Fodor, L.; Thamó-Bozsó, E.; Müller, P.; Braucher, R. Wind erosion under cold climate: A Pleistocene periglacial megayardang system in Central Eurpe (Western Pannonian Basin, Hungary). Geomorphology 2011, 134, 470–482. [Google Scholar] [CrossRef]
Kiss, T.; Lóki, J. Parabolic dunes of the Southern Nyírség. Landscapes and landforms of Hungary. Parabolic dunes of the Southern Nyírség. In Landscapes and landforms of Hungary; Lóczy, D., Ed.; Springer: Cham, Switzerland, 2015; pp. 255–261. [Google Scholar] [CrossRef]
Sipos, G.; Marković, S.B.; Gavrilov, M.B.; Balla, A.; Filyó, D.; Bartyik, T.; Mészáros, M.; Tóth, O.; van Leeuwen, B.; Lukić, T.; et al. Late Pleistocene and Holocene aeolian activity in the Deliblato Sands, Serbia. Quat. Res. 2021, 107, 113–124. [Google Scholar] [CrossRef]
Holuša, J.; Nývlt, D.; Woronko, B.; Matějka, M.; Stuchlík, R. Environmental factors controlling the Last Glacial multi-phase development of the Moravian Sahara dune field, Lower Moravian Basin, Central Europe. Geomorphology 2022, 413, 108355. [Google Scholar] [CrossRef]
Rizzetto, F.; Mycielska Dowgiallo, E.; Castiglioni, G.B. Some aeolian features in the Po plain near Este (North Italy). Geogr. Fis. Dinam. Quat. 1998, 21, 245–253. [Google Scholar]
Taramelli, T. Dei terreni morenici e alluvionali del Friuli. Ann. Reg. Ist. Tecn. Udine 1875, 8, 1–91. [Google Scholar]
Feruglio, E. La Zona Delle Risorgive Del Basso Friuli Fra Tagliamento e Torre: Studio Geologico, Idrologico e Agronomico; Stab. Tip. Friulano: Udine, Italy, 1925; 343p. [Google Scholar]
Comel, A. Ricerche pedologiche sulle dune di Belvedere (Aquileia) e di Grado. N. Ann. Ist. Chim. Agr. Sperim. Gorizia 1951, 2, 18. [Google Scholar]
Lenardon, G.; Marocco, R. Le dune di Belvedere-San Marco. Una antica linea di riva? (2) Considerazioni sedimentologiche. Atti Museo Friul. Storia Nat. 1995, 16, 5–24. [Google Scholar]
Marocco, R. Prima ricostruzione paleo-idrografica del territorio della Bassa Pianura Friulano-Isontina e della Laguna di Grado nell’Olocene. Gortania Geol. Paleontol. Paletnol. 2010, 31, 69–86. [Google Scholar]
Marocco, R. Le dune di Belvedere-San Marco. Una antica linea di riva? (1) Considerazioni geomorfologiche. Atti Museo Friul. Storia Nat. 1992, 13, 57–76. [Google Scholar]
Arnaud-Fassetta, G.; Carre, M.B.; Marocco, R.; Maselli Scotti, F.; Pugliese, N.; Zaccaria, C.; Bandelli, A.; Bresson, V.; Manzoni, G.; Montenegro, M.E.; et al. The site of Aquileia (Northeastern Italy): Example of fluvial geoarchaeology in a Mediterranean deltaic plain. Géomorphologie 2003, 4, 227–246. [Google Scholar] [CrossRef]
Siché, I.; Arnaud-Fassetta, G. Anthropogenic activities since the end of the Little Ice Age: A critical factor driving fluvial changes on the Isonzo River (Italy, Slovenia). Méditerranée 2014, 122, 183–199. [Google Scholar] [CrossRef]
Cimolino, A.; Della Vedova, B.; Nicolich, R.; Barison, E.; Brancatelli, G. New evidence of the outer Dinaric deformation front in the Grado area (NE-Italy). Rend. Fis. Acc. Lincei 2010, 21, 167–179. [Google Scholar] [CrossRef]
Fontana, A.; Mozzi, P.; Bondesan, A. Late Pleistocene evolution of the Venetian-Friulian Plain. Rend. Lincei Sci. Fis. Nat. 2010, 21, 181–196. [Google Scholar] [CrossRef]
Fontana, A.; Monegato, G.; Devoto, S.; Zavagno, E.; Burla, I.; Cucchi, F. Evolution of an Alpine fluvioglacial system at the LGM decay: The Cormor type megafan (NE Italy). Geomorphology 2014, 204, 136–153. [Google Scholar] [CrossRef]
Dal Gobbo, C.; Colucci, R.R.; Monegato, G.; Žebre, M.; Giorgi, F. Atmosphere-cryosphere interactions during the last phase of the Last Glacial Maximum (21 ka) in the European Alps. Clim. Past 2023, 19, 1805–1823. [Google Scholar] [CrossRef]
Bavec, M.; Tulaczyk, S.M.; Mahan, S.A.; Stock, G.M. Late Quaternary glaciation of the Upper Soča River Region (Southern Julian Alps, NW Slovenia). Sediment. Geol. 2004, 165, 265–283. [Google Scholar] [CrossRef]
Ronchi, L.; Fontana, A.; Cohen, K.M.; Stouthamer, E. Late Quaternary landscape evolution of the buried incised valley of Concordia Sagittaria (Tagliamento River, NE Italy): A reconstruction of incision and transgression. Geomorphology 2021, 373, 107509. [Google Scholar] [CrossRef]
Marocco, R. Evoluzione Quaternaria della Laguna di Marano (Friuli-Venezia Giulia). Quaternario 1989, 2, 125–137. [Google Scholar]
Marocco, R. Evoluzione Tardopleistocenica-Olocenica del Delta del F. Tagliamento e delle Lagune di Marano e Grado (Golfo di Trieste) Quaternario 1991, 4, 223–232. [Google Scholar]
Amorosi, A.; Fontana, A.; Antonioli, F.; Primon, S.; Bondesan, A. Post-LGM sedimentation and Holocene shoreline evolution in the NW Adriatic coastal area. GeoActa 2008, 7, 41–67. [Google Scholar]
Vacchi, M.; Marriner, N.; Morhange, C.; Spada, G.; Fontana, A.; Rovere, A. Multiproxy assessment of Holocene relative sea-level changes in the western Mediterranean: Sea-level variability and improvements in the definition of the isostatic signal. Earth Sci. Rev. 2016, 155, 172–197. [Google Scholar] [CrossRef]
Fontana, A.; Vinci, G.; Tasca, G.; Mozzi, P.; Vacchi, M.; Bivi, G.; Salvador, S.; Rossato, S.; Antonioli, F.; Asioli, A.; et al. Lagoonal settlements and relative sea level during Bronze Age in Northern Adriatic: Geoarchaeological evidence and paleogeographic constraints. Quat. Int. 2017, 439, 17–36. [Google Scholar] [CrossRef]
Cerrone, C.; Vacchi, M.; Fontana, A.; Rovere, A. Last Interglacial sea-level proxies in the western Mediterranean. Earth Syst. Sci. Data 2021, 13, 4485–4527. [Google Scholar] [CrossRef]
Geoprodotti—IGM E-Commerce Site. Available online: https://www.igmi.org/geoprodotti (accessed on 6 October 2023).
Planet Team. Planet Application Program Interface: In Space for Life on Earth. San Francisco, CA, USA. 2017. Available online: https://api.planet.com (accessed on 6 October 2023).
Tukey, J.W. Exploratory Data Analysis; Addison-Wesley Publishing Compagny: Reading, UK, 1977; 677p. [Google Scholar]
Shennan, S. Quantifying Archaeology; Edimburgh University Press: Edimburgh, UK, 1988; 364p. [Google Scholar]
Eagle FVG. Available online: https://eaglefvg.regione.fvg.it (accessed on 6 October 2023).
Zakšek, K.; Oštir, K.; Kokalj, Ž. Sky-View Factor as a Relief Visualization Technique. Remote Sens. 2011, 3, 398–415. [Google Scholar] [CrossRef]
Kokalj, Ž.; Somrak, M. Why Not a Single Image? Combining Visualizations to Facilitate Fieldwork and On-Screen Mapping. Remote Sens. 2019, 11, 747. [Google Scholar] [CrossRef]
Hesse, R. Visualisierung hochauflösender Digitaler Geländemodelle mit LiVT. In Computeranwendungen und Quantitative Methoden in der Archäologie. 4. Workshop der AG CAA 2013; Lieberwirth, U., Herzog, I., Eds.; Edition Topoi: Berlin, Germany, 2016; pp. 109–128. [Google Scholar]
Kokalj, Ž.; Hesse, R. Airborne Laser Scanning Raster Data Visualization: A Guide to Good Practice; Založba ZRC: Ljubljana, Slovenia, 2017; 88p. [Google Scholar]
Orengo, H.A.; Petrie, C.A. Multi-scale relief model (MSRM): A new algorithm for the visualization of subtle topographic change of variable size in digital elevation models. Earth Surf. Process. Landf. 2018, 43, 1361–1369. [Google Scholar] [CrossRef] [PubMed]
Abdi, H.; Williams, L.J. Principal Component Analysis. Wiley Interdiscip. Rev. Comput. Stat. 2010, 2, 433–459. [Google Scholar] [CrossRef]
Arcanum Maps—The Historical Map Portal. Available online: https://maps.arcanum.com (accessed on 30 October 2023).
Rossi, M. Kriegskarte, 1798–1805: Il Ducato di Venezia Nella Carta di Anton von Zach; Fondazione Benetton Studi Ricerche-Grafiche V. Bernardi: Treviso, Italy, 2005; 788p. [Google Scholar]
Busetti, M.; Baradello, L.; Caburlotto, A.; Gordini, E.; Zanolla, C.; Accettella, D.; Antonioli, F.; Paganini, P. Indagini geofisiche per lo studio dell’evoluzione Tardo-Pleistocenica ed Olocenica della Laguna di Grado e Marano (Adriatico Settentrionale). In Proceedings of the Workshop finale Project VECTOR (Vulnerability of the Italian Coastal Area and Marine Ecosystems to CLIMATIC Changes and Their Role in the Mediterranean Carbon Cycles), Rome, Italy, 18–19 October 2010. [Google Scholar]
Reimer, P.; Austin, W.E.N.; Bard, E.; Bayliss, A.; Blackwell, P.G.; Bronk Ramsey, C.; Butzin, M.; Edwards, R.L.; Friedrich, M.; Grootes, P.M.; et al. The IntCal20 Northern Hemisphere radiocarbon age calibration curve (0–55 cal kB). Radiocarbon 2020, 62, 725–757. [Google Scholar] [CrossRef]
Daniels, D.J. Surface-Penetrating Radar; Institut. Electric. Engineers: London, UK, 1996. [Google Scholar]
Francke, J. A review of selected ground penetrating radar applications to mineral resource evaluations. J. Appl. Geophy. 2012, 81, 29–37. [Google Scholar] [CrossRef]
Michelutti, G.; Zanolla, S.; Barbieri, S. Suoli e Paesaggi del Friuli Venezia Giulia: 1. Pianura e Colline del Pordenonese. Note Illustrative; Agenzia regionale per lo sviluppo rurale—Servizio della sperimentazione agraria—Ufficio del suolo: Udine, Italy, 2003; p. 512. [Google Scholar]
Bagalini, M. Geomorphological and Stratigraphical Evolution of the Post-LGM Incised Valleys of the Eastern Low Friulian Plain (NE Italy). Master’s Thesis, University of Padova, Padova, Italy, 2022. [Google Scholar]
Schenk, C.D.; Gautier, D.; Olhoeft, G.R.; Lucius, J.E. Internal Structure of an Aeolian Dune using Ground-Penetrating Radar. In Aeolian Sediments: Ancient and Modern; Pye, K., Lancaster, N., Eds.; Blackwell Scientific Publications: Oxford, UK, 2009; pp. 61–71. [Google Scholar] [CrossRef]
Wolfe, S.A.; Moorman, B.J.; Hugenholtz, C.H. Effects of sand supply on the morphodynamics and stratigraphy of active parabolic dunes, Bigstick Sand Hills, southwestern Saskatchewan. Can. J. Earth Sci. 2008, 45, 321–335. [Google Scholar] [CrossRef]
Pye, K. Morphological Development of coastal dunes in a humid tropical environment, Cape Bedford and Cape Flattery, North Queensland. Geogr. Ann. 1982, A64, 213–227. [Google Scholar] [CrossRef]
Pye, K.; Lancaster, N. Aeolian Sediment, Ancient and Modern; Blackwell Scientific Publications: Oxford, UK, 1993; 169p. [Google Scholar]
Livingstone, I.; Warren, A. Aeolian Geomorphology: A New Introduction; Wiley-Blackwell: Oxford, UK, 2019; 318p. [Google Scholar]
Yan, N.; Baas, A.C.W. Parabolic dunes and their transformations under environmental and climatic changes: Towards a conceptual framework for understanding and prediction. Glob. Planet. Chang. 2015, 124, 123–148. [Google Scholar] [CrossRef]
Durán, O.; Herrmann, H.J. Vegetation Against Dune Mobility. Phys. Rev. Lett. 2006, 97, 188001. [Google Scholar] [CrossRef]
Filion, L.; Morriset, P. Eolian landforms along the eastern coast of Hudson Bay, Northern Québec. Collect. Nord. 1983, 47, 73–94. [Google Scholar]
Dal Moro, G.; Stemberk, J. Tools for the efficient analysis of surface waves from active and passive seismic data: Exploring an NE-Italy perilagoon area with significant lateral variations. Earth Planets Space 2022, 74, 140. [Google Scholar] [CrossRef]
Schaffernicht, E.J.; Ludwig, P.; Shao, Y. Linkage between dust cycle and loess of the last Glacial Maximum in Europe. Atmos. Chem. Phys. 2020, 20, 4969–4986. [Google Scholar] [CrossRef]
Miola, A.; Bondesan, A.; Corain, L.; Favaretto, S.; Mozzi, P.; Piovan, S.; Sostizzo, I. Wetlands in the Venetian Po Plain (northeastern Italy) during the Last Glacial Maximum: Interplay between vegetation, hydrology and sedimentary environment. Rev. Palaeobot. Palynol. 2006, 141, 53–81. [Google Scholar] [CrossRef]
Donegana, M.; Fontana, A.; Paiero, G.; Ravazzi, C. Le torbe tardiglaciali di Bannia—Palazzine di Sopra: Palinologia e ricostruzione ambientale. Quad. Del Mus. Archaeol. Del Friuli Occident. 2005, 5, 15–18. [Google Scholar]
Vescovi, E.; Ravazzi, C.; Arpenti, E.; Finsinger, W.; Pini, R.; Valsecchi, V.; Wick, L.; Ammann, B.; Tinner, W. Interactions between climate and vegetation during the Lateglacial period as recorded by lake and mire sediment archives in Northern Italy and Southern Switzerland. Quat. Sci. Rev. 2007, 26, 1650–1669. [Google Scholar] [CrossRef]
Vandenberghe, J. Changing conditions of aeolian sand deposition during the last deglaciation period. Z. Für Geomorphol. 1991, 90, 193–207. [Google Scholar]
Fontana, A.; Mozzi, P.; Marchetti, M. Alluvial fans and megafans along the southern side of the Alps. Sediment. Geol. 2014, 301, 150–171. [Google Scholar] [CrossRef]
Haase, D.; Fink, J.; Haase, G.; Ruske, R.; Pećsi, M.; Richter, H.; Altermann, M.; Jäger, K.D. Loess in Europe—Its spatial distribution based on a European Loess Map, scale 1:2,500,000. Quat. Sci. Rev. 2007, 26, 1301–1312. [Google Scholar] [CrossRef]
Lehmkuhl, F.; Nett, J.J.; Pötter, S.; Schulte, P.; Sprafke, T.; Jary, Z.; Antoine, P.; Wacha, L.; Wolf, D.; Zerboni, A.; et al. Loess landscapes of Europe—Mapping, geomorphology, and zonal differentiation Loess landscapes of Europe—Mapping, geomorphology, and zonal differentiation. Earth Sci. Rev. 2021, 215, 103496. [Google Scholar] [CrossRef]
Berendsen, H.J.A.; Stouthamer, E. Late Weichselian and Holocene palaeogeography of the Rhine-Meuse delta, The Netherlands. Palaeogeogr. Palaeoclimatol. Palaeoecol. 2000, 161, 311–335. [Google Scholar] [CrossRef]
Kasse, C.K. Sandy aeolian deposits and environments and their relation to climate during the Last Glacial Maximum and Lateglacial in northwest and Central Europe. Prog. Phys. Geogr. 2002, 24, 507–532. [Google Scholar] [CrossRef]
Badea, L.; Niculescu, G.; Sencu, V. Atlasul Republicii Socialiste Romania: Harta Geomorfologica, III-1; Institutul de Geographie: București, Romania, 1976. [Google Scholar]
Borsy, Z. Blown sand territories in Hungary. Z. Für Geomorphol. 1991, 90, 1–14. [Google Scholar]
Brezsnyánszky, K.; Síkhegyi, F. Geology. Map at 1:2,350,000 scale. In Hungary in Maps; Kocsis, K., Schweitzer, F., Eds.; Geographical Research Institute, Hungarian Academy of Sciences: Budapest, Hungary, 2009; 35p. [Google Scholar]

Figure 2. DTM of the main study area (a), and San Giuliano Isle (b). Profile a–a’’’ is parallel to the dunes’ ridges, profile b–b’’ is transverse. Name of localities: 1: Beligna; 2: Muson; 3: San Marco; 4: Belvedere; 5: Domine; 6: Viola.

Figure 3. Elaborations of DTMs in the area of Belvedere Pineta Camp Ground (Centenera). (a) DTM and hillshading; (b) slope; (c) simple local relief model (SLRM); (d) local dominance; (e) multi-scale relief model (MSRM); and (f) principal component analysis of previous elaborations in RGB = pc1–pc2–pc3 composite.

Figure 4. Examples of soilmarks and cropmarks evidencing the dunes in visible colors satellite images. (a) Planar morphology of the dunes near Isola Domine in a Planet image; (b) frontal lobes and ridges of dunes in the areas of Beligna and Muson in a satellite image available in Google Earth (Google/Maxar Technology); (c) levelled ridges of dunes in the Muson area in a satellite image from ESRI basemap; (d) satellite image from Google Earth (Google/Maxar Technology) that reveals a dense pattern of parallel traces (as indicated by arrows) forming an angle between 15 and 35°, with respect to the direction of the dune ridge, which corresponds to the intersection of the clinostratified layers of the levelled dunes with the flat ground surface.

Figure 5. Box plot of the spectral signatures of the dunes and the surrounding distal plain (lagoon and alluvial deposits) on satellite imagery of Aster (left) and Sentinel-2 (right). For Aster, we used red (band 2), near infrared (band 3) and short-wave infrared (band 9). For Sentinel-2, we used green (band 3), vegetation red edge (band 6) and short-wave infrared (band 11).

Figure 6. Area of San Marco, an example of integration of different data. (a) Historical cartography with red arrows indicating the slope of the dune ridges; (b) DTM; (c) satellite imagery (Planet: visible bands); (d) reconstruction of the planar morphology of the dunes (in yellow).

Figure 7. Area of Belvedere, an example of integration of different data. (a) Historical cartography with red arrows indicating the scarps of the dunes’ ridges; (b) DTM; (c) satellite imagery (Planet: visible bands); (d) reconstruction of the planar morphology of the dunes (in yellow).

Figure 8. Area of Beligna, an example of integration of different data. (a) Historical cartography with red arrows indicating the slopes of the dunes’ ridges; (b) DTM; (c) aerial photos (1993) (d) reconstruction of the planar morphology of the dunes (in yellow). The red star indicates the location of BLG1 borehole.

Figure 9. (a) Cross-section location in Due Leoni estate. (b) Closer look at the cross-section with visible clinostratified bed. In (c) 1a–1f are the name of the foreset layers which were sampled during the field activity.

Figure 10. Examples of cemented sands (“maltoni”) collected in the plough layer, but formerly formed as carbonate concretions in deeper pedogenetic horizons (Bk/Ck). The presence of a clinostratification is detectable in (a,b). In (b), the rectangle corresponds to the thin sections reported in Figure 11. (c) Sample DueLeoni2, with evident lumpy morphology related to bioturbation features.

Figure 11. Thin sections of sample Muson27 in crossed-polarized light. See Figure 10 for orientation of the sample. (a) The darker layers consist of finer grains; (b) particular of the transition from a layer of fine sands (bottom) to coarser sand and fine gravel (top).

Figure 12. Historical photos of outcrops exposed in sand quarries, where the clinostratification of the cemented sands is evident. (a) Area of Belvedere at the end of the 19th century; in the foreground, a well of the Roman period was also exposed by the excavations (Archive of the National Archaeological Museum of Aquileia); (b) Gorgo Island in the Grado Lagoon [24]; (c) area of Belvedere Centenera [21].

Figure 13. Stratigraphic log of the first 10 m of the core BLG1. Radiocarbon dates are calibrated (cal BP), as reported in Table 4.

Figure 14. (a) BEL1 and (b) MARC5 radargrams with highlighted clinostratification. On the bottom, the location in Belvedere (c) and San Marco (d) of the highlighted radargrams (in red) and the other GPR lines (in white).

Figure 15. Map of the reconstructed planar morphology of the Belvedere–San Marco dunes. The planar morphology was reconstructed by combining historical cartography, topographic data obtained from LiDAR and optical traces visible in aerial photos and satellite images. Examples of the integrated analyses are shown in Figure 6, Figure 7 and Figure 8. (a) Main study area; (b) San Giuliano Isle area.

Figure 16. Map of the continental dunes of Belvedere–San Marco of Aquileia, with indications of the different morphological units and schematic reconstruction of the hypothetical depositional setting that favored their formation. Different colors of the dunes reflect a difference in the parabolic morphology. The source of the sediment blown by the Bora wind (light blue arrows) was likely corresponding to a braided paleochannel of Isonzo River flowing eastward of the dunes’ field and slightly entrenched with respect to the alluvial plain. (a) Main study area; (b) San Giuliano Isle area.

Table 1. Characteristics of aerial photo datasets used in this work.

Year	Frames	Flight Strip	Scale	Format	Format (cm)	Channels	Flight Season
1938	89c–91c	5	1:12,500	Analog	13 × 18	Black/White	April
1945	85	2	1:48,000	Analog	20 × 20	Black/White	September
1954	212,662,663	13A, 12	1:29,000	Analog	23 × 23	Black/White	April–May
1984	225	17	1:32,000	Analog	23 × 23	Black/White	June
1993	1058–1061, 1087–1090	77, 75	1:33,000	Analog	23 × 23	Black/White	March–April

Table 2. Characteristics of satellite datasets used in this work. Spectral resolution are defined as: VNIR (visible and near-infrared), SWIR (short-wave infrared), and TIR (thermal infrared).

Year	Date	Resolution (m)	Bands	Spectral Resolution
Aster	15 March 2003 14 April 2007	15, 30, 90	14	VNIR, SWIR, TIR
Sentinel-2	1 March 2021 24 March 2021	10, 20, 30	13	VNIR, SWIR
Planet	29 April 2017 6 March 2021	3	4	VNIR

Table 3. Characteristics of historical maps used in this work.

Map	Year	Scale	Format	Source
Innerösterreich— 1st Military Survey	1784–1785	1:28,800	Colours	Maps Arcanum
Von Zach	1809	1:28,800	Colours	Rossi
2nd Military Survey	1821–1824	1:28,800	Colours	Maps Arcanum
3rd Military Survey—Austroungarico	1869–1887	1:25,000	Colours	Maps Arcanum
Tavoletta 25V—Grado	1917	1:25,000	Black/White	IGM
Tavoletta 25V—Grado	1927	1:25,000	Black/White	IGM
Tavoletta 25V—Grado	1949	1:25,000	Colours	IGM

Table 4. Radiocarbon ages were measured on samples from the core BLG1.

Depth (m)	Lab Number	Description	Age ¹⁴C BP	Calibrated Age 2σ cal-BP
5.50	ETH-114380	Plant fragment	17,867 ± 88	22,119–21,255
9.73	ETH-55393	Wood	17,912 ± 70	21,915–21,464
9.93	ETH-55394	Wood	17,877 ± 69	21,881–21,427
28.62–28.64	ETH-55384	Peat	22,589 ± 108	27,113–26,729
32.51–32.52	ETH-55385	Peat	28,333 ± 176	32,561–31,868
33.58–33.59	ETH-55388	Peat	43,980 ± 1032	48,394–46,180

Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

© 2024 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https://creativecommons.org/licenses/by/4.0/).

Share and Cite

MDPI and ACS Style

Vanzani, F.; Fontana, A.; Ronchi, L.; Boaga, J.; Chiarini, V.; Hajdas, I. The Dunes of Belvedere–San Marco of Aquileia: Integrating High-Resolution Digital Terrain Models and Multispectral Images with Ground-Penetrating Radar Survey to Map the Largest System of Continental Dunes of Northern Italy. Remote Sens. 2024, 16, 765. https://doi.org/10.3390/rs16050765

AMA Style

Vanzani F, Fontana A, Ronchi L, Boaga J, Chiarini V, Hajdas I. The Dunes of Belvedere–San Marco of Aquileia: Integrating High-Resolution Digital Terrain Models and Multispectral Images with Ground-Penetrating Radar Survey to Map the Largest System of Continental Dunes of Northern Italy. Remote Sensing. 2024; 16(5):765. https://doi.org/10.3390/rs16050765

Chicago/Turabian Style

Vanzani, Federica, Alessandro Fontana, Livio Ronchi, Jacopo Boaga, Veronica Chiarini, and Irka Hajdas. 2024. "The Dunes of Belvedere–San Marco of Aquileia: Integrating High-Resolution Digital Terrain Models and Multispectral Images with Ground-Penetrating Radar Survey to Map the Largest System of Continental Dunes of Northern Italy" Remote Sensing 16, no. 5: 765. https://doi.org/10.3390/rs16050765

APA Style

Vanzani, F., Fontana, A., Ronchi, L., Boaga, J., Chiarini, V., & Hajdas, I. (2024). The Dunes of Belvedere–San Marco of Aquileia: Integrating High-Resolution Digital Terrain Models and Multispectral Images with Ground-Penetrating Radar Survey to Map the Largest System of Continental Dunes of Northern Italy. Remote Sensing, 16(5), 765. https://doi.org/10.3390/rs16050765

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Menu