Correcting for Static Shift of Magnetotelluric Data with Airborne Electromagnetic Measurements: A Case Study from Rathlin Basin, Northern Ireland

. Galvanic distortions of magnetotelluric (MT) data, such as the static shift effect, are a known problem that can lead to incorrect estimation of resistivities and erroneous modelling of geometries with resulting misinterpretation of subsurface electrical resistivity structure. A wide variety of approaches have been proposed to account for these galvanic distortions, some depending on the target area, with varying degrees of success. The natural laboratory for our study is a hydraulically permeable volume of conductive sediment at depth, the internal resistivity structure of which can be used to estimate reservoir viability 5 for geothermal purposes, however static shift correction is required in order to ensure robust and precise modelling accuracy. We propose a method employing frequency–domain electromagnetic data for static shift correction, which in our case are regionally available with high spatial density. The spatial distributions of the derived static shift corrections are analysed and applied to the uncorrected MT data prior to inversion. Two comparative inversion models are derived, one with and one without static shift corrections, with instructive results. As expected from the one–dimensional analogy of static shift correction, at 10 shallow model depths, where the structure is controlled by a single local MT site, the correction of static shift effects leads to vertical scaling of resistivity-thickness products in the model, with the corrected model showing improved correlation to existing borehole wireline resistivity data. In turn, as these vertical scalings are effectively independent of adjacent sites, lateral resistivity distributions are also affected, with up to half a decade of resistivity variation between the models estimated at

Abstract.Galvanic distortions of magnetotelluric (MT) data, such as the static-shift effect, are a known problem that can lead to incorrect estimation of resistivities and erroneous modelling of geometries with resulting misinterpretation of subsurface electrical resistivity structure.A wide variety of approaches have been proposed to account for these galvanic distortions, some depending on the target area, with varying degrees of success.The natural laboratory for our study is a hydraulically permeable volume of conductive sediment at depth, the internal resistivity structure of which can be used to estimate reservoir viability for geothermal purposes; however, static-shift correction is required in order to ensure robust and precise modelling accuracy.
We present here a possible method to employ frequencydomain electromagnetic data in order to correct static-shift effects, illustrated by a case study from Northern Ireland.In our survey area, airborne frequency domain electromagnetic (FDEM) data are regionally available with high spatial density.The spatial distributions of the derived static-shift corrections are analysed and applied to the uncorrected MT data prior to inversion.Two comparative inversion models are derived, one with and one without static-shift corrections, with instructive results.As expected from the one-dimensional analogy of static-shift correction, at shallow model depths, where the structure is controlled by a single local MT site, the correction of static-shift effects leads to vertical scaling of resistivity-thickness products in the model, with the corrected model showing improved correlation to existing borehole wireline resistivity data.In turn, as these vertical scalings are effectively independent of adjacent sites, lateral resistivity distributions are also affected, with up to half a decade of resistivity variation between the models estimated at depths down to 2000 m.Simple estimation of differences in bulk porosity, derived using Archie's Law, between the two models reinforces our conclusion that the suborder of magnitude resistivity contrasts induced by the correction of static shifts correspond to similar contrasts in estimated porosities, and hence, for purposes of reservoir investigation or similar cases requiring accurate absolute resistivity estimates, galvanic distortion correction, especially static-shift correction, is essential.

