WO2021222797A1 - System and method for determining high resolution subsurface velocity estimates - Google Patents

Abstract

It is known that the outer product covariance eigenvectors of a matrix of seismic data are needed to obtain a covariance solution equivalent to the conventional semblance solution. Shown herein, the outer product covariance solution can be obtained in a more computationally efficiently manner using the inner product covariance solution as an intermediate step. More particularly, the inner product eigenvalues and be used to obtain the outer product eigenvectors which, in turn, can be used to construct a high resolution semblance solution at decreased cost. Further, various embodiments provide an improved method of obtaining subsurface velocity estimates from seismic data than has been available previously.

Description

SYSTEM AND METHOD FOR DETERMINING HIGH RESOLUTION SUBSURFACE VELOCITY ESTIMATES CROSS REFERENCE TO RELATED APPLICATIONS [0001] This application claims the benefit of U.S. Provisional Patent Application serial number 63/018,466 filed on April 30, 2020 and incorporates said provisional application by reference into this document as if fully set out at this point. TECHNICAL FIELD [0002] The instant invention relates generally to systems and methods of seismic processing and, more particularly, to systems and methods of calculating an improved seismic velocity analysis. BACKGROUND [0003] It is hard to overestimate the importance of accurate velocity estimates to conventional seismic processing. Without accurate subsurface velocity estimates well focused seismic images of the subsurface cannot be obtained. Of course, absent well focused images the risk associated with drilling wells and other seismic exploration efforts will be greatly increased. In a typical processing sequence, velocity estimation is performed early on as a separate step. It should be noted that velocity estimation actually involves determining a time-varying velocity function at several different points along a seismic line or within a 3D volume. One thing that complicates this determination is that the seismic data itself must often be used to estimate the velocities that will produce the best image from that same data set. [0004] One popular approach to estimating subsurface velocities involves computation of a semblance analysis, also called a semblance panel or a semblance display. This approach has a number of well understood problems, but it still is used by many seismic interpreters to obtain a time versus velocity function which is then provided as input to subsequent seismic processing steps (e.g., stacking, time migration, etc.). [0005] The conventional semblance coefficient (Neidell and Taner, 1971) is widely used today to compute the seismic velocity spectra. Several algorithms have also been developed for computing the semblance coefficient in specialized cases, as well as suggested improvements for conventional semblance (Abbad and Ursin, 2012). Of the latter, Luo and Hale (2012) is especially noted as developers of a method which significantly improves velocity resolution, but at great computational cost. Over the years, alternatives to the conventional semblance coefficient have been proposed, but none have received general acceptance. One noted deficiency of the conventional semblance coefficient is the inability to properly separate noise and signal subspaces (Shan et al, 1985; Kirlin, 1992). [0006] One improvement to a conventional semblance approach came about through the recognition that semblance analysis is inherently noisy due to the inclusion of noise in the recorded data of the seismic data set. As might be suspected, if the noise and signal could be separated better velocity estimates would be obtained. One approach that has had some success involves calculating the eigenvalues of a matrix that is used in the semblance calculation and then using only the eigenvector(s) that correspond to the largest eigenvalues to calculate the semblance. This approach (i.e., a subspace or eigenstructure method) has been shown to reduce the contribution of noise which, in turn, results in a semblance display that has greater resolving power than was available with a conventional approach. [0007] For nearly three decades subspace or eigenstructure, techniques have been investigated as alternatives to conventional semblance for computing seismic velocity spectra, yet still receive little use in seismic applications today due to their high computational costs. Seismic velocity analysis using eigenstructure techniques has developed along two main lines. First, these approaches include the so-called hybrid approach in which the conventional semblance coefficient is reformulated using eigen analysis. In a second more direct approach, eigen analysis and MUSIC techniques are used to define a coherence measure. Those of ordinary skill in the art will recognize that MUSIC is an acronym for MUltiple SIgnal Classification which is a commonly used high-resolution method which is based on some properties of the eigen-decomposition of seismic data. [0008] The hybrid and direct approaches both use eigen decomposition based on computing the spatial data covariance matrix, or the so-called matrix outer product. In a typical application for seismic velocity analysis, in which the spatial dimension of the problem is significantly greater than the temporal domain, the solution for the seismic velocity spectra is computationally intensive. Although smoothing techniques and partial stacking have been used to lower the computational cost, the hybrid and direct eigenstructure approaches have still seen limited use. One feature of an eigen decomposition approach is that by using the largest eigenvalue(s) in the generation of velocity spectra is that a higher a higher resolution result can be obtained which is less bothered by the noise associated with the higher order modes. [0009] The hybrid and direct approaches both use eigen decomposition based on computing the spatial data covariance matrix, or so-called matrix outer product. In a typical application for seismic velocity analysis, in which the spatial dimension of the problem (e.g., the number of traces in a gather) is significantly greater than the temporal domain (e.g., the length of the analysis window), obtaining the seismic velocity spectra is computationally intensive. Though smoothing techniques, partial stacking, and specialized algorithms for finding the lowest order eigenvalues can be used to lower the computational cost, the hybrid and direct approaches have still seen limited use. Despite expectations of improved velocity resolution, time domain direct eigenstructure approaches have not been widely used, due largely to their increased computational times. [0010] Hence, what is needed is a more computationally efficient and effective method of determining seismic velocities using eigenstructure analysis than has been available heretofore. [0011] Before proceeding to a description of the present invention, however, it should be noted and remembered that the description of the invention which follows, together with the accompanying drawings, should not be construed as limiting the invention to the examples (or embodiments) shown and described. This is so because those skilled in the art to which the invention pertains will be able to devise other forms of this invention within the ambit of the appended claims. SUMMARY OF THE INVENTION [0012] Shown herein, the outer product covariance solution of the semblance calculation can be obtained in a more computationally efficiently manner using the inner product covariance solution as an intermediate step. More particularly, the inner product eigenvalues and be used to obtain the outer product eigenvectors which, in turn, can be used to construct a high-resolution semblance solution at decreased computational cost. Further, various embodiments provide an improved method of obtaining subsurface velocity estimates from seismic data as comparted with approaches that have been available previously. [0013] Contained herein is a disclosure of an embodiment which utilizes an eigen decomposition approach to calculating seismic velocity spectra that calculates the temporal covariance data matrix as an intermediate step and uses its eigenvalues to calculate the eigenvectors of the spatial covariance data matrix. This approach appropriately recognizes covariance analysis is a rank reduction process. Whether an outer or inner product is used, the full information content is preserved in identical sets of reduced-rank eigenvalues. Using non normalized definitions for the inner and outer products, a showing is made that the two sets of reduced eigenvalues are identical. A fast-spatial covariance data matrix version of the hybrid algorithm is applied and compared with the results using both conventional semblance and the original hybrid algorithm directly calculating the spatial covariance data matrix. A solution with the resolution of the spatial data covariance matrix is obtained with a much lower computational cost. This fast-spatial covariance solution is significantly more competitive to conventional semblance in computational cost [0014] The foregoing has outlined in broad terms some of the more important features of the invention disclosed herein so that the detailed description that follows may be more clearly understood, and so that the contribution of the instant inventors to the art may be better appreciated. The instant invention is not to be limited in its application to the details of the construction and to the arrangements of the components set forth in the following description or illustrated in the drawings. Rather, the invention is capable of other embodiments and of being practiced and carried out in various other ways not specifically enumerated herein. Finally, it should be understood that the phraseology and terminology employed herein are for the purpose of description and should not be regarded as limiting, unless the specification specifically so limits the invention. BRIEF DESCRIPTION OF THE DRAWINGS [0015] These and further aspects of the invention are described in detail in the following examples and accompanying drawings. [0016] Figs.1A-1D contain an example of multi-channel marine seismic data, tying events of the near offset and far offset through the shot gather: A) seismic section created from far offset channel; B) shot gather 1800; C) shot gather 1800 (same as B); D) seismic section created from near offset channel. [0017] Figs.2A-2C contain an example of conventional semblance analysis: A) conventional semblance velocity spectra; B) total energy; C) signal energy. [0018] Fig.3 contains an embodiment of a Spatial eigen covariance semblance velocity spectra analysis for sum of all eigenvalues. [0019] Fig.4 contains an embodiment of a temporal eigen covariance semblance velocity spectra analysis for sum of all eigenvalues. [0020] Figs.5A-5F contain an embodiment of an eigen analysis for the spatial data covariance matrix: A) scree plot for all 636 eigenvalues; B) scree plot for the 25 largest eigenvalues; C) scree plot for eigenvalues order 26 thru 636; D) semi-log (y) scree plot for the 25 largest eigenvalues; E) scree plot of absolute value for eigenvalues 26 thru 636 [same data as Figure 5 (c)]; F) semi-log (y) scree plot of absolute values for all 636 eigenvalues. In both plots E) and F) here, eigenvalues above 349 are negative as shown previously in Figure 5C. [0021] Figs.6A-6D contain an exemplary comparison of the spatial (A and B) with the temporal eigen analysis velocity spectra (C and D): a) scree plot for all 636 eigenvalues of spatial covariance data matrix; b) eigenvectors for eigenvalues 1 thru 3; c) scree plot for all 25 eigenvalues of temporal covariance data matrix; d) eigenvectors for eigenvalues 1 thru 3. [0022] Figs.7A-7D contain a comparison of the spatial (A and B) with the temporal eigen analysis velocity spectra (C and D): A) scree plot for first 25 eigenvalues of spatial covariance data matrix; B) eigenvectors for eigenvalues 1 thru 3; C) scree plot for all 25 eigenvalues of temporal covariance data matrix; D) translated eigenvectors for eigenvalues 1 thru 3. [0023] Figs.8A-8B contain an exemplary embodiment of a full spectrum semblance eigen analysis: A) spatial covariance semblance velocity spectra; B) temporal covariance semblance velocity spectra. [0024] Figs.9A-9D contain an exemplary embodiment of a Temporal covariance eigen-Semblance velocity spectra results associated with A) eigenvalue 1; B) eigenvalue 2; C) eigenvalue 3; and D) eigenvalue 4. [0025] Fig.10 contains an operating logic suitable for use with an embodiment of the instant invention. DETAILED DESCRIPTION [0026] While this invention is susceptible of embodiment in many different forms, there is shown in the drawings, and will herein be described hereinafter in detail, some specific embodiments of the instant invention. It should be understood, however, that the present disclosure is to be considered an exemplification of the principles of the invention and is not intended to limit the invention to the specific embodiments or algorithms so described. [0027] Hereinafter, a bold upper case X represents a two-dimensional seismic data matrix. Time is represented by the row dimension i, and columns j represent the spatial dimension (e.g., different channels or seismic traces that are spaced apart on the surface of the ground). Individual components of the two dimensional matrix are x_{i ,j} . Row and column vectors can then be represented respectively, as

) where M is preferably an odd number of samples in the time gate (i.e., the length of

the time window) and N is the number of data channels or seismic traces. Note that x_i is a row vector of data samples that are taken from each trace in the data matrix (i.e., it is a spatial vector), whereas x_j is column vector of time samples taken from the same channel or trace (i.e., it is a time vector). The vector xi has dimension 1 by N (i.e., one row by N columns), whereas the vector x_j has dimension M by 1. [0028] The semblance coefficient is defined to be the normalized output/input energy ratio. Under an assumption that the noise sum over all channels at any time is zero, the semblance coefficient is equal to the ratio of the signal energy over the total energy (Neidell and Taner, 1971). The equation for the conventional semblance coefficient from Neidell and Taner (1971, eqn.11) becomes:

