Parsimonious refraction interferometry and tomography

S U M M A R Y We present the theory of parsimonious refraction interferometry and tomography where a densely populated refraction data set can be obtained from two reciprocal and several infill shot gathers. The assumptions are that the refraction arrivals are head waves, and a pair of reciprocal shot gathers and several infill shot gathers are recorded over the line of interest. Refraction traveltimes from these shot gathers are picked and spawned into O(N2) virtual refraction traveltimes generated by N virtual sources, where N is the number of geophones in the 2-D survey. The virtual traveltimes can be inverted to give the velocity tomogram. This enormous increase in the number of traveltime picks and associated rays, compared to the many fewer traveltimes from the reciprocal and infill shot gathers, allows for increased model resolution and a better condition number with the system of normal equations. A significant benefit is that the parsimonious survey and the associated traveltime picking can be an orderof-magnitude less time consuming than that for a standard refraction survey with a dense distribution of sources.


I N T RO D U C T I O N
Refraction tomography is an important imaging tool in earthquake studies (Stein & Wyession 2003), crustal-mantle imaging (Prodehl & Mooney 2012), exploration geophysics (Bishop et al. 1985;Yilmaz 2001) and engineering seismology (Yilmaz 2015).For 2-D engineering seismology, receivers are deployed along a line and common shot gathers (CSGs) are recorded for sources at selected positions along the line.To save costs, a parsimonious survey is carried out where sources are located only at each end of the receiver line as shown in Fig. 1.The first breaks are picked from each shot gather and the traveltimes are inverted by a simple formula that assumes a layered medium with an unknown dip angle for each interface.These two shot gathers are denoted as a pair of reciprocal shot gathers.We often require several short-offset infill shot gathers that can be used to find the near-offset and direct traveltimes.Therefore, parsimonious surveys can significantly reduce acquisition costs and survey efforts, but at the cost of less slowness resolution and certainty in the estimate of the subsurface model.
For a 2-D survey, we now show that two reciprocal shot gathers and several short-offset infill shot gathers can be decomposed into many virtual shot gathers that have many more traveltime picks available for inversion.In some cases, the virtual shot gathers give as much refraction information as a full survey with N shots, where a shot is located at each of the N geophone locations.The result is a tomogram with a much denser ray coverage and better slowness resolution compared to that from the original data.We call this procedure parsimonious refraction interferometry because it uses a stationary-phase principle to decompose the reciprocal traveltimes into many virtual traveltimes associated with shorter ray paths.Unlike the original pair of reciprocal shots, the virtual shot locations are at all of the interior geophone locations.
The next section describes the theory of parsimonious refraction interferometry.It is an application of Fermat's interferometric traveltime principle (Schuster 2005) and closure phase (Schuster et al. 2014) that allows for the decomposition of long ray paths and traveltimes into, respectively, shorter ray paths and traveltimes.Instead of body-wave traveltimes we now apply it to head-wave traveltimes.Section 3 presents the results of applying parsimonious refraction interferometry to both synthetic data and field data.The field data are from 2-D surveys conducted next to the Gulf of Aqaba and over a basin in East Africa.The final section presents a summary and conclusions.

THEORY
Assume two reciprocal sources (A and D) and the irregularly layered medium in Fig. 1, where head waves propagate along the interface between the upper and lower layers.There can be lateral velocity variations in the upper medium and there are N evenly spaced geophones on the recording surface between the two sources.The head-wave traveltime from the source at A to the geophone at C is defined as and the reciprocal traveltime from D to B is where τ xx is the traveltime from x to x along the refraction ray, and reciprocity says τ xx = τ x x .The traveltimes τ AC and τ DB are denoted as a pair of reciprocal traveltimes.