Introduction
The electrical resistivity of a volume of rock is highly sensitive to the presence of laterally and vertically varying amounts of electrically conductive fluids connected via pore spaces or fluid conduits.Due to these potentially strong resistivity contrasts between competent host rock and fluid penetrated rock, electromagnetic (EM) methods, and in particular magnetotellurics (MT), have been used with considerable success to image conductive volumes at depth (Chave and Jones, 2012;Simpson and Bahr, 2005).
Published by Copernicus Publications on behalf of the European Geosciences Union.
As with all EM methods, MT data are highly sensitive to rock fluid content and distribution (i.e.porosity and hydraulic permeability) and can be related to other properties relevant to fluid movement.This has made the method particularly interesting for the exploration of geothermal resources.Indeed, geothermal research was the first commercial application of MT in the late 1950s, though the interpretation of the corresponding conductive structures is not always straightforward (Muñoz, 2014).
The MT data set studied here was acquired in the context of a multidisciplinary geothermal research program (IRETHERM), the overarching aim of which is to identify and evaluate low-enthalpy geothermal resources within Ireland.One such resource (Goodman et al., 2004) is the thick, porous, and permeable succession of Permian and Triassic sandstones found within several concealed sedimentary basins in Northern Ireland, with the Rathlin Basin in particular having significantly elevated estimated geothermal gradients in comparison to the remainder of Ireland (Reay and Kelly, 2010).
The island of Ireland was formed during the Caledonian orogeny by the complex accretion of several continental and island arc fragments during the closure of the Iapetus Ocean between the Early Ordovician (485 Ma) and late Silurian (423 Ma), resulting in seven identifiable terranes that comprise the present-day basement across both Ireland and Great Britain (Mitchell, 2004;Hepworth and Sanders, 2009).Our survey area (Fig. 1) lies within the Central Highlands Terrane of Laurentia, the basement of which comprises mainly mid-to late-Neoproterozoic (1000-635 Ma) metasedimentary rocks classified as the Dalradian supergroup, a metasedimentary and igneous rock succession that was deposited on the eastern margin of Laurentia between late Neoproterozoic ( ≈ 800 Ma) and early Cambrian (≈ 510 Ma) times.Specifically, the basement across the test area is assumed to consist of the latest-Proterozoic (Ediacaran, 635-541 Ma) Argyll Group of psammites and semipelites.
Regional shear and stress during the subsequent late-Paleozoic (350-250 Ma) Variscan orogeny reactivated the Caledonian (490-390 Ma) Tow Valley Fault (TVF), and the ensuing normal and dextral strike-slip faulting resulted in the formation of a rift basin later filled by a succession of sediments to form the Rathlin Basin.Although drilling in the adjacent Magilligan Basin encountered Carboniferous formations at 1347 m total depth, the most basal formations confirmed within the Rathlin Basin are the Permian Enler Group (EG) sandstones and Early-Triassic Sherwood Sandstone Group sandstones (SSG).Both formations are hydrocarbon reservoirs in the Irish Sea to the east (Naylor and Shannon, 2011).The Sherwood Sandstone Group is overlain by the Late-Triassic Mercia Mudstone Group (MMG), which is itself generally overlain by late Jurassic Lower Lias Group (LLG) mudstones.However, in many places significant dolerite and basalt sills (up to approx.100 m in combined thickness) have been encountered, with poorly known spatial extent.The final and youngest successions in the basin are chalks of the Cretaceous Ulster White Limestone Formation (UWLF), with the Antrim Lava Group (ALG) concealing the basin entirely.
To date, two deep boreholes have been completed onshore in the Rathlin Basin, namely the Port More 1 (PM1) and Ballinlea 1 (B1) boreholes, drilled in 1967 and 2008 respectively; however, only data from the former are available as information from the latter is not yet in the public domain.The PM1 borehole was drilled to a total depth of 1897 m and terminated in the EG sandstones, with wireline log data acquired in two separate sections due to technical difficulties (Wilson and Manning, 1978).The upper portion of normal resistivity data covers the uppermost 250 m of the hole, including the Antrim Lava Group, Ulster White Limestone Formation, and the upper portion of the Lower Lias Group.These data provide relatively consistent resistivity estimates of ≈ 80 and ≈ 5 m for the UWLF and LLG sedimentary formations respectively, whereas estimates for the ALG vary from ≈ 5-80 m as it comprises a succession of tuffs and basalts.The lower portion of resistivity data spans the depth interval of ≈ 1050-1450 m, covering the boundary between the MMG and SSG at 1320 m, and provides resistivity estimates of ≈ 4 and ≈ 6 m for the respective groups.The estimates within the SSG also show a higher variance, which may be due to the presence of conglomerates and breccias in the upper portion of the group.The stratigraphic column encountered in the PM1 borehole is displayed in Fig. 2, alongside a plot of the borehole resistivity data.Modelling of regional gravity and magnetic data has been undertaken, with results presented in Mitchell (2004) and Gibson (2004).The density model of Mitchell (2004) (Fig. 3) shows a relatively homogeneous structure along the basin, particularly of the Permo-Triassic section, with a maximum depth to Dalradian basement of approximately 3 km modelled for the Rathlin Basin (located at 33 km distance along profile).This modelling adopts a density of 2.35 Mg m −3 for the assumed basal Carboniferous rocks; this value comes from borehole samples in the adjacent Foyle Basin, and Mitchell (2004) advises that this value may be insufficiently dense to represent Rathlin Basin conditions.If Carboniferous sediments in the Rathlin Basin are of a higher density, a greater thickness of overlying lighter sediments (i.e. the target Permian and Triassic sandstones) would be required to be consistent with the observed gravity anomaly.Magnetic and gravity modelling by Gibson (2004) suggests that the Tow Valley Fault zone consists of a series of major fault segments with varying dip of 20-50 • to the north-west.
Core samples of the EG and SSG sandstones successions gathered from the Port More 1 borehole show promising reservoir properties, with fractional porosities and hydraulic permeabilities ranging from 0.10 to 0.22 and 1-1000 mD respectively (Mitchell, 2004).Equilibrated temperatures taken from both the PM1 and B1 boreholes have previously been used to estimate geothermal gradients.A temperature of 35.4 • C was observed at 582 m depth in the PM1 borehole, whereas a temperature of 99 • C was observed at 2650 m in B1 (Reay and Kelly, 2010).Assuming a surface temperature of 10 • C, simple linear estimation gives calculated geothermal gradients of 43.6 (PM1) and 33.6 K km −1 (B1), both of which are elevated above the typical estimates of ≈ 20-30 K km −1 measured elsewhere across Ireland (Goodman et al., 2004).In conjunction with the promising reservoir properties and basin depth expected from gravity modelling, it has been proposed previously (Goodman et al., 2004;Reay and Kelly, 2010;Pasquali et al., 2010) that the Rathlin Basin may be favourable for geothermal exploitation.As the reservoir potential depends strongly on the intra-basin structure, variations in modelled resistivity may be taken as an excellent proxy for images of the presence of fluids, their distribution, and interconnection within the basin.
The imaging of sub-basalt structures poses difficulties to other commonly employed geophysical methods, particularly seismics (e.g.Martini et al., 2005;Bean and Martini, 2010) due to the negative acoustic impedance contrast at the base of the basalt, and previous reflection experiments struggled to clearly image the sediments through the overlying ALG (Naylor and Shannon, 2011).As the MT method has been successfully applied in sub-basalt investigations both onshore and offshore (Hautot et al., 2007;Jegen et al., 2009;Colombo et al., 2011;Heincke et al., 2014), the method was proposed to study the three-dimensional electrical resistivity distribution of the sedimentary fill of the onshore portion of the Rathlin Basin.
Due to the expected elevated hydraulic properties and saline pore fluids (both factors that increase conductivity) of the proposed hydrothermal aquifer within the basin, it was expected that MT data could be carefully modelled to im- age the properties and distribution of the aquifer formations.The increase in resistivity observed in wireline data from the MMG to the underlying target sediments implies that, depending upon the thickness of units beneath the MMG, MT may not be able to accurately estimate the units' resistivities, as MT is primarily sensitive to a layer's conductance (i.e. the ratio of a layer's thickness to resistivity), and thinner or less conductive layers may be shielded by overlying conductors (Jones, 1992).However, due to the proven thicknesses and similar lithologies of the SSG and EG sandstones, it is still expected that the SSG and EG sediment fill will cause a sufficiently high resistivity contrast against the resistive metasedimentary country bedrock and provide a viable target for MT.The geometry of this interpreted aquifer structure is expected to be compatible with the gravity model presented in Fig. 3.The MT method samples the impedance transfer functions that relate the electric and magnetic field components of EM plane waves that propagate into the Earth.As these EM waves attenuate with dependency on the Earth's lateral and vertical resistivity distribution, the observed MT responses can be employed for estimating the underlying 3-D resistivity distribution (Chave and Jones, 2012).The electric and magnetic field components of these EM source waves are each acquired in (preferentially) orthogonal horizontal directions, allowing the definition of four magnetotelluric transfer functions at the measuring location.These four elements carry information on the value and dimensionality of the subsurface resistivity structure at a range of periods.Many decompositions and analyses have been employed to expose this information (see, e.g., Chave and Jones, 2012, for an overview), with the aim of improving estimation or justifying 1-D (i.e.resistivity varying with depth only) or 2-D (i.e.resistivity varying with depth and one horizontal dimension only) modelling.
Though sensitive to conductive structures at depth, MT data are prone to distortion, primarily of the electric field, due to the presence of galvanic charges on the boundaries of shallow conductivity structures that are unresolvable at the frequency range of the recorded data.One simple form of this galvanic distortion is often easily identified by vertical offsets of the logarithmic apparent resistivity curves and is referred to as the static-shift effect (Jones, 1988), following a similar effect in seismology named "statics".These galvanic signatures are related to inescapable issues in observation of the electric field, wherein point-wise electric field observations, assumed during modelling and inversion, are replaced during MT surveys with voltage difference measurements along fi-nite length dipoles (Poll et al., 1989;Pellerin and Hohmann, 1990), and by issues related to insufficient gridding resolution to describe the lateral variability of the near surface that affect both the electric and magnetic fields (Chave and Smith, 1994;Chave and Jones, 1997).Although the former of these may be handled by appropriate post-processing when modelling and inverting the field measurements, any near-surface inhomogeneity not parameterised in the modelling or inversion process, even at the electrode-scale size, will contribute to the galvanic signatures by distorting the local (primarily) electric fields.For land-based MT data, the magnetic effects of galvanic distortion only occur for a short frequency range of the order of half a decade at most (Chave and Smith, 1994;Chave and Jones, 1997).
Various methods have been proposed to quantify and correct for these static-shift effects, including continuous sampling and filtering of the electric channels (Torres-Verdin and Bostick, 1992), spatial filtering based on mapping of MT data (Berdichevsky, 1989), modelling of parametric homogeneous layers at depth (Jones, 1988), estimation of distortionrelated parameters as unknowns during inversion (Sasaki and Meju, 2006;Miensopust, 2010;Avdeeva et al., 2015;De Groot-Hedlin, 1991;DeGroot-Hedlin, 1995), and finally the use of complementary EM geophysical methods (Sternberg et al., 1988;Pellerin and Hohmann, 1990;Miensopust et al., 2014).These methods can be broadly divided into methods that use intrinsic information from MT data and those that use extrinsic information from other geoscientific data.Whereas both families of methods can account for static shifts between MT modes at a single site and improve interstation shifts, intrinsic information may not yield a correct resistivity in the case of both modes being distorted -as stated in Sternberg et al. (1988), "there is no reason to expect that either of the two MT polarisations will provide the correct resistivity".In contrast, extrinsic methods using purely magnetic measurements (i.e. with no instruments using the subsurface as a component of their circuitry) by definition are unaffected by the electric effects of galvanic distortion and thereby provide a more correct resistivity estimate, albeit generally for much shallower depths than magnetotelluric measurements.The use of extrinsic information may, in some cases, be limited by the requirement that both the extrinsic information and MT data must illuminate some common depth volume within the Earth with a common current system in order to allow reconciliation of the resistivity structure.For example, if the near-surface is three-dimensional (3-D), then the current system from regionally induced currents observed in MT is very different from the current system from locally induced currents from a small EM transmitter-receiver array, and any resistivities estimated by MT would correspond to a different volume than that sampled by the smaller array.
Fortunately, in the case of the Rathlin Basin, broad-scale extrinsic resistivity information is available.As part of the regional Tellus ground and airborne geoscience mapping programme across both Northern Ireland and the Republic of Ireland, airborne frequency domain electromagnetic (FDEM) data were gathered over the target area (Beamish, 2013) in 2005and 2006. Sternberg et al. (1988) and Pellerin and Hohmann (1990) describe how time domain electromagnetic (TDEM) data can be used to estimate and correct static-shift distortion of MT data, assuming that the MT data are onedimensional (1-D) or two-dimensional (2-D) and in the correct geoelectric strike coordinate system, and we propose adapting their method to use the Tellus frequency domain airborne electromagnetic measurement (AEM) data for this purpose.Due to the high density of measurements often available, AEM data allow estimation of near-surface conductivity distributions at resolutions exceeding those of MT, making them an excellent choice for the correction of static-shift effects.The use of airborne EM data for static-shift correction has precedence, as outlined with airborne TDEM data by Crowe et al. (2013).
As three-dimensional inversion of MT data is becoming a more common practice, the effects of static-shift correction on resulting resistivity distributions must be considered.This article focuses on the implementation of this correction scheme for MT data by comparing models found by independent 3-D MT inversion of the observed MT data and the static-shift-corrected MT data.In addition to comparing the model results, the statistical and spatial distribution of calculated static-shift corrections are examined and compared to previous works to verify their validity.Both the absolute resistivity and resistivity gradients are used to evaluate the differences between the two models, with the implications of the differences discussed in the context of the possible geothermal aquifer.Note that the approach presented here is primarily methodological in nature; the geological interpretation is brief and will be investigated at length in a future publication.

Electromagnetic methods
EM methods include a wide variety of techniques that observe electromagnetic induction in the Earth and are most commonly used to image the subsurface distribution of electrical resistivity ρ ( m) or its inverse, electrical conductivity σ (Sm −1 ), and, at high frequencies, electrical permittivity ε (Fm −1 ).Subsurface materials are rarely homogeneous and can generally be described as a mixture of materials with strongly differing properties, as detailed in Nover (2005) and Chave and Jones (2012).The bulk resistivity of a rock is commonly determined by very few of its constituents, with electrically the most important candidates at crustal depths being metallic conductors (sulfides, graphite, iron oxides), clays, and conductive fluids in pore spaces (both saline fluids and partial melts).Note that we use pore space in a general sense for primary and secondary porosity, including pores, fractures, and conduits that may dominate the rock type under consideration.For low-enthalpy geothermal exploration www.solid-earth.net/8/637/2017/Solid Earth, 8, 637-660, 2017 in the upper crust, the most relevant property is the influence of electrolytic conduction by fluids in porous rocks.The relationship between the observed effective resistivity ρ, the pore fluid resistivity ρ i , and the formation porosity φ is classically described for sandstones by Archie's Law (Archie, 1947), with generalisations discussed in Glover (2010).As Archie's Law assumes a clean sandstone matrix with well-established relationships between porosity and hydraulic permeability, its application may not always be appropriate, particularly if clay minerals are present (Mavko et al., 2009;Guéguen and Palciauskas, 1994;Zinszner and Pellerin, 2007).
In a broader sense, electrical conductivity is a proxy measure of hydraulic permeability rather than porosity, as the interconnection of conducting pathways facilitates electric current flow.Due to this strong dependence on the geometry of flow paths on the scale of interest, the relationship between permeability and porosity is highly nonlinear (see, amongst others, Raffensperger, 1996;Pape et al., 1999Pape et al., , 2000;;Luijendijk and Gleeson, 2015).These dependencies have a close relationship with the type and geologic history of the rock considered (see, e.g., Bernabe et al., 2003), and the often complex development of geological units can lead to heterogeneities and preferential flow pathways at all scale lengths (Bjørlykke, 2010).

Magnetotelluric method
The MT method uses impedance transfer functions relating the electric and magnetic field components of vertically propagating EM source field plane waves to image the lateral and vertical resistivity distribution within the Earth.MT signal waves are generated by two sources, namely atmospheric electricity (generating signals of frequency > 10 Hz) and interactions of the Earth's magnetosphere with solar wind (generating signals of frequency < 10 Hz).Recent detailed reviews of MT methods, the underlying assumptions, and their application include Simpson and Bahr (2005), Berdichevsky and Dmitriev (2008), and Chave and Jones (2012).
Resistivity information at a range of depths is inferred by considering planar EM waves in the Earth at a range of frequencies, as their attenuation at a given frequency is a function of the resistivity ρ and magnetic permeability µ (Hm −1 ) of the subsurface material, where the latter is generally assumed not to vary from that of free space, µ 0 .In a uniform half-space (i.e. a space with no lateral or vertical resistivity variation) of resistivity ρ, the scale length δ in metres describing this attenuation, is termed the electromagnetic skin depth and describes the characteristic length over which the amplitude of an EM wave of frequency ω = 2πf decays by a factor of e −1 .This quantity is commonly used as a simple measure for the depth of investigation and radius of influence, although one must beware of its overuse in situations that depart from a uniform half-space, especially in the case of a multidimensional Earth (Jones, 2006).
The resistivity of a select volume of the Earth, as sampled by an EM wave of frequency ω, is determined from complex transfer functions that relate the amplitudes of the horizontal electric E i (Vm −1 ) and magnetic H j (Am −1 ) field components that constitute the wave, defined as the complex MT impedance tensor Z ( ).This is generally formulated in the frequency domain, where the transfer function can be defined by the ratio of the fields: By considering both orthogonal and parallel pairs of fields, an impedance tensor Z with four components (Z xx , Z xy , Z yx , Z yy ) can be defined: These impedances can be restated in more familiar magnitudes and units as an apparent resistivity ρ a (i.e.equivalent half-space resistivity for a wave of that specific frequency for the orthogonal pairs Z xy and Z yx ) and phase lead of the electric field over the magnetic field φ (which is π/4 for a uniform half-space for the orthogonal pairs), defined by In MT surveying, the electrical field E is measured as a voltage difference over an appropriate distance rather than at a point.For the sampling of the magnetic field components, the electromagnetic properties of the volume sampled within a sensor are known and homogeneous (being the internal properties of the sensor itself), and the averaged field sampled through this sensor volume accurately represents a point magnetic field value (to within a length scale of that of the sensor, typically 1.5 m for coil sensors used in broadband MT).However, the volumes sampled for the electric field components are of the order of 100 m (the typical bipole length in broadband MT (BBMT) surveys to acquire data in the frequency range of 100 to 0.01 Hz), with resistivity variations present within that length scale that may bias observations of the electric field.This distortion of primarily (but not exclusively) the electric field is one form of galvanic distortion and often manifests as frequency-independent multiplicative vertical offsets of MT apparent resistivity data for 1-D or 2-D cases (i.e.where impedance tensor diagonal components Z xx and Z yy are 0), when data are plotted as apparent resistivities on a log-log scale, hence the nomenclature of "static shift".As explained in, e.g., Jones (2011Jones ( , 2012)), whereas real multipliers applied to an impedance tensor with 1-D or 2-D form affect only the magnitudes of the impedances, this is not the case if the impedance tensor has a 3-D form with non-zero diagonal elements.Instead, the applied distortions cause mixing between the diagonal and off-diagonal components, affecting both the magnitudes and phases of the impedances.

Frequency domain AEM method
The basic theory for AEM can be found in Ward and Hohmann (1988).The FDEM method, as implemented for the Tellus AEM surveys, uses a pair of small coils as the transmitter-receiver (Tx-Rx) pair.The transmitter can be treated as a magnetic dipole that induces eddy currents in the subsurface at discrete frequencies, allowing the resistivity structure to be characterised by comparing the primary and secondary magnetic fields (H p and H s respectively).In the Tellus surveys Tx and Rx are oriented in a vertical, coplanar configuration, i.e. the coil axes are horizontal (Leväniemi et al., 2009), with magnetic dipole moments parallel to the flight direction.AEM data take the form of ratios of the secondary magnetic fields (i.e.formed by current systems in the ground) to the primary magnetic fields (i.e.emitted by the transmitter coil), stated in parts-per-million.Both inphase (i.e.no phase change) and quadrature (i.e. 90 • phase change) fields are considered; generally, the quadrature data are sensitive to the overall ground resistivity, whereas the inphase data are more sensitive to strong conductors.Unlike the magnetotelluric signal, the induced current systems and ensuing secondary magnetic fields of the AEM method are very much in the near-field region.As a result, the volume of Earth interrogated by the AEM signal cannot be easily reduced to a measure of skin-depth-type attenuation.
Although multidimensional modelling and inversion methods are available (see Auken et al., 2014, for a recent review), AEM data are commonly treated as representative of a one-dimensional (1-D) Earth, with spatial smoothing or other constraints along 2-D flight lines possible to improve spatial continuity.One-dimensional modelling is usually performed based on analytical solutions for the layered Earth case, which are well known for most Tx-Rx configurations and can be found in many publications (Keller and Frischknecht, 1966;Ward and Hohmann, 1988;Kaufman et al., 2014).The particular analytical solution in a layered half-space for vertical, coplanar configuration of the transmitter and receiver coils as used here (i.e.parallel, horizontal magnetic dipoles) is found in Minsley (2011).With the forward solution known, an inversion algorithm can be applied to determine a suitable resistivity model.The inversion of AEM data for this work uses the standard damped leastsquares technique of Jupp and Vozoff (1975) as implemented in Airbeo (Raiche et al., 1985).
3 Rathlin Basin survey MT data were collected at 56 sites across part of the onshore Rathlin Basin (site locations shown in Fig. 1) in May and June 2012.Seven parallel profiles were aligned perpendicular to the bounding Tow Valley Fault to the southeast (thick blue line in Fig. 1), with profile and site separations each of 2 km in order to obtain a near-regular array of site locations, facilitating three-dimensional modelling and inversion.Both BBMT (i.e. from ≈ 300 to 0.001 Hz) and audio-magnetotelluric (AMT) (i.e. from ≈ 10 000 to 10 Hz) data were acquired at each site.Data were recorded with Phoenix Geophysics MTU-5A receivers, with either MTC-50 (for BBMT) or AMTC-30 (for AMT) induction coils used to sample the horizontal magnetic field components (H x and H y respectively).Vertical magnetic field components (H z ) were measured at almost all sites using either the appropriate induction coil or an AL-100 airloop, as deemed appropriate given the local ground conditions.The horizontal electric field components E x and E y were sampled by non-polarising lead-lead chloride (Pb-PbCl) PE4 Phoenix Geophysics electrodes arranged at each site in a cross configuration with electrode separations of typically 80 m.BBMT measurements were taken over a period of three nights at each site followed by an overnight measurement of AMT data.
Robust estimates of the frequency domain MT transfer functions were determined from the MT time series using commercial processing software from Phoenix Geophysics that implements the technique described in Jones and Jödicke (1984) and Jones et al. (1989).This approach is based on cascade decimation (Wight and Bostick, 1980) and uses a least trimmed squares algorithm (Rousseeuw, 1984;Rousseeuw and Leroy, 2003) to achieve the robustness of the estimate.Whilst a dedicated distant remote reference site was not used, the principle of remote referencing, as described in Gamble (1979) and Gamble et al. (1979), was applied by using the horizontal magnetic field components of each simultaneously acquired site as reference fields and selecting the best available site as reference (typically five sites were recorded simultaneously).Finally, the AMT and BBMT data at each site were merged into a single response, spanning from 10 000 to 0.001 Hz at most sites.
Although data were acquired at an array of sites with the intent of 3-D inversion, as the data were expected to be predominantly 2-D in nature due to the expected strong lateral contrast across the bounding Tow Valley Fault, multisite, multi-frequency Groom and Bailey distortion analysis was applied to the data on a cross-basinal profile basis using the "strike" analysis tool (McNeice and Jones, 2001).The results of the analysis are presented in Fig. 4 for four depth bands.Analysis of the data with respect to depth in "strike" shows that the data are predominantly 1-D or 2-D to depths of ≈ 2 km, with increasing rms misfits beyond these depths indicative of either a change in strike direction (i.e.still a 2-D structure, but with a different geoelectric strike direction) or www.solid-earth.net/8/637/2017/Solid Earth, 8, 637-660, 2017 The orthogonal vectors at each MT site location indicate the azimuth of best-fitting geoelectric strike direction, coloured by the rms misfit between the observed MT response and the Groom-Bailey model response for the best-fit strike direction.The orthogonal pair of vectors is required as geoelectric strike estimates have a 90 • ambiguity.The size of orthogonal vectors is classified by the phase difference, with larger vectors corresponding to larger phase differences.Small phase differences are associated with 1-D resistivity structure, whereas larger phase differences occur with a 2-and 3-D structure.Larger rms misfits indicate that the decomposition to a 2-D structure is potentially invalid and can be caused by either significant noise contamination and distortion of the data or a 3-D structure.
a fully 3-D structure.An estimate of the regional geoelectric strike azimuth was computed by arithmetically averaging the strike estimates from the deepest depth band (1780-3000 m) in the south-west half of the model, as these estimates represent portions of the model less affected by the coastal effect.The mean strike azimuth of the south-western part of the model is ≈ 43 • E, with a standard deviation of 10 • .The mean geoelectric strike direction is coherent with the strike of the major structural feature, the Tow Valley Fault; however, it should be noted that as geoelectric strike directions inherently possess an ambiguity of ±90 • , the mean azimuth could also be interpreted as 47 • W (i.e. a bearing of 313 • ).As the site profile azimuths were aligned perpendicular (≈ 55 • W) to the mapped strike of the TVF in the area (≈ 35 • E), both the inversion mesh and input data were rotated to a bearing of 315 • (i.e.midway between the mean geoelectric strike direction and dip direction of the TVF).Note that the data were rotated without decomposition, to retain information that does not conform to the 2-D assumption in Groom-Bailey decomposition.For a similar reason, the correction of static shifts was performed before rotation.Note that inversion of nonrotated data and meshes is presented in the Supplement of this article.The MT responses were inverted for three-dimensional models using the ModEM 3-D MT inversion program (Kelbert et al., 2014;Egbert and Kelbert, 2012), with all four impedance tensor elements and vertical transfer functions as input data.MT data were downsampled to a subset of 28 frequencies, spanning from 1000 to 0.001 Hz (displayed in Fig. 10), with an increased sampling of frequencies in the range most sensitive to the target sediment depths (1-0.01Hz).In order to avoid leverage bias in the search for optimum solutions that minimise rms misfit, poor-quality data (typically located near cultural noise sources) were manually identified and removed from the input data.The mesh used for inversion was 59 × 68 × 82 cells in size (X, Y, Z) with the region of interest (the portion covered by MT sites) populated by cells of lateral extent 400 m by 400 m, with layer thicknesses logarithmically increasing beyond the depths required to accurately model bathymetry.Bathymetry was modelled by spanning the first 50 m depth interval with layers of 5 m, increasing to 25 m for the more distant (generally greater than 5000 m from coastline), deeper bathymetry to a total depth of 300 m.Below the bathymetry, layers increased in thickness at a rate of increase of 1.01, increasing to a rate of 1.5 for depths beyond the volume of interest (i.e.beyond 4 km depth) to a total depth exceeding 1500 km (i.e. at least 10 skin depths, given the initial half-space resistivity of 30 m and lowest frequency of 0.001 Hz).The highly efficient but approximate coast-effect forward modelling approach of Booker (described in Burd et al., 2014) was not used as some of our sites are located very close to the coast (well within one skin depth for moderate frequencies), requiring accurate modelling that could not be guaranteed with the approximate approach.
Inversion algorithms determine an appropriate model by iteratively adjusting a resistivity model, computing its forward MT responses, and comparing these responses to the observed data.Whereas the model steps vary depending on the precise algorithm implemented, there are several key parameters that influence an algorithm's behaviour, such as the data errors, type and degree of regularisation, and initial starting and prior models.In particular, the distance between data and model responses, i.e. the sum of the squared residuals is used to measure the distance between the calculated and observed responses.To meet the assumptions of least-squares theory, these residuals must be standardised, i.e. scaled by the variance of the measurements in order to make them normally distributed over N (0, 1).In order to facilitate convergence of the inversion process, an error floor is commonly applied to input data for MT inversions, wherein the error provided for inversion is defined as the greater of either the observed error or some function of the magnitude of the datum.Separate error floors were used for off-diagonal and diagonal impedance tensor components in this study of 5 and 20 % of the mean magnitudes of (Z xy , Z yx ) and (Z xx , Z yy ) respectively.As the diagonal components observed were typically an order of magnitude smaller than the off-diagonal components due to generally 2-D resistivity structures, the signal-to-noise ratio was significantly poorer, and applying a greater error floor here reduces the leveraging of the modelling process by noise-contaminated data.
The regularisation of an inversion describes the weighting between minimising the data residuals and some penalty function, commonly a roughness penalty that enforces smoothness in order to stabilise the resulting model.ModEM allows for the specification of separate regularisation parameters for the x, y, and z directions, with higher values placing greater weight on the penalty function.Several values of these parameters were tested; however, varying regularisation in the z direction had a negligible effect on the model.Laterally isotropic values of 0.15, 0.3, and 0.45 were tested for the x and y directions, and whereas the overall model misfit did not vary significantly, lateral resistivity structure was strongly affected, with a regularisation of 0.15 leading to discrete features with extreme resistivities (i.e.either very high or very low resistivity) beneath MT stations and poor continuity between sites.Resistivity structures in the models obtained with lateral regularisations of 0.3 and 0.45 correlated very well, with a slightly compressed range of resistivities present when regularisation was set to 0.45.In order to reduce over-smoothing of structural boundaries, final inversions were performed with lateral and vertical regularisation parameters of 0.3.
The starting model for each inversion was a 30 m halfspace, with seawater to depths defined by coastal bathymetry fixed (i.e.invariant) at 0.3 m (marine sediments were not included).Final resistivity models were obtained by two consecutive inversion runs: a first model was determined by inversion from a uniform half-space starting model; then a second starting model was constructed by logarithmically averaging the resistivities in this first model with those of the starting half-space.The inversion of the averaged starting model was found to improve model fits significantly and result in resistivity distributions of greater range and contrast across interpreted structural boundaries.This work flow was found to preserve broad structural outlines in order to guide the inversion algorithm whilst suppressing finer features (typically associated with local minima in model space).
Frequency domain airborne electromagnetic data in our MT survey area were acquired as part of the regional Tellus survey described in Beamish (2013).The AEM data were obtained in 2005 and 2006 across Northern Ireland by the AEM-05 system described by Leväniemi et al. (2009), giving observations of both inphase and quadrature data at four frequencies (24 510,11 962,3005,and 912 Hz).The entirety of Northern Ireland (barring high-flight regions above urban areas and steep topography) is covered by the Tellus AEM data set, which comprises parallel flight lines spaced 200 m apart on a bearing of 345 • at a nominal altitude of 56 m, with a spatial sampling rate of one sample approximately every 15 m along the flight lines.The pre-processing work flow of the AEM data is detailed at length in Beamish et al. (2006) and Leväniemi et al. (2009).The inversion of the FDEM data was performed using the one-dimensional Airbeo code from Amira International (Raiche, 1999) that implements the approach of Jupp and Vozoff (1975).

Correction of static shifts
The method proposed here for the correction of staticshift effects follows the approach of Pellerin and Hohmann (1990), adapted to using airborne FDEM data as extrinsic information.Earlier approaches took advantage of the downwards-propagating transient signal of the TDEM method to directly calculate an empirical, quasi-MT 1-D response (Sternberg et al., 1988), and later work by Pellerin and Hohmann (1990) developed this approach by explicitly modelling the TDEM data to obtain a 1-D resistivity model.The MT forward problem can be solved for this resistivity model and factors δ E x , δ E y that correct static shift between the calculated responses and observed data determined by taking the ratio of apparent resistivities.The resulting set of static-shift-corrected MT data may then be modelled as desired.We propose that the approach of Pellerin and Hohmann (1990) can be equally applied with FDEM data, subject to the same constraints, namely that there should be an overlap between the minimum depth of penetration of the MT sounding and the maximum depth of penetration of the extrinsic (to MT) information and that the dimensionality of the two methods should agree.The first constraint mandates the use of high-frequency AMT data rather than BBMT data to ensure overlapping volumes of sensitivity.A flowchart describing the steps taken in correcting the Rathlin Basin MT data with the Tellus AEM data is presented in Fig. 5.
The approach used here is predicated upon certain key assumptions about the near-surface geology and the induced galvanic distortion, and these assumptions clearly show the limits of the approach.Firstly, as mentioned we assume that the near-surface geology is 1-D in structure, i.e. we can treat the off-diagonal impedance tensor elements Z xy and Z yx independently without rotating or otherwise preparing the data.Secondly, we assume that the galvanic distortion affects only the electric field; as the total electric field is represented by E x and E y , we require only two corrective factors δ E x and δ E y to be applied to the two impedance tensor element pairs corresponding to E x and E y , namely (Z xx , Z xy ) and (Z yx , Z yy ).Finally, whereas the cause of static-shifttype distortion of the MT data is the estimation of the electric field, the airborne FDEM data are observations of the magnetic field (directly proportional to the electric field).Hence, the FDEM data are unaffected by static-shift-type distortion, and any resistivity estimates computed from these data are closer to the true values.For clarity about the effect of the static-shift corrective factors, δ E x and δ E y are often presented here on a logarithmic scale; as the corrective factors transform to additive (or subtractive) changes of resistivity on a logarithm scale, their value is intuitively related to whether a correction to more resistive or more con-ductive true data is required.For example, a (multiplicative) static-shift corrective factor of δ E = 3.16, indicating that the true resistivity ρ T = δ E ρ obs = 3.16ρ obs corresponds in logarithmic scale to an additive static-shift correction, i.e. log 10 (ρ T ) = log 10 (ρ obs ) + 0.5.
The first step of the implemented procedure is the inversion of the uncorrected MT data to obtain a baseline resistivity model for comparison M o .The remaining steps describe how static-shift-corrected MT data Z c were obtained by solving the MT forward problem for a model of the AEM data, and in turn inverted to obtain a corrected resistivity model M c .
-Step 1: 3-D MT inversion of observed MT data Z o with the ModEM code to obtain model M o .Both the inversion mesh and input data were rotated to a bearing of 315 • for inversion.
-Step 2: Modelling of each four-frequency AEM sounding within the survey area as a single-layer structure (i.e.half-space) with Airbeo (Raiche, 1999), resulting in an apparent resistivity value at each location best reproducing the observed AEM data.
-Step 3: Interpolation of AEM half-space models by inverse-distance-weighted (IDW) averaging of logresistivity values to populate the uppermost 200 m of an MT forward modelling mesh with cells of 170 × 170 m.Below the uppermost 200 m, the model reverts to a 100 m half-space.
-Step 4: 3-D solution of MT forward problem for the resistivity model found in Step 3 with ModEM, resulting in a set of high-frequency synthetic MT responses at six frequencies from 10 000 to 1000 Hz for each MT sounding location.The frequencies chosen coincide with those of the downsampled MT data.
-Step 5: Multiplicative static-shift corrective factors δ E x and δ E y found by taking the ratio of the apparent resistivities of the up to six high-frequency MT responses found in Step 4 to those of the observed data (i.e.δ E x = ρ a,xy (synthetic)/ρ a,xy (observed)) over the 10 000-1000 Hz band.Due to either noise contamination or non-1-D behaviour (i.e.violating our assumptions), not all data in the compared frequency band were used; typically, comparison was made using three to four of the six data points that parallel the synthetic responses.The corrective factors are applied to the entire bandwidth of the unrotated observed data to obtain ρ a (corrected), with δ E x applied to all elements dependent on E x , i.e. ρ a,xx and ρ a,xy .δ E y is treated analogously.
-Step 6: The static-shift-corrected data are used as input for 3-D MT inversion with ModEM to obtain an improved resistivity model M c .Both the inversion mesh

Tellus AEM data
Step 2

1-D AEM inversion
Ensemble of halfspace AEM resistivity models Step 3 IDW interpolation and input data were rotated to a bearing of 315 • for inversion.
Whereas multi-layered models that better reproduce the AEM data can also be determined using Airbeo, they were not used in favour of the half-space apparent resistivities for two principal reasons.Firstly, the interpolation of the apparent resistivities to a 3-D MT mesh can be directly computed, whereas multi-layered models require more advanced approaches to reconcile variation in layer thicknesses unless these are explicitly set in the 1-D inversion to facilitate interpolation.Secondly, the depth of sensitivity of the lowest frequency of the AEM data (912 Hz ≈ 60 m, from forward modelling) has a moderately narrow overlap with the skin depth of the highest MT frequencies (10 800 Hz, skin depth ≈ 50 m for a typical near-surface resistivity of ≈ 100 m).Above this overlapping volume, the MT data remain sensitive to but poorly resolve resistivity contrasts, and for the purposes of our work, the added complexity in multi-layered AEM modelling does not significantly improve our results.A map of the half-space models found from the AEM data in Step 2 is shown in Fig. 6.As this work is reliant upon the depth of sensitivity of the FDEM data, we assume that the causes of static shift encountered are locally one-or twodimensional anomalies that perturb estimates of the electric field estimates E i , causing frequency-independent static shifts of the MT impedance data.Additionally, as each electric field component is used to compute two forms of data (i.e.{Z ii , Z ij } ∝ E i ), the same correction will be applied to data sharing a common electric field component.Groom and Bailey distortion analysis of the high-frequency MT data showed that the regional near-surface resistivity structure was 1-D or 2-D at most, validating this approach to correct-ing the static-shift effects.We reiterate that for regions with regionally 3-D structure, this approach would be invalid.
The statistical and spatial distributions of the static-shift multiplicative corrective factors δ E x and δ E y were examined for spatially coherent correlation and features that may be indicative of known regional-scale geological structures.The spatial distribution and magnitudes of the corrective factors across the survey area are shown in Fig. 7. Whereas some spatial correlation between static shifts and geology is apparent, static-shift variation mostly does not coincide with mapped surface geological boundaries, which invalidates the use of regional geological units as predictors of static shift.The histograms in Fig. 8 show the distributions of the logarithmic transforms of δ E x and δ E y , with statistical quantities shown in Table 1 (the quantities of both log δ En and δ En are tabulated, as are those of the mean static-shift correction at each site δ E ).The distributions of both δ E x and δ E y appear close to log-normal, with longer tails towards conductive corrections that likely indicate natural bias introduced by the sampling, such as the geographic location of the MT sites.The bivariate distribution of log δ E x and log δ E y from each site is shown in Fig. 9, where the strong 45 • (i.e.δ E x ≈ δ E y ) clustering indicates that the two electric field components at most sites are similarly affected by static-shift-type effects, with no evident regional preference for static shifts to one polarisation over the other.
From examination of the spatial and statistical distributions of δ E x and δ E y , the galvanic distortions present in these data show no consistent anisotropic behaviour and weak spatial correlation with surface geology.As the accepted theory is that galvanic distortion is typically caused by irresolvable near-surface resistivity inhomogeneities below the level of www.solid-earth.net/8/637/2017/Solid Earth, 8, 637-660, 2017  We propose that this conductor primarily represents the conductive MMG, with a clearly delineated upper boundary.The interpretation of the lower boundary against the SSG (and equally, the boundary between the SSG and the EG sandstones) is hindered by both the smoothing effects of the regularised inversion approach used and the fact that inductive EM responses are intrinsically sensitive to the tops of conductive units (and their integrated conductivity) rather than to the bottoms of conductive units (i.e. the tops of resistive units).The TVF, which forms the south-eastern boundary of the basin, is clearly defined, although the angles of dip  (blue) and δ E y (yellow) in decades (i.e. as additive factors, where +0.5 corresponds to a multiplicative factor of 10 0.5 , ≈ 3).Statistical quantities related to these distributions are tabulated in Table 1.sponds to a multiplicative factor of 10 0.5 , ≈ 3), paired by site.The strong clustering and 45 • trend (i.e.δ E x ≈ δ E y ) indicate that the orthogonal static shifts are generally similar sized across the survey area, with no observable regional trend of greater static shift of one polarisation.
ity distribution of the basin portion of M c as a set of layer-bylayer histograms, each histogram normalised by the number of counts for the mode of that layer.As the model becomes much smoother below depths of approximately 2500 m, several tests were carried out to examine the sensitivity of the model responses to resistivity contrasts beyond these depths.
The insertion of synthetic resistive bodies showed that the model responses are sensitive to them but unable to resolve resistivity contrasts below the conductive basin.It remains possible that older sediments with higher resistivity exist below this depth; however, the resistivity contrasts present within the MT data are insufficient to clarify their existence or extent.Two diagnostic measures were used to assess the changes between M c and the original model M o , caused by static-shift correction.The first diagnostic measure is the logarithmic resistivity difference between the two models, and the second is the normalised cross-gradient (NCG) of the two models, being the cross-product of each model's gradient vectors at each cell, The normalised cross-gradient was introduced in Gallardo and Meju (2003) in the context of joint inversion based on structural similarity and has been used by, e.g., Schnaidt and Heinson (2015) and Rosenkjaer et al. (2015) to highlight structural similarities and differences between resistivity models.
As the logarithmic resistivity difference between M c and M o highlights discrepancies in the absolute resistivity values of the two models, it is more useful for observing contrasts in formation resistivities (i.e.locations where one model has a relative minimum or maximum) rather than comparing structural boundary locations.The figures presented show that the conductive features in M c tend to be between 0.25 and 0.50 of a decade more conductive than their equivalent volumes in M o , whereas the resistive features show similar tendencies towards greater resistivity in M c .The largest vertical variations in (i.e. from a local maximum to a local minimum or vice versa) are typically seen at the upper boundary of the MMG, making the interpretation of this boundary far less subjective in model M c .The resistivity difference also shows greater definition of the thin conductors interpreted as the LLG, which are difficult to delineate on the plots of resistivity alone due to weaker resistivity contrasts.The metasedimentary horst shows a range of both conductive and resistive contrasts; however, as this region shows greater geoelectrical heterogeneity and greater static-shift effects in the observed data, the cause of these differences is variations in the uncorrected data used to model M o rather than the relatively uniform horst seen in M c .The spatial distribution of is insightful, particularly the distribution of points where = 0 (i.e. points where M o and M c have the same resistivity) which generally do not coincide with the resistivity minima or maxima of M c nor with the identifiable resistivity boundaries.Instead, these points in places indicate regions in which the structural geometries differ between models, such as thicker or thinner conductive bodies (e.g. at ≈ 800 m depth, 10 000 m along Profile A, Fig. 12b, where the = 0 line occurs in the middle of a conductor).These differences in geometry are most noticeable as horizontal features in the uppermost 1000 m, where the models are constrained locally by single sites and the static-shift correction affects the model in a generally 1-D manner.The lateral distribution of at these shallow depths is a result of the pseudo-1-D local scaling of the model differing from site to site; as such, we do not expect the distribution of to replicate those of the static-shift corrections themselves (except for the trivial case of surficial and immediately subsurficial layers).
The NCG is largest where the gradient vectors (i.e.resistivity changes) of the two models are orthogonal, such as differences in structures and locations or similar structure but differing magnitude of resistivity change (i.e.difference in the curvature of the model).Both situations of elevated NCG are evident in the models, especially on the zero-lines (i.e.where the resistivity models intersect and reverses polarity, requiring significantly different gradient vectors), and at resistivity minima and maxima of M c (and assumingly M o , although not shown), indicating a difference in curvature due to slight offsets of critical resistivity points or large resistivity differences in co-located critical points.We note again that the significant variation in NCG values computed in the Dalradian horst likely reflects the heterogeneity of this block.Whereas the NCG was originally proposed as a means of determining structural similarity between models, our use of the NCG in this study suggests that it may also be used to  Each layer represents the distribution of resistivities present in that layer of the model, normalised by the number of counts for the mode of that layer (black).Whereas a wide spread of resistivities is present from 300 to 2500 m depth, the column of modal resistivities can be viewed as a "typical" 1-D structure of the basin, with a distinct minimum of ≈ 10 m evident at 1100 m depth, the depth interval of the MMG observed in the Port More 1 borehole.
assess the coincidence of similar structure with different resistivities.
The closeness of MT model responses to the data is commonly judged by the normalised root mean square error (nRMS), defined as where r is the residual between calculated and observed data, normalised by the square root of its variance.In the case of a sufficient number of data N with approximately Gaussian and independent misfits, nRMS error should approach 1 for a model fitting to within 1 standard error; however, the application of an error floor to the observed misfits results in artificially lowered nRMS errors.Additionally, due to the global averaging of residuals, an nRMS error is uninformative about which portion of data are poorly fit.Regardless, the nRMS error remains a useful metric for comparing the relative goodness of fit of a succession of model responses to the same data set.The two models presented here reproduce the observed and corrected MT responses to similar degrees with overall nRMS of 1.77 and 1.96 for M c , the corrected model (66 iterations total), and M o , the original model (52 iterations total), respectively.1) and inversion frequency number (as shown in Fig. 10) for each element of the impedance tensor (a-d).The normalised residual is determined by taking the absolute value of the difference between model response and datum and normalised by the error used in inversion (i.e. the greater of the observed errors or the applied error floor).Normalised residuals of less than 1 indicate a residual smaller than the error.Note that misfits of apparent resistivity use the logarithmic difference to determine the residual.Values in grey represent the lower-quality data that were masked and not included in the inversion.
be drawn from a single nRMS value for the entire model.The measure of misfit displayed is the magnitude of the residual (the absolute value of the difference between the model response and observed datum), normalised by the error used for inversion (the larger of the experimental errors or the applied error floor), similar to the overall nRMS defined in Eq. ( 8), with values of 1 or less indicating a residual smaller than the error (i.e.fit to within 1σ ).In general, the off-diagonal components are well fit, with a handful of sites showing poorer fits with respect to phase components.
As the model responses are noise-free and have internally consistent magnitude and phase behaviours, worse misfits to impedance phases in comparison to magnitudes imply that the affected portions of the data themselves may not entirely satisfy the assumption of being in the far-field region of the signal, likely due to the proximity of a noise source.Interference from a local noise source contaminates both the magnitude and the phase of the data, and inverting such data can lead to a spurious model structure.Data that showed significant noise contamination (identified by phases trending towards 0 • ) were removed prior to inversion.The misfits of the diagonal impedance components show similar behaviour in that the magnitude is better fit than the phases; however, the differences in fit between magnitudes and phases are greater for the diagonals.Some increase in misfit is expected at the high-frequency limits due to the relocation of MT sites to the centre of their respective cells for computation; as the cells are 400 m wide within the survey area, the model is effec-tively 1-D for forward computation of MT responses at such high frequencies, and diagonal impedances are 0 for such resistivity structure.As the observed diagonal MT data can be non-zero due to resistivity variations at scale lengths below those used for modelling purposes, the fit of these data can be improved by the use of a finer mesh, with a corresponding increase in computation time.At longer periods the misfits of diagonal phases are increased for the same reasons as the off-diagonal phases, namely noise contamination in the data; however, due to the generally smaller magnitudes of the diagonal MT data the signal-to-noise ratio can be considerably worse.It is possible that the error floor used for the diagonal components is perhaps too stringent, and lower misfits would be obtained with a higher error floor; however, as the observed data imply a predominantly 2-D structure (as shown in Fig. 4), such a change is unlikely to drastically change the final model.Although examining the misfits of a single model's responses to the data set used in its determination in such a granular fashion is useful in determining which data components are poorly fit and possibly why, taking a similar approach in order to compare two models is not valid in this case.Due to the application of individual static-shift corrections to each site's MT data, the gradients of the respective data spaces of the models are significantly altered.As most inversion algorithms (the nonlinear conjugate gradient method implemented in ModEM included) rely upon gradients within the data space to determine the direction of line searches as part of their optimisation, the static-shift corrections applied all but guarantee that the two inversions presented in this work are the products of different paths through their data spaces.Hence, although the overall mean nRMS estimates are similar, we cannot categorically state that any differences in misfit at the granular, individual datum level are not simply due to the different gradient progressions.
The key test of the two models lies in comparing the measured resistivity values from the PM1 borehole to the vertical resistivity columns from M c and M o corresponding to the location of the borehole (Fig. 17).The input data for the two inversions differ solely in the application of static-shift corrections, assumed to only affect the magnitude of the data; considered as a purely 1-D problem, such a shift results in the rescaling of both layer resistivities and thickness.The nRMS misfits of the two models' responses for this location are 0.97 (M c ) and 1.05 (M o ).The effect in 3-D is more complicated, as the variation occurs in both lateral and vertical directions; however, some geometrical correlation and similar-shaped structures can still be expected in each of the two models.Note that we argue this based upon the perturbation of the data magnitudes rather than the similar nRMS -as mentioned, such an assumption based on misfit is unreasonable.Additionally, for the uppermost extents of the model (i.e. at depths within the inductive volume of only a single MT site), it can be seen that the application of staticshift corrections results in similar effects to a 1-D case, i.e. a model from data that have been statically shifted upwards to higher apparent resistivities will have increased resistivitythickness products (thicker, more resistive layers).With this behaviour in mind, the two models M c and M o have very similar geometries at shallow depths (the uppermost 200 m), with the resistivities of M c significantly closer to the mean of the borehole measured resistivities.Due to the regularisation of the inversion process, neither model can adequately reproduce the highly variable near-surface variations measured through the ALG in the borehole.Comparison of the structures at deeper depths shows that whereas M c has conductive layers at 300 and 1300 m depth, with resistivities of ≈ 3 and 15 m, interpreted as the middle of the LLG sediments and the combined conductances of the MMG and SSG respectively, the equivalent conductors (assuming a scaled geometry) in M o occur at 450 and 1500 m, with resistivities of ≈ 20 and ≈ 8 m.Similarly, the models show a resistor of ≈ 60 and ≈ 30 m at 600 and 750 m depth in M c and M o respectively.
The resistivity columns from the two models at the PM1 borehole site can also be compared on the basis of integrated conductances, i.e. the ratio of a layer's thickness and resistivity.In such a comparison, the LLG sediments are represented as conductances of 31 (M c ) and 13 (M o ) S, the interpreted dolerite sill sequence as conductances of 34 (M c ) and 23 (M o ) S, and the MMG or SSG sediments as conductances of 103 (M c ) and 100 (M o ) S. The greatest contrast in conductance is observed in the LLG, wherein M c shows a conductance almost half an order of magnitude greater than M o .The resistivity measurements in the PM1 borehole do not span the entire LLG interval; however, given the measured values the real conductance of the interval is likely greater still than the 31 S recovered in M c .Hence, although M c does not fully recover such a conductance, the elevated conductance in comparison to M o , coupled with the near-identical lower conductances, indicates that M c is closer to the real structure at this location.
Given the knowledge of the lithology and the measured borehole resistivities, we conclude that M c is a categorically superior resistivity model in that, of the two models, it more accurately and correctly resolves the central depths and resistivity values of the lithological units encountered.Critically, the limited borehole resistivity data in the MMG and SSG suggest that resistivity increases slightly with depth between the two units, and whereas the resistivity column of M o more closely matches the absolute resistivity values, the column of M c better approximates the trend of the observed borehole resistivities.Without further external information to verify the models with, the smaller differences between the M o model resistivity column and the borehole resistivity measurements of the target sediments cannot lead to the acceptance of the entire model, as the shallower structure is poorly recovered.Conversely, M c , the static-shift-corrected model, correlates well with the borehole lithological boundaries from the PM1 borehole and shows a trend in resistivity www.solid-earth.net/8/637/2017/Solid Earth, 8, 637-660, 2017 with depth that matches the trend evident in the wireline resistivity measurements.
Changes in resistivity at depth from static-shift correction have strong implications for the interpretation of the reservoir potential of the area.Assuming Archie's Law holds and depending on the cementation exponent m (formation dependent, but typically ≈ 2 for reservoir formations), a relative change in resistivity f ρ between two models ρ c and ρ o (i.e.ρ c = f ρ ρ o ) corresponds to a change in estimated porosity of φ c = f φ φ o , where f φ ∝ m f ρ .The relationship between porosity and hydraulic permeability is known to be complex and highly nonlinear; it would be highly speculative to quantify how the resistivity perturbations from the correction of static-shift effects would affect any attempt to estimate hydraulic permeability.At this stage, given the paucity of borehole information and unfavourable geoelectric setting (in which the conductive MMG directly overlies the moderately less conductive SSG and EG reservoir targets), the modelled resistivity distribution cannot be used to extend the knowledge of the hydraulic properties of the reservoir targets beyond what is reported from the boreholes within the Rathlin Basin with any useful level of confidence.
The quasi-1-D correction of static-shift effects applied here affects the resulting three-dimensional models to significantly greater depths than expected, with differences between the corrected and uncorrected resistivity models of up to half an order of magnitude present at depths of up to 2000 m (see, for example, the constant depth slice at 2100 m depth in Fig. 11) that do not correlate with the spatial distribution of applied static-shift corrections (Fig. 7).Whereas the structural geometries are very similar between the two models (especially in the shallow regions constrained by single MT sites), interpretations of the target sediment locations in the models differ in significant ways, with both depth to sediments and absolute resistivity affected.As such, it is clear that when seeking accurate resistivity and depth estimates in three-dimensional modelling of MT data, galvanic distortion must be accounted for as its effects are subtle but pervasive.

