Skip to main content

Estimation of temporal and spatial variation of sound speed in ocean from GNSS-A measurements for observation using moored buoy

Abstract

Using Global Navigation Satellite System–Acoustic (GNSS-A) technique, we have been developing observation system on a moored buoy for continuous monitoring of seafloor crustal deformation. The sound speed structure near a warm current has heterogeneity, which is the main cause of a seafloor positioning error. Assuming a sloping structure, previous studies proposed sound speed model to reduce positioning error. We examined the validity of the model by comparing the estimated structure with the actual structure measured at multiple points around our observation site. The result shows that the gradient parameter estimated from GNSS-A data acquired by vessel is appropriate. The numerical examination indicates that modeling error caused by the misinterpretation of the depth of gradient layer occurs, and it can be suppressed by performing acoustic ranging at the point near the centroid of units. From the calculation of estimation error of sound speed variation, the predicted acoustic ranging error observed using the moored buoy staying near the centroid is 9.0 cm or below. Therefore, seafloor displacement can be detected with centimeter class via moored buoy in the basin of a warm current.

Introduction

Spiess et al. (1998) originally presented observation method for seafloor crustal deformation by combining Global Navigation Satellite System (GNSS) positioning with acoustic ranging. It is difficult to estimate the absolute position on the seafloor using GNSS only because the electromagnetic wave propagates poorly in seawater. GNSS–Acoustic (GNSS-A) technique uses a vessel or buoy as a relay point and estimates the relative position on the seafloor by acoustic ranging. Although the estimation accuracy for stationary displacement has been enhanced by the previous studies (Fujita et al. 2006; Honsho et al. 2019; Ikuta et al. 2008; Kido et al. 2008; Yasuda et al. 2017; Yokota et al. 2019), it is still difficult to detect that for non-stationary displacement or separate with high-accuracy coseismic from postseismic one.

Yokota and Ishikawa (2019b) conducted an observation at frequent intervals and were successful to detect the displacement due to the slow slip event. From the results of this study, continuous monitoring is expected to increase the estimation accuracy. Recently, a buoy has been used as the continuous observation platform. Kido et al. (2018) developed the GNSS-A observation system on the moored buoy and had been operating it for almost 1 year. Kato et al. (2005) has utilized a buoy, named GNSS buoy, to observe tsunami. GNSS buoy was successful to detect tsunamis generated by the 2001 Peru earthquake, the 2004 off the Kii peninsula earthquake in Japan, and so on. Our research group then conducted experiments farther from the coast and extended the function of the GNSS buoy by introducing GNSS-A observation equipment. We previously reported the development of a comprehensive observation system for tsunami, ionospheric total electron content, precipitable water, and seafloor crustal deformation (Kato et al. 2018).

The moored buoy we used in this research is owned by Kochi Prefecture, Japan, which is located approximately 30 km off the coast of cape Ashizuri (Fig. 1). We installed three acoustic units on the seafloor around the buoy in 2017 and finished installing observation equipment, such as GNSS antenna and receiver, gyro sensor for attitude measuring, acoustic ranging equipment, computers for data processing, and satellite communication unit for data transfer, on the deck in 2018. The measurements of GNSS and gyro are acquired with 1 Hz sampling. The acoustic ranging at 3 min intervals for units one by one in sequence is continuously operated.

Fig. 1
figure 1

Location of our observation site (pink square) with schematic diagram of the Kuroshio basin

As shown by Yokota et al. (2019), from 2014 to 2016, water temperature at the surface layer around this region had gradient toward the southeast. The temperature gradient was generated by warm current called Kuroshio that flows northward in the Pacific Ocean along the eastern coast of Japan. Yokota end Ishikawa (2019a) also pointed out that sound speed structure has strong gradient when the observation point is located at the edge of the Kuroshio. Kido (2007) presented sound speed model considering horizontal gradient because the heterogeneity in sound speed structure is the main cause of positioning error for GNSS-A technique. Assuming a sloping structure for sound speed to estimate the variation from reference sound speed profile, Yasuda et al. (2017), Honsho et al. (2019), and Yokota et al. (2019) developed an analytical method. They eventually demonstrated that the estimation for the time series of unit position was improved.