Virtual-traveltime equation
To create virtual sources and receivers within the array in Fig. 1, we assume that there is a post-critical distance between the geophone positions B and C, and both traveltimes τ AC and τ DB are associated with the same refractor.Subtracting the reciprocal traveltime τ AD = τ Ax + τ x x + τ x D from the sum τ AC + τ DB gives the virtual traveltime δτ CB : (3)   Here, δτ CB is denoted as an interferometric stationary traveltime because the reciprocal ray path Ax x D, marked by the dashed red ray in Fig. 1, cancels the phase associated with the common ray paths of the purple Ax xC and green Dx x B rays.The result is the virtual traveltime δτ CB associated with the much shorter ray path C x x B denoted by the dashed blue ray.Equivalently, δτ CB is the traveltime of a virtual refraction excited by a virtual source at C and recorded at B. This natural redatuming operation is the key principle underlying seismic interferometry (Snieder 2004;Wapenaar 2004).Eq. ( 3) can be used to generate O(N) virtual shot gathers of head-wave traveltimes, where the number of geophone pairs with post-critical separation is nearly equal to the number N of geophones in the survey, assuming that a few geophones are recording direct waves as first arrivals.Each virtual shot gather will, on average, contain O(N) virtual traveltimes to create O(N 2 ) virtual refraction traveltimes from the 2N traveltimes in two reciprocal shot gathers.This abundance of new traveltimes can be used to invert for the subsurface velocity model with much greater ray density and better model resolution compared to that from the original data set.Eq. (3) satisfies Fermat's interferometric principle (Schuster 2005) because the subtraction of τ AD (red dashed ray) from τ AC + τ DB (solid green and purple rays) gives the same value of δτ CB for all post-critical, i.e. stationary, locations of the reciprocal sources.Eq. ( 3) is also a special case of the closure phase condition discussed in Schuster et al. (2014).
The reciprocal shots are typically too far apart from one another, so the first-arrivals will not be direct waves for far-offset receivers.These direct waves and their traveltimes are needed to estimate the velocities at the near-surface.Therefore, several infill shots are required to get the near-surface velocities from the direct-wave traveltimes.The direct-arrival traveltimes are picked from these infill shot gathers, and then interpolated to other source locations to fill in the missing direct-wave traveltimes in the other CSGs.The number of infill shot gathers depends on the complexity of the subsurface.We suggest that the offset between the shots of the infill shot gathers should be less than the cross-over distance offset of the direct-wave/refraction arrivals in the recorded data.However, in the case of a simple near-surface velocity distribution, this rule can be relaxed and the offset between infill shot gathers can be larger than the cross-over distance.

Condition number, illumination area and slowness uncertainty
The resolution properties associated with the virtual-traveltime equations are now compared to those from the reciprocal traveltime equations.These properties are (1) the condition number of the normal equations, (2) the enhancement of illumination with virtual-traveltime equations and (3) the slowness uncertainty associated with errors in traveltime picking.

Condition number
As Appendix A shows, the extra traveltime equations generated by eq. ( 3) act as a preconditioner that can reduce the condition number of the normal matrix [L T L] by a factor of 3 in some cases, where L is the ray path matrix associated with the traveltime equations.The dense set of virtual traveltimes does not create new information about the subsurface because it is created by adding and subtracting traveltime equations from the reciprocal data.However, the synthetic examples in the next section suggest that the enhanced subsurface illumination and smaller condition number in the virtual [L T L] matrix can significantly reduce the number of velocity models that can approximately satisfy the same data.For refractions, these gaps can be caused by the inability to pick far-offset traveltimes as indicated by the truncated green lines in Fig. 2(a).This problem can be partly remedied if the later-arriving direct waves and refractions can be picked, as illustrated by the dashed lines in Fig. 2(b).In practice, however, such picking is both tedious and time consuming.

Illumination zones and spatial resolution
A large illumination zone does not necessarily mean a good spatial resolution limit dx.Here, dx is the minimum width of a slowness box whose velocity can be resolved from the data.For virtual data, dx depends on the midpoint position and minimum offset between a reciprocal pair of geophones.As an example, the pair of dashed vertical lines at A and B in Fig. 2 arrival refraction times can be picked almost everywhere between the two reciprocal shots, then the spatial resolution dx is equal to the spacing between adjacent geophones.This assumes that the typical wavelength of the arrivals is much shorter than a geophone spacing.

Slowness variance
Traveltime picking is prone to random errors, and so the variation of these data errors give rise to an uncertainty in the estimated slowness model, denoted by the slowness variance σ 2 i in the ith cell.Appendix B shows that the slowness variance associated with the reciprocal data is proportional to 1 l 2 (l is the width of a slowness cell in Fig. A1, see Appendix A for more details), while that for the virtual data is proportional to 2 Nl 2 .Thus, the variance of the virtual slowness can be much less than that for the reciprocal data if N is large and there is only weak correlation in the data errors.

Consistency test
A consistency test is now described which assesses the validity of the assumption that the recorded refraction arrivals are those of head waves.
Assume head wave arrivals in Fig. 3(a) so that is the traveltime of a virtual head wave excited at B and recorded as a 1st-arrival at C. However, this might not be a head-wave traveltime if the source at A is far from the receivers and there is a veloc-ity gradient linearly increasing in depth.As an example, Fig. 3(b) shows a diving wave recorded as a 1st-arrival at B and C when there is a strong velocity gradient that linearly increases in depth.
Diving waves violate the head-wave assumption underlying eq.( 4) so the resulting traveltime τBC is that of a non-physical diving wave between B and C (see the solid black ray in Fig. 3b).Therefore, traveltimes from recorded shot gathers must be tested for consistency with the head-wave assumption.We propose a consistency test identical to that given in Dong et al. (2006) where the traveltimes in a common geophone pair gather (CPG) must equal one another if they are for head waves.As an example, the headwave traveltimes in Figs 4(a)-(c) are transformed by into the traveltimes τBC for the traces in Fig. 4(d).If the difference between these traveltimes is greater than the picking uncertainty for different positions of A on the left side of B then they should be rejected as head-wave traveltimes associated with the same refractor.Some of the head-wave traveltimes associated with the same refractor are not physically related to refraction times because the receiver-receiver offset is not at a post-critical offset (Fig. 5).Such non-physical traveltimes are identified by a sudden jump in the CPG.These non-physical traveltimes are easily identified1 and eliminated, but they leave a gap in the arrival times of the direct waves.The interpolated direct-wave traveltimes from the infill CSGs will automatically fill this gap.