The notation j(i) indicates that the sample on the j^th trace is a function of the several variables including the chosen analysis window location, the length of the data window (“M”), the data values on each trace in the gather, and an implicit velocity term. For example, in some embodiments the selected sample on the i^th trace will be determined via the well-known hyperbolic trajectory (or moveout) expression,

Where T_x is two-way arrival time, T₀ is the two-way normal incidence time, x is the offset distance, and

is dip-weighted rms velocity or stacking velocity. Thus, in equations 1a and 1b the choice of the sample in a given channel will be a function of a selected subsurface velocity and the distance of the trace from the source of the signal, e.g., the shot. Such expressions are well known in the art of semblance calculation and may be found, for example, in Neidell and Tanner (1979), the disclosure of which is incorporated herein by reference as if fully set out at this point. Then, referring to the semblance coefficient as S_NT = A ∕ B , it follows that:

where the index k is a global index in the temporal domain that indexes the center of the time gate, “N” is the number of traces or channels in the gather, “M” is the window length in samples, and j(i) is the sample on the trace j adjusted by an assumed or trial velocity. Note that to the extent that the quantity SNT is a relatively large value, that indicates that the chosen trial velocity is close to the actual subsurface velocity for the current window center k. [0029] Kirlin (1992) defines the covariance data matrix as the outer product of the data matrix. Yet it is also possible to define a covariance data matrix as the inner product of the data matrix. To distinguish between these approaches, spatial and temporal covariance data matrices are defined as follows. [0030] Forming the outer product, the spatial covariance data matrices can be defined as:

where, on the right-hand-side, the vector products of the row vectors defined above are summed over the number of samples in the temporal window. In this equation, is the Hermitian (i.e., conjugate transpose for complex data items or the vector or mat

rix transpose for real data) of the trace / column vector x_i. This spatial data covariance matrix is N × N in size. Recognizing that it is preferred to define the length of the time gate M as an odd number (although that is not required), in the first summation above the local temporal index ^^ is referenced to the global temporal index k ; in the second summation, the local index is independent of the global index. [0031] Then, to form the inner product to define the temporal covariance data matrix, write

