Multiscale characterisation of chimneys/pipes: Fluid escape structures within sedimentary basins

Abstract Evaluation of seismic reflection data has identified the presence of fluid escape structures cross-cutting overburden stratigraphy within sedimentary basins globally. Seismically-imaged chimneys/pipes are considered to be possible pathways for fluid flow, which may hydraulically connect deeper strata to the seabed. The properties of fluid migration pathways through the overburden must be constrained to enable secure, long-term subsurface carbon dioxide (CO2) storage. We have investigated a site of natural active fluid escape in the North Sea, the Scanner pockmark complex, to determine the physical characteristics of focused fluid conduits, and how they control fluid flow. Here we show that a multi-scale, multi-disciplinary experimental approach is required for complete characterisation of fluid escape structures. Geophysical techniques are necessary to resolve fracture geometry and subsurface structure (e.g., multi-frequency seismics) and physical parameters of sediments (e.g., controlled source electromagnetics) across a wide range of length scales (m to km). At smaller (mm to cm) scales, sediment cores were sampled directly and their physical and chemical properties assessed using laboratory-based methods. Numerical modelling approaches bridge the resolution gap, though their validity is dependent on calibration and constraint from field and laboratory experimental data. Further, time-lapse seismic and acoustic methods capable of resolving temporal changes are key for determining fluid flux. Future optimisation of experiment resource use may be facilitated by the installation of permanent seabed infrastructure, and replacement of manual data processing with automated workflows. This study can be used to inform measurement, monitoring and verification workflows that will assist policymaking, regulation, and best practice for CO2 subsurface storage operations.


Overview
Carbon dioxide (CO 2 ) capture and subsurface storage (CCS) within sedimentary basins has been identified as an effective solution for reducing anthropogenic CO 2 emissions in the atmosphere (IPCC, 2005; The Global CCS Institute, 2019). CCS must form a key component of present and future global climate policy, in order to meet anthropogenic greenhouse gas emission reductions of 80-95 % by 2050, and limit model predictions of likely warming of <2 • C relative to pre-industrial levels (IPCC, 2014). Requirements for wide-scale implementation of CCS include: 1) cost-effective CCS technologies; 2) government policy/incentives for negative emissions technologies; and 3) the need for public acceptance/confidence, all of which are intrinsically linked. The primary technological requirement to widespread implementation of CO 2 storage is subsurface site characterisation and containment assurance. Legal regulations governing CCS in Europe exist in the form of the EU CCS Directive on Geological Storage of Carbon Dioxide, 2009/31/EC (2009, which defines requirements for CO 2 storage across the lifetime of a storage site, including closure and post-closure obligations. Factors that must be considered for the characterisation and assessment of potential CO 2 storage complexes and their surroundings include the role and impacts of potential fluid migration pathways causing loss of containment, and the potential flux rates through these pathways. CO 2 can be sequestered into porous and permeable subsurface sandstone reservoirs (Bachu, 2000;Benson and Cole, 2008), such as depleted oil and gas reservoirs and saline aquifers. Sandstone reservoirs of this type are commonly overlain by impermeable cap rocks and overburden stratigraphy, which together provide an effective seal that prevents the upward migration of CO 2 , ensuring safe and permanent storage. Offshore CO 2 storage in sandstone reservoirs has been successfully demonstrated in Europe and globally (e.g., Sleipner, North Sea; Tomakomi, Japan). Several other commercial-scale offshore CO 2 storage projects are also in planning or development stages, such as the Northern Lights project (e.g., Furre et al., 2019a;Global CCS Institute, 2019).
The location and potential intensity of any possible loss of CO 2 containment from a storage reservoir is dependent on the distribution of fluid pathways in the cap rock and overburden, and the ability of these pathways to transmit dissolved, liquid, and/or gaseous CO 2, depending on the pressure-temperature conditions and presence of other fluids. Potential pathways may include anthropogenic sources, such as abandoned wells (e.g., Watson and Bachu, 2009), formation level inherent structures, including natural migration up-dip along permeable stratigraphic horizons (e.g., Tóth, 1980;Hindle, 1997), and the formation or reactivation of fluid escape structures (e.g., Nichols et al., 1994;Frey et al., 2009). Such fluid escape, or seal bypass, structures permit pressure-driven, focused fluid flow, which hydraulically connects deeper strata with the seafloor through inter-connected faults, fractures, and porous-permeable sediment layers . Fluid escape can take place as single blow-out events, episodic/pulsed flow, or as continuous seepage flow. The type of flow can vary depending on the subsurface pressure, stress, and lithological conditions. Therefore, the activity of a fluid escape structure may exhibit temporal variability, which may be cyclical over both short timescales such as tidal cycles (e. g., Boles et al., 2001) or longer-term sea level changes (e.g., Plaza-Faverola et al., 2011;Riboulot et al., 2014). Therefore, the combined understanding of the presence of fluid pathways and their fluid flow regime is critical for the risk assessment of potential subsurface CO 2 escape.

Aims & objectives
In this contribution we provide a broad overview of the integrated geophysical, geological, and geochemical methods which can be applied to the characterisation of focused fluid conduits. In order to achieve this, we use as our context an exemplar study of the fluid conduit beneath the Scanner pockmark, which we undertook as part of the European Union Horizon 2020 project Strategies for the Environmental Monitoring of Marine Carbon Capture and Storage (STEMM-CCS; http://www.ste mm-ccs.eu), together with a partner project CHIMNEY (Characterisation of major overburden leakage pathways above sub-seafloor CO 2 storage reservoirs in the North Sea; Bull et al., 2018). STEMM-CCS and CHIMNEY focussed on determining the permeability of subsurface fluid pathways, and developing better techniques to locate fluid escape structures, so that they can be better quantified and constrained, with relevance to potential fluid flow at CO 2 storage complexes.
The paper has three main aims: (1) Firstly, we describe the various methods which may be applied to the characterisation of focused fluid conduits, and which allow us to resolve the physical parameters of interest (Section 3). For each individual method, we describe their capabilities, consider their benefits, address uncertainties, and deduce any further developments that may be required. In order to contexualise the various techniques which may be used, we describe their application during our investigation of the Scanner pockmark.
(2) Secondly, we discuss the scales of imaging and resolution of the various methods described here, and the co-dependencies which exist between them. This permits integration into a multi-scale approach, and demonstrates how the different techniques may be employed in combination, to ensure appropriate constraint and calibration of the different methods, for complete characterisation of the fluid escape structures (Section 4). (3) Thirdly, based on our findings from this study, we describe a framework which can be used to determine the approaches that are needed to understand potential fluid flow structures in marine environments, in the context of the risk assessment of potential future CO 2 geological storage sites (Section 5).

Focused fluid conduits in seismic data
Seismic chimneys (e.g., Hustoft et al., 2010) or pipes (e.g., Moss and Cartwright, 2010a), referred to hereafter only as chimneys, are observed in seismic reflection data as vertical to sub-vertical anomalies with circular or elliptical planforms, displaying seismic blanking and discontinuous or chaotic reflections (e.g., Løseth et al., 2011). Where free gas is present in the chimney, high amplitude seismic reflections, known as bright spots, with polarity reversals may be observed at discrete intervals, indicating gas accumulation during migration in layers of porous sediments (e.g., Ostrander, 1984). Pull-up of reflectors may also be observed, caused by high seismic velocities, which are commonly attributed to authigenic carbonate accumulations, or where located in the gas hydrate stability zone, to the presence of gas hydrate (e.g., Plaza-Faverola et al., 2010). If CO 2 migrates from a sub-seafloor storage reservoir and reaches the base of these chimneys, and if their permeability is sufficiently high, they could act as CO 2 pathways towards the seafloor and overlying water column. To provide a reliable prediction of potential seafloor seep sites, the degree to which these pathways are able to transmit fluids (i.e. permeability) needs to be better understood.
Chimney-like features can also be the result of seismic imaging artefacts. Seismic imaging artefacts may arise both as a result of data acquisition and/or processing (e.g., Tucker and Yorston, 1973). In particular, for locations where gas is present in the subsurface, the effect this has on seismic velocity determination can have a significant impact on both time and depth migration. An example locality where a seismic artefact was interpreted as a chimney is the Goldeneye field, a prospective CO 2 storage site in the Central North Sea. Following high resolution 3D seismic processing a feature previously interpreted as a fluid escape conduit was later reinterpreted as a seismic imaging artefact caused by a glacial tunnel valley (Dean et al., 2015;Karstens and Berndt, 2015).

Formation of chimneys
Chimney formation was observed on a small scale during a controlled subsurface CO 2 release experiment known as QICS (Taylor et al., 2015;Cevatoglu et al., 2015). In this experiment, CO 2 was released into sediments at an increasing rate of 20-210 kg/day, at 12 m depth below the seabed in shallow water (5− 30 m) in Ardmucknish Bay, Scotland. Repeated seismic reflection data acquisition prior to, during, and after the gas release showed the temporal development of a chimney, formed by gas propagating upwards by fracture generation and reactivation in fine grained sediment . Conditions for hydraulic fracture generation are favourable in shallow (low effective stress) unconsolidated, fine grained sediments, and may be considered a primary mechanism for chimney initiation (Fauria and Rempel, 2011). The upward propagation of fluids may also be facilitated by capillary driven invasion, most prevalent in conditions of high effective stress (Cathles et al., 2010). Further mechanisms for chimney genesis include: erosive fluidisation, localised subsurface volume loss, and syn-sedimentary formation (Lowe, 1975;Sun et al., 2013;Cartwright and Santamarina, 2015).
Based on the observations from QICS, and other experiments (e.g., Fauria and Rempel, 2011;Räss et al., 2018), a conceptual model for the formation mechanism and structure of chimneys in the shallow overburden has been developed . In this model, the first phase of formation is the hydraulic fracturing of low permeability sediments due to high fluid overpressure (Arntsen et al., 2007;Cartwright et al., 2007;Løseth et al., 2009). Seal breaching occurs due to reduced effective stress and leads to either opening of new fractures and/or reactivation of pre-existing fractures, generating a localised connected fracture system. These pathways may permit vertical buoyancy-driven migration of gas-rich pore fluids through the fracture network. Large-scale chimneys (~100− 1000 m wide) are therefore hypothesised to represent a series of interconnected sub-vertical or radial fractures, which allow the vertical flow of gas in the shallow subsurface  due to the elevated permeability relative to the normal 'background' permeability of the host sediment . A transition from fracture to capillary dominant flow behaviour may be observed with increasing depth, due to increased overburden thickness (higher effective stress).
In addition to the role of the geometric structures of chimneys in governing fluid flow, coupled physical and chemical processes act within fractures and pores, resulting in complex feedback mechanisms between porosity-permeability and CO 2 /CH 4 reactivity, which may affect the hydro-mechanical response of the system. From a quantitative perspective, little is known regarding the impact that chimneys have on the upwards migration of CO 2 and CH 4 to the seabed Marín-Moreno et al., 2019).

Study area -Scanner pockmark
We have investigated an exemplar natural fluid escape system located near the centre of the Witch Ground Basin, located 190 km off the north-east coast of Scotland (Fig. 1a), within licence block 15/25 of the North Sea. Here, chimneys are observed, which underlie active natural methane venting sites at the seabed, within ~150 m water depth (Fig. 1b) -the Scanner, Challenger and Scotia pockmark complexes (Gafeira and Long, 2015). Pockmarks are seabed depressions, created by release of over-pressured pore-water and gas from the subsurface (Hovland et al., 2010). In this area of the North Sea, large pockmarks (>6 m deep, >250 m long, and >75 m wide; class 1, Böttner et al., 2019) are continuously active. Evidence for active methane venting at the Scanner pockmark complex is provided by water column imaging and the presence of methane derived authigenic carbonates (MDACs) at the seabed (Judd et al., 1994;Judd and Hovland, 2009). The West Scanner pockmark (Fig. 1c) releases methane at 1600− 2600 kg/day (Li et al., 2020), derived from a combination of biogenic and thermogenic sources (Clayton and Dando, 1996). Smaller pockmarks (class 2) are also distributed across the area (Fig. 1b) with a dominant NNE/SSW orientation (>1500 across 225 km 2 ) and are interpreted as dewatering features attributed to localised pressure changes . The Scanner pockmark complex overlies an area that has been appraised for CO 2 storage potential, the East Mey Storage Site (ACT Acorn Consortium, 2018).

Stratigraphy and seismostratigraphic framework
The lithostratigraphy and seismostratigraphic framework of the ~600 m-thick Quaternary sediment succession containing the Scanner pockmark complex was described by Böttner et al. (2019) and Stoker et al. (2011) (Fig. 2). Deposited within the Witch Ground Basin (Andrews et al., 1990), this complex is underlain by the Hordland and Nordland Groups, of Palaeogene and Neogene age respectively, and which are composed of claystone with limestone and sandstone interbeds (Judd et al., 1994).
The basal Quaternary unit, the Aberdeen Ground Formation (Fm.; unit S1 in Fig. 2), is composed of clay-rich sediments deposited in the Early Pleistocene (up to Marine Isotope Stage, MIS, 13), and displays a laterally continuous, layered seismic character (Stoker et al., 2011;Ottesen et al., 2014). The overlying Ling Bank Fm. (S2) erodes into the top of the Aberdeen Ground Fm., representing a regional glacial unconformity, with deposition of glacial tills in the Early to Middle Pleistocene (~1.2− 0.5 Ma, MIS 12 to 10; Stewart and Lonergan, 2011;Reinardy et al., 2017;Böttner et al., 2019). The glacial tunnel valleys of the Ling Bank Fm. display both a layered and non-layered seismic character. The Coal Pit Fm. (S3) overlies the Ling Bank and Aberdeen Ground Fms., and comprises glacial tills (which include pebbly and muddy sands) deposited in the upper Mid to Late Pleistocene (MIS 6-3; Andrews et al., 1990;Stoker et al., 2011). The Coal Pit Fm. is conformably overlain by Last Glacial Maximum (LGM) deposits, which comprises silty-sandy clays with rare pebbles, deposited during MIS 3− 2 (S4). Unit S4 extends upwards to the base of Scanner pockmark. The Coal Pit Fm. and LGM deposits have a dim and chaotic seismic character and are conformably overlain by the Witch Ground Fm. (S5), composed of silty clay deposited during MIS 2 to 1 (Stoker et al., 2011). The Witch Ground Fm. has two main units: the lower (S5.1) and upper (S5.2) Witch Ground Member respectively. The lower Witch Ground Member has an apparent interbedded seismic character. The upper Witch Ground Member conformably overlies the lower Witch Ground Member, thinning or pinching out towards the northeast. The upper Witch Ground Member is composed of sediments of Holocene age (MIS 1) and has a uniform seismic character (Stoker et al., 2011). The Scanner pockmark depression erodes down to the base of the lower Witch Ground Member.

Overview & data acquisition
The methods used to investigate chimneys and their associated fluid flow can be divided into several types (Table 1). These include: seismic reflection imaging, ocean-bottom seismic methods, including seismic tomography and anisotropy analysis, controlled source electromagnetic surveying (CSEM), active acoustics, passive seismic monitoring, sediment sampling of both the target site and onshore analogues, laboratory rock physics experiments, and process-based numerical modelling and fluid dynamic modelling approaches.
In support of the STEMM-CCS and CHIMNEY projects, four research expeditions were undertaken for data collection at the Scanner pockmark complex. 2D seismic reflection and refraction data, for use in travel time tomography, were acquired using a GI gun source and 18 ocean bottom seismographs (OBS) by RV Maria S Merian cruise MSM63 (Fig. 3a,d;Berndt et al., 2017), in addition to multi-beam bathymetry data, controlled source electromagnetic (CSEM) data (Fig. 5), and Parasound sub-bottom profiling data (Fig. 2a). RRS James Cook cruise JC152 conducted a wide-ranging seismic experiment over the Scanner and Challenger pockmarks (Figs. 3 & 4). Five different seismic sources were used (Bolt and GI airguns, Squid and Duraspark surface sparkers, and a deep towed sparker; Fig. 3), which were recorded by arrays of 25 and 7 OBSs, at the Scanner and Challenger pockmarks respectively ( Fig. 3a; Bull, 2017). The seismic sources were also used to generate multichannel (GI guns, surface sparkers; Fig. 2b) and single channel (deep tow sparker) seismic reflection profiles. During cruises MSM78  and RV Poseidon POS518-2 (Linke and Haeckel, 2018), sediment cores for geological and geochemical analysis were taken from within the Scanner pockmark and a reference site (Fig. 1b,c) using a gravity corer and rock drill (RD2), which acquired core to depths of ~6 and ~33 m below seafloor (mbsf) respectively. Additionally, data from a 3D seismic survey conducted by PGS ( Fig. 2c; CNS Mega Survey Plus) were used to further support the study.

Introduction
Multichannel seismic reflection data record the wavefield reflected from physical discontinuities in the subsurface at a range of distances from an active source of acoustic energy. By doing so, they carry information about seismic wave propagation velocity and anomalies in visco-elastic properties. So long as the geological boundaries correspond to changes in physical properties controlling seismic wave propagation, reflection data imaging provides us with a representation of subsurface sedimentary and tectonic structures with a resolution equal to a fraction of the propagated wavelength (Kallweit and Wood, 1982). Sub-vertical fluid escape structures and sediment deformation (e.g., sediment mobilisation and polygonal faulting) can therefore be imaged, potentially also allowing for a relative dating of geological events by interpreting cross-cutting structures and stratal stacking patterns.
Changes in pore fluid type, especially partial gas saturation, have a strong influence on seismic velocities and the absorption of seismic wave energy (White, 1975;Domenico, 1977). The bulk effect of gas accumulation within sediment pores is a reduction in the sediment compressibility and, therefore, acoustic impedance (Tóth et al., 2014). Typically this results in strong acoustic impedance contrasts, visible in the seismic data as negative polarity reflections known as bright spots . More generally, pore gas manifests itself as local increments in the subsurface reflectivity, as a function of the properties of the encasing medium and the local partial gas saturation (Berndt, 2005;Cartwright, 2007;Løseth et al., 2009;Andresen, 2012;Karstens and Berndt, 2015). Frequency-dependent attenuation and velocity dispersion can also be observed in ultra-high-frequency (in the order of kHz) data, as the seismic frequency approaches the gas bubble characteristic frequency, allowing the detection of gas migration irrespective of the presence of reflective anomalies (Tóth et al., 2015). Therefore, not only do reflection data allow us to image potential gas migration pathways in the subsurface, but also to detect gas accumulation pockets and potentially quantify their volume.
However, as a consequence of the strong energy partitioning at the interface with gas-charged sediments, acoustic blanking may be observed within chimneys (Fig. 2b,c; Fader, 1997), which hampers the Fig. 1. STEMM-CCS and CHIMNEY study area. a) Regional map of northern North Sea plotted over GEBCO (2019)  effectiveness of seismic reflection imaging and inversion in such contexts. As a result, wide-angle transmission data recorded at the seabed, for example on OBS, may also be required, in order to characterise seismic velocities below the gas layer. Wide-angle transmission also helps resolve the issue of ambiguity between the position of the reflector and the true interval velocity along the wave-path, since the short offset range and the limited bandwidth of seismic reflection data result in a limited sensitivity to interval velocities (Jannane et al., 1989). Therefore, integration with diving wave travel time tomography, in addition to reflection waveform inversion techniques (e.g., Brossier et al., 2015) and seismic-to-well calibration are necessary to compensate for the lack of sensitivity to the kinematic properties of the medium, and provide an accurate depth representation.