N U M E R I C A L R E S U LT S
Parsimonious refraction interferometry is now tested for traveltimes generated for two synthetic models and two seismic field surveys.One of the synthetic models is a three-layered velocity model with a velocity gradient and the second one has a strong velocity gradient that increases with depth to demonstrate the diving wave problem.In each synthetic example the traveltimes were generated using a finite-difference solution to the eikonal equation (Qin et al. 1992) and inverted by the multiscale gradient method described in Nemeth et al. (1997).The field examples are for refraction data recorded near the Gulf of Aqaba and over a deep basin in East Africa.Two more data examples are presented in the Supporting Information.

Three-layer synthetic model
The first example is a three-layer model (Fig. 6a), where the first layer has a constant velocity of 400 m s −1 , the second layer has a horizontal velocity gradient of 0.2 m s −1 m −1 , while the third layer has a vertical velocity gradient of 4 m s −1 m −1 with a minimum velocity value of 2600 m s −1 .A finite-difference solution (Qin et al. 1992) to the 2-D eikonal equation is used to compute 240 shot gathers of first-arrival traveltimes (Fig. 7a), where the source and receiver points are located every 5 m on the surface.
We first inverted the first-arrival traveltime values of two shot gathers, with shots located at X = 0 m and X = 1200 m, as well as the near-offset traveltimes picked from five infill shot gathers (Fig. 7b) to get the tomogram shown in Fig. 6(b).In this case there is a poor correspondence between the actual velocity model and the tomogram constructed from the reciprocal and infill data.For Eq. ( 3) is then used to compute the virtual traveltimes from the 480 traveltimes associated with the two reciprocal shot gathers.To account for the fact that the direct-wave traveltimes could not be generated at long offsets, the direct-wave traveltimes were picked from five infill CSGs with the shots evenly located at 200 m intervals between the endpoints of the receiver line.These direct-wave traveltimes were then interpolated to other source locations to fill in the missing direct-wave traveltimes in the other CSGs.The result is the creation of O(57 000) virtual traveltimes computed for virtual shots at each of the geophones (Fig. 7c).These virtual traveltimes are inverted to give the tomogram shown in Fig. 6(d).This tomogram agrees with the standard tomogram, which is not surprising because (1) the virtual traveltimes are nearly the same as the true traveltimes shown in Figs 7(d) and 8(a), (2) there are no strong velocity gradients in the model and (3) there is a dense illumination in the subsurface from the virtual rays.
The validity of the head-wave assumption could has been predicted by determining whether or not the traveltimes in the firstarrival traveltimes in the CPGs are flat.Fig. 9(a) shows that the CPG traveltimes mostly agree with one another for the same CPG.For specified locations of B and C (see ray diagram in upper left of Fig. 9a), the thin black line represents the traveltimes τ BC for all actual shot locations to the left of B. Here, the location C is always at location #240 and B is progressively moved to its left to give rise to larger traveltimes τ BC .

Linear velocity gradient model
Eq. ( 3) is not valid if the recorded traveltimes are mainly associated with diving waves and not with head waves.To demonstrate this point, we used the linear-velocity-gradient model of Fig. 10(a) to calculate the first-arrival traveltimes for 126 shot gathers with a shot interval of 8 m.The calculated traveltimes are shown in Fig. 10(b).Only the first-arrival traveltimes from the two reciprocal shots at x = 0 and x = 2 km (Fig. 10c) are converted by eq. ( 3) to virtual traveltimes for virtual shots at every grid point along the surface (Fig. 10d).The differences between the actual and virtual traveltimes range between 0 and −1.4 s as shown in Figs 10(e) and 8(b).This is consistent with the disagreement between the traveltime matrices in Figs 10(b) and (d).There are two reasons for this disagreement, (1) the diving waves violate the head-wave assumption and (2) in the case of diving waves, sources are not at