In our case, we needed to use the measured in the past for analysis because our buoy was not equipped with the function of measuring sound speed. Therefore, the estimation accuracy for sound speed variation has more direct effects on the positioning accuracy. Accordingly, this study aimed to evaluate estimation accuracy for its temporal and spatial variations. We then measured sound speed at multiple points around seafloor units in parallel with acoustic ranging and obtained its three-dimensional structure. We applied the sloping structure model to GNSS-A data acquired using vessel and obtained seafloor unit position and sound speed variations. Then, we compared the estimated horizontal gradient with the measured structure. We also calculated the anticipated ranging error observed using buoy to examine the capability of the buoy to detect the seafloor displacement.

First, we presented sound speed model and analytical method considering its sloping structure. Then, we showed the analytical results of array shape and horizontal gradient of sound speed by using GNSS-A data acquired using vessel. We discussed the estimated acoustic ranging error observed using buoy. The modeling error caused by misinterpretation of sound speed structure is also discussed.

Observation and data

GNSS-A observation using vessel

We conducted observation using research vessel Yugemaru from 2017 to 2019 to determine the initial position of seafloor units and survey the sound speed structure at the site. Figure 2a shows the track maps for six epochs. First, the vessel ran clockwise and counterclockwise around the unit array on the circular route in either case. The vessel also navigated along straight orbits passing both the vicinity of the center of the circular route and the point right above each seafloor unit. The acoustic ranging was performed at intervals of 8 s for units one by one in sequence for about 5–6 h for each epoch. During acoustic ranging, the onboard GNSS receiver recorded satellite signal and the gyrocompass measured attitude of the vessel. The position of onboard transducer can be calculated from the measured attitude and the position of GNSS antenna placed on the vessel estimated using kinematic positioning. We also measured conductivity, temperature, and depth (CTD) of sea water in parallel with acoustic ranging.

Fig. 2
figure 2

Outline of conducted observation using vessel for six observational epochs. a Two-dimensional ship track maps. Pink squares and star symbols represent seafloor units and the point where observations for conductivity, temperature, and depth (CTD) were conducted. b The variation of sound speed in ocean along to the depth. c The sound speed perturbation calculated by Eq. (2)

Sound speed structure obtained from CTD measurements

We conducted underway CTD observations at the points as indicated by the star symbols on the circular route in Fig. 2a. We performed at multiple points to understand sound speed structure at our observation site. Figure 2b shows sound speed profiles at 1 m intervals of water depth for each epoch, which is calculated using the equation proposed by Del Grosso (1974). The similar variation patterns according to the depth were measured in 2017 and 2018, whereas the vertical gradient at the depth 200–300 m measured in 2019 were extremely large.

The horizontal structure of sound speed has been generally considered to be sloping in the region near the Kuroshio basin (Kido 2007; Yasuda et al. 2017; Yokota et al. 2019; Honsho et al. 2019). We applied this model to detect horizontal gradient of sound speed from CTD measurements. Then, sound speed at the point (x, y, z) is modeled by

$$ V\left(x,y,z\right)={V}_0(z)+\delta {V}_x(z)\cdot x+\delta {V}_y(z)\cdot y, $$
(1)

where V0(z) (m/s) is sound speed at the origin and δVx and δVy (m/s/m) are the constant gradients of east and south components, respectively. Using the least squares method, the parameters V0(z), δVx(z), and δVy(z) at 1 m intervals of the depth from multiple profiles (Fig. 2b) can be obtained. Figure 2c shows the perturbation from V0(z) calculated by

$$ \Delta V\left(z,\phi \right)=\delta {V}_x(z)\cdot R\;\sin \phi +\delta {V}_y(z)\cdot R\;\cos \phi, $$
(2)

where R (m) is radius of circular route and ϕ (radian) is the direction of gradient axis from the true north. We considered that if the estimation errors of δVx(z) and δVy(z) are extremely large, the horizontal structure at a particular depth is not sloping. Then, the perturbation was not painted at the depth where the radius of error ellipse for ΔV(z, ϕ) with 95% confidence probability was higher than its estimation value.

The results in 2017 and 2019 showed that a strong gradient mainly from southeast to south is generated at the depth 200–400 m. Contrarily, the structures in 2018 have weak gradient or complicated gradient field. Yokota and Ishikawa (2019a) categorized a structure in accordance with the relative position between observation site and the Kuroshio basin. The structure in shallow depth has weak gradient in the area inside the Kuroshio. On the edge of the Kuroshio, there is strong thick gradient layer. The unstable disturbance is likely to be generated because there is no strong current flow in the area outside the Kuroshio. According to this categorization, our observation site was considered to be on the edge for the epochs in 2017 and 2019 and outside the Kuroshio for the epochs in 2018.