2D vs 3D seismic reflection
Two-dimensional (2D) seismic imaging assumes that subsurface properties are invariant with respect to the direction normal to the survey line. In the case of inherently three-dimensional structures, such as pockmarks and chimneys, this assumption represents a potentially significant source of error. Whilst a chimney can be assumed, in a simple case, to be radially symmetrical around the depth axis, there potentially may be out-of-plane reflections that affect imaging. Three-dimensional (3D) seismic data have substantially advanced the knowledge of subsurface fluid migration features, which were previously often discarded as poorly imaged zones and seismic artefacts in 2D seismic data . The analysis of 3D seismic data is an effective method to map fluid accumulations in the subsurface, identify permeability barriers, and constrain subsurface geometries of entire fluid flow systems Løseth et al., 2009;Andresen, 2012;Karstens and Berndt, 2015). Where only 2D seismic data are available, application of 3D processing of 2D acquisition geometries (e.g., Lin et al., 2019) has proved extremely useful for improving the quality of the final results, provided that variations of the streamer position during acquisition are monitored within the desired accuracy (Whiteside et al., 2013) or estimated from the data (Clay and Vardy, 2018).
However, the resolution of conventional seismic data, laterally (often above 12.5 m) and vertically (~10 m for a dominant frequency of 40 Hz), is often not sufficient to image the seismic expression of fluid flow systems in detail. Recent developments in 3D high-resolution seismic techniques allow imaging of the shallow subsurface in much greater detail than previously (e.g., Planke et al., 2009;Petersen et al., 2010;Plaza-Faverola et al., 2015). An example technology is the P-Cable (e.g., Planke et al., 2009), which has been used to map shallow gas accumulations, gas hydrate systems, and fluid flow structures, such as chimneys (e.g., Plaza-Faverola et al., 2011;Bünz et al., 2012) down to a resolution of 3 m. In contrast to conventional seismic frequencies of ~5-120 Hz (suitable for monitoring deep reservoirs), high resolution P-Cable uses frequencies of up to 350 Hz, surpassing the resolution of conventional 3D seismic.
Collectively, the seismic experiments conducted at the Scanner pockmark complex utilise different imaging resolutions and depths, by applying a multiple-frequency 2D surveying approach using a number of seismic sources ( Fig. 3b; Bull et al., 2018). Progressive extension of the seismic bandwidth to higher frequencies and shorter shot intervals results in higher vertical and lateral resolution. Surface sparker, deep towed sparker, and sub-bottom profiler data (the latter is single-channel, and may be considered a hydro-acoustic technique, but we include it here with the other seismic approaches as it is used for the same purpose) with frequencies up to several kHz reveal near surface structural features previously not discernable using lower frequency sources, such as sediment slumping/flank collapse within the pockmark, more detailed characterisation of gas accumulation at the top of the Aberdeen Ground Fm., and direct observation of fluid migration towards the base of the pockmark (Fig. 3a,b). Therefore, multi-frequency surveys with different sources (e.g., Bolt Gun, GI Gun, Chirp, Boomer, and Sparker) and high-resolution 2D or P-cable streamers can complement conventional 3D seismic data (e.g., Fig. 3c) and lead to a better understanding of the nature and internal architecture of chimneys, particularly in the shallow overburden stratigraphy.

