Multiscale porosity changes along the pro- and retrograde deformation path: an example from Alpine slates

. Estimating the porosity of slates is of great interest for the industries dealing with sub-surface areas such as CO 2 sequestration, nuclear waste disposal and shale gas but also for engineering purposes in terms of mechanical stabil-ity for underground or surface constructions. In this study, we aim at understanding estimates of the porosity of slates from the Infrahelvetic ﬂysch units (IFUs) in the Glarus Alps (east-ern Switzerland). Surface and sub-surface samples were collected along a temperature gradient from 200 to 320 (cid:14) C and therefore give the opportunity to link pore types along this temperature and deformation path. In addition, we indicate which porosity is the effect of surface processes and indicate the contribution of artiﬁcially induced porosity. The developed workﬂow consists of a combination of bulk rock measurements including helium pycnometry (He pycnometry) and mercury intrusion porosimetry (MIP) with image analysis. Image analysis was performed with high-resolution scanning electron microscopy (SEM) on broad ion beam (BIB) prepared cross sections (BIB-SEM). Different vein generations provide evidence of porosity formation at depth, as they present “paleo-porosity”. Towards peak metamorphic conditions (prograde path), porosity reduces to < 1 vol%, indicated by matrix porosity detected by BIB-SEM. During exhumation (retrograde path) porosity increases due to the formation of microfractures interpreted as the effect of unloading (open fractures). At the surface, porosity is further increased due to the formation of macro-fractures (fracture apertures up to 1 mm), which are interpreted as being either due to the effect of weathering processes such as freeze and thaw cycles or artiﬁcially induced by sample preparation. Additionally, porosity and pore morphology are strongly dependent on mineralogy, sample homogeneity and strain, which change dynami-cally in time and space.

Abstract. Estimating the porosity of slates is of great interest for the industries dealing with sub-surface areas such as CO 2 sequestration, nuclear waste disposal and shale gas but also for engineering purposes in terms of mechanical stability for underground or surface constructions. In this study, we aim at understanding estimates of the porosity of slates from the Infrahelvetic flysch units (IFUs) in the Glarus Alps (eastern Switzerland). Surface and sub-surface samples were collected along a temperature gradient from 200 to 320 • C and therefore give the opportunity to link pore types along this temperature and deformation path. In addition, we indicate which porosity is the effect of surface processes and indicate the contribution of artificially induced porosity. The developed workflow consists of a combination of bulk rock measurements including helium pycnometry (He pycnometry) and mercury intrusion porosimetry (MIP) with image analysis. Image analysis was performed with high-resolution scanning electron microscopy (SEM) on broad ion beam (BIB) prepared cross sections (BIB-SEM). Different vein generations provide evidence of porosity formation at depth, as they present "paleo-porosity". Towards peak metamorphic conditions (prograde path), porosity reduces to < 1 vol%, indicated by matrix porosity detected by BIB-SEM. During exhumation (retrograde path) porosity increases due to the formation of microfractures interpreted as the effect of unloading (open fractures). At the surface, porosity is further increased due to the formation of macro-fractures (fracture apertures up to 1 mm), which are interpreted as being either due to the effect of weathering processes such as freeze and thaw cycles or artificially induced by sample preparation. Additionally, porosity and pore morphology are strongly dependent on mineralogy, sample homogeneity and strain, which change dynamically in time and space.

