Interdisciplinary geoscientific approach to radioactive waste repository site selection: the Březový potok site, southwestern Czech Republic

ABSTRACT Interpretation of the deep geological and tectonic setting is a key step for safety assessment of the potential sites for a deep geological repository of radioactive waste. However, the conventional geological mapping is difficult in the agricultural landscape with lacking rock outcrops. To address this task, we combined geological and hydrogeological mapping, morphostructural analysis and geophysical methods to construct a new geological map and a 3D model of the Březový potok site in the southwestern Czech Republic. They provide unique constraints for the understanding of the subsurface geological and tectonic architecture. The results are compared with the archive geological data and the usefulness of individual geophysical methods to detect the lithological contacts and tectonic zones is discussed. We present this multi-disciplinary research strategy as a workflow for detailed geological characterization, which can be successfully applied on other sites with similar geological environment.


Introduction
The site selection for the deep geological repository is a subject of extensive research in many disciplines worldwide (e.g. Behlau and Mingerzahn, 2001;Bukovská et al., 2019;Liebscher et al., 2020;Tirén et al., 2001;Wang, 2014). The geological and tectonic characteristics at depths of up to 1 km are crucial for the safety assessment of the potential sites. However, little is known about such a subsurface setting and its characterization is a challenging task confronting geoscience community. Geological maps are fundamental for the visualization of geological data and for reconstruction of the subsurface 3D configuration.
Conventional geological mapping of basement rocks is strongly limited by sedimentary cover and vegetation and often resulting in subjective and variably realistic outcomes. Modern geological maps require a more complex approach by the acquisition of data from both the surface and depth. However, only limited possibilities of indirect and non-invasive techniques exist. A combination of geological and hydrogeological mapping, morphostructural analysis and geophysical surveys is a promising approach. Both the morphostructural analysis of relief (e.g. Burbank and Anderson, 2001;Jelínek, 2008;Ritter et al., 2002) and geophysical measurements (e.g. Demanet et al., 2001;Ganerød et al., 2006;IAEA, 2003) are often used to trace the unexposed geological contacts and tectonic features, and can complement to traditional surface mapping. Moreover, such integrated survey strategy can provide sufficient data to understand the subsurface architecture of the study area and to construct reliable 3D geological modelling (e.g. Caumon et al., 2009;Thornton et al., 2018;Zanchi et al., 2009), which is a powerful tool to approximate and to visualize the deep geological and tectonic features. It can further serve as a geological background for subsequent scientific activities on the site.
The objective of this study is to present new geological, hydrogeological, morphostructural and geophysical data used to construct geological map and 3D model of the Březový potok potential site ( Figure  1). The paper is based on the results of the 'Geophysical works to describe the geological structures of potential DGR sites in the Czech Republic' project (Mixa et al., 2020) realized as a part of the long-term national strategy of the Czech Republic for the siting process of radioactive waste repository (Radioactive Waste Repository Authority; www.surao.cz). In the project, a total of nine sites (Figure 1(b)) were investigated by the identical methodology. The Březový potok site is used as an example to present results and introduce an applied multi-disciplinary approach. It is also a unique opportunity to test the suitability of various geophysical methods for the mapping of subsurface geological structures.

Location
The studied area is located c. 100 km southwest of Prague and belongs to the southwestern part of the Bohemian Massif (Figure 1; e.g. Matte, 1990;Schulmann et al., 2009;Soejono et al., 2020). The locality covers the contact of two important lithotectonic units (Figure 1(c), Main Map): (1) the high-grade metasediments of Moldanubicum in the north (Košler et al., 2014;Vrána et al., 1995) and (2) the igneous rocks of the Central Bohemian Plutonic Complex in the south (CBPC; Holub et al., 1997).

