Sonic Log Measurements Research Articles

Experimental study of seismic dispersion: influence of clay mineral content

SUMMARY We conducted an extensive study on the elastic properties of Opalinus Clay and the overlying and underlying rock formations, which range in the overall clay mineral content from nearly 0 to 60 wt.%. Our laboratory experiments focused on seismic and ultrasonic frequencies to determine the extent to which seismic dispersion affects elastic parameters and seismic wave velocities. The results comprise the static stiffness from undrained triaxial cycles (axial–confining stress: 8–10 MPa), the dynamic stiffness at seismic frequencies (0.5–143 Hz), intrinsic attenuation (0.5–20 Hz), compressional velocity measurements (0.5–2 Hz) and ultrasonic velocity measurements (250 and 500 kHz). We compared these laboratory results to in situ sonic logging measurements to assess the role of frequency in measured elastic parameters. The results suggest a notable correlation between clay mineral content and dispersion. Specifically, high clay mineral content leads to increased dispersion, even within the seismic frequency band. The overall dispersion of P-wave velocity in the frequency range from 1 Hz to 500 kHz is up to 16%. This frequency dependency is crucial when establishing a relationship between sonic well log data and static stiffness for geomechanical modelling. The results are discussed with respect to possible dispersion mechanisms, including the role of bound water in clay.

Prediction of Compressional Slowness from Conventional Well Log Data using the Gradient Boosting Algorithm

Compressional slowness (DTCO) is the most basic parameter in geophysics, petrophysics, and geomechanics. These parameters can be obtained through the sonic log tool. However, equipment constraints, relatively new technology, and high cost of measurement make the parameters generated by sonic logs unavailable in old wells or wells being developed. Therefore, it is essential to predict sonic logs, especially in the case of compressional slowness prediction. Using machine learning, predictions can be generated by studying data on existing log wells. One of the algorithms that can produce predictions on continuous data, such as log values, is gradient boosting. MAPE and RMSE were used as evaluation metrics. The inputs used are gamma ray log data (GR), density (RHOB), porosity (NPHI), and shear slowness (DTSM). MAPE results show an error value of 12.28% with an RMSE of 10.74, indicating that the predictive model obtained has good results and performance. Using hyperparameter tuning in machine learning can reduce the error rate by 2.29% with faster processing times. In addition, it was found that the quantity of training wells can affect the resulting error value. The existence of this research can help a petrophysicist, geologist, and geophysicist characterize a reservoir with limited data. The use of this method also has the potential to be an alternative solution when sonic log measurements are expensive.

Acoustic Dispersion in Low Permeability Unconventional Reservoir Rocks and Shales at In Situ Stress Conditions

Utilizing laboratory measurements on dry and fully brine-saturated well-lithified shale reservoir rocks and comparing them to log data and theoretical modeling, we find no statistically significant intrinsic dispersion from seismic to sonic and laboratory measurement frequencies due to fluid effects. Under in situ stress conditions, the Gassmann zero-frequency P-wave velocity prediction for a Permian basin sample was within 0.3% of the measured velocity on the brine-saturated sample at an ultrasonic frequency. At deviatoric stresses ranging from 1.7 to 70 MPa on the same sample, the percent error in the Gassmann P-wave velocity prediction ranges from 0.3 to 2.2%. These results are consistent with the lack of significant dispersion in the reported velocities and theoretical predictions in the Cotton Valley shale. Based on the Biot–Gassmann model, the characteristic frequency in both formations occurs at ~1010 Hz. Applying a squirt flow model also predicts transition to the high-frequency regime occurring at ~109 Hz for both formations. The comparison of the sonic log measurements for four different shale/tight rock formations to ultrasonic measurements taken from core samples under in situ stress conditions confirms these findings. To our knowledge, this is the first extensive comparison of the sonic log to ultrasonic core data measurements. For these rocks, there is no clearly observable or predicted significant dispersion due to fluid effects at in situ stress conditions within the frequency range for which measurements are made.

Application of the differential effective medium theory to determine fluid saturation and crack density from measured P- and S-wave velocities

