© ESO 2018

1. Introduction

The detection of g modes remains a major objective of helioseismology. The benefit of detecting these modes would be to obtain the structure and dynamics of the very inner core of the Sun. There have been several claims of g-mode detection (see Appourchaux et al. 2010, for a review). Recently, Fossat et al. (2017) using the propagation time of the p-mode wave packet claimed to have detected the signature of g modes. In order to test that detection claim, we made longer data sets using a new calibration strategy for the GOLF data.

Since the beginning of helioseismology, solar radial velocities have always been measured using solar spectral lines. The intensities are typically measured in the blue and red wings (I_b and I_r) of the line and the displacement is deduced by calibrating the ratio $\frac{I_{\mathrm{r}}-I_{\mathrm{b}}}{I_{\mathrm{r}}+I_{\mathrm{b}}}$ with respect to known radial velocities (see Elsworth et al. 1995, and references therein). The purpose of the ratio is mainly to remove the effect of the Earth’s atmospheric variations. Concerning space-based instruments, the ratio is used for reducing the effect of the change of transmission in the course of the lifetime of the instrument. For ground-based instruments, the ratio is measured at a very high cadence (from tens of Hz to a few kHz) while for space-based instruments, a slow cadence can be used (slower than 1 Hz). The signal inferred from the ratio is somewhat affected by signals not related to radial velocities, but due to the effect of radiative transfer across the line (Ulrich et al. 2000).

The calibration of the GOLF data was rendered more complicated by the fact that the two-wing measurement could not be made after 31 March 1996 (see Gelly et al. 2002). This introduced an additional effect on the p-mode velocity signal since on either wing the intensities are intrinsically modulated by intensity fluctuations due to the p modes themselves and the granulation background. The fraction of the intensity fluctuations to velocity fluctuations due to p modes is about 0.12 (Renaud et al. 1999). The measurement of the GOLF instrument was then done using only one wing (red or blue) using three different methods (Ulrich et al. 2000; Gelly et al. 2002; García et al. 2005). Ulrich et al. (2000) relied on a detailed modelling of the line profile for inferring the residual velocities. Gelly et al. (2002) tried to minimise the yearly modulation of p-mode power to optimise the calibration. The change of p-mode power from 1996 to 2002 due to solar activity was 10%, typical of what is observed by the Birmingham Solar Oscillation Network¹ (see Howe et al. 2015). García et al. (2005) relied on instrumental calibrations and a non-linear calibration method developed by Pallé et al. (1993).

The recent g-mode detection claim of Fossat et al. (2017) was performed using the calibration done by García et al. (2005). There are two types of data available on the GOLF website²: those used by Fossat et al. (2017) sampled at 80 s, and a time series sampled at 60 s. Both data have the same length and are produced using the same time series sampled at 20 s but binning over 4 and 3 samples, respectively (García 2018 priv. comm.). In the course of reproducing the findings of Fossat et al. (2017), we used both time series to produce our own version of their Figs. 10 and 16. Figure 1 shows the output of the procedure used by Fossat et al. (2017) on two different times series. It is clear that Figs. 10 and 16 of Fossat et al. (2017) cannot be reproduced with the time series sampled at 60 s. Very recently, Schunker et al. (2018) reproduced Fig. 10 of Fossat et al. (2017). They showed that using different fitting procedures, the prominent peaks at 210 nHz and its acolytes would smear out or even disappear, as seen for the time series sampled at 60 s. This leads us to investigate how the GOLF data are calibrated and whether a different calibration might or might not reproduce the results of Fossat et al. (2017).

Fig. 1. Left: correlation of the power spectrum as obtained by Fossat et al. (2017) as a function of frequency lag for two different time series sampled at 80 s (top) and 60 s (bottom); this is comparable to Fig. 10 of Fossat et al. (2017). The green vertical lines correspond to frequencies at 210 nHz, 630 nHz, and 1260 nHz. Right: sum of the correlation for l = 1 and l = 2 modes as obtained by Fossat et al. (2017) as a function of rotation frequency for two different time series sampled at 80 s (top) and 60 s (bottom); this is comparable to Fig. 16 of Fossat et al. (2017).