Methods/experimental

Model definition for travel time of acoustic signal

Travel time between onboard transducer xs = (xs, ys, zs) and seafloor unit xu= (xu, yu, zu) is defined as the integration of slowness S, which is an inverse of sound speed, along the propagation path:

$$ T={\int}_0^LS\left(x,y,z,t\right)\mathrm{dl}, $$
(3)

where L = |xs− xu| is path length (m). To express the equation of travel time simply, we used S instead of sound speed. As shown in Eq. (1), several previous works considered that horizontal sound speed structure is sloping. We also assumed that horizontal structure consists of N layers, the thickness of each layer is H (m), and the gradient is uniform within each layer. Then, the slowness at the depth in nth layer at the time t is modeled by

$$ S\left(x,y,z,t\right)=a(t){S}_0(z)+\delta {S}_{x,n}\cdot x+\delta {S}_{y,n}\cdot y, $$
(4)

where S0(z) is slowness (s/m) at the origin andSn = (δSx,n, δSy, n) are the east and north components of horizontal gradient (s/m/m) in nth layer, respectively. According to Ikuta et al. (2008), the temporal variation of slowness a(t) is represented by a cubic B-spline curve:

$$ a(t)=\sum \limits_{i=0}^3{\alpha}_i{B}_{i,4}(t), $$
(5)

where αi is control points and Bi,4(t) is basis function. We set 3 h for the knot interval in this study to detect long-term variation caused by the ocean tide. Assuming that acoustic signal propagates linearly as described in Fig. 3a, we obtained dl = L/D dz (D = zs − zu). By substituting Eq. (4) into (3), the theoretical travel time can be written as

Fig. 3
figure 3

Schematic diagram drawing a propagation path and b horizontal slowness model having constant gradient from seafloor to sea surface. It is assumed that acoustic signal transmitted from vessel or buoy xs (green circle) propagates linearly to seafloor unit xu (pink square)

$$ {T}_{\mathrm{cal}}=L\cdot a(t){\overline{S}}_0\left({z}_{\mathrm{s}},{z}_{\mathrm{u}}\right)+\frac{L}{D}\Delta {T}_{\operatorname{grad},} $$
(6)

where the average sound speed at the origin is expressed as

$$ {\overline{S}}_0\left({z}_{\mathrm{s}},{z}_{\mathrm{u}}\right)=\frac{1}{D}{\int}_{z_{\mathrm{u}}}^{z_{\mathrm{s}}}{S}_0(z)\mathrm{dz}, $$
(7)

and the contribution of horizontal gradient to travel time is

$$ \Delta {T}_{\mathrm{grad}}=\sum \limits_{n=1}^N{\int}_{z_{\mathrm{u}}+\left(n-1\right)H}^{z_{\mathrm{u}}+ nH}\left[\delta {S}_{x,n}\cdot x+\delta {S}_{y,n}\cdot y\right]\kern0.5em \mathrm{dz}. $$
(8)

By substituting the coordinate of the point on the path \( x={x}_{\mathrm{u}}+\frac{z-{z}_{\mathrm{u}}}{z_{\mathrm{s}}-{z}_{\mathrm{u}}}\left({x}_{\mathrm{s}}-{x}_{\mathrm{u}}\right) \) into Eq. (8), the east component of ΔTgrad is given by

$$ {Hx}_{\mathrm{u}}\sum \limits_{n=1}^N\delta {S}_{x,n}+\frac{H^2}{D}\left({x}_{\mathrm{s}}-{x}_{\mathrm{u}}\right)\sum \limits_{n=1}^N\left(n-\frac{1}{2}\right)\cdot \delta {S}_{x,n} $$
(9)

The contributions of gradient to travel time are determined by its strength and thickness and depth of gradient layer. As suggested in Eq. (9), the parameters of strength δSx,n and thickness H are inseparable. Therefore, the previous works generally set H to 500–1000 m (Yasuda et al. (2017), Honsho et al. (2019)). Since the depth of our observation site is below 800 m (Fig. 1), this study assumed that strength and direction of gradient are constant from sea bottom to sea surface (Fig. 3b). Then, we obtained