SUMMARY The differential effective medium (DEM) theory studied in this paper describes elastic moduli of a fractured medium with help of differential equations, where crack density is the independent variable and fluid saturation is a parameter. The effective medium is isotropic for randomly oriented flat ellipsoidal cracks and thus fully characterized by two elastic constants. In this article we derive an analytical solution of the equation for Poisson’s ratio and we transform the differential equation for Young’s modulus into a non-linear algebraic equation. Fluid saturation and crack density can then be determined from measured wave propagation velocities by a simple algorithm. We also derive approximate solutions for elastic moduli as a function of crack density and saturation, which allows to quantify the uncertainty of the result due to measurement errors. The DEM theory leads to higher crack densities than the self-consistent method and to lower crack densities than the non-interacting theory for measured elastic moduli, while all three methods give similar fluid saturation fractions. As an example of application of our theoretical results, we study weathered granite in the Strengbach water catchment in the Vosges mountains in France. We have performed full waveform sonic logging measurements in an 86-m-deep borehole located at an altitude of 1130 m above sea level, which is used for hydrogeophysical and geochemical studies of a granitic aquifer. The logging data allows us to investigate P and S waves in the depth range between 40 and 80 m. The P- and S-wave propagation velocities take average values of 5.0 and 2.7 km s–1, respectively, with the highest values of 5.8 and 3.2 km s–1 at 75–80 m depth. From these velocities we obtain a water saturation of 75 ± 25 per cent. The crack density describes the degree of weathering of the granite, which generally decreases with depth, but takes high values near layers of strongly weathered granite. Crack density is on average 0.5, with the highest value of 1.0 at 65 m and the lowest value of 0.2 at 75 m depth. The analysis of the full waveform logging data by the DEM method supports results from previous geochemical and hydrological studies in the Strengbach catchment which concluded that water is stored in deeper layers of the granitic aquifer.

Assessing the potential of passive seismic receiver functions for ore body exploration

ABSTRACTWe use a passive seismic high‐frequency receiver function method to image the shallow structure around a mine site. This is a relatively new application of a standard method for mapping major discontinuities in the crust and upper mantle at a scale small enough to be relevant in an exploration context. Data collected in a 21‐instrument array is inverted for isotropic velocity structure. The retrieved velocities in the south‐eastern part of the array match very well with available sonic log measurements in the top 800 m. Based on the differences in receiver function behaviour for stations across the array, recovered velocity profiles and their similarity to sonic log measurements, the target area is split into two separate regions along a north‐northeast trending line that correlates well to the strike of the dominant structure in the area – the Navan Fault. Our results demonstrate the ability of receiver functions to provide both qualitative and quantitative information in an exploration environment.

Data-driven control over short-period internal multiples in media with a horizontally layered overburden

SUMMARY Short-period internal multiples, resulting from closely spaced interfaces, may interfere with their generating (bandlimited) primaries, and hence they pose a long-standing challenge in their prediction and removal. A recently proposed method based on the Marchenko equation enables removal of the entire overburden-related scattering by means of calculating an inverse transmission response. However, the method relies on time windowing and can thus be inexact in the presence of short-period internal scattering. In this work, we present a detailed analysis of the impact of band-limitation on the Marchenko method. We show the influence of an incorrect first guess, and that adding multidimensional energy conservation and a minimum phase principle may be used to correctly account for both long- and short-period internal multiple scattering. The proposed method can currently only be solved for media with a laterally invariant overburden, since a multidimensional minimum phase condition is not well understood for truly 2-D and 3-D media. We demonstrate the virtue of the proposed scheme with a complex acoustic numerical model that is based on sonic log measurements in the Middle East. The results suggest not only that the conventional scheme can be robust in this setting, but that the ‘augmented’ Marchenko method is superior, as the latter produces a structural image identical to one where the finely layered overburden is missing. This is the first demonstration of a data-driven method to account for short-period internal multiples beyond 1-D.

A new crustal shear-velocity model in Southwest China from joint seismological inversion and its implications for regional crustal dynamics

SUMMARY Southwest (SW) China is located in a transition site from the active Tibetan Plateau to the stable Yangtze craton, which has complicated tectonic deformation and severe seismic hazards. We combine data from ambient noise, teleseismic body and surface waves, and petroleum wells to better constrain the crustal shear-velocity structure in SW China. We jointly invert the Rayleigh wave dispersion (5–40 s period), Rayleigh wave ZH ratio (20–60 s period), and P-wave receiver function for 114 permanent stations with a stepwise linearized joint inversion method. Compared to previous tomography results, we observe higher shear velocity in the sedimentary rocks within the Sichuan Basin, which is consistent with sonic logging measurements. Our model reveals widespread low-velocity zones in the mid-lower crust, and their boundaries correlate well with major fault systems. Between two main mid-crustal low-velocity channels, a prominent high-velocity region surrounded by earthquakes is observed in the inner zone of the Emeishan large igneous province (ELIP) and around the Anninghe-Zemuhe fault zone. These observations are comparable to regional tomography results using very dense arrays. Based on the results, we suggest that mid-lower crustal ductile flow and upper-crustal rigid fault movement play equally important roles in controlling the regional deformation styles and earthquake distribution in SW China. Our results also resolve thick crust-mantle transition zones beneath the eastern Tibetan Plateau and the inner zone of the ELIP due to ‘top-down’ and ‘bottom-up’ crust-mantle interactions, respectively. Our new model can serve as a reference crustal model of future high resolution model construction in SW China.

