Research ArticleGEOPHYSICS

Far-field super-resolution imaging of resonant multiples

See allHide authors and affiliations

Science Advances  20 May 2016:
Vol. 2, no. 5, e1501439
DOI: 10.1126/sciadv.1501439


We demonstrate for the first time that seismic resonant multiples, usually considered as noise, can be used for super-resolution imaging in the far-field region of sources and receivers. Tests with both synthetic data and field data show that resonant multiples can image reflector boundaries with resolutions more than twice the classical resolution limit. Resolution increases with the order of the resonant multiples. This procedure has important applications in earthquake and exploration seismology, radar, sonar, LIDAR (light detection and ranging), and ultrasound imaging, where the multiples can be used to make high-resolution images.

  • seismic resonant multiples
  • super-resolution imaging


In 1873, Abbe recognized a fundamental diffraction limit of optics: whenever an object is imaged by a classical optical lens, features of the object smaller than λ/2 become unrecoverable. Such fine, detailed information is lost because light emerges with these fine features as evanescent energy, which decays exponentially away from the object and is not carried by the propagating waves. Evanescent energy can be used for imaging to achieve a resolution that is higher than the Abbe limit (1, 2). As an example, Lerosey et al. (3) demonstrated that subwavelength (that is, super-resolution) imaging of the source location could be achieved by refocusing scattered microwave energy from near-field scatterers. They denoted this refocusing operation as a time-reversal mirror (4, 5), which is similar to reverse time migration by exploration geophysicists (68) except that no velocity model is required because the Green’s functions for refocusing are recorded in the data.

Super-resolution imaging is important for many applications in optical imaging, including nanolithography and imaging of live biological specimens (911). For nanolithography, subwavelength focusing is used with ultraviolet light to etch nanoscale transistors, but such light cannot be used with live biological specimens without harming the subject. Recent developments show that super-resolution can be achieved in the far field of both sources and receivers. For example, Huang and Zheludev (12) demonstrated that the construction of a mask at the aperture of an optical or plasmonic device will create super-oscillatory constructive interference in the object plane. The super-oscillatory characteristic means that bandlimited functions can oscillate arbitrarily faster than the highest Fourier components they contain (13). Numerical simulations show that optical super-resolution imaging can be achieved with a resolution higher than 1/10 λ in the far-field region. However, the imaged spot is characterized by a significantly diminished amplitude and is surrounded by severe sidebands. The energy of the focused spot is often several orders of magnitude less than that of the side lobes, making it difficult to apply this technology to seismic waves.

Analogs of practical super-resolution imaging in seismology are currently limited to a few applications where the objects of interest are in the near field of the sources or receivers (14). However, most of the scattering bodies of interest in seismology are in the far field of both the sources and receivers; hence, near-field scattering methods are not applicable. We now show both theoretically and experimentally that subwavelength imaging of reflectors with arbitrary dip angles is possible with resonant multiples. This largely overcomes the limitation of resonant super-resolution to mostly horizontal layers (15).

Resolution limits of resonant multiples

A resonant multiple is defined as a multiple where the upgoing and downgoing ray paths coincide (see Fig. 1). For example, the offset between the source and the receiver in Fig. 1A is zero, such that the nth-order free-surface multiple has a travel time described byEmbedded Image(1)where v is the velocity of the topmost layer and N = 0 for primary reflection. The vertical resolution in locating the depth of the interface can be estimated by first taking the derivative of Eq. 1 with respect to z to getEmbedded Image(2)

Fig. 1 Ray diagrams for resonant multiples.

(A and B) Ray diagrams for a flat layer first-order resonant multiple (A) and for a tilted layer third-order resonant multiple (B). The topmost layer represents the free-surface boundary. S represents the position of the coincident source and receiver, which is just below the free surface. x marks the point that is revisited by the multiples. The rays in the tilted layer are described by Embedded Image, and its mirrored ray path Embedded Image has the same path length Embedded Image. (B) is a third-order multiple because the ray reflects on the free surface three times: twice at A1 and once at A2. The ray for a primary reflection is described by Embedded Image. (C) Pure resonant multiples where the rays at the source point S and at A2 are perpendicular to the free surface. Unlike the illustration in (B), higher-order rays in (C) will completely coincide with Embedded Image.

Defining the dominant period by T0, substituting δTT0/2 and λ = vT0 into Eq. 1, and rearranging give the vertical resolution limit δz → Δz (16)Embedded Image(3)

