Ecogeomorphological Transformations of Aeolian Form—The Case of a Parabolic Dune, Poland

: The range of natural environmental degradation caused by anthropogenic activity may include geomorphological forms such as dunes resulting from the build-up activity of the wind. In effect, such environmental transformation affects changes connected not only with their relief, but also with the presence and health of diverse plant and animal inhabitants. The subject of the survey was a parabolic dune with asymmetric shape, the sand of which was subjected to exploitation over many years. Terrain data acquired by means of GNSS (Global Navigation Satellite Systems) served to elaborate the present relief of the surveyed dune and to reconstruct its primary relief. These were mainly places where the impacts of human activities were recorded. For this purpose, ordinary kriging (OK) estimation was performed. Simultaneously, satellite data and UAV (Unmanned Aerial Vehicle) imaging were acquired, and subjected to image fusion in order to acquire near infrared bands (NIR), red, green, blue in high spatial resolution. These in turn were applied so as to estimate the condition of the vegetation overplanting the dune and surrounding terrain. The correctness of the modelling was veriﬁed by cross-validation (CV), which disclosed low error values. Such values in present and primary relief were, respectively, mean error (ME) at − 0.009 and − 0.014, root mean square error (RMSE) at 0.564 and 0.304 and root mean square standardised error (RMSSE) at 0.999 and 1.077. Image fusion, with use of pansharpening allowed a colour-infrared composition (CIR) and a Modiﬁed Chlorophyll Absorption in Reﬂectance Index 1 (MCARI1) to be obtained. Their analysis disclosed that vegetation on the dune is characterised by worse health condition as compared with the surrounding area thereof. The proposed approach enabled the environmental condition of the surveyed dune to be analysed, and thereby it allows for a determination of the consequences of further uncontrolled sand recovery without taking into account the historical cartographic materials customarily considered to be the main source of information.

Increasingly, with regard to mineral raw material grabbing, such operations are often illegal or unmonitored as they take place in underdeveloped or developing countries. Such grabbing includes sand which, together with gravel, constitutes the important natural aggregates used in the building trades. The problem concerns areas where uncontrolled and irrational raw material management contributes to their significant depletion [4,17,18]. However, rational exploitation of these raw materials may have a benefit, e.g., new habitats

Object of the Survey
The object of the survey was a parabolic dune together with its neighboring terrains. The parabolic dune is an accumulative aeolian form with its elongated arms directed to the windward direction that is usually stabilised by vegetation, and is mainly found in areas with a semi-dry climate or along seacoasts [40,49]. The object of this survey was formed in the central part of Poland and is situated in the village of Dudki, in the Kobiele Wielkie commune, Radomsko county, Łódź voivodship. The site of the surveyed geomorphological form arose due to the morpho-shaping action of wind, and is determined by the following geographical coordinates: 50.994875 • N, 19.648821 • E (Figure 1). Over the years, the surveyed dune has been subjected to continuous transformations by the forces of human activity (Figure 2), directed at the exploitation of the Quaternary aeolian sands that had formed on Pleistocene sands, and of the fluvio-glacial gravels from the mid-Poland maximum glaciation stadial. North and northeast from the surveyed object, Holocene peats can be found, as well as peaty silts on the sands and muds of the fluvial floodplain terraces [50]. Available cartographic material analysis, including reference to historical maps (among other items), show [51] that the dune's degradation began most probably about 100 years ago. The recovered sand was and still is utilised mainly for building purposes.  (e,f) portion B of the dune-exhibiting low-level human activity in the last year; (i,j) portion C of the dune-evidencing low-level human activity in the last year and the progressive succession of vegetation; (c,d) blowout temporarily filled with water; (g,h) pine monoculture overgrowing the dune; (k,l) slope processes on the dune caused by human activity.