$$ \Delta {T}_{\mathrm{grad}}=\frac{1}{2}D\left[\delta {S}_x\left({x}_{\mathrm{s}}+{x}_{\mathrm{u}}\right)+\delta {S}_y\left({y}_{\mathrm{s}}+{y}_{\mathrm{u}}\right)\right], $$
(10)

by substituting N = 1, H = D, and Sn =  S = (δSx, δSy) into Eqs. (8) and (9). ΔTgrad indicates a sum of product of gradient value and horizontal distance from the origin to the middle of path. Accordingly, theoretical travel time model applied in this study is

$$ {T}_{\mathrm{cal}}\left({\boldsymbol{x}}_{\mathrm{u}},{\boldsymbol{x}}_{\mathrm{s}},{S}_0,a,\nabla \boldsymbol{S},t\right)=L\left\{a(t){\overline{S}}_0\left({z}_{\mathrm{s}},{z}_{\mathrm{u}}\right)+\frac{1}{2}\left[\delta {S}_x\left({x}_{\mathrm{s}}+{x}_{\mathrm{u}}\right)+\delta {S}_y\left({y}_{\mathrm{s}}+{y}_{\mathrm{u}}\right)\right]\right\}, $$
(11)

obtained from Eqs. (6) and (10).

Analytical process for detecting seafloor displacement

We need to determine the shape of seafloor unit array in advance for detecting precisely seafloor displacement observed using buoy. Then, first, we determine the array shape using data acquired using vessel, and second, we estimate the displacement of array from the initial position with assumptions that every single unit moves in the same direction and also at the same speed due to the crustal deformation.

The theoretical travel time from the transducer of the vessel to the unit i for kth shot at time tk is given by Tcal(xu,i+ δxm, xs,k, S0,m, a, Sm, tk) for mth epoch. The position of unit i for mth epoch is expressed as xu,i+ δxm, where xu,i is its initial position and δxm is the vector of array displacement. To determine the array shape, the parameters to be estimated are xu,i, δxm, and temporal variation a(tk) and horizontal gradient Sm from the slowness profile S0,m. The average of their profiles was applied for S0,m, since this analysis is performed without the information of horizontal sound speed structure obtained from multiple profiles (Fig. 2c).

As shown in Fig. 2c, the horizontal structure changes day by day. Although the time variation of S should be considered, we treat S as constant value during a certain time in accordance with Honsho et al. (2019). We assumed that S is constant during observation for single epoch since the observation using vessel takes approximately 6 h for this site. For the observation using buoy, we were planning to divide observation period by a prescribed time interval and to consider S constant within interval.

The theoretical travel time from the buoy is also given by Tcal(xu,i+ δxp, xs,k, S0, a, Sq, tk). δxp and Sq are those of the pth and qth divisions when observation period is divided by a prescribed time interval for each. For the determination of array displacement, the unknown parameters are δxp,a(tk), and Sq. We need to use the profile measured by vessel in the past because our observation system on the buoy has no function for CTD measuring. Both analyses use the least squares method and determine the unknown parameters that minimize the sum of residuals between observed and theoretical travel time.

Results

Array shape of seafloor units

First, we examined the effect of considering horizontal gradient of sound speed for array shape determination. For this purpose, we compared array shape determined using sound speed model considering and ignoring gradient. Here, the model considering gradient (Eq. (4)) is called constant gradient model, and model ignoring gradient is called stratification model. By conducting analysis for each single epoch, a set of six shapes was obtained. Figure 4 compares the obtained array shapes for each sound speed model. The centroid position of each shape is arranged to the origin. For both horizontal and vertical shapes, the dispersion was suppressed by considering gradient. Table 1 also compares averaged distance between units with the standard deviation. The maximum dispersion of horizontal one between unit 1 and unit 3 was decreased from approximately 10 to 7 cm by considering gradient. The dispersions of vertical one were also decreased.

Fig. 4
figure 4

Determined horizontal and vertical array shape by using sound speed models: a constant gradient model and b stratification model (ignoring the gradient). The shape for epoch (1) is displayed with the scale of axes. The differences in horizontal and vertical shapes from which of epoch (1) are magnified by × 2000 and × 15, respectively, to make the dispersion easier to be seen. The epoch number corresponds the number indicated in Fig. 2a

Table 1 Averaged distances between seafloor units determined by using different sound speed models

