Evaluation of spectral induced polarization field measurements in time and frequency domain

Abstract Spectral induced polarization (SIP) measurements have been demonstrated to correlate with important parameters in hydrogeological and environmental investigations. Although SIP measurements were often collected in the frequency domain (FDIP), recent developments have demonstrated the capabilities to solve for the frequency-dependence of the complex conductivity through measurements collected in the time domain (TDIP). Therefore, the aim of our field investigations is a comparison of the measured frequency-dependence at a broad frequency range resolved through FDIP and TDIP. In contrast to previous studies, we conducted measurements with different instruments and measuring technologies for both FDIP and TDIP. This allows for investigating the robustness of different measurements and assessing various sources of errors, for the assessment of the advantages and drawbacks from different measuring techniques. Our results demonstrate that data collected through different instruments are consistent. Apparent resistivity measurements as well as the inversion results revealed quantitatively the same values for all instruments. The measurements of the IP effect are also comparable, particularly FDIP readings in the low frequencies ( 0.1 s). However, data quality for higher frequencies in FDIP (i.e., early times in TDIP) show larger variations, which reflects the differences between the instruments to deal with the electromagnetic contamination of the IP data. Concluded in general, the different instruments and measuring techniques can provide consistent responses for varying signal-to-noise ratio and measuring configurations.