Approach and methods
Geological and hydrogeological mapping and geophysical surveys were performed along the 14 profiles of a total length of 73.5 km ( Figure 2). The profiles were positioned to determine the location of tectonic zones and lithological contacts in the site area. The surveys were applied to detect and verify the rock types and especially the tectonic zones. In addition, the morphostructural analysis of the digital elevation model (DEM) was used to extract natural relief features (morpholineaments) genetically linked to the lithology and tectonics. Geophysical surveys included the following methods: resistivity profiling, electric resistivity tomography, gravity method, seismic refraction and reflection and gamma-  ray spectrometry. Moreover, electromagnetic, vertical electric sounding and magnetic data were acquired along with some parts of the profiles but these methods did not provide applicable results and will not be further discussed. Full methodological details of the geological and hydrogeological mapping, used geophysical methods, morphostructural analysis and 3D modelling are listed in Appendix 1. Spatial distribution of geological and hydrogeological indicators, morphostructural lineaments and geophysical anomalies are shown in Figures 3, 4 and 5 of Major Map. All the geophysical data are shown in the Appendix 2 and accessible on the SÚRAO's website (www.surao.cz), the 3D PDF of the 3D model is given in the Appendix 3.

Geological and hydrogeological mapping
In the north, the Moldanubicum is composed of the E-W trending and steeply dipping migmatized paragneiss and quartzitic gneiss containing lenses of calcsilicate and amphibolite. Few irregular bodies of paragneiss and biotite orthogneiss up to km-scale isolated within the granodiorite were detected in the southern margin of the locality. The CBPC is represented by the Blatná granodiorite (medium-grained amphibolebearing biotite granodiorite) dominating in the central part and by the Červená granodiorite (weakly porphyritic fineto medium-grained amphibole─biotite granodiorite) prevailing in the southern part. The direct contact between these two types is not exposed. However, the field observations suggest that their contacts are gradational and diffuse. This is in agreement with previous geochronological and geochemical studies, which interpreted the Červená type as a marginal mafic variety of the Blatná type (Holub et al., 1997;Janoušek et al., 2010). Moreover, several smallscale bodies of the Polánec granite (medium-to coarse-grained biotite granite) intruded into the metasediments of the Moldanubicum in the northeastern part of the area. The contacts of granodiorites with metasediments have not been directly observed, however, it was well localized by fragment mapping and their intrusive nature was reported out of the locality (Žák et al., 2012). High-grade rocks of the Moldanubicum underwent two-stage tectonic evolution recorded by two generations of metamorphic fabrics, subvertical E-W trending metamorphic foliation heterogeneously refolded into the younger subhorizontal deformational fabric. Such a structural pattern is typical for most of the Moldanubicum (e.g. Bukovská et al., 2019). In contrast, both the granodiorite types of the CBPC locally show only weak magmatic to solid-state fabric (Žák et al., 2012).
Furthermore, both the granodiorites and metasediments are cross-cut by numerous E-W trending dykes of aplite, lamprophyre, muscovite and biotite granites. The site has a strongly denuded character and the Quarternary deposits consist mostly of fluvial sand and gravel and alluvial sandy clay and gravel infilling shallow valleys. Several bodies of colluvial sandy clay with rounded to angular clasts and blocks cover the gentle hillslopes in the northern part. The youngest sequence is represented by the fishpond organic-rich sandy loam and marshy sediments, widespread in the northern and central parts. No map-scale fault was detected on the surface due to the poor exposure conditions and only few mesoscopic faults have been observed (Figure 3(a)). Fault zones are indirectly indicated by fragments of hydrothermal quartz fault infill (Figure 3(b)) and the rock alteration. Occurrence of these features is interpreted as spatially associated with assumed tectonic zones (Figure 3 of Major Map).
Hydrogeological mapping along the profiles revealed a total of ninety features including springs, natural outflows, seeps, wetland meadows and wells ( Figure 3(c,d)) that are relatively homogeneously distributed throughout the study area ( Figure 3 of Major Map). Such objects are considered as surface manifestations of geological structures such as rock contacts or tectonic zones (e.g. Bense et al., 2013;Singhal and Gupta, 2010).