This work has been divided in two parts. The present paper (Part I) explains how a new calibration of the GOLF data has been obtained. The accompanying paper (Part II) investigates the results obtained by Fossat et al. (2017) and compares it with different time series. The first section of this present paper explains the calibration procedure, the extraction of the velocity, and the time corrections applied. The second section compares the calibration results obtained with other data sets. Subsequently, we present our conclusions.

2. A new calibration of the GOLF data

2.1. Extraction of the velocity

The data used for the new calibration start on 11 April 1996 and end on 10 April 2018. The technique used for the calibration is derived from that of García et al. (2005) using the so-called X method. Since GOLF operates only in one wing of the Sodium lines, the following ratio is computed as a proxy to the two-wing signal:

$$X=\frac{1}{2}\frac{P^{+}+P^{-}}{\langle P^{+}-P^{-}\rangle },$$(1)

where P⁺ and P⁻ are the signals measured in the same wing (blue or red) using a weak magnetic field modulation that induces a wavelength change. The bracket denotes a low-frequency filtering of the slope proxy (i.e. the derivative). The division by 〈P⁺ − P⁻〉 is needed since it provides a way to normalise to the slope of the line. The procedure for computing the radial velocity starts by computing the photomultiplier (PM) signals as follows:

1. Read daily data for both PMs and magnetic modulation (4 signals as $P_{1,2}^{+}$, $P_{1,2}^{-}$; see Gabriel et al. 1995, for details).

2. Compute daily median of the slope proxy ($P_{1,2}^{+}-P_{1,2}^{-}$) for low frequency filtering.

3. Either no time correction is performed or if a time correction is measured, shift daily data by integer multiples of 10 s (no interpolation is used to avoid introducing additional noise).

4. Concatenate all of the daily data and save the the four signals and the slope proxy.

We note that the time correction in Step 3 is only made in a second iteration after the time shift measurement is performed on the unshifted time series. The procedure for the time shift measurement is explained in Sect. 2.2. In a second stage, we compute the X_1,2 ratio as follows:

1. Bin the data originally sampled at 10 s over two samples resulting in a 20 s cadence.

2. Using a spline, interpolate the daily slope proxy over the 20 s sample to a uniform temporal sampling.

3. Compute the X₁ and X₂ ratios according to Eq. (1).

4. Cut the time series of the X₁ and X₂ ratios in sub-series of 20 days.

5. Interpolate to remove one-sample spikes greater than 11% of the median of the sub-series.

6. Compute residuals of the sub-series with respect to a two-day smoothed version of the sub-series.

7. Detect variations greater than 1% in the residuals then replace the remaining outliers by interpolation and remove the slow changes using a sixth-order polynomial.

The cut-offs used for filtering the outliers are a compromise between the level of the outliers and the signals to be kept; this could be perceived as arbitrary but in fact this procedure was chosen according to our experience with the data. The calibration of the X₁ and X₂ ratios with respect to the solar radial velocity is done using the technique described in García et al. (2005), in which they assume that the ratio X can be modelled as follows:

$$X=f(V)=\frac{a+bV}{1+\alpha V^2+\beta V^4},$$(2)

where a, b, α and β are parameters, and V = V_grav + V_orb is the radial velocity including the gravitational redshift (V_grav = 636 m s⁻¹) and the contribution generated by the spacecraft orbiting the Sun (V_orb). The residual velocity (V_X, made of the oscillation signals and other contributions) is excluded from the model since this is what we want to extract. The orbital spacecraft velocities are obtained from the VIRGO Data Center, which provides computation of the solar radial velocity on a 60 s cadence. The velocity used here is the predictive, not the measured, velocity since the SoHO operation team stopped providing the reconstructed velocity after 1998 because the predictive one was good enough. The spacecraft velocity is then interpolated onto the sample series (20 s cadence). In theory, one would invert Eq. (2) to obtain the residual velocity (V_X) as V_X = f⁻¹(X) − V but this is not done in practice. The residual velocity is inferred as in García et al. (2005) using the following computation of the residuals:

$$V_{\mathrm{X}}=\frac{1}{b}\left(X(1+\alpha V^2+\beta V^4)-a\right)-V.$$(3)

This so-called inversion performs better than any other true inversion obtained by solving Eq. (2). The true inversion provides yearly modulations of the p-mode amplitude that are not physical. We also tried using the calibration procedure of Elsworth et al. (1995) consisting in fitting a polynomial to the ratio X as a function of V, then using the derivative to obtain the velocity residuals, but this scheme also produces the same yearly modulation. The reason for the exceptional performance of Eqs. (2) and (3) is not yet fully understood. The calibration of the velocity is done for three time segments, each corresponding to a different operational mode of GOLF, as in García et al. (2005). The ratio is fitted in two passes; on the second pass, residuals larger than 150 m s⁻¹ are excluded from the fit (about 0.2% of the points are excluded). A fitted ratio is shown for the blue wing in Fig. 2 for PM1. We compute the residual using Eq. (3). Then on that residual, we detect spikes greater than 12 m s⁻¹ by applying a high pass filter based on a two-day triangular smoothing; about 0.15% of the points are rejected. Table 1 gives the fitted parameters according to Eq. (2) for the ratios X₁ and X₂. The final velocity is the average of the velocity residuals of the times series of PM1 and PM2. Figure 3 gives the final residual velocity for PM1. There is still a yearly modulation of the residual velocity that is related to the variation of the temperature along the orbit which is not properly taken into account. Nevertheless this residual modulation has no impact on the p-mode amplitude (see Sect. 3).

Fig. 2. X ratio as a function of the velocity for PM1: one point per day. The fit of the X ratio is in green. The fit is for the blue wing after February 2002.

Fig. 3. Velocity residuals for PM1 as a function of time: one point per day. The two vertical lines indicate when the modes of operation changed from blue wing to red wing and vice versa.

2.2. Correction of time shifts

As outlined in García et al. (2005), the GOLF data may lose synchronisation with respect to the basic temporal cadence. In principle, ancillary data allow us to quantify the resulting induced shifts thanks to a daily pulse (DP) generated on board SoHO. The precise time datation according to the Temps Atomique International (TAI) is available in the header of the raw daily data files. However, in a few cases, the GOLF clock could not be synchronised on the DP, and in order to check for potential unknown time jumps, we used data from other instruments such as the Global Oscillation Network Group³, the Birmingham Solar Oscillation Network⁴, and the Sun Photometer Monitor of the Variability of Irradiance and Gravity Oscillations instrument⁵. The measurement of the time shift is done using the p modes themselves with the following procedure:

1. Select a sub-series from two instruments (e.g. GOLF, BiSON, GONG, SPM) with a typical duration of one day. When the cadence differs, resample the sub-series to a common temporal cadence of 60 s using linear interpolation (e.g. for BiSON). The integration is one day as well as the sampling (for display purposes, other integration times of 10 days and 45 days have been used with one-day sampling).

2. Compute the backward difference filter (BDF)⁶, for each instrument but GONG, for which it has already been applied.

3. Compute the cross correlation C for a range of ±13 mins by steps of 1 min.

4. Compute the cross correlation envelope C_env using the Hilbert transform $\mathscr{H}(C_{\mathrm{env}}=\sqrt{C^2+\mathscr{H}{(C)}^2}$)⁷.

5. Fit the cross correlation envelope using least squares with a Gaussian to obtain an estimate of the time shift between the two data sets. This estimation is used for the following step.

6. Fit the cross correlation using least squares with a time-shifted cosine function modulated by a Gaussian envelope. The Gaussian envelope and the cosine function have the same time shift.

The BDF is used to reduce the low-frequency solar noise especially when correlations are made between intensity and velocity signals. Since GONG uses the BDF by default⁸, this filter is not applied to GONG. The precision of the time shift obtained by fitting the correlation of the p modes is far more precise than the fitting of the correlation provided by the envelope. The former is related to the p-mode phase velocity while the latter is related to the p-mode group velocity.