Introduction
The method of the spectral induced polarization (SIP) gathers the low frequency electrical properties of the subsurface and can be measured in the time domain (TDIP) and frequency domain (FDIP). SIP has emerged in the last decade as an important technique for hydrogeological and environmental investigations (e.g., Kemna et al., 2012). Based on strong polarization effects observed in electronic conductors (e.g., Wong, 1979), the IP method was initially developed for the prospection of mineral deposits like sulphides (e.g., Seigel et al., 2007 and references therein). Developments in the instrumentation over the last decades permit detection of weaker signatures associated with other materials, such as soils (Vanhala and Soininen, 1995) or hydrocarbon contaminants (e.g., Vanhala et al., 1992). Further laboratory investigations in biofilms (e.g., Davis et al., 2006), and other biological materials such as roots (e.g., Weigand and Kemna, 2017) and wood (e.g., Martin, 2012) have paved the way for the application of FDIP to characterise biological processes in the emerging discipline of biogeophysics (e.g., Atekwana and Slater, 2009). Another benefit of the method is the strong correlation between SIP parameters and the hydraulic conductivity k observed in laboratory investigations (e.g., Binley et al., 2005;Zisser et al., 2010). Particularly for samples with high clay content, the FDIP method has demonstrated an improved estimation of k in comparison with previous approaches based only on measurements of the bulk conductivity (e.g., Slater, 2007, and references therein]. Although laboratory studies have promised a broad range of processes and materials to be characterized through FDIP, to date, fieldscale imaging investigations are still rare. Hördt et al. (2007) reported the first FDIP imaging at the field scale to characterise hydraulic conductivity. Yet, they neglect the frequency dependence of the complex conductivity and assume the so-called constant phase model. Based on the comparison of data at two different frequencies (1 and 10 Hz), Williams et al. (2009) demonstrated the applicability of the FDIP method to monitor processes accompanying the bioremediation of a uranium contaminant plume. Using a broader frequency bandwidth, Flores Orozco et al. (2011) demonstrated the sensitivity of the method to characterise changes in the redox conditions of the subsurface due to microbial activity. Also the accumulation of metallic minerals could be observed (Flores Orozco et al., 2013). Further applications in biogeophysics revealed the possibility to detect fungi infections in oak trees (Martin and Günther, 2013). Regarding site characterization, imaging data collected between 0.06 Hz and 1000 Hz in a hydrocarbonimpacted site (Flores Orozco et al., 2012a) revealed the possibility to use FDIP images to delineate the geometry of the source-zone and the plume of dissolved contaminants. In a similar frequency bandwidth the SIP method, together with inversion strategies using frequencyregularization schemes, permitted the characterization of hydraulic conductivity (Kemna et al., 2012).
SIP measurements have been traditionally performed in the frequency domain. In recent years, time domain measurements have also emerged as a suitable technique to gain information about the frequency-dependence (e.g., Fiandaca et al., 2013). The derivation of SIP parameters based on TDIP measurements has been demonstrated for the estimation of hydraulic parameters (Maurya et al., 2018a) and the evaluation of contaminated sites (Johansson et al., 2014(Johansson et al., , 2015. Furthermore, the robustness of the existing TDIP instruments has permitted to open new areas of research, such as the monitoring of CO2 injections (e.g., Doetsch et al., 2015).
The comparison of IP measurements collected in time and frequency domain was a topic of research during an early development of the method (e.g., Zonge et al., 1972;Voorhis et al., 1973). Zonge et al. (1972) demonstrated similar but not identical results between time and frequency and phase measurements due to the inherent differences in the respective parameter definition. Recent developments in IP data collection, data processing and modeling open the door to further investigations aiming at defining potential advantages and drawbacks in TDIP and FDIP surveys. Based on the comparison of imaging results, Flores Orozco et al. (2012a) demonstrated that TDIP, single-and multi-frequency FDIP provide consistent results regarding the quantification of hydrocarbon contaminants in groundwater. Yet, that study used only integral chargeability measurements from TDIP and did not record the full waveform for the estimation of SIP parameters. A recent study by Maurya et al. (2018b) demonstrates the advantages of the TDIP measurements to provide SIP information in a broader frequency bandwidth with a reduced acquisition time and increased depth penetration compared to their FDIP measurements. Yet, most of the FDIP observations from that study are not consistent with previous results observed at the field scale, e.g., the removal of large electrode separation (deep) measurements or higher frequencies (> 10 Hz) (Flores Orozco et al., 2012a;. Such apparent inconsistency asks for further research to better understand the actual limitations of the FDIP method, and whether those are associated to the measuring instruments, the processing of the data or to the measurement technology itself. The objective of this study is the investigation of the quality of the IP readings and their frequency bandwidth using two different instruments for both TDIP and FDIP, thus extending previous studies. We used four, for us available and well-known, commercial instruments and techniques to permit the collection of dipole-dipole arrays with the same settings as possible to avoid favouring particular devices. The measurements presented here were conducted primarily at a profile with short electrode spacing (a = 1 m) in an area characterized by a shallow polarizable geological unit, which permitted readings with high signal-to-noise ratio (S/N). We also discus results from a profile with larger electrode spacing (a = 5 m) and lower S/N for deeper investigations. After introducing the IP methods and the section about the settings used in the geological background, we investigate the spatial consistency and the frequency-dependence observed. After presenting the results from both profiles, we discuss the data evaluation in detail, also considering the signal strength. Furthermore, notes about optimal field procedures follows the comparison about acquisition time and field effort of the different instruments.

The induced polarization (IP) method
The IP technique is an extension of the direct-current (DC) resistivity method using four-electrode arrays, i.e., two electrodes are used to inject a current and a second pair of electrodes measure the resulting voltages. In FDIP, a phase-shifted voltage relative to an alternating current is measured, with the recorded values given in terms of the magnitude (i.e., transfer resistance given by the voltage-to-current ratio) and phase shift (between current and voltage) of the electrical impedance. The injection of a current at different frequencies (commonly below 1 kHz) results in an impedance spectrum (e.g., Ward, 1990). In contrast, most of the commercially available instruments work in the time domain, which is based on injections of a square wave. Accordingly, the voltage-to-current ratio during the current injection provides a measure of the transfer resistances. Readings of the voltage decay after switching off the current is used to quantify the IP effect.
Inversion of TDIP measurements have until recently generally been carried out by inverting for resistivity and integral chargeability, where often the resistivity inversion is done first and used as a base for the chargeability inversion (Oldenburg and Li, 1994). Another approach is to treat it as complex conductivity data with real and imaginary components . Not accounting for the IP adequately in the data acquisition and inversion does not only affect the IP model, but can also lead to under-estimation of the DC resistivity as the measured potentials may have not reached the full magnitude when integrating the signal while there is still an IP build-up in progress (e.g., Olsson et al., 2019). If TDIP data are acquired in a sufficiently broad range of time windows (IP decay integration intervals), reflecting a corresponding spectral range, it is also possible to invert for spectral properties using e.g. a Cole-Cole model (Fiandaca et al., 2013) or maximum phase angle model . It is however challenging to acquire TDIP data with a sufficiently broad range of time windows to allow meaningful spectral inversion in practical applications, and it is only quite recently that instrumentation with suitable properties and matching signal processing algorithms have become available (Olsson et al., 2016).
Inversion of IP measurements solve for the distribution of the complex conductivity (σ * = σ′ − iσ′′), or its inverse, the complex resistivity (ρ * = ρ′ + iρ′′) with the imaginary unit i ¼ ffiffiffiffiffiffiffi ffi −1 p . The real component of the complex conductivity (σ′) represents the electrical conduction (i.e., energy loss) properties of the materials, which -in case of negligible metallic conductors -is controlled by ionic conduction (i.e., dependent on porosity, saturation and fluid conductivity) and by surface conduction processes taking place at the grain-fluid interface, which may dominate in case of fine grained materials. The imaginary component (σ′′) represents the capacitive (i.e., energy storage) properties of the materials. For further details we refer to Ward (1990) and Kemna et al. (2012) and references therein. Moreover, FDIP measurements also provide information about the frequency dependence of the complex conductivity σ * (ω), with ω being the angular frequency. Modern FDIP instruments are able to record data in a broad frequency bandwidth (e.g., Kemna et al., 2012) deploying tens to hundreds of electrodes for the collection of thousands of measurements in SIP imaging (i.e., tomographic) data sets.
Regardless of whether the measurements are collected in time or frequency domain, the acquisition of SIP field data in a broad frequency/time bandwidth is particularly challenging due to particular sources of errors such as polarization of electrodes and electromagnetic (EM) coupling, which may affect the SIP data quality. Polarization of the electrodes used previously for current injection can result in distortions of the voltage readings, which especially affect measurements at low frequencies (< 10 Hz) (e.g., Dahlin, 2000;LaBrecque and Daily, 2008). On the other hand, measurements at high frequencies or early times (FDIP: > 100 Hz, TDIP: > 10 ms) can be affected by unwanted EM inductive or capacitive effects. The former arises due to the self-induction of EM fields in the near surface during the flow of the current along the cables connecting the electrodes and the instrument, and it is proportional to the conductivity of the ground, the frequency and the square of the cable length (i.e., Hallof, 1974). Capacitive EM effects originates from parasitic currents caused by potential differences between the electrodes and between the electrodes and the subsurface (Zimmermann et al., 2008;Flores Orozco et al., 2013). Another source of capacitive coupling is the inductance between the cables (or the electronics) used for current injection and potential readings (e.g., Zimmermann et al., 2008;Zhao et al., 2013Zhao et al., , 2015. Different approaches have been proposed to minimize EM contamination of SIP measurements, such as the deployment of separated cables for current injections and potential readings (Dahlin and Leroux, 2012), the use of a single treat of shielded cables (e.g., Flores Orozco et al., 2013), or the application of processing techniques to remove the EM effects (e.g., Pelton et al., 1978;Kemna, 2000;Fiandaca, 2018).

Instruments
In this study, we collected SIP measurements deploying four different instruments, two for FDIP and two for TDIP. In case of FDIP measurements, we deployed a Radic Research instrument SIP256C and a Data-Acquisition System 1 (DAS-1). For TDIP we used an ABEM Terrameter LS 2 and a Syscal Pro Switch 72. A picture of the used instruments can be seen in Fig. 1.
The SIP256C instrument is a multi-channel instrument measuring the electrical impedance (amplitude and phase) over 6 decades of frequencies (between 1 mHz and 1000 Hz, sinusoidal). The current injection is done by a cable which interconnects the electrodes. To avoid capacitive cable effects, the SIP256C uses multiple receivers (RU -remote unit) that digitize the measured potential difference at each electrode and transfer these data via fibre-optic cables. Each RU can be used for current injection, whereas potential differences can be measured between two adjacent (or further distant) electrodes. The SIP256C has an internal 50 W high-voltage source capable of transmitting currents up to 250 mA using maximum voltage of 400 V. For this survey 24 RUs were available, therefore it was necessary to use a roll-along scheme for the shown profiles. The cable length between the electrodes can be adjusted to 1 m or 10 m.
The DAS-1 is an eight-channel acquisition system designed for direct current (DC) resistivity, being able to work in both TDIP and FDIP. The DAS-1 can acquire the frequency domain in two modes: the FDIP mode, where the system is performing a real-time Fourier Transform of data values collected at even intervals performed during the current injection and in the spectral IP mode where data can be collected along 17 selected frequencies between 0.015 and 225 Hz. In contrast to the SIP256C the DAS-1 use only squared waveform. We used the spectral IP mode and connected 64 electrodes to the instrument using a multicore cable. The DAS-1 has a 250 W built-in transmitter capable of transmitting a constant current (max. 2.5 A) or a constant voltage (max. 475 V). In our measurements we requested injections to the maximum current possible, yet automatic adjustment of the transmitted current is performed internally based on the voltage measurements of the first potential dipole for each current dipole to avoid overload. For the measurements we used the same multi-core cable to connect both DAS-1 and Syscal Pro Switch 72.
The Terrameter LS 2 (Terrameter) is designed for DC resistivity and TDIP data acquisition. It supports data acquisition with the traditional 50% duty cycle, which we used in this study, and a 100% duty cycle, where a square wave with instant polarity change is the transmitted signal. The instrument has a built-in constant current transmitter capable of transmitting a current 2.5 A, limited by 600 V voltage or 250 W power. It records induced potentials in up to 12 channels simultaneously and supports up to 64 electrodes. Full waveform data were saved at a data rate of 3750 samples/s, for all input and current transmission monitoring channels. The input channels have an analogue low-pass filter with around 1500 Hz cut-off frequency. For the measurements presented in this study we used two multi-core cables, each one manufactured with 21 outputs (2 m separation) for the connection with the electrodes for the 1 m-electrode separation profile and four multicore cables, each one manufactured with 21 outputs (5 m separation) for the 5 m-electrode separation profile.
The Syscal Pro Switch 72 (Syscal) is a ten-channel acquisition system designed for resistivity and TDIP measurements. TDIP measurements were possible in a 50% duty cycle, with 20 voltage gates (defined by the operator) of the decay curve. Full-wave form data were collected at a sampling rate of 10 ms, yet only for the voltages. A low-pass filter (10 Hz cut-off frequency) is incorporated at the analogto-digital (A/D) converter. The 250 W built-in transmitter can transmit a constant voltage for a maximum current of 2.5 A. We deployed 72 electrodes connected by standard multi-core cables and requested injections with a voltage of 800 V in the current dipole, which will result in the current possible. Automatic adjustment of the transmitted current is however performed by the instrument based on the voltage measurements adjacent to the current dipole to avoid overload of channels. For the measurements presented in this study we used two multi-core cables, each manufactured with 36 take-outs (5 m separation).

Investigation area and experimental setup
The measurements were carried out in the Eastern Thuringian Shale Mountains, Germany. Airborne electromagnetic surveys conducted in the region revealed the occurrence of a shallow low resistivity (<5-10 Ωm) anomaly, corresponding to Graptolite blackshale, which is embedded between diabase in the NW and quartzite in the SE (Fig. 2  left). The blackshale is known to contain pyrite, which is commonly related to high IP effects (e.g., Pelton et al., 1978;Flores Orozco et al., 2011. In this study we present IP data collected along two profiles, where TDIP and FDIP measurements were performed with electrode separations of a = 1 m (IP 1) and a = 5 m (IP 5). The position, length and orientation of the profiles were selected to i) investigate changes in the IP response associated to different signal-to-noise ratio, which are expected to decrease with increasing the separation between electrodes and ii) assess possible variations in the IP response at different position over the Graptolite shale. We used the same electrodes with contact resistances of maximum 2 kΩ for all measurements and devices. Stainless steel electrodes were chosen due to the expected high polarization effects at this blackshale site (and therefore neglect electrode polarization, which we tested in the beginning of the measurements with unpolarizable Ag-AgCl electrodes) and the easier handling of the electrodes. The position of all electrodes was recorded by a DGPS instrument. Whereas the short profile IP 1 is only evenly dipping, the longer profile IP 5 shows significant topography.
The length of profile IP 1 was 41 m, using 42 electrodes. For FDIP we used 14 frequencies in the range of 0.156 to 1000 Hz (SIP256C) and 0.125 to 225 Hz (DAS-1). In TDIP we used a 50% duty cycle with a 4 s pulse each (ON-time, OFF-time) for both instruments. The Terrameter measured the decay curve in 13 time gate windows, whereas the Syscal used 20 time gate windows. To compare the data directly, the Syscal data were re-gated to the Terrameter time gate windows. Additionally, we processed the full waveform data from the Terrameter instrument to achieve a broader time range covered by a larger number of gates. The number of current injections are the same for all instruments. Hence, the numbers of quadrupols varies due to the different number of available and automatically used channels. In particular, the SIP256C is varying due to the roll-along scheme (two segments with an overlap of 7 electrodes). The number of stacks (TDIP) or periods (FDIP) differs for the different instruments (between 2 and 4) because we have chosen the smallest possible number for the specific conditions. An complete overview of the settings can be also found in Table 2.
The length of profile IP 5 was 315 m, using 64 electrodes. This profile was measured as a single profile using Terrameter and Syscal and, in roll-along mode, with DAS-1 and SIP256C. To this end, the profile was divided into three segments with an overlap of 14 electrodes (DAS-1), and 6 electrodes for the SIP256C due to the limited amount of available RUs. The measurement settings were similar to the settings for profile IP 1 and can be found in Table 3.
Considering, the DAS-1 permits only eight potential measurements (channels) for a given current dipole, the protocols were limited to eight dipole separations. Only for SIP256C measurements, we measured skip-0 data for all possible 22 dipole separations due to the automatic recording of all RUs in parallel. Skip-3 measurements were collected with all possible current dipoles for the DAS-1, Terrameter and Syscal instrument to get a full spatial distribution. Due to limitations in the setup of the RUs we could measure only every fourth dipole (e.g., C1-C2 = 1-5, P1-P2 = 5-9, P1-P2 = 9-13, ..) with the SIP256C.

Data processing
Even though all instruments are able to save the raw time series of the registered voltages, both FDIP instruments directly provide internally processed amplitude and phase values along with standard deviation from stacking. For the SIP256C instrument the time series are calculated by means of an FFT after some basic processing (50 Hz and drift filters) is done. The DAS-1 also uses a set of filter functions to reject low frequency noise (≤ 5 Hz) and/or power-line noise. These filters require two full waveforms. Data collected at higher frequencies (≥ 7.5 Hz) using the average over four full waveforms.
In FDIP, an automatic processing protocol was applied, which results in the removal of very few obvious outlier spectra which do not match the neighboring data points. For the SIP256C data set we removed 365 of 526 data points for profile IP 1 and 349 of 759 data points for profile IP 5. For the DAS-1 data set 644 of 835 data points were removed for profile IP 1.
In time domain, the waveforms are gated by the instrument to reduce the amount of data to a limited number of time windows, typically multiples of the power line and train periods of 20 and 60 ms, respectively. The time lengths of the windows is increasing, starting from 0.06 s for the earliest time until 0.84 s for the latest. Here, we used 13 gates which were originally designed for avoiding train disturbances, for processing the time-series internally. For the reasons of comparison, we used the Syscal full waveform data computed to the same time gates as were measured by the Terrameter.
2.4.1. Full-waveform processing Fig. 3 provides an impression on the full waveform processing for Terrameter data of profile IP 5 for a selected current dipole (26-27) and three different dipole separations (n = 1/4/7) exhibiting different voltages (400 mV/ 1 mV/ 0.1 mV) and thus different S/N ratio. Besides the raw data, a low-pass filtered curve (30 Hz corner frequency) is plotted. It shows two 16 s long periods with a lead-in of 6 s, 2 s delay between the signals, and 2 s lead-out. Whereas for n = 1 the voltages are large and very clean, the amplitude decreases with increasing dipole separation and therefore the S/N ratio is reduced significantly. There is a significant offset and low-frequency (background) drift that can be caused by spontaneous potential and electrode charge-up effects. Whereas for n = 4 the decay is still visible, this is almost impossible for n = 7, where some erratic behaviour is observed.
Modern post-processing techniques can improve the signals and extend the spectral range of the TDIP data by the procedures described by Olsson et al. (2016). That includes (1) drift removal, (2) harmonic noise cancellation, (3) spike removal and (4) tapered gating. So far, we could apply the described algorithms only to the full waveform data from the Terrameter instrument mainly due to the much lower sampling rate of the Syscal (100 Hz) and the internal low pass filter. Furthermore, the exact current waveform is not known for the Syscal instrument.
Although the post-processing by Olsson et al. (2016) used to increase the bandwidth is sophisticated, it will not deal with the above described type of erratic low frequency noise. The origin of it is not known, but it could possibly be caused by nearby wind power plants (<1 km away) or electrically powered railway traffic (5-6 km away). Therefore, the processing routine is usually followed by manual inspection of the decays resulting in the deactivation of erratic data points in particular for the very early times (compare Fig. 8 with respect to Fig. 12, gray coloured points) or (rarely) obvious strange decays of single quadrupoles for the inversion. In parts, also negative data, which do not follow a systematic trend and couldn't be explained by geometrical effects were removed. In the following, we show all available data and use only manual processed data for the TDIP inversion of the Terrameter. Here, we removed only very minor quadrupols for profile IP 1 (6/835 points) but many more for profile IP 5 (144/460 points).

Raw data visualization and inversion
For the data processing, visualization and inversion we used the freely available open-source package pyBERT, based on the Python modeling framework pyGIMLi (Rücker et al., 2017). It provides specialized modules for FDIP data  and now also TDIP data. Both modules allow importing of different instrument formats, data pre-processing (filtering, outlier removal and combination of subsets), different kinds of inversion and visualization of data and models as spectra/decays or (pseudo)sections.
In our study, inversion was performed using a Gauss-Newton minimization using first-order smoothness regularization on irregular triangular meshes (Günther et al., 2006). The FDIP inversion procedure has been described in more detail by Martin and Günther (2013) for inverting single-frequency, and by  for inverting multi-frequency data simultaneously. Unlike other authors, they do not assume a certain behaviour (like constant phase or Cole-Cole) beforehand, but use spectral constraints to regularize the individual models along the frequency axis.
For TDIP, a similar approach is followed, only that resistivity inversion is done as the first step. While the apparent chargeability is inverted into subsurface chargeability using the approach of Seigel (1959), the models for neighboring time gates are constrained along the time axes so that the model decays are just assumed smooth, without requiring any model.
After inversion, the model spectra, or decays, can be fitted using a Cole-Cole model in order to reduce the number of unknowns to four for each model cell. Both data and model spectra (FDIP) as well as decays (TDIP) can be transferred into a Debye distribution by a Debye decomposition Nordsiek and Weller (2008). For a fixed number (here 100) of predefined relaxation times (logarithmically spaced between 1 ms and 4 s) a partial chargeability distribution is computed using smoothness constraints between neighboring relaxation times. This allows comparison between frequency and time domain data or models.

Results
3.1. Profile IP 1 3.1.1. Apparent resistivity The apparent resistivity for data collected with the four instruments along profile IP 1 is shown in Fig. 4. Note, the upper part (DD n) denotes dipole-dipole skip-0 measurements, while the bottom part corresponds to skip-3 readings (DD4-n). Whereas the measuring arrays of DAS-1, Terrameter and Syscal are identical (skip-0, n max = 8 and full skip-3), the SIP256C provides the full skip-0 readings (automatically) but only a very limited number of skip-3 data, meaning only the upper eight rows can be compared directly between all instruments.
The pseudosections of the data are practically identical, which also holds for the skip-3 data of the three instruments using multi-core cables. In general, higher apparent resistivity values (>200 Ωm) can be found in the upper layer whereas low apparent resistivity values are observed in greater depth, disconnected by slightly higher values in the middle. The sections also show that towards the end of the profile, the high apparent resistivity values decrease.

Analysis of spatial consistency
To compare individual frequencies from the FDIP measurements, Fig. 5 shows the apparent phase pseudosections of skip-0 data with a maximum electrode separation of n = 8 for the selected frequencies f = 0.156 Hz/0.125 Hz, f = 1.2 Hz/1.0 Hz and f = 80 Hz/75 Hz (SIP256C/ DAS-1). Negative values are indicated light gray and out-ofthe scale values by dark gray. White values indicate missing data (e.g. due to roll-along of the SIP256C). The white rectangular boxes indicate selected data points discussed in the following section.
The apparent phases of the lowest frequency (first row) look nearly identical and altogether very smooth with higher values towards the end of the profile and higher separations. Whereas the pseudosections for data collected with the SIP256C data remain smooth for all frequencies, the DAS-1 data lose spatial consistency and show to be more affected by outlier for measurements collected at higher frequencies (> 75 Hz).
For TDIP, Fig. 6 shows the apparent chargeability pseudosections of the skip-0 data (n max = 8) for selected time gates (at t = 0.03-0.04 s (first row, only to the right), t = 0.07-0.13 s (second row), t = 0.31-0.43 s (third row) and t = 2.47-3.13 s (forth row)) in analogy to Fig. 5. In contrast to the left (Syscal) and middle (Terrameter) column, where the originally exported data from the instrument were used; the pseudosections from the right column show the processed pseudosections for the Terrameter full waveform data (according to Olsson et al., 2016). Therefore, even pseudosections for early times (0.03-0.04 s, first row) could be realized. Here, only very few outlier needed to be removed. Similar to the FDIP data (Fig. 5), all IP pseudosections show smooth patterns up to n = 8 in both, early and late, times. The early times seem to be a bit smoother for the Syscal (compared to the Terrameter), most likely due to the internal lowpass filter of the Syscal. Fig. 7 shows the apparent phase spectra for measurements collected using electrodes 26 and 27 for current injection and electrical impedance along three dipoles with electrode separation factors of n = 1, 4, and 7 (M-N electrodes: 28-29, 31-32, and 34-35) representing spectra with lower signal strength. The solid lines indicate SIP256C measurements, the dashed lines DAS-1 data.

Analysis of the frequency-dependence and decay curve
It can be seen that the SIP256C spectra are quite flat over the first three frequency decades (0.1-100 Hz) for all three selected quadrupoles with a slightly phase maximum at around 40 Hz for smaller separation factors. At higher frequencies, the data are still smooth but show a indication of EM coupling which starts earlier for the higher electrode separation factors (n = 7) at around 50 Hz.
The DAS-1 spectra are comparable for smaller electrode separation factors (n = 1, 4) at lower frequencies < 30 Hz. Above that frequency, the data are much earlier affected by EM coupling and start to become noisily. Furthermore, for measurements with larger separation between current and potential dipoles (pink line), the spectra collected by the two instruments show some systematic shift and a systematic deviation starting at about 3 Hz.
The TD decay curves recorded by both instruments are presented in Fig. 8, using the same quadrupoles as in Fig. 7. The dashed line with triangles represents the unprocessed data from the Terrameter, the solid line with triangles the processed data. Strong colours with circles indicate original Syscal data. Gray symbols/lines indicate processed but out-of-scale data.
The unprocessed decay curves show a consistent shape, particularly for later times (> 0.1 s), but reveal small shifts in the magnitude between the Syscal and Terrameter that might be caused by the slightly different waveform (transmitted waveform, delays etc.). Also, the first gates of the Syscal are systematically disturbed due to the (above described) low-pass filter. After processing of the Terrameter data, a significantly improvement of the spectral range can be observed for early times. In particular for small electrode separations (n = 1, blue line) consistent data can be obtained from approx. 0.005 s. However, for higher electrode separations/lower S/N ratios (n = 4, 7) the processing can reveal only very small or no plausible data for earlier times (> 0.05 s).
To compare the FDIP and TDIP measurements directly, Debye decompositions for the TD full waveform data (Terrameter) and the SIP256C data were carried out. In Fig. 9 the spectral chargeability distribution as a function of relaxation times are plotted for the above shown quadrupoles (TDIP: dashed lines, FDIP: solid lines).
In general, the FDIP and TDIP curves from the matching electrode separations (n = 1/4/7) are in the same range, but the shape of the curve varies. For n = 7 (red/pink lines), both curves (TD and FD) show a maximum at a similar relaxation time (around 0.5 s). Also for n = 1 (blue lines) both curves start and stop at almost the same values. For n = 4 (green lines) the magnitude of the chargeability is similar for the middle part of the relaxation times but varies a lot for the minimum and maximum relaxation times. The reason for the difference is yet to speculate and might be due to he slightly different spectral content, e.g. caused by the measuring scheme. Furthermore, due to the smoothness constraints in the inversion, ambiguities in the Debye decomposition can arise.

Inversion results
The results of the spectral constrained inversion for profile IP 1 are displayed in Fig. 10 for all instruments, combining skip-0 and skip-3 data. For the Terrameter we used the processed data. Only a few values needed to be removed manually for all data sets.
All apparent resistivity spectra (both TD and FD) could be fitted within a range of 2% to 4% for the relative root-mean square between measured and modelled data. The fit of the IP data inversion is different for both measuring scheme and cannot easily be compared. Note that the FDIP inversion directly inverts for a phase image, whereas for the TD data the spectral chargeability is inverted and constrained and later transferred into a maximum phase model using another fitting process, i.e. inversion. The obtained root-mean square for the SIP256C data (5.7 mrad) was lower than for the DAS-1 (11.2 mrad). On the other hand, the data fit for both TDIP instruments was comparable (8.2% and 8.8% for Syscal and Terrameter, respectively).
On the left-hand side, the resistivity sections are shown. The variation between the individual instruments is small, only SIP256C is marked with slightly higher values at depth due to the less amount of large-separation data. In general, we observe a resistive (200-500 Ωm) overburden with a thickness of 2-4 m that is thinning from x = 30 m on. Below, more conductive (< 50 Ωm) material is found, forming almost a layer with some lateral discontinuities (e.g., around 10 m and between 25 and 30 m). Partly very low values of about 5-10 Ωm show up between the discontinuities.
On the right-hand side, the maximum phase shift for all instruments are shown. Note that for this is a simple maximum function for FD, whereas for TD it represents a fit with a constant-phase model. In accordance with the resistive overburden layer, the first meters are characterized by low phases. Below three meters depth, the phase increases significantly up to 120 mrad. Interestingly, the high phase layer is continuously for both FDIP instruments, but only pronounced in the right part of the profile for the TDIP measurements. Even though we would expect the same spectral information, it is differently contained in the measured decays and spectra and therefore projected into the subsurface by the inversion process. The different model parameterization (spectral phase and chargeability) might be one reason for the different imaging results that could be equivalent considering the measuring accuracy.
We interpret the overburden layer as Quaternary deposits of sand and loam. The good conductor with the highly polarizable zones is interpreted as the Silurian graptolite blackshale (compare Fig. 2). According to core-drilling information in that area, these shales contain iron sulphides (e.g., pyrites) which explains the high IP values due to electrode polarization mechanisms, associated to metallic minerals.

Profile IP 5
For the discussion of the influence to different signal strengths, the data from profile IP 5 (electrode spacing of 5 m) are also discussed. Here, we focus only on one FDIP instrument (SIP256C) and one TDIP instrument (Terrameter). Reasons for that are the better performance at high frequencies from the SIP256C device (in contrast to the DAS-1) and the possibility of the processing for the Terrameter full waveform data.

Analysis of spatial consistency
In Fig. 11 the pseudosections for readings collected with eight separation levels (n max = 8) are shown. On the left hand side the results from the SIP256C instrument for the apparent resistivity (top row) and the phase for three selected frequencies (f = 0.156 Hz/2.5 Hz/ 40 Hz -row 2/3/4) are plotted. At the right hand side the apparent resistivity for the Terrameter is shown together with the pseudosections of three different time gates (at t = 0.09-0.13 s (second row), t = 0.51-0.71 s (third row), t = 2.7-3.7 s (fourth row)). Note the missing values in the SIP256C data due to the three roll-along parts (triangle shape at 100 and 200 m).
For both resistivity pseudosections the data show consistent patterns. But due to the full removal of data points with erratic decays, the pseudosection for the Terrameter instruments is reduced. Also, in FDIP occasional readings needed to be removed due to erratic behaviour at some electrodes (white lines in pseudosection).
In TDIP, the different time gate data are much noisier for measurements along this profile compared to profile IP 1, revealing a decrease in data consistency with increasing electrode separation and also to earlier times. The FDIP phase measurements are also considerably noisier compared to profile IP 1, but still consistent up to n = 8. However, the difference between the individual frequencies indicate only minor spectral content at this profile.

Analysis of the frequency-dependence and decay curve
In analogy to profile IP 1 the spectra (left) and decay curves (right) for selected quadrupoles from profile IP 5 are presented in Fig. 12. Gray colours indicate bad data points that are not used for subsequent analysis and inversion (see below). Especially higher S/N ratio (n = 1, 4) show very smooth curves over the entire frequency range in FD. Only higher electrode separations (n = 7) seems to be noisier but still usable up to 50 Hz. In TD, the processed data show good data up to three decades for high S/N ratio (n = 1, blue lines). With decreasing S/N ratio (n = 4 -green line, n = 7 -red line) early time data can be calculated but may not considered for further spectral analysis/interpretation due to unrealistic and/or bad negative data (black circles) indicating inductive effects.

Inversion results
For reasons of comparison we show the inversion results for the lowest frequency, i.e. around 0.1 Hz (Fig. 13). To this end, the TDIP data have been simultaneously inverted and the resulting model decays have been fitted using a Cole-Cole term assuming a low c value that reflects the more or less constant phase behaviour seen in Fig. 12. Note that the root-mean-square errors were significantly higher compared to IP 1 with 5.8% (SIP256C) and 7.4% (Terrameter) for the apparent resistivities, 14.7 mrad (SIP256C) for the phases, and 26.4 mV/V (Terrameter) for the spectral chargeabilities.
Analogue to profile IP 1 (which is a small part of that profile and located between 180 and 221 m, see Fig. 2), the first 2-3 m are characterized by higher resistivities. Below the higher resistivities a highly conductive body (< 5 Ωm -indicating the graptolite shale) shows up in the middle of the profile (between 100 and 250 m) and reaches greater depths. Also, a near surface conductive zone at approx. x = 55 m can be recognized. The differences in resistivity are small for both domains. Only the interruption of the body at 170 m in the Terrameter data marks some differences.
In contrast, the results in phase from both domains (TD and FD) show higher variation due to the limited data quality for high dipole distances. Below a near surface low phase layer, higher phases occur almost along the entire profile. But whereas the layer is thin (5 m) in FDIP, the visible high phase anomalies in TDIP reflect three zones with extended dimensions at x = 30 m, 100 m and 230 m. The highly Fig. 6. Apparent chargeability pseudosections (skip-0) from TDIP instruments (left: Syscal, middle: Terrameter unprocessed, right: Terrameter processed) at different time gates: t = 0.03-0.04 (first row), t = 0.07-0.13 s (second row), t = 0.31-0.43 s (third row) and t = 2.47-3.13 s (forth row). White colour indicate removed data (e.g. due to outlier), light gray are negative values and dark gray out-of-the scale values. The white rectangles indicate the data discussed later in Fig. 8. polarisable thin layer in FDIP correspond well with the transition zone in resistivities between high and low values. In contrast, the high phase zones in TDIP are mainly in agreement with the low resistivity areas, except the part at x = 150 m.
According to the interpretation of profile IP 1, the large highly conductive zone in the middle of the profile can be interpreted as the Silurian Graptolite blackshale (compare Fig. 2). The zones of high polarization are not exactly in agreement with this body. In fact, the highest phase shifts occur in the transition areas to both sides of the blackshale body. Obviously, most of the iron sulphides were formed here.

Discussion
The apparent resistivity values measured with all four instruments used in this study are practically identical and reveal a broad low resistivity zone in depth. Also, the IP values are comparable and show partly high values. The low resistivity and high IP data are most likely caused by the presence of iron sulphides in the blackshale (Dill, 1986) and lead to high signal-to-noise (S/N) ratios (e.g., Pelton et al., 1978;Flores Orozco et al., 2011. This favours the IP comparison, resulting in pseudosections with consistent patterns, for both TDIP and FDIP.
Furthermore, the inversion results show almost identical results for the resistivity sections for both profiles IP 1 and IP 5. In general, the maximum phase angle sections show similar results as well but revealing also some differences between the IP domains by showing not so pronounced high polarisable layer in TD for profile IP 1 and partly stronger polarization zones in IP 5. Reason for that might be the limited data quality for larger dipole distances and therefore lower S/N ratios.

Detailed data evaluation
Besides resistivity data, FDIP measurements at low frequencies are quantitatively comparable between both instruments (SIP256C and DAS-1). However, the DAS-1 seems to be limited in the quality of the readings at frequencies above 10 Hz. In contrast, consistent data can be observed with the SIP256C with regard to both spatial and frequency variations for the first three frequency decades (0.1 Hz -100 Hz). To higher frequencies (until 1000 Hz) the curves are still smooth but highly affected by EM coupling. If this coupling is taken into account and modelled, the useful frequency range can even be wider. The high consistency in the SIP256C readings in a broad frequency bandwidth, compared to the DAS-1, can be explained due to the digitization of the potential differences directly at the RUs, placed next to the electrodes avoiding coupling between the transmitter and the receiver. Furthermore, the RU technique minimizes cross-talk between cables, which significantly improve readings at higher frequencies (Radic, 2008;Schmutz et al., 2014). Measurements with the DAS-1 performed within this study were conducted with common multi-core cables, as those used for DC resistivity, without shielding between the different wires. Hence, measurements will be primary affected by cross-talk effects within the cables, explaining the lack of spatial consistency in the DAS-1 data at frequencies above 10 Hz. This effect has been addressed before (e.g., Flores Orozco et al., 2011), yet in the present study we decided upon this setting to clearly present the contamination of FDIP data at frequencies above 10 Hz. As a result, low-frequency data collected with the DAS-1 are quantitatively comparable with those obtained with the SIP256C clearly revealing the possibility to perform FDIP measurements with high quality using easy-handling multi-core cable as those used in DC resistivity.
The two used TDIP instruments (Syscal and Terrameter) reveal pseudosections with the same patterns for measurements collected at different times within the decay-curves, over a broad time range (>0.1 s). However, the well-coinciding pseudosections (Fig. 6) may mask small systematic shifts, visible only after a detailed analysis of the actual recorded decay-curves (Fig. 8). A direct comparison of the decays between both instruments is also affected by the low-pass filter integrated in the Syscal device.
The time range of TDIP measurements can be extended by carry out an enhanced processing using the full waveform data. Here, we applied the procedure to the full waveform data from the Terrameter instrument as published in (Olsson et al., 2016). Obviously, the usable spectral range could be extended significantly, in particular for high signal-tonoise ratios. Even for profile IP 5, early times of almost 5 ms could be reached. With decreasing S/N early time data can be calculated but may not considered for further spectral analysis/interpretation due to unrealistic and/or unreliable negative data (black circles) indicating inductive effects. Note, that processing the Syscal full waveform data the same way would be more complicated due to the fact that the input filter needs to be considered, and less efficient since the data sampling rate is more than a magnitude lower.  Comparing the bandwidth of the TDIP and FDIP data, at least three decades spectral range can be reached in the frequency-domain (Fig. 7), even for lower S/N ratios. Also, for the highest measured frequency (f = 1000 Hz), pseudosection with high spatial consistency and with only minimal occurrence of outlier could be registered (with the SIP256C, Appendix Fig. 14). Nevertheless, the EM effects, which occur at high frequencies, needs to be modelled and taken into account for modeling/inversion which is behind the scope of this study.
The TDIP instruments achieved a bandwidth of two decades from the decay data (Fig. 8), but the low frequency bandwidth could have been increased by using longer transmitted pulses and measuring the decays during correspondingly longer time (max. half a decade with the limitations set by the GUI in the instrument software). Increasing the bandwidth in the upper frequency/early time range in TDIP is not possible directly from the instruments decay data but requires processing of the full waveform (e.g., Fiandaca et al., 2013;Olsson et al., 2016). Doing so, the spectral range could be increased significantly by almost one decade for high S/N ratios. For lower S/N (n ≥ 4) the data for the early times getting implausible and could therefore not be considered in the inversion process.
Our data are not in full agreement with the observations reported by Maurya et al. (2018b) with the same instruments (SIP256C and Terrameter), where only two decades of spectral range could be observed after processing in FDIP, and four decades for TDIP. As could be seen, in TD significant improvements could be done by the processing but we do not find any reason for removing up to two decades of high frequency FDIP data because of inductive coupling. In our data we were able to gather also apparent induction free data above 10 Hz. Furthermore, in our understanding, even induction-contaminated data could be used by considering and modeling these effects.

Signal strength
Comparing the signal strength for both profiles in general (a = 1 m for profile IP 1 and a = 5 m for profile IP 5), one can see that larger electrode distances significantly worsen the S/N ratio for both, TDIP and FDIP, measurements. As shown in Fig. 3 for the IP 5 Terrameter data, the S/N ratio for the shortest n-factor (n = 1) is very favourable and decreases with increasing n-factor, resulting in that the measured raw data signal drowns into noise for high n-factor (n = 7). In FD, higher electrodes distance also affect the IP measurements but mainly for the DAS-1 instrument and minor for the SIP256C. Whereas the DAS-1    provides similar results and bandwidth as both TD instruments, the differences for higher frequencies (> 30 Hz) compared to the SIP256C are significant. Even for large dipole separation and lower S/N ratio the SIP256C data are still consistent over a broader frequency range in contrasts to the decay curves/time gates in TD (Figs. 11 and 12, right side). Again, the main reason for that is the above mentioned remote unit technique used by the SIP256C instrument which leads to a minimization of the EM coupling within the data acquisition cables. As a result, good data even at high frequencies can be measured but, on the other side, this is associated with a greater time effort in the field.

Comparison of acquisition time and field settings
Due to the slightly different settings used for different instruments (Table 2), we calculated the acquisition time with identical settings for profile IP 1 to compare the data dispassionate. As basis we used only skip-0 data for n max = 8 with a stack/period of 2, and a pulse length of 4 s resp. a frequency range between 0.25 and 225 Hz (DAS-1)/ 0.3-266 Hz (SIP256C). We also took the roll-along into account. In Table 1 the respective acquisition times can be seen. It is evident that the times for the DAS-1 is not even half of the times from SIP256C instrument. One reason for that is the time series transmission in the SIP256C that takes even longer as more RUs are connected (here: 24 RU). This means also, that automatically more quadrupols are measured (n max = 24 instead of n max = 8) for each current injection. Furthermore, the use of real sinusoidal waves instead of Fourier transformations is slower but also enhances accuracy.
Comparing acquisition time in TDIP, the Terrameter needs almost twice the time compared to the Syscal. The discrepancy is mainly caused by limitations in the hardware design of the built-in relay switch in relation to the receiver channels (even for a 12 channel system), which is unfavourable for dipole-dipole array with n-factor 1-8. It leads to 7 measurements for one current transmission in most cases, so that one extra current transmission is needed for each current electrode pair for the skip-0 measurements. On the other hand, the acquisition time for the skip-3 measurement is almost identical for both instruments, indicating that this measurement sequence is more favourable (see Table 2). Some differences may also be explained by the different waveform (length of the total signal). In contrast to the Syscal, the Terrameter has a constant current transmitter and checks the optimal settings for current injection by sending test pulses before each measurement to gather the used current pulse sequence exactly. Additionally, it measures background levels (2 s, see signals in Fig. 3) before and after each measurement cycle to calculate and re-trend the background drift.
Obviously, measuring in FDIP and TDIP, similar acquisition times can be reached. FDIP field measurements are not per se more timeconsuming than measurements in TDIP. As example, measurements collected with the DAS-1 (FDIP) and the Terrameter (TDIP) required similar acquisition times and field procedures (Table 1) and they reveal consistent pseudosections at low frequencies/late times (up to 10 Hz resp. 0.1 s). The use of multi-core cables are easy to deploy in the field and save acquisition time in both domains. However, the data collected at high frequencies resp. early times reveal a poor spatial consistency, reflecting the contamination of the data by cross-talk and electromagnetic coupling (Dahlin and Zhou, 2006;Flores Orozco et al., 2011Kemna et al., 2012). The use of separated multi-core cables for current injection and potential measurements can minimize the EM effects (Dahlin and Leroux, 2012). Shielded cables can also improve the data quality (Flores Orozco et al., 2013). Furthermore, measurements could be performed using algorithms to take into account the geometry of the individual wires in the multi-core cables to enhance data quality in measurements at high frequencies (Zhao et al., 2013). Further acquisition time saving in TDIP measurements can be achieved by using a 100% duty cycle instead of the standard 50% duty cycle which was used in this study to realize a comparison between the Terrameter and the Syscal instrument. Olsson et al. (2015) demonstrated that the measuring time can be decreased by almost a factor of two.
In this study, we used the dipole-dipole array as this performs potential readings with electrodes which were not used as current electrodes before, minimizing the contamination of the data by means of polarization effects at the electrodes. Other arrays, such as Wenner alpha, leads to polarization effects taking place at the potential electrodes. However, the signal strength decreases rapidly with increasing dipole separation whereas arrays as the multiple gradient array are able to provide enhanced S/N ratio. Furthermore, these arrays may also provide lower acquisition time (Dahlin and Zhou, 2006), yet some studies have also revealed similar results for data collected with both multiple gradient and dipole-dipole configurations (Flores Orozco et al., 2018b. Using adequate configuration systems can also avoid (unintentional) charge-up effects due to ground contact test currents (up to  Table 1 Comparable acquisition times for profile IP 1 (dipole-dipole, skip 0, n max = 8, stacks = 2, pulse length = 4 s (Syscal/Terrameter) resp. frequency range f = 0.25-225 Hz (DAS-1)/ 0.3-266 Hz (SIP256C)).

SIP256C
DAS-1 Terrameter Syscal 20 mA), measuring the potential just after using the electrodes for current injection or due to previous IP measurements (Dahlin, 2000). Such charge-up effects can lead to large errors in late time decay TDIP data due to the linear background trend removal typically applied in such instruments, but can be handled using the higher order functions for trend removal as suggested by Olsson et al. (2016). Based on our observations we suggest that field procedures can benefit of extensive mapping at lower frequencies using multi-core cables, and the selection of particular profiles for the acquisition in a broader frequency range, e.g., using an RU system or sophisticated cable setups to minimize cross-talk. In general, the quality of IP field data benefit very much from the effort done in the field, e.g. good ground contact (e.g. watering/using gel for the electrodes), optimized arrays and the measurement of normal and reciprocal data to retrieve error models by statistical analysis of the reciprocity (Flores Orozco et al., 2012b;Flores Orozco et al., 2018;Zhou and Dahlin, 2003). Using recent instrument generations together with the newest data processing and inversion routines can significantly increase the spectral range and thereby the use of the IP methods in general for further applications.

Conclusions
Both, FDIP and TDIP data provide consistent results for measurements collected with different instruments. The DAS-1 data have shown that measurements in FDIP can be used routinely in the same way as in TDIP, however with similar limitations in data quality for high frequencies. In general, only small differences could be observed by comparing the bandwidth in frequency-domain and time-domain by the use of the latest processing algorithm for TD full waveform data. In our study we were able to measure up to three decades without apparent EM noise in FDIP (SIP256C) even for higher electrode separations, resp. low S/N ratios. In TDIP two decades, and after data processing up to three decades, for high S/N ratios could be reached with the field settings presented here. By using longer pulses it could have been extended even more (up until half a decade).
To increase the quality of high-frequency/early-time IP data with multi-core cable instruments (DAS-1, Syscal, Terrameter), more effort is needed, e.g., by using separated and/or shielded cables for current and potential measurements. The SIP256C data demonstrates that the collection of broad bandwidth SIP spectra is possible under field conditions, also for investigations with larger electrode spacing, resp. low S/N ratios. However, it goes along with larger acquisition times and advanced field settings which can be reduced using the newest generation of instruments.
For standard IP applications (e.g., mineral exploration) the use of state-of-the-art instruments and processing might not be always needed but for more detailed investigations over particular lowpolarisable targets (e.g. in hydrogeological applications) it can be crucial.

Declaration of Competing Interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.
In Fig. 14 is shown the apparent resistivity and phase for the highest measured frequency at f = 1000 Hz (profile IP 1). Apart from the consistent apparent resistivity data, also the phase pseudosection show data with only minor outlier. Especially for higher frequencies, these data are affected by unavoidable induction EM effects which needs to be modelled and considered in the inversion. Even when we don't aim for the modeling in this study (which is behind our scope) the high frequency data can be used and are therefore classified as good data.