Advanced seismic wavefield analysis
In addition to providing reflectivity images of the subsurface, multichannel seismic data attributes, such as amplitude as a function of offset, can also be exploited to infer kinematic and dynamic properties of the subsurface (Ostrander, 1984). Pre-stack waveform inversion techniques (Virieux and Operto, 2009) can be applied to obtain quantitative characterisation of the effect of partial gas saturation on the elastic properties (P-wave velocities and Poisson's ratio). In particular, partial gas saturation has a strong effect on bulk modulus, associated with a relatively low influence on shear modulus, which corresponds to a lowering in P-wave impedance and a reduction of V p /V s ratio (Ostrander, 1984;Tóth et al., 2014;Provenzano et al., 2018). These contrasts can in turn be exploited to quantify gas accumulation using appropriately calibrated rock-physics models (Tóth et al., 2014;Provenzano et al., 2018).

Table 1
Summary of the methods used for the characterisation of a chimney, displaying the assessed parameters and the method specific traits. The parameters used for complete characterisation of a chimney include: 1) fracture geometry, orientation, connectivity, and subsurface structure; 2) physical propertiese.g., porosity, permeability, and resistivity; and 3) fluid presence, distribution, composition, and flux/flow rate. A method can (green) or cannot (red) be used to assess a given parameter. The traits of each method include: a) temporal applicability -whether the method can be applied at different times to observe temporal variability at the site; b) co-dependencies -whether input from another method is required for calibration and/or constraint; c) time and resource intensiveness -whether data can be produced with a given amount of resources on the timescale of weeks (low), months (medium) or greater than six months (high); and d) cost -a relative economic cost scale for completion (that can be qualitatively described as a method requiring desk time (low), laboratory time (medium) and/or ship time (high)).

Fig. 3. Layout of the Scanner pockmark seismic experiments. a)
OBS locations for JC152 seismic tomography, anisotropy, and passive seismic experiments (white triangles), and MSM63 seismic tomography experiment (orange triangles). b) Far-field source frequency spectra of the seismic sources used in the JC152 tomography and anisotropy experiments. c) Seismic acquisition tracks for JC152, with sources labelled, showing the multiazimuthal geometry of source coverage. Only the airgun (Bolt and GI) and surface sparker (Duraspark and Squid) sources are shown. Blue box shows extent of area shown in a). Red line indicates location of the profiles shown in Fig. 4. d) Acquisition track for the MSM63 acquisition.

Time-lapse/4D seismic reflection
Time-lapse seismic imaging refers to the acquisition of seismic data at the same location at different points in time, in order to assess temporal changes in the subsurface. Where both seismic datasets are 3D, this may be referred to as 4D seismic imaging. This technique has been applied extensively to subsurface reservoir monitoring, including conventional and unconventional hydrocarbon production (e.g., Watts et al., 1996;Landrø et al., 1999;Barkved et al., 2003). Changes in the subsurface due to fluid flow may change properties such as fluid saturation, temperature, porosity, and pressure, and, hence, the elastic and stress-strain properties and behaviour of the rock, which cause a change in seismic response (Johnston, 2013). Such temporal changes can be used to better understand the formation and development of fluid flow features, as well as providing constraints for multiphase thermo-hydro-mechanical simulations. Recognising temporal changes in subsurface fluid flow systems is integral for monitoring CO 2 storage reservoirs and potential leakage in the overburden.
A principal challenge in time-lapse/4D seismic surveying is accounting for the repeatability, such that the produced images represent true temporal changes and not seismic artefacts associated with acquisition and/or processing. Attributes such as the normalised root-mean square noise are used to measure quantitatively the quality of the survey repeatability. Generally, repeatability is excellent where sediments are well stratified and undisturbed (e.g., Waage et al., 2019). However, since chimneys are often much more chaotic seismic structures, repeatability can be poor, and detection of changes in fluid flow between individual time-lapse surveys requires careful interpretation (Waage et al., 2019).
Time-lapse and/or 4D seismic data were not acquired at the Scanner pockmark complex. This technique has, however, been applied to the studies of actively-seeping chimneys, such as at Lomvi pockmark on the Vestnesa Ridge, offshore W Svalbard (Bünz et al., 2012;Smith et al., 2014;Panieri et al., 2017), and for experimental fluid injections, including the Sleipner storage site (e.g., Arts et al., 2004;Chadwick et al., 2004Chadwick et al., , 2019Boait et al., 2012;Eiken, 2019), the QICS experiment , and the STEMM-CCS Goldeneye release experiment (Flohr et al., 2021;Roche et al., in review). Time-lapse seismic data therefore play a key role in the understanding of both naturally occurring fluid flow, and CCS monitoring in pre-, active and post-CO 2 injection phases (Lumley, 2010).

Seismic tomography and full waveform inversion
The P-and S-wave velocity structure of a chimney and its surroundings can be determined by applying travel time tomography and full waveform seismic inversion. These techniques are applied with the principal aims of: locating gas-bearing zones; delineating the shape of the chimney; determining the presence and features of fractures (open or cemented, size and connectivity); and characterising the sediment properties within and outside the chimney.
Seismic travel time tomography is an inversion technique in which observed travel times (e.g., Fig. 4b) are compared to those computed through a discretised and parameterised model representation of the subsurface, with defined parameters that control the balance between minimising the data misfit and generating a model with the minimum required structure to fit the data. Travel time tomography approaches may be isotropic (e.g., Zelt and Barton, 1998) or anisotropic (e.g., Dunn et al., 2005), where in the latter the direction and magnitude of velocity anisotropy is fit and calculated directly, whereas in the former the anisotropy manifests as travel time residuals which can be analysed to determine the model anisotropic properties (Dunn and Toomey, 2001). Travel time tomography is a relatively low-resolution (normally hundreds of m to km laterally, hundreds of m vertically) imaging technique, with the resolution ultimately controlled by the shot and receiver spacings and further being dictated by the discretisation and parameterisation of the model.
The final OBS tomography model can be used as the initial model for full waveform inversion of OBS data (FWI; e.g., Virieux and Operto, 2009). FWI uses the entire recorded seismograms, and, therefore, includes amplitude and phase information instead of solely using the travel times. As a result, FWI represents a higher resolution technique, which is now widely used in crustal-tectonic (e.g., Morgan et al., 2013;Górszczyk et al., 2017;Davy et al., 2018) and engineering (e.g., Smithyman et al., 2009) contexts, in addition to being widely adopted across the hydrocarbon sector (e.g., Sirgue et al., 2010;Prieux et al., 2013). The effectiveness of FWI is highly sensitive to the starting model being capable of predicting the travel time within half a dominant period, in order to avoid cycle-skipping (Virieux and Operto, 2009). This approach is, therefore, commonly preceded by travel time tomography.
Two phases of data acquisition ( Fig. 3a) were conducted at the Scanner pockmark complex to acquire data suitable for the application of first arrival seismic travel time tomography and full waveform inversion. The first survey, conducted during cruise MSM63, utilised an array of 18 OBSs, with shots generated using two GI airguns. The OBS and profile locations are shown in Fig. 3a & d respectively. A second, more extensive experiment was conducted during cruise JC152, coincident with the anisotropy experiment (described in section 2.4, below), and utilising an array of 18 OBSs located in and around Scanner pockmark and 7 at a nearby reference site, located ~1 km southeast of the pockmark, where industry 3D seismic data showed no evidence of the presence of gas (Fig. 3a). These instruments recorded shots from all of the seismic sources used in the experiment (Fig. 3c). Fig. 4c-f shows the shots recorded on an OBS located within Scanner pockmark for an ~E-W oriented profile for each of the airgun and surface sparker sources.
While travel time tomography approaches are not necessarily reliant on a priori information, it is important to consider additional constraints, where available, since a poorly constrained inversion may produce a biased result, based on its parameterisation and/or starting model (Zelt et al., 2003). Interpreted multi-channel seismic profiles may be used to inform the initial P-wave velocity model. In turn, the results of travel time tomographic modelling can be employed to improve the processing and interpretation of multi-channel seismic profiles, by providing independent velocity information. The wider-angle transmission regime involved in OBS methods also assists in resolving the issue of acoustic blanking beneath gas accumulations, such as those observed beneath Scanner pockmark (Fig. 2c), and provides a further tool for the verification of whether observed chimneys represent actual subsurface geological features, or reflection imaging artefacts.
The two seismic approaches, multi-channel seismic reflection imaging and tomography, are highly complementary to one another and may be utilised most effectively when applied synchronously. The results of seismic tomography are also complementary to laboratory scale rock physics experiments, since they place constraints on larger scale in situ . Data in c-e) are shown for profiles and OBEM instruments coloured in yellow (c-d) and red (e), respectively. b) Sketch of survey setup, with towed DASI source, towed Vulcan receivers and three component OBEM seafloor receivers. c-d) Example data at 1 Hz from N to S trending profile (yellow line) for the closest towed receiver c) and the OBEM d) (yellow circle on a). The x axis shows the source position (black antenna with white field lines in b)) along profile, with the location of Scanner pockmark set at 0 m. The OBEM in d) is located at ~700 m along profile. The observed data (black dots), the vertical electric field amplitude (E z ) and the total horizontal electric field amplitude (P max ) respectively, are compared to 2D (c), towed receiver) and 1D (d), OBEM) forward modelled data (coloured lines) for a subsurface with half-space resistivities of 1 to 2 Ωm. The towed receiver data, sensitive to the top ~50-100 m, agree well with predicted data for a 1.2 Ωm subsurface resistivity, while the OBEM data at larger offset to the source agree better with higher resistivities, indicating an increase in resistivity with depth. In c), the location of the grey box corresponds to observed E z values when the transmitter is above the pockmark. In d), a slight increase of P max at ~-200-400 m (highlighted with grey box), however, may indicate a localised increase of resistivities in the subsurface. e) OBEM instrument P max data example from SW to NE trending profile (red line). At about 100-400 m from the Scanner pockmark position, a localised increase of P max is observed, potentially related to an increase in resistivity in the subsurface in the vicinity of the seismic chimney. physical properties, and, thus, fracture orientations and fracture density (Amalokwu et al., 2017;Jin et al., 2018), which then can be used to complement the geo-mechanical models.

Seismic anisotropy
Fractures in sedimentary settings play a crucial role in defining the physical properties of subsurface reservoirs, as they enhance porosity and permeability, or conversely may contribute to reservoir compartmentalisation. Thus, fracture orientation, size, volumetric density, and connectivity are of interest to the understanding of subsurface fluid flow. Attribute analysis of stacked seismic images, which includes techniques such as coherence analysis (e.g., Bahorich and Farmer, 1995), can detect larger fractures. However, these techniques are unable to image smaller fractures, below the spatial resolution of the seismic image. Therefore, in order to determine fracture properties at sub-seismic scale, seismic anisotropy analysis can be applied, which utilises the directional dependence of transmitted seismic signals. A range of theories have been developed to describe the elastic response of fractured rocks (e.g., Hudson, 1981;Thomsen, 1995). While these theories generally agree for dry rock, they differ considerably where fluids and fluid flow between cracks and pores are present (Liu et al., 2000).
Estimating fracture sizes from narrow-band observations of seismic anisotropy may lead to misinterpretations, as there is an ambiguity whereby a medium containing a small number of large fractures will generate the same response as a medium containing a larger number of smaller fractures (e.g., Maultzsch et al., 2003). Therefore, a frequency-dependent approach to studying anisotropy, which is sensitive to the length-scale of the causative mechanism for the anisotropy, is required. Properties such as fracture scale length and fluid saturation can then be inferred from the frequency-dependence of anisotropic attributes, as predicted by theoretical work (Chapman, 2003;Jakobsen and Chapman, 2009). The impact of partial saturation on anisotropy and attenuation in materials of known fracture density and orientation has been studied using models that can link laboratory and field datasets (Amalokwu, 2016), which are required to make robust determinations of permeability.
The Scanner pockmark anisotropy experiment used several seismic sources with different frequencies, in the range ~10 Hz to 2 kHz, as described above (Fig. 3b). This broadband dataset was specifically designed for the measurement of frequency-dependent anisotropy, to permit enhanced fracture characterisation. In order to achieve maximal azimuthal coverage, which is necessary for determining the directionality of anisotropy, profiles were acquired at multiple orientations through the OBS array (Fig. 3c). Several approaches are available for investigation of frequency-dependent anisotropy associated with chimneys. Details of the two active source and one passive source methods that were used at Scanner pockmark are provided below.