On the right-hand-side of the previous equation, the vector products of the column vectors are summed over the number of spatial samples or traces, N. This temporal data covariance matrix is M × M in size where M is the length of the time gate and where the indices are i = k

with M preferably odd using the local

index temporal frame. [0032] As first shown in Kirlin (1992), and recently in Barros et al (2015), the conventional semblance coefficient, rewritten in matrix form using the spatial covariance data matrix, is

where the identity column vector is ; and where r_{epresents the trace of the sp}

a_{tial covariance data matrix. After eigen decomposition}

to find the eigenvalues of the numerator can be written as the sum over all of the matrix eigenvalues / eigenvector

s:

where u_j is the jth eigenvector of and λ_j is the associated jth eigenvalue.

[0033] The right-hand-side of the previous equation can be separated into signal and noise subspaces by partitioning the eigenvectors / values into two groups. The first group (j = 1, M) corresponds to the largest eigenvalues of (signal) and a second group (j=M+1, N)

corresponds to noise. This equation may then be rewritten as

[0034] Here, the columns of the matrix U are the eigenvectors of the spatial data covariance matrix As will become clear below, the magnitudes of the eigenvalues

associated with the signal separate well from the eigenvalues of the noise. [0035] As suggested in Barros et al (2015, Eq. 24), a temporal covariance semblance coefficient, based on the temporal covariance matrix, is

[0036] Referring back to the fundamental definition of the conventional semblance coefficient, Eq.2, it becomes readily apparent that S_{cov -S} of Eq.5 expands mathematically to the same form. But as S_{cov -T} is expanded for Eq.9, it is just as apparent that this is not the conventional semblance coefficient S_NT . Eq.9 is not equivalent to the conventional semblance definition as multichannel coherence. Instead, it represents coherence among the arrays. [0037] Semblance analysis formed using the spatial data covariance matrix has shown that the information content is primarily contained within the first M eigenvectors. Yet reproducing the analysis using the temporal covariance data matrix is insufficient to correctly produce the conventional semblance map. The eigenvalues of the temporal covariance data matrix should contain all the information necessary to obtain the semblance coefficient, yet there is no explicit equation for S _c as a function of

However, recognize that covariance analysis is a rank reduction process. Irrespective of whether the outer or inner product is calculated, the full information content is preserved in identical sets of reduced- rank eigenvalues. This is verified using non normalized definitions for the inner and outer products. [0038] In an embodiment of the instant fast-spatial covariance, the temporal covariance matrix is used to find the reduced rank eigenvalues, then the spatial covariance matrix is used to recover the corresponding spatial eigenvectors. If an eigen decomposition is performed on the temporal covariance data matrix r, M eigenvalues k_j with M respective eigenvectors v_j are obtained, the latter of size 1 × M . With respect to and the

eigenvalues have the special relationship k_j = λ_j for j = 1, 2, 3, … , M ; however, the eigenvectors are significantly different in size. Therefore, v_j ≠ u_j . However, the set of eigenvectors U can be recovered or calculated through the relationship

where X^T is the transpose of the matrix data matrix X. Note that the previous expression includes an implicit trial velocity term and a window location in time, as do equations (11) through (16) that follow. For the individual eigenvectors u_j, the following expression holds:

Computationally, this is of value in many calculations for seismic velocity analysis if M ≪ N , i.e., if the size of the window in time is much less that the number of traces in the gather which it usually is. By way of example, CMP fold for a 2D line might range from 60 to 240 traces, whereas the length of the analysis window in a velocity analysis might be 15 to 50 or so. [0039] By inspection, the semblance equation for the temporal data covariance matrix can be written as follows. For the numerator, taking j = 1, 2, 3, … , M and calculating u_j from the relationship above

[0040] For the denominator, using the definitions of supra, it follows

therefore,

Therefore,

This equation (14), if implemented will produce a numerical result identical to that of the traditional semblance coefficient calculation. However, it is important to note that the instant approach uses fewer than all, i.e., L < M, of the eigenvectors to form the fast semblance coefficient, S_FS, i.e.,

Preferably, L will be a small integer and more preferably it will be only 1 so that a fast semblance is calculated using a single eigenvector and eigenvalue. In the event that only a single fast semblance coefficient is needed (e.g., if the goal is to produce one of the panels in Fig.9), the previous equation reduces to:

where p ≤ M. Thus, equations (15) and (16) provide a much more efficient method of calculating subspace solutions to obtain a high-resolution semblance based on one or more eigenvalues / eigenvectors at the computational cost of Barros’ (2015) low resolution semblance coefficient. [0041] Characterization of the outer and inner product solutions, respectively, as high resolution and low resolution has led to two misunderstandings. The different approaches to the outer and inner product solutions does not really permit such a direct comparison. The outer product solution presented in Barros et al (2015) benefits from the noise cancellation benefits of spatial averaging. The inner product solution, as presented by Barros, is incomplete and does not acknowledge the distinction of computing an outer product (spatial) semblance versus computing an inner product (temporal) semblance. Semblance in seismic analysis is defined as multichannel coherence, and the inner product covariance does not meet this definition. Mathematically, the expansion of the outer product solution is identical to conventional semblance, as discussed previously. Such is not the case for the inner product solution. By observation the inner product solution, cross correlation across spatial arrays, is different, although the velocity profiles may in fact be very similar and may be the motivation for Barros et al (2015) to characterize the temporal solution as low resolution. Barros et al (2015), without showing the inner product solution, apply the coherence measure with a good result. [0042] Barros et al (2015) uses the properties of the outer and inner products to argue that the largest eigenvalues should be equivalent. It is important to note that the instant embodiment takes a different approach and uses the eigenvectors of the inner product solution to compute the eigenvectors for the outer product solution. The reformulated inner product velocity spectra for the largest eigenvalue as lead to a shaper velocity profile than has been available in the past. This of course is a motivation for wanting to lessen the computational costs of a subspace solution. [0043] This approach shows that by computing the temporal covariance data matrix as an intermediate step (fast-spatial) the spatial covariance eigenvectors can be computed without the computational cost of performing the spatial covariance eigen analysis. Instead the computational cost is greatly reduced and competitive with computing the conventional semblance coefficient. MATLAB run times for the examples here for the conventional semblance coefficient were just over 20 minutes. The run times for the hybrid algorithm to calculate the eigenstructure semblance coefficient took approximately 80 hours. Using the alternative fast-spatial hybrid approach took in the neighborhood of 70 to 90 minutes. This fast- spatial covariance approach should be easily adaptable to the direct approach of Key and Smithson (1990). [0044] There are important lessons with implications for calculating the coherence measure in Barros et al (2015). Barros et al (2015) define an analytical expression for a temporal hybrid approach, but do not show their results. The instant results for a temporal hybrid suggest their result is not a temporal covariance semblance calculation. The relative success of their temporal coherence measure may be a result of the superiority of their coherence measure algorithm despite an incorrect definition for their temporal covariance semblance analysis. The instant results using the fast-spatial algorithm can be applied to the temporal coherence measure; the result should have the resolution of spatial coherence measure at the cost of their current temporal coherence measure. [0045] As compared with the in Barros et al (2015), the instant embodiment uses a shot gather from a long array rather than a CDP gather. Note, that a person of ordinary skill in the would readily be able to convert the above to operations on a CMP gather. Also, the instant approach is less rigorous to properly align our normal moveout. While a phase alignment is possible, a MUSIC algorithm is not applied here. Also, for purposes of discussion only the real form of the data is used herein, while Barros et al (2015) likely used the analytic (or complex) form of the seismic data. [0046] Turning now to the figures, Figs.1A-1D contain an example of a multi- channel marine seismic data set where reflections from the same seismic events as recorded by longer (1A) and shorter (1D) offsets are compared. Fig.1A contains a seismic section created by extracting a far offset channel from multiple shot gathers. Fig.1B contains an example of one such shot gather. Fig.1C contains the same shot gather as 1B and Fig.1D contains a seismic section created from a near offset channel from multiple shots. [0047] Figs.2A-2C contain an example of a conventional semblance velocity analysis. Fig 2A is a conventional semblance velocity spectra, 2B is the total energy, and 2C contains the signal energy, all computed by conventional means. Fig.2A presents the semblance velocity spectra in an image plot of Two-Way Travel Time versus Average Velocity. Recall from the conventional semblance equation that the semblance coefficient is just the nondimensional ratio of the signal energy divided by the total energy which is shown Fig.2B and 2C, respectively. The semblance velocity spectra of Fig.2A is the baseline that will used to compare against the eigen analysis results of the instant invention. [0048] Fig.3 contains an embodiment of a Spatial eigen covariance semblance velocity spectra analysis for sum of all eigenvalues. It should be compared with Fig.2A. [0049] Fig.4 contains an embodiment of a temporal eigen covariance semblance velocity spectra analysis for sum of all eigenvalues. It should be compared with Figs.2A and 3, the conventional and spatial covariance velocity spectra maps respectively. As can be seen, the results of the temporal covariance are a significant departure from the other solutions. This difference can be seen in the definition of semblance, which is multichannel coherence or correlation across the seismic traces. Temporal covariance, however, is array coherence, or correlation across the arrays. [0050] Figs.5A-5D provide the scree plots of the eigenvalue distribution in both linear and semi-log y plots. For display the eigenvalues are rank ordered from largest to smallest, then displayed with the value on the y-axis and the rank-ordered index on the x- axis. Figure 5A, displaying all 636 eigenvalues, shows that almost all of the values are close to zero. With Figure 5B a closer view is provided of only the leading 25 eigenvalues. With some experimentation it was determined that eigenvalues above index 25 were several orders of magnitude smaller, and so they were displayed separately in Figure 5C. The eigenvalues above index 25 are more than 15 orders of magnitude smaller than the largest eigenvalue as shown in Figure 5A. In Figure 5D the higher order eigenvalues are displayed with a semi-log y plot to reveal a linear drop off in magnitude. [0051] To further investigate the results in Fig.5C, additional plots are provided in Figs.5E and 5F. To more closely assess symmetry, the eigenvalues of Figure 5C were displayed by taking the absolute values. This result is shown in Figure 5E. In looking at Figure 5E recognize that eigenvalues above index 349 are actually negative as can be seen from Figure 5C. Finally, in Figure 5F, again taking the absolute values, all eigenvalues are displayed on a semi-log y plot. [0052] Figs.6A-6D compare the spatial and temporal eigenvectors before translation of the temporal eigenvectors. In Fig.6, the eigenvalues from the spatial-covariance eigen analysis, Fig.6A, can be compared with the eigenvalues from the fast-spatial covariance eigen analysis, Fig.6C. The spatial-covariance eigen analysis provides the 636 eigenvalues described earlier, whereas the temporal covariance eigen analysis provides only 25 eigenvalues. Although this may not be readily apparent from these subplots, the 25 eigenvalues from the temporal analysis are identical with the first 25 eigenvalues from the spatial analysis. The spatial information contained within the eigenvectors from the spatial analysis is much higher resolution than the spatial information contained within the eigenvectors for the temporal analysis. This is readily apparent by comparing the first three eigenvalues in Fig.6B with those in Fig.6D. [0053] Figs.7A-7D compare the spatial and temporal eigenvectors after translation. After the translation relationship was applied to the temporal eigenvectors the results are displayed again. The scree plots for eigenvalues 1 thru 25 for both the spatial and temporal analysis are compared in Figs.7A and 7B, respectively. Likewise, the eigenvectors from the spatial and temporal analysis in Figs.7B and 7D, respectively, are compared. From these displays it can be seen the results are essentially identical. [0054] The eigenvectors of the spatial covariance analysis have been successfully generated, for the reduced-rank eigenvalues, by using the data matrix to translate the eigenvectors from the temporal covariance analysis. This was accomplished at a much lower computational cost than performing the full spatial covariance analysis. This is referred to herein as fast-spatial covariance. [0055] Figs.8A-8B contain an exemplary embodiment of a full spectrum semblance eigen analysis showing a spatial covariance semblance velocity spectra (Fig.8A) and a temporal covariance semblance velocity spectra (Fig.8B). The results were computed by summing the results of the 25 individual eigen components. [0056] Figs.9A-9D contain an exemplary embodiment of a Temporal covariance eigen-Semblance velocity spectra analysis results associated with A) eigenvalue 1; B) eigenvalue 2; C) eigenvalue 3; and D) eigenvalue 4. In Fig.9A the semblance velocity spectra map from the largest eigenvalue is displayed. Figs.9B, 9C, and 9D show the eigen- Semblance velocity spectra associated with eigenvalues 2, 3, and 4, respectively. Careful inspection shows sharper velocity profiles associated with the first two. Since the spatial and fast-spatial analysis results are essentially identical only the results are present for the fast- spatial analysis. [0057] Barros et al (2015) shows how spatial covariance analysis can be applied to compute the conventional semblance coefficient with the hybrid approach from Kirlin (1992). They also define a coherence measure for a direct application of spatial covariance analysis. An important contribution of Barros et al (2015) is the attempt to apply temporal covariance analysis to the hybrid approach and to define the modified coherence measure. Shown here as Equation 9, the analytic equation shown in Barros et al (2015) does not successfully define a temporal covariance equation for semblance. Here one embodiment of the instant approach takes advantage of the properties of rank reduction with covariance analysis to apply a reformulated temporal covariance data analysis (Equation 14) to the hybrid approach. This analysis has important implications for the modified coherence measure of Barros et al (2015). Barros et al (2015) show the value of the modified coherence measure. Even when using a temporal covariance analysis that is deficient, and likely incorrect, the velocity spectra look reasonable though degraded from their spatial covariance coherence measure. By adopting the reformulated temporal covariance approach taught herein as an intermediate step to efficiently calculate the spatial covariance eigenvectors, both the efficiency of the temporal covariance algorithm and the resolution of the spatial covariance algorithm can of the temporal covariance algorithm and the resolution of the spatial covariance algorithm can be realized. There are several important things to note: 1) The spatial covariance semblance algorithm summed for all reduced-rank eigenvalues accurately reproduces the results of the conventional semblance algorithm. No data smoothing, no partial stacking, and no special algorithms to recover the lowest order eigenvalues (much less preconditioning of data) are needed. 2) The temporal covariance semblance algorithm is not a true semblance since it measures coherence across arrays rather than coherence across channels and may inappropriately be seen as a low resolution semblance measure and does not constitute a semblance calculation. 3) The reduced-rank eigenvectors computed from the temporal covariance matrix can be multiplied by the data matrix to generate the eigenvectors that would be for the spatial covariance semblance without the high computational cost. 4) The largest eigen-component can be used to generate a seismic velocity spectra that is of a higher resolution than the full reduced-rank spectra spatial covariance semblance equivalents or the conventional semblance whether the spatial covariance or conventional semblance are used in its construction and does so at a much lower computational cost that is referred to herein as fast-spatial covariance. 5) That information can be recovered, with minimal computational cost, using the matrix relationships between the spatial and temporal domains suggested by Barros et al (2015). Thus, and according to the present embodiment, it has been shown herein that the properties of rank reduction with covariance analysis allows for the successful use of a reformulated temporal (fast- spatial) covariance analysis to compute the semblance with eigen analysis at much greater efficiency. It is applicable to both hybrid and direct approaches. Improvements in the normal moveout alignment are necessary to realize the full value of the instant approach. [0058] Barros et al (2015) have made an important contribution with their modified spatial coherence measure for calculating the seismic velocity spectral map, but at great computational cost. Their attempt to develop the modified temporal coherence measure seems successful, particularly since the computational costs are more competitive with the computational cost of the conventional semblance coefficient but results in a low-resolution solution. There are two short comings with their results. It is not appropriate to compare their results to refer to the modified spatial coherence as high resolution and the modified temporal coherence as low resolution since the spatial approach uses smoothing and spatial averaging to improve the results. Secondarily, their temporal covariance is not a true semblance because the cross correlations are not performed between traces. Our solution gives a high resolution results with their temporal times. [0059] The fast-spatial analysis taught herein is an intermediate step and improvement in the approach of Barros et al (2015). This should achieve the competitive computational cost of the temporal-covariance analysis; after applying the reformulation to calculate the spatial eigenvectors the results should have the resolution of the modified spatial coherence measure. [0060] Using the fast-spatial covariance analysis that is competitive with the conventional semblance coefficient, it is now feasible to further investigate the inherent value of the coherence measure. Spatial averaging can be used to further improve the result. [0061] Note that the mathematical algorithm taught herein doesn't change for a shot gather versus a CDP gather. The main practical use will be for CDP gathers but the algorithm does not change. Both Kirlin (1992) and Barros et al (2015) show the spatial covariance semblance provides a higher resolution result than the conventional semblance. And if the sum up all spectral components of the spatial covariance semblance result is calculated a conventional semblance result is obtained. [0062] Turning now to Fig.10 and the expanded textual representation of that figure below, below is provided a generalized flow chart suitable for use with an embodiment. Note that in this embodiment of the calculation of a fast spatial semblance coefficient only L eigenvalues and eigenvalues will be used in the calculation with L preferably much smaller than M. If all M eigenvalues are calculated and used, the resulting semblance coefficient will be equal to the conventionally calculated value. Collect seismic survey Chose minimum and maximum velocity for analysis, number of velocities, NV Choose window (gate) size M DO {Number of gathers (CMP or shot) to analyze} Select next gather

