2022 | 75/1 | 51–66 | 14 Figs. | 3 Tabs. | www.geologia-croatica.hr Journal of the Croatian Geological Survey and the Croatian Geological Society 1. INTRODUCTION Landslides are the most widespread geological event affecting 4.8 million people and caused more than 18.000 deaths worldwide in the period from 1998 to 2017 (CRED, 2018). In the same period, landslides, along with volcanoes and other mass movements caused damage of US $ 61 billion. Given that world population growth at the moment shows a linear trend (ROSER et al., 2013), an increase of negative effects related to landslides is to be ex­ pected and these expectations according to SCHUSTER (1996) are based on the following reasons: increased urbanization and development in landslide-prone areas, continued deforestation of landslide-prone areas, and increased regional precipitation caused by changing climatic patterns. The systematic monitoring of landslide occurrences and the damage they cause is not established in Croatia on a national level. However, by inspecting the Register of Damages from Natu- ral Disasters at the Ministry of Finance, it is possible to extract some data for the period from 2015 to 2019. According to the source, the total damage from landslides in that period was US $ 22 million, counting only reported landslides and their casualties in those municipalities or cities that have declared a natural di­ saster in the specified period. This data indicates that landslides in Croatia clearly present a harmful threat to humans, their pro- perties, and the environment. Still, raising public awareness of landslides as geohazard events in Croatia remains a challenging task. As the Earth’s surface is the only human habitat, it is logical to preserve it, which can be achieved especially through preven­ tive activities, primarily through effective planning and manage­ ment (DAI et al., 2002). This approach includes restriction of de­ velopment in the landslide-prone area; use of excavation, grading, Preliminary analysis of a LiDAR-based landslide inventory in the area of Samobor, Croatia Vlatko Gulam*, Iris Bostjančić, Nina Hećej, Marina Filipović and Radovan Filjak Croatian Geological Survey, Sachsova 2, 10000 Zagreb, Croatia * Corresponding author: vlatko.gulam@hgi-cgs.hr doi: 10.4154/gc.2022.12 Abstract The paper presents an analysis of the LiDAR-based landslide inventory for the area near Samo- bor, in northwestern Croatia with two main objectives: i) to define the geological units (obtained from Basic Geological Map of SFRY) most susceptible to landslides, and ii) to analyse the limi- tations of the Basic Geological Map and its applicability in landslide susceptibility map design. Within the study area of 63.8 km2, 874 landslide polygons were manually outlined, covering an area of 2.15 km2. The landslide outline confidence level, landslide index and the relief energy map were used to analyse the landslide susceptibility of a particular geological unit. By that, units in the same state of stress, i.e., in the same relief energy group were compared. This preliminary analysis has shown that the geological units Pl,Q, M3 1,2, and 1M3 1 are the most susceptible to landslides and that older geological units, Pc and K1,2, are also prone to landslides. Still, land- slides within those older units can be considered as old and inactive. As for the limitations of the Basic Geological Map of SFRY, three things emerged, namely scale, the geological unit defining approach, and the neglect of regolith. Despite the limitations presented, the usability of the Ba- sic Geological Map of SFRY in the development of small-scale landslide susceptibility maps is emphasized. However, instructions that should attribute engineering geological features to the geological units outlined in the Basic Geological Map should be prepared in the near future. landscaping, and construction codes; use of physical measures (drainage, slope- geometry modification, and structures) to pre­ vent or control landslides; and development of warning systems. Preventative directions focused on the first step aimed to de­ termine the spatial distribution of landslide-prone areas, assess the probability of landslide occurrence, and estimate the threat of elements within the area affected by landslides. This can be achieved by making several types of maps for the area of inter­ est, covering susceptibility, hazard, and risk (CHACÓN et al., 2006; FELL et al., 2008). However, in order to be able to produce such maps by accepting the general principle that “the past and the present are key to the future” (VARNES & IAEG, 1984), it is necessary to create detailed landslide inventories for a particu­ lar area (GUZZETTI et al., 2012). From that perspective, airborne light detection and ranging (LiDAR) data improved the creation of landslide inventories dur­ ing the last decade, especially when considering their precision and detail. In this study, a detailed LiDAR-based landslide inven­ tory is prepared for the Samobor area located in the northwestern part of Croatia. This area is characterized by high geological di­ versity which within the study area enables as many as 18 geo­ logical units (GUs) to be distinguished within Basic Geological map "of SFRY 1:100.000, Zagreb sheet (ŠIKIĆ et al., 1977), herein after referred to as BGM" (BGM). Based on the developed landslide inventory and available geological data, an extensive spatial analysis is performed, with the main purpose of gathering knowledge of landslide occurrence within each GU. It provides the basis for accomplishing the main objectives of this paper: i) to define the GUs within the study area most susceptible to landslides, and ii) to analyse the limitations of the available geological maps and their applicability in the landslide susceptibility map (LSM) design. Article history: Received October 01, 2021 Revised manuscript accepted November 01, 2021 Available online February 28, 2022 Keywords: landslide inventory, LiDAR, Basic Geological Map, geological units, landslide susceptibility G eo lo gi a C ro at ic a 52 Geologia Croatica 75/1 2. STUDY AREA The boundaries of the study area are defined within the safEarth project co-financed by the Interreg IPA CBC Programme Croatia - Bosnia and Herzegovina - Montenegro. One of the main goals of this project was to create LSMs at different scales and in dif­ ferent areas. In order to achieve this goal, an airborne LiDAR survey was performed on six pilot polygons, for which detailed landslide inventories were created. The selection of pilot polygons was based on the geological analysis, attempting to cover as wide a variability of geological environments as possible. The special geological features of the Samobor polygon, which is presented within this paper, set it apart from other safEarth polygons exactly due to the large number of GUs and the fact that a part of the polygon is composed of carbonate rocks of Triassic age. The study area covers approximately 63.8 square kilometres and is located 19 kilometres southwest of the centre of Zagreb. More precisely, it is situated between Sveta Nedelja, Samobor, and Rude in the north, and Drežnik Podokićki and Rakov potok in the south (Fig. 1). Most of the area belongs to the Sava river basin (38.38 km2), and a smaller part to the Kupa river basin (25.46 km2). 2.1. Geological settings Landslides as exogenetic processes are the result of the interaction between several Earth systems; the lithosphere, atmosphere, and hydrosphere. Professionals and scientists often analyse different natural and human factors that can trigger the landslides (e.g., LI et al., 2012; ZEZERE et al., 1999; LUKIĆ et al., 2018). However, to bring the rock material into a state of motion, the rock material first must be in a state of labile equilibrium, i.e., a large percenta ge of the rock material strength must be “consumed”. This primarily depends on two factors, namely the strength of the rock material, and then the stress (i.e., potential difference) that the rock material suffers. As this paper for the first time com­ prehensively analyses the landslides within the study area, it is logical to present that area according to these two primary fac­ tors as follows: – rock material strength is described through geological characteristics; – potential difference is described through the endogenetic movements. 2.1.1. Geological units The sheets of Basic Geological Map at scale 1:100.000 that cover the area of Croatia were prepared, which for the most part lasted through the second half of the 20th century. It represents the most comprehensive overview of geology in the entire Croatian territory and as such represents an unavoidable basis in any type of research that requires knowledge of the geological conditions of a particular area. According to BARNES & LISLE (2004), there are four main groups of geological maps, and the BGM falls into the category of Regional geological maps, and as such it is produced by col­ lecting field data together with the analysis of photogeological documentation. The geological settings of the study area (Fig. 1) are extensively described within the BGM sheet of SFRY 1:100.000 (ŠIKIĆ et al., 1977) and an accompanying guide (ŠIKIĆ et al., 1979), according to which an overview of GUs pre­ sented within the Samobor area is given. The study area is composed of Palaeozoic, Mesozoic-Palaeo­ gene, Neogene and Quaternary deposits. The oldest formation of Permian age (P2,3), is represented by brownish sandstones with interlayers of shale and siltstone with sporadic transitions to quartz conglomerate. Sediments of the Lower Triassic (T1) appear in the form of purple and reddish mica sandstone and siltstone and in a smaller amount of ooid limestone. In the upper part, ooid limestones with marl intercalations predominate. The sediments of the Middle Triassic (T2) are represented by the transition from diagenetic dolomitized stromatolite, dolo­ pelmicrite and fenestral dolopelmicrite. Also, sporadically, thick- layered and massive crystalline limestone and dark mudstones are observed. The massive dolomite with weakly expressed stratification of the Upper Triassic (T3) are characterized by frequent interca­ lation of dolopelmicrite, fenestral dolopelmicrite and characteri- stic thin laminated dolostromatolite. The deposits of the Cretaceous period (K1,2) are quite litho­ logically heterogeneous and are represented by medium­grained greenish sandstone, with shale, radiolarite and chert alternation. Within the sedimentary rocks a larger mass of basic igneous rocks like diabase, spilitized diabase and spilite (bb) are observed. The Palaeocene (Pc) deposits lie transgressively on older Tri­ assic or Cretaceous deposits and they are mostly represented by polymictic conglomerate with imbrication of marl, siltite and sand­ stone, as well as bioclastic limestone and biolithite are observed. The sediments of Ottnangian age (1M2 1) are lithologically heterogeneous with conglomerate, breccia, and gravel dominate while siltstones are subordinate. The Badenian(2M2 2) sediments are characterised by biocal­ carenite, conglomerate and biolithite limestone developed in the ridge and coastal facies, while silty marl and clayey limestone developed in deeper marine sedimentary environments. Sarmatian deposits (1M3 1) disconformably overlie Badenian deposits in the form of biocalcarenite and sandstone deposits in the marine environment with reduced salinity, as well as charac­ teristic deep­ water laminated marl and silt with sporadic inter­ calations of sand. The shallow lake deposits of the Lower Pannonian (M3 1,2) better known as the „Croatica layers“ are characterized by strati- fied clayey limestones of irregular strata that may be partly inter­ bedded with marl. In the upper Pannonian, marls are deposited in a deeper lacustrine environment with interlayers of sand and sandstone. Concordant on the marl, are pro-delta and delta deposits of the Upper Pontian (Pl1,2) in the form of sand, calcite-rich siltstone and subordinate marl. Sediments of Pliocene-Pleistocene age (Pl, Q) were formed in the transitional, fluvial to lacustrine environment and are rep­ resented by interbeds of gravel with sand lenses, sand with gravel lenses, clayey to sandy silt and clay with sporadically strong li­ monitization. Loess (l) and loess-like sediment (lb) are mostly silty clay and clayey silt, while proluvial deposits (pr) are in the form of coarse-grained, slightly rounded gravel mixed with sand and clay. Gravel, sand and clay are deposits of alluvium (a), and deposits of the second Sava terrace (a2). 2.1.2. Endogenetic movements within the study area The ascent of internal energy originating in the Earth’s core drives a complicated set of geological processes (HUGGETT, 2017). This energy induces endogenetic forces that continuously G eologia C roatica 53Gulam et al.: Preliminary analysis of a LiDAR-based landslide inventory in the area of Samobor, Croatia elevate or build up parts of the earth’s surface and hence the exo­ genetic processes fail to even­out the relief variations of the sur­ face (UPPAL, 2006). In that way, they give rise to the fundamen­ tal relief forms and condition the manifestation of exogenetic processes (ERSHOV et al., 1988). From the aforementioned fundamentals, it is clear that before any pioneering research of the exogenetic processes, it is neces­ sary to determine whether the studied area is in a state of uplift or subsidence. Thus, if the studied area is characterized by a pro­ nounced uplift concerning the surrounding areas, the resulting Figure 1. Study area: Location map (on the bottom left) and BGM of SFRY 1:100,000, Zagreb sheet, according to ŠIKIĆ et al. (1977). G eo lo gi a C ro at ic a 54 Geologia Croatica 75/1 height difference can be considered as the accumulated energy of that area. Such accumulated energy will cause a whole range of exogenetic processes, e.g., erosion (initiation of erosive surface flows), gullying (a set of exogenetic processes that take place along the erosive concentrated surface water flow), landslides, and rockfalls. Since the intensity of exogenetic processes depends, among other things, on the rate of uplift, it is necessary to quantify that rate. One of the first maps dealing with the neotectonic activity in Croatia is the neotectonic map at a scale of 1:1,500.000 (Fig. 2) produced by the Federal geological institute from Belgrade (ĆIRIĆ, 1967). The map shows that the studied area is characterized by ne­ otectonic uplift, while a good part of the southwestern part is characterized by neotectonic subsidence. Thereby, the scale of the presented map should be kept in mind. While analysing the endogenetic movements in the study area, a proposal for future research has arisen. The correlation between neotectonic activity and the intensity of exogenetic pro­ cesses would present a valuable study, the results of which could be used as one of the crucial data sets in landslide susceptibility assessment on a national level. 3. DATA AND METHODS 3.1. Remote sensing data sets As predicted by JABOYEDOFF et al. (2012), the LiDAR sensor has become a standard tool for landslide analysis in this modern time. In the framework of the safEarth project, an airborne sur­ vey was performed by Flycom Technologies d.o.o. during the spring of 2018. A LiDAR survey, with requested point density of 20 points per m2 and an accuracy of ±10 cm in each direction for each collected point, enabled the derivation of a high-resolution DTM with a cell size of 0.5 m. A topographic survey resulted with digital orthophotos with a 10 cm resolution. The DTM was used to create several crucial derivatives (e.g., ARDIZZONE et al., 2007; VAN DEN EECKHAUT et al., 2007; JABOYEDOFF et al., 2012) for the production of landslide in­ ventories: – hillshade grid with a cell size of 0.5 m, derived by a sun azimuth angle of 315° and altitude angle of 45°; – slope grid with a cell size of 0.5 m; – contour lines with one m elevation equidistance. 3.2. Landslide inventory creation PIKE (1988) emphasises that a geometric signature is a subset of the geomorphic signatures, a broader and largely undeveloped concept that includes much more than topography. By studying the works that cite the mentioned author, it is possible to interpret the “signature” that exclusively defines the morphometry of a geomorphological phenomenon (GUZZETTI et al., 2012), e.g., a landslide. Even though technically this possible interpretation is not entirely correct, one can conclude that a distinct signature of landslide morphometry should for certain be noticed on high­ resolution elevation models and as such can be used for landslide investigations (JABOYEDOFF et al., 2012). To create a landslide inventory, landslide polygons in this study were outlined manually within the GIS environment using visual analysis and interpretation of the topographic surface pro­ vided by high-resolution DTM derivatives. Primarily, the com­ bination of two data sets was used, slope grid as a top layer with 50 % transparency and the hillshade grid as a bottom layer. In some cases where landslide boundaries were not completely clear and sharp, the use of 3D visualisation software (Esri ArcSceene) was necessary to obtain the best possible results (Fig. 3). Each polygon within the landslide inventory is accompanied by a landslide outline confidence level (LOCL) aimed to describe the level of visual clarity of the individual landslides. Specifically, the LOCL of each landslide was assessed according to how clearly the boundary of the landslide can be traced and features within the landslide body can be observed (e.g., cracks, surface depressions, hummocky topography). The LOCL is expressed by a grade from one to ten, where the highest grades are ascribed to those landslides with the most visible features. A very similar ap­ proach has been presented by BURNS & MADIN (2009). To conduct the analyses, the LOCL grades were clustered here into four groups (Fig. 3). 3.3. Field prospection The field prospection was performed as a part of this work, pri­ marily to gain new knowledge regarding the engineering geo­ logical properties of materials present in the study area. A field prospection at 190 observation points gave a detailed overview of: – granulometric composition of Quaternary deposits which are outlined on the BGM; – lithological composition and durability of the bedrock ma­ terial; – types and characteristics of regolith, primarily thickness and granulometric composition. It should be emphasized that the data from the field prospect­ ing are insufficient for extensive statistical analyses. However, the experience gained on that occasion greatly helped in the ar­ ticulation of certain, for this paper, important conclusions. Figure 2. A fragment of the Neotectonic map of the SFRY (ĆIRIĆ, 1967). G eologia C roatica 55Gulam et al.: Preliminary analysis of a LiDAR-based landslide inventory in the area of Samobor, Croatia 3.4. GIS analysis 3.4.1. Landslide index Landslide index is used to express the spatial relationship between delineated landslide polygons and GUs. According to BARTEL­ LETTI et al. (2017), it is calculated as the ratio between the area of landslides occurring in a given GU and the area of that particular GU. In the same way, the landslide index of the study area is cal­ culated as the ratio between the total landslide area and the entire study area. Since landslide indices are expressed in percentages, it was possible to compare the landslide index of each GU with the landslide index of the study area and draw the first conclusions about the susceptibility of each GU to landslide occurrence. 3.4.2. Landslide density To get a better insight into the spatial distribution of landslides over the entire investigated area, a landslide density map was created. Landslide inventory is transformed to a grid format, Figure 3. Examples of landslides with different LOCL grades presented on DTM derivatives used in the process of landslide inventory development. G eo lo gi a C ro at ic a 56 Geologia Croatica 75/1 where the value “one” is added to each cell that contains a land­ slide, and value “zero” to each cell without a landslide. Using the Focal statistics ArcGIS tool, a landslide density is calculated for each cell within a circle neighbourhood with a radius of 150 m. A landslide density map is used to delineate the zones for detailed analysis, in order to analyse the spatial distribution of landslides within specific GUs. 3.4.3. Relief energy A relief energy map is derived from the high-resolution DTM (resolution 0.5 cm). Using the Focal statistics ArcGIS Tool, the elevation range is calculated by subtracting the maximum and minimum elevation within the circle neighbourhood with a radius of 150 m. The relief energy map is used as a representation of the stress state within the study area.Figure 4. Distribution of the geological unit area within the study polygon. Figure 5. Landslide inventory of the study area. G eologia C roatica 57Gulam et al.: Preliminary analysis of a LiDAR-based landslide inventory in the area of Samobor, Croatia 4. RESULTS 4.1. Statistical sample analysis Statistical sample analysis is performed to see if the area of the individual GU is statistically significant within the study area. To determine the statistical significance, a simple histogram chart is analysed, comparing the area of each GU with the average area of GUs, which is taken as a reference line (Fig. 4). The presented histogram shows that GUs a2, 1M2 1, and ββ do not present a sta­ tistically significant sample, since the histogram bars of those units are well below the reference line. Accordingly, the area of the GUs pr, Pl12, and Pc have questionable statistical significance with histogram bars a little below the reference line. This should be kept in mind when discussing those GUs, hence the results reached could be misinterpreted. For the other 12 GUs, it can be pointed out that the area they cover determines a statistically sig­ nificant sample, i.e., the conclusions derived from the analysis performed are based on relevant and sufficient data. 4.2. Landslide inventory map A total of 874 landslides were outlined within the investigated area. The delineated landslides cover an area of approximately 2.15 square kilometres, which is 3.4 % of the total study area. The average values of the landslide area, length, and width are 2.462 m2, 70 m, and 37 m, respectively. The landslide inventory map shows that the vast majority of landslide polygons are concentrated in the central part of the study area (Fig. 5). 4.3. Landslide index After the visual analysis of the landslide inventory, the landslide share for the individual GUs was calculated. Some basic statistics related to the GUs and accompanying the landslide inventory are summarised in Table 1. The landslide share analysis was con­ ducted for all 18 GUs. In the analysis of landslide number and landslide index dis­ tribution within GUs (Fig. 6), all landslides (LOCL grades from 1 to 10) were taken into account. Considering the landslide number (Fig. 6a), it is shown that the bars for the GUs Pl12, M3 1,2, 1M3 1, 1M2 1, Pc, K1,2 and T1 rise above the reference line which should classify them as suscepti­ ble or prone to landsliding. This analysis served as a starting point in the study aimed at identifying the “problematic” GUs, regard­ ing the statistical sample size (Fig. 4). The problems are related to the GUs Pl,Q and 1M2 1, as there are 293 delineated landslide polygons within the unit Pl,Q and only 1 within the unite 1M2 1 (Table 1). From the above, it is obvious that the representation of an individual landslide as a point loses the data of its surface area, which equalizes the importance of small and large landslides, masking landslide susceptibility characteristics for a particular GU. However, GUs with large areas and a high number of land­ slides may have a landslide density per km2 below the average of the whole study polygon (e.g., Pl,Q). The analysis that involves the surface area of the landslides within the study polygon is presented through landslide indices (Fig. 6b). This analysis shows somewhat different results to the previous one, as the bar for the GU Pl,Q stands above the refer­ ence line, and for GU 1M2 1 it is well below the reference line. The reason for that could be large landslide polygons delineated within that GU and this assumption is confirmed in Table 1, where the GU Pl,Q is presented as the second one in the rank list of average landslide area for all GUs. Table 1. Basic statistics regarding geological units and delineated landslides within the research area. G eo lo gi ca l u ni t U ni t s ha re w ith in th e st ud y ar ea [% ] N um be r o f l an ds lid es To ta l l an ds lid e ar ea in th e ge ol og ic al u ni t [ km 2 ] Av er ag e la nd sl id e ar ea in th e ge ol og ic al u ni t [ m 2 ] Av er ag e la nd sl id e le ng th in th e ge ol og ic al u ni t [ m ] Av er ag e la nd sl id e w id th in th e ge ol og ic al u ni t [ m ] a2 0.1% 0 - - - - a 5.7% 25 0.02653 1061 37 30 pr 1.3% 0 - - - - lb 2.7% 1 0.00289 2894 122 32 l 5.8% 4 0.00369 923 72 16 Pl,Q 35.9% 293 0.93064 3176 62 45 Pl12 2.3% 35 0.07266 2076 59 33 M3 1,2 11.4% 215 0.42321 1968 82 28 1M3 1 4.9% 84 0.17129 2039 81 30 2M2 2 6.7% 25 0.05746 2298 97 32 1M2 1 0.1% 1 0.00037 366 23 21 Pc 2.5% 33 0.10969 3324 91 44 K1,2 5.9% 71 0.18799 2648 71 41 bb 0.0% 0 - - - - T3 4.0% 8 0.00954 1193 64 26 T2 4.4% 20 0.02634 1317 62 27 T1 3.3% 41 0.09073 2213 58 45 P2,3 3.1% 18 0.03836 2131 59 47 Figure 6. a) Number of landslides per km2 within geological units and number of landslides per km2 within the study area; b) Landslide index distribution for the geological units and landslide index for the study area. G eo lo gi a C ro at ic a 58 Geologia Croatica 75/1 Other GUs are in a similar relationship to the reference line as in the previous analysis. Grey bars of the units Pl,Q, Pl12, M3 1,2, 1M3 1, Pc, K1,2, and T1 (Fig. 6b) are above the reference line which characterises the listed GUs as susceptible to landsliding. Such results did not entirely match the engineering impression that was formed during the field prospection within the study polygon and during the creation of the landslide inventory. This refers to the older GUs Pc, K1,2, and T1, and especially to unit Pc. Although landslides can be observed within the mentioned GUs, they share somewhat different characteristics than landslides within the units Pl,Q, Pl12, M3 1,2, and 1M3 1. Such observations encourage further landslide index analy­ sis in which landslides of certain LOCL were taken into the ac­ count. Landslides were grouped into four groups, based on the following LOCL grades: 1) from one to three (Fig. 7a), 2) four to five (Fig. 7b), 3) six to seven (Fig. 7c), and 4) eight to ten (Fig. 7d). Part of the LOCL grade can very likely be attributed to the age of a landslide (MCCALPIN, 1984; BELL et al., 2012) since the grade depends only on the clarity of landslide feature expression, which gradually vanishes on the terrain surface as the landslide ages. In the group of younger GUs (Pl,Q, Pl12, M3 1,2, and 1M3 1), there is no obvious trend as the bars sporadically stand above or below the reference line in all the presented charts (Fig. 7). Still, all of these units have been represented with bars that are well above or very close below the reference line. According to the LOCL, grades attributed to the landslides, both, old and recent landslides are present within those units. That could imply that these young GUs are as prone to landslides nowadays as they were in recent geological history. On the other hand, the bars for older GUs (Pc, K1,2, and T1) show some trend from the first (Fig. 7a) to the last chart (Fig. 7d). It can be observed that bars related to older GUs are mostly above the reference line in the charts in Figure 7a and 7b, and well below the reference line in Figure 7c and 7d. This means that landslides within the older GUs are mostly rated with lower LOCL grades which could led to the con­ clusion that landslides are characterized by “blurry” features and are probably older than the landslides within the younger GUs. According to the accompanying lithological characteristics (a mixture of gravel, sand, silt, and clay), the alluvial unit should be classified as a unit very susceptible to landsliding. Still, be­ cause these are the youngest deposits (Quaternary age) the endo­ genetic movements have not had the opportunity to bring them into a state of labile equilibrium. Therefore, a high landslide share within this unit (Fig. 6a and 6b; Fig. 7a) is unexpected and is a result of the scale of the geological map used, which is elaborated within section 5.1.1. 4.4. Landslide density map The landslide density map (Fig. 8), as a derivative of the landslide inventory, highlights very similar zones as the landslide inven­ tory map (Fig. 5). Yet the landslide density map more clearly sin­ gles out subzones of similar characteristics in terms of landslide density, around which a more detailed discussion is presented in section 5. The subzones on the landslide density map marked as poly­ gons of detailed analyses are separated in the area of two GUs, namely Pl,Q and 2M2 2. The polygons one, two, and three are placed within GU Pl,Q and were outlined using a landslide den­ sity map and an energy relief map (presented within the next chapter). Although located in the same GU, these polygons indi­ Figure 7. Landslide index distribution for the geological units and landslide in- dex for the study area, considering: a) landslides with LOCL graded between one and three; b) landslides with LOCL graded between four and five; c) land- slides with LOCL graded between six and seven; d) landslides with LOCL graded between eight and ten. G eologia C roatica 59Gulam et al.: Preliminary analysis of a LiDAR-based landslide inventory in the area of Samobor, Croatia cate very contrasting environments in terms of landslide density and relief energy, discussed in more detail in section 5.2.2. Polygons four, five, six, seven, and eight presented on the landslide density map are taken from the BGM and represent en­ tirely GU 2M2 2. During the analyses, an interesting difference in the share of landslides was noticed within those polygons, which is also further discussed in section 5.2.2. 4.5. Relief energy map In general, a relief energy map can be used to determine the zones in which rocks are found in a similar state of stress. Just by visual observation (Fig. 9), there is an obvious trend going from the northwest part of the study polygon (the highest values of relief energy) to the southeast part of the study polygon (the lowest values of relief energy). In addition to the presented map, the analysis of relief energy for GUs was performed and the histogram chart was made (Fig. 10). This simple chart in a way follows the mentioned relief energy zones as the oldest geological units are located in the northwest part of the polygon with the highest average relief energy values. This obvious spatial trend can be followed in the histogram chart shown in Fig. 10, as GUs presented through average relief energy can be divided into the three following groups: Figure 8. Landslide density map of the study area. G eo lo gi a C ro at ic a 60 Geologia Croatica 75/1 a) Quaternary group – units: a2, a, pr, lb, and l; b) Cenozoic and Cretaceous group – units: Pl,Q, Pl12, M3 1,2, 1M3 1, 2M2 2, 1M2 1, Pc, K1,2, and ββ; c) Triassic and Permian group - units: T3, T2, T1, and P2,3. 5. DISCUSSION 5.1. Geological units prone to landslides Numerous studies in the professional and scientific literature dis­ cuss the conditions under which landslides may occur. In fact, the complete procedures for obtaining landslide susceptibility are based on the analyses of various conditioning factors, i.e., geo- environmental variables. Still, there are no established rules on which conditioning factors to select within the study, thus those selections are mostly based on expert opinion (JEBUR et al., 2014). According to REICHENBACH et al. (2018), based on re­ view of 565 articles, geo-environmental variables are classified into 23 classes within five clusters: i) morphological (slope, as­ pect, curvature, elevation, other morphometric, geomorphologi­ cal), ii) geological (geo-lithological, distance to fault, geo-struc­ tural), iii) land cover (land use/cover, soil, distance to road, tree, other anthropic), iv) hydrological (river/catchment, distance to Figure 9. Relief energy map of the study area. G eologia C roatica 61Gulam et al.: Preliminary analysis of a LiDAR-based landslide inventory in the area of Samobor, Croatia river), and v) other variables (precipitation, earth observation, geotechnical, seismic, other, landslide related, other climatic). All of the factors mentioned above certainly affect the land­ slide processes, and some of them are relatively easy to measure and the related data are often available. Therefore, their use in the process of creating the susceptibility, hazard, and risk maps is logical and above all justified. However, when assessing the sus­ ceptibility of geological or EG units, or better yet when compared to each other, it seems necessary from an engineering point of view to return to the foundations of mechanical settings. In other words, to consider two factors that are crucial for landslide oc­ currence which are the stress conditions within certain parts of the relief, and the strength (internal resistance force) of the ma­ terial within a particular GU. This means that conclusions about whether one GU is more susceptible to landslides than another can be drawn only when the analysed units are in a state of equal stress state, i.e., when their strength is consumed in equal or similar amounts. For such an estimate of landslide susceptibility characteristics for a par­ ticular GU, it is logical to use the combination of the two para- meters presented in this paper, namely landslide index and relief energy. Relief energy is defined as the height difference within the particular terrain surface, and as such, it can indicate a certain state of stress within the particular surface area. By accepting this approach, every energy relief group (Fig. 10) could be con­ sidered as a group that suffers similar stress. Such an approach may also imply that a comparison of landslide susceptibility levels of particular GUs can be performed only within a single energy relief group. For estimation of the strength as a starting point, the land­ slide index can be used. However, for deeper analyses of landslide susceptibility within the one relief energy group, one should cer­ tainly consider the geological composition of the particular GU. 5.1.1. Relief energy group A This group consists of Quaternary deposits the genesis of which is related to water and wind erosion. Namely, these young alluvial, proluvial, and loess sediments have been deposited on floodplains in a wider area and as such are certainly prone to landsliding. Still, endogenetic movements have not yet had the opportunity to uplift them. Consequently, this group is characterized by the lowest level of relief energy, in a range from one to fifteen m (Fig. 10). The low levels of the landslide index are therefore logical and expected. It is clear that in all the graphs (Fig. 6 and Fig. 7) the landslide index in the GUs of the relief energy group A is much lower than the average level of the landslide index for the entire study area. This observation is somewhat disturbed by a slightly larger landslide index for the GU a (alluvium - Fig. 6a and 6b; Fig. 7a). This can be related to the scale of the geological map used, which is elaborated on in more detail in section 5.2.1. 5.1.2. Relief energy group B Relief energy group B has a moderate level of energy in relation to other groups formed within this paper with the range of aver­ age relief energy (Fig. 10) from 25 m (GU Pl,Q) to 43 m (GU Pc). Values of the landslide indices are very scattered within relief energy group B (Fig. 6 and Fig. 7), and that is especially empha­ sized in Fig. 7. This scattering indicates that the GUs are very diverse in terms of susceptibility to landsliding, but it can also indicate different bedrock material, different regolith material, and finally different types and mechanisms of landsliding. The following brief analysis of the landslide index related to the GUs of this group can serve as a starting point for future rec­ ommendations associated with the reclassification of the GUs defined on the BGM into the EG units. This can be especially important for the creation of the EG maps at a small scale as a crucial input data in LSM design. Without going into deeper analyses, one could conclude that within the relief energy group B according to the landslide index (Fig. 6b), there are two groups of GUs. Namely, those with bars standing above the reference line (landslide susceptible units – Pl,Q, Pl12, M3 1,2, 1M3 1, Pc, and K1,2) and those with bars standing below (landslide unsusceptible units – 2M2 2, 1M2 1, and ββ). Such conclusions would be very hasty for several reasons, which will be briefly elaborated below. Statistical sample The low landslide index for the 1M2 1 and ββ units presented in the scope of this research can be attributed to an insufficient statistical sample (Fig. 4). Although from an engineering point of view it would be justified to classify these GUs as unsuscep­ tible to landsliding according to their lithology, it would be wise to support such conclusions in future research by analysing the landslide inventory in areas where these GUs are better repre­ sented. Although the high landslide index is defined for GUs Pl12 and Pc (Fig. 6b), it should be repeated that the statistical sample analysis for those units confirmed a questionable statistical sig­ nificance (Fig. 4). GU Pl,Q covers a very large part of the study area (35.9 %) and according to the landslide index (Fig. 6b) could be labelled as a high landslide susceptible unit. Similarly, GU 2M2 2 covers a statistically significant area (6.7 %) and according to the landslide index (Fig. 6b) could be labelled as a low landslide susceptible unit. Still, the lithology of these units is very diverse so the land­ slide susceptibility level should be defined with caution. The litho logical heterogeneity of the parent rock which strongly influ­ ences the landslide density within those units is discussed as a part of one of the BGM limitations in section 5.2.2. Landslide outline confidence level (LOCL) As already described, the LOCL generally defines the visibility and sharpness of a landslide on DTM derivatives used in land­ slide inventory creation. By accepting the assumption that the LOCL grade can be assigned to the landslide age, it is possible to conclude that GUs M3 1,2 and 1M3 1 are still somewhat more sus­ ceptible to landslides in recent climatic conditions than GU K1,2, Figure 10. Average relief energy distribution for the geological units within the study area. G eo lo gi a C ro at ic a 62 Geologia Croatica 75/1 due to the higher LOCL grades. The answer to the question why unit K1,2 is less susceptible to landslides in recent climatic condi­ tions than units M3 1,2 and 1M3 1, should be sought in the EG char­ acteristics of the near-surface parts of the above mentioned GUs. GUs M3 1,2, 1M3 1, and K1,2 are characterized by sliding sur­ faces located in the regolith formed above the bedrock. Thus, several important features of both regolith and bedrock need to be considered to assess their susceptibility to landslides, and those are: thickness and granulometric composition of the rego­ lith; susceptibility of the bedrock to the mechanical weathering, i.e. formation of the regolith; and the permeability of the bedrock. The granulometric composition, especially the proportion of the coarse-grained fraction, strongly influences the physico-me­ chanical properties of the regolith-like materials (PARK et al., 2018; INDRAWAN et al., 2006; SHAKOOR & COOK, 1990). With the increase of the coarse-grained fraction, the water per­ meability also increases, which prevents the growth of pore pres­ sure which strongly affects the slope stability. Similarly, the wa­ ter permeability of the bedrock strongly contributes to the oscillation of the pore pressure at the contact with the upper rego­ lith. Therefore, if the bedrock is characterised by low permeabi- lity (e.g., marl) it is very likely that the high pore pressure zone will be formed at the contact surface of the regolith and bedrock, presenting a potential sliding surface. By taking all this into account in combination with the knowledge gained during the fieldwork and by studying the lite- rature, it is possible to distinguish two different models. For GUs M3 1,2 and 1M3 1 the regolith can be defined as 1 to 5 m thick, with a very small proportion of a coarse-grained frac­ tion, and low permeability. The bedrock is susceptible to the for­ mation of a thicker regolith and is poorly water-permeable. With such characteristics, this EG model is defined as a high landslide susceptible one, due to the possible rapid disturbance of the me­ chanical stability on slopes. For GU K1,2 the regolith can be defined as 1 to 3 m thick, with a variable proportion of the coarse-grained fraction, and variable levels of water permeability. The bedrock is susceptible to the formation of a medium­thick regolith and with variable level of water permeability. With such characteristics this EG model is defined as being of low landslide susceptibility espe­ cially in the current climatic conditions. The approach presented seems to be a good direction that can be used in the reclassification of the GUs defined within the BGM into the EGUs required for LSM development. This ap­ proach has already shown good results in landslide susceptibility mapping of the Sisak-Moslavina County (BOSTJANČIĆ et al., 2021), which certainly encourages future more serious attempts in a similar direction. Figure 11. a) Relief energy for the alluvial GU within the boundary outlined on the map at a scale of 1:100.000; b) Relief energy for the alluvial GU within the bound- ary outlined on the map at a scale of 1:5,000. G eologia C roatica 63Gulam et al.: Preliminary analysis of a LiDAR-based landslide inventory in the area of Samobor, Croatia 5.1.3. Relief energy group C This group has the highest average relief energy within the study polygon, ranging from 53 m to 64 m, which may be a conse­ quence not only of more pronounced endogenetic uplift, but also the strength of the bedrock material. Landslide indices of GUs within this group indicate extremely low (T2 and T3) and low (T1 and P2,3) landslide susceptibilities (Fig. 6b). Slightly higher land­ slide indices for GUs T1 and P2,3 are very likely due to a thicker regolith superficial layer with characteristics very similar to GU K1,2. Geological units T2 and T3 are characterised by a very thin and permeable regolith superficial layer (up to 1 m), and perme­ able bedrock (dolomite and limestone) which is why they can be classified as having an extremely low susceptibility to landslides. 5.2. Limitations of the BGM in LSM creation One of the pre-defined goals of this paper was to test the usabi- lity of the BGM in the process of LSM creation. Interesting observations elaborated later in this section try to present the well­known limitations of this map in a simple and illustrative way. However, here, based on logical and numerically measurable details, closely related to the occurrence of landslides, an attempt was made to show the level of the limitations of the geological maps that cover the whole Croatian territory. 5.2.1. Scale The scale of any type of map affects the level of detail that can be shown on the map and the spatial accuracy of the cartographic elements (PERKINS & PARRY, 1990), e.g., geological bounda­ ries and faults. How the scale affects the precision of the line ele­ ment (in this case the geological boundary) is shown by a simple analysis of the boundary of the alluvial deposits. During the analysis of the Quaternary geological group, for an alluvial GU, a slightly higher landslide index was observed (Fig. 6 and Fig. 7a). Besides, a somewhat higher value of the ave- rage relief energy was also noticed (Fig. 10). As these Quaternary sediments are deposited on floodplains, such values of landslide index and average relief energy are not logical. This inaccuracy is attributed to the scale of the BGM, which was used to deline­ ate the GUs polygons. In order to confirm the stated claim, anal­ ysis of the polygon outlining part of the alluvial deposits in the study area is presented. In Figure 11 two different geological boundaries are shown. The first is taken from the BGM at a scale of 1:100,000 that was used in this study (Fig. 11a) and the second one is taken from the geological map at a scale of 1:5,000 which was made for the part of study area within the safEarth project (Fig. 11b). Presented polygons were used to clip the relief energy map. The cells that were left within the polygons were analysed and compared. The polygon outlined from the small-scale map contains parts of high energy relief, while the polygon outlined from the large-scale map contains mild energy relief. The basic statistics of these differences are given in Table 2. The age and the genesis of the alluvial material imply that these deposits are placed in plains where endogenetic movements have not yet uplifted the terrain. In that sense, the low average relief energy is largely to be expected. As the lower value of this parameter is obtained using the large-scale geological map, it cer­ tainly seems more accurate and thus more acceptable. Accordingly, the BGM should be used with caution in the process of LSM creation, i.e., it can be used for the production of LSM at the same or similar scale, and in no circumstances at scales larger than 1:25,000. 5.2.2. Geological units defining approach Instructions for BGM design define the approach for outlining the GUs within the different types of rocks (igneous, metamorphic and sedimentary; SAVEZNI GEOLOŠKI ZAVOD, 1962). Thus, for sedimentary rocks, the first criterion is age, while for igneous and metamorphic rocks it is lithological composition. The pre­ sented approach in separating sedimentary GUs from the aspect of engineering geology causes quite a problem, as lithology largely conditions the engineering geological properties. With such an approach, GUs with very heterogeneous lithological com­ position have been outlined, which makes the assessment of their susceptibility to landslides very difficult, and in some cases even impossible. This problem is described through the analysis of two GUs outlined within the study area, Pl,Q and 2M2 2. Pl,Q This unit consists of gravelly, sandy, silty and clayey sediments, mostly deposited in the form of lenses. Solid rocks are present in the form of breccias and conglomerates but their prevalence is almost negligible. Only a few scientific papers analyse Pl,Q deposits (NOVO- SEL-ŠKORIĆ et al., 1986; ŠIKIĆ, 1995; HEĆIMOVIĆ, 2009), but unfortunately, none of them seem to analyse in more detail the lithological heterogeneity of these deposits. In terms of land­ slide activation, the importance of detecting the spatial distribu­ tion of individual fine and coarse-grained lithological members within Pl,Q deposits is presented. This GU is characterised by landslides with a deeper failure surface and the main triggering factor is the increase in pore pres­ sure due to the existence of lenses of coarse fraction material in the fine-grained soil mass. It should certainly be emphasized that the explanation given below regarding the landslide triggering mechanism within the GU Pl,Q came from insufficiently detailed research, which sug­ gests that the given conclusions should be confirmed in the future by conducting research that is even more detailed. Three zones within GU Pl,Q based on different landslide density were separated (Fig. 8), in which the main criterion was a contrasting spatial distribution of landslides. Since all three polygons were isolated in the GU Pl,Q with the primary assump­ tion that they share similar lithological characteristics (with the same resistance force, i.e., strength), the cause of the different landslides density had to be sought in the state of stress within each outlined polygon. Thus, the polygons were added to the en­ ergy relief map (Fig. 9) to determine their average relief energy. Finally, for all three zones, the landslide share and average relief energy were overlapped and analysed (Fig. 12). Zone one, which has the highest landslide share, somewhat unexpectedly does not have the highest average relief energy. Zone two with the lowest landslide share, unexpectedly, has the Table 2. Basic statistics regarding relief energy for the alluvial geological unit within the maps of different scales. Statistics Geological map [1:100000] Geological map [1:5000] Min 1 1 Max 56 44 Average 14.0 11.8 Variance 95.2 51.1 Stan. deviation 9.8 7.1 G eo lo gi a C ro at ic a 64 Geologia Croatica 75/1 highest average relief energy. Even the data related to zone three have somewhat unexpected results because that zone has the low­ est average relief energy and a relatively high landslide share. The analysis presented led to re-examination of the first as­ sumption, which is that the entire GU Pl,Q is characterized by a similar lithology. During preliminary field profiling in all three zones, the field determination of soil samples was conducted. It suggested that zones one and three are characterized by a very heterogeneous lithological environment, i.e., an environment of very diverse granulometric composition in which the fine-grained fraction is somewhat predominant. Zone 2 is more likely charac­ terized by coarse-grained material (sand and gravel). This conclusion made it possible to design the two charac­ teristic EG models for the GU Pl,Q (Fig. 13). Model A (Fig. 13a) is related to zones one and three (Fig. 12), and model B (Fig. 13b) to zone two (Fig. 12). Model A describes landslide activation caused by increased pore pressure in the fine-grained material and erosion of the landslide toe. In model B, landslide activation in the current state of stress is not possible because the present lithology drains very easily, i.e., an increase in pore pressure is not possible. Landslide foot erosion is also not possible because surface erosive flows cannot develop due to the very rapid infil­ tration of surface water into deeper horizons The GU Pl,Q has been indicated in many previous works as a unit susceptible to landslides (JURAK et al., 1998; MIHALIĆ et al., 2008). The graphs in Fig. 7 certainly confirm these claims, however, when designing future LSMs (especially at a smaller scale for which new LiDAR substrates will not be made), this unit should be considered as highly heterogeneous. The recommen­ dations (guidelines) that should follow such maps must indicate such properties of these deposits, which would certainly require very detailed research necessary for the production of LSMs at a larger scale. 2M2 2 The Badenian geological complex stands out with a lower land­ slide index compared to other GUs within the relief energy group b (Fig. 6 and Fig. 7). Although the BGM (ŠIKIĆ et al., 1977) con­ sider this complex as a single GU, from the lithological point of view (VRSALJKO et al., 2006; VRSALJKO et al., 2007) these Figure 13. Characteristic models for GU Pl,Q; MODEL A - Very heterogeneous lithology from fine-grained to coarse-grained fraction; MODEL B - Homogeneous coarse-grained (sand and gravel with very little coarse-grained silt) lithology. Figure 12. Average relief energy and landslide share for the three zones within the GU Pl,Q outlined on the landslide density map. G eologia C roatica 65Gulam et al.: Preliminary analysis of a LiDAR-based landslide inventory in the area of Samobor, Croatia are very diverse deposits, which can be divided into at least two EG units: 1) EG unit 1 - biocalcarenite, conglomerate, and biolithite limestone; 2) EG unit 2 - silty marl and clayey limestone. These lithological units differ in water permeability and level of weathering durability. As such, they create regolith complexes of different thicknesses and granulometric features with very di­ verse water-permeable values (Table 3). From all the above, it is obvious that the units described in this way also differ in their landslide susceptibility level, which was attempted to visualize through the following very simple analysis. The polygons of Badenian deposits (2M2 2) from the BGM (Fig. 1) have been grouped according to the visual criterion of landslide density (Fig. 8). The polygons are labelled on the map from four to eight. For each polygon, the corresponding landslide index was calculated and the obtained results are presented in Figure. 14. As expected, the outlined polygons have very different land­ slide indices, which very likely confirms the thesis of diverse li­ thology. The assumption, which has been verified to some extent by fieldwork, suggests that polygons five, six, and eight are part of EG unit 1 due to their very low landslide index. EG units 1 and 2 equally build polygon seven, due to the landslide index close to the landslide index of the entire study area. Polygon four is proba- bly related to EG unit 2 due to the high landslide index. The pres­ ence of limestone within polygon four can also be confirmed by the large number of sinkholes which are characteristic of this li­ thology. 5.2.3. Neglecting the regolith Regolith is a general term for the layer of unconsolidated or un- cemented, weathered material which overlies bedrock over much of the land surface. It includes rock fragments, mineral grains, and all other superficial deposits such as unlithified in situ sap­ rolite, ash, colluvium, alluvium or drift (KEAREY, 2001; ALLABY & ALLABY, 2003). Superficial deposits that can be characterized as regolith, according to the presented definition, are mostly of Quaternary age, and as such, following the BGM instructions are ignored in most cases (e.g., they are not outlined if their thickness is less than 5 m; SAVEZNI GEOLOŠKI ZA­ VOD, 1962). During the fieldwork in the study area, it has been confirmed for older GUs (M3 1,2; 1M3 1; 2M2 2; 1M2 1; Pc; K1,2, ββ, T3, T2, T1, and P2,3) that landslides occur almost exclusively in the regolith layer. Thereby, it can be concluded that for research re­ lated to landslides, this is the most important limitation of the BGM. 6. CONCLUSION Here, a LiDAR-based landslide inventory is presented and ana­ lysed for the first time, for 63.8 km2 of the Samobor area, (located in the northwestern part of Croatia). Using high-resolution DTM LiDAR derivatives, 874 landslide polygons were outlined, cover­ ing an area of 2.15 km2 (3.4 %). According to the BGM at a scale of 1:100.000, the study area is characterised by strong geological diversity, where 18 GUs were distinguished. This enabled extensive spatial analysis to correlate the GUs and landslide occurrence and to define the units most susceptible to landslides. Thereby, two factors crucial for landslide occurrence were considered. Namely, the stress state expressed as relief energy, and the strength of the material described by the landslide index. The analysis is based on the approach that the same average relief energy implies the same average stress state. For this reason, landslide indices were compared for GUs within the same relief energy group. The results of statistical sample analysis were also taken into account, which clearly showed the representation of each unit in the observed area. As a final result, GUs Pl,Q (gravel, sand, clay, and silt), M3 1,2 (clayey limestone and -marl), and 1M3 1 (biocalcarenite, sand­ stone, marl, silt with intercalation of sand) stood out as being the most susceptible to landslides. The analysis has also shown that older GUs, Pc (conglomerate, marl, siltite, and sandstone) and K1,2 (sandstone, shale, radiolarite and chert), are prone to landslides. Still, according to lower LOCL grades, which are given to each landslide based on the visibility and sharpness of landslide fea­ tures, the landslides within those units can be considered as old and inactive. Analysis of landslide inventory together with the field prospection revealed some significant limitations of the BGM when considering their usage in LSM design. In the paper, spe­ cial attention was given to the limitations related to the scale, geological unit defining approach, and neglecting the regolith. However, the analysis also confirmed that the BGM has an im­ portant application in small­scale landslide susceptibility map­ ping, especially as a valuable set of maps that cover the whole territory of Croatia based on equal methodology. It is very diffi­ cult to expect that in the near future a superficial map for Croatia will be developed. Therefore, to try to reduce the negative im­ pacts of the BGM limitations, it seems logical to propose the in­ structions to attribute EG features to GUs outlined on the BGM. Using such instructions in the process of small-scale LSM deri­ vation will increase the accuracy of assessed landslide suscepti­ bility, which will encourage the use of LSMs, especially in the spatial planning sector. Table 3. Basic engineering geological features of the Badenian lithofacies forma- tion gained from preliminary fieldwork. EG unit Bedrock Regolith EG unit 1 Highly permeable with low potential for regolith formation Up to 1 m thick, large proportion of coarse-grained fraction, highly water-permeable EG unit 2 Very poorly permeable with very high regolith formation potential Up to 5 m thick, very low coarse-grained fraction, low water permeability Figure 14. Landslide indices for the Badenian polygons (2M2 2) shown within Fig. 8. G eo lo gi a C ro at ic a 66 Geologia Croatica 75/1 ACKNOWLEDGMENT LiDAR-based digital terrain model used within this study was obtained in the framework of the safEarth project, co-financed by the Interreg IPA Cross-border Cooperation Programme Croa­ tia - Bosnia and Herzegovina - Montenegro 2014-2020. REFERENCES ALLABY, A. & ALLABY, M. (2003): A Dictionary of Earth Sciences.– Oxford univer­ sity press, Oxford, 619 p. ARDIZZONE, F., CARDINALI, M., GALLI, M., GUZZETTI, F. & REICHENBACH, P. (2007): Identification and mapping of recent rainfall-induced landslides using elevation data collected by airbone Lidar.– Natural Hazards and Earth System Sci­ ences, 7, 637–650. doi: 10.5194/nhess-7-637-2007 BARNES, J.W. & LISLE, R.J. (2004): Basic Geological Mapping.– John Wiley & Sons Ltd, Chichester, 184 p. BARTELLETTI, C., GIANNECCHINI, R., D’AMATO AVANZI, G., GALANTI, Y. & MAZZALI, A. (2017): The influence of geological–morphological and land use settings on shallow landslides in the Pogliaschina T. basin (northern Apennines, Italy).– Journal of Maps, 13/2, 142–152. doi: 10.1080/17445647.2017.1279082 BELL, R., PETSCHKO, H., RÖHRS, M. & DIX, A. (2012): Assessment of landslide age, landslide persistence and human impact using airborne laser scanning digital terrain models.– Geografiska Annaler: Series A, Physical Geography, 94, 135–156. doi: 10.1111/j.1468‐0459.2012.00454.x BOSTJANČIĆ, I., FILIPOVIĆ, M., GULAM, V & POLLAK, D. (2021): Regional-Scale Landslide Susceptibility Mapping Using Limited LiDAR-Based Landslide Inven­ tories for Sisak-Moslavina County, Croatia.– Sustainability, 13, 4543. doi: 10.3390/ su13084543 BURNS, W.J. & MADIN, I.P. (2009): Protocol for inventory mapping of landslide de­ posits from Light Detection and Ranging (LIDAR) imagery.– Special paper 42, Oregon Department of Geology and Mineral Industries, 30 p. CHACON, J., IRIGARAY, C., FERNANDEZ, T. & EL HAMDOUNI, R. (2006): Engi­ neering geology maps: landslides and geographical information systems.– Bulletin of Engineering Geology and the Environment, 65, 341–411. doi: 10.1007/s10064- 006-0064-z CRED, UNISDR (2018): Economic Losses, Poverty & Disasters 1998-2017.– Available at https://www.preventionweb.net/files/61119_credeconomiclosses.pdf [Accesed 20.07.2021.] ĆIRIĆ, B. (1967): Neotektonska karta SFR Jugoslavije 1:1500000 [Neotectonic map of SFR of Yugoslavia 1:1500000 – in Croatian].– Federal geological institute, Belgrade. DAI, F.C., LEE, C.F. & NGAI, Y.Y. (2002): Landslide risk assessment and management: an overview.– Engineering Geology 64, 65–87. doi: 10.1016/S0013- 7952(01)00093-X ERSHOV, V.V., NOVIKOV, A.A. & POPOVA, G.B. (1988): Fundamentals of geology.– Mir Publishers, Moscow, 358 p. FELL, R., COROMINAS, J., BONNARD, C., CASCINI, L., LEROI, E. & SAVAGE, W.Z. (2008): Guidelines for landslide susceptibility, hazard and risk zoning for land use planning.– Engineering geology, 102/3–4, 85–98. doi: 10.1016/j.eng­ geo.2008.03.022 GUZZETTI, F., CESARE MONDINI, A., CARDINALI, M., FIORUCCI, F., SANTAN­ GELO, M. & CHANG, K.T. (2012): Landslide inventory maps: New tools for an old problem.– Earth-Science Reviews, 112/1–2, 42–66. doi: 10.1016/j.earsci­ rev.2012.02.001 HEĆIMOVIĆ, I. (2009): Klastične naslage (pliokvartar - Pl, Q) [Clastic deposits (Plio Quaternary - Pl, Q), In Croatian]. In: VELIĆ, I. & VLAHOVIĆ, I. (eds.): Tumač Geološke karte Republike Hrvatske 1:300000, Croatian geological survey, Zagreb, 141, 94–95. HUGGETT, R.J. (2017): Fundamentals of Geomorphology.– Routledge, Abingdon, 579 p. INDRAWAN, I.G.B., RAHARDJO, H. & LEONG, E.C. (2006): Effects of coarse- grained materials on properties of residual soil.– Engineering Geology, 82, 154– 164. doi: 10.1016/j.enggeo.2005.10.003 JABOYEDOFF, M., OPPIKOFER, T., ABELLÁN, A., DERRON, M.-H., LOYE, A., METZGER, R. & PEDRAZZINI, A. (2012): Use of LIDAR in landslide investi­ gations: a review.– Natural Hazards, 61, 5–28. doi: 10.1007/s11069-010-9634-2 JEBUR, M.N., PRADHAN, B. & TEHRANY, M.S. (2014): Optimization of landslide conditioning factors using very high-resolution airborne laser scanning (LiDAR) data at catchment scale.– Remote Sensing of Environment, 152, 150–165. doi: 10.1016/j.rse.2014.05.013 JURAK, V., MATKOVIĆ, I., MIKLIN, Ž. & CVIJANOVIĆ, D. (1998): Landslide haz­ ard in the Medvednica submountain area under dynamic conditions.– In: MARIĆ, B., LISAC, Z. & SZAVITS-NOSSAN, A. (eds.): Geotechnical hazards, XIth Dan­ ube-European conference on soil mechanics and geotechnical engineering, Poreč, Croatia. KEAREY, P. (2001): The new penguin dictionary of geology.– Penguin group, London, 327 p. LI, Y., CHEN, G., TANG, C., ZHOU, G. & ZHENG, L. (2012): Rainfall and earthquake- induced landslide susceptibility assessment using GIS and Artificial Neural Net­ work.– Natural Hazards and Earth System Sciences, 12, 2719–2729. doi: 10.5194/ nhess-12-2719-2012 LUKIĆ, T., BJELAJAC, D., FITZSIMMONS, K.E., MARKOVIĆ, S.B., BASARIN, B., MLAĐAN, D., MICIĆ, T., SCHAETZL, R.J., GAVRILOV, M.B., MILANOVIĆ, M., SIPOS, G., MEZŐSI, G., KNEŽEVIĆ-LUKIĆ, N., MILINČIĆ, M., LÉTAL, A. & SAMARDŽIĆ, I. (2018): Factors triggering landslide occurrence on the Ze­ mun loess plateau, Belgrade area, Serbia.– Environmental Earth Sciences, 77, 519. doi: 10.1007/s12665-018-7712-z MCCALPIN, J. (1984): Preliminary age classification of landslides for inventory map­ ping.– Proceedings of the Annual Symposium on Engineering Geology and Soil Engineering, 21, Moscow, Idaho, 99–120. MIHALIĆ, S., OŠTRIĆ, M. & VUJNOVIĆ, T. (2008): Landslide susceptibility mapping in the Starca Basin (Croatia).– Proceedings of 2nd European Conference of Inter­ national Association for Engineering Geology, 2008, Madrid, Spain. NOVOSEL-ŠKORIĆ, S., MAGDALNEIĆ, Z. & KOROLIJA, B. (1986): Pliokvartarni sedimeneti jugoistočnog dijela Krndije – Sjeverna Hrvatska [Plioquaternary sedi­ ments of the southeastern part of Krndija - Northern Croatia - in Croatian].– Geološki vijesnik, 39, 11–18. PARK, T.-W., KIM, H.-J., TANVIR, M.T., LEE, J.-B. & MOON, S.-G. (2018): Influence of coarse particles on the physical properties and quic undrained shear strength of fine-grained soils.– Geomechanics and Engineering, 14/1, 99–105. doi: 10.12989/ gae.2018.14.1.099 PERKINS, C.R. & PARRY, R.B. (1990): Information Sources in Cartography - Guides to Information Sources.– Bowker-Sour, London. PIKE (1988): The geometric signature: quantifying landslide - terrain types from digital elevation models.– Mathematical geology, 20/5, 491–511. REICHENBACH, P., ROSSI, M., MALAMUD, B.D., MIHIR, M. & GUZZETI, F. (2018): A review of statistically-based landslide susceptibility models.– Earth- Science Reviews, 180, 60–91. doi: 10.1016/j.earscirev.2018.03.001 ROSER, M., RITCHIE, H. & ORTIZ-OSPINA, E. (2013): World Population Growth. –Published online at OurWorldInData.org. Retrieved from: ‘https://ourworldin­ data.org/world-population-growth’ [Accessed 13.07.2021.] SAVEZNI GEOLOŠKI ZAVOD (1962): Uputstvo za izradu osnovne geološke karte FNRJ [Instructions for making a basic geological map – in Croatian].– Savezni geološki zavod, Beograd, 99 p. SCHUSTER, R.L. (1996): Socioeconomic significance of landslides.– In: TURNER, A.K. & SCHUSTER, R.L. (eds.): Landslides: Investigation and Mitigation.– Spe­ cial Report 247, Transportation Research Board, National Research Council, Na­ tional Academy Press, Washington, DC, 12–35. SHAKOOR, A. & COOK, B.D. (1990): The effect of stone content, size, and shape on the engineering properties of a compacted silty clay.– Environmental and Engi­ neering Geoscience, XXVII/2, 245–253. doi: 10.2113/gseegeosci.xxvii.2.245 SHAN, J. & TOTH, C.K. (2009): Topographic laser ranging and scanning - Principles and processing.– Taylor & Francis Group, LLC, New York, 654 p. ŠIKIĆ, K. (1995): Prikaz geološke građe Medvednice.– In: ŠIKIĆ, K. (ed.): Geološki vodič Medvednice. Institut za geološka istraživanja, INA-Industrija nafte, Zagreb, 7–30. ŠIKIĆ, K., BASCH, O. & ŠIMUNIĆ, A. (1977): Osnovna geološka karta RH 1:100000, list Zagreb L38−80 [Basic Geological Map of SFRY 1:100000, Zagreb sheet – in Croatian].– Geološki zavod, Zagreb, Savezni geološki zavod, Beograd. ŠIKIĆ, K., BASCH, O. & ŠIMUNIĆ, A. (1979): Osnovna geološka karta SFRJ 1:100000. Tumač za list Zagreb L38−80 [Basic Geological Map of SFRY 1:100000, Geology of the Zagreb sheet – in Croatian].– Institut za geološka istraživanja, Zagreb, Savezni geološki zavod, Beograd, 75 p. UPPAL, S. (2006): Fundamentals of Physical Geography.– National Council of Educa­ tional Research and Training, New Delhi, 136 p. VAN DEN EECKHAUT, M., POESEN, J., VERSTRAETEN, G., VANACKER, V., NYSSEN, J., MOEYERSONS, J., VAN BEEK, L.P.H. & VANDEKERCKHOVE, L. (2007): Use of LIDAR-derived images for mapping old landslides under for­ est.– Earth Surface Processes and Landforms, 32, 754–769. doi: 10.1002/esp.1417 VARNES, D.J. & IAEG (1984): Landslide hazard zonation: a review of principles and practice.– IAEG Commission on Landslides and Other Mass Movements, Natural Hazards 3, UNESCO, Paris, 63 p. VRSALJKO, D., PAVELIĆ, D., MIKNIĆ, M., BRKIĆ, M., KOVAČIĆ, M., HEĆIMOVIĆ, I., HAJEK-TADESSE, V., AVANIĆ, R. & KURTANJEK, N. (2006): Middle Miocene (Upper Badenian/Sarmatian) Palaeoecology and Evolu­ tion of the Environments in the Area of Medvednica Mt. (North Croatia).– Geolo­ gia Croatica, 59/1, 51–63. VRSALJKO, D., HEĆIMOVIĆ, I. & AVANIĆ, R. (2007): Miocene deposits of Northern Croatia.– In: GRGASOVIĆ, T. & VLAHOVIĆ, I. (eds.): 9th International Sympo­ sium on Fossil Algae Croatia, Field Trip Guidebook and Abstracts, Croatian Geo­ logical Survey, Zagreb, 143–153. ZEZERE, J.L., DE BRUM FERREIRA, A. & RODRIGUES, M.L. (1999): The role of conditioning and triggering factors in the occurrence of landslides: a case study in the area north of Lisbon (Portugal).– Geomorphology, 30, 133–146. doi: 10.1016/ S0169-555X(99)00050-1