The limit states that Δz is the minimum separation in depth between two reflectors that are distinguishable from one another in the data. The greater the order N of the resonance, the higher the resolution in depth. If the bed is tilted with respect to the horizontal recording surface, then this resolution limit is calculated along the direction perpendicular to the interface in Fig. 1B.

The resolution limit for seismically imaging a dipping reflector is now derived for a resonant multiple with arbitrary order N. An example of a third-order (N = 3) resonant multiple is depicted in Fig. 1B, where the resonant multiple makes N bounces on the free surface, for instance, at A1 and A2. Let the distance from S perpendicular to the reflector be Embedded Image. Then, we have Embedded Image and Embedded Image.

The total two-way travel time T is given byEmbedded Image(4)

Here, v is the homogeneous velocity in the upper layer and Embedded Image allows Eq. 4 to subsume the case of a primary reflection (N = 0), where T = 2d/v. Perturbing the travel time δT with respect to the perturbation in δd givesEmbedded Image(5)

Setting vδT → λ/2 gives the nth-order resolution limit along the direction perpendicular to the interfaceEmbedded Image(6)where Embedded Image is the resolution coefficient. Setting α = 0 and N = 1 gives the vertical resolution limit Δd → Δz = λ/8 for a horizontal reflector, which agrees with Eq. 3 for N = 1. Equation 6 also reveals that the resolution limit becomes higher with decreasing interface dip angle.

The resolution limit in Eq. 6 is that for bulk resolution where the entire layer interface in Fig. 1B is shifted by δd along the perpendicular to the interface. For a localized perturbation at B1 where only the point B1 is perturbed by the amount δd, the depth resolution limit is inversely proportional to the number of bounces at B1 that is, Δd = 0.25λ/(M cos γ), where γ is the angle of incidence at B1 and M is the number of bounces at B1. In Fig. 1B, for example, M = 2, even though the ray diagram is that for a third-order multiple. If the rays at the source point S and at A2 are normal to the free surface as in Fig. 1C, then this forms a pure resonance system. Unlike the illustration in Fig. 1B, the higher-order pure resonance rays in Fig. 1C will completely coincide with Embedded Image and M will grow with the order N of the multiple.

Imaging resonant free-surface multiples

Relocating a reflection event in a trace to its reflecting point on the interface is known in the geophysics community as migration. For example, if the zero-offset trace d(g, t) only contains a diffraction from a point scatterer at xo, then it can be represented byEmbedded Image,(7)where w(t) is the impulsive-like source wavelet (with dominant period T0) associated with the point source and geophone at g and Embedded Image is the travel time for waves to propagate from g to xo. Migration answers the following question: Given the trace d(g, t), what are the possible locations of point scatterers that can explain an event that arrives at g with travel time Embedded Image? For a homogeneous medium, the answer is that the scatterers could have been located anywhere along the annulus in Fig. 2A. Here, the radius of the inner annulus centered at g is equal to Embedded Image, and the outer annulus radius is Embedded Image. Thus, the thickness of the annulus is Δrprim = λ/2 = 0.5Tov, which is the difference between primary reflections that reflect off the inner and outer borders of the annulus. This answer, also known as a single-trace migration image, is rather vague because the scatterer could have been located anywhere along the thick annulus. To reduce uncertainty in the answer, the migration images from other traces at gB can be added together to obtain the stacked migration image m(x)Embedded Image,(8)where B is the set of recording positions along the horizontal recording line; Embedded Image is a preconditioning factor that compensates for geometrical spreading, obliquity, and wavelet stretch; Embedded Image is the Fourier transform of the recorded data d(g,t); and x is the trial image point.

Fig. 2

(A and B) Single-trace migration images computed from a trace containing a primary (A) and a first-order resonant multiple (B). The model is a point scatterer (yellow ball) in a homogeneous medium, and the thickness (that is, the resolution limit) of the resonant annulus (defined by the dashed yellow lines) is half of that of the primary annulus.

Similarly, the single-trace migration image for a first-order resonant event is shown in Fig. 2B. For example, the thickness of the migration annulus for a flat reflector directly below the geophone position is Δrmult = λ/4 = 0.25Toc, where Δrmult/c is the travel time difference between the first-order resonant reflections that scatter from the outer and inner parts of the resonant annulus. Hence, the spatial resolution associated with the resonance migration operator is twice that of the primary migration operator.

The general migration operator for migrating first-order resonant multiples associated with the Fig. 2B model is given byEmbedded Image(9)where Embedded Image is the Fourier transform of the source wavelet. Here, τgx and τxA are the propagation times from the trial image point x to the geophone at g and the normal incidence point at A on the free surface, respectively (17).