Shear-wave splitting
The measurement of seismic travel time anisotropy using shear-wave splitting (SWS) is an established technique for determining the orientation and density of fracture networks (e.g., Crampin, 1985). Shear-wave splitting occurs when a polarised shear wave enters an anisotropic medium. When this occurs, the shear wave is split into two orthogonal components, oriented perpendicular and parallel to the fracture normal direction. The two split shear waves travel at different speeds through the fractured medium, resulting in differences in their travel time. SWS analysis uses converted P-to-S waves to determine the orientation of a symmetry axis associated with the fracture normal direction. SWS distinguishes between different anisotropic symmetry systems, which produce characteristic patterns that can be observed in transformed radial and transverse components of the horizontal geophone records (Fig. 4b). This approach can be used to differentiate between dominant horizontal transverse isotropy, which may be associated with vertically aligned fractures, and vertical transverse isotropy, which may arise in the presence of concentric fractures.
The orientation of these symmetry planes and the directions of the fast and slow S-wave arrivals can be used to determine the fracture orientations. The measured delay between the fast and slow S-wave arrivals is used to determine the intensity of anisotropy, which is related to the fracture size and/or density (e.g., Crampin, 1985;Mueller, 1992;Li, 1997). Raw estimates of SWS give only depth-averaged estimates of fracture properties. Therefore, to determine the depth variation of both fracture orientation and density, a layer-stripping approach is required (Haacke et al., 2009), which recursively compensates and removes the anisotropy measured in shallow layers. Shear-wave splitting analysis is applied to the study of the chimney beneath Scanner pockmark in order to map the geometry and extent of the fracture network. This approach allows us to distinguish between different hypotheses for the structure of chimneys.

Seismic attenuation
While S-wave velocity anisotropy, observed through SWS analysis, is sensitive to both open and closed fracture networks, attenuation anisotropy is primarily sensitive to open fracture networks and the fluids that may be present within these (e.g., Worthington and Hudson, 2000;Chapman, 2009). Attenuation may occur through mechanisms including scattering and wave-induced fluid flow (Baird et al., 2013). Scattering of seismic waves due to aligned heterogeneities has long been recognised to be frequency dependent (e.g., Shapiro and Hubral, 1995;Werner and Shapiro, 1999).
Field observations of attenuation anisotropy have been reported from walkaround vertical seismic profile data (Varela et al., 2006;Maultzsch et al., 2007;Bouchaala et al., 2019), with the observed effects often being consistent with the predictions of poroelastic models accounting for wave-induced fluid flow (Chapman, 2003). If wave-induced fluid flow is the dominant mechanism for the observed attenuation anisotropy, then the phenomenon may provide a link to fracture induced permeability. Chapman (2009) considered the case of multiple sets of fractures having different fluid connectivity.
Application of seismic attenuation and attenuation anisotropy analysis at Scanner pockmark was directed at investigating both the geometry (cf. SWS) and contents of subsurface fracture networks, where these are present. As with travel time anisotropy, the broad-band frequency nature of the various seismic sources used is required to discriminate between the properties of fractures of different sizes, which may otherwise not be distinguished.

Ambient noise anisotropy
Methods for analysing surface waves from ambient noise utilise lower frequencies, typically <1− 2 Hz, than the active seismic sources. This approach further extends the frequency range of frequencydependent analyses to determine fracture properties. Rayleigh wave velocities can be measured on the vertical and hydrophone components of the OBS. By taking spectrograms of cross-correlations between OBS pairs (processed according to Bensen et al., 2007) located on opposite sides of the chimney, Rayleigh wave phase-velocity dispersion can be measured (Yao et al., 2006) and used to invert for vertical shear-wave velocity structure. Applying Radon transforms to a set of cross-correlations between all possible OBS pairs within an array and summing these over time, a technique known as 2D beamforming (Lacoss et al., 1969), shows the slowness at which Rayleigh waves cross the array at all azimuths. This method facilitates full-azimuth observations of anisotropy (Alvizuri and Tanimoto, 2011).
However, there are limitations on the frequencies that can be observed using passive methods, which arise due to the OBS array geometry. The largest spacings within the array, termed the aperture, used at Scanner pockmark result in a lowest usable frequency of ~0.5 Hz (Rost and Thomas, 2009). For beamforming there are also challenges arising from the array geometry. While the OBS array located at the Scanner pockmark is relatively symmetrical, the array at the reference site is not, and so has a highly asymmetric response in the slowness domain, which causes aliasing effects that can cause artificial anisotropy to be observed. For analysis of azimuthal dependence to be effective, any artificial anisotropy must be removed (Picozzi et al., 2010).
Cross-correlation analysis and 2D beamforming can be used to observe the dispersion and azimuthal dependence of Rayleigh waves and can be inverted for 1D shear wave velocity structure at Scanner pockmark in the 0.5-1.5 Hz range. This technique contributes both to constraining the structure of the chimney and extending the range of frequencies which can be used for considering frequency-dependent anisotropy characteristics.

Controlled source electromagnetic surveying
The marine controlled source electromagnetic (CSEM) method is used for mapping variations in the electrical resistivity of the subsurface (e.g., Cheeseman et al., 1987;Edwards, 2005). The CSEM method used in this study involves towing a horizontal electric dipole close to the seabed that transmits an alternating electromagnetic field. The electromagnetic energy diffuses through the seawater and seafloor and is recorded by electric field receivers that are towed in-line behind the transmitter and others that are stationary on the seafloor (Constable, 2013). The measured field amplitude and phase lag between source and receivers relate to the electrical resistivity of the seafloor. When more resistive Earth materials are encountered (e.g., low porosity sediments), elevated electric fields are observed at the receivers and the phase lag is reduced.
The CSEM method is sensitive to changes in the bulk electric resistivity structure of the Earth, a property that depends on lithology and mineral composition but is particularly sensitive to connected porosity and pore fluid content (Palacky, 1988). Due to this sensitivity to sediment pore fluids, the marine CSEM method has been widely used for hydrocarbon exploration (e.g., Ellingsrud et al., 2002;Constable, 2010). Electrical resistivity is also an indicator for whether pathways for fluids to the seabed exist, as the presence of aligned permeable and conductive fractures may lead to electrical anisotropy (e.g., Naif et al., 2015).
The CSEM survey at Scanner pockmark was conducted during cruise MSM63, using the University of Southampton deep-towed active source instrument (DASI; Sinha et al., 1990), towed at 20− 40 m above the seafloor. The source signal comprised a ~100 A, 1 Hz square wave. Twelve profiles were acquired, oriented in four different azimuths to assess the electrical anisotropy. Two different types of receiver were utilised to record the induced field: an array of 14 ocean bottom electromagnetic (OBEM) field receivers (6 three-component receivers and 8 horizontal-component receivers; Fig. 5a), here measuring the electric field only, and two three-axis Vulcan electric field receivers towed behind the source ( Fig. 5b; Constable et al., 2016). The CSEM data from both receiver types are processed with a Fourier transform over 1 s-long windows to obtain amplitude and phase data in the frequency domain, and stacked over 60 s-long windows to improve the signal-to-noise ratio (Myer et al., 2011). Estimating the electrical resistivity distribution from the data requires the use of modelling procedures. For the towed receivers (Fig. 5c), 2D forward modelling for a homogeneous subsurface with one resistivity value (half-space) is performed to compare the modelled data to the observed vertical electric field amplitude (E z ) data and estimate the subsurface resistivity. For the seafloor receivers (Fig. 5d,e), the magnitude of the major axis of the polarisation ellipse traced by the electric field, P max , is compared to 1D forward modelling results for half-space resistivities. P max is independent of the receiver orientation. Resistivities appear to increase with depth. Inversion algorithms (e.g., Constable et al., 1987) which optimise the fit between the observed and predicted data are required to estimate a heterogeneous resistivity structure (analysis for the towed receivers is presented in Gehrmann et al., in review).
There are trade-offs between the thickness and resistivity of anomalous features (e.g., Edwards, 1997) that can be overcome using prior geological knowledge or constraints from complementary data. For example, seismic reflection data can provide high resolution imaging of geological stratigraphy and structures, while CSEM is able to detect the presence and type of fluids in the pore space, such as hydrocarbons, CO 2 , and variations in salinity. The difference between changes in pore fluid or lithological variations may be distinguished by interpreting CSEM and seismic data in conjunction (e.g., Hoversten et al., 2006). Geological interpretations benefit from combined analysis with, for example, wide-angle seismic data (e.g., Goswami et al., 2015), seismic reflection data (e.g., Weitemeyer et al., 2011;Attias et al., 2016;Berndt et al., 2019), magnetic data (Gehrmann et al., 2019), and well logs (e.g., Harris et al., 2009).
The CSEM method is, overall, a relatively low-resolution method (10− 1000 m), due to the diffusive nature of EM fields. The specific resolution and penetration depth depend on the source-receiver geometry, the number of receivers at different offsets to the source, and the frequency spectrum of the source signal (Edwards, 1988). The dual-receiver approach used here includes towed instruments sensitive to the shallow subsurface and OBEM receivers sensitive to the deeper subsurface. For the purposes of understanding potential subsurface fluid escape structures, the CSEM method provides a tool to map variations in resistivity related to connected porosity and fluid/gas content of seafloor sediments.

Passive seismic methods
Passive seismic monitoring involves the use of seismometers and sound recording devices, such as OBSs, located at the seafloor. These instruments provide continuous recordings, and, therefore, capture events during seismic acquisition interludes. In the Scanner pockmark seismic experiment, passive seismicity was recorded using the same OBS array as for the anisotropy experiment (Fig. 3a).
Various types of events were detected using this approach, including: short (<1 s) duration events, similar to those previously detected in gas seepage areas (Tary et al., 2012;Bayrakci et al., 2014;Batsi et al., 2019), and which may represent collapses and/or the formation of small pockmarks; medium (few seconds) duration events, comparable to a volcano-tectonic tremor (e.g., Latter, 1981;Harrington and Brodsky, 2007), rarely observed in this tectonic context, and which may represent the movements of fluids in conduits in the shallow subsurface; and events related to bubbles escaping at the seafloor. The low frequency content displayed by these types of events (<30 Hz) means it may also be possible to identify these events during periods of, for example, sparker seismic acquisition, due to the much higher frequency content of the source signal (≥200 Hz).
For recorded seismicity associated with gas bubble escape at the seafloor, the Minnaert (1933) equation can used to estimate the bubbles' radii from the frequencies of the events. However, this approach may not be able to distinguish between the acoustic signature of many small bubbles and one big bubble, which highlights the necessity to calibrate this approach using additional methods (e.g., active acoustics).
Passive monitoring of seismicity associated with methane venting and/or fluid movement in the subsurface at Scanner pockmark aims to quantify the gas flux through the chimney to the seabed and into the water column. Unlike active acoustic methods, this approach represents a continuous monitoring tool. The microseismic events detected through passive monitoring may be associated with the local tidal cycle and, hence, may indicate a possible correlation between tides and the movement of fluids from a subsurface reservoir to the seabed. Furthermore, if the seismicity is a result of fluid movement at depth, it may be applicable to studies of fluid flow beneath the cap rock or overburden, where the fluids do not directly reach the seabed.