We also showed the positional shift of unit array between two models for each epoch in Fig. 5. Violet and blue vectors in Fig. 5a denote the estimated horizontal gradient and the centroid shift, respectively. Two arrows point in almost the same direction. Since violet arrow points in the direction of increasing sound speed, the estimated unit position is shifted toward faster sound speed by ignoring gradient. Moreover, the shift amount seems to be proportional to the quantity of gradient parameter. Vertical array shift in Fig. 5b suggests that the depth of unit is estimated shallower than the true position for the unit where the sound speed is faster.

Fig. 5
figure 5

The positional shift of determined array shape for each epoch using sound speed models: constant gradient (blue triangle) and stratification model (yellow triangle). The shape for constant gradient model is displayed with the scale of axes. a Horizontal shift amount of the shape for stratification model is magnified by × 2000. Violet arrow denotes the estimated horizontal gradient vector and points in the direction of increasing sound speed. Blue arrow indicates the centroid shift. b Vertical shift amount is magnified by × 15. Horizontal axis indicates gradient axis for each epoch

Second, we conducted analysis to determine common array shape for six epochs by using constant gradient model. Table 2 displays the determined position of seafloor units expressed as the east–north–up coordinates centering the centroid of the array (latitude 32.48662 north, longitude 133.20850 east, and at ellipsoidal height − 756.15 m). From the calculated major axial radius of error ellipse with 95% confidence probability, the unit position was determined with 1.5 cm accuracy.

Table 2 Position of seafloor units estimated by using multiple campaign data

Horizontal sound speed structure

We also obtained horizontal gradient from the analysis to determine common array shape in the above section. The estimation results are summarized in the left column in Table 3. The magnitude of horizontal gradient ΔV (m/s/km) was converted from S by \( \Delta V\approx \left|\mathbf{\nabla}\boldsymbol{S}\right|/{\overline{S}}_0^2\times {10}^3 \). The estimated direction of gradient axis ϕ is also displayed. The estimated ϕs of epochs (1) and (2) in 2017 are southeast, which corresponds with the measured structure (Fig. 2c). The estimated ΔV value of epoch (3) is about half those of epoch (1). The effect of gradient on the observed travel time was canceled because the measured structure of epoch (3) is complicated. The estimated ΔV is smallest for epoch (4) because the structure has weak gradient. Although the structure of epoch (6) has large gradient at the depth 200–300 m, ΔV is smaller than those of epochs (1) and (2). It suggests that the effect of gradient is suppressed because the thickness of gradient layer is thin.

Table 3 Horizontal gradient of sound speed estimated from different profiles

Estimation error of sound speed variations

A comparison with actual structure demonstrated in the previous section that horizontal gradient estimated from GNSS-A data with sound speed profile obtained in parallel with acoustic ranging was appropriate. Although our buoy has no function for measuring sound speed, the profile measured in the past needs to be used. Accordingly, we examined the effect of substituting the profile obtained on other days on the estimation of sound speed variations. We applied the profiles measured on 6 June 2017, and 5 June 2019, to the calculation for other epochs and obtained temporal variation and horizontal gradient from those reference profiles.

Table 3 compares the horizontal gradient estimated using different profiles. The root-mean-square (rms) errors of residual between observed and calculated travel time are also shown. The rms error was 0.04–0.06 ms in the case of using profile obtained in each epoch. We then considered estimation using other profile to be validated when the rms error is less than 0.06 ms. Based on this index, the profile of epoch (1) should be selected for epochs (2)–(4) and that of epoch (6) should be selected for epoch (5). As shown in Fig. 2b, observed variation patterns along the depth for two epochs in 2019 were quite different from others. Thus, applying the profile markedly different from the actual pattern has adverse effect on the estimation. Even if the suitable one was selected as the reference, the estimation errors for ΔV and ϕ are maximum of 0.03 m/s/km and 3°, respectively.

Figure 6 also compares the estimated temporal variation of averaged sound speed (\( 1/a(t){\overline{S}}_0 \)). As mentioned above, the profile of epoch (6) should not be selected as the reference for epoch (1)–(4). However, regardless of the used profile, the temporal variation patterns were successfully identified. The offset in the figures represents the average of the residual of estimation values between used profiles. The results of all epochs indicate that using suitable profile can suppress the offset and the maximum value is 0.11 m/s.

Fig. 6
figure 6

