# Three-dimensional geoelectrical imaging beyond 1 km depth for mineral exploration: framework of deep electrical resistivity tomography and induced polarization with advanced strategies

Bitnarae Kim, Jacques Deparis, François Bretaudeau, Simon Védrine, Jochen

Kamm, Uula Autio, Cedric Patzer, Suvi Heinonen, Julien Gance, Mathieu

Darnet

► To cite this version:

Bitnarae Kim, Jacques Deparis, François Bretaudeau, Simon Védrine, Jochen Kamm, et al.  Three-dimensional geoelectrical imaging beyond 1 km depth for mineral exploration: framework of deep electrical resistivity tomography and induced polarization with advanced strategies. Geophysical Journal International, 2026, 244 (3), 10.1093/gji/ggaf460. hal-05481366

HAL Id: hal-05481366

https://brgm.hal.science/hal-05481366v1

Submitted on 28 Jan 2026

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
Geophys. J. Int. (2026) **244**, 1–22
Advance Access publication 2025 November 17
Research Paper

https://doi.org/10.1093/gji/ggaf460

# Three-dimensional geoelectrical imaging beyond 1 km depth for mineral exploration: framework of deep electrical resistivity tomography and induced polarization with advanced strategies

Bitnarae Kim<sup>1</sup>, Jacques Deparis,<sup>1</sup> François Bretadeau<sup>1</sup>,<sup>1</sup> Simon Vedrine,<sup>1</sup> Jochen Kamm<sup>1</sup>,<sup>2</sup> Uula Autio,<sup>3</sup> Cedric Patzer,<sup>2</sup> Suvi Heinonen<sup>1</sup>,<sup>4</sup> Julien Gance,<sup>5</sup> Mathieu Darnet<sup>1</sup> and SEEMSDEEP working group

Accepted 2025 November 5. Received 2025 October 27; in original form 2025 May 7

## SUMMARY

Direct current electrical resistivity tomography (ERT) is a widely used geophysical method for near-surface investigations, offering high-resolution imaging for geological, engineering and environmental applications. While traditional ERT surveys typically target depths of 0–200 m, technological advancements have enabled deeper investigations, commonly referred to as Deep ERT. In this study, we explore the practical challenges and methodological improvements associated with Deep ERT, particularly when combined with induced polarization (IP) measurements. Rather than other electromagnetic methods, ERT offers a more straightforward framework for analysing IP effects, which can potentially correlate with the volume fraction of ore minerals. Nevertheless, deep IP investigations are often challenged by weak signal strength and various sources of electromagnetic interference. To address these challenges, we evaluated key strategies including survey planning, high-power current injection, unconventional electrode configurations and advanced signal processing techniques. The adoption of nodal geophysical recording systems eliminates the logistical constraints of cabled multi-electrode setups, improving flexibility and data acquisition efficiency. Additionally, continuous full time-series recording allows for enhanced noise filtering and signal stacking, ultimately increasing the signal-to-noise ratio and extending the effective exploration depth. We demonstrate this methodology through a comprehensive case study conducted at the Koillismaa Linear Intrusion Complex in Finland, where a 3-D Deep ERT-IP survey successfully delineated conductive and chargeable anomalies at depths exceeding 1.5 km. These anomalies closely align with independent gravity and borehole logging data, consistent with the mafic-ultramafic intrusion structures. Our results emphasize the importance of balancing data quality, survey efficiency and spatial resolution in survey design. This work not only provides a robust workflow for the implementation of Deep ERT-IP surveys but also represents the first documented successful acquisition of high-quality IP data at these substantial depths, significantly advancing the state of deep geoelectrical exploration.

**Key words:** Electrical properties; Europe; Electrical resistivity tomography (ERT).

## 1. INTRODUCTION

Direct current electrical resistivity tomography (ERT) is one of the most well-established geoelectrical methods for near-surface investigations, with a wide spectrum of applications in geological, engineering and environmental sciences. While magnetotellurics (MT),

audio magnetotellurics (AMT) and controlled source magnetotellurics (CSAMT) are the standard methods to image resistivity in 2-D or 3-D at depths larger than a few hundred metres, and even up to several kilometres (Vozoff 1990; Spichak & Manzella 2009; Balasco et al. 2015; Darnet et al. 2020b), airborne electromagnetics (AEM) is usually the method of choice to map efficiently very

<footnote n="1"><sup>1</sup>BRGM, 3 Av. Claude Guillemin, 45100 Orléans, France. E-mail: b.kim@brgm.fr</footnote>

<footnote n="2"><sup>2</sup>Geological Survey of Finland (GTK), Vuorimiehentie 5, 02151 Espoo, Finland</footnote>

<footnote n="3"><sup>3</sup>Geological Survey of Finland (GTK), Walleniuksenväylä 10, 96100 Rovaniemi, Finland</footnote>

<footnote n="4"><sup>4</sup>Institute of seismology, University of Helsinki, Pietari Kalminkatu 5, 00014 University of Helsinki, Finland</footnote>

<footnote n="5"><sup>5</sup>IRIS Instruments, 1 Av. Buffon, 45100 Orléans, France</footnote>
large areas over the first 500 m (Wynn 2002; Darnet et al. 2020a; Reninger et al. 2020).

On the other hand, ERT is usually the ideal choice for achieving high-resolution imaging in near surface application while maintaining the best sensitivity through the use of multiple efficient galvanic sources and electric field measurements (Asch & Morrison 1989; Wilkinson et al. 2010; Everett 2013). Conventional ERT surveys typically rely on cabled multi-electrode measurement systems. While there are no fundamental limitations preventing the application of ERT at greater depths, increasing logistical challenges and costs impose practical constraints. Although most ERT applications target depths in the range of 0–200 m, several studies have demonstrated its potential for investigating deeper geological structures (Alfano 1984; Storz et al. 2000; Tamburriello et al. 2008; Pucci et al. 2016; Troiano et al. 2019; Balasco et al. 2022, 2015; Rizzo et al. 2022).

Driven by technological advancements, deeper ERT applications have gained increased attention, leading to improvements in efficiency and exploration capabilities. However, the term ‘Deep ERT’ has often been used inconsistently to describe surveys with varying depths ranges. In this study, we adopt the definition proposed by Balasco et al. (2022), which defines Deep ERT investigations as those achieving exploration depths of at least 500 m.

One of the key advantages of ERT compared to other resistivity methods like MT or time-domain EM, is its ability to effectively capture induced polarization (IP) properties from the same data set (Kemna et al. 2012; Revil et al. 2022). While electrical resistivity characterizes a rock’s ability to impede the flow of electric charges, which is influenced by factors such as mineral composition, fluid content and porosity, IP refers to the delayed response of a material in storing and releasing electrical charge. This phenomenon is often attributed to electrochemical interactions at mineral–fluid interfaces, providing valuable insights into the composition and texture of subsurface materials.

The IP effect was first recognized by Conrad and Marcel Schlumberger in the 1920s and 1930s, who observed voltage decays associated with ‘provoked’ polarization currents in the subsurface (Seigel et al. 2007; Kemna et al. 2012). These field-based time-domain observations represent the earliest detection of IP effects in geophysics. Initially developed for mineral exploration in the 1960s and 1970s, the IP method soon proved to be a powerful tool for identifying disseminated mineralization within porous or fractured host rocks, where electrochemical polarization effects are significant. Disseminated sulfides, commonly associated with economically important ore deposits of gold, silver, nickel and copper, produce strong IP responses (Jones et al. 1982; Seigel et al. 2007; Revil et al. 2022). Because resistivity values alone can yield ambiguous geological interpretations, IP measurements provide complementary information that enhances subsurface characterization.

Later, Pelton et al. (1978) demonstrated that the IP phenomenon could be quantitatively modelled using the Cole–Cole relaxation model, originally developed for describing complex dielectric behaviour (Cole & Cole 1941). In ERT, when measurements are performed at sufficiently low frequencies (typically with periods longer than a few seconds), resistivity and IP responses can be effectively separated within the same data set. The amplitude ratio between measured voltage and injected current reflects resistivity, while the IP effect is contained in the phase shift between the two signals. Although there are several study cases to extract IP data from other types of electromagnetic such as in controlled source electromagnetics (CSEM) (Hohmann et al. 1970; Marchant et al. 2013; Porté et al. 2023), in AEM (Kratzer & Macnae 2012; Macnae 2016; Viezzoli & Manca 2020), and in MT (Gasperikova & Morrison 2001; Gasperikova et al. 2005), ERT method has obviously used normally to extract IP information with its facilities (Seigel et al. 2007).

Recent technological advances such as the nodal geophysical recording systems now common for seismic, EM and ERT (Lajau-nie et al. 2019; Dean 2024), have led to the use of networks of independent receiving and injecting dipoles that can be easily installed on the surface without any geometrical limitations. The dipoles are therefore decoupled and separately distributed, eliminating the need for long cables, and removing cross-talk between the cables that limit the performance of ERT-IP. Additionally, those nodal systems perform continuous recording of the full time-series for both voltages and current, which opens possibilities for flexible and adaptive stacking and usage of advanced signal post-processing tools. On the other hand, flexible high-voltage power transmitters capable of injecting DC into the ground with current intensities reaching 40A are commercially available (e.g. TXU-30, Phoenix; GGT-30B, Zonge). Those technological advancements allow improvement of signal-to-noise ratio and then increased distance between emitting and receiving systems, enabling the investigation of subsurface resistivity patterns at greater depths. Due to these improvements, Deep ERT has been more actively studied in the last decade (Lajau-nie et al. 2019; Brettaudeau et al. 2021; Balasco et al. 2022; Porté et al. 2023).

While ERT methods have shown their potential for deep resistivity imaging, IP methods for deep depths have very few reported case studies. IP signal is always weaker than DC signal and it makes that IP exploration can be disturbed much by a number of different source, including telluric currents, industrial current interference, major power lines (Liu et al. 2016; Revil et al. 2022). To resolve the weak IP signal problem, the principal ways are utilizing survey planning, such as employing high-power source current transmitter and unconventional array (Liu et al. 2017). Another way is adopting advanced data processing procedure to suppress the electromagnetic interference from IP signal. Digital filter based on the Discrete Fourier Transform is widely used in a variety of FDIP instruments with designing high-pass, low-pass, band-pass and notch filter (Liu et al. 2016).