Conclusions
An approach for the correction of static-shift-type galvanic distortion in MT data utilising airborne FDEM data has been tested that follows the use of TDEM data in previous methods.The new approach was tested on an MT data set from Northern Ireland, using a publicly available regional data set of airborne frequency domain electromagnetic data to create a set of corrected MT data.Three-dimensional inversion of each magnetotelluric data set recovers structures with similar geometries; however, structures in the near-surface show scaling of resistivity-thickness products proportional to the static-shift correction applied.When compared to geophysical borehole logs it is clear that the model from static-shiftcorrected data reproduces the observed resistivity with significantly greater fidelity.Significant suborders of magnitude variations in resistivity are caused in the model by the correction of the static-shift components of galvanic distortion, not only in the near-surface but extending down to the target sediment depths (≈1500 m).As the test area is a prospective geothermal exploitation site (a sedimentary aquifer with elevated temperatures), electrical resistivity could be used to infer the heterogeneous distribution of hydraulic properties within the reservoir; however, the suborder of magnitude perturbations caused by inadequate consideration of galvanic distortion would lead to a gross misestimation of these physical properties.
The model determined by the inversion of static-shiftcorrected data was found to better recover the resistivity structure observed in a nearby borehole in comparison to the model from observed data.Based on these observations, we conclude that airborne FDEM data provide sufficiently accurate resistivity estimates to allow the correction of staticshift effects in MT data.Note that this approach as discussed here is valid only for locales where the near-surface resistivity distribution is approximately one-dimensional.Given the often regional acquisition and open availability of such AEM data, it is hoped that the approach demonstrated here could be further tested with other MT surveys.Pending further case studies, FDEM could in future be considered as another alternative method to evaluate and correct static-shifttype distortion.Additionally, whereas our approach assumes one-dimensional, single-layer models for the AEM data in deriving the static-shift corrections, future advances could investigate what effect more advanced AEM modelling (i.e.multiple layers or, where applicable due to AEM acquisition specifications, full 3-D modelling) would have on the computed forward MT responses and associated static-shift corrections.