P-to-S-wave velocity ratio in organic shales

In situ P- and S-wave velocity measurements in a variety of organic-rich shales exhibit P-to-S-wave velocity ratios that are significantly lower than lithologically similar fully brine-saturated shales having low organic content. It has been hypothesized that this drop could be explained by the direct influence of kerogen on the rock frame and/or by the presence of free hydrocarbons in the pore space. The correlation of hydrocarbon saturation with total organic content in situ makes it difficult to separate these possible mechanisms using log data alone. Theoretical bounding equations, using pure kerogen as an end-member component without associated gas, indicate that kerogen reduces the P- and S-wave velocities but does not in general reduce their ratio enough to explain the observed low velocity ratio. The theoretical modeling is consistent with ultrasonic measurements on organic shale core samples that indicate no dependence of velocity ratios on the kerogen volume alone. Sonic log measurements of P- and S-wave velocities in seven organic-rich shale formations deviate significantly (typically more than 5%) from the Greenberg-Castagna empirical brine-saturated shale trend toward lower velocity ratios. In these formations, and on core measurements, Gassmann fluid substitution to 100% brine saturation yields velocity ratios consistent with the Greenberg-Castagna velocity trend for fully brine-saturated shales, despite the high organic content. These sonic and ultrasonic measurements, as well as theoretical modeling, suggest that the velocity ratio reduction in organic shales is best explained by the presence of free hydrocarbons.

A comparative analysis of time–depth relationships derived from scientific ocean drilling expeditions

Time–depth relationships (TDR) are required for correlating geological information from drill sites with seismic reflection profiles. Conventional time–depth domain conversion is implemented using P-wave velocity data, derived from downhole sonic logs, calibrated with vertical seismic check-shots. During scientific ocean drilling expeditions, immediate seismic correlation is carried out using laboratory velocities measured on recovered core material. As these three velocity measurements vary significantly in signal frequency, resolution and acoustic pathways, they carry potential for substantial TDR differences and consequent miscorrelation to seismic profiles. Our analytical work uses the comprehensive scientific ocean drilling dataset to quantify these differences in core-seismic integration. TDRs are calculated and compared at sites where check-shot, sonic log, and laboratory velocity measurements cover the same depth segments of the drill hole. We find that the maximum differences between the TDRs (TDR diffmax ) reach up to 55%, which can cause fundamental errors in the seismic correlation. No direct relationship to porosity and bulk density of the cored material is observed. Instead, higher TDR variability is found at sites with carbonate content > 70%, particularly with coarser grain texture. Sites containing primarily igneous and siliciclastic sequences show less than 10% TDR diffmax . This semi-quantitative criterion indicates that downhole logging should be conducted during drilling expeditions, especially at sites with carbonate sequences, or low core recovery, to ensure accurate core-seismic integrations.

Assessing uncertainties in high-resolution, multifrequency receiver-function inversion: A comparison with borehole data

We use teleseismic P-to-S converted waves from a permanent station to estimate the uncertainties in a 1D elastic model of the shallow crust (0–7 km depth) obtained from the inversion of receiver function (RF) data. Our earth model consists of layers with a constant S-wave velocity [Formula: see text] and P- to S-wave velocity ratio ([Formula: see text]). We apply a Bayesian formulation and transdimensional Monte Carlo sampling to compute the posterior uncertainties of the earth model. The model uncertainties rely on a realistic representation of the data uncertainties, and we estimate directly from the stacking of the teleseismic data, a full-error covariance matrix. To explore the effect of the number of teleseismic events and the RF frequency content, we compare the results of inverting a single RF computed for a cut-off filter frequency of 4 Hz with the joint inversion of four RFs computed from independent ensembles in a larger pool of events for cut-off frequencies of 0.5, 1, 2, and 4 Hz. The inversion results are compared with the lithostratigraphy and sonic-log measurements from a 7 km deep borehole drilled near the seismic station. The inversion of a single RF results in larger uncertainties in the recovered [Formula: see text] profile and in the depth to seismic discontinuities compared with the multifrequency inversion. Moreover, the multifrequency inversion predicts more accurately the depth to a velocity inversion at approximately 6 km below the surface and matches more closely the borehole sonic-log data. Our results indicate that RF data can be used to map shallow (3–5 km depth) crustal interfaces with uncertainties in the order of 300–500 m, whereas uncertainties are consistently smaller (<300 m) for interfaces in the top kilometer.