The workflow for implementing Eq. 9 is to (i) select any trial image point at x in Fig. 2B and compute the imaging time 2(τgx + τxA) in Eq. 9 and (ii) use Eq. 9 to migrate the resonant multiple. The above procedure is valid for heterogeneous velocity models where the travel times are computed by ray tracing. Another way of implementing Eq. 9 is to use wavepath migration (18, 19). Its workflow is similar to the previous one, except that the trial image point at x has to be on the wavepath from g to the image.


We now present the results of resonant imaging for both synthetic data and field data from two different experiments. The first field data experiment is sonar imaging of the seafloor with 120- and 200-kHz ultrasonic sources. The second experiment is a conventional exploration survey where the airgun source has a frequency content from 5 to 60 Hz.

Synthetic test

The synthetic experiment is designed to illustrate super-resolution imaging achieved by migrating resonant multiples. Eight point scatterers are placed in a homogeneous background model, and diffractions are generated by Kirchhoff modeling (20). The data are recorded on the surface and then migrated by wavepath migration to obtain the images in Fig. 3. The resonant multiples are migrated using Eq. 9 to produce the image in Fig. 3B, which has a finer resolution than the primary image in Fig. 3A. The top four point scatterers are both vertically and horizontally distinguishable from one another in the resonant multiple migration image compared to the blurred primary image.

Fig. 3

(A and B) Migration images of point scatterers by migrating primary reflections (A) and resonant multiples (B). The red dots are the positions of the point scatterers. λ represents the dominant wavelength (1 km) associated with the source wavelet.

Field data tests

The resonant multiples in two field data sets will be tested in their ability to achieve super-resolution imaging. The first data set consists of 120- and 200-kHz echo sounder data that can be used to map out the topography of the seafloor and to detect organisms in the water column such as schools of fish. The experiments were carried out in shallow waters close to the King Abdullah University of Science and Technology (KAUST) by the survey ship shown in Fig. 4A. The second data set is a low-frequency marine survey used to map out deep oil deposits in a sedimentary basin.

Fig. 4

(A) The survey ship with the echo sounder. (B and C) Data at 120 kHz (B) and 200 kHz (C) after band-pass filtering. (D) Reflectivity images of the seafloor obtained by migrating (B) and (C).

Sonar survey. Figure 4 (B and C) shows the recorded 120- and 200-kHz sonar profiles. The seafloor mapped out by the 200-kHz echo sounder is used as the ground truth for the 120-kHz echo sounder. That is, the 120-kHz image computed from the resonant multiples in the data in Fig. 4B should achieve about the same resolution as the 200-kHz data in Fig. 4C. Sonar readings were taken 10 times per second as the boat moved at a speed of approximately 1.5 m/s over coral reefs. The primary reflections and first-order resonant multiples were migrated to produce the reflectivity images in Fig. 4D. As predicted from our theory, the image obtained from the 120-kHz resonant multiples has about the same resolution as the image computed from 200-kHz primaries. Also, the resonant multiple image has almost twice the resolution of the 120-kHz primary image. This verifies the theoretical prediction that migrating first-order resonant multiples achieves almost twice the resolution that imaging the primaries does. This assumes no significant loss of high frequencies, which is true in this example.

Exploration seismic survey. Figure 5 (A and B) shows two migration images, one from primary reflections and the other from resonant multiples. The data set is collected by a geophysical survey company. The primary reflections are migrated using a conventional Kirchhoff method (20), whereas the resonant multiples are moveout-corrected and stacked for a higher signal-to-noise ratio (17) and then migrated using Eq. 9. Here, moveout correction aligns the non–zero-offset multiples with zero-offset resonant multiples. As expected, the resonant multiple image in Fig. 5B has about twice the resolution of the primary image in Fig. 5A.

Fig. 5

(A and B) Migration image computed from first-order resonant multiples (B) has almost twice the resolution of the primary migration image (A). The blue arrow indicates migration artifacts caused by nonresonant multiple events, which bounce on different reflectors. The seafloor and top of salt are indicated by red arrows.