Introduction
Slates as representatives of sheet-silicate-rich rocks delineate a fine-grained anisotropic microstructure and are common low-grade metamorphic rocks. These rock types play an important role in the industry of underground storage, such as nuclear waste disposal and geological carbon sequestration (e.g. Loon, 2008;Thury and Bossart, 1999) but are also of importance in the surface and sub-surface building industry because of the low mechanical strength and associated geotechnical problems (e.g. Blümling et al., 2007). Moreover, slates are crucial in mountain building processes, since they can act as a mechanically weak phase and localize large amounts of strain in low-to very low-grade metamorphic domains (e.g. Labaume et al., 1997;Milliken and Reed, 2010;Warr et al., 2014). Slates also act as important fluid sources, liberating water with progressive compaction and metamor-I. V. Akker et al.: Multiscale porosity changes along the pro-and retrograde deformation path phic reaction (e.g. Dielforder et al., 2015;Rieke and Chilingarian, 1974). The pore network defines the major fluid pathways, and a detailed investigation of microstructural porosity is a key to better understanding the fluid flux and circulation of fluids in collisional orogens and orogenic wedges.
The porosity of a rock develops along the deformation path and with changing pressure and temperature (PT) conditions. A possible evolution starts with early compaction and diagenesis on the prograde path, evolving to maximum burial, followed by metamorphism and possibly tectonic deformation on the retrograde path, eventually resulting in surface exposure (Horseman et al., 1996). Porosity is produced by processes such as grain dilation, grain boundary sliding, cavity formation and dissolution (Herwegh and Jenni, 2001;Fusseis et al., 2009;Gilgannon et al., 2017), and porosity is sealed by grain size reduction, compaction and precipitation, as observed for example in clay fault gauges (Holland and Urai, 2006;Laurich et al., 2014).
The first aim of this study is to give an estimate of the total present-day porosity of a slate which underwent this full deformation-metamorphic cycle. The second goal is to establish what causes the present-day total porosity and which evolutionary processes account for which type of porosity.
For this study, we focus on surface-collected slates from the Infrahelvetic flysch units (IFUs) in the Glarus Alps (eastern Switzerland), which are affected by large out-ofsequence thrusts, such as the Glarus thrust (e.g. Badertscher et al., 2002;Burkhard et al., 1992;Ebert et al., 2007;Herwegh et al., 2008;Pfiffner et al., 2011;Poulet et al., 2014;Rahn and Grasemann, 1999;von Daeniken and Frehner, 2017). This specific area is of interest because it covers a temperature gradient, which positively correlates with background strain. Tectonic processes change along this gradient from soft sediment behaviour including particulate flow, at the lowest temperatures, to pressure solution at higher temperatures (Dielforder et al., 2016). Investigating porosity over this gradient would give an opportunity to link different pore types to such specific tectonic processes.
A variety of methods have been developed over the past decades to determine porosity, pore morphology and pore size distribution, each having a characteristic pore size range (Anovitz and Cole, 2015;Busch et al., 2017). Some common direct approaches include gas expansion techniques such as helium pycnometry (He pycnometry) and mercury intrusion porosimetry (MIP) (Anovitz and Cole, 2015;Tiab and Donaldson, 2015). These techniques both yield information on the bulk-interconnected porosity in 3-D volumes. Recent developments of surface-preparation-based ion beam milling (focused ion beam, FIB, and broad ion beam, BIB) allow imaging porosity down to nanometre scale on an almost perfectly flat surface with scanning electron microscopy (SEM). These images are suitable for high-quality image segmentation and the subsequent quantification of porosity in 2-D (e.g. Keller et al., 2013b;Klaver et al., 2012).
Porosity estimates obtained with bulk rock measurements yield information on large sample volumes but remain unresolved in terms of the contribution of different pore morphologies. Therefore, most studies use a combination of bulk rock measurements with image analysis (using He pycnometry: e.g. Houben et al., 2016;Yang et al., 2016;MIP: e.g. Hemes et al., 2013;Klaver et al., 2012Klaver et al., , 2015; MIP, spectrally induced polarization (SIP), nuclear magnetic resonance (NMR), micro-computed tomography (µ-CT): Zhang et al., 2018). Most of this research is focused on shales or clay(stones). In this study, we applied these two complimentary techniques to slates. In contrast to shales, slates show an enhanced cohesion and higher rock strength and exhibit a more anisotropic microstructure due to increased metamorphism, deformation and dehydration. Particularly because of the latter, (micro)fractures play an important role as fluid pathways in these slate samples. We aim to find out if such fractures result only from surface processes or form already at depth. In addition, there is a strong contrast in porosity between samples collected from sub-surface drill cores and surface outcrops in terms of the effect of unloading due to stress changes during exhumation. Finally, it is necessary to distinguish original porosity from artificially induced porosity, as the porosity of the delicate slate samples is likely affected by sample collection and preparation.