Integrated petrophysics and rock physics modeling for well log interpretation of elastic, electrical, and petrophysical properties

Rock and fluid volumetric properties, such as porosity, saturation, and mineral volumes, are generally estimated from petrophysical measurements such as density, resistivity, neutron porosity and gamma ray, through petrophysical equations. The computed petrophysical properties and sonic log measurements are generally used to estimate the petro-elastic relationship between elastic and rock and fluid volumetric properties used in reservoir characterization. In this paper, we present a unified workflow that includes petrophysical relations and rock physics models for the estimation of rock and fluid properties from elastic, electrical, and petrophysical (porosity, density, and lithology) measurements. The multi-physics model we propose has the advantage of accounting for the coupled effect of rock and fluid properties in the joint petro-elastic and electrical domains, and potentially reduce the uncertainty in the well log interpretation. Furthermore, the presented workflow can be eventually extended to three-dimensional reservoir characterization problems, where seismic and electromagnetic data are available. To demonstrate the validity of the methodology, we show the application of this multi-physics model to both laboratory measurements and well log data.

Distributed temperature sensing applied during diaphragm wall construction

Distributed temperature sensing (DTS) can be used to monitor the production process of diaphragm walls. DTS is able to differentiate between already present or fresh bentonite suspensions during refreshing of the bentonite slurry in the trench. During concrete casting, DTS is able to differentiate between the bentonite suspension and concrete. As a result, the continuity of the casting process and the arrival of good grade concrete at crucial locations in the trench can be monitored. Tests conducted on laboratory models provided reference information for interpretation of field data. Field experiences have shown the benefits of DTS tests and the predictive value of the reference measurements. Results are compared with crosshole sonic logging measurements at the same location.

Critical investigation of calculation methods for the elastic velocities in anisotropic ice polycrystals

Abstract. Crystallographic texture (or fabric) evolution with depth along ice cores can be evaluated using borehole sonic logging measurements. These measurements provide the velocities of elastic waves that depend on the ice polycrystal anisotropy, and they can further be related to the ice texture. To do so, elastic velocities need to be inverted from a modeling approach that relate elastic velocities to ice texture. So far, two different approaches can be found. A classical model is based on the effective medium theory; the velocities are derived from elastic wave propagation in a homogeneous medium characterized by an average elasticity tensor. Alternatively, a velocity averaging approach was used in the glaciology community that averages the velocities from a given population of single crystals with different orientations. In this paper, we show that the velocity averaging method is erroneous in the present context. This is demonstrated for the case of waves propagating along the clustering direction of a highly textured polycrystal, characterized by crystallographic c axes oriented along a single maximum (cluster). In this case, two different shear wave velocities are obtained while a unique velocity is theoretically expected. While making use of this velocity averaging method, reference work by Bennett (1968) does not end with such an unphysical result. We show that this is due to the use of erroneous expressions for the shear wave velocities in a single crystal, as the starting point of the averaging process. Because of the weak elastic anisotropy of ice single crystal, the inversion of the measured velocities requires accurate modeling approaches. We demonstrate here that the inversion method based on the effective medium theory provides physically based results and should therefore be favored.

Mechanical characteristics of laminated sand–shale sequences identified from sonic velocity and density correlations

In the past, numerous relations between sonic velocity (i.e. an inverse of sonic travel time DT) and rock density were developed that were suitable for certain fields or certain rocks only. In this paper, the correlations between sonic and density log measurements were investigated again using the data information from multiple fields of Gulf of Mexico (GOM) and North Sea (NS). For sandstone-shale sequences, the dominant linear relationships between sonic travel time (DT) and formation density (RHOB) were observed, and compared with the popular Gardner’s method (Gardner et al. in Geographics 39(6):770–780, 1974) along with other known methods. The implications of data clusters distinguished by formation lithologies and rock mechanical strengths were revealed from cross-plot analysis. The probabilities and uncertainties of the developed correlations were determined using the actual histograms from the collected data of GOM and NS fields coupled with Monte Carlo simulations.