ENDDO {Number of gathers (CMP or shot) to analyze} [0063] As can be seen, the starting point is the acquisition of a 2D or 3D seismic survey (box 1005). The methods by which such are obtained are well known to those of ordinary skill in the art. Next, range of velocities will need to be specified as well as the number of velocities that are to be used (box 1008). When a semblance plot is calculated that diagram covers a range of different velocities (e.g., see, Fig.2A, where the range is 1200 m/s to 1800 m/s). Additionally, an analysis window size in time will need to be specified, preferably in terms of samples. The window size might be any number of samples but in some cases, it might be chosen to be approximately the same as the pulse length of the seismic source (e.g., 50ms or 25 samples) but there is no particular reason it needs to be that length. In some embodiments the window length might be between about 15 and 50 samples, although other window sizes are certainly possible and well within the ability of one of ordinary skill in the art to select. [0064] Box 1010 is reached if the current CMP gather has been fully processed. Of course, those of ordinary skill in the art will recognize that in a typical velocity analysis many gathers along a seismic line or within a 3D volume will be selected so that an estimate of the velocity variation laterally can be obtained. Fig.10 illustrates only how one such gather might be processed but it should be clear how the instant approach could be generalized to any number of gathers. Note that when the term “gather” is used herein, that term should be understood to include a collection of seismic traces such as conventional CMP or shot gathers, as well as gathers that are constructed by using traces from multiple conventional gathers, gathers that might be created by summing together or removing traces from such a gather, etc. [0065] As is indicated in box 1013, during the first pass a starting trial velocity will be selected. This would typically be the lowest velocity in a chosen velocity range under consideration. For simplicity in the discussion that follows and for purposes of illustration only, it will be assumed that the chosen gather will have a constant velocity moveout correction applied (box 1015) at the starting trial velocity. Note that this formulation is discussed only purposes of illustrating more clearly the claimed method and simplifying the discussion that follows. It is not necessary that the gather have moveout pre-applied. As an example, rather than applying a constant velocity moveout to the entire gather, a hyperbolic path may be calculated through traces to determine the window start or center and the time shifts associated with each trace as a function of offset from the shot. Then data values from each trace can be extracted given the calculated path. However, for purposes of the claims that follow, when it is said that a trial velocity is selected and applied to a gather that statement should be construed in the broadest sense to mean that a trial velocity is used to guide the selection of the analysis window on each trace, whether moveout is applied to each trace, time shifts for each trace are calculated using the hyperbolic equation, or any other scheme by which the appropriate analysis windows are identified for each trace in the gather. In Fig.10, the matrix X in box 1015 is the windowed data matrix (M x N) after moveout at the current trial velocity. [0066] Subsequent passes through box 1013 will preferably select the next higher (or lower) trial velocity value (e.g., the starting velocity will be incremented or decremented) and the current trial velocity will be applied to the original gather by, for example, moveout at that velocity so that a new X matrix can be identified and/or filled. This will be repeated until the highest (or lowest) velocity in the specified range is reached. Once all of the velocities in the designated range have been processed, typically another gather will be selected (box 1010) elsewhere along the line or within the 3D volume via decision item 1035. The next gather will then be analyzed the same way by resetting the trial velocity to its first or initial velocity value and resetting the window location to its first or initial value. [0067] In box 1015, the starting sample of the first analysis center will be selected during the first pass. In subsequent passes, other windows will be selected, and this process will continue until the selected time range has been covered. Note that in some cases, k might be increased by the window size each iteration, e.g., if the window is 25 samples in length k might be 1000, 1025, 1050, etc., until the maximum window position is analyzed. In other cases, there might be some overlap between the successive windows, e.g., window starts of 1000, 1020, 1040, etc., with a window 25 samples long. As noted previously, successive analysis windows will be selected until the maximum analysis time is reached. For example, in Fig.2A the range of windows covers 5000ms to 7000ms. Selection of the length of each window, the number of windows, the amount of window overlap, and filtering (or not) of the window contents is well within the skill of one of ordinary skill in the art and such selections are routinely made by seismic processors. As is indicated by decision item 1033, once the last time window has been processed in this example another trial velocity will be selected and the window location will be reset to its starting value. [0068] Next, and preferably, the spatial vectors x_i will be built (box 1018). These vectors are comprised of one sample from each of the traces (or channels) in the gather as has been discussed previously, where the sample chosen from each trace is determined as a function of the selected trial moveout velocity and current window location with each trace being shifted in time by the moveout step as compared with the others. Additionally, as is indicated in box 1020 the temporal vectors will also be constructed, i.e., a time window will be extracted from each trace in the gather. Note that each vector is a time window from one of the traces and the actual window that is chosen will be a function moveout velocity that pass. So, note that when the terms “build” or “construct” or “form” or “assemble” are used those terms should be broadly construed to include instances where the data values are moved in memory from the gather to a separate vector as well as instances where the data values are not moved but accessed in place within each seismic trace. That being said, creation of separate arrays tends to be more efficient in many implementations. [0069] Box 1023 describes the calculation of the and matrices from the