Estimated temporal variation of averaged sound speed. Black lines indicate estimated variation from the profile observed on each day. Pink broken lines indicate estimated variation from the profile observed on a 6 June 2017, and b 5 June 2019, respectively. The offset was calculated by averaging the residual of value estimated from different profiles

Discussion

Separability between array displacement and sound speed gradient

Honsho et al. (2019) discussed the requirements to detect array displacement separately from horizontal gradient. In two-dimensional example, they theoretically demonstrated that the measurements obtained from the point survey at the center of an equilateral triangle array have little information to separate them. The buoy we applied is in a similar situation as shown in Fig. 7 showing its movement for a month observed through GNSS. The position of buoy had been stable within a square of approximately 400 m × 400 m. Therefore, we also considered the separability with three-dimensional example for observation using buoy.

Fig. 7
figure 7

Location of seafloor units and buoy observed in March 2019

According to Eq. (11), the observed travel time for unit i for kth shot is given by

$$ {T}_{i,k}={L}_{i,k}\left\{{\overline{S}}_0+\frac{1}{2}\left[\delta {S}_x\left({x}_{\mathrm{u},i}+\delta x+{x}_{\mathrm{s},k}\right)+\delta {S}_y\left({y}_{\mathrm{u},i}+\delta y+{y}_{\mathrm{s},k}\right)\right]\right\}, $$
(12)

where Li,k = |xs,k − (xu,i + δx)|. The quantity of travel time changes caused by the array displacement δx becomes

$$ \delta {T}_{i,k}^{\mathrm{disp}}=\left({L}_{i,k}-{L}_{i,k}^{\prime}\right){\overline{S}}_0\approx \frac{{\overline{S}}_0}{L_{i,k}^{\prime }}\left[\delta x\left({x}_{\mathrm{u},i}-{x}_{\mathrm{s},k}\right)+\delta y\left({y}_{\mathrm{u},i}-{y}_{\mathrm{s},k}\right)\right], $$
(13)

where \( {L}_{i,k}^{\prime }=\left|{\boldsymbol{x}}_{\mathrm{s},k}-{\boldsymbol{x}}_{\mathrm{u},i}\right| \). The temporal variation of \( {\overline{S}}_0 \) is ignored here for brevity. The contribution of horizontal gradient to travel time is

$$ \delta {T}_{i,k}^{\mathrm{grad}}\approx \frac{1}{2}{L}_{i,k}^{\prime}\left[\delta {S}_x\left({x}_{\mathrm{u},i}+{x}_{\mathrm{s},k}\right)+\delta {S}_y\left({y}_{\mathrm{u},i}+{y}_{\mathrm{s},k}\right)\right]. $$
(14)

The position of onboard transducer xs,k is fixed at the origin if the buoy stops at the centroid of unit array. Assuming that the array shape is equilateral triangle and all units are at the same depth, the path length L and the incident angle of acoustic signal θ are independent of unit and shot (Fig. 8c). Then, the quantities of travel time changes in Eqs. (13) and (14) are rewritten as

$$ \delta {T}_i^{\mathrm{disp}}\approx {\overline{S}}_0\mid \boldsymbol{\delta} \boldsymbol{x}\mid \sin \theta \cos \left(\gamma -{\psi}_i\right), $$
(15)
Fig. 8
figure 8

Schematic diagrams drawing a azimuth angle ϕ of horizontal gradient of slowness; b azimuth angle ψi of unit position xu,i; c incident angle θi,k of acoustic signal from buoy xs,k to seafloor unit xu,i; d azimuth angle γ of displacement vector δx

and

$$ \delta {T}_i^{\mathrm{grad}}\approx \frac{1}{2}{L}^2\mid \nabla \boldsymbol{S}\mid \sin \theta \cos \left({\phi}_S-{\psi}_i\right), $$
(16)

where the directions of array displacement and slowness gradient are γ and ϕS and xu,i is expressed by azimuth ψi (Figs. 8a, b, d). From these equations, the separation of two quantities is theoretically impossible when the buoy stops at the centroid.