Morphostructural analysis
The morphostructural analysis of DEM shows the complex geomorphological pattern defined by three groups of lineaments (green lines in Figure 4 of Main Map): the dominant NW-SE oriented lineaments accompanied by minor NNW-SSE and E-W trending lineaments. The granodiorites of the central part display slightly higher density of lineaments compared to the northern metasedimentary counterpart. This difference is probably related to a more regular fracture system typical for granitoids (Scheidegger, 2001). The morphostructural analysis revealed two categories of localized lineament zones (blue stripes in Figure 4 of Main Map) consisting of subparallel interconnected lineaments. Such zones are classically considered as expressions of large-scale faults or important mechanical boundaries (Burbank and Anderson, 2001;Ganas et al., 2005;Jelínek, 2008). The major zones form two broad NW-SE trending belts running across the entire area in its central part and three NW-SE oriented features in the southeast. In the granodiorites, the NNW-SSE trending zones interconnect the major NW-SE zones while the E-W and NNE-SSW-oriented zones are less important. On the contrary, such E-W and NNE-SSW trending features are dominant in the northern metasediments. The resulting morphostructural network reflects lithological, structural and surface elements and provides one of the data sources for the geological map construction ( Figure 4).

Geophysical surveys
In case of lack of basement outcrops and homogenous lithologies, the appropriate geophysical methods can be used to identify the subsurface discontinuities. Thus, resistivity, seismic and gravity data help to distinguish small variations of the physical properties of rocks and determine the contacts between the geological units, tectonic zones and thickness of the sedimentary cover (Appendix 2). The resistivity profiling was performed all along each profile and combined with  the electric resistivity tomography in some key locations (8 km of data acquisition). These methods are appropriate to differentiate the basement from the sedimentary cover, as solid rocks are characterized by a higher resistivity than sedimentary and weathered rocks. Thus, the variation of resistivity reveals up to 20 m thick sedimentary cover and the terminations of the tectonic zones in the underlying basement rocks but the apparent resistivity cannot detect the contact between the rocks with similar mineralogy (Figures 5(a) and 6(b)). The seismic surveys consist of 34 km refraction profiles and 8 km reflection profiles. It refines the interpretations based on resistivity data as most of the discontinuities in the seismic refraction sections correspond to the tectonic zones detected by the resistivity methods ( Figures 5 and 6). In some cases, the orientation of the tectonic zones can be inferred from the dip of zones of decreased seismic velocities (Figures 5(b) and 6(b); Figure 8 of Main Map). Seismic refraction profiles also confirm the thickness variation of the sedimentary cover. Seismic reflection reveals subvertical and subhorizontal discontinuities at deeper levels than resistivity methods and seismic refraction (Appendix 2). The subhorizontal discontinuities (at a depth of c. 100-300 m below surface) can be interpreted as an abrupt increase in rock-mass quality related to the decrease of fracturing and weathering effects with depth (e.g. Barton, 2006). Gravity data were acquired at the length of 30 km and combined with the existing complete Bouguer anomaly map. The NNW-SSE gradient in the north of this map (Figure 7 of Major Map) highlights the major lithological boundary of the site, contact of the CBPC granodiorites with the high-grade metamorphic rocks of the Moldanubicum. This gradient is also observed in the new acquired gravity profiles when crossing this contact. Otherwise, the slight localized minima on the gravity profiles often correspond to the discontinuities delineated by resistivity and seismic methods ( Figure 6). These decreases in the gravity signal are interpreted as tectonic zones located in the basement or sometimes as dykes. The detection of aplite and pegmatite dykes could not be achieved by the applied geophysical methods, due to their similar petrophysical properties with the hostrocks. The lamprophyre dykes have been well localized by the gamma-ray spectrometry based on high U and Th and low K contents.