temporal and spatial vectors. In box 1025, the eigenvalues and eigenvectors of will be calculated using any appropriate numerical method. Note that if the end goal is

to calculate a single high-resolution semblance, only a single eigenvalue / eigenvector needs to be obtained. If that is the case, such would substantially speed up the calculation. On the other hand, since the M largest eigenvalues of are equal to the eigenvalues of the eigenvalues of can ed to calculate any of the

be us eigenvectors of using a well-known relationship from

linear algebra.

[0070] However, preferably as many of the eigenvectors of as are needed will be calculated using equation (11) supra (box 1030). Note as is indicat

ed in box 1030 that in most embodiments only a single L, which is less than M, eigenvalues / eigenvectors will be used to calculate the semblance. For example, Fig.9A indicates the calculated semblance that results if only one eigenvector / eigenvalue is used. In this case the single eigenvalue is the largest one of the M eigenvalues. In other embodiments, only the second largest, third largest, etc., eigenvalues / vectors might be used. See, e.g., Figs.9A – 9D. In still other cases, the semblance associated with two or more arbitrary eigenvalues / eigenvectors might be summed together, e.g., the first and third, the second and fifth, the first, third, and sixth, etc. Those of ordinary skill in the art will recognize that the instant method provides a rapid way to calculate single-eigenvector / eigenvalue semblances and those semblances might be used individually or combined in arbitrary ways. This embodiment one example of the instant fast (or high resolution) semblance calculation. However, it should be noted that if all M of the eigenvectors/eigenvalues are used in equation (15) the resulting semblance will be the same as the conventionally calculated display (e.g., compare Figs.2A and 3), albeit at a greatly increased computational cost. [0071] It should be clear that many different time windows will need to be extracted and the calculations performed for a given trail velocity in order to form a single column of the display in, for example, Fig.3. Further, many different trail velocities will need to be used in order to form a full panel of semblance coefficients. [0072] After the semblance plot is obtained, in most applications the time versus velocity curve associated with the highest data values will be “picked” (chosen) to obtain a velocity profile for the selected gather. Then, a seismic processor can use this time-varying velocity information in subsequent processing steps to obtain a more accurate picture of the subsurface of the earth that is beneath the seismic survey. For example, such velocity information might be used to improve CMP stacks or create a velocity field for use in seismic migration to assist in the location of subsurface hydrocarbons. A potential well location could be selected using the improved seismic image and, based on the seismic and other measurements, a well might be drilled at the selected well location. [0073] By way of summary, the instant methods make it possible to replace calculation of semblance with an equivalent higher resolution covariance approach at a comparable computational cost that involves calculating the largest eigenvalue of the temporal covariance matrix, then using the data matrix to calculate the corresponding spatial eigenvector. This can also be done using any number of smaller eigenvalues with the resulting semblance estimates optionally being combined via a weighted average or some other means to produce a more accurate calculation, a more informative velocity display, and, ultimately, better velocity picks which will result in a better representation of the subsurface of the earth that lies beneath the seismic survey. Better seismic images will be of value to individuals who engage in the exploration for hydrocarbons in the subsurface and can potentially reduce the risk associated with drilling a well. [0074] The foregoing has outlined in broad terms some of the more important features of the invention disclosed herein so that the detailed description that follows may be more clearly understood, and so that the contribution of the instant inventors to the art may be better appreciated. The instant invention is not to be limited in its application to the details of the construction and to the arrangements of the components set forth in the following description or illustrated in the drawings. Rather, the invention is capable of other embodiments and of being practiced and carried out in various other ways not specifically enumerated herein. Finally, it should be understood that the phraseology and terminology employed herein are for the purpose of description and should not be regarded as limiting, unless the specification specifically so limits the invention. [0075] One important insight of provided by the instant approach is the recognize that covariance analysis is a rank-reduction process that preserves all information. It doesn't matter whether the temporal or spatial covariance approach is used until the physics are considered. Then you have to find the spatial eigenvectors that is meaningful for semblance. That is one important aspect of this invention. [0076] It is to be understood that the terms "including", "comprising", "consisting" and grammatical variants thereof do not preclude the addition of one or more components, features, steps, or integers or groups thereof and that the terms are to be construed as specifying components, features, steps, or integers. [0077] The singular shall include the plural and vice versa unless the context in which the term appears indicates otherwise. [0078] If the specification or claims refer to "an additional" element, that does not preclude there being more than one of the additional element. [0079] It is to be understood that where the claims or specification refer to "a" or "an" element, such reference is not to be construed that there is only one of that element. [0080] It is to be understood that where the specification states that a component, feature, structure, or characteristic "may", "might", "can" or "could" be included, that particular component, feature, structure, or characteristic is not required to be included. [0081] Where applicable, although state diagrams, flow diagrams or both may be used to describe embodiments, the invention is not limited to those diagrams or to the corresponding descriptions. For example, flow need not move through each illustrated box or state, or in exactly the same order as illustrated and described. [0082] Methods of the present invention may be implemented by performing or completing manually, automatically, or a combination thereof, selected steps or tasks. [0083] The term "method" may refer to manners, means, techniques and procedures for accomplishing a given task including, but not limited to, those manners, means, techniques and procedures either known to, or readily developed from known manners, means, techniques and procedures by practitioners of the art to which the invention belongs. [0084] For purposes of the instant disclosure, the term “at least” followed by a number is used herein to denote the start of a range beginning with that number (which may be a ranger having an upper limit or no upper limit, depending on the variable being defined). For example, “at least 1” means 1 or more than 1. The term “at most” followed by a number is used herein to denote the end of a range ending with that number (which may be a range having 1 or 0 as its lower limit, or a range having no lower limit, depending upon the variable being defined). For example, “at most 4” means 4 or less than 4, and “at most 40%” means 40% or less than 40%. Terms of approximation (e.g., “about”, “substantially”, “approximately”, etc.) should be interpreted according to their ordinary and customary meanings as used in the associated art unless indicated otherwise. Absent a specific definition and absent ordinary and customary usage in the associated art, such terms should be interpreted to be ± 10% of the base value. [0085] When, in this document, a range is given as “(a first number) to (a second number)” or “(a first number) – (a second number)”, this means a range whose lower limit is the first number and whose upper limit is the second number. For example, 25 to 100 should be interpreted to mean a range whose lower limit is 25 and whose upper limit is 100. Additionally, it should be noted that where a range is given, every possible subrange or interval within that range is also specifically intended unless the context indicates to the contrary. For example, if the specification indicates a range of 25 to 100 such range is also intended to include subranges such as 26 -100, 27-100, etc., 25-99, 25-98, etc., as well as any other possible combination of lower and upper values within the stated range, e.g., 33-47, 60- 97, 41-45, 28-96, etc. Note that integer range values have been used in this paragraph for purposes of illustration only and decimal and fractional values (e.g., 46.7 – 91.3) should also be understood to be intended as possible subrange endpoints unless specifically excluded. [0086] It should be noted that where reference is made herein to a method comprising two or more defined steps, the defined steps can be carried out in any order or simultaneously (except where context excludes that possibility), and the method can also include one or more other steps which are carried out before any of the defined steps, between two of the defined steps, or after all of the defined steps (except where context excludes that possibility). [0087] Further, it should be noted that terms of approximation (e.g., “about”, “substantially”, “approximately”, etc.) are to be interpreted according to their ordinary and customary meanings as used in the associated art unless indicated otherwise herein. Absent a specific definition within this disclosure, and absent ordinary and customary usage in the associated art, such terms should be interpreted to be plus or minus 10% of the base value. [0088] The invention is described herein with a certain degree of particularity, but is should be understood that the invention is not limited to the embodiment(s) set forth herein for purposes of exemplification. [0089] Thus, the present invention is well adapted to carry out the objects and attain the ends and advantages mentioned above as well as those inherent therein. While the inventive device has been described and illustrated herein by reference to certain preferred embodiments in relation to the drawings attached thereto, various changes and further modifications, apart from those shown or suggested herein, may be made therein by those of ordinary skill in the art, without departing from the spirit of the inventive concept the scope of which is to be determined by the following claims.