Our buoy stays stable but actually moves 50–100 m a day. To show the separability in the actual condition, we performed the numerical calculation assuming that the buoy moves around the centroid. Synthetic travel time for 1 year was generated from Eq. (12) by setting array displacement (δx, δy) = (−44.2, 35.7) mm/year and horizontal gradient (|S| = 0.2 ns/m/m, ϕS = γ = 309°). We examined the effect of the movement range of xs,k by giving the random value within squares of 10 m × 10 m and 50 m × 50 m according to the uniform distribution. Figure 9 compares the results of δx and S estimated once a week and every 12 h, respectively. A comparison between two ranges of movement shows that the estimation accuracy becomes higher with wider movement range. Therefore, the actual condition for our buoy meets the requirement for determining array displacement correctly.

Fig. 9
figure 9

Results of the numerical experiment with different movement ranges of buoy: estimated a array displacement and b horizontal gradient of slowness. Pink lines indicate the regression lines derived from the estimated δx and δy.

Acoustic ranging error caused by estimation error of sound speed

We examined the effect of substituting profile obtained on other days for the estimation of sound speed variation (Table 3 and Fig. 6). To examine the capability of buoy to detect seafloor displacement, we discussed the estimated acoustic ranging error for observation using buoy in this section. From the positional relationship between buoy and seafloor units (Fig. 7), the average of observed travel time can be calculated and summarized in Table 4.

Table 4 Observed travel time from buoy to each unit and horizontal distance

First, we consider ranging error due to the misinterpretation of sound speed average. As indicated in Fig. 6, the estimation error of averaged sound speed ϵV was maximum of 0.11 m/s. Since the maximum of Tobs is 0.82 s for unit 3 (Table 4), the maximum ranging error given by ϵVTobs is 9.0 cm. As indicated in Table 3, estimation error of horizontal gradient ϵΔV was less than 0.03 m/s/km. ϵV caused by ϵΔV is also given by ϵV = xgradϵΔV, where xgrad is horizontal distance from the origin to the middle of path along the gradient axis. Since the maximum of xgrad is approximately 521 m for unit 1 (Table 4), the ranging error due to the estimation error of horizontal gradient becomes 1.0 cm as a maximum.

We demonstrated that acoustic ranging from buoy without CTD measuring parallelly will probably be performed with under 10 cm accuracy. Consequently, proposed analytical method can achieve seafloor positioning with accuracy of several centimeters class, which is sufficient for detecting seafloor displacement.

Validity of proposed sound speed model

The quantities of contribution of horizontal gradient to travel time is determined by its strength, thickness, and depth of gradient layer, as well as signal route passing through the layer, as indicated in Eq. (8). Particularly, the depth of layer determines the ratio of the contributions between xu and xs into the change of travel time. Assuming that sound speed structure has constant gradient from seafloor to sea surface, the ratio is fixed to 1:1, as shown by the coefficients of xu and xs in Eq. (10). However, the actual structures are not so simple (Fig. 2c). Then, modeling error can occur, which can cause seafloor positioning error.

In this section, we demonstrated the contribution of horizontal gradient to observed travel time by numerical calculation and verify the validity of the assumption for sound speed model. We examined five types of structure shown in the left column in Fig. 10. The contribution of gradient ΔTgrad was calculated for each structure by using Eq. (8) assuming that acoustic ranging is performed from the point in the displayed area. The distribution of ΔTgrad for each unit is described by the color scale.

Fig. 10
figure 10

Numerical calculation for the contribution of horizontal gradient of sound speed for five types of structure to observed travel time, which is denoted by ΔTgrad in Eq. (10). The magnitude ΔV and azimuth ϕV of sound speed gradient were given by the displayed value. The color scale in three right figures shows the distributions of ΔTgrad calculated for each unit by setting ranging position (xs, ys) within a range of − 1000 to 1000 m

Type A is the structure having constant gradient from sea surface to seafloor, which is the same as the introduced model in this study. The distribution of ΔTgrad suggests that if the gradient is ignored, the estimated position for every single unit can be shifted into the direction of gradient axis. Additionally, the unit where sound speed is faster can be estimated to be shallower. This result completely corresponds to that in Fig. 5.

The structure of type B has gradient layer in the middle depth and the thickness of layer is set to 100 m. The distribution of ΔTgrad is nearly the same as type A. Accordingly, even if the thickness of actual gradient layer is thinner than the assumed, the effect of gradient can be removed when the middle depth of actual layer corresponds to that of the assumed.

Type C has gradient layer in the depth of 200–300 m. In this case, the signal propagates through the layer near from the ranging position xs. Then, the contribution ratio of xu becomes small. Contrarily, if the structure has gradient layer near the seafloor, ΔTgrad is almost independent of xs.