Aqaba field data
A seismic survey was carried out near the Gulf of Aqaba (Hanafy et al. 2014), where 120 vertical-component geophones were deployed at 2.5 m intervals along a line.A 90-kg accelerated-weight drop was used for a source at every geophone position to record 120 common shot gathers.A typical shot gather is shown in Fig. 11(a), and more than 14 000 first-arrival traveltimes were picked from the 14 400 traces in the CSGs (Fig. 12c).The signal-to-noise ratio (SNR) is excellent so traveltimes could be picked at almost every trace, as verified by the densely populated traveltime matrix in Fig. 12(c).
A source at each end of the line was used to generate two reciprocal shot gathers.The first-arrival traveltimes were picked and used to compute more than 14 000 virtual traveltimes from eq. (3).As in the synthetic example, some of these traveltimes are not physically related to refraction times because the receiver-receiver offset is not at a post-critical offset.These non-physical traveltimes are easily identified and eliminated, but they leave a gap in the arrival times of the direct waves.To fill in this gap, we included the directwave traveltimes from five infill CSGs with a shot interval of 50 m (Fig. 12a).Inverting the traveltimes by a multiscale tomography method (Nemeth et al. 1997) gives the tomograms shown in Figs 12(b),  (d) and (f).It is clear that there is a mostly excellent agreement between the virtual and standard tomograms shown in Figs 12(d) and (f), respectively.In the standard tomogram, the shooting effort required 1.5 d, and more than a day of effort was needed for picking (14 400) traveltimes.For the virtual tomogram, it is estimated that we would need less than 1 hr of shooting time and less than 30 min of picking time for the reciprocal and infill traveltimes.Moreover, the five shot gathers recorded for the direct-wave traveltimes would not require much more effort to record, because only one strike of the weight-drop source is needed to record short-offset traces with adequate SNR.The subsequent generation of the O(14 000) virtual traveltimes is automatically computed by the computer in less than a second.Fig. 13(a) shows that the average differences between the virtual and actual traveltimes are well below the estimated picking error of 4 ms in the data.
The consistency test results for the recorded Aqaba data are depicted in Fig. 14(a).For specified locations of B and C (see ray diagram in upper right), the thin black line represents the traveltimes τ BC for all actual shot locations to the left of B. Here, C is always at location #120 and B is progressively moved to its left to give rise to larger traveltimes τ BC .
Most of the recorded traveltimes τ BC to the left are flat to within a quarter of the dominant period ( 1 4 T, where T = 32 ms), and therefore indicate passage of the consistency test.However, the traveltimes associated with the bumps are rejected.This plot is for one location of C, and so the traveltimes for other locations of C should also be tested for consistency with the head-wave assumption.
It is tempting to apply the consistency test to the virtual traveltimes τij , not just the recorded traveltimes τ ij .However, this will not be useful as shown in the following exercise.Assume a source at s to the left of A (Fig. 15).The receiver locations are at B and C. Recorded traveltimes, from two reciprocal shots at s and D, are used to compute the virtual traveltimes τ ij for all interior locations associated with the i and j indices.These virtual traveltimes are used to generate the CPG traveltime τBC in eq. ( 5) (see Fig. 15a).Plugging eq. ( 4) into eq.( 5) we get which is true for any post-critical location s to the left of A (Fig. 15b).
Even if the recorded reciprocal traveltimes are those for diving waves, the virtual CPG time τBC will appear to be a head-wave time because it is the same for any virtual source location s to the left of A. Thus, the consistency test cannot be applied to virtual traveltimes to distinguish head-wave arrivals from diving waves.

East Africa field data
Parsimonious refraction interferometry is now applied to a data set consisting of strong diving waves.In this case, the survey is conducted over a deep basin where the two reciprocal shots at the end of the line are separated by 690 m, and there are 139 verticalcomponent geophones evenly distributed along the line at 5 m intervals.The standard survey is conducted with the accelerated weight drop source located at 10 m intervals along the line.A typical shot gather is shown in Fig. 11(b).
The reciprocal traveltimes are picked and used to generate the virtual traveltimes.The comparison between the reciprocal, virtual and standard traveltimes is depicted along the left column of Fig. 16.Unlike the Aqaba data, there are notable discrepancies between the virtual and standard traveltimes (Fig. 13b).These discrepancies are manifested in the notable dissimilarities between the standard and virtual tomograms in the right column.The reason for these discrepancies is revealed by the consistency test results illustrated in Fig. 14(b), which shows strong violations of the head-wave assumption.