REFERENCES Abbad B. and Ursin B.2012. High-resolution bootstrapped differential semblance. Geophysics, 77 (3), U39-U47. Asgedom E. G., Gelius L. J. and Tygel M.2011. High-resolution determination of zero-offset common-reflection-surface stack parameters. International Journal of Geophysics., doi:10.1155/2011/819831, 10 pp. Barros T., Lopes R. and Tygel M.2015. Implementation aspects of eigendecomposition- based high-resolution velocity spectra. Geophysical Prospecting, 63, p.99-115. Biondi B. L. and Kostov C.1989. High-resolution velocity spectra using eigenstructure methods. Geophysics, 54 (7), 832-842. Key S. C. and Smithson S. B.1990. New approach to seismic-reflection event detection and velocity determination. Geophysics, 55, (8), 1057-1069. Kirlin R. L.1992. The relationship between semblance and eigenstructure velocity estimators. Geophysics, 57 (8), 1027-1033. Li F. and Liu H.1999. Subspace based seismic velocity analysis (Chapter 7). In: Kirlin, R. L.; Done, W. J.; Editors, 1999. Covariance Analysis for Seismic Signal Processing. Geophysical Developments Series, No.8, Society of Exploration Geophysicists, pp.109140. Luo S. and Hale D.2012. Velocity analysis using weighted semblance. Geophysics, 77 (2), U15U22. Neidell N. S. and Taner M. T.1971. Semblance and other coherence measures for multichannel data. Geophysics, 36 (3), 482-497. Ng D. and Martinson D.2011. RV Langseth Data Reduction Summary, MGL1111, Dutch Harbor, AK – Dutch Harbor, AK, Not Final, V1.2, 2011-09-02. Lamont Doherty Earth Observatory, Columbia University, 64 pp, +5 Appendices. Scholl D. W., Wood W. T., Barth G. A. and Childs J. R.2012. The Bering Sea Basin: New drilling and geophysical observations and evidence for an important contribution of thermogenic methane to interstitial deposits of methane hydrate. Eleventh International Conference on Gas in Marine Sediments, Nice, France, 4-7 September 2012, pp.74.-77 Shan, T.-J., Wax M. and Kailath, T.1985. On spatial smoothing for direction-of-arrival estimation of coherent signals. IEEE Transactions on Acoustics, Speech, and Signal Processing, ASSP-33 (4), 806-811. Wang Y. Y., Chen J.-T. and Fang W. H.2001. TST-MUSIC for joint DOA-delay estimation. IEEE Transactions on Signal Processing, 49 (4), 721-729. Williams R. T, Prasad S. and Mahalanabis A. K.1988. An improved spatial smoothing technique for bearing estimation in a multipath environment. IEEE Transactions on Acoustics, Speech, and Signal Processing, ASSP-36 (4), 425-432. IEEE Transactions on Acoustics, Speech, and Signal Processing, ASSP-36 (4), 425-432.

Claims

CLAIMS What is claimed is:

1. A method of seismic exploration, comprising the steps of:

(a) accessing a seismic data set, said seismic dataset comprising at least one seismic gather of traces, each of said seismic gather of traces comprising a plurality of seismic traces;

(b) selecting a seismic gather of traces from said seismic data set and a selected plurality of seismic traces associated therewith;

(c) selecting a window length, M

(d) selecting a trial velocity value;

(e) applying the trial velocity value to said selected gather of seismic traces, thereby forming a velocity adjusted seismic gather containing a plurality of velocity adjusted data values x_i,j;

(f) selecting a window starting location k;

(g) forming a data matrix X from said velocity adjusted gather of seismic traces and said plurality of velocity adjusted data values, where

where A is a number of said plurality of seismic traces associated with said selected seismic gather;

(h) building M spatial vectors x_i, where

(i) building N temporal vectors where

(j) calculating a temporal covariance data matrix from said temporal vectors according to the equation

where is said temporal covariance data matrix;

(k) calculating a spatial covariance data matrix from said spatial vectors according to the equation

where is a transpose conjugate of said spatial vector x_i, and

where

is said spatial covanance data matrix;

(l) calculating one or more eigenvalues k_h of for n=1, L:

where L is a number of said one or more eigenvalues if k_n , and where L <M;

(m) calculating L eigenvectors v_n of

where n=1, L;

(n) calculating L eigenvectors, U_j of

according to the equation:

X^T is a transpose of said data matrix X;

(o) calculating a fast semblance coefficient value, SFS, according to the equation:

where

SFS is said fast semblance coefficient value, is a trace of said temporal covariance matrix,

1_M is a row vector of ones of length and is a row vector of ones of length N;

(p) performing at least steps (f) through (o) a plurality of times for a plurality of different window starting locations to produce a plurality of fast semblance coefficient values at said selected tnal velocity;

(q) performing at least steps (d) through (p) for a plurality of different selected trial velocities thereby producing said plurality of fast semblance coefficient values for each of said selected velocities;

(r) forming said plurality of fast semblance coefficient values associated with for each of said selected velocities into a semblance display;

(s) using said semblance display to select a plurality of subsurface velocity estimates;

(t) using said plurality of subsurface velocity estimates to process at least a portion of said seismic data set to obtain an image representative of a subsurface of the earth;

(u) using said image representative of the subsurface of the earth in a seismic exploration program to locate subsurface hydrocarbon; and (v) determining a well location based at least in part on said image representative of the subsurface of the earth.

2. The method of claim 1, wherein L is equal to 1.

3. The method according to claim 1, wherein Mis between 15 and 50 samples in length.

4. The method of claim 1, further comprising the step of:

(w) drilling a well at said determined well location.

5. A method of seismic exploration, comprising the steps of:

(a) accessing a seismic data set, said seismic dataset comprising at least one seismic gather of traces, each of said seismic gather of traces comprising a plurality of seismic traces;

(b) selecting a seismic gather of traces from said seismic data set and a selected plurality of seismic traces associated therewith;

(c) selecting a window length, M;

(d) selecting a trial velocity value;

(e) selecting a window starting location k,

(f) using at least said window starting location, said selected trial velocity value, and said window length to form a temporal covariance data matrix, where is

said temporal covariance data matrix;

(g) using at least said window starting location, said selected velocity value, and said window length to form a spatial covariance data matrix, where is said

spatial covariance data matrix;

(h) calculating one or more eigenvalues of

(i) calculating one or more eigenvectors of

(j) using said one or more eigenvectors of and said one or more eigenvalues of

to calculate one or more eigenvectors of

(k) using said calculated one or more eigenvectors of and said one or more

eigenvalues of to calculate a fast semblance coefficient value;

(l) performing at least steps (e) through (k) a plurality of times for a plurality of different window starting locations k to produce a plurality of fast semblance coefficient values at said selected velocity; (m) performing at least steps (d) through (1) for a plurality of different selected trial velocities thereby producing a different plurality of fast semblance coefficient values at each of said selected trial velocities;

(n) forming said different plurality of fast semblance coefficient values at each of said selected trial velocities into a semblance display;

(o) using said semblance display to select a plurality of subsurface velocity estimates;

(p) using said plurality of subsurface velocity estimates to process at least a portion of said seismic data set to obtain an image representative of a subsurface of the earth; and

(q) using said image representative of the subsurface of the earth in a seismic exploration program to locate subsurface hydrocarbons.

6. The method of claim 5, further comprising the step of:

(r) determining a well location based at least in part on said image representative of the subsurface of the earth.

7. The method of claim 5, further comprising the step of:

(s) drilling a well at said determined well location.

8. The method of claim 5, wherein the step of selecting a window starting location k comprises the step of:

(1) forming a data matrix X from said velocity adjusted gather of seismic traces and said plurality of velocity adjusted data, where

where X is said velocity adjusted gather of seismic traces and N is a number of said plurality of seismic traces associated with said selected seismic gather;

(2) building M spatial vectors x;, where

(3) building N temporal vectors x

where

9. The method of claim 8, wherein the step of using at least said window starting location, said selected trial velocity value, and said window length to form a temporal covariance data matrix, comprises the step of calculating said temporal covariance data matrix from said temporal vectors according to the equation

where is said temporal covariance data matrix;

10 The method of claim 9, wherein the step of using at least said window starting location, said selected trial velocity value, and said window length to form a spatial covariance data matrix, comprises the step of calculating said spatial covariance data matrix from said spatial vectors according to the equation

where

is a transpose conjugate of said spatial vector x_i. and where

is said spatial covanance data matrix;

11 The method of claim 10, wherein the step of calculating one or more eigenvalues of comprises the step of calculating one or more eigenvalues k_h of for n=1, L:

where L is a number of said one or more eigenvalues if k_h , and where L ≤ M;

12 The method of claim 11, wherein the step of calculating one or more eigenvectors of comprises the step of calculating one or more eigenvectors v_n of

where n=1, L;

13. The method of claim 12, wherein the step of using said one or more eigenvectors of and said one or more eigenvalues of to calculate one or more eigenvectors of

comprises the steps of: calculating L eigenvectors, u_j of

according to the equation:

X^T is a transpose of said data matrix X;

14. The method of claim 13, wherein the step of using said one or more eigenvectors of and said one or more eigenvalues of

to calculate a fast semblance coefficient value comprises the step of calculating a fast semblance coefficient value, S_FS, according to the equation:

where

S_FS is said fast semblance coefficient value, is a trace of said temporal covariance matrix,

1_M is a row vector of ones of length M, and

is a row vector of ones of length N;

15. The method of claim 11, wherein L is equal to 1.

16. The method according to claim 5, wherein Mis between 15 and 50 samples in length.

17. The method of claim 6, further comprising the step of:

(s) drilling a well at said determined well location.