Cartographic and Remote Sensing Data
The analysis of source materials included historical and contemporary cartographic materials about the study area [51,[53][54][55]. Due to the lack of historical remote sensing data (aerial imagery) covering the study area, contemporary remote sensing materials (orthophotomap and satellite imagery) were analysed and evaluated. These provided an important additional source of data for acquired historical maps, mainly in terms of the scale of transformations in the study area. Three Tactical Maps of Poland at various scales were analysed: 1:100,000 [53,54], 1:25,000 [55] and the Topographic Map of Poland 1:10,000 [51] (Figure 3). Landforms are depicted on these maps as contour lines. Due to the inappropriate scale of the historical maps, none were used to obtain the spatial data necessary for the reconstruction of the primary relief of the studied dune. Appropriately processed spatial data acquired in the field with a GNSS receiver (Section 3.2), were used to reconstruct the primary relief. The reference standard was the Topographic Map of Poland 1:10,000 from 1986 [51] (Figure 3e) because of its appropriate scale.

GNSS Measurements
In June 2020, topographic mapping of the surveyed aeolian form and the surrounding terrain was performed utilising geodetic equipment that included a GNSS receiver (Spectra Precision SP20) with access to real time kinematic correction (RTK) from reference stations ASG-EUPOS, as well as a precision single frequency antenna (Ashtech ASH111660). Application of the precision antenna allowed the performed measurement accuracy to be increased, especially under the tree canopies of the closed forest area where satellite access was obstructed. GNSS measurements were performed in the coordinate system ETRS89/Poland CS92, with code EPSG: 2180, including quasigeoid model PL-geoid-2011 and height system PL-KRON86-NH. Satellite levelling of the GNSS measurements to the currently valid geoid was performed by means of the following Formula (1): where H orth is the orthometric height (the height over the defined geoid determined through Earth gravitation and mean sea level (MSL)), H ellipse is the GPS height above the reference ellipsoid, and N is the geoid height. As a result, 4045 measurement points were obtained ( Figure 4a). These provided the source of the necessary input data for further advanced geostatistical analysis. The points included all changes in dune morphology, both natural and human activity. This approach made it possible to conduct an accurate reconstruction of the present state of the studied dune, particularly in places of sand mining activity. Due to the dry period in which the measurement was conducted, the blowout was not filled with water. Therefore, its partial mapping was performed so as to provide a near-complete view of the present study area. On the basis of the collected terrain data, a hypsometric map and a digital elevation model (DEM) with a resolution of 1 m were generated. To conduct relief reconstruction from the period before sand exploitation began, the points referring to the fragments transformed anthropogenically, and the blowout located to the west, were removed from the acquired data set. These indicated weak spatial autocorrelation (Section 4.1). Thus, the total number of points used in the reconstruction was equal to 1999 ( Figure 4b).

UAV Mission
A polygon mission for 2D mapping using an UAV (Unmanned Aerial Vehicle) was performed on 21 June 2020. This time was also optimal for collecting field data corresponding to the summer period of the growing season (Section 3.6). In order to obtain the current orthophotomap, the flight of a DJI Mavic 2 Pro unmanned aerial vehicle was undertaken at the height of 100 m above the terrain within the limits of the surveyed dune. During this flight, 160 photos with approximated Ground Sampling Distance (GSD) of 2.34 cm/pixel were registered ( Figure 5). The generation of an orthomosaic of the surveyed object and adjacent terrains was the selected task. Table 1 presents the detailed mission parameters.

Geostatistic Analysis
The present and primary relief of the surveyed parabolic dune and terrain situated in its surrounding area was elaborated by applying the ordinary kriging (OK) method according to the protocol of Zarychta et al. [16]. An experimental semivariogram was used to describe the spatial pattern of the variable under study (elevation). This was defined as half of the variance between two locations z(u α ) and z(u α + h), distant by vector h [16,56,57] (2): where N(h) is the number of the pairs of points. In order to describe the spatial variability of the present and primary relief, a Gaussian model with a practical range a [57] was used as a theoretical model. This was adjusted to the experimental semivariogram obtained from terrain data (3): where c 0 is the spatial uncorrelated variance ("nugget effect"), c is the spatial correlated variance ("partial sill").
The OK method provides an estimate of the variable at unmeasured locations with known values of the measured variable [57] (4): where λ OK α is the weight assigned to n(u) random variables Z(u α ). The correctness of each of the conducted modellings was verified by cross-validation (CV). It consisted of removing the value of each sampling location in the dataset one by one and predicting the removed values to compare theirs with the measured value. The following errors were calculated: mean error (ME), root mean square error (RMSE), root mean square standardised error (RMSSE).
The geostatistical modelling was performed in commercial software ArcGIS 10.7.1 (ESRI, Redlands, CA, USA) and Surfer 20 (Golden Software, Golden, CO, USA).

Data Processing from Sentinel-2A and UAV Image Fusion
In order to complete the data and increase analysis range, satellite images from Sentinel-2A, level L1C from 1 July 2020 were acquired (data were obtained from the EarthExplorer website through the courtesy of the U.S. Geological Survey). This was subjected to radiometric correction in order to eliminate the disadvantageous impact of aerosols (among other issues) on the value of the electromagnetic radiation reaching the Earth's surface. The products of the L2A level were thereby obtained. The lack of accessible atmospheric window caused by cloudy weather was the reason of divergence in terms of imaging acquisition from Sentinel-2A in relation to the UAV mission. The scene acquired from Sentinel-2A was then used in the condition analysis of the vegetation planted on the dune and throughout its surroundings. Due to the insignificant size of the terrain surface, 10 m band resolution: 8 (NIR-near infrared), 4 (RED), 3 (GREEN), 2 (BLUE), proved to be too low, and precise analysis execution was impossible. Therefore, Sentinel-2A and UAV image fusion was undertaken using pansharpening. In effect, resolution of the selected 8, 4, 3, 2 bands from Sentinel-2A scene was increased up to 1 m on the basis of the data acquired from the UAV. The procedure was conducted in accordance with the assumptions adopted by Zhao et al. [58] and Jenerowicz and Woroszkiewicz [59] with minor modification. It relied on the selection of only four bands (8, 4, 3, 2) from Sentinel-2A.
Hence, there was no necessity to conduct the Gram-Schmidt transformation for Sentinel-2A-derived data. In addition, the imaging resolution from UAV was downscaled to 1 m in order to eliminate artefacts that occurred during the flight execution and further stages of the photogrammetric processing. The UAV panchromatic band (PAN UAV ) was calculated by averaging the red, green, and blue bands [60]. For Sentinel-2A, the panchromatic band (PAN SEN ) was calculated as the average of channels 8, 4, 3, 2 [61]. The mean and standard deviation were then calculated for PAN UAV and selected bands (8, 4, 3, 2) from Sentinel-2A. PAN UAV was then stretched so that its mean and standard deviation matched PAN SEN . The PAN UAV band was swapped with PAN SEN using a Gram-Schmidt transformation where the data were transformed into higher resolution multispectral bands [58].

Analysis of the Vegetation Condition
In 2020, during the growing season, botanical data about the occurrence (presence/absence) of trees, shrubs, herbaceous plants, mosses and lichens were also acquired in the study area. Analysis of the condition of the vegetation currently growing upon the dune and in its surroundings was performed on the base of a created composition of false colours with the use of near infrared (CIR-colour-infrared), bands 8 (NIR), 4 (RED), 3 (GREEN) and Modified Chlorophyll Absorption in Reflectance Index 1 (MCARI1) calculated according to the following Formula (5) [62]: where N IR-near infrared band, RED-red band, GREEN-green band. The MCARI1 was selected due to its resistance to pigment changes and reactivity to changes in leaf area index (LAI) [62]. Calculation of CIR and MCARI1 was preceded by satellite and UAV image fusion for enlarging the spatial and spectral range of the data obtained from the UAV. This allowed the earlier acquired orthophotomap and satellite imaging from Sentinel-2A to be combined.

Present and Primary Relief
The use of advanced mathematical modelling, combined with GNSS data, enabled a complete relief reconstruction of the parabolic dune under study. The present relief of the surveyed aeolian form was established, taking into account the surrounding terrain. In contrast, the primary relief included only the dune due to the lack of detailed spatial information on the relief of the surrounding area. This resulted from the weak spatial autocorrelation for the GNSS-derived dataset prepared to perform the primary relief reconstruction of the parabolic dune. Considering that the description of the spatial variations would be incomplete for the surrounding area, it was limited. Hypsometric maps and DEMs presenting relief of the selected dune were elaborated with the terrain data obtained from GNSS (Section 3.2) using ordinary kriging (OK) (Section 3.4). The OK procedure was performed by modelling the spatial variability described with the help of the semivariogram. The matching of the Gaussian model to the experimental semivariogram calculated from the GNSS data ( Figure 6) was checked by cross-validation. The model parameters are presented in Table 2, and the cross-validation errors are presented in Table 3.
Using the ordinary kriging estimator, its standard deviation (SD OK ) was calculated. The following SD OK values were obtained for the present relief: 0.251 m, and the primary relief: 0.196 m, respectively. Hence, the SD OK values were taken into account for the analysed estimated elevation values.   On the basis of the results (Figure 7a,b), it was concluded that the current crest of the dune reaches approximately 240.4-251.1 m a.s.l. and that the maximum relative height in its uppermost portion reaches 10.6 m. In addition, the crest is the most deformed fragment of the dune, within which the three degraded portions are indicated clearly: the eastern one (A) of surface area 794.9 m 2 , and two southern ones: (B) of surface area 285.8 m 2 , and (C) of surface area 2779.5 m 2 ( Table 4). Portions A and C are not being exploited, in contrast to portion B from which sand for building purposes is currently being recovered ( Figure 1). In every portion of the depleted dune portion, in its upper portions, the presence of slightly concave, almost vertical fragments of the youth slope was noted ( Figure 2). The surveyed parabolic dune possesses two arms, a well-developed northern arm and a poorly preserved southern arm, both of which are smoothly bent in a northwest and southwest direction.  The primary relief (Figure 7c,d) analysis of the surveyed dune allows one to conclude that it was relatively well-developed before the beginning of the sand exploitation for building purposes and distinguished by the regular shape characteristic of a typical aeolian form. The exception was the southern arm in its embryonic stage of development. The crest portion was well-developed and rose to the height of 251.2 m a.s.l. (Table 4). It dropped rather smoothly in the southern-eastern direction where it reached the height of 241.7 m a.s.l. More clearly so than is the case currently, the shorter arm was well-distinguished, situated at the height of 234.3 m a.s.l. and had reached the width of 38.4 m. The northern arm was also well-developed and had the shape of a regular levee and rose to the height of 248.1-245.8 m a.s.l. and had reached the width of 55.4 m. In turn, the ridge portion had been slightly fragmented and within its range, five trough-like depressions at an early stage of development were indicated. The biggest of these was situated at the junction of the ridge and northern arm portions and had reached a surface area of approximately 1035.7 m 2 (Figure 7c,d). The reconstructed primary relief of the surveyed parabolic dune is similar to that shown on the three Tactical Maps of Poland [53][54][55] and Topographic Map of Poland [51] ( Figure 3). The arms of the dune-north and south-are clearly visible only on the Topographic Map of Poland 1:10,000 [51], due to the large scale. The crest and the ridge portion are similarly so and the direction of this is consistent with that obtained from the reconstruction. Although the map does not show the relief from before the start of mining, rather only that from the period of its duration, traces of its degradation related to sand exploitation are noticeable on the eastern and southern sides (Figure 3e). The other maps are too small in scale to make a detailed comparative analysis of the different parts of the surveyed dune. On these, only the crest and the ridge of the studied object are distinguishable [53][54][55] (Figure 3a-c).

Vegetation Condition on the Dune and Its Surrounding
During image fusion, the mean and standard deviation were calculated for the acquired panchromatic band from the UAV (PAN UAV ) and for the panchromatic band from Sentinal-2A (PAN SEN ) ( Table 5). The calculated statistics allowed the PAN UAV histogram to be appropriately stretched to PAN SEN (Figure 8). Consequently, the high-resolution multispectral images were achieved by applying the Gram-Schmidt transformation process. It was also possible to create colour-infrared (CIR) compositions and calculate Modified Chlorophyll Absorption in Reflectance Index 1 (MCARI1).  The colour-infrared (CIR) composition indicated the vegetation occurrence with differentiated habitat requirements within the dune range, including its depleted portions and the adjacent terrains. In this range, the exception was a blowout periodically filled with water situated to the west of the surveyed object where species composition was connected especially with the hydrophytes (Figures 2 and 9a).
The calculated Modified Chlorophyll Absorption in Reflectance Index 1 (MCARI1) points to low chlorophyll content (0-0.42) in the leaves of the plants overgrowing the depleted dune fragments (A, B, C). This is vegetation that is in the initial growth stadium. Higher MCARI1 values were noted in the ridge dune portions where pine monoculture occurs. In the blowout to the west, the MCARI1 took its lowest values (0-0.22) where a small number of plant species were found, including, among others, species from the genera Juncus and Typha. In turn, the highest index values (>0.81) were noted for the grassland situated to the northeast from the surveyed object, which confirmed the very good health condition of the vegetation in this region (Figures 2 and 9b). Results from CIR and MCARI1 were confirmed by the field observations that were performed in June 2020. Currently, the facial portion of the dune, including its arms, has been overplanted by evenly aged sections of Scots pine (Pinus sylvestris) production monocultures. The mean age of the trees is approximately 35-40 years. Moreover, in these dry and poor habitats, natural Scots pine renewals, young specimens of the silver birch (Betula pendula) from natural seeding, and a predominating share of common heather (Calluna vulgaris) and lingonberry (Vaccinium vitis-idaea) were found, whereas the slopes were covered primarily with grey hair grass (Corynephorus canescens), as well as wood smallreed (Calamagrostis epigejos). Operation measures leading to the stabilisation of volatile sands through pine monoculture planting probably took place around the year 1980. The role given to the artificial pine forest stand was the prevention of the wind-blown dispersal of aeolian material onto the nearby grassland in areas to the northeast of the dune. In turn, in the blowout periodically filled with water because of low groundwater level, among other plant life noted was a marked predominance of compact rush (Juncus conglomeratus), sweet flag (Acorus calamus), common reed (Phragmites communis), and white water-lily (Nymphaea alba) (Figure 2). Within non-forested parts of the dune, the presence of mosses and lichens growing on the ground or the litter layer was also recorded.

Discussion
A geostatistical approach was used in the performed survey. This allowed the reconstruction of the most probable shape and relief of the surveyed aeolian form, particularly the exploited parts, and to elaborate its present relief. In both cases, the Gaussian model was applied to describe the spatial variability. A theoretical semivariogram of this type effectively describes the topographic elevation of gently undulating hills [57]. Therefore, in the present case it was deemed to be appropriate for use on the surveyed aeolian form. The results of the cross-validation (CV), where mean error (ME), root mean square error (RMSE), and root mean square standardised error (RMSSE) were calculated (Table 3), confirm the degree of the appropriate adjustment of the model.
For the present relief and primary relief, the values of the errors obtained were appropriate. The values of ME and RMSE were close to zero and were −0.009 and 0.564 for the present relief and −0.014 and 0.304 for the primary relief, respectively. Moreover, the values of RMSSE were close to one and were 0.999 for the present relief and 1.077 for the primary relief, respectively. The same theoretical Gaussian model was used by Zarychta et al. [16], who obtained similar results of the CV in the survey relating to the reconstruction of the Liban Quarry in Cracow. In this case, the error values reached ME 0.004, RMSE 0.28 and RMSSE 0.98, respectively. Bzdęga et al. [63] obtained similar results of the CV for spherical models and exponential model in modelling the occurrence of the invasive plant species Sosnowsky's hogweed (Heracleum sosnowskyi) and knotweeds (Fallopia spp.) These are recognised as the most troublesome alien invasive species in Poland [64] and pose a threat to native biodiversity. Their spread also affects the economy and human health. In this work, the obtained error values were ME 0.002-0.013, RMSSE 0.213-0.509, and RMSSE 0.991-1.164, respectively. As a result, the correct variogram structural analysis allowed for obtaining predicted probabilities of occurrence with quantified adequacy of the underpinning models for the targeted species. Fubelli et al. [65] achieved slightly worse results than those found in our research. In the geostatistical reconstruction of the PaleoFarfa River alluvial plain, they found higher CV error values with a higher density of sampling points. In the case of reconstructing top surfaces of the Gelasian fluvial to deltaic deposits, the values for ME of 0.085, RMSE of 4.085, and RMSSE of 0.595 were achieved, while the values in the reconstruction of basal surfaces of the Gelasian fluvial to deltaic deposits reached ME −0.006, RMSE 3.433, and RMSSE 0.421, respectively. Simultaneously the obtained RMSSE values were different from one, indicating a poor fit of the used semivariograms and underestimating of the variability in predictions. Furthermore, RMSE values were 3.521 higher for the reconstruction of top surfaces of the Gelasian fluvial to deltaic deposits and 3.129 higher for the reconstruction of basal surfaces of the Gelasian fluvial to deltaic deposits in relation to those obtained in the reconstruction of the parabolic dune. Similarly, higher RMSE values were found by Cvijetinović et al. [66] by comparing different interpolation methods in terrain reconstruction. For the ordinary kriging method, the obtained RMSE value of 2.35 was 2.046 higher than in the parabolic dune relief reconstruction.
The geostatistical modelling performed allowed the use of ordinary kriging (OK), which in the case of the analysed dune, best represented elevation where point data were missing. Overall, according to Cvijetinović et al. [66], better results may be obtained from applying the OK and multiquadratic methods. However, as the authors noted, geostatistical methods allowed the possibility of calculating the most probable terrain surfaces. This is also confirmed by the studies conducted by Zeng et al. [67], who used ordinary kriging for modelling the correlation structures of column-averaged CO 2 dry air mole fractions (Xco 2 ), considering that the selected estimator was relevant.
In our study, the obtained kriging standard deviation (SD ok ) is prediction dependent and represents the error predicted by the model in the location of the point measurement [68]. Hence, the calculated SD ok values were similar to those obtained by Zeng et al. [67], but much lower than those achieved by Cvijetinović et al. [66]. The SD ok values herein were 0.251 m for the present relief and 0.196 m for the primary relief, respectively, while in the study of Zeng et al. [67], SD ok values for the Xco 2 variable ranged from 1.87 to 2.36 ppm. In contrast, Cvijetinović et al. [66] obtained SD ok of 2.28 m. This was 2.029-2.084 higher than that found for the parabolic dune reconstruction. The obtained SD ok values are considered to be prediction dependent and to indicate a slight change in elevation due to the use of the ordinary kriging (OK) estimator compared to performed predictions.
The reconstruction of the primary relief of the dune surveyed in the form of a hypsometric map and digital elevation model (DEM) was obtained as a result of processing points from field measurements. Their credibility, apart from the CV, was also confirmed by the analysed cartographic materials (historical), i.e., the Topographic Map of Poland 1:10,000 (1986) [51]. Its comparison with the developed hypsometric map and digital elevation model (DEM) revealed that these correctly predict elevation values to a very large extent. Although there are no measurement points within the exploited parts of the dune, the obtained contour lines are analogous to those contained in the 1986 map [51] due to the determination of the spatial variability of the examined variable (elevation) and the application of ordinary kriging. The results of the research on the surveyed dune relief were completed with a health condition analysis of the vegetation occurring within its range and on the surrounding terrain. This was conducted on the basis of the field observations and on that of the colour-infrared (CIR) analysis and Modified Chlorophyll Absorption in Reflectance Index 1 (MCARI1) that were preceded by satellite Sentinel-2A and UAV image fusion. Application of the image fusion allowed imaging of significantly higher spatial resolution to be obtained, with a simultaneous preservation of the spectral resolution [59]. The proposed approach emerged as the better solution during the CIR composition and the MCARI1 calculation as compared to the use of only satellite images with high resolution. Such an outcome was confirmed by the research results by Zhao et al. [58], who used image fusion for improving cultivation classification quality in the finer spatial scale. The MCARI1 values depend on the relative abundance of chlorophyll. Thus, the MCARI1 and CIR obtained in our survey disclosed that definitely better habitat conditions for vegetation prevail on the terrains surrounding the dune than within the dune area. Furthermore, many factors affect this, including dune topography and the availability of water and plant nutrients. The analogical impact of topography on the condition and growth of vegetation on an inland dune in Toruń was also shown by Sewerniak and Jankowski [21].
On the basis of the performed analyses, it was concluded that the surveyed aeolian form is located on flatland approximately 15 m from a no-outflow terrain depression that is periodically filled with water. This shallow depression is a blowout that consists primarily of a flat surface built from loose, fine-grain material prone to wind movement and transportation. Hence, it was eroded and transported by winds from the west and deposited onto the nearby terrains, and was accumulated there. Consequently, a small aeolian form was built with a clearly extended crest portion and two declining arms that are faced against the wind, which was gradually stabilised by the vegetation that took root there. The windward slope was smooth as compared to the leeward slope that constituted the zone of the prevailing sand material deposition. In geomorphological terms, the surveyed object is an inactive inland parabolic dune, and this is suggested by the bending arc of its ridge portion, which in turn can be related with the direction of the prevailing winds. These results were confirmed by Yan and Baas [69] and Hesp [70], who presented and explained the exact course of the mentioned dune processes covering the formation of aeolian forms, dynamics, and sediment transport, and their evolution. Dunes of this type were formed in Central Poland, especially in the Łódź voivodship where terrains prone to aeolian processes are found [48]. The surveyed parabolic dune primarily had a well-formed northern arm, a more poorly formed southern arm, and asymmetric slopes. It was an active form, which is suggested by its crest portion displacement. In order to stop the aeolian processes and limit the spread of the sand onto the nearby grassland, pine monocultures were introduced on the ridge portions of the dune, as occurs very frequently on inland dunes in Europe [21]. In Poland, similar actions were undertaken in the case, among others, of dunes in Poręby and Ogrodzisko [48].
The sand exploitation from portions A, B, and C of the surveyed dune (Figures 2a,b,e,f,i-l and 7) contributed to its structure and relief disturbance. In the mining process, steep walls were created that are characteristic of open pit minings (among otherssand pits). These activated the slope processes and simultaneously initiated the process of slope retreat. The standing non-exploited dune portion A is distinguished by its slope. At the foot of this are slump features in the form of soil heaps. In turn, in the remaining dune portions B and C, regular fragments of the crest break off from top portions of the slopes, which cause trees and bushes to fall over and tumble down the slope, transporting them toward the foot of the dune.
Inappropriate and uncontrolled natural resource exploitation precipitates subsequent geomorphological processes that involve irreversible changes in the relief of the surveyed aeolian form which was also emphasised by Goudie [40] in the case of the parabolic dunes. The lack of action connected with sand recovery abandonment will, in a short time, lead to the complete destruction of the dune, as was observed in the case of a majority of the Holocene dunes formed in the Toruń area [71].
The terrain within which the surveyed aeolian was formed is subject to intrinsic regeneration due to ecological succession as takes place, among others, in the Mojave and Sonoran Deserts (United States of America, USA) [72], the Tengger Desert (China) [73], refuse dumps at the Heidaigou opencast coal mine (China) [74], in gravel-sand pits (Czech Republic) [75], and at Pustynia Błędowska (Poland) [76]. However, if on the terrain covered by the survey, anthropogenic activity continues apace and intensifies in regard to the uncontrolled exploitation of sand and the excessive uprooting of trees from the dune terrain and its surroundings, the aeolian processes will re-activate, producing a natural imbalance on this terrain, entailing considerable economic and material losses.
Our research is innovative because in the proposed solution, the detailed reconstruction of the primary relief and significantly depleted sites was performed without taking into account the historical cartographic materials customarily considered to be the main source of information. Generally, the analysis of geomorphological forms with their application may be difficult or almost impossible, especially because they often present the historical topographic condition in a highly simplified manner.
The proposed solution in the shape of reconstructing the aeolian form with the use of data derived from a GNSS receiver is a low-cost solution and may find application, among others, in the construction of artificial coastal dunes that would constitute a protective barrier against waves and storm surges [45,77,78]. The issue of the reconstruction of the primary dune relief may be useful in practices relating to the consolidation of dunes [21], as well as in environmental management and the formulation of recultivation plans.

Conclusions
Human economic activity contributes to great changes in the landscape. Uncontrolled mineral raw material exploitation is one of the many aspects of such activity. This may affect the various geomorphological forms that are frequently the potential source of the raw material recovered for various economic purposes. The parabolic dunes found in Central Poland are an example of such forms, the crest portions of which have been subjected to significant transformations due to excessive and uncontrolled sand exploitation for building purposes. The performed modelling of the present and primary dune relief allowed the undertaking of the assessment of the changes that have occurred over the years due to human intervention. The surveys confirmed that anthropogenic activity has contributed to the activation of the aeolian and slope processes. Within the surveyed form, vegetation of varied habitat requirements was grown. Currently, the ridge portion of the dune is overplanted with pine monoculture that was introduced to prevent the sand from being wind-blown on the nearby grassland.
The significant achievement of our study, despite insufficient data on the former morphology of the area, but with the help of appropriate and advanced geostatistical tools, was reconstructing the relief of the studied area, including the parabolic dune with its adjacent areas, as hypsometric maps and as digital elevation models (DEMs) based on a set of point data obtained from GNSS receiver.
In conclusion, it should be emphasised that the present elaboration is the first research in which methodical solutions in reference to the relief changes of the accumulative aeolian forms was proposed. This utilises geostatistical modelling for full reconstruction of an accumulative aeolian form relief in an anthropogenically transformed area, and includes assessment of vegetation condition. It constitutes a valuable source of information for institutional bodies involved in the management and protection of the objects of this type that are situated not only on inland terrains, but on coastal terrains as well.