Active acoustic methods
Hydro-acoustic surveying involves the use of either single-or multi-beam echosounders, which transmit one or more beams, respectively, of monochromatic frequency beneath the ship. The delay associated with the returning pulse is typically used to map the underlying bathymetry, while the strength or amplitude of the returning pulse can provide information on the physical properties of the reflecting medium. Given the strong impedance contrast between water and gas, hydro-acoustic methods are sensitive to the presence of gas bubbles in the water column (Fig. 6). The gas bubbles rise due to buoyancy, drawing surrounding fluid into and upwards with the plume, with lateral dispersion as a function of bubble rising velocity and water current velocity. The gas plume comprises bubbles of various sizes, with larger bubbles possessing relatively higher rising velocity and smaller bubbles possessing relatively lower rising velocity. The plume can be observed to rise to the thermocline at ~40 m below sea surface ( Fig. 6; Böttner et al., 2019;Li et al., 2020).
Multi-beam data can be used to map the three-dimensional extent of gas plumes escaping from pockmarks and, across multiple surveys, observe variations over time. Seabed currents may cause the plume to separate based on differing bubble size distributions, which can be mapped using multi-beam data (Li et al., 2020). Similarly, observing the target strength of gas plumes in single-beam data can be used to estimate the bubble size distribution and gas flux being released from the pockmark, provided a range of frequencies covering the bubble resonance frequency are used (Li et al., 2020).
Active hydro-acoustic surveying was applied in the study of the active venting at the Scanner pockmark in order to quantify the gas flux from the seabed and to study the gas plume evolution in the water column. Repeat surveys can also be used to determine if there is any temporal variation, for example due to tidal cycles (e.g., Boles et al., 2001;Schneider von Deimling et al., 2010). Temporal patterns may also be compared with observations from passive seismics, which provide a continuous record of gas venting at the seabed, or of the migration of fluids at depth.

Drilling and coring
Direct sampling of seabed sediments via coring allows for a local characterisation of lithology, porosity, and other geophysical and geochemical attributes at centimetre-scale or less resolution, providing a verification of geophysical properties and parameterisation for numerical models. Detailed description of sediments in terms of their lithology, physical, and chemical properties allows determination of the geological origin of the sediments, the extent and characteristics of specific stratigraphic units, the possible origin of the pore fluids, and a very high-resolution, visual observation of changes in sediment texture and fabric. Changes in sediment colour, e.g., the presence of black layers that form when methane is consumed by anaerobic methane oxidation with sulphate, and the structure of the sediment, e.g., the presence of outgassing and fluidisation structures, provide evidence for the presence of gas in the sediment, as well as evidence for palaeo-fluid migration. Further, non-invasive 3D analysis of sediment structure is possible using XCT imaging techniques (Fig. 7b-e), which can highlight the spatial distribution of dense materials such as iron sulphides that can form during anaerobic methane oxidation (Fig. 7c,e). Geochemical analyses of pore fluids can be used to determine their origin and, hence, the presence of migration pathways and the connectivity of different layers. Geochemical analyses of the solid phase allow understanding of the chemical reactions that are occurring or have occurred in the sediments and can explain physical properties, such as changes in permeability caused by carbonate or iron sulphide precipitation. Methane derived authigenic carbonates in the sediments and on the seafloor are indicators for prolonged presence of fluid flow, and can be dated to derive information about fluid flow history.
Sediment cores from the Scanner pockmark complex and a reference site located ~7 km to the east (Fig. 1b,c) were recovered using the BGS Rockdrill (RD2) and a gravity corer during cruises MSM78 and POS518-2 (Fig. 7b). RD2 is a remotely operated multi-barrel wireline subsea robotic sampling system that continuously cores 1.72 m-long sections, with a diameter of 6.1 cm. The maximum drilling depth was 33 mbsf at the reference site. Gravity cores were collected from the uppermost 5.75 m of the sediment column using a 12.5 cm-diameter tube driven into the seabed under weight-assisted free-fall (Linke and Haeckel, 2018).
Immediately following collection of the sediment cores, pore waters and sediment samples were extracted from both RD2 and gravity cores.
Pore waters were collected with 0.2 μm-pore diameter Rhizons (Seeberg-Elverfeldt et al., 2005;Dickens et al., 2007) at approximately 30 cm depth intervals, and preserved for analyses of cations (e.g., dissolved Ca, Fe, Mn, Si, B), anions (e.g., Cl, SO 4 , total alkalinity), nutrients (NH 4 , PO 4 , Si, NO x ), and dissolved sulphides. Samples of the sediments were collected with cut-off 5 cm syringes and analysed for their methane concentrations and methane carbon isotopic composition, porosity, grain size, and organic and inorganic carbon content. The cores were then scanned using a GEOTEK multi-sensor core logger (MSCL; Fig. 7a), to measure acoustic wave velocity, density, resistivity, and magnetic susceptibility, before being sectioned horizontally, photographed, and logged. The chemostratigraphy of the sediment was then measured in high resolution (1 mm) with the ITRAX XRF core scanner.
Some caution is required during core logging and interpretation, as artefacts such as small-scale fractures and other disturbances can be Fig. 6. Hydro-acoustic imaging of the Scanner pockmark gas plume, using EM712 echosounder data at a frequency range of 40-100 kHz. Highest backscatter is displayed in red (seabed). The gas plume is clearly visible as a strong (green) zone of backscatter emanating from the pockmark, rising to the thermocline at ~40 m below the sea surface. Figure  introduced during the drilling procedure. This is particularly important for soft, non-lithified sediments such as those found in our study area, as the coring process may reduce porosity, due to dewatering compaction and shear straining (Bouma and Boerma, 1968;Blomqvist, 1985). Concentrations of some reactive geochemical parameters (e.g., dissolved Fe, dissolved inorganic carbon, H 2 S and CH 4 ) may further be affected by outgassing or chemical oxidation during the retrieval and subsampling of the cores. These effects can be minimised by subsampling immediately after core recovery, or by the use of in-situ measurement techniques such as lab-on-chip sensors (e.g., Whitesides, 2006).
While coring and drilling approaches allow very high (cm-scale) vertical resolution, the horizontal resolution is limited by the number of sediment cores that can be taken within the time frame of an expedition. Therefore, appropriate sampling strategies are required to effectively integrate the information from sediment cores with information gathered using other approaches. For the purpose of our study, previously acquired seismic, real-time fast-track processed seismic lines, and prior information about the shallow geology were used to inform our drilling strategy. However, at the Scanner pockmark drilling locations, RD2 settings were configured for recovery of authigenic carbonates, rather than unconsolidated sediments, and, as a result, sediment recovery was rather low (~50 %), due to washout of loosely consolidated sand intervals.
The integration of sediment core and geophysical data maximises their inherent value, by providing localised sub-seismic resolution sitespecific constraints (Vardy, 2015;Vardy et al., 2017;Provenzano et al., 2018). Multi-sensor core logging of sediment cores permits the measurement of geophysical and geotechnical parameters, providing a way to: calibrate the two-way-time seismic reflection data to depth; constrain seismic inversion and interpretation within a range of plausible values (e.g., Provenzano et al., 2018); interpret and cross-validate the results of seismic tomography and CSEM methods; and interpret physical properties in terms of changes in sediment properties, such as porosity (e.g., Vardy, 2015;Baasch et al., 2017). Appropriate calibration of the instruments is necessary for a meaningful comparison between samples acquired using different coring techniques and core liner materials, and the need to account for laboratory temperature and humidity conditions.

Geomechanical testing and sediment core analysis
Rock physics experiments permit the direct examination of the changes in geophysical, geomechanical, and geochemical sample properties under controlled conditions in the laboratory. Two rock physics experimental setups were used in combination for the study of coupled thermo-hydro-mechano-chemical phenomena occurring in a rock-brine-CO 2 system: (1) multi-flow CCS and (2) acoustic pulse tube rigs.
Firstly, the multi-flow CCS rig permits the analysis of geophysical signatures (ultrasonic P-and S-wave attributes, and electrical resistivity) of reservoir rocks subjected to brine and CO 2 flow-through tests at controlled reservoir pressure and temperature conditions, resulting from the hydromechanical changes occurring in the rock sample due to CO 2 -brine fluid substitution effects and rock-fluid reactivity. Secondly, the acoustic pulse tube rig allows the measurement of seismic velocity and attenuation of P-waves in sediment and rock samples over a frequency range of 1− 10 kHz and S-waves at 500 Hz. This frequency band complements the ultrasonic frequencies obtained by the multi-flow CCS rig and helps to constrain improved rock physics models accounting for fluid mobility and frequency-dependent seismic effects (Batzle et al., 2006;Mavko et al., 2009). In addition, the P-and S-wave frequency ranges of the acoustic pulse tube measurements partially overlap with the sparker sources used for the seismic anisotropy assessment of the Scanner pockmark (Fig. 3b), permitting the calibration of the seismic data under controlled conditions in the laboratory.
The experimental approach used provides a geophysical and hydromechanical assessment of fluid flow in (i) CO 2 reservoir rocks and (ii) the shallowest section of the marine sediment overburden. In the former case, the reservoir rock experiments simulate variable spatial-temporal states during (drainage) and post-(imbibition) CO 2 injection, considering homogeneous (background) and fractured (slip-related features) reservoir samples. In the latter case, the tests focus on the hydrodynamic response of near-seabed sediments subjected to variable CO 2 flow rates.
A study of reservoir rocks was conducted using natural reservoir samples (Utsira sand formation, Sleipner; Falcon-Suarez et al., 2018), and quartzose synthetic sandstones manufactured in the lab following the procedure described in Falcon-Suarez et al. (2019), including homogeneous (Falcon-Suarez et al., 2016, 2017 and fractured samples (Muñoz-Ibáñez et al., 2019). In all tests, an increasing CO 2 -to-brine partial flow condition was imposed through the originally brine-saturated samples, in order to progressively saturate the sample in CO 2 . The confining stress, pore pressure, and temperature conditions slightly varied between tests, within the ranges 16− 40 MPa, 7− 12 MPa, and 20− 35 • C, respectively. The fluid substitution process and related hydromechanical response of the rock were assessed though the continuous monitoring of ultrasonic P-and S-wave velocity and attenuation, electrical resistivity, axial and radial strains, and hydraulic conductivity. Fig. 9 shows a comparison between the results obtained using a synthetic (Falcon-Suarez et al., 2016) and a natural  sandstone, with similar composition and porosity. The arrival of CO 2 in the rock samples lead to increased resistivity, together with a Fig. 7. Gravity cores of sediments recovered from directly beneath East Scanner Pockmark (GC17; Fig. 1c) during the MSM78 expedition. The cores were analysed using a) Multi-Sensor Core Logging (MSCL), which includes density data and b-e) X-ray micro-CT analysis. b) and d) Analysis of a 25 cm length x 12.5 cm diameter core section composed of fine silt. Coarser grained horizons (lighter grey) show evidence of fluidisation structures and layering, picked out by the dotted lines. c) Same image as b) but with lower density material removed, showing the presence of disseminated iron sulphide (bright spots, close-up shown in e) that has precipitated within the coarser grained horizons. The presence of iron sulphide is indicative of transport of methane-rich fluids through the coarser grained sediment horizons.
sudden drop (>5%) in P-wave velocity and increase in the attenuation (more significant in the synthetic sample), while subsequently increasing the CO 2 content had only limited further effects. Forced imbibition is used to study the trapping efficiency of the rock through evaluation of the residual CO 2 content and the geophysical signatures associated with this value. In this regard, Muñoz-Ibáñez et al. (2019) studied the hydrodynamic behaviour of a fractured reservoir analogue sample, and found that non-connected fractures may act as larger cavities than background porosity, and these cavities are energetically favourable for CO 2 storage. Therefore, non-connected fractures have the potential to enhance CO 2 trapping and the overall sealing efficiency of CO 2 storage reservoirs (Muñoz-Ibáñez et al., 2019), in contrast to the role played by a connected fracture network (i.e., the hypothesised structure of a chimney).
In addition, gravity core samples from cruise POS527 (Achterberg and Esposito, 2018) were used to conduct flow-through experiments, to analyse the coupled hydro-mechanical-chemical response of seafloor sediments subjected to sudden gas release episodes (e.g., pockmarks). These experiments included (i) a set of permeability tests to monitor the sediment reactivity to CO 2 (i.e. dissolution/precipitation), and (ii) a simulation of the STEMM-CCS CO 2 release experiment (Flohr et al., 2021) to collect linked geophysical-hydromechanical information about the target site under controlled conditions in the laboratory. The data collected during these experiments indicate permeabilities of ~10 − 16 m 2 in the upper part of the sediment column. Preliminary analyses of reactive transport models  indicate that the sediment is largely unreactive to the exposure of CO 2 saturated brine, while the grain size may be conditioning the CO 2 distribution within the sediment; from clay-to sand-rich layers the distribution becomes more homogeneous (Roche et al., in review).
For further understanding of hydraulic fracturing and gas migration mechanisms, flow-through experiments were carried out on gravity core samples recovered from locations adjacent to the STEMM-CCS CO 2 release site during the RV Poseidon cruise POS534 (Schmidt, 2019). Combining geomechanical testing with XCT imaging (Fig. 10d; Deusner et al., 2016) allows the visualisation of fracture formation in response to pressure accumulation and fluid injection. Fracture formation in the poorly consolidated clay sediments started at a low critical pressure of ~5− 15 kPa above the confining stress (Fig. 10a). Fracture propagation was essentially perpendicular to the major principal stress direction, likely resulting from the formation of disk-shaped gas bubbles with high aspect ratios, as expected in clay sediments (Boudreau et al., 2005). Fracture formation was first observed from rapid gas pressure dissipation (Fig. 10b, event 4) and the appearance of a small gas accumulation in the XCT image projections (Fig. 10c). Mass balancing of the injected gas and changes in confining volumes after fracture opening until the end of the first injection interval indicate a high apparent fracture permeability (Fig. 10b). Flow-through experiments and coupled thermo-hydro-chemo-mechanical studies can provide important insight into fracture formation and propagation in natural sediments, and are an important tool to develop and validate numerical models of chimney evolution. However, mechanical and hydraulic constraints in small-scale laboratory experiments will not provide an exact replication of the in-situ field environment, which must be carefully considered during experimental planning and data interpretation.