Super-resolution imaging with resonant multiples has possible applications in many fields of science. Figure 6 depicts four possible applications. For 35-GHz radar imaging of cloud structure (, the internal multiples in Fig. 6A might provide higher resolution of cloud structure than the primary reflections. Figure 6B suggests that resonant reflections in full-waveform (21) bathymetric LIDAR (light detection and ranging) data ( could be used to get better estimates of the water-floor topography or underwater objects than possible with standard bathymetric LIDAR. This might also be true for topographic LIDAR where resonant reflections between the ground and trees can be used to better characterize the underside features of the tree canopy ( In an urban environment, it might be possible to image the underside of a bridge with resonant multiples. If successful, this opens the possibility of achieving both ground-based and topography LIDAR using just the aerial LIDAR unit. If a mirror (see the pink mirror in Fig. 6B) is attached in front of a ground-based LIDAR, then at least twice the resolution can be achieved with full-waveform LIDAR imaging. This assumes that the topography of the land surface is known. In the ultrasound imaging of the fetus in Fig. 6C, the resonant multiples should provide higher resolution than the primaries if they can be separated from one another. For nondestructive testing and engineering, the surface features of an object can be imaged with super-resolution, as depicted in Fig. 6D.

Fig. 6

(A to D) Resonance rays for imaging using a 35-GHz radar (A), LIDAR (B), ultrasound (C), and nondestructive testing (D).

Increasing the resolution of reflectivity images by migrating multiples is similar to the principle of a Fabry-Perot interferometer, where a light beam is allowed to bounce multiple times between two mirrors in an arm of the interferometer before it interferes with the beam from the other arm. If the arms are slightly different in length, then multiple bounces amplify the slight phase difference between the two different beams, as markedly illustrated with the recent detection of gravitational waves by a Fabry-Perot interferometer (22). Increasing the number of resonant bounces in a medium increases the sensitivity of the data to small changes in the medium.


Sonar survey

Acoustic measurements were made using a Simrad EK60 echo sounder (7° beam width) operating at two frequencies (120 and 200 kHz). The pulse duration was 64 μs. Raw data (as shown in fig. S1, A and B) were demeaned and band-pass–filtered such that the peak frequency of each trace in the 200-kHz data was approximately twice that of the 120-kHz data (peak frequencies are shown in fig. S1C). After filtering, the primary reflections and resonant multiples were separated on the basis of their arrival time and normalized individually (as shown in Fig. 4, B and C). This normalization was for the enhancement of the multiple signal. Last, the primary reflections and resonant multiples were imaged using a conventional Kirchhoff method (20) and based on Eq. 9 to produce the images of the seafloor in Fig. 4D.

Exploration seismic survey

The marine data set was collected by a geophysical survey company. There are 334 shots with a shot interval of about 80 m. Each shot has 1126 hydrophones as receivers. The primaries and multiples were separated on the basis of their different arrival times. Then, moveout corrections were applied to primaries and multiples, respectively, to align the non–zero-offset data with the zero-offset data (as shown in fig. S2). After the alignment, stacking across different offsets yields a higher signal-to-noise ratio, zero-offset primary reflections, and resonant multiples (17, 20). Then, the data were filtered by a low-pass filter with a high cutoff frequency of 22 Hz (as shown in fig. S3). Finally, the stacked primaries and resonant multiples were imaged separately to obtain images shown in Fig. 5 (A and B).


We have presented a theory for transforming the coherent noise of resonant multiples into signals that can be used for super-resolution imaging. The theory is validated with both synthetic data and field data tests. The limitations to this method are that the resonant multiples must be separated from the primaries, and there must be no significant attenuation in the model; otherwise, the multiples will lose high frequencies compared to the primaries. This procedure has implications for exploration and earthquake seismology, especially if resolution can be more than doubled compared to previous methods. We also believe that it has important applications in radar, sonar, LIDAR, nondestructive testing, and ultrasound imaging, where the noise from multiples can be transformed into high-resolution images.


Supplementary material for this article is available at

fig. S1. Raw sonar data and their peak frequencies after processing.

fig. S2. Moveout corrections in the common midpoint gathers align the non–zero-offset multiples (covered by yellow shadows) with the zero-offset resonant multiples.

fig. S3. Near-offset marine data before and after moveout correction and stacking.

This is an open-access article distributed under the terms of the Creative Commons Attribution-NonCommercial license, which permits use, distribution, and reproduction in any medium, so long as the resultant use is not for commercial advantage and provided the original work is properly cited.


Acknowledgments: Funding: This work was supported by funding from KAUST. Author contributions: G.S., B.G., and Y.H. conceived the theory. A.R. and B.G. collected the sonar data. B.G. implemented the synthetic and field data results. G.S., B.G., and Y.H. wrote the paper. Competing interests: The authors declare that they have no competing interests. Data and materials availability: All data needed to evaluate the conclusions in the paper are available at or are present in the paper and/or the Supplementary Materials. Additional data related to this paper may be requested from the authors.
View Abstract

Stay Connected to Science Advances

Navigate This Article