Figure 4 gives an example of what is fitted. We found that the GONG data for l = 0 and the VIRGO/SPM data do not present any temporal jumps as shown in Fig. 5. As it is unlikely that both time series would have time jumps at the same date, the blue SPM and GONG then provide two references for a constant time base.The GONG data were then used as a time reference to check and correct GOLF datation. On close inspection, Fig. 5 shows that the time delay in this time series is modulated with a periodicity of almost 6 months. This is due to the SoHO halo orbit⁹, which has a periodicity of 178 days resulting in a time modulation of 1.54 s. The periodicity of 178 days is primarily constrained by the solar radiation pressure (Chidambararaj & Sharma 2016).

Fig. 4. Cross correlation between the first difference of the GOLF radial velocity and the first difference of the l = 0 GONG radial velocity as a function of lag in minutes (black line), with the envelope shown (orange line) and the fitted cross correlation (green line).

Fig. 5. Time difference between the SPM Blue/VIRGO data and the l = 0 GONG data as a function of time computed using a 10-day window (black dots) or a 45-day window (green dots).

Figure 6 shows the measured time delay of the unshifted GOLF time series together with the time delay of the corrected time series using the shifts listed in Table 2. The measured time shifts as given in Table 2 are used for correcting the time series after the listed dates of the table. We used an integration time of one day to find the exact dates of the jumps. After correction, the time series are checked again for potential errors in the date of the jumps. We note that the intrinsic time shift of about 12 s measured by Renaud et al. (1999) between the blue wing and the red wing is automatically corrected by our procedure. We also outline that the TAI time retrieved from the header of the file does not always reflect the presence of a time jump, and even sometimes the time jump is detected before the exact date provided by the change deduced from the TAI time. As for SPM, we can see the six-month modulation in the corrected time series even more clearly in Fig. 6.

Fig. 6. Time difference between the newly calibrated GOLF data and the l = 0 GONG data as a function of time; (top) unshifted GOLF data; (bottom) shifted GOLF data according to Table 2 computed using a 10-day integration window (black dots) or a 45-day integration window (green dots). The median value of the time shift for the corrected time series is 47.3 s. About 99.8% of the corrected values are within 5 s of this median value.

The typical rms precision obtained is of the order of 0.5 and 2.5 s for 10 days of integration, for GOLF and SPM, respectively. Figure 7 shows the time delay of the time series used by Fossat et al. (2017) with respect to GONG and to the newly corrected GOLF time series. There are seven time shifts not corrected in the time series used by Fossat et al. (2017). When the time shift agrees between the two versions of the GOLF time series, the resulting mean time difference is 0.4 ± 42 ms.

Fig. 7. Top: time difference between the GOLF time series sampled at 60 s and the GONG l = 0 data as a function of time. Bottom: time difference between the GOLF time series used by Fossat et al. (2017) and the newly calibrated GOLF data as a function of time, both sampled at 80 s.

3. Comparison of the calibration with other data sets

In order to double check the velocity calibration and time correction, we used data from GONG l = 0 and from BiSON¹⁰ and compared them with the previous calibration of the time series used by Fossat et al. (2017). Figure 8 shows the comparison of the rms p-mode amplitude and the rms p-mode noise for these different time series. The computation of the rms amplitude and noise was done as follows.

Fig. 8. Top: rms mode amplitude as a function of time for the newly calibrated GOLF data (black), for the GOLF data as used in Fossat et al. (2017; orange), for the BiSON data (red), and for the GONG data (green). Bottom: rms noise as a function of time. The sub-series are 30 days in length.

1. Select a sub-series of the data set, then compute the power spectrum. Typical duration of the sub-series is 30 days.

2. Compute the integrated power (P_T) between 2500 μHz and 3500 μHz normalised by the integration factor of the sampling window (i.e. divide the power spectrum by [sinc(πντ)]², where τ is the integration time of a sample). This provides a proxy of the power of the sum of the noise power and the p-mode power.