Onshore outcrop analogues
In addition to direct sampling of an actively venting system, further characterisation and the long-term development of fluid escape structures can be determined from analysis of onshore outcrop analogue systems. Field outcrops permit characterisation of physical structure, geometry, and physical properties, and understanding of past physical, chemical, and hydro-mechanical processes, through direct observation and sampling. However, fluid escape outcrop analogues are challenging features for comparative analysis with active subsurface systems (e.g., Roberts et al., 2017). Following the cessation of active fluid flow, a number of diagenetic and surface processes may alter the physical properties of the rocks, which can include post-dewatering compaction, cementation, groundwater-induced geochemical alteration, weathering, and erosion. XCT techniques can be utilised to acquire 3D digital rock images of the outcrop samples, whereby the cement can be selectively removed using image processing techniques, allowing the study of the physical and hydraulic properties of the rock formation prior to cementation (Fig. 8d; e.g.. Böttner et al., in review).
There is often a discrepancy between the scales of fluid flow features observed in onshore outcrops and their seismic manifestations, whereby chimneys observed in seismic data can extend several hundreds of metres in diameter and up to 1000 m in height (Berndt, 2005;Cartwright et al., 2007;Løseth et al., 2009;Andresen, 2012;Karstens and Berndt, 2015). Outcrops of sand-filled fractures in California (Fig. 8c), which are sand intrusions that form due to exceedance of lithostatic pressure and subsequent remobilisation and fluidisation of sediments (Huuse et al., 2010;Hurst et al., 2011), have been observed to reach similar dimensions to their corresponding seismic discernible features Ross et al., 2014). However, the largest known field outcrops of hydrocarbon-derived carbonate conduits, located in Varna, Bulgaria, are only up to 10 m high and ~1 m in diameter, and are spread over a comparatively small area (2.4 × 10 5 m 2 ; De Boever et al., 2006aBoever et al., , 2006b. The lack of outcrop evidence for fracture systems that extend up to 1000 m in vertical height and with large spatial extent is likely a consequence of the preferential preservation of hard parts, as suggested by the observed sand-filled fractures and carbonate pipes in California and Varna respectively (De Boever et al., 2006a;Hurst et al., 2011). However, a benefit of field outcrop analogues is that sub-vertical structures can be mapped at mm-scale resolutions, which is not possible using seismic reflection imaging. Consequently, more accurate quantification of fault/fracture properties, such as measurement of aperture or cement infill type, is achievable using outcrop based studies. The study of outcrops can be complemented with the mapping of a target area, using unmanned aerial vehicles (UAVs) that allow mapping over large distances (km-scale) at very high accuracy (cm-scale). These data can then be used to calculate digital elevation models, image mosaics, and high-resolution point clouds which, once geo-referenced using ground control points, can help close the observational gap between seismic data and geological field sampling.

Process-based numerical modelling
Process-based numerical modelling, for simplicity called hereafter numerical modelling, can be used to provide quantitative insights into: how chimneys form and propagate through time and space; the distribution and transport of mass (water, CO 2 , CH 4 , higher hydrocarbons, salts, grains etc.) and energy (heat) components among different phases (solid, liquid, gas) through pores and fractures; the complex coupled physical and chemical processes between the above components and properties; and the changes in petrophysical (e.g., porosity and permeability) and chemical properties (e.g., dissolved gas distribution in the formation fluid brine). Therefore, numerical modelling approaches are useful in helping to quantitatively understand the governing processes that can explain observations, but they are, in turn, reliant on appropriate calibration relating to the target of interest.
As part of the STEMM-CCS and CHIMNEY projects, reactive transport modelling has been conducted (by Marín-Moreno et al., 2019) to examine the migration of injected CO 2 through a generic chimney, using Scanner pockmark as a proxy for pressure, temperature, geological, mineral, and chemical conditions. In the model, the chimney extends upwards in the subsurface reservoir through the cap rock to a shallowest depth of 400 mbsf, where it meets the overburden (Fig. 11a). Modelling was conducted in a radially symmetric geometry approximation, with injection of CO 2 taking place at 500 mbsf and 0 m radial distance. The time scale of this analysis is 100 years, longer than potential geological CO 2 injection activities.
Modelling indicates that, over these timescales, the porosity (Fig. 11b) and permeability (Fig. 11c) changes due to dissolution and precipitation of minerals are negligible for the representative chimney being modelled (Marín-Moreno et al., 2019). Hence, for this particular location and over the time scale of this analysis, the modelling of chemical reactions (Fig. 11d,e) can be dismissed. For injected CO 2 saturations of ≥30 %, and with overburden isotropic matrix permeabilities and chimney effective vertical permeabilities (average permeability of the matrix and the fractures in the chimney) of ≥10 − 14 m 2 , the models indicate that CO 2 reaches the seabed, at 500 m above the injection point, in less than 100 years (Fig. 11a,f).
The modelling conducted assumes that the injection of CO 2 occurs far away from the chimney, in isothermal conditions, and that there are no mechanical changes . However, the non-isothermal upwards migration of CO 2 may create stress-and thermal-induced deformations that could affect the aperture of existing fractures, or create new fractures, and so, in turn, affect the migration processes. Thus, assessment of the impact of thermal, hydraulic, and mechanical changes of fractures within chimneys may require further consideration.
The dynamics of fluids migrating within the shallow overburden sediments can also be assessed using a Navier-Stokes-Darcy multi-fluid flow model (AnsdMF), which simulates the mass and momentum exchanges of multi-fluid phases (Behzadi et al., 2004;Bannari et al., 2008) and couples the turbulent flows in the ocean with the flows in sediments (Fig. 12a). The AnsdMF modelling was conducted in 2D, and comprises a sediment layer with a defined thickness, porosity, and fracture properties (length and aperture; Fig. 12b). These subsurface characteristics may be further constrained using other methods discussed in this study. Above the seabed is a turbulent bottom boundary layer. The location and rate of gas release beneath the seabed are then specified and the gas is tracked as it migrates and disperses through the sediment, and is expelled from the seabed (Fig. 12c).
The AnsdMF model is robust and efficient for simulations of multifluid flow in sediments with a complex structure and variable fracture properties, and fluid pathways to a turbulent ocean. In addition, AnsdMF can be applied over a wide range of scales, ranging from micro/ pore-to field-scale (≤10 − 2 to ≥10 3 m). The observations from computational fluid dynamic models can be calibrated against observations derived from active seismics, which provide a measure of fluids expelled from the seabed, and passive seismics, which can additionally detect flow within the subsurface, to compare flux rates of subsurface fluid migration with observed discharge from the seabed. Further, bioecosystem models could be coupled with AnsdMF and advanced oceanic hydrodynamic models to assess the impacts of fluid migration on benthic ecosystems and in the water column (Blackford et al., 2020).

Parameters
In this paper, a number of methods have been identified which can be utilised and integrated to generate a comprehensive characterisation of subsurface fluid pathways associated with chimneys. The parameters influencing fluid flow within chimneys can be broadly divided into three categories: 1) structural properties including the subsurface stratigraphic structure, the presence, geometry, distribution and connectivity of fractures, and principal stress magnitude and orientation; 2) physical properties of the sediments, including thermal, hydraulic, mechanical, elastic, and electrical properties; and 3) fluids, which includes their presence, phase, and distribution in the subsurface, chemical composition, and volume/mass flux calculations (Table 1; Fig. 13). In general, at least one approach from each of the method groupings (Table 1 left-hand column) is required to reach a sufficiently detailed and holistic understanding of the structure and processes involved in the venting of gas via a chimney.
We have also identified five method traits: resolution, length scale, temporal applicability, co-dependency, and resource use. The method parameters and traits are summarised in Table 1 and are further discussed below. Lithostratigraphy of the Panoche Hills which includes the Panoche Fm., comprised of up to 6.5 km of turbiditic sandstone, overlain by greater than 600 m of mudstone and siltstone cap rock/overburden sealing units. c) Map view of sand filled fractures intruding sub-vertically through the sealing units of the Moreno Fm. The map view is analogous to the spatial extent mappable with UAVs. d-f) Direct sampling and observation of a sand-filled fracture at different length scales, which includes d) pore-scale XCT analysis and e) core-scale laboratory analysis. Table 2 shows (a) the order of magnitude resolutions and (b) length scales (area and/or depth) achievable for each method applied to the study of a chimney. Table 2a demonstrates that there is a resolution gap between the geophysical (seismic, acoustic, CSEM) and geological/ geochemical (rock sampling, laboratory) methods around the m-scale. Numerical modelling approaches allow bridging of the resolution gap, by the integration of data from the geophysical, geological and geochemical methods.