Designing an effective Deep ERT-IP surveys necessitates a balanced approach that accounts for both data quality and cost-efficiency. Critical factors for successful data acquisition include validating the dc—direct current assumption to maintain data accuracy and numerical modelling validity, achieving a high signal-to-noise ratio for reliable measurements, ensuring sufficient sensitivity to detect deep targets and obtaining an adequate spatial resolution at the desired target depth. In addition to technical considerations, economic factors such as cost-effectiveness and survey duration are essential components of survey planning. As the target depth increases, the survey size and costs grow, while the target detectability and spatial resolution decreases. These elements should be carefully analysed quantitatively and integrated systematically, tailored to the specific environmental and logistical conditions of each survey site. However, these preliminary procedures before field surveys are easily neglected in many field cases.

In this paper, we look into the requirements as conducting Deep ERT-IP surveys and illustrated entire workflow to address practical challenges encountered in real field surveys. First, we define what is practical challenges for survey designing for Deep ERT-IP surveys and demonstrated pre-test procedure to show how to decide survey geometries based on quantitative analysis. These strategies
<figure><img src="image_2.png" /><figcaption>**Figure 1.** Schematic illustration defining the required aperture of a survey configuration ($d_{aperture}$) to achieve a maximum investigation depth ($z_{max}$) with full lateral coverage ($d_{full}$). L represents the distance between the two outer electrodes in a symmetrical, collinear four-electrode arrays.</figcaption></figure>

are demonstrated with a real case study in the Koillismaa Linear Intrusion Complex in Finland where we conducted a 3-D Deep ERT-IP survey dedicated to mineral exploration at depths greater than 1 km. To efficiently cover a large area and economically investigate deep targets in 3-D, we employed nodal systems which record full time-series data at each station, thereby avoiding the use of excessively long cables, and applied an optimized unconventional array specifically designed for deep exploration. Finally, we will show recovered 3-D resistivity and IP images aligning existed other data for targeting very deep depths (> 1.5 km), which is realized by adopting our reliable workflow.

# 2. SURVEY DESIGN FOR DEEP ERT-IP

## 2.1. Practical challenges in survey design

Designing ERT and IP surveys for deep targets (e.g. >1 km depth) presents several technical and logistical challenges that must be carefully considered in survey planning. These challenges intensify with increasing target depth and include the following:

1. **Resource Allocation:** Mobilizing the necessary geophysical instruments and logistical support, while controlling costs and survey duration.

2. **DC Assumption:** Ensuring that the dc—direct current assumption remains valid.

3. **Signal-to-Noise Ratio (SNR):** Achieving a sufficient SNR for effective data acquisition.

4. **Sensitivity:** Ensuring that the system is sensitive enough to detect deep targets.

5. **Spatial Resolution:** Maintaining an adequate spatial resolution at target depths.