3. Compute the integrated power between 1000 μHz and 1500 μHz normalised by the integration factor of the sampling window. For providing a proxy of the noise power in the p-mode range, we scaled this power by a factor taking account of the frequency bandwidth (factor 2) and extrapolation factor (factor $1/\sqrt{2}$); we then multiply this power by $\sqrt{2}$ to obtain the proxy (P_noise).

4. The final p-mode power is obtained as P_p-mode = P_T − P_noise.

5. The rms amplitude is then $A_{p\mathrm{-mode}}=\sqrt{P_{p\mathrm{-mode}}}$, while the rms noise is $A_{\mathrm{noise}}=\sqrt{P_{\mathrm{noise}.}}$

6. Apply an ad-hoc formation height correction for the amplitudes (if required).

We apply this procedure to all data but GONG. For GONG, additional corrections of the power spectrum were required because of the BDF used, which requires us to divide by the filtering factor T(ν) = 4 sin²(πνΔt), where Δt is the sampling cadence. For the p-mode frequency range and a 60 s cadence, the filtering factor in power ranges from 0.8 to 1.5, while for the noise it ranges from 0.14 to 0.31.

For the last step, the formation height was corrected for GOLF during the red-wing mode by multiplying the amplitudes by 1.11. For GONG, the resulting velocities were multiplied by 1.17 in order to take into account the different formation height of GONG; the factor was adjusted post-facto. No other correction for formation height was applied to the remaining data sets.

We note two different amplitude regimes in the time series used by Fossat et al. (2017), occurring at the beginning of the time series and after 2002 for about 2 yr, both in the blue wing. These two differences are related to the way that the temperature of the cathode of the photomultipliers is used to correct the calibration (Garcia, 2018, priv. comm.). We also note the variation of p-mode amplitude with the solar activity cycle with a typical variation of 9% between solar maximum and solar minimum, as already measured by Howe et al. (2015) with the BiSON data. Figure 8 indeed shows that a calibration that does not take into account the calibrated instrumental characteristics performs better in terms of a more consistent p-mode amplitude. The most likely reason for this better performance may lie in the stability of the thermal environment provided after November 2002; thereby not requiring additional instrumental corrections. Closer inspection of the p-mode amplitude also reveals that there is no yearly modulation present when GOLF observes in the blue wing. On the other hand, a modulation is clearly seen when the observation is done in the red wing, which is due to the fact that the p-mode temperature fluctuations have a much larger contribution in that wing.

It is also interesting to compare the noise performance of the different instruments. Again there is a clear difference for the GOLF noise between the two time series at the beginning of the time series (in the red wing) and after 2002 for about 2 yr. We also note that the GOLF noise does indeed increase with time primarily due to the photon noise increasing by a factor of five due to the number of photons dropping by a factor of 25. Here, contrary to the p-mode amplitude, there is a clear yearly modulation of the noise in the GOLF data and in all data sets.

4. Conclusion

The advantage of the new calibration is that it requires only very limited knowledge of the instrument. Although it may seem contradictory at first, the calibration is more in line with other instruments that are very well calibrated. On the other hand, some aspects remain to be fully explained. For instance, the yearly modulation in the velocity apparent in Fig. 3, which is possibly related to thermal effects induced by the variation of the SoHO–Sun distance throughout the orbit. The current time series is also corrected for time shifts to keep 99.8% of the corrected time series within ±5 s. There is a residual periodic time modulation due to the halo orbit, which is about 178 days. In theory, the g-mode detection by Fossat et al. (2017) is not affected by either the changing amplitude of the p modes or the measured time shifts. In practice, this new calibration approach may remove any doubts in the results of Fossat et al. (2017) related to these changes and time shifts. This work will serve as a basis for an analysis of longer time series for testing the g-mode discoveries made by Fossat et al. (2017).

The new calibrated series are available as FITS files at the CDS and online¹¹.

Acknowledgments