Samples and geological background
Slates were selected from the IFU in the Glarus Alps in eastern Switzerland (Fig. 1a). The flysch units are deposited in the underfilled Northern Alpine Foreland Basin and accreted in sequence to the accretionary wedge during the subduction of the European margin (Lihou, 1996). The IFU consists of the North Helvetic Flysch (NHF), which is overthrusted by the South Helvetic Flysch and Ultrahelvetic Flysch, the latter including the Sardona and Blattengrat units. The names of the flysch units refer directly to the paleo-graphic realm in which they were deposited (Lihou, 1996). The Ultrahelvetic Flysch is itself overthrusted by the Subhelvetic units and the Helvetic nappes along the Glarus thrust (Pfiffner, 1986;Schmid, 1975;Trümpy, 1969). The lower greenschistfacies metamorphic conditions in the thrust wedge relate to an overburden of 8-12 km (Frey et al., 1980;Lihou, 1996;Rahn et al., 1995).
The flysch units include heterogeneous endmembers from volcano-clastics, turbiditic sandstones and dark slates (Fig. 1b;Siegenthaler, 1974). These rocks show indications of soft sediment deformation (Dielforder et al., 2016), such as folding of bedding-parallel calcite veins. Late Cretaceous limestones and marls are indicative of the Sardona and Blattengrat units, which occur towards the south of the study area (Fig. 1c;Bisig, 1957;Leupold, 1942;Lihou 1995;Trümpy, 1969  The compositions of the slates from the flysch units are heterogeneous, comprising layered mineralogical variations owing to turbiditic sedimentation processes. In this study, we focus on the clay-(mica)-rich layers, excluding sandy layers from our analysis. We collected eight samples from the flysch between the village Weesen in the north and the locality Flims in the south (Fig. 1a). All samples are from surface outcrops, except for samples D and E, which are from an abandoned sub-surface slate mine (Landesplattenberg). The N-S sample profile across the Flysch units covers a change in peak metamorphic conditions from 200 • C in the north to 320 • C in the south determined by calcite-dolomite thermometry and Raman spectroscopy of carbonaceous material (Ebert et al., 2007;Lahfid et al., 2010). This gradient cor-relates to a strain gradient, as the most southern units were subducted to the deepest levels prior to in-sequence imbricate stacking. He pycnometry is a gas expansion technique based on Boyle's law to measure volume (Anovitz and Cole, 2015). He pycnometry was performed on oven-dried (115 • C) plugs with a 2.5 cm diameter and geometrical volumes > 5 cm 3 (Fig. 2). Plugs were drilled in the hand specimens (Fig. 3a), and the flat surfaces of the plugs were polished to obtain a perfectly cylindrical shape. The geometrical volume was determined by an average of six calliper measurements of the height and four measurements of the diameter of the plug. Pycnometric density measurements were obtained using an AccuPyc 1330 pycnometer. The accuracy of the measurements is related to sample quantity and cell size. The accuracy of the porosity measurements of the machine falls between 0.02 and 0.1 vol% (Viana et al., 2002). Helium was used as a specific gas as its small atom radius was assumed to percolate the sample volume down to the smallest pore sizes. The drilling of plugs resulted in some cases in the splitting of the rock along the foliation planes, resulting in several thin cylindrical disks ( Fig. 3b and c). In the cases where the volume of one disk was < 10 % of the sample chamber, several disks of adjacent sample material were analysed simultaneously. This was done to increase the total sample quantity and reduce the analytical error. When the plug volume exceeded 10 vol% of the sample volume chamber, an average of two adjacent plugs was taken. The amount of purges was set to 10 and the connected porosity was calculated using Eq. (1). connected porosity = 1 − bulk density grain density (1) For He pycnometry, the sample volume has shown to effect the measurements, with larger sample volumes resulting in lower porosity (Table S1 in the Supplement). Duplicate measurements yield absolute differences of about 0.03-0.07 vol%. The difference between measuring different layers and several cylinders simultaneously, compared to an average of separate layers, yields a maximum difference of about 0.5 %. Therefore, the error, used as the standard deviation of the different measurements, was set to 0.5 vol%.

Mercury intrusion porosimetry
MIP was applied to oven-dried (115 • C) sample volumes of clay(mica)-rich layers taken from the close vicinity of the plugs used for He-pycnometry measurements. Sample material was cut in wet conditions, using water, by means of a  diamond saw so that the measured pieces had a 4-5 mm edge length. In one measurement, three to five pieces were analysed. The analyses were performed using Thermoelectron Pascal 140/440 equipment. In a first step, mercury was intruded to a maximum pressure of 200 kPa. Then the sample was moved to the Pascal 440 machine for high-pressure analysis up to 200 MPa. Pore sizes were calculated using Washburn's equation (Washburn, 1921). For the analyses, a cylindric pore shape with a contact angle of 140 • and a surface tension of mercury of 0.48 N m −1 porosimetry was assumed, leading to intruding a minimal pore radius of 3.7 nm.
MIP was performed in a multi-cycle experiment, including two full injection and extrusion cycles. It is assumed that "ink-bottle type pores" remain filled after the first intrusion and that in a second intrusion the pores are successively filled with mercury (e.g. big pores followed by smaller pores). The first cycle indicates the complete pore volume and the second cycle only the connected pores, which include the small pore fraction (Fig. 4). Subtraction of the second cycle from the first cycle gives the size of the neck pore entrances at which the ink-bottle type pores become filled (Kaufmann, 2010;Kaufmann et al., 2009). All errors of this technique fall into the 5 vol% range.

High-resolution imaging techniques
A regular petrographic light microscope was used to examine the petrography and microstructure of the samples. In addition, a Virtual Petrography (ViP) dataset was created from a subset of thin sections at the Structural geology, Tectonics and Geomechanics (GED) institute, RWTH Aachen, using an automated petrographic microscope, which under various polarization and illumination conditions (parallel polarized and crossed polarized) digitizes complete thin sections at high resolution (e.g. 45 000 × 30 000 pixel maps for a thin section of 3 × 2 cm; Virgo et al., 2016). The Petroscan Tile-Viewer software, which allows changing the polarization angle digitally from the interpolated layers, was used to interrogate the dataset. Such an environment works as a "digital microscope" in which it is possible to directly select microstruc-  Figure 4. Conceptual model of first and second MIP cycle. Numbers indicate the succession in pore filling, which is related to a stepwise increase in pressure, linked to a certain radius of a pore neck entrance. The second cycle indicates the remaining mercury from the first cycle.
tural domains of interest. For this study, we selected only areas from clay-(mica)-rich intervals and excluded large-scale fractures ( Fig. 3d and e).

Ion polishing sample preparation
To examine the differences between mechanical polishing and ion polishing of thin sections in terms of induced artificial porosity, such as fractures, first mechanically polished thin sections were produced and examined. This was followed by surface low-angle ion polishing performed on the 2-D surfaces of such thin sections (circular areas with a diameter of 2.6 cm). Ion polishing was performed with a Leica EM TIC 3X argon ion stand-alone polisher (6 kV, 2 mA, 3 • , 8 h) at the Institute of Geological Sciences, University of Bern.
Additionally, slope cut ion polishing was performed on the sample material. Samples were treated by two different polishing techniques: 2-D surfaces on thin sections and slope cutting on sample material. These samples were cut from clay-(mica)-rich layers with a diamond saw under wet conditions to dimensions of about 5 × 7 × 3 mm. Cross sections were made in the 7 × 3 mm side and were 1-2 mm 2 in size. Slope cutting was performed with a JEOL SM-09010 BIB argon cross-section polisher (6 kV, 90 • , 8 h) at the GED, RWTH Aachen. The BIB cross sections were coated with tungsten. The slope cutting polishes smaller areas but does so with very high quality, which allows microstructural imaging down to the nanometre scale (e.g. Hemes et al., 2015).
SEM imaging using samples prepared with different techniques is shown in Fig. 5. The mechanically polished thin section shows many polishing damage effects such as breakouts of pores. These effects are not seen in either of the ionpolished samples (ion-polished 2-D thin section and slope cutting polishing (BIB-SEM)). The advantage of the ion- cal polishing because there is (1) no smearing of clay, (2) no breakage of mineral grains and (3) no polishing dust, which could fill pores.

Microstructure imaging
A Zeiss EVO 50 SEM with backscatter and secondary electron detectors at the Institute of Geological Sciences, University of Bern, was used to obtain images from the mechanically polished thin sections. Image mosaics (0.4-0.6 mm 2 ) consisting of 10 to 15 images were acquired using a magnification of 400×. Images from BIB-prepared samples (see Sect. 3.2.1) were taken with a Zeiss Supra 55 SEM with backscatter and secondary electron detectors (GED, RWTH Aachen). Both overview image mosaics with a magnification of 125× and 1250× were acquired and in addition high magnification images (up to 20 000×) were collected. For testing the representative element area (REA) relative to the pore distribution the box counting method was used (Kameda et al., 2006;Klaver et al., 2012), where a stepwise growing box was applied to the segmented pore images in which the total area% of pores within the box was measured. At a certain box area the total area% of pores does not define a significant variation anymore. This box area was taken as the REA. We checked this for two different samples using BIB-SEM images (Fig. 6). Sample A shows a minimum REA of about 900 µm 2 , whereas sample C requires a minimum REA of at least 25 000 µm 2 . The microstructure reveals the differences between those REAs (Fig. S2). In samples where large heterogeneities in the microstructure (e.g. the local presence of fossils) exist, the calculated REA is larger than in samples with a relatively more homogeneous microstructure. In the following, image analysis is presented from areas between 1750 and 31 500 µm 2 .

Image analysis
Images were first cleaned by manually removing artefacts from dust or charging effects. Images were segmented using manual thresholding in ImageJ version 1.51n (Schindelin et al., 2012;Schneider et al., 2012). This software was also used for the pore analysis. The outlines of the analysed particles were projected as overlay on the original images in Photo-shop CS 5.1 to optimize the thresholding values. Pores were classified by their morphology, based on the ratio of the minor axis to the major axis (b / a). Pores with ratios of < 0.2 were considered as fracture. This criterion comes from the b / a vs. area % plot shown in Fig. S1, which shows that for values b / a < 0.3, the porosity steeply increases as the effect of fractures. Manual pore segmentation yields a human error of about 10 % . Data from two to three mosaics taken on several different places within the same sample were averaged, and the standard deviation was taken as mean error. The highest given standard deviation is then applied to all samples. The standard deviation of the porosity obtained by image analysis from mechanically polished thin sections is 0.5 area%, and it is 0.15 area% for BIB-SEM prepared samples.

Microstructure, petrography and sample heterogeneity
The slates are defined by alternating intervals of silty and clay-(mica)-rich layers (millimetre to centimetre scale) owing to turbiditic sedimentation processes. The clay-(mica)rich layers include four dominant mineral phases: mica, chlorite, quartz and calcite. In minor amounts there are also clay, dolomite, albite, pyrite and rutile. Mica occurs in two distinct phases.
(1) There are relatively large mica aggregates up to 20-50 µm in size, often deformed. These aggregates occur both aligned and randomly oriented with respect to the foliation plane (Fig. 7b, d). (2) There are relatively small mica grains (1-5 µm) mostly aligned in the foliation plane (Fig. 7d). The alignment of the sheet silicates in the foliation plane forms the main planar fabric. The low-temperature endmember (sample A: T max = 200 • C; Fig. 7a and b) exhibits no preferred alignment of minerals. This is a low strain endmember as inferred from the occurrence of un-deformed fossils in a clay-rich matrix. Sample E (T max = 250 • C; Fig. 7c and d) shows the alignment of mica and an increase in microfractures along the grain boundaries of mica. The hightemperature endmember (Sample H: T max = 320 • C; Fig. 7e and f) is characterized by a strong and defined spaced cleavage with pressure solution seams and shows no fossils. This sample also represents the high strain endmember, as indicated by the elongation of detrital quartz and calcite crystals. An overview of the microstructures of all the samples is given in Fig. S2.
Vertical and horizontal sample heterogeneity is determined by the approach shown in Fig. 3. We analysed two adjacent drill cores (core axis is bedding-perpendicular), which were measured as sub-layers along the foliation planes. He-pycnometry results show that the bedding-perpendicular porosity shows substantial variations between 1.4 and 5.2 vol% (sample plug G1: Fig. 3b). Only small variations < 0.8 vol% are seen parallel to the bedding, meaning that within the sample volume the sub-layers of each plug correlate laterally very well and are therefore relatively homogeneous. Despite the bedding-perpendicular porosity variation, image analysis of the same sample material demonstrates a homogeneous lithological composition (Fig. 3c). The porosity variation across these layers is primarily due to the local occurrence of fractures.

Pore morphology
Porosity consists of two main groups classified by their pore morphology. (1) There is macro-fracture porosity, which includes macro-fractures with apertures up to several millimetres and a length of a few millimetres to centimetres. Such large-scale fractures are immediately visible in thin sections and greatly influence bulk rock measurements. (2) The matrix porosity consists of matrix pores, including inter-and intra-particle pores and microfractures (apertures between 1-5 µm). Matrix porosity is measured by image analysis on areas selected in between macro-fractures.
Image analysis on mechanically polished thin sections has shown that porosity estimates from these thin sections are unreliable. The only reliable information that these analyses yield is an insight into the contribution of microfractures. In contrast, BIB-prepared samples give a good insight into the complete matrix porosity. -(1) Open fractures: in most cases these fractures are bedding-parallel and occur on the foliation planes as seen in the preparation plugs and by UV light detection of the fracture network of a sub-sample (Fig. 3e).

Macro-fracture porosity
In thin sections such fractures are recognized by the relatively large size in both aperture (> 10 µm) and length, as these run through the complete section. Sometimes material is smeared into the fractures during mechanical polishing (Fig. 8). They show no precipitation of secondary minerals after fracturing.

Matrix porosity
Image analysis on BIB-prepared samples allow us to divide the matrix porosity into two types of pores.
(1) Microfractures are defined by an aspect ratio of b / a < 0.2, have apertures of about 1-5 µm and are not continuous (length up to 30 µm). In Fig. 8, a special type of fracture porosity is seen in fine-grained mica-rich layers, which accommodate brittle strain. Inside the fractures pyrite precipitated, demonstrating a geological origin of these fractures. In such samples, the fracture porosity is the largest contribution to the total porosity. These fractures are often located along the grain boundaries of micas and in some cases along quartz (Figs. 8, 9d  and 10a). (2) The second pore type comprises non-fracturerelated pores; mostly these are interlayer pores associated with mica and clay minerals or due to dissolution along grain boundaries of quartz ( Fig. 9c and d). Intra-particle pores occur in minor quantities in calcite, dolomite (Figs. 9c, d and 10b) or quartz ( Fig. 10e and f). Such pores are often related to quartz and calcite cement, which occurs in some cases as new infill of fossils (Fig. 10c) or fluid inclusions.

Quantification of porosity
All porosity estimates obtained by the different techniques applied in this study are given in Table 1 and Fig. 11. Note that the techniques each have their own characteristic resolution (Fig. 2a). Moreover, He pycnometry and MIP both obtain a connected porosity. Image analysis yields microfractures and inter-and intra-particle pores. Microfractures and inter-particle pores have the greatest contribution to the total porosity in our sample set as seen by microstructural imaging (Figs. 8 and 9). The porosity of the bulk rock of the slates, obtained by He pycnometry, shows that the samples have a connected porosity in the range between 0.7 and 7.2 vol% (Fig. 11). Moreover, He-pycnometry results show a bulk density of 2.55-2.73 g cm −3 and grain density of 2.73-2.78 g cm −3 .
The MIP results show variations in porosity from 0.4 to 2.7 vol%, which does not directly correlate with the Hepycnometry data in all samples. Pore volumes with a pore throat radius and fracture aperture between 5 and 20 nm contribute most to the total porosity ( Fig. 12a and b). In addition, a slightly increased frequency of pores with radius and fracture aperture between 100 and 200 nm is seen. In Fig. 13 the neck pore entrances are plotted vs. the cumulative pore volume. Samples C, G and H have the highest bulk porosity (Fig. 13). The second cycle contains much less remaining mercury in the pores, as the intrusion and extrusion curve approach each other (Fig. S3). The correction which is applied comes from the assumption that for extrusion another contact angle applies (110 • ; Kaufmann, 2010). Hence, the correction is only for extrusion data. The hysteresis between the second intrusion and extrusion almost disappears by applying the above explained pore size correction. This is more evident for sample C than for sample H, which indicates a slight increase in hysteresis.
Image analysis from thin sections shows fracture porosities between 0.04 and 1.8 area%. Image analysis from BIB-SEM reveals total (fracture plus matrix) porosity estimates between 0.3 and 1.4 area%, much lower compared to the    above presented bulk approaches. Pore radii and fracture aperture distributions from image analysis (BIB-SEM) show a similar trend as the MIP data (Fig. 12) for pore sizes > 30 nm. Where there is a peak in the 100-200 nm pore radius in the BIB-SEM data, there is also an elevation in the MIP data, although not with the same magnitude. In addition, the MIP data indicate a major pore radius of 10 nm, which cannot be measured by image analysis because of resolution limits of the SEM. Figure 11 shows an exponential decrease with asymptotic behaviour for decreasing sample input from bulk rock measurements to image analysis. There exists no direct trend along the temperature and indirect strain gradient from sample A to H for any of the measurement techniques (see varying colour-coded symbols in Fig. 11 for the different approaches). Nonetheless, knowing the microstructure enables us to interpret the outcome of the porosity results. Sample C, with a relatively high amount of large mica grains, yields the highest porosity values (for all techniques). The highestgrade samples G and H, also show relative high porosity for MIP, due to bedding-parallel microfractures, whereas sample A, which is the weakly deformed homogeneous clay sample, plots at the lower end of the diagram, indicating lower porosity (for MIP and BIB-SEM). The sub-surface samples (samples D and E), which are from the close vicinity of sample C, in general show a lower porosity in He pycnometry and MIP than for the mechanically polished thin section of sample C.

Discussion
In this work, we use a workflow which includes validating bulk rock measurements with high-resolution imaging using BIB-SEM to differentiate different types of porosity of slates. In a second step, we link these estimates for the present-day porosity to processes along the pressure-temperature-time (PTt) deformation path.

Contribution to porosity estimates
Artificially induced porosity is often due to sample collection and preparation. During sample collection, artefacts like fracturing could arise from hammering, transport or storage (drying). Sample preparation, like cutting, and mechanical polishing of thin sections could induce significant additional artificial porosity, such as breakouts, resulting in apparent porosity values without geological meaning (Fig. 5). Therefore, for image analysis procedures on slate samples, sample preparation like ion polishing techniques, either on 2-D thin sections or on cross sections of sub-samples, are necessary (e.g. Hemes et al., 2013;Keller et al., 2013a;Klaver et al., 2012Klaver et al., , 2015. The bulk rock measurements are strongly influenced by the presence of natural and possible artificial fractures and sample heterogeneity. This results in non-correlating Hepycnometry and MIP data. Differences in porosity estimates from bulk rock measurements and image analysis techniques are related to the different principles underlying the techniques, which are the sample volumes (Ougier-Simonin et al., 2016) and the limits in spatial resolution (e.g. Hemes et al., 2013). In addition, Kaufhold et al. (2016) acknowledge the resolution issue by a comparison of MIP with FIB-SEM data from shales and show that only 20 vol% of the total porosity can be resolved with FIB-SEM; the remaining porosity content should be measured with a bulk rock measurement technique such as He pycnometry. Additionally, large microfractures are not incorporated into the representative domain of the BIB-SEM but do play an important role in the MIP data (Hemes et al., 2013). This correlates with our observations: larger cracks are measured in the bulk rock methods, whereas the BIB-SEM data only provide matrix porosity as we avoid imaging the large fractures. Hence, using a combination of such methods, the REA of the image analysis does not necessarily need to cover the full range of pore sizes. However, the BIB-SEM investigation links the bulk He pycnometry and MIP to each other. MIP and He pycnometry both measure connected porosity. MIP includes pores with a > 3.7 nm eq pore throat radius and He pycnometry > 1 nm eq pore throat radius. Both methods do not provide information on the pore type while the BIB-SEM does.
Porosity types from different scales seem to be related to each other (Fig. 11). A sample with many small-scaled microfractures has also the highest amount of macro-fractures, as seen by their relatively high contribution to the bulk rock porosity. There is a correlation of pore radii and fracture apertures between the BIB-SEM and MIP datasets in the overlapping measured scale (Fig. 12), even though BIB-SEM yields sections of pore bodies in 2-D images and MIP pore volume (3-D). For the overlapping size range of these two methods, these methods are comparable to each other (see also , indicating that 2-D data are sufficient to estimate average values. Each measurement technique has its characteristic pore size range it can analyse, which means that combining different methods results in overlapping scales, allowing us to resolve the porosity from nanometre to centimetre scale. Our analysis shows that sample inhomogeneity between the internal slate layers can be large, comprising primary sedimentological variations, such as clay-(mica)-rich layers alternating with silt-rich layers. For bulk rock measurements the representative element volume (REV) is dependent on the clay content and grain size distribution of sand grains for the different sediment layers (e.g. carbonates and quartz): the higher the clay content relative to the non-clay content, the smaller the REV (Keller, 2015).
In light of a comparison of our slate results with previous studies performed on clay stones as precursor rocks of slates, MIP-measured porosity ranges reported for the poorly compacted Boom clay (local name of Oligocene clay in Belgium) are about 35-40 % (Boisson, 2005). Boom clay has a high pore connectivity, and it is shown that in samples with a higher clay content the pore connectivity is controlled by pores in the clay matrix, whereas in samples with a smaller clay content, large inter-particle pores have the largest contribution to the total connectivity (e.g. Hemes et al., 2013). A second clay stone type, the Opalinus clay is more compacted and cemented than the Boom clay (about 10 % porosity from MIP: Houben and Urai, 2013 and references therein). It shows smaller values for connected porosity and like in the Boom clay, the clay matrix controls the permeability . In contrast to these clay stones, the slates from our study show 0.4-2.7 vol% porosity from MIP measurements, and the pores in the clay matrix cannot be resolved with a SEM. This means the sheet-silicate ma- trix is much tighter in these slates, which experienced higher PT conditions and in contrast to the Boom clay and Opalinus clay, show pressure solutions seams, elongated grains, and a much smaller number of fossils due to recrystallization and diagenesis. In these higher-grade samples, the clay minerals are replaced by micas; and due to their strong anisotropy, (micro)fractures, as also illustrated by the BIB-SEM images of sample H, are the main contributor to pore connectivity. As a result, the porosity and pore morphology from slates in this study are not directly comparable to the pore types reported for diagenetic clays (e.g. Boom clay; Opalinus clay) but show an advanced stage of porosity reduction owing to an increased metamorphic overprint.

Processes influencing porosity
To understand what part of the total connected porosity is fracture porosity and what part is matrix porosity, image analysis can be used to establish a threshold. Figure 11 shows that image analysis gives an estimate of the matrix porosity. Subtracting the BIB-SEM matrix porosity from the MIP bulk rock measurements yields the macro-fracture porosity. Subtracting the MIP porosity from the He pycnometry gives insight into the smallest matrix pores.
Applying theoretical considerations, Dielforder et al. (2016) estimated that porosity decreases from 60 % to 10 % for convergent plate boundaries within the upper 5 km of the accretionary wedge along the prograde metamorphic gradient (see references in Dielforder et al., 2016). Along this prograde evolution, the tectonic processes change from soft sediment deformation (such as particulate flow) to pressure solution (Dielforder et al., 2015(Dielforder et al., , 2016. The porosity evolution found in this study is shown in Fig. 14 and indicates different vein generations on both the prograde and retrograde evolution. Veins are indicative of brittle failure and can therefore be used as indicators of ancient fracture porosity at the stage prior to fracture sealing by precipitation of the vein minerals. The earliest formed fracture porosity is shown by deformed shear veins overprinted by the main foliation (stage I). From initial sedimentation to maximum burial the rock texture changes from an un-deformed matrix with randomly oriented sheet silicates to a slate with a well-developed spaced cleavage and strong structural anisotropy (Fig. 7). The porosity decreases from 60 vol% (Dielforder et al., 2016) to values < 0.9 vol% at maximum burial (peak metamorphism). This is reflected by the matrix porosity obtained by image analysis on BIB-SEM prepared sub-samples from the underground mine of Landesplattenberg. At these highest temperatures reached, dynamic recrystallization of bedding-parallel veins takes place (Fig. 14: stage II). On the prograde path, the sediments are strongly influenced by compaction together with dehydration and diagenesis as shown by diagenetic quartz and/or calcite cement precipitation. Diagenetic cementation processes play an important role as a porosity-reducing mechanism in these slates accompanying dissolution by pressure solution (e.g. Katsube and Williamson, 1994). While pressure solution generally reduces porosity, pores are formed by dissolution in and around dolomite (Fig. 10b). This newly created porosity, however, is far less than the amount of porosity destroyed by pressure solution.
Along the temperature and correlated strain gradient from Weesen (200 • C) to Flims (320 • C), there seems to be no relation between porosity and temperature/strain either in the bulk rock measurements or in the images analysis results. However, local variations in mineralogy and strain influence the porosity: a high sheet-silicate content leads to a strong anisotropic fabric, resulting in an increase in microfractures and therefore in a higher porosity (sample C). An only weakly deformed sample, with a less defined fabric, exhibits a low porosity (sample A). This suggests that porosity and pore morphology are closely related to the mineralogical composition and spatial heterogeneity of the phase distributions and strain (e.g. Janssen et al., 2011;Keller et al., 2013a). Together with aforementioned temporal porosity changes by fracturing and vein mineral precipitation, the prograde local porosity evolution can therefore be fairly variable and dynamic in time but results in a general decrease in the average porosity with increasing burial.
The porosity on the retrograde path is mainly characterized by fracturing and the formation of bedding-parallel, sometimes partly recrystallized veins, veins oblique to the bedding and abundant microfractures ( Fig. 14; stages III, IV and V). In a first stage, these microfractures evolve in an isolated but spatially dispersed manner along clay-rich layers. Sometimes pyrite and/or calcite precipitated in such microfractures clearly evidences a synmetamorphic process (Fig. 8). The respective porosity values are represented by the fracture porosity obtained from bulk rock measurements of subsurface samples (Fig. 14). Such microfractures are likely the effect of unloading.  Dielforder et al. (2016) link the pressure solution processes and the illite-muscovite transformation to an embrittlement during the prograde path. Such embrittlement allows the formation of microfractures (e.g. Zeng et al., 2013). Along the retrograde path, this effect is enhanced because temperatures decrease again, allowing the further increase in size and amount of fractures. With progressive exhumation, the microfractures grow and eventually interconnect, transferring microfractures to fractures (Fig. 14).
In the light of near-surface porosity formation, plotting fracture porosity from bulk rock measurements of surfacecollected samples and comparing them with the underground sample reveals the contribution to porosity formation by surface processes. In the surface-collected samples, which are exposed to strongly varying climatic conditions, such fractures could be due to freeze and thaw cycles (Cárdenes et al., 2012) but also to unloading during deglaciation resulting in exfoliation jointing (Jahns, 1943;Ziegler et al., 2016). Particularly in the case of the mechanically highly anisotropic slates, the preexisting foliations are prone to the reactivation of these structural planes by fracturing during the youngest climate-and weathering-related processes. Additionally, macro-fractures can in many cases also be artificially induced, requiring discrimination between artificial and young geological fracturing processes. Principally, we can list two discrimination criteria. Near-surface fractures experience infiltration of meteoric water resulting in oxidation of Fe-bearing mineral phases such as iron sulfides (e.g. pyrite). Consequently thin layers of iron hydroxides form, appearing as faint reddish staining along the fracture planes. This staining provides clear evidence of an in situ crack formation in the rock within its natural environment. In contrast, preparation-induced fractures will not show such staining phenomena. Moreover, these fractures mostly contain small amounts of polishing material, which is smeared into the newly created fracture during this mechanical treatment (Fig. 8). This accumulation is not possible in the case of resin, used during thin-section preparation, has already filled in the natural fractures. Despite these two unequivocal criteria, it might be difficult to discriminate between the two late fracture types in cases where none of the two types of evidence is present.

Conclusions
In this study, we use a combination of He pycnometry and MIP as bulk rock measurements and image analysis to obtain estimates of the porosity of slates. Bulk rock measurements obtain the total connected porosity, which can be subdivided into fracture and matrix porosity determined by image analysis. Image analysis from BIB-SEM images yields information about the matrix porosity including intra-and inter-particle pores and microfractures. This matrix porosity reflects the porosity at maximum burial (peak metamorphism) and is strongly influenced by diagenetic processes on the prograde path. Subtracting the matrix porosity from the total porosity obtained by bulk rock measurements yields the fracture porosity and matrix porosity below BIB-SEM resolution. Most of the microfractures are bedding-parallel or occur along the grain boundaries of mica and are ascribed to the effect of unloading on the retrograde path. Macro-fractures (fracture apertures up to 1 mm) are in many cases related to surface processes such as freeze and thaw cycles or are artificially induced by sample preparation. Different vein generations show that the formation of porosity is not restricted to unloading or surface processes but fracturing took place along the entire PTt-deformation path. From a temporal point of view, all porosity, excluding porosity formed by surface or artificial processes must be fairly cyclic, given the stages of pore opening owing to fracturing and subsequent sealing by mineral precipitation. Moreover, porosity and pore morphology are strongly dependent on mineralogy, homogeneity and strain. The multiscale approach that was developed on shales, which links microporosity to macroporosity by combining bulk rock measurements with image analysis, is, in this study, successfully applied to slates.
Data availability. All data underlying this manuscript are available in the Supplement.
Author contributions. MH and AB designed the project. Together with IVA, MH and AB collected the samples and IVA carried out the He pycnometry measurements, ion 2-D polishing and SEM work on thin sections. JK carried out the MIP measurements. JLU, JK and GD provided and carried out, together with IVA, the BIB-SEM work. IVA analysed all data and prepared the manuscript with contributions from all co-authors.
Competing interests. The authors declare that they have no conflict of interest.