Resolution and scale
There is an inherent trade-off between resolution and length scale, whereby methods capable of achieving higher spatial resolution can typically only be applied at reduced length scales. This is particularly an issue in the vertical dimension. This trade-off arises in seismic and CSEM studies due to an intrinsic relationship between the imaging resolution and depth scale, whereby lower frequency sources have greater depth penetration, due to the stronger attenuation of higher frequency signals over lower frequencies (Resnick, 1990). However, vertical imaging resolution is also related to bandwidth, with lower bandwidth signals resulting in a lower vertical resolution than higher bandwidth signals. Taking the seismic approaches utilised in this study as an example, the lower frequency airgun sources generate the best depth penetration (≥100 s of m), but are limited to high m to 10 s of m or greater vertical resolution, while sparker sources may only penetrate a few 100 s of m, but can generate m-scale resolution, and sub-bottom profiler recording (which we include here with seismic approaches as they are applied for the same purpose) can achieve resolutions on the dm-scale, but with only very limited (10 s of m) depth penetration. Therefore, where frequency-dependent attenuation is a factor in measurements made using a particular method, Table 2 must be considered in this context.
An optimal subsurface imaging approach would assess parameters at the highest resolution possible over the areal and depth extent of the chimney. However, where these features are located at depth, beneath overburden and sealing layers, attenuation of higher frequencies results in a depth-resolution trade-off, with increasingly only larger features visible. This represents a key motivation both for the application of a multi-frequency surveying approach in the Scanner pockmark seismic experiment, utilising a range of seismic sources (e.g., airguns, sparkers), and the study location chosen. Together, these allow the development of a cross-scale, sufficiently resolved image of a chimney at shallow depth, where it can be studied in detail, with other techniques (laboratory experiments, numerical modelling) then able to extend the observations to the deeper subsurface, where they may be detected by geophysical techniques at lower resolution. The use of a broad seismic source frequency range for the anisotropy experiment is also critical, due to the frequency dependence of anisotropy detection and interpretation (Chapman, 2003).
In addition, some of the methods display differing sensitivities to the magnitudes of the properties being measured (e.g., Constable, 2010). For example, while the elastic properties of sediments are sensitive to both lithological and pore fluid variations, their electrical resistivity is dominated by the pore fluid component. The differing response of these two techniques as gas saturation varies has been shown by Constable (2010) for a 50 % porosity sand. The largest effect on P-wave velocity occurs for the first few percent of gas fraction and changes very little above this (White, 1977;Lee, 2004). In contrast, the resistivity response increases exponentially with increases in gas saturation (Archie, 1952). As a result, gas saturations determined using seismic and CSEM methods may differ (e.g., Goswami et al., 2015;Attias et al., 2020;Bialas et al., 2020), with CSEM being more sensitive to the wider saturation-range response measurable.

Temporal applicability
The integrated methods described in this paper also influence the temporal scales over which fluid flow through chimneys can be understood (Table 1). Only a few of the methods can be used to determine temporal information related to pore fluid changes in the subsurface (measured at discrete time intervals using time-lapse/4D seismic), fluxes through the subsurface (measured discretely by active acoustics, and measured/modelled continuously using passive seismics and computational fluid dynamics), and fluid-induced chemical and physical changes (geochemical methods and numerical modelling).
Repeat geochemical data collection and seismic surveying can be conducted, commonly on the time scale of every year or few years. However, measurement of true temporal changes using seismic data, rather than seismic artefacts, is limited by the repeatability of the surveying technique. This makes time-lapse seismic studies most suitable for application to areas of well stratified and undisturbed sediments, rather than chaotic features such as chimneys (e.g., Waage et al., 2019). Developments in time-lapse seismic repeatability, such as continued use of true-4D imaging, or the installation of permanent seabed sensing infrastructure such as OBS or fibre optic cables, have the potential to improve this.
Passive seismic and active acoustics can resolve temporal flow variations on the order of seconds to weeks, providing information on the ascent of individual bubbles, but also how flow is modulated by diurnal tides. Of these, however, only passive seismic provides continuous temporal information over extended periods of time. These two methods have, therefore, different applicabilities to short-and longer-term monitoring activities of fluid escape associated with chimneys. Both methods also play an important role in the calibration of computational fluid dynamic models of fluid migration through the chimney and out of the seabed. These models can then be adapted to provide further Fig. 9. Laboratory flow-through experimental results. Brine− CO 2 flow-through tests were conducted for two samples, a 38 % porosity synthetic sandstone (Syn) and an intact core sample from Utsira Sand formation (Sleipner field). The tests simulate cyclic overpressure scenarios (ΔP P ) in a Sleipner-like CO 2 storage reservoir during pre-injection, CO 2 injection (drainage) and post-injection (imbibition) stages. Electrical resistivity, P-wave velocity, V P , and P-wave attenuation, Q P information about the temporal evolution of fluid flow processes through chimneys in different geological settings.
Of the methods that can be used to assess changes in the properties of chimneys over time, only numerical modelling approaches can be extended to much longer timescales (centuries to millennia), compared to experimental, field or laboratory techniques. Field-based outcrop analysis may also be considered to provide additional temporal information, allowing understanding of flow over geological timescales, long after active flow has ceased. However, the integration of these techniques enables the development of new multiscale coupled seismichydro-mechano-chemical modelling approaches, with the potential of enhancing predictions of chimney behaviour over time. Table 1 shows that a range of geophysical methods can be used to determine the site parameters of interest to the study of fluid flow structures. However, the successful application of individual methods relies on a network of co-dependencies between different methods for calibration and constraint (Table 1). In particular, inverse methods, such as travel time tomography and CSEM, are subject to inherent trade-offs between model parameters, such as velocity or resistivity respectively, and layer thickness (Edwards, 1997). Travel time tomography is also problematic when trying to resolve low velocity zones, as the energy is diverted away from them. As a result, use of prior geological knowledge, or constraints from other data, is important in the modelling process (e. g., Harris et al., 2009;Goswami et al., 2015;Attias et al., 2016).

Co-dependencies
Some co-dependencies are more significant than others. Seismic reflection imaging provides constraints on several of the methods utilised in this study, which explains its widespread use in subsurface exploration and monitoring programmes. Seismic reflection imaging and seismic inversion are co-dependent, although it is possible for these methods to be deployed fully independently. It is advantageous for travel time tomography to have input from reflection imaging of the Fig. 10. Fracture formation due to gas migration in poorly consolidated clay sediments. Flow through experiments were carried out on gravity core samples recovered during RV Poseidon cruise POS534 (Schmidt, 2019). Injection experiments were carried out with N 2 because of its low solubility compared to CO 2 . a) Gas pressure and confining stress. b) Pressure-volume relations during and after fracture formation. c) X-ray projection images at different time points. d) Setup of the X-ray micro-CT imaging system. Gas was injected in intervals (dark blue arrows) followed by pressure dissipation periods. Fracture formation was rapid after gas pressure increased above the confining stress (1). After reaching a critical pressure for fracture formation, an increase in confining stress (2) triggered confining fluid removal (3). The gas pressure dissipated rapidly after fracture opening (4), which subsequently resulted in the accumulation of gas at the sample surface. The offset between confining stress and gas pressure during ongoing gas injection (5) likely resulted from measurement inaccuracies because all system components were rated for high-pressure operation. Fig. 11. Reactive transport modelling of CO 2 migration through a generic chimney using the North Sea, Scanner pockmark complex as a proxy for pressure, temperature, geological, mineral, and chemical conditions. The model geometry comprises CO 2 injection at a radial distance of 0 m and depth below seabed of 500 m. The chimney extends vertically to the top of the cap rock (400 mbsf). Results are shown after 100 yr, with a constant injected CO 2 saturation of 50 %, an overburden isotropic permeability of 10 − 13 m 2 and chimney vertical and horizontal permeabilities of 10 − 13 m 2 , and 10 -15 m 2 . Variation of a) CO 2 saturation, b) porosity and c) permeability with depth and radial distance. Inset in a) shows the temporal evolution over 100 yr of CO 2 saturation with depth at a distance of 1.25 m from the injection point. Red numbers are CO 2 saturations. Changes in d) mineralogy and e) concentration of liquid species at a distance of 1.25 m from the injection point. f) Models of CO 2 gas reaching the seabed in less than 100 yrs for tested injected CO 2 saturations of 10 %, 30 %, and 50 %; overburden isotropic permeabilities of 10 − 13 m 2 , and 10 -14 m 2 ; and chimney vertical permeabilities ranging from 10 -11 m 2 to 10 -16 m 2 (ratios of chimney vertical to horizontal permeability from 1 to 1000 were tested). Figure   subsurface in order to constrain the solution of the non-unique tomographic problem towards realistic subsurface structures. However, it is also important that modelling is conducted in an objective manner that does not overly bias the result, and suitable modelling approaches can be applied which do not require/depend on this additional starting constraint. The velocity information which is obtained from travel time tomography or FWI can be used to improve the quality of seismic reflection imaging, including assisting with migration and depth conversion, in which accurate interval velocity information plays a central role. Furthermore, it can provide velocity information in areas of seismic blanking or potential artefacts, such as those associated with the presence of gas in the chimney. Converted-wave seismic anisotropy may also be partially reliant on seismic reflection data, to identify seismic horizons where mode conversion may occur, and an independent shearwave velocity constraint is also required to interpret from which layers these signals may be originating, and for layer-stripping (Haacke and Westbrook, 2006).
Laboratory experimental approaches can provide insights into the dominant processes governing the response of the system. The experimental strategy and the ability to approximate the real behaviour of subsurface chimneys, or other subsurface fluid migration pathways, depends on suitable calibration and constraints from field geophysical and geochemical experiments. Laboratory experiments are directly dependent on the quality of direct sampling. One way to minimise this dependency is to use synthetic samples to replicate subsurface material. In the study of Scanner pockmark, the main limitation of direct sampling was dewatering of the water saturated, poorly consolidated samples that Fig. 13. Summary of relationships between methods used and parameters investigated in the study of the chimney beneath Scanner pockmark. Triangle plot highlighting the three main parameters (Table 1) used for characterisation of a chimney/fluid escape system, where each method is located within the regions of the triangle corresponding to the parameters it can be used to resolve.

Table 2
Summary of the resolution (a) and the length scale (b) over which each geophysical method can be used to assess parameters, for complete characterisation of chimneys. Note that for methods subject to frequency-dependent attenuation with increasing depth penetration (e.g., seismic, CSEM), we show here the full range of resolutions and length scales. Thus, as depth penetration increases, the smallest achievable resolution and length scale will also increase.
likely occurred during the core extraction process, leading to alterations in the physical properties of the core material. Future optimisation of the core extraction process is required, to ensure improved preservation of poorly consolidated cored material. Autoclave pressure coring tools could be used, such as the IODP Pressure Core Sampler (Graber et al., 2002) or MeBo-MDP (Pape et al., 2017), that are able to retrieve sections of a drill core under in-situ pressure, thereby avoiding disturbances of the cored sediment due to degassing that occur during conventional coring. Hence, pressure coring permits the investigation of cored materials whilst retaining the in-situ physical, structural, and chemical properties (e.g., permeability, gas content, resistivity, salinity). Furthermore, the acquisition of borehole geophysical measurements (e. g., cross borehole surveys using wireline logging tools) is another common method for obtaining in-situ measurements, provided that the drilled borehole is stabilised by a liner (e.g., Crain, 2010). The results from rock physics experiments can be used to derive and validate the results of travel time tomography, FWI, and seismic anisotropy studies.
Multiphase and multicomponent numerical modelling of the thermal, hydrological, mechanical, and chemical response of marine geological systems with chimneys is currently limited by the robustness of available in-situ geophysical and chemical characterisation methods, especially for fracture geometrical, physical, and chemical properties, and the difficulty in accurately recreating these in the models. How well individual processes are represented in the numerical modelling results depends strongly on the in-situ geophysical and geochemical data available to constrain the model input parameters. This further highlights that an appropriate direct rock sampling strategy is an essential part of the multi-disciplinary approach for chimney characterisation.
Active acoustic and passive seismic methods for detecting fluid flow can be applied independently of other approaches (Fig. 13). A watercolumn sound velocity profile is required for accurate calibration of the shipboard echo-sounders, and, depending on the subsurface geological structure, a detailed velocity model may be necessary for accurate relocation of micro-seismic events. Cross-calibration of these two approaches is important, in order to extend the temporal observations of more punctual active acoustic surveys, and to resolve ambiguities in the apparent response of different bubble size distributions. These approaches also provide an important calibration for computational fluid dynamic modelling, either by indicating the potential flux entering the modelled system or measuring the expulsion from the seabed.