C O N C L U S I O N S
The theory of parsimonious refraction interferometry and tomography is presented where a dense set of virtual refraction traveltimes is computed from refraction traveltimes picked from a pair of reciprocal and several infill shot gathers.In theory, a virtual shot gather of refraction traveltimes can be computed for a virtual shot at each of the N geophones in the 2-D reciprocal survey.This means that O(N 2 ) virtual traveltimes can be created from virtual shots placed at each of the N geophones.An advantage of virtual data is that they can reduce the condition number of [L T L] by a factor of three compared to that of the reciprocal data.The accompanying increase in ray coverage reduces the number of slowness models that satisfy the traveltime equations.In addition the slowness variance associated with virtual data can be less than that obtained from reciprocal data.
Tests with synthetic data and field data verify that virtual refraction traveltimes can be inverted to give tomograms that closely resemble those computed from traveltimes recorded in a dense survey.If the first arrivals are mostly head-wave arrivals, not strong diving waves, then dense 2-D refraction surveys might be replaced by inexpensive reciprocal surveys with as few as two shots placed at each end of the line and several infill shot gathers.
Some opportunities presented by parsimonious refraction interferometry include the following.
(i) The virtual traveltimes can be used to design narrow picking windows for an autopicker to pick first-arrival traveltimes from the dense survey data.For example, the Aqaba data only used 240 traveltime picks to automatically generate O(14 000) virtual traveltimes.
(ii) More efficient surveys can be carried out because a reciprocal survey that includes several infill shot gathers is much less time consuming than a standard survey.
(iii) Traveltimes from refractions can be used to estimate statics (Yilmaz 2001).The problem is that land surveys have coarse spacing between recording stations and so prevent an estimate of the shallow velocity model.Can an inexpensive source truck towing a shortoffset land streamer be used to quickly obtain short-offset reciprocal gathers?If so, the virtual traveltimes can be inverted to give a detailed estimate of the shallow velocity distribution.A source is required at both the front and back of the land streamer.
(iv) Refraction surveys in environmentally sensitive sites, such as archaeological surveys, might use a parsimonious survey to minimize the number of sources.Some limitations of parsimonious refraction interferometry are the following.
(i) The virtual traveltimes can have up to three times the variance of the recorded traveltimes for random uncorrelated picking errors.
(ii) The first arrivals are assumed to be mostly head waves, which is strictly not true for velocities that strongly increase with depth.Surprisingly, our tests suggest that parsimonious tomography has some tolerance to traveltimes associated with diving waves.We do not always expect this type of tolerance for wide-offset refractions from the deep part of a basin or crust.
(iii) First-arrivals with low SNRs will contain large picking errors at far-offset traces.These large picking errors are expected to propagate to all source-receiver offsets in the virtual data.
(iv) If the reciprocal pair of shots is too far apart then this will lead to poor lateral resolution and gaps in the illumination of the subsurface.This problem can be partly mitigated by picking later direct and refraction arrivals, or by conducting additional reciprocal surveys with shorter offsets across the survey area.
(v) Out-of-the-plane refractions in the reciprocal survey will degrade the accuracy of the virtual refractions.
(vi) The recorded first-arrival traveltimes must satisfy the consistency assumption, which is difficult if only two reciprocal shot gathers are recorded.One possibility is to examine these two reciprocal gathers for evidence of head waves, where first arrivals have the same slope for a large range of source-receiver offsets.Diving waves might be recognized by the continuous decrease of the slope dt/dx with increasing offset.Alternatively, one can record some extra shots in the interior over a portion of the line and test the traveltimes from these shot gathers for head-wave consistency.Please note: Oxford University Press is not responsible for the content or functionality of any supporting materials supplied by the authors.Any queries (other than missing material) should be directed to the corresponding author for the paper.

A P P E N D I X A : C O N D I T I O N N U M B E R O F [ L T L]
We will derive the formulas for the condition number κ of the [L T L] associated with the reciprocal and virtual-traveltime equations.The traveltimes are for the refractions in the Fig. A1 slowness model.By definition, the 1-norm condition number is

A1 κ recip for Reciprocal Data
Assume the two-layer model in Fig. A1, where the second layer has lateral heterogeneities with an unknown slowness s i in each of the square boxes.The slowness of the top layer is homogeneous and is assumed to be known, the thickness of this top layer is also known, and each box has a width of l.The goal is to estimate the condition number of [L T L] associated with the reciprocal refraction traveltimes where a source is located at each end of the model.

Parsimonious refraction interferometry 707
The traveltime equations for a shot at position 0 in Fig. A1(a) and one at X = 5 can be represented by Ls = t, where l 0 0 0 0 l l 0 0 0 l l l 0 0 l l l l 0 l l l l l 0 0 0 0 l 0 0 0 l l 0 0 l l l where t i (t i ) is the traveltime for a rightgoing (leftgoing) refraction along a horizontal ray path.The traveltimes t i and ti are considered to be the recorded traveltimes after subtracting the upgoing and downgoing traveltimes in the top layer, which are assumed to be known.
The L matrix can be decomposed into where ) where the inverse to Q is the backward difference matrix (Schuster 1988) and its inverse is In general, the normal matrix [L T L] for 2N reciprocal equations is the Toeplitz matrix ) which is defined as Parsimonious refraction interferometry 709 Therefore, for even values of N, the 1-norm of the normal matrix for reciprocal data is The inverse matrix for eq.( A8) is where α = 0.5N +2 N +3 and β = − α N 3 + 0.5(N +1)

3
for N ≥ 3. The 1-norm of this inverse is Substituting eqs (A10) and (A12) into eq.(A1) gives the condition number for the reciprocal traveltimes equations: ) so that κ(L T L) recip grows quadratically with the number of unknowns N.