Figure 1 .
Figure 1.Geological map of Rathlin Basin area, overlain with magnetotelluric acquisition sites numbered and denoted by black crosses.Note the large areal extent of the Antrim Lava Group basalts, with minimal surface expression of the underlying sedimentary basin.Profiles A, B, and C denote the locations corresponding to Figs. 12-14.Inset map shows the northern half of Ireland, with the location of the magnetotelluric survey area shown by the red rectangle.The yellow line indicates the location of the density profile shown in Fig. 3.

Figure 3 .
Figure3.Density model fromMitchell (2004), showing the broad geological structure predicted prior to this study.Density profile locality is shown in Fig.1(inset).Although the profile crosses the Rathlin Basin south-west of the MT survey area, a similar lithological sequence is expected elsewhere in the basin, including the survey area.

Figure 4 .
Figure 4. Visualisation of dimensionality of MT data, decomposed using the "strike" analysis program McNeice and Jones (2001) over a range of depth bands (in metres below sea level).The orthogonal vectors at each MT site location indicate the azimuth of best-fitting geoelectric strike direction, coloured by the rms misfit between the observed MT response and the Groom-Bailey model response for the best-fit strike direction.The orthogonal pair of vectors is required as geoelectric strike estimates have a 90 • ambiguity.The size of orthogonal vectors is classified by the phase difference, with larger vectors corresponding to larger phase differences.Small phase differences are associated with 1-D resistivity structure, whereas larger phase differences occur with a 2-and 3-D structure.Larger rms misfits indicate that the decomposition to a 2-D structure is potentially invalid and can be caused by either significant noise contamination and distortion of the data or a 3-D structure.