Resource use, challenges and development
It is also important to consider resource demand in the study of focused fluid flow conduits, such as the chimneys at Scanner pockmark, and how this may impact future studies relevant to CCS operations. Resources may include both direct monetary costs associated with data acquisition, in addition to less-well defined time-based costs, which cover the data processing and interpretation. Ultimately, there is a tradeoff between the field experimental methods which can be performed using a given set of resources (including time), and the total amount of information which can be determined about the subsurface structure.
Offshore site characterisations are generally spatially discrete snapshots in time. Undertaking studies of this nature requires field expeditions, including ship time, which results in relatively high monetary costs. Furthermore, the spatially discrete nature of offshore site investigations means that subsurface information must be gathered from all potential future sites of interest to CCS, to understand the geometry of the structures that may be present.
Depending on the availability of prior information about the specific imaging target, such as its spatial and depth scale, and the practicalities of implementing different imaging techniques, it may be possible to conduct experiments where multiple data sets are acquired concurrently, reducing the time-cost of the operation. For example, during the cruise JC152 seismic experiments at Scanner pockmark both the travel time tomography and seismic anisotropy datasets were collected using the same sets of shots, while the GI, Duraspark and Squid sparker sources were also recorded using a towed seismic streamer. Passive seismicity, for the study of fluid migration through the subsurface and ambient noise, was recorded using the same OBS array as for the JC152 tomography and anisotropy experiments (Fig. 3a). However, the parameters of a single seismic experiment may not be able to ideally meet the needs of all the different types of analysis that may be required for full and detailed characterisation of the subsurface feature(s) of interest. For example, for ambient noise studies, larger OBS station spacings would allow lower frequencies to be utilised in the analysis (e.g., Rost and Thomas, 2009). However, this conflicts with the smaller spacing required to achieve the maximum shallow subsurface resolution in tomographic models.
In a number of cases, the time-based cost of the techniques described in this study is relatively low, due to the relative simplicity of data processing using well-established methods (e.g., seismic reflection data, active seismics, core/field sample analysis). However, while many of the methods discussed in this study are well-established and extensively used in both research and industrial settings, there is potential for future reductions in resource costs. Methods such as travel time tomography, CSEM, passive seismics, and seismic anisotropy involve significant time in manual data picking and/or modelling, resulting in high relative time costs. Future developments in data processing may permit the replacement of manual processing stages with automated workflows, which may be facilitated by machine learning and artificial intelligence (e.g., Kim and Nakata, 2018;Kong et al., 2019), resulting in reduced time costs. Further, from the geophysical experiments at Scanner pockmark, it can be determined that significant processing time could be saved if OBS and OBEM spatial coordinates, degree of tilt, and orientations are more accurately constrained during data acquisition (e.g., Haacke and Westbrook, 2006). Installing permanent OBSs or fibre optic cables would allow easier and more cost-effective implementation of repeat seismic surveying, where multiple cruise legs are required. Compared to offshore exploration costs, laboratory experiments can be run at much lower expense, especially the measurement of basic parameters, while the investigation of the effect of parameter changes (e.g., temperature, pressure) on the system may also be tackled using carefully designed experiments.
Ultimately, the STEMM-CCS and CHIMNEY studies of the chimneys at Scanner pockmark provided significant information on the structure of these features, which may be applied in future studies. The multi-scale and multi-parameter design of these investigations demonstrates that the integration of subsurface geophysical characterisation of a chimney with laboratory experiments, process-based numerical modelling approaches, and onshore analogues is possible. We also anticipate that our findings will be applicable for future targets, in particular those located at depth, beneath the overburden and sealing cap rocks, where higher resolution subsurface imaging is much more challenging. Therefore, the findings of the STEMM-CCS and CHIMNEY projects have the potential to reduce future resource intensiveness, whereby the pre-existing knowledge of fluid flow systems gained from studying the Scanner pockmark complex can be applied to new sites, reducing the need for additional site specific data acquisition.

Relevance to geological CO 2 storage
The study of chimneys as potential fluid flow systems is highly relevant to the risk assessment of future CCS operations. For the assessment of a prospective CO 2 storage site, a site characterisation is first conducted to understand the potential containment risks, based on pre-existing seismic, well, and environmental baseline data (e.g., Furre et al., 2019b). As the project progresses in maturity, additional data may be acquired based on identified risks or uncertainties. Containment risks are commonly identified using a bowtie risk assessment, where monitoring and modelling activities are implemented as barriers to reduce A.H. Robinson et al. the likelihood of a top event (Dean and Tucker, 2017). A top event may be described as the unintended migration of CO 2 from a primary or a secondary sealing layer. In this context, chimneys could be considered as one such risk element. Fluid escape through abandoned/legacy wells and faults are other common risk elements that require consideration (IPCC, 2005).
If chimneys are identified as a potential threat from the initial site characterisation (based on existing data), stakeholders would be required to follow a process to determine unknown parameters which may contribute to the containment risk. Fig. 14 demonstrates a potential approach to this, where the findings of our multi-scale approach to determining chimney parameters (Section 4; Tables 1 and 2) are integrated into an existing Measurement, Monitoring and Verification (MMV) risk-based decision framework (Dean and Tucker, 2017). The outstanding risk question is defined as specifically as possible, and the unknown factors which contribute to this are separated from the parameters which are sufficiently well resolved from pre-existing resources. The remaining unknown factors, which may include fracture connectivity, permeability, and/or fluid volume flux and composition, are then targeted using appropriate methods drawn from Section 3, in order to mitigate the outstanding risks. The selection of which method(s) to use is based on each methods' suitability for determining the parameters of interest (Table 1; Fig. 13), at the appropriate resolution and scale for the problem (Table 2), and any potential dependencies it may have (Table 1). Tailored deployment of new resources is then conducted, in a monitoring step. This may, for example, involve further data acquisition and/or modelling studies. Alternatively, this could also include the re-processing or re-appraisal of existing data or models. In order to optimise the risk assessment process, this step should factor both the additional cost and time of further acquisition against the benefits of obtaining the new information using these approaches. The staging of resource deployment is also important at this stage, such that the maximum amount of new information is resolved using the minimum necessary additional resource cost/time. It would not, therefore, necessarily be either feasible or desirable to use all of the methods discussed in this study, highlighting that stakeholders need to efficiently determine the parameter-specific traits they require for appropriate de-risking of the chimney. The outcomes of the monitoring step are then used to re-appraise and verify the containment risk-level attributed to the chimney, before a decision is made. If this verification stage indicates that all of the outstanding risks have not yet been constrained sufficiently, the process should be repeated until they have. In this respect, the risk assessment and mitigation process can be considered to be a cycle (Fig. 14).
Each prospective CO 2 storage site will have unique characteristics and, therefore, will be different. However, potential risks need to be assessed holistically, which will result in different resource deployment choices for different sites. Based on the outline workflow we have described, and using Tables 1 and 2 as a guide, we can consider a number of simple examples which could be resolved using this approach. If a fracture network at a spatial scale of 10 s of m was deemed to be the greatest threat/uncertainty for prospective CCS operations, then seismic reflection or seismic tomography would be the most suitable method to use, while CSEM data may also assist in determining fracture connectivity and fluid content. If the sediment matrix porositypermeability or the effects of microfractures were the principal cause of concern, laboratory experiments and/or modelling, using appropriate field analogues to define the input parameters, would be the optimal approach. As a final example, if CO 2 emission and flux rate from the seabed over a spatial scale of >1 km were the greatest threat/uncertainty, active acoustics would be the most appropriate method to use.
Alternatively, further data acquisition may not be required. In the case of the Goldeneye prospective storage site, located in the Outer Moray Firth Basin, Central North Sea, a chimney-like feature was observed during the risk assessment phase. However, rather than further data acquisition, the first monitoring step applied was to conduct high resolution re-processing of existing data, alongside pre-stack interpretation. Together, these approaches resolved the issue sufficiently, finding that the observed feature instead represented a seismic imaging artefact caused by a glacial tunnel valley, without the need to acquire additional data (Dean et al., 2015).
Active natural methane escape is observed directly within the shallow subsurface of the exemplar chimney at the Scanner pockmark complex in the North Sea. However, prospective CO 2 storage reservoirs are more commonly located at greater depth and are capped by several sealing layers, for example the storage units at Sleipner are located at a depth of ~800− 1000 mbsf. Therefore, understanding the role of chimneys in applied CCS environments may also represent a deeper subsurface problem. Compared to shallow chimneys, deeper chimneys may exhibit a number of differences, such as 1) increased stress -fluid pressure -temperature conditions; 2) the fluid of interest may be in a different phase (liquid compared to gas); and 3) changes in the sediment Fig. 14. Outline of a decision-making cycle illustrating the steps in the risk-based assessment of a site where chimneys may pose a risk to potential CCS operations. The multi-parameter, multi-scale approach discussed in Sections 3 and 4 of this paper should be used in steps 3 and 4, the mitigating and monitoring steps, in order to optimise the decision making process. and pore fluid composition. At greater depths (higher effective stresses) there is likely to be more consolidated material, which may result in reduced porosity-permeability of the sediment matrix. Moreover, fracture apertures tend to decrease with increasing effective stress (increasing depth). Therefore, if a chimney existed at greater depths, it is more likely to have lower effective permeability, although this is dependent on the rock material properties and microscopic fracture roughness (e.g., McDermott and Kolditz, 2006;Rutqvist, 2015). At greater depths the chimney would need to contain considerable stress-persistent fracture apertures to be deemed a risk. Further research is required to more accurately constrain the differences expected between deeper and shallower chimneys. The study of deeper chimneys may also alter the experimental approach. For example, a number of the applied geophysical methods lose resolution with depth, due to the increasing attenuation of higher frequency signals that are required to image structures at shorter length scales (Resnick, 1990). However, by applying the principles discussed here and in Fig. 14, a multi-component, appropriate target-scale risk appraisal approach can be developed for a chimney structure located at greater depth.
Overall, our integrated multi-scale experimental approach can be incorporated into pre-existing measurement, monitoring and verification approaches (e.g., Dean and Tucker, 2017). Future monitoring approaches must also be supported by robust baseline surveys and detailed understanding of the risks of CO 2 emissions that extend into the marine environment (Blackford et al., 2014(Blackford et al., , 2015Flohr et al., 2021). Our understanding of fluid flow through chimneys must be improved to decrease containment risks at prospective storage sites and thereby realise the potential of CCS for mitigating global warming.

Summary
Chimneys (also referred to as seismic chimneys and/or pipes) may represent an important class of subsurface fluid escape pathway. The study and understanding of chimneys are highly important to risk assessments for potential CCS applications. To this end, the STEMM-CCS and CHIMNEY research programmes have investigated chimneys in the North Sea. The results of these investigations will contribute to the further development of models for the understanding of fluid flow through such structures.
In this paper we have described an approach for the multi-scale, multi-method geophysical characterisation of chimneys within sedimentary basins, using the Scanner pockmark as a case study. In the development of this multi-method approach, we note the following: • While individual experimental approaches, conducted in the field or in the lab, are able to determine one or more of the physical properties of chimneys and their role in fluid flow, integration between different methods is required to fully characterise the system. • Many of the described approaches rely on having suitable constraints and calibration from other methods. It is necessary, therefore, to ground-truth the methods, both in the field and in the laboratory, wherever possible, in order to ensure that conclusions determined are robust. • The method approaches available are applicable at different physical length scales and resolution, with a common trade-off between these two factors. Hence, it is also important that the techniques used to study fluid flow structures are applied across scale and resolution ranges, as well as different temporal scales, in order to build a holistic and sufficiently detailed understanding. • We propose that the integrated methodology presented to study the chimney beneath Scanner pockmark is applicable beyond this specific study area, and thus our work contributes to a wider understanding of this class of feature and its potential importance in fluid migration.

Declaration of Competing Interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.