A2 κ virt for virtual Data
We will now derive the condition number κ virt = N 2 /2 + N for the L T L from virtual data, which is about three times smaller than that for reciprocal data.This suggests that there will be many fewer slowness models that satisfy the noisy virtual data compared to the reciprocal data.
As discussed previously, virtual refraction traveltimes can be created from the reciprocal traveltimes.Thus, the virtual-traveltime equations for the ray paths depicted in Fig. A1(b) for a virtual source at position 1 are given by where t 1 is the vector of virtual traveltimes.In this case L 1 can be decomposed into the product of the summation Q matrix in eq. ( A4) and the column-mask matrix M 1 : 0 0 0 0 0 0 1 0 0 0 0 0 1 0 0 0 0 0 1 0 The role of the column-mask matrix M i is to mask columns (1, 2, . . .i) of Q.
In general, for N unknown slowness cells, N virtual shot gathers can be created similar to those shown in Fig. A1.The picked traveltimes give rise to the overdetermined system of equations (A16) The above system of equations lead to the N × N normal equations where [L T L] is the symmetric matrix which is also classified as a bisymmetric matrix (Tao & Yasuda 2002), where one symmetry axis is along the matrix diagonal and the other is along the axis that intersects the bottom left and top right of the matrix.In this case [L T L] i j = i(N − j + 1) for i ≤ j; [L T L] i j = j(N − i + 1) for i > j. (A19) Therefore, the formula for the 1-norm of [L T L] is where the maximum column sum is at j = N/2 for N an even number.

Figure 1 .
Figure 1.Two-layer model where the grey medium is faster than the top layer; the reciprocal sources are at A and D and are associated with the dashed red ray.The dashed blue ray path is associated with the virtual refraction ray that is excited by the virtual source (blue star) at C and terminates at B. The refraction traveltimes associated with the reciprocal shots (green and red stars) can be decomposed into the virtual refraction traveltime generated by the blue star.Here, x and x are two points along the interface between the upper and lower layers.

Figure 2 .
Figure 2. The red and green zones approximate the areas of illumination associated with traveltimes picked from (a) first arrivals and (b) both first and the dashed later arrivals.The red (green) lines indicate direct-wave (refraction) arrivals.The horizontal resolution limit dx is denoted by two-sided arrows for different portions of the subsurface.The resolution limit dx cannot be less than the distance between adjacent geophones.

Figure 3 .
Figure 3. Virtual traveltimes created with (a) short-offset T BC = T B C + T BC − T B C and (b) long-offset TBC = T AC + T DB − T AC virtual traveltimes.The long-offset (short-offset) sources generate diving waves (head waves) as first arrivals so that T BC = TBC .The black solid ray corresponds to the virtual ray path associated with either T BC or TBC .

Figure 5 .
Figure 5. Ray diagram showing the dashed non-physical rays such that τ BC = τ DB + τ AC − τ D A is a non-physical virtual traveltime which does not correspond to a physical head wave.In comparison, τ BC = τ DB + τ AC − τ DA corresponds to the virtual traveltime of a physical head wave along the dotted black ray.The zone of non-physical virtual traveltimes is between the origin at B and the sudden change in the traveltime slope.

Figure 4 .
Figure 4. Virtual traveltimes τ BC created by computing (a) τA 1 C − τ A 1 B , (b) τ A 2 C − τ A 2 B and (c) τ A 3 C − τ A 3 B .If the input traveltimes are those of head waves, then the output (d) traveltime τ BC will be that of a virtual head wave recorded at C. The common geophone pair gather (CPG) of virtual traces in a-c share the same pair of geophones at B and C, but have different source locations.

Figure 6 .
Figure 6.(a) Three-layer model, (b) tomogram inverted from the two reciprocal and five infill shot gathers, (c) standard tomogram inverted from the 57 600 actual traveltimes in 240 shot gathers, with a shot at each geophone location.(d) Tomogram inverted from O(57 000) virtual traveltimes created from the first-arrival traveltimes picked from two reciprocal and five infill shot gathers.

Figure 7 .
Figure 7. (a) The actual traveltime picks, (b) traveltimes for two reciprocal and five infill shot gathers used as input to the proposed method, (c) the virtual traveltimes calculated from eq. (3) and (d) the difference between the actual and virtual traveltime picks.

Figure 8 .
Figure 8.The histogram of the difference between the true and virtual traveltimes of the (a) three-layer and (b) gradient synthetic models.The average differences in (a) are less than 5 ms, while in (b) it is as high as 1.4 s.

Figure 9 .
Figure 9. Consistency test results for the (a) three-layer and (b) gradient synthetic models.

Fig. 2
Fig. 2(a)  shows that a reciprocal pair of shot gathers can give rise to gaps in the illumination zones of the direct arrivals and refractions.

Figure 10 .
Figure 10.(a) Model with velocity linearly increasing in depth.(b) True, (c) reciprocal and (d) virtual traveltimes.(e) The difference between the true and virtual traveltimes.
(a) indicates the reciprocal refraction traveltimes with minimum offset dx = |x A − x B | that are pickable from the leftgoing and rightgoing refractions.This means that picking later refraction arrivals in Fig. 2(b) will lead to wider subsurface illumination, but the estimated velocity model can still suffer from poor resolution as demonstrated by the wide double-sided black arrow.If both the first-arrival and later-

Figure 11 .
Figure 11.Typical common shot gathers recorded at the (a) Aqaba and (b) East Africa field sites.