Figure 5 .
Figure 5. Flowchart illustrating the steps required to (a) find an original, uncorrected resistivity model M o from MT data (Step 1 only) and (b) find an improved, static-shift-corrected resistivity model M c from MT data, with coincident airborne FDEM data (data from the Tellus project used for this study).For clarification, the MT data used in Step 5 are in their original coordinate system, whereas the MT data in Steps 1 and 6 are rotated to a bearing of 315 • for inversion.IDW in Step 3 refers to inverse distance weighting.

Figure 6 .
Figure6.Map of airborne FDEM half-space resistivity models, overlain with geological boundaries (solid black lines, as in Fig.1) and MT site locations (black cross symbols).The model exhibits good correlation with the mapped boundaries, as well as variation within each geological unit.Note that data from one flight line are missing, resulting in a pale dashed line in the figure.These absences were interpolated across in the synthetic model.

Table 1 .
Statistical measures of the log 10 of the corrective factors δ E x and δ E y (i.e.additive in logarithmic domain), as displayed in Fig.8.In addition to the means, medians, and standard deviations of log 10 δ E , the means and medians of δ E (i.e.multiplicative) are also shown.Note that δ E is the site-wise mean static shift -i.e. the mean of δ E x and δ E y on a site-by-site basis.AMT data, the weak spatial correlation with surficial geology is unlikely to be random and instead reflects the variation in heterogeneity of formations4 Model evaluation and discussionInitial comparison of the two models, one with (M c ) and one without (M o ) static-shift correction, shows greater resistivity contrasts and sharper delineation between resistive and conductive volumes in the corrected model M c , and for these reasons M c is preferred for interpretation and evaluation of structures in the Rathlin Basin.Due to the static-shift correction, M c also more likely better represents the correct resistivity distribution of the real Earth.Figures 11-14 present the resistivity distribution ρ of horizontal and vertical slices taken through the preferred model M c (profile locations are marked in Fig.1).Measures of model comparison are also shown on these figures and are discussed below (note that no images of the resistivity distribution of M o are presented, as the resistivity differences with respect to M c are subtle and difficult to perceive given the dynamic range of resistivities in the visualisation).The two models generally show good correlation of geoelectric structural geometries but differ in the absolute resistivity estimates and the exact extent of these structures.The resistivity structure of M c is generally quite simple, with the most prominent feature being the extensive conductor, of a resistivity less than 10 m, that lies between 900 and 2000 m depth and extends from the Tow Valley Fault in the south-east for a distance of up to 10 km towards the northwest.