Ultrasonic measurements of attenuation and velocity of compressional and shear waves in partially frozen unconsolidated sediment and synthetic porous rock

The presence of partially frozen liquid in the pore spaces of porous materials has significant effects on elastic wave propagation. Although the characterization of partially frozen systems using velocity information has been well developed, application of the attenuation information is limited because the attenuation mechanisms in partially frozen systems are poorly understood. We have conducted ultrasonic wave transmission measurements with changing temperatures from 0°C to [Formula: see text] to estimate the effect of partially frozen liquids grown in unconsolidated (unconsolidated sediment) and consolidated (synthetic porous rock) materials on the velocity and attenuation of P- and S-waves. Our experimental results determined that the existence of partially frozen liquid in the unconsolidated and consolidated materials increases the velocity and attenuation for temperatures of 0°C to around the freezing point (i.e., [Formula: see text]), thus experimentally validating the unintuitive observations of high velocity and high attenuation. We interpreted the differences in velocity-versus-temperature curves as both the difference in inherent stiffness between the matrix of the consolidated material and the ice frame of the partially frozen unconsolidated material and the microscale ice distribution in pore spaces. We have also attributed the difference in the attenuation-versus-temperature curves in the unconsolidated and consolidated materials between the P- and S-waves to the difference of attenuation mechanism between the P- and S-waves. Our findings can be used for interpreting the velocity and attenuation results from the sonic logging measurements.

Spatial sensitivity functions for rapid numerical simulation of borehole sonic measurements in vertical wells

Borehole sonic measurements are routinely used to measure dynamic elastic properties of rock formations surrounding a wellbore. However, shoulder-bed effects, mud-filtrate invasion, and near-wellbore damage often influence the measurements and bias the interpretations. Inversion-based methods can reduce the influence of complex geometric conditions in the estimation of formation properties from sonic logs, but they require fast-forward modeling algorithms. We have developed new spatial sensitivity functions for rapid modeling of borehole sonic measurements, which were equivalent to the Green’s functions of borehole modal wave measurements. An adaptive finite-element (FE) method was used to calculate 1D axial (vertical) and 2D (axial-radial) sensitivity functions, which combined the dependence of sonic logs on the wave mode, frequency, formation properties, and logging-instrument geometry. Semianalytical formulations were also used to calculate the sensitivity functions, which involved modeling the wave propagation through an effective layered medium. The 1D axial sensitivity functions were also extended to quantify propagation properties of nondispersive modes in a borehole. Calculated semianalytical sensitivity coefficients were efficient for rapid simulation of modal frequency dispersion in the presence of a logging tool. Simulated sonic logging measurements across synthetic thinly bedded and invaded formations were compared with numerical simulations obtained with the finite-element method. Our results confirmed the efficiency, reliability, and accuracy of the approximate sonic simulation method. The maximum relative error of the rapid simulation method was consistently less than 4%, with only 2% of the central processing unit time and 2% of the memory usage compared with FE simulations.

Propagation of elastic waves through textured polycrystals: application to ice.

The propagation of elastic waves in polycrystals is revisited, with an emphasis on configurations relevant to the study of ice. Randomly oriented hexagonal single crystals are considered with specific, non-uniform, probability distributions for their major axis. Three typical textures or fabrics (i.e. preferred grain orientations) are studied in detail: one cluster fabric and two girdle fabrics, as found in ice recovered from deep ice cores. After computing the averaged elasticity tensor for the considered textures, wave propagation is studied using a wave equation with elastic constants c=〈c〉+δc that are equal to an average plus deviations, presumed small, from that average. This allows for the use of the Voigt average in the wave equation, and velocities are obtained solving the appropriate Christoffel equation. The velocity for vertical propagation, as appropriate to interpret sonic logging measurements, is analysed in more details. Our formulae are shown to be accurate at the 0.5% level and they provide a rationale for previous empirical fits to wave propagation velocities with a quantitative agreement at the 0.07–0.7% level. We conclude that, within the formalism presented here, it is appropriate to use, with confidence, velocity measurements to characterize ice fabrics.

Analytic Simulation of the Elastic Waves Propagation in the Neighborhood of Fluid Filled Wells with Monopole Sources