Figure 12 .
Figure 12.Gulf of Aqaba traveltime picks for the (a) reciprocal, (c) standard and (e) virtual field data.The corresponding tomograms are on the right.

Figure 13 .
Figure 13.The histogram of the differences between the true and virtual traveltimes obtained from the (a) Aqaba and (b) East Africa data.The average differences in (a) are less than 5 ms, while in (b) it is as high as 120 ms.

Figure 14 .
Figure 14.Consistency test results for the (a) Aqaba and (b) East Africa data.The East Africa data in (b) fail the consistency test while it is passed for the Aqaba data.In all of the figures, the input data were traveltimes picked from the standard survey where there is a dense distribution of shots.

Figure 15 .
Figure 15.A sketch shows the calculation of common pair gather τBC from (a) virtual ( τAC and τAB ) and (b) recorded data (τ sC and τ sB ).

Figure 16 .
Figure 16.East Africa traveltime picks for the (a) reciprocal, (c) standard and (e) virtual data.The corresponding tomograms are on the right.

Figure 1 .
Figure 1.(a) Aqaba velocity model, (b) tomogram inverted from the two reciprocal and some additional direct-wave traveltimes, (c) standard tomogram inverted from the 14 400 actual traveltimes in 240 shot gathers, with a shot at each geophone location, and (d) tomogram inverted from virtual and interpolated direct-wave traveltimes.The Aqaba model is based on results from a seismic survey near the Gulf of Aqaba.Figure 2. (a) The actual traveltimes, (b) two reciprocal and five infill shot gathers used as input to the proposed method, (c) the virtual traveltimes calculated from eq. (3) and (d) the difference between the actual and the virtual traveltimes.Figure 3. (a) The traveltime differences and (b) the corresponding histogram of the synthetic model.Figure 4. The CPG of the (a) synthetic model and (b) Wadi Qadid field example.Figure 5.A sample CSG from Wadi Qadid field data.Figure 6. Wadi Qadid traveltime picks for the (a) reciprocal, (c) standard and (e) virtual traveltimes.The corresponding tomograms are on the right-hand column of figures.Figure 7. (a) The traveltime differences and (b) the corresponding histogram of Wadi Qadid example.

Figure 2 .
Figure 1.(a) Aqaba velocity model, (b) tomogram inverted from the two reciprocal and some additional direct-wave traveltimes, (c) standard tomogram inverted from the 14 400 actual traveltimes in 240 shot gathers, with a shot at each geophone location, and (d) tomogram inverted from virtual and interpolated direct-wave traveltimes.The Aqaba model is based on results from a seismic survey near the Gulf of Aqaba.Figure 2. (a) The actual traveltimes, (b) two reciprocal and five infill shot gathers used as input to the proposed method, (c) the virtual traveltimes calculated from eq. (3) and (d) the difference between the actual and the virtual traveltimes.Figure 3. (a) The traveltime differences and (b) the corresponding histogram of the synthetic model.Figure 4. The CPG of the (a) synthetic model and (b) Wadi Qadid field example.Figure 5.A sample CSG from Wadi Qadid field data.Figure 6. Wadi Qadid traveltime picks for the (a) reciprocal, (c) standard and (e) virtual traveltimes.The corresponding tomograms are on the right-hand column of figures.Figure 7. (a) The traveltime differences and (b) the corresponding histogram of Wadi Qadid example.

Figure 3 .
Figure 1.(a) Aqaba velocity model, (b) tomogram inverted from the two reciprocal and some additional direct-wave traveltimes, (c) standard tomogram inverted from the 14 400 actual traveltimes in 240 shot gathers, with a shot at each geophone location, and (d) tomogram inverted from virtual and interpolated direct-wave traveltimes.The Aqaba model is based on results from a seismic survey near the Gulf of Aqaba.Figure 2. (a) The actual traveltimes, (b) two reciprocal and five infill shot gathers used as input to the proposed method, (c) the virtual traveltimes calculated from eq. (3) and (d) the difference between the actual and the virtual traveltimes.Figure 3. (a) The traveltime differences and (b) the corresponding histogram of the synthetic model.Figure 4. The CPG of the (a) synthetic model and (b) Wadi Qadid field example.Figure 5.A sample CSG from Wadi Qadid field data.Figure 6. Wadi Qadid traveltime picks for the (a) reciprocal, (c) standard and (e) virtual traveltimes.The corresponding tomograms are on the right-hand column of figures.Figure 7. (a) The traveltime differences and (b) the corresponding histogram of Wadi Qadid example.