Integration of the results into the tectonic zone network, geological map and 3D model
We combined all the results and constructed a new geological map at a scale of 1:25,000 (Main Map) with a special focus on the tectonic zone network. The map shows the graphical interpretation of geological units, rock bodies, the location and characteristics of their contacts and trends, and the relationships between major fault zones. Four lithologically defined units were established: (1) the high-grade metamorphic complexes of the Moldanubicum, (2) the magmatic rocks of the CBPC, (3) the dykes and (4) the sedimentary cover, containing sixteen distinct rock types. The geological map with tectonic zone network (Main Map) has been compiled based on the analysis of the geophysical anomalies, surface mapping and spatial distribution of the tectonic zone markers such as fault infill fragments, hydrogeological features and DEM lineaments. All the interpreted tectonic zones were subsequently classified by supposed width and length after Andersson et al. (2000).
Using the new map, an interpretative 3D geological model was constructed to a depth of 1.5 km below the surface (Figure 7; Figure 9 of Main Map; Appendix 3). The 2D surface orientation and shape of the rock bodies and faults were extrapolated based on the geophysical data and structural data interpretation. At the depth exceeding the resolution of the geophysical methods, the continuation of individual rock bodies and tectonic zones has been modelled continuously with the last known orientation. Uncertainties in the geological contacts and tectonic zones occurrence and geometry increase with depth. Hence, the resulting 3D model shows updated knowledge of the deep geological setting.
Both the presented geological map and 3D model represent our interpretation of all new and archive data available today from the Březový potok site. However, further detailed subsurface geological and geophysical studies are required to verify the main geological and tectonic structures of the site.

Comparison with previous maps
Comparison of the new results with archive geological maps ( Figure 8) revealed several significant differences. The original map scheme (Figure 8(a)) was compiled from varied documentations at different scales (Franěk et al., 2018). Our multi-disciplinary approach resulted in a new map (Figure 8(b)) with improved lithology contacts that include simplification of the Moldanubian metasedimentary rocks, diffusion-type intrusive contacts of two granodiorites of the CBPC, and precise location of dyke bodies. The contacts of two major regional units remain more or less identical, however, the location, shape and the nature of the contact between individual rock bodies changed significantly. An occurrence of biotite orthogneiss in the SE part of the area was newly recognized whereas the presence of muscovite-biotite orthogneiss in the NE corner was not confirmed (Figure 6). The extent and shape of the Polánec granite bodies were also specified. The dykes are more frequent and homogeneously distributed than expected. In general, the spatial distribution of faults is more uniform in comparison with older maps. Importantly, the tectonic zone network is now better balanced with respect to the brittle deformation within the rock massif and important fault zones are well

Conclusions
The presented geological map and the 3D model synthesize the newest advancements in the understanding the geological environment of the Březový potok site. We used conventional geological mapping enhanced by hydrogeological survey, morphostructural and geophysical analyses. This multi-disciplinary approach provides a significant progress in site investigation and is proposed as a complex strategy that can be easily applied to other sites with only minor adjustments.  The case study of the Březový potok demonstrates the complementarity between the resistivity, seismic and gravity methods. They respectively provided: (1) the distinction between the sedimentary cover and the basement rocks, (2) the boundary between the CBPC and Moldanubicum and (3) the tectonic zone network. Only c. 12% of the hydrogeological indicators spatially correlate with the low-resistivity anomalies. Thus, the resistivity methods are generally useful for the identification of underground flows without surface expressions (e.g. during dry periods) and can be used to verify the hydrogeological function of assumed tectonic zones. In contrast, the vertical electric sounding, electromagnetic and magnetic methods were not suitable for research of geological and tectonic structures in a given geological conditions.
The new map and the 3D model will provide detailed geological and tectonic constraints for the planning of further geological and geophysical campaign, technical works (e.g. deep-boreholes drilling), geotechnical stability models and hydrogeological flow models. All these results will serve as a crucial dataset for safety assessment and final selection of the site for the deep geological repository of radioactive waste.

Software
The ArcGIS software (ESRI) was used for editing the spatial data and map construction. The resistivity data were processed in Grapher (Golden Software) and RES2DINV, gravity data in Geotools-GravModeler (LaCoste and Romberg) seismic refraction data in Reyfract (Intelligent Resources), seismic reflection data in ReflexW (Sandmeier Software) and gammaray spectrometry data in Surfer (Golden Software) softwares. The morphostructural analysis and 3D model were performed using ArcGIS, MOVE (Petroleum Experts) and Surfer (Golden Software) softwares. All figures were created in CorelDRAW.