An analytic formulation to understand the scattering, diffraction and attenuation of elastic waves at the neighborhood of fluid filled wells is presented. An important, and not widely exploited, technique to carefully investigate the wave propagation in exploration wells is the logging of sonic waveforms. Fundamental decisions and production planning in petroleum reservoirs are made by interpretation of such recordings. Nowadays, geophysicists and engineers face problems related to the acquisition and interpretation under complex conditions associated with conducting open-hole measurements. A crucial problem that directly affects the response of sonic logs is the eccentricity of the measuring tool with respect to the center of the borehole. Even with the employment of centralizers, this simple variation, dramatically changes the physical conditions on the wave propagation around the well. Recent works in the numerical field reported advanced studies in modeling and simulation of acoustic wave propagation around wells, including complex heterogeneities and anisotropy. However, no analytical efforts have been made to formally understand the wireline sonic logging measurements acquired with borehole-eccentered tools. In this paper, the Graf’s addition theorem was used to describe monopole sources in terms of solutions of the wave equation. The formulation was developed from the three-dimensional discrete wave-number method in the frequency domain. The cylindrical Bessel functions of the third kind and order zero were re-derived to obtain a simplified set of equations projected into a bi-dimensional plane-space for displacements and stresses. This new and condensed analytic formulation allows the straightforward calculation of all converted modes and their visualization in the time domain via Fourier synthesis. The main aim was to obtain spectral surfaces of transfer functions and synthetic seismograms that might be useful to understand the wave motion produced by the eccentricity of the source and explain in detail the new arising borehole propagation modes. Finally, time histories and amplitude spectra for relevant examples are presented and the validation of time traces using the spectral element method is reported.

Frequency-domain simulation of logging-while-drilling borehole sonic waveforms

Numerical simulation of sonic logging-while-drilling (LWD) borehole measurements is challenging because of significant wave propagation effects due to the massive drilling collar occupying a large portion of the borehole. In addition, the internal structure of the LWD tool can have a significant impact on the measured dispersions of Stoneley and quadrupole modes. The collar is typically constructed with a set of inner periodic grooves, which act as a mechanical filter to attenuate undesirable collar modes. Reliable numerical simulation and interpretation of LWD sonic waveforms requires that all features and dimensions of the drilling collar be included in the simulation model. Furthermore, the presence of the drilling collar can prompt numerical instabilities due to backward propagating modes in the perfectly matched layer (PML) commonly used to truncate the computational domain. This problem can be circumvented with the implementation of artificial viscoelastic attenuation in the collar whenever the simulations are intended to reproduce only wave propagation within the surrounding rock formations. In addition, reliable modeling of borehole wave propagation in the presence of high contrasts in material properties and the internal structure of the LWD collar requires a numerical method capable of accurately and stably resolving all spectral scales present in the model. We implemented an automatic [Formula: see text]-adaptive finite-element method in the frequency domain combined with a PML technique to simulate LWD sonic logging measurements. Examples of the application verified the accuracy and reliability of the simulated borehole and formation propagation modes in the presence of casing and internal structures in the LWD collar. The presence of steel casing and quality of casing/formation bond significantly influence the propagation modes excited in a borehole. However, it is still possible to estimate the formation shear slowness using monopole and quadrupole sources regardless of the quality of cement bond in fast formations. Assessment of the formation compressional slowness was significantly impeded by the strong pipe mode. Estimation of formation shear slowness in slow formations is practically impossible due to the presence of casing and a strong annulus mode when the quality of casing bond is poor.

Quantifying uncertainties in attenuation estimation at methane-hydrate-bearing zones using sonic waveform logs

Application of seismic attenuation estimation using sonic waveform data is limited because the estimation methods have not yet been fully developed. Although the median frequency shift method is considered to be effective and robust compared to the conventional spectral ratio method, we demonstrated that the median frequency shift methods strongly depend on reference data under lower signal-to-noise ratios. We modified an existing median frequency shift method not to depend on arbitrarily choosing a reference value and to quantify the uncertainties in attenuation estimation. Furthermore, we implied the optimum selection of receiver pairs used for more stable attenuation analysis. Our numerical experiments supported the advantage of the proposed method. Although our main findings by applying the proposed methods in methane hydrate-bearing sediments are almost consistent with past field sonic logging measurements, we find some differences in the magnitude of attenuation values compared to existing sonic attenuation measurements and discuss various possible factors. We believe that more stable and reliable attenuation results can lead to clarifying various factors affecting attenuation estimation, such as the effect of scattering, near-field effects, and source-coupling effects. Furthermore, we emphasized the importance of scattering effect caused by the heterogeneity of the formation and demonstrated the limitation of characterizing the 1D heterogeneity using the sonic logging data spatially sampled at 0.15 m to adequately estimate the effect of scattering attenuation.