https://orcid.org/0000-0002-3663-182X
SUMMARY
Although controlled-source electromagnetic (CSEM) methods have higher sensitivity to thin resistive bodies than the magnetotelluric (MT) method, their delineation by the inversion requires CSEM data with high signal-to-noise ratio (SNR). This study aims to enhance the SNR of CSEM data by increasing the number of stacks. To efficiently stack long-term data, we use an EM-accurately controlled, routinely operated signal system (ACROSS), which can transmit accurately controlled waveforms by synchronizing the transmitting waveforms with a 10 MHz Global Positioning System signal. We conducted a CSEM survey using the EM-ACROSS in the Kusatsu-Shirane Volcano to demonstrate that the SNR can be improved by extensive observation data and the CSEM inversion can delineate hydrothermal systems, including resistive bodies of vapour-rich reservoirs. Our EM-ACROSS simultaneously transmitted waveforms from two dipoles during a 192-h of the survey; five-component receivers located 4–6 km away from the transmitter captured EM-ACROSS signals ranging between 146 and 192 h. By stacking extensive observation data using a weighted method, the CSEM responses show minimal error levels, with standard errors <2 per cent for most frequencies. The SNR roughly followed the square root of the stacking times. 3-D inversion of the collected CSEM data delineated a relatively resistive body, interpreted as a vapour-dominated reservoir below a cap-rock layer, while the MT inversion failed to recover the same. This highlights the ability of an EM-ACROSS-based CSEM survey to delineate hydrothermal systems including vapour-dominated reservoirs, and provides a compelling rationale for establishing CSEM as a standard methodology in hydrothermal imaging. Furthermore, this study suggests that the enhanced imaging capabilities of CSEM data can be further improved when integrated with MT data.
1 INTRODUCTION
Electromagnetic (EM) methods can determine subsurface and subseafloor electrical resistivity structures and are useful for exploring hydrothermal systems (Bertrand et al. 2022; Ishizu et al. 2022a), ore deposits (Smith 2014; Mittet & Avdeeva 2023; Ishizu et al. 2024), freshwater aquifers (Blatter et al. 2019; Attias et al. 2020; Haroon et al. 2021) and hydrocarbon reservoirs (Eidesmo et al. 2002; Constable 2010; Hesthammer et al. 2010). The targets are detected by exploiting their resistivity properties that are higher or lower than those of the surrounding rock.
Controlled-source EM (CSEM) methods, which use artificial sources, have been employed to detect resistive bodies, such as oil and gas (Myer et al. 2015; Schaller et al. 2017), gas hydrate (Weitemeyer et al. 2006; Schwalenberg et al. 2017) and freshwater (Gustafson et al. 2019; Micallef et al. 2020), due to their superior sensitivity to thin resistive bodies compared to the magnetotelluric (MT) method using natural plane-wave sources (Ziolkowski et al. 2007; Constable 2010; Strack 2014; Streich 2016). CSEM methods also facilitate the acquisition of high signal-to-noise ratio (SNR) data in noisy environments by adjusting the source against a particular noise (Streich et al. 2013; Grayver et al. 2014; Védrine et al. 2023). Due to the higher sensitivity to thin resistive bodies and noise resistance, CSEM methods have also piqued the interest of researchers as a monitoring technique (Wirianto et al. 2010; Tietze et al. 2015; Minami et al. 2018; Darnet et al. 2020; Bretaudeau et al. 2021; Commer et al. 2022).
Even with their above-mentioned advantages over the MT method, delineating thin resistive bodies by CSEM inversions requires high-SNR data, which can be obtained by either increasing the source dipole moment (dipole length × output current) or extending the observation period for data stacking. In onshore areas, limitations on cable lengths and output currents of source dipoles are imposed due to poor accessibility and high-earth resistance. Therefore, we aim to enhance the SNR by extending the observation duration for data stacking. Given that noise is Gaussian with zero mean, SNR will increase as the square root of the stacking times. To efficiently stack long-term data, we use an EM-accurately controlled, routinely operated signal system (ACROSS), which can transmit accurately controlled waveforms over long measurement intervals, such as several weeks by synchronizing its waveforms with a Global Positioning System (GPS) signal (Ogawa & Kumazawa 1996; Nakajima et al. 2000). Accurately controlling repetitive waveforms helps efficiently stack long-term data as the spectra of transmitting currents only shows peak values at the signal frequencies without dispersion to other frequencies. The EM-ACROSS waveforms consist of superposed sinusoids with a limited number of frequencies. Thus, our EM-ACROSS survey can be classified as a frequency-domain CSEM approach. Receiver time is also based on GPS clocks, as EM-ACROSS requires time synchronization between the transmitter and receivers.
We conduct a CSEM survey using the EM-ACROSS in a volcanic area. The EM-ACROSS simultaneously transmit regulated waveforms from two dipoles using synchronization with a 10 MHz GPS signal. We deploy five-component receivers with offsets of 4–6 km from the transmitter. The range of transmitting frequencies is set to 0.1–100 Hz to include both near and far-field data. Near-field data produce vertical currents yielding high sensitivity to resistive bodies compared to that of the MT method. Plane-wave assumption is valid for far-field data, resulting in comparable sensitivity to the MT method. The study area is the Kusatsu-Shirane Volcano in Japan, where phreatic eruptions have repeatedly occurred at the Yugama crater (Yaguchi et al. 2021). The MT method has been routinely used for investigating various volcanic areas (Ingham et al. 2009; Bertrand et al. 2022) including the Kusatsu-Shirane Volcano (Tseng et al. 2020). However, due to its lower sensitivity to resistive bodies, the MT method can fail to detect resistive bodies of vapour reservoirs, which are necessary for the in-depth understanding of hydrothermal systems in volcanic areas and forecasting phreatic eruptions (Stix & de Moor 2018).
The objectives of this study are to demonstrate that stacking can improve SNR of the CSEM data obtained using the EM-ACROSS and inversion of the measured CSEM data can delineate hydrothermal systems up to 2 km in depth, including resistive bodies of vapour reservoirs. We first describe our EM-ACROSS instruments and a weighted stacking method for data processing. Then, we present the result of the EM-ACROSS application to the Kusatsu-Shirane Volcano, focusing on relation between SNR and stacking time and comparison of CSEM and MT data inversion results.
2 METHODS
2.1 CSEM instrument
The developed EM-ACROSS transmitter (Fig. 1) is an improved version from the prototype by Tseng et al. (2019). The time base of our waveform generator (WF1974; NF Corporation, Yokohama, Japan) was synchronized with a 10 MHz signal from the GPS clock (Syncserver S650, Microchip Technology Inc., Chandler, AZ, USA) to repeatedly transmit accurately controlled waveforms. Two power amplifiers (AA2000XG2; Takasago International Corp., Ltd, Tokyo, Japan) were used to amplify the accurately controlled waveforms generated by the waveform generator by a factor of eighty in voltage-controlled mode. The amplified waveforms were transmitted from north–south (NS) and east–west (EW) oriented dipoles. A data logger (SC-ADH10K; Scimolex Corp., Inc., Yamanashi, Japan) connected to the amplifiers recorded the output currents at sampling rates of 1 and 10 kHz. Receivers (MTU-5C; Phoenix Geophysics Ltd, Scarborough, Canada) recorded electric fields in the NS (Ex) and EW (Ey) directions as well as magnetic fields in the NS (Hx), EW (Hy) and vertical directions (Hz). The x-, y-, and z-directions represent the geographic north, east and vertical depths, respectively. The MTU-5C receivers used the same gain and resolution settings as typical MT measurements. The SC-ADH10K data logger and MTU-5C receivers also used GPS clocks for time synchronization.
Schematic diagram of the EM-ACROSS transmitter. The time base of the waveform generator is synchronized with a 10 MHz GPS signal from the GPS clock, enabling signal control with high precision. Two power amplifiers boost the accurately controlled waveforms from the waveform generator 80 times. The amplified waveforms are injected from the NS and EW dipoles to the ground. A data logger connected to the amplifiers records the output currents at a sampling rate of 1 and 10 kHz.
Our EM-ACROSS waveforms consisted of superposed sinusoids with eight different frequencies (Tseng 2020):
where I is source current (A) for time t (s), |${{A}_j}$| is amplitude (A), |${{f}_j}$| is frequency (Hz) and |${{\phi }_j}$| is phase shift (°) for index number j. Table 1 lists the coefficients of EM-ACROSS waveforms. The period of the waveforms for both the NS and EW dipoles was 10 s. The frequency ranges for the NS and the EW dipoles were 0.2–46.1 and 0.3–46.3 Hz, respectively. The frequency of the NS waveform differed slightly from that of the EW waveform, allowing for simultaneous receiver recording without interference. Furthermore, as EM-ACROSS transmission only affects the MT signals at a restricted number of frequencies (i.e. 16 frequency number), MT data can be acquired while the EM-ACROSS transmission is running. The amplitude of |${{A}_j}$| was set to higher values at lower frequencies (Tseng 2020) to overcome natural noise, which is often higher at lower frequencies (Macnae et al. 1984). Different phase shifts of |${{\phi }_j}$| were used to lower the peak amplitude value (Tseng 2020).
| j | 1 | 2 | 3 | 4 | 5 | 6 | 7 | 8 |
|---|---|---|---|---|---|---|---|---|
| |${{A}_j}$| | 2.3 | 1.7 | 1.3 | 1.1 | 0.7 | 0.5 | 0.3 | 0.2 |
| |${{f}_j}$| (Hz) for NS | 0.2 | 0.5 | 1.1 | 2.9 | 5.9 | 10.1 | 22.7 | 46.1 |
| |${{f}_j}$| (Hz) for EW | 0.3 | 0.7 | 1.3 | 3.1 | 6.1 | 10.3 | 22.9 | 46.3 |
| |${{\phi }_j}$| (°) for NS | –15.3 | 24.5 | –41.8 | 151.7 | 115.6 | –28.9 | –36.8 | 81.4 |
| |${{\phi }_j}$| (°) for EW | –28.5 | 53.7 | 161.2 | –16.8 | 149.9 | 69.5 | 66.0 | –26.9 |
| j | 1 | 2 | 3 | 4 | 5 | 6 | 7 | 8 |
|---|---|---|---|---|---|---|---|---|
| |${{A}_j}$| | 2.3 | 1.7 | 1.3 | 1.1 | 0.7 | 0.5 | 0.3 | 0.2 |
| |${{f}_j}$| (Hz) for NS | 0.2 | 0.5 | 1.1 | 2.9 | 5.9 | 10.1 | 22.7 | 46.1 |
| |${{f}_j}$| (Hz) for EW | 0.3 | 0.7 | 1.3 | 3.1 | 6.1 | 10.3 | 22.9 | 46.3 |
| |${{\phi }_j}$| (°) for NS | –15.3 | 24.5 | –41.8 | 151.7 | 115.6 | –28.9 | –36.8 | 81.4 |
| |${{\phi }_j}$| (°) for EW | –28.5 | 53.7 | 161.2 | –16.8 | 149.9 | 69.5 | 66.0 | –26.9 |
| j | 1 | 2 | 3 | 4 | 5 | 6 | 7 | 8 |
|---|---|---|---|---|---|---|---|---|
| |${{A}_j}$| | 2.3 | 1.7 | 1.3 | 1.1 | 0.7 | 0.5 | 0.3 | 0.2 |
| |${{f}_j}$| (Hz) for NS | 0.2 | 0.5 | 1.1 | 2.9 | 5.9 | 10.1 | 22.7 | 46.1 |
| |${{f}_j}$| (Hz) for EW | 0.3 | 0.7 | 1.3 | 3.1 | 6.1 | 10.3 | 22.9 | 46.3 |
| |${{\phi }_j}$| (°) for NS | –15.3 | 24.5 | –41.8 | 151.7 | 115.6 | –28.9 | –36.8 | 81.4 |
| |${{\phi }_j}$| (°) for EW | –28.5 | 53.7 | 161.2 | –16.8 | 149.9 | 69.5 | 66.0 | –26.9 |
| j | 1 | 2 | 3 | 4 | 5 | 6 | 7 | 8 |
|---|---|---|---|---|---|---|---|---|
| |${{A}_j}$| | 2.3 | 1.7 | 1.3 | 1.1 | 0.7 | 0.5 | 0.3 | 0.2 |
| |${{f}_j}$| (Hz) for NS | 0.2 | 0.5 | 1.1 | 2.9 | 5.9 | 10.1 | 22.7 | 46.1 |
| |${{f}_j}$| (Hz) for EW | 0.3 | 0.7 | 1.3 | 3.1 | 6.1 | 10.3 | 22.9 | 46.3 |
| |${{\phi }_j}$| (°) for NS | –15.3 | 24.5 | –41.8 | 151.7 | 115.6 | –28.9 | –36.8 | 81.4 |
| |${{\phi }_j}$| (°) for EW | –28.5 | 53.7 | 161.2 | –16.8 | 149.9 | 69.5 | 66.0 | –26.9 |
2.2 Weighted stacking method
The weighted stacking method (Nagao et al. 2003, 2010) was used to enhance the SNR of the EM-ACROSS data. It incorporated noise-level information into stacking weights. We acquired EM-ACROSS data over two weeks, during which the noise level substantially changed. Stacking noisy data without down-weighting can decrease the SNR of the EM-ACROSS response. Therefore, noise-level information should be integrated into the data processing.
The weighted stacking method divided receiver time-series data into segments, conducted a Fourier transformation, performed noise estimations, computed weights based on the estimated noise levels, and calculated the CSEM response with weights (Fig. 2). First, we divided the entire data by the duration of each segment Ts (s), yielding M segments, with each segment represented by m (1 to M). We removed the linear trend from the time-series data in each segment of m, subsequently transforming the time-series data into frequency domain spectra X(f) using a fast Fourier transform (FFT), where f represents the frequency. The spectra X(fj) at signal frequency fj are referred to as signal channels, and the other spectra, excluding the signal channels, are referred to as noise channels. The noise level |${{n}_{mj}}$| for segment m and signal channel j was estimated using noise channels as follows (Nagao et al. 2010):
where K is the one-side search number for noise channels around the signal channel j, k is the integer index for the noise estimation (|$- K \le k \le K$|), |$\Delta f$| is the frequency interval of the spectra, calculated as 1/Ts, and L is the total number of noise channels used in the noise estimation. If |$X( {{{f}_j} + k\Delta f} )$| for |$- K \le k \le K$| excluding k = 0 does not include signal channels, L becomes 2K. The weight |${{w}_{mj}}$| for segment m and signal channel j was computed as follows (Nagao et al. 2010):
Flow of data processing using the weighted stacking method.
Eq. (3) indicates that when the noise level |${{n}_{mj}}$| is low, the weight |${{w}_{mj}}\ $|is large.
We also separated the time-series data of the injected current of the transmitter for each duration Ts, which were transformed into a spectrum S(f) in the frequency domain using FFT with a time window of Ts. The duration of each segment Ts was same for the receiver and transmitter data. CSEM transfer function|$\ {{G}_{mj}}$| for segment m and signal channel j was computed as follows:
where |${{X}_{mj}}$| and |${{S}_{mj}}$| represent the spectra of X and S for segment m and signal channel j, respectively. The stacked transfer function Gj for signal channel j was computed using the weight:
and its standard error |$\mathrm{ S}{{\mathrm{ E}}_j}\ $|of Gj was estimated as
The error bar is defined as |$\mathrm{ S}{{\mathrm{ E}}_j}$|. Using the weighted stacking method described above, we computed the transfer function of five-component CSEM data (Ex, Ey, Hx, Hy and Hz).
3 RESULTS AND DISCUSSION
3.1 Measured CSEM data
We deployed our EM-ACROSS transmitter with a 486 m NS dipole and a 961 m EW dipole at the Ishizu sulfur mine ruins (Fig. 3). The currents injected from the transmitter were ∼5 and ∼3 A for the NS and EW dipoles, respectively (Fig. 4a). This difference was due to the use of the constant voltage mode in the amplifiers and lower earth resistance for the NS dipole although |${{A}_j}$| were set as same for both dipoles (Table 1). The transmission time used in our analysis was 192 h, beginning at 17:00 on 2000 September 25, and ending at 17:00 on 2020 October 3 (Table 2). Note that the times in this study are indicated by Japan Standard Time (JST). We divided the time-series data of the transmitted current recorded at a sampling rate of 1 kHz into M segments of Ts = 600 s each. The time-series data of each segment were transformed to the frequency domain using FFT with a time window of 600 s. The spectra of the injected current S(f) exhibit peaks only in the signal channels (Fig. 4b). Moreover, during the 192 h, the largest phase changes relative to the phase values at their starting points were only 0.05° and 0.1° for the NS and EW dipoles, respectively (Figs 5a and b). The phase change is defined as |${{R}_j} - {{R}_1}$| where |${{R}_j}$| and |${{R}_1}$| are the phases at segment j and the first segment, respectively. These findings illustrate that our EM-ACROSS can repeatedly transmit accurately controlled signals for 192 h. The diurnal phase variation suggests that temperature effects in the transmitter electronics probably caused the phase variation. In contrast, the amplitudes decreased by 4 per cent relative to their initial levels (Figs 5c and d). The amplitude change is defined as |$100 \times ( {{{R}_j} - {{R}_1}} )/{{R}_1}$| where |${{R}_j}$| and |${{R}_1}$| are the amplitudes at segment j and the first segment, respectively. The automated meteorological data acquisition system reported that the accumulated precipitation for the first and second days (September 25 and 26) of transmission was 52 mm, while the accumulated precipitation for the remaining days from the third day (September 27) to the final day of transmission (October 3) was 1 mm. The amplitude changes at the transmitter were aligned with precipitation data and consistent across all frequencies and indicating that the near-surface resistivity changes caused by wetting and drying in the near-surface soils were responsible for the amplitude changes. Additionally, the amplitude changes might be contributed by the deterioration of electrodes due to corrosion over the 192 hr of continuous transmission.
Map of the study area. (a) Location of Kusatsu-Shirane Volcano, 200 km NW of Tokyo, marked by a yellow star. (b) Enlarged map of Kusatsu-Shirane Volcano with transmitter dipole positions (blue and red lines) and receiver positions (white triangles). Blue and red lines indicate NS and EW dipoles, respectively. The figures are plotted using GMT (Wessel et al. 2019).
Waveform of EM-ACROSS transmitter in the Kusatsu-Shirane Volcano. (a) Measured transmitting current as a time-series injected from NS and EW dipoles, 17:00 (JST) on 2020 September 25. The period of both waveforms is 10 s. (b) Amplitude spectra of the transmitting currents for 10 min from 17:00 (JST) on 2020 September 25. Blue and red represent values for NS and EW dipoles, respectively. NS current values is higher than EW values due to lower earth resistance for NS dipole and use of constant voltage mode of the amplifiers although |${{A}_j}$| in Table 1 are set same for both dipoles.
Stability of EM-ACROSS transmitting current over time. (a) Phase, (b) phase difference, (c) amplitude and (d) amplitude difference of the transmitted current for NS dipole at frequencies of 0.2 Hz (light blue) and 46.1 Hz (blue) and EW dipole at frequencies of 0.3 Hz (magenta) and 46.3 Hz (red).
Start and end times of transmitter and receiver data used for the analysis.
| Start time (JST) | End time (JST) | Duration (hr) | |
|---|---|---|---|
| Transmitter | 17:00 on Sep 25 | 17:00 on Oct 03 | 192 |
| Rx1 | 17:00 on Sep 26 | 14:00 on Oct 03 | 165 |
| Rx2 | 17:00 on Sep 25 | 14:00 on Oct 03 | 189 |
| Rx3 | 17:00 on Sep 25 | 14:00 on Oct 03 | 189 |
| Rx4 | 17:00 on Sep 25 | 14:00 on Oct 03 | 189 |
| Rx5 | 17:00 on Sep 25 | 17:00 on Oct 03 | 192 |
| Rx6 | 17:00 on Sep 25 | 17:00 on Oct 03 | 192 |
| Rx7 | 13:00 on Sep 27 | 15:00 on Oct 03 | 146 |
| Rx8 | 14:00 on Sep 26 | 17:00 on Oct 03 | 171 |
| Start time (JST) | End time (JST) | Duration (hr) | |
|---|---|---|---|
| Transmitter | 17:00 on Sep 25 | 17:00 on Oct 03 | 192 |
| Rx1 | 17:00 on Sep 26 | 14:00 on Oct 03 | 165 |
| Rx2 | 17:00 on Sep 25 | 14:00 on Oct 03 | 189 |
| Rx3 | 17:00 on Sep 25 | 14:00 on Oct 03 | 189 |
| Rx4 | 17:00 on Sep 25 | 14:00 on Oct 03 | 189 |
| Rx5 | 17:00 on Sep 25 | 17:00 on Oct 03 | 192 |
| Rx6 | 17:00 on Sep 25 | 17:00 on Oct 03 | 192 |
| Rx7 | 13:00 on Sep 27 | 15:00 on Oct 03 | 146 |
| Rx8 | 14:00 on Sep 26 | 17:00 on Oct 03 | 171 |
Start and end times of transmitter and receiver data used for the analysis.
| Start time (JST) | End time (JST) | Duration (hr) | |
|---|---|---|---|
| Transmitter | 17:00 on Sep 25 | 17:00 on Oct 03 | 192 |
| Rx1 | 17:00 on Sep 26 | 14:00 on Oct 03 | 165 |
| Rx2 | 17:00 on Sep 25 | 14:00 on Oct 03 | 189 |
| Rx3 | 17:00 on Sep 25 | 14:00 on Oct 03 | 189 |
| Rx4 | 17:00 on Sep 25 | 14:00 on Oct 03 | 189 |
| Rx5 | 17:00 on Sep 25 | 17:00 on Oct 03 | 192 |
| Rx6 | 17:00 on Sep 25 | 17:00 on Oct 03 | 192 |
| Rx7 | 13:00 on Sep 27 | 15:00 on Oct 03 | 146 |
| Rx8 | 14:00 on Sep 26 | 17:00 on Oct 03 | 171 |
| Start time (JST) | End time (JST) | Duration (hr) | |
|---|---|---|---|
| Transmitter | 17:00 on Sep 25 | 17:00 on Oct 03 | 192 |
| Rx1 | 17:00 on Sep 26 | 14:00 on Oct 03 | 165 |
| Rx2 | 17:00 on Sep 25 | 14:00 on Oct 03 | 189 |
| Rx3 | 17:00 on Sep 25 | 14:00 on Oct 03 | 189 |
| Rx4 | 17:00 on Sep 25 | 14:00 on Oct 03 | 189 |
| Rx5 | 17:00 on Sep 25 | 17:00 on Oct 03 | 192 |
| Rx6 | 17:00 on Sep 25 | 17:00 on Oct 03 | 192 |
| Rx7 | 13:00 on Sep 27 | 15:00 on Oct 03 | 146 |
| Rx8 | 14:00 on Sep 26 | 17:00 on Oct 03 | 171 |
Eight receivers were placed around the Yugama crater and at 4.5–6 km from the transmitter (Fig. 3). The times used for the data analysis vary among the receivers because of their installation and collection periods (Table 2). The shortest and longest times were 146 and 192 hr, respectively. The receiver time-series data, sampled at 150 Hz, were processed using the weighted stacking method. First, the time-series data were separated for each Ts = 600 s segment. The time-series data from each segment was then transformed to the frequency domain using FFT with a time window of 600 s. Consequently, a frequency resolution of the spectra |$\Delta f$| was 1/600 Hz. The noise levels in the signal channels were estimated using K = 50.
The noise levels in the observed data changed substantially over time. Fig. 6 shows fluctuations in the noise level for Ex at a frequency of 0.2 Hz at Rx6, which is 4.5 km away from the transmitter. The changes were periodic, with noise levels at night (18:00–04:00) being 5–10 times higher than those during the day. Numerous cabbage farms in the vicinity of the study site use electric fences as a deterrent to wild animals, typically activating them only during nocturnal hours. This nocturnal activation pattern likely contributes to the elevated noise levels observed at night. In addition to the periodic nighttime noise, strong pulse noise was detected on September 29 from 06:00 to 07:00. These fluctuations in the noise level over time imply that a weighted stacking method is required to address the noise level changes. Note that the CSEM signals before stacking were still within the noise at the frequency, thus the signal amplitudes were coherent with the noise level.
Time variability of Ex data and its noise level from Rx6. Black lines and magenta circles show the estimated noise level and observed Ex data (i.e. signal plus noise) at a frequency of 0.2 Hz, respectively. Electric field data are not normalized by the source dipole moment.
The SNR of the CSEM response obtained using the weighted stacking method improved with increasing stacking time. For the stacking periods of 3 and 20 hr, the noise levels were greater than the signal levels at frequencies of 0.2 and 0.3 Hz of Ex at Rx6 (Fig. 7a). In contrast, stacking for 192 hr significantly lowered the noise level and yielded the SNR of ∼6 for data at 0.2 and 0.3 Hz. SNRs for data at 46.1 and 46.3 Hz also improved with increasing stacking time (Fig. 7b). The stacking period of 192 hr yielded a high SNR of >50 at frequencies of 46.1 and 46.3 Hz. This result demonstrates that our EM-ACROSS could significantly increase the SNR by stacking data with a sufficiently long observation time, as expected based on the system design.
Increasing SNR during the stacking period. Amplitude spectra of Ex at Rx6 obtained by stacking for 3 (left), 20 (centre), and 192 h (right); (a) and (b) show different frequency ranges of 0.1–0.4 and 46.0–46.4 Hz, respectively. Blue squares indicate signals at 0.2 and 46.1 Hz from the NS dipole; red squares show signals at 0.3 and 46.3 Hz from the EW dipole. Results are obtained using the weighted stacking method.
Given that noise is Gaussian with zero mean, SNR increases as the square root of the stacking number. Fig. 8 shows the relationship between SNR and stacking number for Rx6 at frequencies of 0.2 and 0.3 Hz. SNR of the five-component CSEM data has improved with increasing stacking numbers for both conventional stacking and weighted stacking methods. The conventional stacking method refers to using the same weights for all segment data. The improvement in SNR is more significant with weighted stacking compared to conventional stacking. The improvement is attributed to appropriately setting weights based on noise levels. The SNR improvement roughly followed the square root of the stacking times. This agreement suggests that if the noise level is measured at the beginning of the observation, the user can optimize the stacking time based on the required SNR level.
SNR improvement through stacking for Rx6. Blue solid and dashed lines represent the SNR at a frequency of 0.2 Hz for the NS transmitting dipole obtained through weighted and conventional stacking, respectively. Similarly, the red solid and dashed lines represent the SNR at a frequency of 0.3 Hz for the EW transmission. The black dotted line indicates the values scaled by the square root of the stacking times from the initial SNR. Vertically, from top to bottom, the panels represent those of for Ex, Ey, Hx, Hy and Hz.
As the weighted stacking method obtained higher SNR than the conventional stacking method (Fig. 8), the weighted stacking method yielded CSEM responses with fewer errors compared to the conventional stacking method (Figs 9 and 10). We compared Ex data at Rx6 computed using both the weighted and conventional stacking methods (Fig. 9). Both methods utilized a stacking time of 192 hr. Error bars for Ex amplitude at 46.1 Hz were 2.6 per cent for the weighted method and 21 per cent for the conventional method (Fig. 9a). Additionally, the weighted stacking method improved response estimations for phases (Fig. 9b). Enhancements achieved by the application of the weighted stacking method are discernible across all five components and receivers (Fig. 10). These results suggest that utilizing EM-ACROSS with a long stacking duration in conjunction with the weighted stacking method enables accurate estimation of the CSEM responses.
Comparison of the weighted and conventional stacking methods for Ex at Rx6. (a) Amplitude and (b) phase obtained using the weighted stacking and conventional stacking methods. Blue and red circles indicate responses obtained using the weighted stacking for NS and EW dipoles, respectively. Black inverse and normal triangle indicate responses obtained using the conventional stacking for NS and EW dipoles, respectively.
Comparison of phase errors obtained by weighted stacking and conventional stacking methods. Horizontally, from left to right, the panels display error bars of phase data for Ex, Ey, Hx, Hy and Hz. Vertically, from top to bottom, the panels represent those of Rx1 through Rx8.
We compared the MT responses at Rx4 acquired with and without EM-ACROSS transmission running to demonstrate that MT data can be collected concurrently with EM-ACROSS data. Rx4 was deployed three days before starting EM-ACROSS transmission. Therefore, its recordings include data both during and before EM-ACROSS transmissions. The gain and resolution settings same as for typical MT measurements were applied to MTU-5C receivers for both data acquisition with and without EM-ACROSS transmission running. The processing time for both MT data without (06:00 on September 23 to 06:00 on September 24) and with EM-ACROSS transmissions (06:00 on September 29 to 06:00 on September 30) was 24 hr. We employed the EMpower software (Phoenix Geophysics Ltd) to process the MT data, whose frequencies were selected to avoid the EM-ACROSS signal frequencies. Note that MT measurements include burst cycles of higher acquisition rates in addition to 150 Hz sampled data. Fig. 11 shows the consistency of the MT data acquired with and without EM-ACROSS transmission, indicating that the MT data can be collected concurrently. This gives EM-ACROSS an advantage over other CSEM systems that interfere with simultaneous MT data collection. The apparent resistivity data significantly differed between frequencies ranging from 0.1 to 3 Hz, corresponding to the MT dead band with a weak MT source (Egbert & Livelybrooks 1996). We infer that the deviation arises from low SNR for either the observation time with and without EM-ACROSS transmission or both. In contrast, stacking for 192 hr of the EM-ACROSS data yielded low-error CSEM responses at frequencies of 0.1–3 Hz (Fig. 12), suggesting that a CSEM survey using the EM-ACROSS is effective in capturing low-error responses in an environment where obtaining low-error MT responses is challenging.
Comparison of the MT data at Rx4 acquired with and without EM-ACROSS transmission running. (a) Apparent resistivity and (b) phase of the MT data without (square) and with (cross) EM-ACROSS transmissions. Red and blue show the xy and yx components of the impedance tensor, respectively.
Data fitting for Rx4. (a) Amplitude and (b) phase of the five-component CSEM data (Ex, Ey, Hx, Hy and Hz). Points with error bars and lines indicate the observed data and the response from the inverted model in Fig. 13(a), respectively. Error bars are defined as standard errors of the data. Blue and red represent data for NS and EW transmitter dipoles, respectively.
3.2 Inverted resistivity structure
To obtain the 3-D resistivity structure, we inverted the CSEM data at the eight receivers using the inversion code with the data-space Occam algorithm developed by Ishizu & Ogawa (2021) and Ishizu et al. (2022b). Forward modelling of the inversion code employed the finite difference method with a primary-secondary field approach. To avoid numerical singularities at the source positions, the total EM fields were divided into primary and secondary fields. The primary fields were analytically computed from a 1-D primary model.
The inversion code received five-component EM data (Ex, Ey, Hx, Hy and Hz) as log10-scaled amplitude and linear-scale phase values. The primary, initial and prior models consisted of a 10 Ωm subsurface layer and a 108 Ωm air layer above the ground surface. The models were divided into 80 × 70 × 81 cells in x-, y- and z-directions. x = 0 and y = 0 m were set at the centre of the receivers, whereas z = 0 m was set at the ground surface. For horizontal cells, 50 m size was used in a survey area of –1.0 km < x < 1.0 km; –0.8 km < y < 0.8 km, and the cell size outside the survey area was exponentially increased to the boundary positions. For vertical cells, the smallest size of 5 m was used for the ground to a depth of 100 m, and the size gradually increased toward the bottom. A numerical singularity develops if the primary and inverted resistivity models vary at the source positions. To avoid this, we fixed the resistivity values at a specific location around the source sites to 10 Ωm; this value was the same as the resistivity of the subsurface layer of the primary model. Other CSEM inversion studies also employed comparable approaches to prevent the occurrence of numerical singularities (Grayver et al. 2014; Peng et al. 2023). As the transmitter and receiver offsets were more than five times the source dipole lengths, we modelled the transmitter dipole as a point source, resulting in minor numerical errors (Streich & Becken 2011). We applied a minimum error setting of 5 per cent to the data to avoid overfitting the inversion to the observed data. The root mean square (rms) misfit of the initial model was 19.3. After 27 iterations, the inversion yielded a model (Fig. 13a) with an rms misfit of 1.5 that fits observed data (Figs 12, 14 and Figs S1–S7, Supporting Information).
Inverted resistivity models from CSEM and MT data of the Kusatsu-Shirane Volcano. (a) CSEM model and (b) MT model. White triangles indicate receiver positions. C1c–C2c and R1c represent conductors and resistors in the CSEM model, respectively; C1m–C3m represent conductors in the MT model
Weighted residuals for each component across all sites and frequencies. Weighted residuals are defined as the difference between measured data and the model response, normalized by the standard error. The panels are arranged in the same order as in Fig. 10.
The resulting resistivity model specified conductive and resistive structures related to volcanic activity (Fig. 13a). A conductor C1c below the Yugama crater is consistent with a conductor, previously specified by an audio-frequency MT (AMT) survey (Nurhasan et al. 2006; Tseng et al. 2020). Tseng et al. (2020) interpreted the conductor as a low-permeability smectite layer that functions as a cap-rock layer for underlying hydrothermal fluids. Here, we follow the interpretation of Tseng et al. (2020) for C1c. The cap-rock layer C1c is thin at the centre of the Yugama crater and thick towards the east. We infer that smectite–illite transformation (∼200° C) takes place deeper towards the east than the centre of the Yugama crater due to lower subsurface temperature. Our model revealed a conductor C2c beneath the Yugama crater and the western edge of C1c. The significantly low resistivities (∼0.5 Ωm) indicate that C2c is a hydrothermal reservoir with a high concentration of chloride. Geochemical study in the Kusatsu-Shirane Volcano (Ohba et al. 2000) suggests that the chloride concentration resulted from the input of magmatic fluids. Tseng et al. (2020) also reported conductive structures (C2 and upper part of C3 underlying C2 in Tseng et al. 2020) at the C2c location. They interpreted C3 as a brine reservoir and C2 as the fluid path upwelling from the C3 reservoir, associating C2 with the microseismic hypocentres. Thus, C2c in our CSEM model is also linked to upwelling fluids observed as microseismic activities. Meanwhile, a seismic tremor study suggested that a two-phase liquid–gas reservoir exists at a depth of 300 m below the Yugama crater (Nakano et al. 2003); the resistive structure of R1c is imaged there. We interpret that R1c is a two-phase liquid–gas reservoir underlying the cap-rock layer. Due to the high-resistivity values, we infer that R1c is a vapour-dominated two-phase liquid–gas reservoir. Of note, gas-phase resistivity is typically much higher than that of liquid phase (Milsch et al. 2010; Watanabe et al. 2021).
The MT data at the eight receivers were inverted to compare the CSEM and MT models. We used WSINV3DMT (Siripunvaraporn & Egbert 2009) for the MT data inversion. The initial and prior models, as well as model discretization for the MT inversion, were the same as those used for the CSEM inversion. The MT data comprised impedance tensors and tippers at eight frequencies ranging from 0.1 to 100 Hz; this frequency range was similar to that of the EM-ACROSS data. We applied the minimum error settings of 5 per cent of |$\left| Z_{xy} Z_{yx} \right|^{1/2}$| for the impedance tensor; 10 per cent of each component for the tipper. |${{Z}_{xy}}$| and |${{Z}_{yx}}$| represent xy and yx components of the impedance tensor. After two iterations, the inversion yielded the model (Fig. 13b) with the rms misfit of 3.5.
The MT model (Fig. 13b) aligns with the CSEM model but displays distinct features. The conductive body C1m in the MT model corresponds to the cap-rock layer C1c identified in the CSEM data. However, the MT model reveals the differing characteristics below the cap-rock layer: resistor R1c is detected in the CSEM model, whereas conductor C3m is identified in the MT model. This discrepancy may be due to the higher CSEM sensitivity to a resistor than the MT data (Ishizu & Ogawa 2021). To demonstrate the sensitivity difference, we present the near-field and far-field zones for the impedance of the observed CSEM data (Fig. 15). In the far-field zone, the CSEM data align with the MT data due to the valid plane-wave assumption. However, the presence of vertical electric current in the near-field zone causes the CSEM data sensitivity to differ from the MT data. The near-field effect is notable at frequencies <3 Hz in the observed data (Fig. 15). This effect can be supported by the synthetic response of the 10 Ωm half-space model at frequencies <3 Hz: apparent resistivities deviate from 10 Ωm; phases of xy and yx components shift from 45° and –135° to 0° and –180°, respectively (Fig. 15). Because of these characteristics, CSEM data are more sensitive to a resistive body of a vapour-dominated reservoir than MT data.
Apparent resistivity and phase of impedance of CSEM data at Rx4 to show near and far field zones. (a) Apparent resistivity and (b) phase of impedance of the observed CSEM data (solid line) and the synthetic response of the 10 Ωm half-space model (dashed line). Red and blue lines show the xy and yx components of the impedance tensor, respectively.
Additionally, while our CSEM inversion identifies the conductor C2c, our MT inversion detects a moderately low-resistivity zone C2m to the west of C2c. A previous study using the inversion of AMT data from 91 sites also reported conductors (C2 and upper part of C3 underlying C2 in Tseng et al. 2020) at the C2c location; and resistivity values of C2 and C3 in Tseng et al. (2020) are more consistent with C2c than C2m. Given that the inversion of AMT data from the 91-sites yielded a reliable model for these conductors, the comparison indicates that the CSEM data recovered the conductor C2c more accurately than the 8-site MT data. We infer that this improved delineation is a result of lower errors in the CSEM data compared to the MT data. However, it is important to note that MT data using frequencies <0.1 Hz can explore deeper depths than CSEM data. The joint inversion of CSEM and MT data has the potential to enhance inversion results by filling in low-sensitivity zones in each data set.
The cap-rock layer and the hydrothermal reservoirs are crucial structures for the hydrothermal system in the Kusatsu-Shirane Volcano. Inversion of the EM-ACROSS data delineate them, indicating its effectiveness for investigating volcanic areas. The cap-rock layer and hydrothermal reservoirs are associated with the occurrence of phreatic eruptions in the Yugama crater (Ohba et al. 2000). Their monitoring is valuable for understanding the mechanism of phreatic eruptions and forecasting their likelihood. The EM-ACROSS can be used for monitoring as it can repeatedly transmit an accurately controlled signal over extended measurement intervals. We plan to use the EM-ACROSS to elucidate the temporal changes in the cap-rock layer and hydrothermal reservoirs.
4 CONCLUSIONS
We conducted a CSEM survey using the EM-ACROSS in the Kusatsu-Shirane Volcano to demonstrate that stacking can improve SNR and CSEM data inversion can delineate hydrothermal system including resistive bodies of vapour-rich reservoir. The EM-ACROSS transmitted waveforms accurately controlled by synchronizing with a 10 MHz signal from the GPS clock from 0.5 km NS and 1.0 km EW dipoles. Accurately controlling waveforms helps efficiently stack long-term data as the spectra of transmitting currents only shows peak values at the signal frequencies without dispersion to other frequencies. The CSEM survey, with a long stacking duration of 146–192 h paired with the weighted stacking method, yielded high SNR > 5 at five-component receivers placed 4.5–6 km from the transmitter. The SNR increased with stacking time, scaling approximately with the square root of stacking time. MT data could be collected concurrently with EM-ACROSS data because the EM-ACROSS transmitted signals with a limited number of frequencies without frequency leakage and our receivers used gain and resolution settings same with typical MT measurements. This is an advantage of the EM-ACROSS over other CSEM systems that interfere with simultaneous MT data. We applied 3-D inversion to the CSEM data from eight receivers. The inverted resistivity model specified a cap-rock layer and hydrothermal reservoirs below the cap-rock layer. The high- and low-resistivity regions estimated in the hydrothermal reservoirs were interpreted as a vapour-dominated two-phase reservoir and a liquid-dominated hydrothermal reservoir with high-chloride contents, respectively. The MT inversion yielded a model that aligned with the CSEM model but failed to recover the resistive body of vapour-dominated two-phase reservoir detected by the CSEM inversion. Although we demonstrate the effectiveness of the EM-ACROSS by applying it to the Kusatsu-Shirane Volcano, the high sensitivity to thin resistive bodies and noise resistance indicates that it is also useful for exploring various targets including hydrocarbons, metal deposits and geothermal resources.
CONFLICT OF INTEREST
The authors declare that this study was conducted in the absence of any commercial or financial relationships that could be construed as conflicts of interest.
ACKNOWLEDGMENTS
This study was supported by the Japan Society for the Promotion of Science, KAKENHI (grant numbers 20H01992, 21KK0081, 22K14104, and 23K17803), ERI JURP 2019-G-11 in Earthquake Research Institute, the University of Tokyo, and the Ministry of Education, Culture, Sports, Science and Technology (MEXT) of Japan, under its Earthquake and Volcano Hazards Observation and Research Program. We would like to thank Dr Amir Haroon and the anonymous reviewer for their detailed and helpful comments, which led to significant improvements in this paper and Dr Weerachai Siripunvaraporn for providing WSINV3DMT.
AUTHOR CONTRIBUTIONS
KI drafted the manuscript. YO supervised the study. KI conducted the CSEM data processing. KI and YO interpreted the data. TK designed the fundamental aspects of the transmitter system, while HT established its basic operation. YO, KI, HT and TK set up the system for the CSEM survey used in this paper. All authors participated in field surveys. All authors have read and approved the final version of the manuscript.
DATA AVAILABILITY
The processed CSEM data and the inverted results are available on the Zenodo open-access repository: 10.5281/zenodo.14245638.
REFERENCES
Author notes
Now at: Department of Earth Resources Engineering, Kyushu University, Fukuoka, 819-0395, Japan.
Now at: Sumiko Resources Exploration & Development Co. Ltd, Tokyo, 105-0001, Japan.