Figure 7 .Figure 8 .
Figure 7. Spatial distribution of static correction factors δ E x (a) and δ E y (b) in decades (i.e. as additive factors, where +0.5 corresponds to a multiplicative factor of 10 0.5 , ≈ 3).Blue shading indicates resistive corrections to the data, whereas red indicates conductive corrections.Note that for the majority of sites, δ E x and δ E y are similar in magnitude and polarity.The magnitudes of static corrections are generally less than 0.5 decades in size, corresponding to apparent resistivity estimates being shifted by factors of less than 3.

Figure 9 .
Figure 9. Bivariate distribution of calculated static-shift corrections δ E x and δ E y in decades (i.e. as additive factors, where +0.5 corre-

Figure 11 .
Figure 11.Top row shows resistivity slices through the static-shift-corrected model M c taken at 850 (a), 1550 (d), and 2100 m (g) below sea level.Middle row shows the resistivity difference between the models in decades ( (M c , M o ) = log 10 (M c /M o )) for the same depths, where red shows M c more conductive than M o and blue more resistive.Bottom row shows the magnitude of the normalised cross-gradient (the cross-product of the gradient vectors of each model ∇M, as defined in Eq. 7), as a diagnostic of structural similarity between the models, with 0 (blue) showing parallel gradient vectors (i.e.very similar structure) and 1 (red) showing orthogonal gradient vectors and structural disagreement.Note that the difference and cross-gradient plots are overlain by the 10 m contour from the corresponding resistivity slice.Magenta lines in panel (g) indicate the location of Profiles A, B, and C.

Figure 12 .Figure 13 .Figure 14 .Figure 15 .
Figure 12.Profile A taken along the axis of the concealed basin through the static-shift-corrected resistivity model M c (location shown in Figs. 1 and 11g).The resistivity is shown in panel (a), the resistivity difference (M c , M o ) is shown in panel (b), and the cross-gradient of M c and M o is shown in panel (c).Contours on the difference and cross-gradient plots show the 10 m contour.For presentation a vertical exaggeration of 1.5 is used.

Figure 16 .
Figure16.Visualisation of the normalised residuals of the magnetotelluric responses of M c compared to observed data at each site (numbered as per Fig.1) and inversion frequency number (as shown in Fig.10) for each element of the impedance tensor (a-d).The normalised residual is determined by taking the absolute value of the difference between model response and datum and normalised by the error used in inversion (i.e. the greater of the observed errors or the applied error floor).Normalised residuals of less than 1 indicate a residual smaller than the error.Note that misfits of apparent resistivity use the logarithmic difference to determine the residual.Values in grey represent the lower-quality data that were masked and not included in the inversion.

Figure 17 .
Figure17.Left-hand panel shows the resistivity columns from both M c (blue) and M o (orange) adjacent to the Port More 1 borehole, plotted with the observed normal resistivity data.Due to the regularisation in the inversion process, neither model can reproduce the high variability of resistivities within the ALG observed in the uppermost 100 m of normal resistivity data.Note that the resistivity column of M c approaches the observed low resistivities of the LLG and reproduces the trend observed in the MMG and SSG resistivities, albeit with resistivities approximately half a decade greater.In contrast, the resistivity column of M o correlates poorly with the observed normal resistivity data, with overestimated resistivities and layer thicknesses.Right-hand panel shows interpreted structure of each resistivity column (patterned, on either side), with the observed borehole lithology in the centre column.The interpreted M c structure shows significantly improved correlation with the observed lithology; in particular, the upper boundary of the MMG against the dolerite sills is recovered much closer to the true depth in M c than in M o .The static-shift correction factors δ E x and δ E y applied at this site were 0.40 and 0.38 (−0.40 and −0.42 in decades) respectively.