Type D simply reproduces the measured structure of epoch (3), which consists of two layers. The distribution of ΔTgrad shows that major effect of gradient is cancelled because the direction of gradient axis is opposite. For this reason, the estimation value of ΔV was small (Table 3). Lastly, type E reproducing the structure of epoch (6) also consists of two layers, and the gradient of upper layer is much stronger. The effect of strong layer is dominant, and the distribution of ΔTgrad is similar to that of type C because the depth of layer is shallower than middle.

The numerical calculation suggests that modeling error occurs if the true depth of gradient layer is shallower or deeper than that of the assumed, such as types C and E. It makes the error increase at the point far from the origin. If the buoy stays stable near the array centroid, the modeling error hardly occurs because the contribution of xs can be ignored. In this case, the estimated value as the gradient parameter includes the information of depth and thickness of layer. Therefore, the assumption for sound speed structure in this study is basically appropriate for observation near the array centroid. If the buoy stays stable at the point far from the centroid, the modeling error could be proven fatal. Yokota and Ishikawa (2019a) proposed estimating the contributions of xu and xs to ΔTgrad separately. This method can reduce the modeling error irrespective of the ranging position. If we intended to introduce their idea for a point survey, we would be faced with another problem: insufficient number of observation data because we have only three units. Accordingly, the depth of gradient layer is required to be determined by using the Monte Carlo method or some other calculation techniques in the case of point survey far from the array centroid.

Conclusions

We have been developing GNSS-A observation system on the moored buoy to monitor seafloor crustal deformation continuously. The sound speed structure in the near a warm current has heterogeneity that is main cause of seafloor positioning error. Some previous works proposed sound speed model assuming a sloping structure. We applied this model to our analysis using data acquired through vessel and obtained horizontal gradient of sound speed with seafloor unit position. To validate the estimated gradient, we measured sound speed at multiple points around the unit array and obtained the three-dimensional structure. A comparison with measured structure showed that the estimated gradient was appropriate.

Since the buoy we applied in this research is not equipped with the function of measuring sound speed, we needed to use the profile measured in the past. Then, we also compared the estimation result for sound speed variations using profile observed on other days. The maximum estimation error for averaged sound speed and magnitude of horizontal gradient were 0.11 m/s and 0.03 m/s/km, respectively. Therefore, the estimated ranging error for observation using buoy will be less than 9.0 cm. This suggests that the presented analytical method can detect seafloor displacement within a required accuracy.

We indicated that it is difficult to separate array displacement from horizontal gradient of sound speed when the buoy completely stops at the unit array centroid. We also suggested that the misinterpretation of the gradient layer depth causes modeling error, and it possibly affects array displacement determination when the moored buoy stays stable at the point far from the array centroid.

Availability of data and materials

The datasets used and analyzed during the current study are available from the corresponding author on reasonable request.

Abbreviations

CTD: conductivity, temperature, depth:

ENU: east–north–up

GNSS: Global navigation satellite system

GNSS-A: GNSS-acoustic

rms: root mean square

STD: standard deviation

References

Download references

Acknowledgements

The authors are grateful to Kochi Prefecture giving permission to use the moored buoy for this survey. We also thank Dr. Akira Futamura and the crew of the research vessel “Yugemaru” owned by the National Institute of technology, Yuge College.

Funding

This study was supported by JSPS KAKENHI Grant Number JP16H06310: “A challenge to develop GNSS buoy system for high-functional tsunami monitoring and continuous observation of ocean-bottom crustal movements.”

Author information

Authors and Affiliations

Authors

Contributions

All authors participated in observation to obtain data used for analysis. NK developed the analytical software and performed analysis. The sound speed structure model used for analysis was designed by KT and NK. TK and YT contributed to the development and improvement of observation system. All authors read and approved the final manuscript.

Corresponding author

Correspondence to Natsuki Kinugasa.

Ethics declarations

Competing interests

The authors declare that they have no competing interests.

Additional information

Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Kinugasa, N., Tadokoro, K., Kato, T. et al. Estimation of temporal and spatial variation of sound speed in ocean from GNSS-A measurements for observation using moored buoy. Prog Earth Planet Sci 7, 21 (2020). https://doi.org/10.1186/s40645-020-00331-5

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://doi.org/10.1186/s40645-020-00331-5

Keywords