An effective survey design must strike a balance between material, cost and time (Challenge #1) while meeting the requirements for DC validity, SNR, sensitivity and spatial resolution (Challenges #2–5). This section provides a general framework for addressing these practical difficulties, guiding survey designs for Deep ERT-IP exploration.

### 2.1.1. Selecting survey geometries

The ideal survey setup would maximize coverage (using large offsets), maintain high density (with closely spaced transmitters and receivers dipoles), use a high frequency of injected current (to enable stacking without losing time and compromising dc—direct current assumptions), and record high-SNR data. However, logistical constraints typically require adjustments to this ideal configuration. An

optimized survey design should aim to minimize field effort while ensuring that the data acquired is sufficient to achieve the survey's imaging goals.

In practice, an ideal Deep-ERT acquisition geometry would resemble conventional arrays such as Dipole–Dipole or Wenner–Schlumberger arrays, with multiple spacings for the transmitter (AB) and receiver (MN) pairs. This approach ideally involves multiplexed multi-electrode resistivity systems arranged on extended cable setups, providing the largest number of quadrupoles for a given number of electrodes and thus ensuring high-spatial resolution.

The required aperture of a survey configuration ($d_{aperture}$), defined as the total horizontal extent of electrode deployment, depends on two main factors: the lateral distance over which full subsurface illumination is required ($d_{full}$), and an additional length proportional to the maximum target depth ($z_{max} = \alpha L$) (Fig. 1). This relationship can be expressed as

$$d_{aperture} = d_{full} + L, \quad (1)$$

where $\alpha$ denotes the proportionality factor between the largest quadrupole spacing ($L$) and the achievable penetration depth ($z_{max}$). The effective maximum imaging depth in inverted resistivity sections is conventionally determined using the median depth principle introduced by Edwards (1977), adjusted by an empirical depth factor. The empirical depth of investigation for collinear arrays such as the Schlumberger, Wenner, or Dipole–Dipole configuration is typically about 17–25 per cent of the outermost electrode spacing (Barker 1989). For practical survey design aimed at reliable deep imaging, a more conservative estimate of 10–15 per cent of the maximum electrode spacing is often recommended (Sean 2019).

For instance, to image a target at a depth of 1 km ($z_{max}$) with a required lateral illumination ($d_{full}$) of approximately 4 km with assuming $\alpha$ is 17 per cent (i.e. $L \approx 5.9$ km), a total cable length of around 10 km would be required to achieve the necessary survey aperture ($d_{aperture}$) and ensure sufficient sensitivity. However, deploying such extensive cable arrays can become impractical, particularly for large-scale 3-D surveys.

An alternative approach is to use independent nodal systems that are GPS-synchronized and recording voltage continuously, combined with a mobile transmitter setup (Carrier et al. 2019; Lajaukie et al. 2019; Mazzini et al. 2021). In this type of configuration, receiver nodes are fixed with short MN distances, while the transmitter (Tx) is moved sequentially. The Pole-Dipole approach (Loke et al. 2022) is often used as it minimizes logistics by conserving a far-field common electrode for all the injection dipoles. The far-field pole for such an array would need to be located 5–10 times
<figure><img src="image_3.png" /><figcaption>**Figure 2.** Voltage normalized amplitude and phase shift due to electromagnetic diffusion in an homogeneous half-space for different resistivity from 1 Ω.m to 100 kΩ.m for three different transmitter–receiver distance : (a) 1 km, (b) 4 km, (c) 12 km. Normalized amplitude must be close to 1 to validate the dc—direct current assumption for resistivity imaging. Phase shift appearing at higher frequency is due to EM induction and must be close to 0° to validate dc—direct current assumption.</figcaption></figure>

the maximum Tx–Rx distance from the survey area, often requiring distances exceeding 50 km, which is usually unfeasible.

A more practical approach for large depth is to employ a Dipole–Dipole array with both Tx and Rx nodal stations, arranged with irregular spacings that more easily fit with the various field constraints. Rx dipole lengths are kept relatively short (usually tens of meters) to keep Rx installation fast and easy, while ensuring adequate voltage measurement dynamics. Because voltages are usually small and highly variable on short dipoles and large distance, high-quality low noise 24bit or even 32bit analogue-to-digital converter (ADC) data loggers with appropriate gain must be preferred. Still, flexibility in the receiver configuration is limited to the number and spacing of Rx sites, directly restricted by the number of instruments available and deployment time. For the transmitters, flexibility lies in varying the dipole lengths (100 m to several km), orientations, Tx position spacing, overlap and recording duration (linked to the injected current period and stacking number). This setup sacrifices some data density relative to the streamer array but is far more feasible in the field, especially as the number of nodal data loggers and fixed AB and MN distances limit deployment options, and as the accessibility in the field is limited (land occupation, anthropic structures and noise, topography, etc). Quantifying this data density reduction is crucial to control its impact and determine the degree to which flexible parameters can be adjusted without compromising data quality.

### 2.1.2. Period of injected current

After defining the target depth and array aperture, selecting an appropriate period for the injected current is critical (Challenge #2). Although dc—direct current theoretically has zero frequency, practical ERT surveys inject alternating currents (AC), usually three-level rectangular wave (+A/0/–A/0, referred as 50 per cent duty cycle) in case of measuring time domain IP or symmetric square

wave (+A/–A, referred as 100 per cent duty cycle). If the injection period is too short (i.e. at higher frequencies), electromagnetic induction effects become significant. The phenomenon, known as the skin effect, causes an exponential reduction in signal amplitude at greater depths.

The skin depth ($δ$), defined the depth at which the amplitude of an electromagnetic wave decreases to approximately 37 per cent (1/e) of its original amplitude in a homogeneous medium (Sheriff 2002), can be expressed as

$$δ = \sqrt{\frac{\rho}{\pi f \mu}} \approx 503 \sqrt{\rho T}, \qquad (2)$$

where $ρ$ is the resistivity of the medium (Ω.m), $μ$ is magnetic permeability (1.2566 × 10⁻⁶ H m⁻¹ for most materials) and $T$ is injection period (s), which is the inverse of the source frequency ($f$). For resistivity imaging, the dc—direct current assumption is generally valid when the targeted investigation depth ($z_{max}$) is less than about one-third of the skin depth ($z_{max} < δ/e$). Therefore, to satisfy the dc—direct current assumption and ensure accurate imaging at the maximum target depth, the period $T$ must meet the following condition

$$T > 3.6 \times 10^{-5} \frac{z_{max}^2}{\rho_{min}}, \qquad (3)$$

where $ρ_{min}$ is the minimum expected resistivity of the survey area.

Fig. 2 depicts for three different in-line Tx–Rx distances (1, 4 and 12 km) the normalized voltage amplitude and phase shift due to electromagnetic diffusion in a homogeneous half-space for different resistivity from 1 Ω.m to 100 kΩ.m for injection frequencies ranging from 1000s to 10 kHz. Note that in the following, for the sake of clarity, we decided to represent the phase of the frequency-domain data in degrees to represent clearly the large phase shifts induced by EM effect, whereas IP as a physical parameter in the models will be represented in mrad, as it is the convention. Amplitude and phase variations with frequency reflects EM wave propagation.
Normalized amplitude at a given frequency must be close to 1 to validate the DC assumption for resistivity imaging.

Phase shift must be close to 0° to separate the IP response from EM induction. For instance, in a 1 kΩ.m medium, DC assumption is respected below 1 s up to Tx-Rx offset of 4 km, but offsets as large as 12 km require periods below 10 s to keep an amplitude error below 1 per cent. The criterion for IP is more restrictive. At 10 s in a 1 kΩ.m, phase due to induction is 0.2° at 4 km distance, but is 1.5° at 12 km. Such a phase shift which appears only for the largest offsets could easily be interpreted as deep IP. In the presence of more conductive sediments or altered formations with resistivities as low 10 Ω.m, the amplitude error at a 10 s period can reach 10 per cent at 4 km and 50 per cent at 12 km, while the phase shift increases to 10° at 4 km and 25° at 12 km, effectively preventing any reliable IP interpretation without full electromagnetic modeling.

This condition is highly sensitive to site-specific factors as it depends on the underground resistivity. For shallower targets or highly resistive media, such as hard rock environments, the DC assumption is often satisfied even at relatively high frequencies. On the other hand, lower frequencies are necessary for more conductive environment. This usual criterion is not valid for IP, however, especially if small IP signature coming from deep chargeable bodies are investigated. In such cases, lower frequencies must be used.

Furthermore, the optimal injection period for IP measurements also depends on the target mineralogy. Grain size significantly influences IP responses: fine-grained minerals such as clays exhibit peak responses at higher frequencies compared to massive sulphides or graphite, which respond better at lower frequencies. Broad-band spectral IP analysis can provide additional insights, such as distinguishing between massive sulphides and graphite, as demonstrated by (Pelton et al. 1978). Therefore, selecting sufficiently low frequencies based on the mineralogical and geological context of the target site is highly advisable.

### 2.1.3. Signal-to-noise ratio and sensitivity

Maintaining a sufficient SNR and sensitivity (challenges #3 and #4) at depth requires addressing two key questions: (1) What is the expected magnitude of the target response, considering the geological feature under investigation? (2) What acquisition parameters ensure that the target response stands out clearly from background noise? The first question can be answered through preliminary forward modelling, using a simplified resistivity model and survey setup, incorporating unitary currents, transmitter-receiver dipole lengths, and the required array aperture. The second question requires field-based preliminary tests to assess local conditions.

A crucial first step is to estimate the optimal ground contact resistance for the transmitters, as this determines the maximum injectable current based on the transmitter's voltage and power constraints. This factor may also influence the selection of transmitter equipment, as different units vary in power capacity and maximum voltage output. The second step involves continuous passive electric recordings over a 1–2 d period at representative test sites. This provides a baseline noise assessment, considering spatial and temporal variations.

To evaluate sensitivity and SNR prior to fieldwork (detailed in Section 3), we performed 1-D EM simulations using a semi analytical solution (Key 2009) for a single transmitter and in-line receiver dipoles at distances ranging from 0 to 10 km. The model was computed for a unitary transmitter dipole moment (1 A–1 m) and receiver dipole length (1 m), across frequencies from 0.01 Hz to

<figure><img src="image_4.png" /><figcaption>**Figure 3.** 1-D resistivity model used for the sensitivity and signal to noise estimation</figcaption></figure>

10 kHz. Frequency-domain EM modelling was preferred, as electrical noise varies significantly with frequency and also provides an estimate of expected inductive effects for the given resistivity model. Of course, more accurate 2-D or 3-D simulations can be used also, especially if specific 2-D or 3-D structures are justified by the geological settings. dc—direct current based simulation could serve as an alternative with collected noise at specific frequency.

For this study, we used a simple 1-D model consisting of a 500 Ω·m layer at depths of 1500–1800 m within a 5 kΩ·m homogeneous background (Fig. 3). The target structure was modelled as a 100 Ω·m layer, 100 m thick, located at 1.5 km depth on top of the 500 Ω·m layer. The computed electric responses were scaled based on Tx and Rx dipole length, and Tx current source, then compared to estimated noise levels for different stacking times.

Fig. 4 presents model responses at Tx-Rx offsets of 1, 5 and 10 km for different Tx configurations (1A–100 m, 1A–1 km, 5A–1 km) and recording durations (1 and 30 min) to evaluate the impact of stacking on the SNR. Every Rx dipole length is set to 50 m. The background response is shown in blue, the total response (background + target) in red and the isolated target response in purple, all plotted on a logarithmic scale. The Rx and Tx dipole lengths, current and stacking time can be adjusted to assess how the target response emerges relative to noise at various Tx-Rx distances. In practice, the target response must be at least an order of magnitude higher than the noise level for reliable detection. To maintain clarity, model responses in Fig. 4 are kept unitary, while noise is scaled according to the current and Rx/Tx dipole lengths, expressed in V A⁻¹ m⁻².

The recorded electrical noise is highly frequency-dependent, with peaks at 50 Hz and higher harmonics from anthropogenic sources. The noise level reaches a minimum around 1 s but increases at lower frequencies due to telluric currents. The first type of noise is site-dependent, while the second varies based on regional conditions, time of day, solar activity and local effects.

The phase shift of the electric field response is analysed for the same model conditions (Fig. 5) as the amplitude responses, with Tx-Rx offsets of 1, 5 and 10 km. The results in Fig. 4 indicate that while amplitude responses at 10 km do not exhibit significant induction effects below 1 s, phase variations become more pronounced. This suggests that a 1 s period is suitable for resistivity imaging but inadequate for ERT or IP, both of which require longer injection periods of at least 10 s. To account for potential local variations in resistivity, a margin should be considered, with periods of approximately 10 s for resistivity imaging and 20 s for IP, ensuring the mitigation of inductive effects at all offsets.
<figure><img src="image_5.png" /><figcaption>**Figure 4.** A comparison of 1-D modelling responses with and without a conductive layer and passive electric noise recorded at the field site. The modelled electrical responses are shown for a case with a conductor, without a conductor and their difference. Each row presents the electrical responses based on the offset between the transmitter and receiver. The modelling was conducted using various configurations: (a) source amplitude of 1 A with a source dipole length of 100 m, (b) source amplitude of 1 A with a dipole length of 1 km, (c) source amplitude of 5 A with a dipole length of 1 km and (d) source amplitude of 5 A with a dipole length of 1 km. Passive signal recording is 1 min in panels (a)–(c), and 30 min in panel (d).</figcaption></figure>

By selecting a 10-s injection period, we observe that the 1 A–100 m scenario (Fig. 4a) fails to detect the target, while the 1 A–1 km configuration (Fig. 4b) provides insufficient signal strength beyond 1 km. However, if contact resistivity estimation and transmitter unit performance indicate that a maximum current of 5 A can be achieved (scenario 5 A–1 km in Fig. 4c), sufficient SNR can be obtained at 1 km, an almost satisfactory response at 5 km, but inadequate detectability at 10 km as the recording time is 1 min.

Increasing recording time is the final approach to improving SNR. The last scenario in Fig. 4 demonstrates that a 30-min recording (180

stacks) with the same configuration significantly enhances SNR and target detectability at all offsets (Fig. 4d). However, this extended acquisition duration has major implications: recording at 10-s injection periods for each transmitter position limits survey efficiency. A standard working day allows for only about a dozen transmission points, which must be factored into survey design considerations. The primary constraint is the reduced number of achievable transmitter positions within a fixed survey duration, necessitating a balance between resolution, depth penetration and acquisition time.
<figure><img src="image_6.png" /><figcaption>**Figure 5.** Phase shift of the modelled electric field response as a function of frequency for different transmitter–receiver offsets (1, 5 and 10 km) from the same 1-D model presented in Fig. 3. The modelled phase responses are shown for a case with a conductor, without a conductor and their difference.</figcaption></figure>

## 2.2. Comparing survey arrays to evaluate spatial resolution and survey efficiency

### 2.2.1. Unconventional array using nodal systems

The traditional cabled dipole–dipole (DPDP) array (Fig. 6a) employs a pair of current electrodes and a separate pair of potential electrodes. One of the advantages of DPDP array is its ability to acquire multiple-spacing data with the same configuration, as electrode separations increase in multiples of a base spacing (e.g. *n* = 1, 2...; Fig. 6b). With its small dipole spacing, DPDP array offers dense data coverage but has a limited depth of investigation. In contrast, large dipole spacing extends the investigation depth but results in coarser spatial resolution and demands longer cables, complicating field logistics. Thus, this configuration is well-suited for dense, linear profiles with providing high-spatial resolution with small, regular spacing but it requires longer survey profiles and increased electrode spacing.

On the other hand, configurations using nodal systems eliminate the need for excessively long cables and offers flexible electrode separation. In this configuration, transmitter stations—functioning as the location of current electrodes—are spaced significantly farther apart (e.g. 1 km) to increase the dipolar moment and enhance the SNR for deep targets. Meanwhile, receiver stations can be spaced similarly to the Tx stations or place at smaller intervals to preserve spatial resolution. Each Rx station is equipped with a nodal data logger, maintaining a short dipole spacing (e.g. 50 m) to facilitate efficient field deployment (Fig. 6c). To avoid singularity issues, Tx and Rx stations should be altered along the survey line.

For 3-D deployment, orthogonal dipole pairs can be established from nodal data loggers (Loke et al. 2022). This configuration optimizes the balance between requiring resources, such as the number of transmitter–receiver pairs and lengths of cable, while maintaining adequate SNR to resolve deep-seated structures (Fig. 6d). Given its specialized application for deep target imaging, this study refers to this configuration as the ‘Deep Focused Array (DFA)’.

### 2.2.2. 2-D synthetic example

Spatial resolution refers the ability to distinguish closely spaced two or more different bodies. For such an unconventional and irregular ERT survey, quantifying resolution can be challenging because it

depends on many parameters, including surface Tx and Rx spacing, dipole lengths, background resistivity and noise level. As a decision-making tool, synthetic modelling and inversion tests under various scenarios are essentially helpful. Here, we demonstrate how the loss in resolution can be quantified through 2-D examples, and how we can make a decision with balancing between an unconventional nodal array, DFA, against a conventional cabled dipole–dipole array with different electrode spacing.

We generated a simple 2-D ERT-IP model (Fig. 7) informed by borehole resistivity logs at the field site (Tirronniemi et al. 2024). The model includes a background resistivity $10^4 \Omega \cdot m$ with a chargeability of 5 mV $V^{-1}$, and a 500 m-thick resistive layer ($10^5 \Omega \cdot m$) at 500 m depth, interrupted by a 4 km lateral gap. Two less resistive anomalies ($10^3 \Omega \cdot m$), measuring $1 \times 1 \text{ km}^2$ and $2 \times 1 \text{ km}^2$, are located at depths of 1 and 1.5 km, respectively, each with a high chargeability of $100 \text{ mV } V^{-1}$.

We tested four conventional dipole–dipole configurations: a multiple-spacing array (250 m, 500 m, 1 km), and mono-spacing array a 250 m array, a 500 m array and a 1 km array. Additionally, we implemented two DFA setups, one with 500 m Tx spacing and one with 1 km Tx spacing. Both DFA configurations used Rx stations at 500 m intervals with 50 m dipole lengths, chosen to keep field deployment economically feasible.

All 2-D ERT-IP forward modelling and inversion were performed using the open-source geophysical software package, pyGIMLi (Rücker et al. 2017). Forward modelling employed a finite element method on unstructured triangular mesh, and a Gauss–Newton approach minimized the cost function to recover resistivity and IP models. The unstructured meshes were generated automatically depending on the survey geometry and model configurations. The model with the DFA required much refining mesh near receiver electrodes due to the smaller dipole spacing (50 m) compared to the conventional dipole–dipole array spacings (250 m, 500 m and 1 km). For 2-D inversion, we added 3 per cent Gaussian random noise and the weight for resistivity data and 1 per cent for IP data, and the weight value of smoothness constraints is 10 for both models.

The inversion results from the conventional DPDP array highlight the effect of electrode spacing on resolution. The multiple-spacing (250 m, 500 m, 1 km) (Fig. 7b) and the dedicated 250 m array (Fig. 7c) provided the most accurate resistivity and phase models, effectively resolving both shallow and deep anomalies. As the spacing increased to 500 m and especially 1 km (Figs 7d and
<figure><img src="image_7.png" /><figcaption>**Figure 6.** (a) Configuration of conventional dipole-dipole (DPDP) array and (b) examples of corresponding pseudo-sections for spacing factors *n* = 1, 2. (c) Configuration of the deep focused array (DFA) and (d) an example of the corresponding pseudo-section, where the transmitter (Tx) spacing is twice the receiver (Rx) spacing. In the DFA setup, Rx dipole electrodes are connected to nodal data loggers at each Rx station and Rx dipole spacing (*α*<sub>Rx</sub>) can be significantly smaller than the Tx dipole spacing (*α*<sub>Tx</sub>).</figcaption></figure>

e), the recovered anomalies became more smoothed, with reduced contrast, making it increasingly difficult to differentiate the separate conductive bodies.

In contrast, the DFA exhibited strong spatial resolution across both tested configurations (Fig. 7f), maintaining reasonable quality imaging even at 1 km Tx spacing (Fig. 7g). The most remarkable feature is the resistivity inversion results with DFA shows better spatial resolution rather than DPDP array, which tends to smooth anomalies horizontally. However, the DFA shows vertically extended resistive layers at the boundary areas of the profile, likely due to resolution limitations in those regions.

Regarding to IP inversion results, both DFA configurations successfully recovered deep IP anomalies. Although the 1 km Tx spacing produced slightly weaker phase values compared to the 500 m Tx spacing, the recovered structures of deep IP anomalies remained sufficiently well-resolved interpretation.

To compare the logistical demands associated with different survey arrays from the synthetic model in Fig. 5, we illustrated requiring resource factors in Table 1 such as total electrode count, dipole configurations, data set size, working time and cable requirements. The working time is set at 30 min per measurement to ensure high data quality with a low signal-to-noise ratio, providing sufficient stacking. While DPDP arrays with small spacing (250 m and mixed-spacing) yield the best spatial resolution, they also demand significant logistical effort, requiring extensive cabling (> 18 km), a large number of electrodes (73), and long working times (2130–3660 min, depending on setup), making them impractical for large-scale deep surveys. In contrast, the 1 km DPDP array greatly reduces survey time (480 min) but sacrifices resolution, particularly in detecting deep chargeability anomalies.

On the other hand, DFA offers a more efficient alternative, reducing field effort while maintaining reasonable resolution. The

required cable length for each Tx and/or Rx station of DFA is dramatically lower (1750 m Rx for Tx 500 m spacing, and 1750 m Rx for Tx 1 km spacing), making field deployment significantly more manageable. The working time for DFA is also reduced (1050 min for Tx 500 m and 480 min for Tx 1 km), making it comparable to the most efficient DPDP setup but with notably better resolution than the 1 km DPDP array.

However, Table 1 illustrates that the trade-off with DFA is a reduction in data density and potential underestimation of anomaly values. The number of data sets collected by DFA (1332 for Tx 500 m, 666 for Tx 1 km) is significantly lower than in high-resolution DPDP arrays (3397 for DPDP multispacing, 2485 for DPDP 250 m), which impacts the accuracy of recovered values. While DFA does not fully match the resolution of high-density DPDP arrays, it still outperforms DPDP with 1 km spacing in both resistivity and phase models while requiring significantly fewer resources. This example highlights the practical trade-offs involved in survey design, demonstrating the balance between data resolution and field efficiency.

### 3. DATA COLLECTION, PROCESSING AND QC

We demonstrate the complete process of the field exploration using Deep ERT-IP for targets depth greater than 1–2 km over the Koillismaa Linear Intrusion Complex, Finland. As discussed in the previous chapter, several challenges needed addressing: (1) covering an area as large as possible while preserving sensitivity, high SNR and spatial resolution for deep depths without losing too much near-surface information, and (2) limited resources such as number of instruments, poles, cables, man power and time. The approach proposed here ensured the depth of investigation for regions deeper
<figure><img src="image_8.png" /><figcaption>**Figure 7.** 2-D ERT-IP synthetic model designed for deep anomaly detection. (a) The true model consists of a disconnected resistive layer and two conductive anomalies with high chargeability. Recovered 2-D resistivity and IP models are shown for different electrode configurations: (b) dipole-dipole array with mixed electrode spacings (250 m, 500 m, 1 km), (c) dipole-dipole array with 250 m spacing, (d) dipole-dipole array with 500 m spacing, (e) dipole-dipole array with 1 km spacing, (f) deep-focused array with 500 m Tx spacing, 500 m Rx station spacing, and 50 m dipole spacing at each Rx station and (g) deep-focused array with 1 km Tx spacing, 500 m Rx station spacing, and 50 m dipole spacing at each Rx station.</figcaption></figure>

<table><caption>**Table 1.** Example of requiring resources depending on survey arrays from the synthetic example in Fig. 7.</caption><thead><tr><th rowspan="2"></th><th colspan="3">Cable system</th><th colspan="3">Nodal system</th></tr><tr><th>DPDP<br>multi-spacing</th><th>DPDP<br>250 m</th><th>DPDP<br>500 m</th><th>DPDP<br>1 km</th><th>DFA<br>(Tx 500 m, Rx 500 m)</th><th>DFA<br>(Tx 1 km, Rx 500 m)</th></tr></thead><tbody><tr><td>No. Total<br>electrode<br>Electrode spacing<br>(Tx/Rx)</td><td>73<br>250 m/250 m</td><td>73<br>250 m/250 m</td><td>37<br>500 m/500 m</td><td>19<br>1 km/1 km</td><td>111<br>500 m/500 m for station<br>50 m for dipole</td><td>93<br>1 km/500 m for<br>station 50 m for<br>dipole</td></tr><tr><td>No. Dipole<br>(Tx/Rx)<br>No. Data<br>Working field time<br>(No. Tx dipole ×<br>30 min)</td><td>122/122<br>3397<br>3660 min</td><td>71/71<br>2485<br>2130 min</td><td>35/35<br>595<br>1050 min</td><td>16/16<br>136<br>480 min</td><td>35/35<br>1332<br>1050 min</td><td>16/35<br>666<br>480 min</td></tr><tr><td>Length of cable</td><td>> 18 km</td><td>> 18 km</td><td>> 18 km</td><td>> 18 km</td><td>500 m for Tx/1750 m for<br>Rx</td><td>1 km for Tx/1750 m<br>for Rx</td></tr></tbody></table>
<figure><img src="image_9.png" /><figcaption>**Figure 8.** Generalized geological map of the northeastern part of the Fennoscandian Shield, showing the locations of early Palaeoproterozoic mafic-ultramafic intrusions (modified from Karinen, 2010).</figcaption></figure>

than 1–2 km with using survey resources economically. Additionally, we utilized recording instruments capable of acquiring full time-series data continuously for 24 hr with battery support, enabling coverage of extensive area approximately 90 km² (18 km by 5 km) without the need for excessively long cables (>2 km).

## 3.1. Background study of the site

The survey area, approximately 150 km northeast of Oulu, is located in Koillismaa Linear Intrusion Complex in northeastern Finland (Fig. 8), which is composed of the Koillismaa as well as Näränkävaara intrusion, connected by an ~50 km long high bouguer anomaly zone. The later has been related to a deep-seated magma flow system connecting these exposed intrusions (Salmirinne & Iijina 2003; Karinen et al. 2010).

In 2021, the drilling of a 1.7 km deep borehole was completed, targeting the central part of the gravity anomaly that links the surface outcrops of the mafic-ultramafic intrusions from Koillismaa in the west to Näränkävaara in the east (Karinen et al., 2021). From the deep borehole, geophysical logging and core data analysis were conducted (Fig. 9). The uppermost section of the borehole until 510 m traversed through Archean basement gneisses. This section is characterized by highly altered gneisses, with occasional intrusions of leucogranite dykes and diabase dykes. From depths of about 1250 m, the borehole intersects ultramafic rocks including peridotite and pyroxenite, together with quartz diorite. Below 1410 m, granites appear interlayered and locally mixed with mafic-ultramafic rocks, which are interpreted as the main source of the pronounced gravity anomaly (Tirronniemi et al. 2024). The drilling generated over 1700 m of diamond drill core intersecting unaltered mafic-ultramafic rocks at depth, which has been dated to be of same age as the exposed intrusions.

The resistivity logging data are measured in frequencies 0.1, 10 and 500 Hz. These log data show very high overall values around 10 kΩ.m, with some smooth variations. While the overall trends are less pronounced compared to other geophysical measurements, a distinct increase in resistivity is observed near a depth of 400 m,

indicating a potential lithological interface. Additionally, a gradual resistivity increase is noted between 1250 and 1500 m. Two small but clear anomalies at 1550 and 1700 m also show a localized decrease in resistivity down to few hundred Ω.m. However, the resistivity values vary significantly, making it difficult to determine clear patterns from depths of 1500 m onwards, which are mainly associated with peridotite. Higher resistivity values are indicative of more resistive minerals such as quartz diorite and peridotite, while relatively conductive zones are correlated with altered Archean basement gneisses. Resistivity data ranges at higher source frequency (10 and 500 Hz) shows similar trends with data at 0.1 Hz but more conductive near approximately a depth of 1500 m.

IP responses were obtained by percent frequency effect (PFE), which quantifies the percent difference in apparent resistivity measured at two different frequencies

$$ \text{PFE} = \frac{100(\rho_{dc} - \rho_{ac})}{\rho_{ac}}, \qquad (4) $$

where $ρ_{dc}$ and $ρ_{ac}$ are the low- and high-frequency apparent resistivities. We obtained two PFE data using two different high-frequency apparent resistivities: $PFE_{10 Hz}$ is calculated from resistivity at 0.1 and 10 Hz source frequency and $PFE_{500 Hz}$ is from resistivity at 0.1 and 500 Hz source frequency, respectively. Near depth between 400 and 700 m, the high $PFE_{500 Hz}$ values in this depth range correspond to zones of altered Archean basement gneiss, diabase and gabbro. These rock types are characterized by varying degrees of alteration and mineral content that contribute to the elevated chargeability. At greater depths, starting from 1250 m, both $PFE_{10 Hz}$ and $PFE_{500 Hz}$ show another set of prominent peaks corresponding to the two low-resistivity zones. Notably, peridotite layers around 1500 m depth demonstrated distinctive magnetic and IP responses in petrophysical laboratory measurements. The potential spatial extend of those anomalies is one of the questions addressed by this Deep ERT-IP survey.
<figure><img src="image_10.png" /><figcaption>**Figure 9.** Stratigraphic column and geophysical logging data across various source frequencies and rock types. The resistivity data presented are measured at source frequencies of 0.1, 10 and 500 Hz. Percent frequency effect (PFE) are calculated based on the percent differences between resistivities at 0.1 and 10 Hz (PFE₁₀ Hz), and between 0.1 and 500 Hz (PFE₅₀₀ Hz), respectively.</figcaption></figure>

## 3.2. Survey design, deployment and data collection

The survey was designed following the recommendations outlined in the previous section (Chapter 2) and the geophysical field work was conducted from 2023 August 18 to September 2. Given the target depth of 1.5 km, a lateral extend on the order of a few kilometers, as well as some field constraints, a 3-D array with a 6 × 16 km² aperture was selected (Figs. 10 and 11). Based on the logging data (Fig. 9), the estimated background resistivity is around 10 kΩ.m. Therefore, a fundamental injection period of 8 s was chosen to ensure that the dc-direct current assumption remains valid for both resistivity and IP imaging at Tx-Rx offsets up to 16 km, resulting in amplitude error below 0.1 per cent and phase errors below 0.3°. Although conductive zones could increase EM-induced phase shifts at larger offsets. Such effects would require corrections based on EM modelling to avoid misinterpretation of subtle IP signatures hidden by EM effects. Hopefully such corrections were unnecessary at this site due to the consistently high-background resistivity.

To optimize sensitivity and SNR, we employed a pool of 114 nodal electric receivers, each equipped with two orthogonal 50 m dipoles using non-polarizable electrodes (Fig. 9). Because the field survey was also mutualized with other geophysical methods (in particular MT and CSEM), we employed three types of nodal data loggers each with a distinct sampling frequencies to mobilize receivers as much as possible: 19 ADU07-e and ADU08-e (Analog Digital Unit, Metronix; also equipped with magnetometers) recorded at 4096 Hz, 77 RAUs (Remote Autonomous Unit, Sercel seismic node + BRGM adapter) recorded at 2000 Hz, and 18 VFullwers (FWVs, IRIS Instruments) recorded at 100 Hz (see Table 2 and Figs 11 a-c). All the recording systems were installed during the first days of the campaign on a coarse grid with approximate spacing of 500 m, and started synchronous data acquisition at the beginning of the source injection.

Preliminary field tests show that injections poles installed in swampy areas minimized grounding resistance, enabling injection current ranging from 3 to 8 A, depending on site conditions and type of transmitter used. Following the procedure described above, we selected transmitter dipoles approximately 1 km length, and injected current using a symmetric square waveform across multiple frequencies, ranging from 512 Hz down to 0.125 Hz, and to

0.03 125 Hz (32 s period) at selected stations. To ensure an adequate number of stacks and improve signal quality, data were acquired at the fundamental frequency of 0.125 Hz for approximately 30 min per transmitter position. In this study, data at 0.125 Hz were primarily interpreted for Deep ERT-IP, while the 0.03 125 Hz recordings were used for data quality control and verification of spectral consistency.

To optimize field operations, adjacent transmitter dipoles shared common poles, minimizing logistical effort. Switching times between transmitter positions were reduced by installing multiple poles and cables simultaneously, allowing for a roll-along deployment. Three different Tx instruments were used sequentially to adapt to site accessibility and grounding conditions: GGT-3 (Zonge), TXU-30 (Phoenix) and TIP6000 (IRIS Instruments) (Figs 11 d-f).

Each Tx system presents distinct advantages and limitations in terms of maximum power, maximum voltage, sampling frequencies and level of portability. Their technical specifications are summarized in Table 2. The average injection current amplitudes obtained were 3–4 A for the GGT-3 and TIP6000, and 7–8 A for the TXU-30. In order to remove each data logger’s signature at the processing (internal calibration, low pass-filter, etc), injected currents were recorded continuously for each transmitter using all three types of data loggers (I-Fullwave, RAU and ADU) equipped with high-precision closed-loop hall effect current sensors.

The main survey activities took place from 2023 August 21 to 30. This period included the installation of nodal recording systems and the recording of natural signals from August 21 to 24, followed by current injection from August 25 to 30. In total, 25 grounded 1 km dipoles were connected by electrical cables energized over four production days (black stars in Fig. 9). Additionally, at night and during non-production periods passive electric field recordings, using the same instruments were utilized for MT processing. However, these MT results are beyond the scope of this paper.

## 3.3. Data processing

In typical ERT methods, the data are usually recorded with a small number of stacks for each transmitter, and then processed efficiently
<figure><img src="image_11.png" /><figcaption>**Figure 10.** Survey geometry at the field site. Small stars mark current injection poles; transmitter dipoles are defined by the line between the stars. Receiver stations are indicated as diamonds for ADU (Analog Digital Unit, Metronix), V-Full Waver (IRIS Instruments) and RAU (Remote Autonomous Unit, Sercel). The shaded area corresponds to the Bouguer anomaly. The location of the deep borehole shown in Fig. 8 is marked with a big star.</figcaption></figure>

in the time-domain. However, for the Deep ERT-IP data presented in this study, the recordings are significantly longer (several tens of minutes) and contain a large number of stacks. As a result, it is more effective to apply robust processing in frequency domain to compute the normalized voltages, as described for instance in Streich *et al.* (2013).

The first step involves data editing, where all measurements are assigned according to their GPS position, dipole length and the corresponding start and stop times for each transmitter. Voltage and current time-series from the various data loggers are then down sampled to a common frequency (100 Hz) and converted into a unified format. For each Tx-Rx pair, the processing algorithm selects

the appropriate time segments, which are further subdivided into overlapping subsegments of 38 periods for 0.125 Hz data. Fourier coefficients are then computed at the frequency of interest, resulting in two vectors of complex-valued Fourier coefficients for current $I(\omega)$ and voltage $V(\omega)$. A least-squares minimization algorithm is then used to estimate the complex transfer function that best fits all available time segments

$$V^{\text{Tx-Rx}}(\omega) = I^{\text{Tx-Rx}}(\omega) I^{\text{Tx}}(\omega), \qquad (5)$$

where $I^{\text{Tx}}(\omega)$ is the vector of all current coefficients for the transmitter, $V^{\text{Tx-Rx}}(\omega)$ corresponds to the synchronized voltage coefficients for the receiver at frequency $\omega$.
<figure><img src="image_12.png" /><figcaption>**Figure 11.** Geophysical instruments deployed in the field: (a–c) recording systems—(a) ADU (Metronix), (b) RAU (Sercel seismic node with BRGM adapter) and (c) V-FullWaver (IRIS Instruments); (d–f) transmitting systems—(d) TXU-30 (Phoenix), (e) GGT-3 (Zonge) and (f) TIP6000 (IRIS Instruments).</figcaption></figure>

The least-square residual error provides an estimate of the measurement uncertainty for the transfer function *T*. Since this approach focuses exclusively on the fundamental harmonic frequency of the injected current, it naturally filters out anthropogenic noise, low-frequency telluric currents and spontaneous polarization effects, eliminating the need for additional low-pass filtering. Additionally, the least-squares approach effectively removes outliers caused by transient noise in the time-series.

At low frequencies, the real part of $T^{\text{Tx-Rx}}$ corresponds to the normalized voltage-to-current ratio, which is used for resistivity imaging, while the phase of $T^{\text{Tx-Rx}}$ represents the IP response. In the absence of polarization effect, this phase should tend toward zero.

The normalized voltage obtained ranges between 2 µV A⁻¹ and 100 mV A⁻¹, with a resulting uncertainty usually below 1 µV A⁻¹. Since the noise in the current measurement is usually low, these values essentially reflect detectable voltage signals for the largest offsets, as low as 5 µV, still above the data logger's instrumental noise (approximately 0.2 µV at 0.125 Hz for all instruments). The error in estimating the ratio *T* is below 0.5 per cent for most dipoles. The injected signal is detectable on nearly all receivers, up to 17 km from the transmitters. The weakest detectable signals (around 5 µV) generally correspond to a weak component on a measurement site located at large distance from the transmitter (more than 10 km). The other component is generally stronger. The error on the phase estimated at 8 s is not well correlated with the transmitter-receiver distance, but those errors range between 0.02° for the strongest signals (~100 mV) to ±0.6° for the weakest signals (5–10 µV) for measured phase between 0 and 1.5°.

An example of data processing from one sounding, using Tx 1 and Rx 51 (marked in Fig. 10), is illustrated in Fig. 12: (a) the

time-series of both current and voltage for a Tx–Rx pair, (b) the corresponding Fourier transform and the averaged Fourier coefficients for the fundamental frequency of 0.125 Hz along with a few harmonics and (c) the final normalized voltage response, $T^{\text{Tx-Rx}}(\omega)$. The final response is presented in amplitude and phase for the entire spectrum, obtained by similar processing multiple fundamental frequencies and their harmonics ranging from 0.03125 Hz (32 s) to 800 Hz. Each colour corresponds to values of *T* extracted from different fundamental frequencies (0.03125, 0.125, 0.5, 2, 8, 32, 128, 512 Hz).

Although processing the full spectrum response is not necessary for subsequent ERT-IP inversion, it provides valuable quality control: First, the continuity of the computed response across different independently processed harmonics serves as an excellent tool for detecting noise or inconsistencies between *V* and *I*. Additionally, analysing amplitude and phase variations with frequency helps distinguish IP effects from EM induction.

Both phase and amplitude should exhibit asymptotic behaviour at low frequencies, tending towards a stable value where the dc—direct current assumption remains valid. In the example of Fig. 12, the amplitude varies slowly and slightly (by approximately 2 per cent) from 1 Hz to almost 50 Hz before increasing significantly at higher frequencies due to induction effects. The phase, which exceeds +5° at high frequencies near 200 Hz, decreases at lower frequencies. Without any IP effects, the phase should converge to 0°. In this case, the phase asymptotically stabilizes around -1° at low frequency (below 10 Hz), indicating the presence of IP effects rather than EM induction. The EM effect at 0.125 Hz can therefore be considered acceptable for dc—direct current processing.

The obtained data are of quite high quality for both voltage and phase measurements. In the time-series data, a square wave shape of the observed voltage is visible even at very long distances between
<table><caption>**Table 2.** Summary of the field survey instruments characteristics.</caption><thead><tr><th>Data logger</th><th>No. Channel</th><th>Sampling frequency [Hz]</th><th>Dynamic range (no gain) [V]</th><th>ADC resolution [bit]</th><th>Signal resolution [µV]</th><th>No. Unit</th><th>Dipole length [m]<br>North<br>East</th></tr></thead><tbody><tr><td>Metronix<br>ADU07/ADU08</td><td>5 (E + H)</td><td>4096</td><td>±1</td><td>24/24(HF), 32(LF)</td><td>470</td><td>19</td><td>50–100<br>50–100</td></tr><tr><td>Sercel RAU + BRGM<br>adaptor</td><td>2</td><td>2000</td><td>±2.5</td><td>24</td><td>0.30</td><td>77</td><td>50<br>50</td></tr><tr><td>IRIS Instrument<br>FullWave</td><td>2</td><td>100</td><td>±15</td><td>16</td><td>1</td><td>18</td><td>50<br>50</td></tr><tr><th>Transmitter Unit</th><th>Min. frequency [Hz]</th><th>Max. frequency [Hz]</th><th>Max. current [A]</th><th>Max. voltage [V]</th><th>Max. power [kW]</th><th>Transmitter<br>Generator</th><th>Weight [kg]</th></tr><tr><td>Zonge GGT3</td><td>DC</td><td>8192</td><td>20</td><td>1000</td><td>3</td><td>30<br>30</td><td>55 (ZMG-3)<br>213 (Europower)</td></tr><tr><td>IRIS Instrument</td><td>0.03125 (32 s)</td><td>2</td><td>15</td><td>3000</td><td>6</td><td>30</td><td>EP20000TE (S))</td></tr><tr><td>TIP6000</td><td>0.003 906 (256 s)</td><td>9600</td><td>40</td><td>1000</td><td>20</td><td>52.5</td><td>996 (QAS 30 KD<br>Atlas Copco)</td></tr><tr><td>Phoenix TXU-30</td><td></td><td></td><td></td><td></td><td></td><td></td><td></td></tr></tbody></table>

<caption>**Figure 12.** Example of signal processing: (a) time-series of current and voltage for a Tx-Rx pair; (b) general Fourier transform of the signal and Fourier coefficient computed for the selected fundamental frequency and its first harmonics; (c) amplitude and phase of the normalized voltage for the whole frequency spectrum. Each colour represents the values extracted from a different fundamental frequency data subset and each dot results from processing a different harmonic.</caption>

Tx and Rx (greater than 10 km), though the IP response is not always distinct in these cases. Data quality can also be assessed using harmonic frequency data after processing: signal amplitudes of the harmonic components decrease exponentially with increasing harmonic order relative to the fundamental frequency. When the full spectrum, including harmonic frequency data, is plotted, the

<figure><img src="image_13.png" /></figure>
<figure><img src="image_14.png" /><figcaption>Figure 13. Examples of processed data from transmitter Tx1 for three receiver stations: (a) Rx70 (~2.8 km from Tx1), (b) Rx51 (~4.8 km) and (c) Rx55 (~6.6 km). Each row presents, from left to right: (1) absolute voltage amplitudes across the fundamental and harmonic frequencies (at a fundamental frequency of 0.03125 Hz, corresponding to a 32 s period), (2) IP decay curves in the time domain and (3) phase responses in the frequency domain including harmonic components.</figcaption></figure>

spectrum exhibits a well-aligned trend across frequencies in both voltage and phase. The high-data quality data can be attributed to the survey being conducted in rural region, resulting in less anthropogenic electromagnetic noise, as well as the sufficiently long injecting time approximately 30 min.

The processed voltage, IP decay and phase data exhibit systematic variations with increasing transmitter-receiver distance (Fig. 13). At shorter offsets (Rx70, ~2.8 km), the voltage amplitudes are high and decrease smoothly with frequency, while the time-domain IP decay curve shows a clear exponential trend, indicating a strong and stable polarization response. The corresponding frequency-domain phase response is well defined, maintaining consistent values around -15 mrad. At intermediate offsets (Rx51, ~4.8 km), voltage amplitudes decrease slightly, and the IP decay becomes less distinct but still coherent, with the phase response remaining interpretable despite minor fluctuations at higher frequencies. At the farthest offset (Rx55, ~6.6 km), the IP decay curve becomes noisier but retains an overall decaying pattern. Importantly, the voltage amplitude at the fundamental frequency remains above 0.1 mV V⁻¹, confirming that the phase data are still reliable for interpretation.

These examples collectively illustrate the expected degradation of signal quality with increasing distance and highlight the advantage of frequency-domain processing for weak-field conditions. While time-domain IP responses were also analysed, frequency-domain interpretation proved more robust against noise contamination, particularly for larger offsets, allowing reliable recovery of phase information even when the time-domain decay was obscured.

Nodal systems usually measure the voltage data in two orthogonal directions (such as Fullwaver system of IRIS; Lajaunie et al. 2019). Occasionally, some data points may present singular low-voltage values (lower than noise level) due to an orthogonal direction between the source dipole and the measuring dipole. Such a

configuration is associated with a singular geometrical factor. To address this issue, we adjusted combining two vectors with two voltage components ($r_a$ and $r_b$) into a single vector ($r_c$) (Loke et al. 2022) as follows

$$r_c = [r_a \sin(\alpha) + r_b \sin(\beta)]i + [r_a \cos(\alpha) + r_b \cos(\beta)]j, \quad (6)$$

where $\alpha$ and $\beta$ are azimuths in the x and y directions, i and j are unit vectors in the x and y directions. The amplitude of the combined vector is then given by

$$ \begin{aligned} r_c &= [(r_a \sin(\alpha) + r_b \sin(\beta))^2 + (r_a \cos(\alpha) + r_b \cos(\beta))^2]^{\frac{1}{2}} \\ &= [r_a^2 + r_b^2 + 2r_a r_b \cos(\beta - \alpha)]^{\frac{1}{2}}. \end{aligned} \quad (7) $$

The combined vector data for the field data is presented in Fig. 14. For systems using three electrodes, the centre electrode remains fixed, and the distance of the new defined electrode can be calculated using eq. (6). On the other hand, the centre electrode is redefined as the midpoint of two dipoles along the x and y-axis as systems using four electrodes. The new electrode position is then calculated by combining the vectors of the two dipoles.

Following standard data processing steps, data quality was generally satisfactory; however, some outliers persisted, primarily due to the weak signals resulting from the large distances between Tx and Rx. In the final stage of data processing, we removed bad data points, including voltage values below 0.1 mV A⁻¹ and outlier phase data less than -100 mrad. For the 3-D inversion, we used only the data corresponding to the 0.125 Hz source frequency, resulting a final data set of 428 filtered data points.

The apparent resistivity values mostly ranged between 10³ and 10⁴ Ω.m, indicating a predominantly resistive subsurface at the field site (Fig. 15a). The 3-D pseudo-section of apparent resistivity (Fig. 15b) highlighted a highly resistive zone (>10⁵ Ω.m) near the surface, particularly on the eastern side, while more conductive values (around 10³ Ω.m) appeared at greater depths. The phase
<figure><img src="image_15.png" /><figcaption>**Figure 14.** Examples of calculating combined data in a data processing for the field data, showing original E-field vectors and vector sum E-field vectors. Adjustments to electrode positions from their original positions address low-voltage issues caused by dipole alignment. (a) Electric field distribution at each voltage dipole when injecting current at Tx 1 and (b) zoomed-in view of a specific area for detailed changes.</figcaption></figure>

<figure><img src="image_16.png" /><figcaption>**Figure 15.** (a) Data histogram and (b) 3-D pseudo-sections of apparent resistivity; (c) data histogram and (d) 3-D pseudo-sections of apparent phase.</figcaption></figure>
<figure><img src="image_17.png" /><figcaption>**Figure 16.** 3-D inversion results of the resistivity model, showing (a) a cross-section of the 3-D volume, (b) depth slices at various depths (500, 1000, 1500 and 2000 m) with comparing gravity data and (c) comparisons between observed data, predicted data from the inversed model, and relative errors.</figcaption></figure>

data showed amplitudes concentrated between 0 and 80 mrad, with a peak around −10 mrad (Fig. 15c). In the 3-D pseudo-section of phase (Fig. 15d), low IP responses (around −7 mrad) were observed at shallow depths. However, at depths exceeding 2 km, a distinct anomalous region emerged, exhibiting a maximum absolute phase value exceeding 40 mrad.

## 4. RESULTS AND DISCUSSION

We inverted the ERT and IP data using the open-source code, py-GIMLi (Rücker et al. 2017), which minimizes the cost function using a Gauss–Newton scheme to obtain a 3-D resistivity and phase model. An unstructured tetrahedral mesh was used, consisting of 206 915 nodes and 986 080 cells. To focus on the variations in the
<figure><img src="image_18.png" /><figcaption>**Figure 17.** Inverted 3-D resistivity models depending on the initial model value ($m_0$), (a) 3 k Ωm, (b) 5 k Ωm and (c) 10 k Ωm, to check the data dependence.</figcaption></figure>

model, we applied L1 norm constraints to ERT and IP inversions for both data and model, with a lambda value of 5 for the smoothness constraints. The pyGIMLi inversion algorithm first inverts a resistivity model and then an IP model based on the recovered resistivity model. The total number of iterations for the 3-D ERT inversion was 8, resulting in a final root-mean-square (RMS) error 3.4 per cent. For the 3-D phase inversion, the total number of iterations was 6, with a final RMS error 5.7 per cent.

The recovered 3-D resistivity model (Fig. 16a) defines a resistivity range approximately from 97.8 Ωm to more than 600 000 Ωm. To ensure visualization of reliable structures, regions with sensitivity values lower than 1e-9 were excluded. Surface resistivity distribution appear heterogeneous due to the swampy conditions at the study site. Below approximately 500 m depth, a highly resistive layer (>30 000 Ωm) extends to about 1500 m. A relatively conductive feature (<5000 Ωm) is clearly visible near the centre of the survey profiles. Depth slices (Fig. 16b) illustrate that this conductive anomaly emerges at approximate 1000 m depth, extending along the primary profile orientation. A comparison between observed and modelled resistivity data demonstrates close agreement (Fig. 16c).

To evaluate the stability and robustness of the 3-D inversion results, we performed ERT inversions using three different homogeneous initial models with resistivity values of 3, 5 and 10 kΩ·m (Fig. 17). The resulting resistivity distributions demonstrate a high degree of consistency in the recovered structures, indicating that the inversion outcomes are not strongly dependent on the choice of starting model. In all three cases, the main conductive zones near the central part of the model and the resistive features towards

the boundaries are well reproduced, with only minor variations in absolute resistivity values. This consistency suggests that the inversion converged toward a stable solution primarily driven by the data rather than the initial conditions. Therefore, the recovered model can be considered reliable for interpreting the subsurface resistivity structure.

The inverted phase model (Figs 18a and b) reveals negligible IP responses from the surface down to approximately 1000 m depth. Below this depth, a phase anomaly characterized by relatively high magnitudes (~ -30 mrad) emerges distinctly. This anomaly spatially coincides with the conductive feature identified in the resistivity model (Fig. 16b). Although predicted IP data indicate weaker magnitudes at greater depths compared to observed data, overall trends align well, with small absolute errors (Fig. 18c).

The near-surface of the study area is characterized by strong heterogeneity due to the presence of swampy zones that locally reduce resistivity values. Although these conductive surface layers are shallow, their influence on the inversion results is non-negligible, as they can distort the electric field distribution and potentially affect the recovered deep resistivity values. However, the use of a true 3-D acquisition geometry mitigates this limitation by providing lateral illumination and cross-line sensitivity that are absent in 2-D surveys. Even when transmitter and receiver spacing is relatively sparse, 3-D coverage enables the inversion to integrate multidirectional constraints, allowing localized near-surface heterogeneities to be better captured and their effects on deeper structures minimized. Quantitatively assessing the gain in shallow resolution achieved by 3-D coverage compared to 2-D remains an open research question in ERT, and, as in seismic studies, it deserves dedicated future
<figure><img src="image_19.png" /><figcaption>**Figure 18.** 3-D inversion results of the phase model, showing (a) a cross-section of the 3-D volume, (b) depth slices at various depths (500, 1000, 1500 and 2000 m) with comparing gravity data and (c) comparisons between observed data, predicted data from the inversed model, and absolute errors.</figcaption></figure>

investigation. Nonetheless, this study demonstrates that 3-D acquisition plays a crucial role in improving the reliability of deep imaging in geologically heterogeneous terrains.

Comparisons between gravity and ERT-IP plane sections indicate that the linear Bouguer gravity anomaly, as resulting from an ultra-mafic intrusion (Tirroniemi et al. 2024), spatially aligns well with

ERT (Fig. 16b) and IP (Fig. 18b) anomalies at depths exceeding 1500 m. These correlated anomalies delineate a conductive zone with elevated IP responses. The phase model shows broader IP anomalous region compared to the conductive anomaly, reflecting lower IP data resolution. While phase anomalies exhibit modest contrast with background values, it shows good coherence with
conductive and bouguer anomaly structures. Therefore, we can expect recovered IP value could be underestimated during inversion but still enough to recognize the IP anomaly structure more than 1.5 km depths.

The drill core visual log indicates that the mafic section does not contain massive sulfide mineralization but includes trace amounts of Cu-Ni-sulfides. Previous studies of the Näränkävaara layered intrusion reported interstitial pyrrhotite-dominated sulfides (approximately 3 vol per cent) within olivine- pyroxenite units at depths near 1750 m (Järvinen et al. 2019). Given the geological continuity between this study area and the Näränkävaara intrusion, it is plausible that alteration processes between pyroxenite and peridotite—predominantly olivine-bearing rocks—enhance chargeability, consistent with the strong IP effects observed in borehole logs (~50 per cent of PFE in Fig. 9). The apparent discrepancy between the 3-D inverted model and borehole IP logs can be explained by differences in data resolution and sampling scale. Borehole measurements directly probe the rocks adjacent to the well, capturing fine-scale heterogeneities. In contrast, surface-based 3-D inversions integrate signals over larger volumes, resulting in smoother resistivity and phase distributions.

Our analysis underscores the challenges of imaging deep structures using ERT and IP, particularly in maintaining sufficient SNR and sensitivity. As target depth increases, achieving a detectable response requires higher transmitter currents, longer dipole lengths and extended acquisition times, which in turn limit the number of feasible transmitter positions and affect overall survey efficiency.

A similar trade-off exists for receiver positioning: deeper targets necessitate a larger array aperture, leading to sparser spatial sampling for a fixed set of instruments. While reduced spatial resolution is a drawback, it remains acceptable if the target is still identifiable. Therefore, survey design must strike a balance between cost, resolution and acquisition efficiency.

As demonstrated in our synthetic tests (Section 2.3), smaller electrode spacing improves imaging resolution when paired with an appropriate depth of investigation. However, this approach significantly increases survey duration and logistical challenges, particularly when using conventional cabled systems. For large-scale deep surveys, nodal systems offer a practical alternative, substantially reducing acquisition time and cable requirements while maintaining adequate resolution.

## 5. CONCLUSION

This study successfully demonstrates the capability of 3-D Deep ERT and IP to image geological structures at depths exceeding 1.5 km. The field survey, conducted at the Koillismaa Linear Intrusion Complex in Finland, resulted in reliable resistivity and IP images clearly delineating conductive and chargeable anomalies. These anomalies align well with independent gravity and borehole logging data and are interpreted to reflect mafic-ultramafic intrusion structures. Significantly, this work represents the first documented case of successfully acquiring IP data in 3-D at these substantial depths.

To achieve these results, we first conducted we thoroughly analysed background electrical noise levels from field measurements, enabling precise quantitative evaluation of the expected SNR and sensitivity. Additionally, comprehensive synthetic tests to quantitatively assess the feasibility of imaging very deep targets using the unconventional configuration, DFA, utilizing nodal data acquisition

systems. These preparatory procedures provided essential decision-making tools for optimized survey design, significantly enhancing survey efficiency and data reliability, and ultimately enabling the successful imaging of deep geological structures.

The use of nodal systems significantly streamlined the large-scale field operations by providing flexibility in array geometry and eliminating the need for excessively long cables. To ensure data quality at very deep targets, we implemented extensive current injection times, reaching up to approximately 30 min per transmitter position at the lowest frequency. Furthermore, processing full time-series data in the frequency domain was essential in identifying and effectively removing inconsistent noise, distinguishing IP effects from EM induction.

This study emphasizes that Deep ERT-IP surveys require careful and quantitative planning to overcome both technical and logistical challenges. The systematic workflow described herein—incorporating rigorous noise-level assessments, optimal transmitter configurations, extended data acquisition times and advanced logistics via nodal instrumentation—enabled acquisition of physically meaningful geoelectrical images. This comprehensive approach is broadly applicable not only in mineral exploration but also for other deep-depth investigations, such as hydrogeological and geothermal studies. Future work will aim at integrating these geoelectrical models with complementary multidisciplinary geophysical data sets and analysing spectral IP data for improved geological interpretation.

## SEEMSDEEP WORKING GROUP

Jochen Kamm, Uula Autio, Suvi Heinonen, Michał Malinowski, Cedric Patzer, Tuomo Karinen, Veera Pajunen, Toni Eerola, Tuija Luhta, Viveka Laakso, Olli Salmela (Geological Survey of Finland, GTK); Mathieu Darnet, Francois Bretaudieu, Simon Vedrine, Bitnarae Kim, Florent Beaubois, Jacques Deparis (French Geological Survey, BRGM); Andrzej Górszczyk, Brij Singh, Szymon Oryński (Institute of Geophysics, Polish Academy of Sciences); Thomas Kalscheuer, Karol Sierszen, Mehrdad Bastani (Uppsala University); Catherine Truffert, Fabrice Vermeersch, Julien Gance (IRIS Instruments); Marek Wojdyła, Łukasz Sito (Geopartner); Antti Kivinen, Janne Kaukolinna (GRM Services); Yuriy Koltun (Laakso Minerals).

## DATA AVAILABILITY

The resistivity and IP data, used for the 3-D inversion in the study, are available with restricted access at [Deep 3D ERT-IP Method from Koillismaa Linear Intrusion Complex in Finland], accessible via https://doi.org/10.18144/41b89187-6a52-43b1-b870-947627e4b5b3. The modelling and inversion were conducted version 1.4.3 of the open-source pyGIMLi software, which is available at https://www.pygimli.org/.

## ACKNOWLEDGMENTS

This work has been conducted as part of the SEEMS DEEP (Seismic and electromagnetic methods for deep mineral exploration) project that was funded through ERA-MIN3 Joint Call 2021 (National funding agencies: Finland: Business Finland, Grant ID: 640/31/2022; France: ANR, Grant ID: ANR-22-MIN3-0006-02). We want to further acknowledge FLEX-EPOS and Geophysical Instrument Pool Potsdam (GIPP) for providing seismic and electromagnetic instruments.
REFERENCES

Alfano, L., 1984. Geoelectrical deep explorations by means of direct cur-
rents, *Problems and Methods for Lithospheric Exploration*, pp. 169–177,
ed., Cassinis, R., Springer-Verlag, Berlin.

Asch, T. & Morrison, H.F., 1989. Mapping and monitoring electrical resistivity with surface and subsurface electrode arrays, *GEOPHYSICS*, **54**, 235–244.

Balasco, M., Giocoli, A., Piscitelli, S., Romano, G., Siniscalchi, A., Stabile, T.A. & Tripaldi, S., 2015. Magnetotelluric investigation in the High Agri Valley (southern Apennine, Italy), *Nat. Hazards Earth Syst. Sci.*, **15**, 843–852.

Balasco, M., Lapenna, V., Rizzo, E. & Telesca, L., 2022. Deep electrical resistivity tomography for geophysical investigations: the state of the art and future directions, *Geosciences*, **12**, 438.

Barker, R.D., 1989. Depth of investigation of collinear symmetrical four-electrode arrays, *GEOPHYSICS*, **54**, 1031–1037.

Bretaudeau, F., Dubois, F., Bissavetsy Kassa, S.G., Coppo, N., Wawrzyniak, P. & Darnet, M., 2021. Time-lapse resistivity imaging: cSEM-data 3-D double-difference inversion and application to the Reykjanes geothermal field, *Geophys. J. Int.*, **226**, 1764–1782.

Carrier, A., Fischanger, F., Gance, J., Cocchiararo, G., Morelli, G. & Lupi, M., 2019. Deep electrical resistivity tomography for the prospection of low-to medium-enthalpy geothermal resources, *Geophys. J. Int.*, **219**, 2056–2072.

Cole, K.S. & Cole, R.H., 1941. Dispersion and absorption in dielectrics I. Alternating current characteristics, *J. Chem. Phys.*, **9**, 341–351.

Darnet, M., Wawrzyniak, P., Coppo, N., Nielsson, S., Schill, E. & Fridleif-
sson, G., 2020a. Monitoring geothermal reservoir developments with the
controlled-source electro-magnetic method—a calibration study on the
Reykjanes geothermal field, *J. Volc. Geotherm. Res.*, **391**, 106–437.

Darnet, M., Wawrzyniak, P., Tarits, P., Hautot, S. & d'Eu, J.-F., 2020b. Map-
ping the geometry of volcanic systems with magnetotelluric soundings:
results from a land and marine magnetotelluric survey performed during
the 2018–2019 Mayotte seismovolcanic crisis, *J. Volc. Geotherm. Res.*
**406**, 107–046.

Dean, T., 2024. A comparison of two land nodal seismic recording systems, *Explor. Geophys.*, **55**, 277–287.

Edwards, L.S., 1977. A modified pseudosection for resistivity and IP, *GEOPHYSICS*, **42**, 1020–1036.

Everett, M.E., 2013. *Near-surface Applied Geophysics*, Cambridge University Press.

Gasperikova, E. & Morrison, H.F., 2001. Mapping of induced polarization using natural fields, *Geophysics*, **66**, 137–147.

Gasperikova, E., Cuevas, N.H. & Frank Morrison, H., 2005. Natural field in-
duced polarization for mapping of deep mineral deposits: a field example
from Arizona, *GEOPHYSICS*, **70**, B61–B66.

Hohmann, G.W., Kintzinger, P.R., Van Voorhis, G.D. & Ward, S.H., 1970. Evaluation of the measurement of induced electrical polarization with an inductive system, *GEOPHYSICS*, **35**, 901–915.

Järvinen, Ville., Halkoaho, Tapio., Konnunaho, Jukka., Heinonen, Jussi S. & Rämö, O. Tapani., 2019. Parental magma, magmatic stratigraphy, and reef-type PGE enrichment of the 2.44-Ga mafic-ultramafic Näränkävaara layered intrusion, Northern Finland, *Mineralium Deposita*, **55** (8), 1535.

Jones, D.H., Bell, B.S. & Hansen, J.H., 1982. Induced polarization survey of sulfide bearing rocks in eastern Tennessee and western North Carolina, *Transp. Res. Rec.*, **892**, 13–19.

Karinen, H.M., Peltonen, J.E., Kähönen, M. & Tikkanen, H.O., 2010. Prediction of acute mountain sickness by monitoring arterial oxygen saturation during ascent, *High Alt. Med. Biol.*, **11**, 325–332.

Karinen, T., 2010. The Koillismaa intrusion, northeastern Finland: evidence for PGE reef forming processes in the layered series, *Geological Survey of Finland, Bulletin*, **404**.

Karinen, T., Heinonen, S., Konnunaho, J., Salmirinne, H., Lahti, I. & Salo, A., 2021. Koillismaa Deep Hole-Solving the mystery of a geophysical anomaly, *Lithosphere 2021, Eleventh Symposium on Structure, Composition and Evolution of the Lithosphere, Programme and Extended Abstracts.* pp. 55–58.

Kemna, A. et al. 2012. An overview of the spectral induced polarization method for near-surface applications, *Surf. Geophys.*, **10**, 453–468.

Key, K., 2009. 1D inversion of multicomponent, multifrequency marine CSEM data: Methodology and synthetic studies for resolving thin resistive layers, *Geophysics*, **74** (2), F9.

Kratzer, T. & Macnae, J.C., 2012. Induced polarization in airborne EM, *GEOPHYSICS*, **77**, E317–E327.

Lajaunie, M., Gance, J., Nevers, P., Malet, J.P., Bertrand, C., Garin, T. & Fer-
hat, G., 2019. Structure of the Séchilienne unstable slope from large-scale
three-dimensional electrical tomography using a Resistivity Distributed
Automated System (R-DAS), *Geophys. J. Int.*, **219**, 129–147.

Liu, W., Chen, R., Cai, H. & Luo, W., 2016. Robust statistical methods for impulse noise suppressing of spread spectrum induced polarization data, with application to a mine site, Gansu province, *China J. Appl. Geophys.*, **135**, 397–407.

Liu, W., Chen, R., Cai, H., Luo, W. & Revil, A., 2017. Correlation analysis for spread-spectrum induced-polarization signal processing in electro-magnetically noisy environments, *GEOPHYSICS*, **82**, E243–E256.

Loke, M.H., Wilkinson, P.B., Gance, J., Malet, J., Truffert, C. & Leite, O., 2022. Measurement and inversion strategies for 3-D resistivity surveys with vector arrays, *Geophys. Prospect*, **70**, 578–592.

Macnae, J., 2016. Quantitative estimation of intrinsic induced polarization and superparamagnetic parameters from airborne electromagnetic data, *Geophysics*, **81**, E433–E446.

Marchant, D., Haber, E. & Oldenburg, D.W., 2013. Inductive source induced polarization, *Geophys. J. Int.*, **192**, 602–612.

Mazzini, A., Carrier, A., Sciarra, A., Fischanger, F., Winarto-Putro, A. & Lupi, M., 2021. 3D deep electrical resistivity tomography of the Lusi Eruption site in East Java, *Geophys. Res. Lett.*, **48**, e2021GL092632.

Pelton, W.H., Ward, S.H., Hallof, P.G., Sill, W.R. & Nelson, P.H., 1978. Mineral Discrimination and removal of inductive coupling with multifrequency IP, *GEOPHYSICS*, **43**, 588–609.

Porté, J., Bretaudieu, F. & Girard, J.F., 2023. 3-D complex resistivity imaging using controlled source electromagnetic data: a multistage procedure using a second order polynomial parametrization, *Geophys. J. Int.*, **233**, 839–860.

Pucci, S. et al. 2016. Deep electrical resistivity tomography along the tec-
tonically active Middle Aterno Valley (2009 L’Aquila earthquake area,
central Italy), *Geophys. J. Int.*, **207**, 967–982.

Reninger, P.-A., Martelet, G., Perrin, J. & Dumont, M., 2020. Processing methodology for regional AEM surveys and local implications, *Explor. Geophys.*, **51**, 143–154.

Revil, A., Vaudelet, P., Su, Z. & Chen, R., 2022. Induced polarization as a tool to assess mineral deposits: a review, *Minerals*, **12**, 571.

Rizzo, E., Giampaolo, V., Capozzoli, L., De Martino, G., Romano, G., Santilano, A. & Manzella, A., 2022. 3D deep geoelectrical exploration in the Larderello geothermal sites (Italy), *Physics of the Earth and Planetary Interiors*, **329-330**, 106906.

Rücker, C., Günther, T. & Wagner, F.M., 2017. pyGIMLi: an open-source library for modelling and inversion in geophysics, *Comput. Geosci.*, **109**, 106–123.

Salmirinne, H. & Iljina, M., 2003. Koillismaan kerrosintrusiokompleksin tulokanavamuodostuman painovoimatulkinta ja alueen malmimahdollisuutet (osa 1), *Geol. Tutkimuskeskus Q*, **21**, 203–212.

Sean, 2019. *Quick Tip: Depth of Investigation for ERI Surveys* [WWW Document], AGIUSA. URLhttps://www.agiusa.com/blog/quick-tip-dept h-investigation-eri-surveys( accessed 4.21.25).

Seigel, H., Nabighian, M., Parasnis, D.S. & Vozoff, K., 2007. The early history of the induced polarization method, *Leading Edge*, **26**, 312–321.

Sheriff, R.E., 2002. *Encyclopedic Dictionary of Applied Geophysics*, 4th edn., Society of Exploration Geophysicists.

Spichak, V. & Manzella, A., 2009. Electromagnetic sounding of geothermal zones, *J. Appl. Geophys.*, **68**, 459–478.

Storz, S.J., 2000. Electrical resistivity tomography to investigate geol-
logical structures of the earth's upper crust, *Geophys. Prospect.*, **48**,
455–471.
Streich, R., Becken, M. & Ritter, O., 2013. Robust processing of noisy land-based controlled-source electromagnetic data, *GEOPHYSICS*, **78**, E237–E247.

Tamburriello, G., Balasco, M., Rizzo, E., Harabaglia, P., Lapenna, V. & Siniscalchi, A., 2008. Deep electrical resistivity tomography and geothermal analysis of Bradano foredeep deposits in Venosa area (Southern Italy): preliminary results, *Ann. Geophys.*, **51**.

Tirronniemi, J. et al. 2024. Koillismaa deep hole-final report, GTK Open File Work Rep., 73pp.

Troiano, A., Isaia, R., Di Giuseppe, M.G., Tramparulo, F.D. & Vitale, S., 2019. Deep electrical resistivity tomography for a 3D picture of the most active sector of Campi Flegrei caldera, *Sci. Rep.*, **9**, 15124.

Viezzoli, A. & Manca, G., 2020. On airborne IP effects in standard AEM systems: tightening model space with data space, *Explor. Geophys.*, **51**, 155–169.

Vozoff, K., 1990. Magnetotellurics: principles and practice, *J. Earth Syst. Sci.*, **99**, 441–471.

Wilkinson, P.B., Meldrum, P.I., Kuras, O., Chambers, J.E., Holyoake, S.J. & Ogilvy, R.D., 2010. High-resolution electrical resistivity tomography monitoring of a tracer test in a confined aquifer, *J. Appl. Geophys.*, **70**, 268–276.

Wynn, J., 2002. Evaluating groundwater in arid lands using airborne magnetic/EM methods: an example in the southwestern US and northern Mexico, *Leading Edge*, **21**, 62–64.