The GOLF instrument onboard SoHO is based on a consortium of institutes (IAS, CEA/Saclay, Nice and Bordeaux Observatories from France, and IAC from Spain) involving a cooperative effort of scientists, engineers, and technicians, to whom we are indebted. GOLF data are available at the MEDOC data and operations centre (CNES/CNRS/Univ. Paris-Sud) at medoc.ias.u-psud.fr. SoHO is a mission of international collaboration between ESA and NASA. This work uses data obtained by the Global Oscillation Network Group (GONG) program, managed by the National Solar Observatory, which is operated by AURA, Inc. under a cooperative agreement with the National Science Foundation. The data were acquired by instruments operated by the Big Bear Solar Observatory, High Altitude Observatory, Learmonth Solar Observatory, Udaipur Solar Observatory, Instituto de Astrofísica de Canarias, and Cerro Tololo Interamerican Observatory. We thank Rafael García, Bernhard Fleck, Antonio Jiménez, Ton von Overbeek, Catherine Renaud and William Thompson for fruitful discussions.

References

Appendix A: Timing of helioseismic instruments

The time of the mid point of the first sample for each helioseismic instruments is given in Table A.1. For GONG, because of the use of the BDF, the first sample x is computed by taking x_GONG(t₀) = v(t₀) – v(t₀ – 1) where t₀ is the time in minutes given in the file (for GONG the file directly provides x_GONG(t₀)) and v(t₀) is the velocity. Because of the use of the BDF for GONG, the time delay is computed by also taking the BDF of the other instruments, that is by computing x_Other(t₁) = v(t₁ + 1) – v(t₁), where t₁ is the time in minutes given in the file. For the other instruments, the file directly provides v(t₁); therefore when t₀ = t₁, the time difference between GONG and the other instruments is −60 s.

All Tables

All Figures

Fig. 1. Left: correlation of the power spectrum as obtained by Fossat et al. (2017) as a function of frequency lag for two different time series sampled at 80 s (top) and 60 s (bottom); this is comparable to Fig. 10 of Fossat et al. (2017). The green vertical lines correspond to frequencies at 210 nHz, 630 nHz, and 1260 nHz. Right: sum of the correlation for l = 1 and l = 2 modes as obtained by Fossat et al. (2017) as a function of rotation frequency for two different time series sampled at 80 s (top) and 60 s (bottom); this is comparable to Fig. 16 of Fossat et al. (2017).

In the text

	Fig. 2. X ratio as a function of the velocity for PM1: one point per day. The fit of the X ratio is in green. The fit is for the blue wing after February 2002.
In the text

	Fig. 3. Velocity residuals for PM1 as a function of time: one point per day. The two vertical lines indicate when the modes of operation changed from blue wing to red wing and vice versa.
In the text

	Fig. 4. Cross correlation between the first difference of the GOLF radial velocity and the first difference of the l = 0 GONG radial velocity as a function of lag in minutes (black line), with the envelope shown (orange line) and the fitted cross correlation (green line).
In the text

	Fig. 5. Time difference between the SPM Blue/VIRGO data and the l = 0 GONG data as a function of time computed using a 10-day window (black dots) or a 45-day window (green dots).
In the text

Fig. 6. Time difference between the newly calibrated GOLF data and the l = 0 GONG data as a function of time; (top) unshifted GOLF data; (bottom) shifted GOLF data according to Table 2 computed using a 10-day integration window (black dots) or a 45-day integration window (green dots). The median value of the time shift for the corrected time series is 47.3 s. About 99.8% of the corrected values are within 5 s of this median value.

In the text

	Fig. 7. Top: time difference between the GOLF time series sampled at 60 s and the GONG l = 0 data as a function of time. Bottom: time difference between the GOLF time series used by Fossat et al. (2017) and the newly calibrated GOLF data as a function of time, both sampled at 80 s.
In the text

	Fig. 8. Top: rms mode amplitude as a function of time for the newly calibrated GOLF data (black), for the GOLF data as used in Fossat et al. (2017; orange), for the BiSON data (red), and for the GONG data (green). Bottom: rms noise as a function of time. The sub-series are 30 days in length.
In the text