Figure 4 .
Figure 1.(a) Aqaba velocity model, (b) tomogram inverted from the two reciprocal and some additional direct-wave traveltimes, (c) standard tomogram inverted from the 14 400 actual traveltimes in 240 shot gathers, with a shot at each geophone location, and (d) tomogram inverted from virtual and interpolated direct-wave traveltimes.The Aqaba model is based on results from a seismic survey near the Gulf of Aqaba.Figure 2. (a) The actual traveltimes, (b) two reciprocal and five infill shot gathers used as input to the proposed method, (c) the virtual traveltimes calculated from eq. (3) and (d) the difference between the actual and the virtual traveltimes.Figure 3. (a) The traveltime differences and (b) the corresponding histogram of the synthetic model.Figure 4. The CPG of the (a) synthetic model and (b) Wadi Qadid field example.Figure 5.A sample CSG from Wadi Qadid field data.Figure 6. Wadi Qadid traveltime picks for the (a) reciprocal, (c) standard and (e) virtual traveltimes.The corresponding tomograms are on the right-hand column of figures.Figure 7. (a) The traveltime differences and (b) the corresponding histogram of Wadi Qadid example.

Figure 5 .
Figure 1.(a) Aqaba velocity model, (b) tomogram inverted from the two reciprocal and some additional direct-wave traveltimes, (c) standard tomogram inverted from the 14 400 actual traveltimes in 240 shot gathers, with a shot at each geophone location, and (d) tomogram inverted from virtual and interpolated direct-wave traveltimes.The Aqaba model is based on results from a seismic survey near the Gulf of Aqaba.Figure 2. (a) The actual traveltimes, (b) two reciprocal and five infill shot gathers used as input to the proposed method, (c) the virtual traveltimes calculated from eq. (3) and (d) the difference between the actual and the virtual traveltimes.Figure 3. (a) The traveltime differences and (b) the corresponding histogram of the synthetic model.Figure 4. The CPG of the (a) synthetic model and (b) Wadi Qadid field example.Figure 5.A sample CSG from Wadi Qadid field data.Figure 6. Wadi Qadid traveltime picks for the (a) reciprocal, (c) standard and (e) virtual traveltimes.The corresponding tomograms are on the right-hand column of figures.Figure 7. (a) The traveltime differences and (b) the corresponding histogram of Wadi Qadid example.

Figure 6 .
Figure 1.(a) Aqaba velocity model, (b) tomogram inverted from the two reciprocal and some additional direct-wave traveltimes, (c) standard tomogram inverted from the 14 400 actual traveltimes in 240 shot gathers, with a shot at each geophone location, and (d) tomogram inverted from virtual and interpolated direct-wave traveltimes.The Aqaba model is based on results from a seismic survey near the Gulf of Aqaba.Figure 2. (a) The actual traveltimes, (b) two reciprocal and five infill shot gathers used as input to the proposed method, (c) the virtual traveltimes calculated from eq. (3) and (d) the difference between the actual and the virtual traveltimes.Figure 3. (a) The traveltime differences and (b) the corresponding histogram of the synthetic model.Figure 4. The CPG of the (a) synthetic model and (b) Wadi Qadid field example.Figure 5.A sample CSG from Wadi Qadid field data.Figure 6. Wadi Qadid traveltime picks for the (a) reciprocal, (c) standard and (e) virtual traveltimes.The corresponding tomograms are on the right-hand column of figures.Figure 7. (a) The traveltime differences and (b) the corresponding histogram of Wadi Qadid example.

Figure 7 .
Figure 1.(a) Aqaba velocity model, (b) tomogram inverted from the two reciprocal and some additional direct-wave traveltimes, (c) standard tomogram inverted from the 14 400 actual traveltimes in 240 shot gathers, with a shot at each geophone location, and (d) tomogram inverted from virtual and interpolated direct-wave traveltimes.The Aqaba model is based on results from a seismic survey near the Gulf of Aqaba.Figure 2. (a) The actual traveltimes, (b) two reciprocal and five infill shot gathers used as input to the proposed method, (c) the virtual traveltimes calculated from eq. (3) and (d) the difference between the actual and the virtual traveltimes.Figure 3. (a) The traveltime differences and (b) the corresponding histogram of the synthetic model.Figure 4. The CPG of the (a) synthetic model and (b) Wadi Qadid field example.Figure 5.A sample CSG from Wadi Qadid field data.Figure 6. Wadi Qadid traveltime picks for the (a) reciprocal, (c) standard and (e) virtual traveltimes.The corresponding tomograms are on the right-hand column of figures.Figure 7. (a) The traveltime differences and (b) the corresponding histogram of Wadi Qadid example.

)Figure A1 .
Figure A1.Ray paths for a common shot gather with the source at positions (a) 0, (b) 1 and (c) 2. Here, it is assumed that the source at position 0 is an actual source and the sources at 1 and 2 are virtual sources generated from a reciprocal pair of shot gathers.

Figure A2 .
Figure A2.Condition number versus N computed from the 1-norm product ||[L T L]|| • ||[L T L] −1 ||, where [L T L] is the N × N normal matrix numerically computed by a MATLAB code.The red (green) stars denote the condition numbers numerically computed for the reciprocal virtual data associated with the checkerboard model in Fig. A1.The black lines, computed by eqs (A24) and (A13), are in perfect agreement with the numerically computed condition numbers.