An introduction to seismic diffraction

Despite its unique properties the diffracted seismic wavefield is still rarely exploited in common practice. Although the first works on seismic diffraction date back at least as far as the 1950s, a first rigorous theoretical framework for diffraction imaging only evolved decades later and many important questions still remain unanswered until the present day. While this comparably slow progression can partly be explained by the lack of densely sampled high quality recordings, recent advances in acquisition and dedicated processing suggest we might be at the door step to a paradigm shift in which seismic diffraction could play an important role. Despite the fact that most major progress – in terms of data acquisition and processing – has been achieved for the reflected wavefield, upon closer inspection it becomes obvious that the concept of diffraction is deeply ingrained in migration-type seismic imaging. With the aim of complementing existing contributions on the topic, this chapter is an attempt to provide an intuitive introduction to the process of seismic diffraction. Discussed are the deep conceptual roots in optics, physical links to the Kirchhoff integral as well as diffraction types and their importance in different contexts of application. By means of controlled synthetic and academic as well as industry-scale field data examples, I suggest a simple integrated framework for non-invasive diffraction separation and high-resolution imaging, which remains computationally affordable and can be reproduced by the reader. Different applications suggest that the faint diffracted background wavefield is surprisingly rich and, once it is given a voice, it announces highly resolved features such as faults, fractures, and erosional unconformities, which remain notoriously hard to image conventionally. Extending the dominant theme of high-resolution seismic imaging, I illustrate how the superior illumination due to the uniform radiation of diffraction carries the additional potential for drastically reduced acquisitions and discuss the possibility of a systematic extraction of inter-scatterer traveltimes from coda waves. Author’s post-print version, please refer to the original publication: B. Schwarz (2019), An introduction to seismic diffraction. Advances in Geophysics 60, Chapter 1, pp. 1-64, Elsevier. DOI: 10.1016/bs.agph.2019.05.001 Distributed under a Attribution-NonCommercial-NoDerivatives 4.0 International License Benjamin Schwarz 2 Seismic diffraction


INTRODUCTION
The first accounts on seismic diffraction date at least back to the 1950s, where the concept was recognized to form a very useful and important ingredient in initial attempts at migration (Hagedoorn, 1954). In addition, it was first identified as a primary carrier of information when faults and fractures are investigated (Krey, 1952). In particular in the 1970s and early 1980s, seismic diffracton phenomena were studied quantitatively by computationally solving the forward problem of modelling seismic waves in complex media (Hilterman, 1970(Hilterman, , 1975Trorey, 1977;Berryhill, 1977;Buchen and Haddon, 1980). Kunz (1960) was arguably the first to emphasize the potential value of diffraction for dedicated seismic fault interpretation. When the field of optics is concerned, however, it becomes apparent that a rigorous theoretical framework for the treatment of diffraction phenomena is even of significantly older origin and can be traced back to the pioneering works of Grimaldi, Huygens, Babinet, Fresnel and Kirchhoff (a concise but compelling account can be found in Born and Wolf, 2013). Aptly named, the process of Kirchhoff migration -be it in the depth or time domain -represents a computationally feasible adaptation of integral representations first formulated by Kirchhoff for the treatment of electromagnetic wave phenomena to seismic imaging (Schneider, 1978;Wiggins, 1984).
Kirchhoff's construction can be viewed as a mathematical manifestation of the Huygens-Fresnel principle, which simply states that any imaginable wavefield and its interaction with changes in medium properties can be described by the interference of infinitely many point sources excited with every progression of the wavefront. The drastic implications become immediately obvious when Young's celebrated slit experiments are considered, in which striking and initially unexpected interference patterns result from wave energy entering geometrical shadows of obstacles, which is not predicted by the older framework of geometrical optics. Interestingly, roughly at the same time as when the importance of seismic diffraction was first investigated, a powerful extension of geometrical optics was suggested, that can likewise accurately describe observable diffraction phenomena (Keller, 1962). Largely based on this far-field approximation, a first theoretical framework for the treatment of diffraction in seismic exploration was formulated (Klem-Musatov andAizenberg, 1984, 1985) culminating in a standard text book, which until the present day must be considered a reference on the topic (Klem-Musatov, 1994).
Triggered by these important theoretical developments, the first impactful applications of seismic imaging with diffractions appeared in the literature, which impressively suggested that the diffracted wavefield bears the potential to facilitate high-resolution imaging of important small-scale geological features (e.g. Landa et al., 1987;Kanasewich and Phadke, 1988). Despite these initial successes, however, the following decade largely saw a stagnation of developments, owing in large part to spectacular developments occurring in the fields of migration and first successful attempts at full-waveform inversion, in which diffraction is principally honored, but arguably does not take a prominent role (e.g. Etgen et al., 2009;Virieux and Operto, 2009). In addition, the introduction of insightful quantitative attributes in seismic interpretation, which, for example, enabled the detection of faults through image processing in the image domain, seemed to make the sometimes tedious investigations into the weak diffracted wavefield superfluous (e.g. Bahorich and Farmer, 1995;Marfurt et al., 1998). It was not until the late 1990s that it was discovered that in particular in the context of seismic monitoring, relevant small-scale structural changes can reliably be detected with dedicated diffraction processing (Landa and Keydar, 1998), which again sparked interesting Seismic diffraction new developments in the field with continued stimulation until the present day.
Without a doubt one of the major stumbling blocks for the success of diffraction imaging is the extreme weakness of the wavefield when compared with reflections and diving waves. Likewise severely complicating quantitative investigations, diffractions are known to heavily interfere with themselves and the rest of the recorded wavefield, which makes them barely or not at all recognizable on individual traces. In addition, the steep inclination of diffracted contributions together with their overall weakness often tends to lead to the misconception that uncorrelated noise rather than coherent signal is observed. Confronting these major challenges, the successful and non-invasive separation of diffractions evolved into a major cornerstone of the field. The important works of Khaidukov et al. (2004), Fomel et al. (2007) and Moser and Howard (2008) arguably were the first to rigorously establish integrated strategies for the successful separation and focusing of diffracted energy to arrive at highly resolved images of the subsurface. Extending the formerly dominant mindset of migration-type imaging, Sava et al. (2005) and Fomel et al. (2007) additionally suggested to use the separated diffracted wavefield for velocity model building. From there onwards, a variety of new separation and imaging workflows were developed, which helped diffraction-based imaging and inversion gain momentum in recent years (e.g. Bansal and Imhof, 2005;Berkovitch et al., 2009;Dell and Gajewski, 2011;Klokov and Fomel, 2012;Santos et al., 2012;Dafni and Symes, 2017;Bauer et al., 2017).
Originating from the works of Berkovitch et al. (2009), Gajewski (2011) along with Bakhtiari Rad et al. (2018), it was recently observed that local coherence measurements can lead to an automated and highly non-invasive separation of weak diffracted contributions from the rest of the wavefield (Schwarz and Gajewski, 2017b;Schwarz, 2019). In this framework, in contrast to previous investigations, the reflected rather than the diffracted wavefield is specifically targeted through coherence measurements and subsequently subtracted to reveal less dominant formerly hidden contributions and make them accessible for further processing. Through its ability to suppress particularly strong interfering reflections, this approach shares strong similarities with the process of plane-wave destruction (Fomel, 2002;Fomel et al., 2007). As the concept of coherence, in addition to targeted wavefield separation, also allows the enhancement and amplification of weak signals and becausejust like diffraction itself -it is deeply rooted and intuitively treated in the field of optics, it underpins essentially all considerations that will be made in this chapter.
Here the general aim is to complement existing theoretical accounts (Klem-Musatov, 1994;Klem-Musatov et al., 2016a,b), with an intuitive easy-to-access introduction to the topic of seismic diffraction. In contrast to more complete accounts, which have an emphasis on the forward problem, i.e. the modeling of wavefields, here I take a more pragmatic, datadriven approach, in which the practical utilization of diffracted signatures is the dominant incentive. A more rigorous and complete treatment of diffraction phenomena, in terms of derivations from first principles can be found in the aforementioned works. Following a brief introduction of diffraction types and definitions and the deep connections to optics, it will be demonstrated by means of coherence arguments and conventional Kirchhoff migration, how seismic diffraction paints a detail-rich complementary picture, backed up by controlled synthetic investigations and complex academic and industry-scale field data examples. Following this treatment of coherent diffraction separation and imaging, additional unique properties of the diffracted wavefield and their potential, as well as current limitations imposed by insufficient acquisition design are illustrated through insightful examples. . While for reflection, equation 1 amounts to Huygens's envelope construction, the Kirchhoff operator is accurate for the diffraction case. R s and R g denote local wavefront curvatures at source and receiver, which can be defined geometrically, when time migration is concerned.

SEISMIC DIFFRACTION
Seismic subsurface imaging has been strongly influenced by optical imaging. This becomes particularly apparent when wavefield-based processing and, in particular, the concept of migration are concerned. While many developments in optics utilized the assumption of monochromatic waves, a rigorous framework for the treatment of pulsed, partially coherent light has also been formulated (Born and Wolf, 2013), which closer resembles the case of typically localized disturbances that seismic waves represent. The process of diffraction was first encountered in the study of light phenomena and found its strongest experimental manifestation in Young's slit experiments. In addition, the wavefield concepts of coherence and interference will prove to be very useful in confronting practical challenges that arise in seismic diffraction imaging. To provide a solid foundation and a deeper intuition for diffraction phenomena in general, the following passages briefly review some important principles of optics that will turn out to be useful in seismic applications. For a more detailed historical overview and a deeper treatment of classical optics, I refer the reader to the classic text book of Born and Wolf (2013) or the recently published work by Klem-Musatov et al. (2016a).

Optical roots
According to historical accounts Leonardo da Vinci was the first to have observed that some portion of the light enters the shadow of an obstacle, which is not explained in the frame of geometrical optics (Meyer, 1934). However, it was only centuries later that Francesco Maria Grimaldi (1665) sought a first systematic treatment of the process and coined the term diffraction to refer to it (e.g. Cecchini and Pelosi, 1990). Christiaan Huygens, one of the first proponents of the wave theory of light then introduced a construction, in which arbitrary wavefronts could be thought of as envelopes of many elementary circular wavefronts originating from imaginary point sources that are locally excited by the incoming energy. A fundamental limitation of the Huygens construction was the fact that only wavefronts -likewise a concept complementary to the ray picture -was concerned. Although ray-theoretical, asymptotic (  led to the development of a first theory of interference, it was Fresnel who in his memoir combined Young's treatment with the wavefront envelope construction to arrive at a first complete (classical) wave theory of light that is principally capable of explaining all observable phenomena (e.g. Born and Wolf, 2013). The resulting Huygens-Fresnel principle simply states, that arbitrary wavefields and their interaction with structural changes can be described by viewing them as a superposition of elementary wavefields emitted by point sources.
It is a unique and very powerful property of the wave equation, that the sum of separate solutions likewise forms a solution. This superposition principle can be physically observable, in the form of noticeably different wavefield components interfering with each other, but can also be of conceptual nature, in that coherent wavefields of arbitrary shape can be thought of as being composed of simple ones. It is exactly this question of observability that turns out to be at the heart of seismic diffraction imaging. Building on the unifying work Seismic diffraction  Figure 2). In the transitional regime (for the case of 50 scatterers) diffraction tails do not fully interfere destructively resulting in a non-negligible contribution to the coda of the main event.
of Fresnel, at the end of the nineteenth century, Kirchhoff provided a solid and practical mathematical foundation that is still influencing a variety of imaging applications across different disciplines (e.g. Lyrintzis, 1994;Daniel et al., 2003;Zhuge et al., 2010). Kirchhoff's integral formulation and its adaptation to seismic problems can be considered one of the major milestones, as even nowadays a majority of seismic depth imaging workflows are still based on it (e.g. Etgen et al., 2009). It can be viewed as a practical approximation to the Huygens-Fresnel principle and can be formulated as a generalized diffraction stack where I denotes the reconstructed image amplitude. In equation 1, ξ = (x 0 , z 0 ) is the image point in depth, whereas ξ = (x 0 , t 0 ) represents its counterpart in time migration (e.g. Hubral, 1977;Schneider, 1978;Schleicher et al., 1993a). In practice, surface observations Seismic diffraction of a scalar wavefield D (representing, for example, pressure) are integrated over a limited spatial aperture, ideally corresponding to the projected first Fresnel zone. To arrive at the correct source pulse after this summation, somewhat not intuitively, the temporal derivative of the wavefield rather than the wavefield itself needs to be considered (Newman, 1990). In equation 1 δ represents Dirac's delta function, w denotes an amplitude weight, and x is a vector parameterizing the lateral positions of sources and receivers located in the limited surface aperture considered in the migration. The Kirchhoff integral thus links the wavefield in depth (or image time, indicated by the subscript 0) with measurements carried out at the acquisition surface. The traveltime t diff = t s + t g corresponds to diffraction at a pointscatterer and t s and t g are the constituent one-way traveltimes linking the image point ξ with the surface positions of the sources and receivers, respectively (compare Figure 1(a)).
By incorporating a weight function w and taking the temporal derivative of the wavefield in equation 1, amplitude-preserving migration with the Kirchhoff integral is enabled (Bleistein, 1987;Schleicher et al., 1993a). However, as was demonstrated in the context of diffraction imaging, if only structural information is desired and absolute amplitudes are not of strict importance in the migrated result, an unweighted diffraction stack (w = 1 and ∂D/∂t = D) leads to drastic simplifications (Khaidukov et al., 2004;Moser and Howard, 2008). Owing to the fact that traveltimes are invoked, which in practice are typically computed by means of ray tracing or solving the eikonal equation, equation 1 is intrinsically of approximative nature (e.g. Etgen et al., 2009). However, it was found by other authors and will become apparent in the following that Kirchhoff migration is well-equipped for the task of diffraction imaging. Figure 1 compares the underlying imaging principle, which is of conceptual nature for reflections, with the case of observable diffraction for which the Kirchhoff operator provides the accurate kinematic description. In the appendix, a detailed discussion and derivation of Kirchhoff's diffraction integral 1 is provided.
Likewise originating from optics, Babinet's principle states that a small localized scattering structure and a small gap, i.e. a small disruption in an otherwise smooth material interface, are equivalent in that they both cause a kinematically equivalent diffraction pattern. Its importance becomes immediately obvious when typical geological structures are considered. Small gaps, faults or intrusions constitute important structural features of the Earth's crust and are likewise honored by Kirchhoff's theory. Figures 2 and 3 illustrate the validity and importance of the Huygens-Fresnel principle for synthetic wavefields diffracted at a varying number of point scatterers, simulated in the Fourier domain with accurate time integration via the Rapid Expansion Method (REM, Tessmer, 2011). Depending on the seismic wavelength, the process of reflection and transmission can be fully reconstructed, if a sufficient number of point scatterers (exciting an equal number of elementary waves) is considered (Figure 2). Emphasizing the special role the process of diffraction plays, it can be observed in Figure 3 that the transition from concept to physical reality likewise translates into the data (time) domain. Circular wavefronts in object space or the geology appear as hyperbolic features in data space or the image, that can be clearly distinguished for a limited amount of scatterers. In full analogy, Babinet's equivalence principle is likewise illustrated in Figure 4 by means of the accurate simulation of diffraction occurring at a single scatterer (Figure 4(a)) and a single slit in an otherwise continuous interface ( Figure  4(b)).

Diffraction types and properties
In the framework of conventional geometrical ray theory, the process of diffraction is not explained and can only be postulated (Keller, 1962). In fact, diffraction can be summarized and defined as the collection of all observations that appear forbidden in the classical ray picture (Born and Wolf, 2013). However, once their occurrence is postulated for smallscale structural features, in the far field diffracted waves, just like reflections or direct arrivals, can often be reliably approximated by rays (Keller, 1962;Klem-Musatov, 1994). Following the systematic classification of Keller (1962), different types of diffraction can be distinguished. Some of the most important diffraction phenomena are illustrated in Figure 5). While the classic text book on diffraction imaging by Klem-Musatov (1994) closely follows the definitions provided by Keller (1962), there sometimes appears to be conflict in how the process of diffraction is defined in the geophysical literature. In the seismological community, for example, diffraction at the Earth's core-mantle boundary can Seismic diffraction  (e)). The special case of surface diffraction (d), occurring when incident wavefields graze a locally (partially) smooth interface, is of particular importance in earthquake seismology (Aki and Richards, 2002). routinely be observed (compare Figure 5(d)). According to Keller (1962), this observation can be identified with the process of surface diffraction (Aki and Richards, 2002). However, it is noted by Keller (1962) that the same waves are also referred to as creeping waves in a different context. So it can be concluded that terminology needs to be treated with care to prevent confusion.
Similar to the example of surface diffraction, depending on the context, it is widely debated whether edge diffraction, i.e. diffraction occurring at sharp edges ( Figure 5(b) or (c)), or point diffraction occurring, for example, at structure-less point-like objects ( Figure  5(a)) can be considered more fundamental. Without intending to add to this debate, I believe it is somewhat helpful to consider earthquakes, which are the most energetic natural sources of seismicity and have arguably been studied much more intensely than seismic diffraction processes. In the far-field, which corresponds to the ray-theoretical realm of validity, the concept of a point source is useful, whereas in the near field, ray densities increase significantly leading to the need for corrections to be applied. The same holds true for diffraction. The notion of a point diffraction, in a sense, can be related to the ray-theoretical view underlying Huygens elementary wavefront construction, which ignores variation in amplitude and phase, while the process of edge diffraction is physically more complete and honors the full wavefield properties, which lets it appear more closely related to Fresnel's combination of interference with the wavefront picture.
For the sake of simplicity, and to illustrate the remarkable robustness and usefulness of Kirchhoff migration, all the examples presented in this chapter will utilize the point diffraction or scattering picture of uniform radiation, which is well describable through Seismic diffraction wavefronts. It is important to note, that in realistic media the concept of omni-directional scattering often represents an idealization. Diffraction at edges and at localized disruptions of geological horizons represent more common causes and the amplitude behavior is often more complicated (e.g. Greenhalgh and Manukyan, 2013). However, from a wavefront perspective, i.e. taking a kinematic viewpoint, it can be argued that more intricate diffraction patterns -following the superposition principle -can be separated into two parts, one honoring Snell's law, the other representing a highly illuminating, uniformly radiating component. It is this illuminating component that encodes very localized (high resolved) information about the scattering geometry and provides superior illumination that can be utilized in seismic imaging.
It can generally be argued that, in contrast to conventionally utilized reflected energy, diffraction can only be fully understood, if wave fields rather than single isolated recordings are considered. Its occurrence is intimately linked to abrupt, very localized changes in medium properties on the order of the predominant seismic wavelength, laterally and vertically. Consequently, in contrast to reflection which occurs when abrupt changes are encountered locally in only one spatial direction, diffraction is caused by discontinuous features that have a radius of curvature that is small compared to the seismic wavelength. Two of the most characteristic properties of diffracted fields originating from point-like features or edges are their overall weakness (due to the rapid geometrical spreading) and the fact that, amplitude and phase considerations aside, their overall shape in space and time is essentially independent of the encountered discontinuity. Although, again for simplicity, examples will be restricted to two dimensions, it is important to appreciate that diffraction is intrinsically three-dimensional and that more types can be distinguished spatially (compare Figure 5(e)). However, for the same reasons as before, the processes of point or edge diffraction remain particularly useful in that the connected wavefronts radiate uniformly and can be viewed as a physical (observable) manifestation of elementary waves originating from secondary sources whose locations are related to Earth structure.
Returning to Figure 2 and Figure 3, it becomes clear that in the Huygens-Fresnel picture of interference, a terminating reflector, up to the point of termination, can be viewed as a collection of a large number of closely spaced point diffractors, whose superposition results in the reflected wavefield. At the termination point, however, the destructively interfering counterpart is missing, resulting in an observable diffracted wavefield that albeit its faintness can be utilized for further processing. In this way diffractions indeed share some crucial characteristics of Green's functions (Dirk Gajewski, personal communication) and it comes as no surprise that migration techniques, dealing with squeezing out the propagation effects and leaving only structural information in the seismic image, have been shown to be naturally equipped constructs for diffraction imaging in time and depth (e.g. Moser and Howard, 2008;Khaidukov et al., 2004;Silvestrov et al., 2016). As will be shown in the next section coherence represents a collective wavefield property that helps to couple the useful wavefront picture to waveform characteristics.

DIFFRACTION SEPARATION AND IMAGING
Research in seismic diffraction has so far mostly been concerned with two overarching goals: (1) the migration-type imaging of small-scale discontinuities and related to that, (2) the successful separation of the often orders of magnitude weaker diffracted contributions from the Seismic diffraction rest of the wavefield. Starting with the pioneering works on wavefield imaging (Hagedoorn, 1954), crustal fault characterization (Krey, 1952;Kunz, 1960), and the development of a more rigorous theoretical framework (Trorey, 1970(Trorey, , 1977Berryhill, 1977;Klem-Musatov andAizenberg, 1984, 1985;Klem-Musatov, 1994), the problem of interference of diffracted with the stronger reflected wavefield components was soon recognized. As a result, first strategies for their extraction and/or practical utilization started to emerge in the 1970s and 1980s (Hubral, 1975;Schilt et al., 1981;Kanasewich and Phadke, 1988;Landa et al., 1987). Generally, two different schools of thought exist -one first aims at separating the diffractions directly in the time domain using some sort of coherence argument (e.g. Fomel, 2002;Bansal and Imhof, 2005;Fomel et al., 2007;Berkovitch et al., 2009;Dell and Gajewski, 2011;Schwarz and Gajewski, 2017b) followed by subsequent focusing; the other supresses the reflected energy directly during migration (e.g. Moser and Howard, 2008;Klokov and Fomel, 2012;Yin and Nakata, 2017;Dafni and Symes, 2017). While both these mindsets have their advantages, it can be argued that a distinct time-domain separation step generally leaves the user with more flexibility. In addition, insufficient knowledge of the velocity structure does not automatically compromise the extraction. On the contrary, time-domain prestack diffraction extraction principally enables dedicated diffraction-driven velocity inversion workflows that utilize the unique illumination properties of these weak wavefields (Sava et al., 2005;Fomel et al., 2007;Santos et al., 2012;Bauer et al., 2017).
One of the main benefits of dedicated diffraction processing lies in the potential for high-resolution imaging of small-scale discontinuities. This chapter on seismic diffraction is aimed at widening this view to other potential impactful applications. However, as the former conventional branch has matured considerably over the past decades and, more importantly, because it without a doubt adds value in the pursuit of the ultimate goal of an accurate characterization of the subsurface, this section will briefly review a straightforward way of achieving goals (1) and (2) in an integrated fashion. Building on the intuitive notion of wavefronts, the presented framework is reasonably simple to implement and should therefore be accessible and reproducible for people in both, industry and academia. After explaining the inner workings of this method, the section will be concluded with a brief discussion of a simple class of time-domain focusing operators and the application to a variety of challenging synthetic and field data examples.

Coherence and wavefronts
Before the recent rise of ambient noise interferometry in earthquake seismology and seismic exploration, a signal was regarded useful for imaging if it is coherent. Even surface-related multiple reflections or the effect of ground roll, both previously considered coherent noise, nowadays serve a purpose and can be utilized for imaging and inversion. Originating from optics, coherence is a fundamental property of lasers, which emit very intense narrow beams of electromagnetic radiation. Owing to their finite duration, seismic excitations are only partially or locally coherent, but similar to lasers, through constructive interference, even weak emerging energy can be amplified through beam forming or stacking (Mayne, 1962). In turn, it was found by Taner and Koehler (1969) and Neidell and Taner (1971) that the coherence of a seismic wavefield can be systematically investigated by means of trial data summations within a predefined aperture. As a discrete and reduced version of the Kirchhoff integral, such a data summation over all traces located at x in the neighborhood Seismic diffraction Figure 6: A simple synthetic example illustrating the typical challenge of interference that has to be confronted for targeted diffraction imaging. Shown are (a) the input data, (b) the reflection reconstruction via stacking, and (c) the adaptive subtraction of (b) from (a).
of a reference location x 0 can be written as where C k denotes the reconstructed coherent wavefield, D is the input data volume and t 0 stands for the reference traveltime for which the reconstruction is performed. The index k is introduced to honor conflicting coherent contributions with different emergence slopes and curvatures that might intersect at the reference data point (x 0 , t 0 ). The trajectory t k (x i ) describes the wavefront-consistent traveltimes of the k-th locally coherent event, as they are observed at neighboring measurement locations x i . Corresponding to w in equation 1 the division by the number of traces entering the data fold in the aperture n represents a normalizing weight of the summation.
As is illustrated in Figure 6(a), when diffractions are concerned, the presence of conflicting event dips is very likely as they usually are not an isolated phenomenon but rather appear in numbers and tend to interfere with each other and more energetic reflected contributions. Consequently, the data can be written as a superposition of interfering coherent wavefields and uncorrelated noise N , Within the framework of automated data enhancement and wavefield reconstruction (e.g. Gelchinsky et al., 1999a,b;Jäger et al., 2001;Baykulov and Gajewski, 2009;Xie and Gajewski, 2017) this was considered a fundamental problem of coherent summation. However, there exist several attempts of honoring multiple event dips, ranging from operator extrapolation (Höcht et al., 2009;Baykulov and Gajewski, 2009) to Fresnel volume inspired superposition strategies Xie and Gajewski, 2018). All of these approaches are phrased as optimization problems, in which, conventionally, the semblance norm is maximized using either local or global solvers (or a combination of the two). The additional summation in expression 4, performed within a time window δt around the reference time t 0 , honors the band-limited nature of the signal under investigation. The semblance norm has gained strong attention over the past decades and is still frequently favored over other coherence measures. This can partly be explained by its relative ease of implementation, its computational efficiency and its overall robustness. From a physical perspective it can be viewed as the normalized beam energy or the ratio of stacked to overall intensity contained in the data. Therefore it intuitively reflects a wavefield's potential to self-interfere.
If sufficiently local apertures are considered, a closed-form relationship between the traveltimes t k (x i ) and the reference t 0 can be found. For the most general 3-D prestack case, paraxial two-way traveltimes for source separations ∆x s i = x s i − x s 0 and receiver separations ∆x where p s k and p g k represent two-dimensional slope vectors and M s k and M g k are 2x2 matrices encoding the observed curvature of the considered k-th coherent event on the source and on the receiver side, respectively (Schleicher et al., 1993b). Equation 5 corresponds to a second-order Taylor series expansion of the squared traveltime. The components of the slope vectors p s k and p g k represent horizontal slownesses and have the dimension s/m. For common sub-configurations, like, for example, the zero-offset section (i.e. the poststack domain), the Seismic diffraction common-source or the common-receiver gather, the hyperbolic traveltime formula 5 reduces to formally equivalent expressions that relate to auxiliary one-way wave propagation and, consequently, can likewise be employed in passive seismic investigations (e.g. Schwarz et al., 2016;Schwarz and Gajewski, 2017a;Diekmann et al., 2018). So as a by-product, the use of equation 5 in conjunction with equation 4 allows for the automated extraction of first and second order attributes of the event's shape, which can be directly related to the slope and curvature of the emerging wavefront causing the observed trace-to-trace traveltime differences.
Adaptive reflection subtraction Figure 7 shows the estimated maximized semblance as well as the absolute wavefront emergence angle and a normalized version of the curvature attribute for the 2D synthetic example introduced in Figure 6(a). All three attributes were estimated automatically without any user involvement aside from the predefinition of the semblance time window δt and the size of the data aperture. As can be observed throughout, the reflected energy, characterized by small curvature and emergence angle estimates, is consistently favored by the coherence analysis. Accordingly, as will be more thoroughly demonstrated later, we often have to a good approximation where C ref denotes the coherent reflected contribution at data point (x 0 , t 0 ) and t ref is the corresponding traveltime trajectory. Rather than aiming at honoring the interfering wavefields' full complexity, which is computationally demanding and challenging to implement successfully, it was recently suggested by Schwarz (2019) to make use of the undesired directional Seismic diffraction Figure 9: Prestack data (a), estimated maximized coherence (c), and diffraction separation result (b) for a collection of common-source gathers of the Sigsbee 2a synthetic dataset. Seismic diffraction filter characteristics imposed by expression 2 and, consequently, the semblance (equation 4) to extract diffractions from pre-and poststack seismic, as well as ground-penetrating radar data. In contrast to previously proposed strategies (e.g. Berkovitch et al., 2009;Dell and Gajewski, 2011;Bakhtiari Rad et al., 2018) which directly target the weak diffracted wavefield, a conventional reflection stack is suggested, followed by its subsequent adaptive subtraction from the full input wavefield, To accommodate data imperfections, local waveform variations and more general deviations from the underlying assumption of local coherence, the subtraction process in equation 7 is formulated in an adaptive fashion through the introduction of perturbative time shifts τ 0 and local amplitude scaling coefficients α 0 . In order to preserve even the faintest diffracted signatures, the adaptation, just like the coherent summation and the semblance estimation, is performed within a data aperture. The residual misfit of this adaptation, as well as the scaling coefficient and time shift fields are displayed in Figure 8. All three quantities show distinct correlations in regions where diffractions and reflections heavily interfere. In cases where diffraction is the dominating contribution, as might, for example, be the case at the high-impedance contrast observed at complex sediment-salt interfaces, the estimated slope and curvature attributes can be used to design diffraction-preserving wavefront filters that can objectively inform the coherent subtraction. For more details on the adaptation step and implementational aspects I refer the reader to the work by Schwarz (2019).
In summary, the results of the automated coherent summation and the subsequent adaptive subtraction alongside the full wavefield data input are presented in Figure 6. To demonstrate that this methodology is reasonably flexible to handle different data configurations, the results of the coherence maximization and the diffraction separation, performed on a set of neighboring common-source gathers simulated for the complex Sigsbee 2a synthetic model, are shown in Figure 9. Again, the coherence analysis consistently favors the reflected contributions, resulting in a successful amplitude-preserving separation of the diffracted wavefield in the prestack domain.

Diffraction focusing
Kirchhoff depth migration uses ray tracing or eikonal solvers to compute diffraction traveltimes for trial subsurface image locations ξ = (x 0 , z 0 ) and known velocity structure. In contrast to that, the discussed multi-dimensional coherence analysis, which can be expressed by a very similar mathematical construct, is based on an adaptive and flexible traveltime kernel (equation 5) that is likewise suited to describe diffracted and reflected contributions. While the former, even when formulated in the time domain, is in need of detailed knowledge of a velocity field for computing diffraction traveltimes, the latter in turn can facilitate the estimation of velocity structure from the estimated wavefront attributes (Duveneck, 2004;Bauer et al., 2017).
In fact, when considering the diffraction case and a zero-offset reference frame, i.e. midpoint and half-offset (∆x s i = ∆x g i ) rather than source and receiver coordinates and confinement to the diffraction case and vanishing wavefront inclination (i.e. p s k = p g k = 0), the number of degrees of freedom reduce drastically and the traveltime operator 5 can be used for Kirchhoff-type time migration (Zhang et al., 2001). However, previous investigations Seismic diffraction suggest that double-square-root-type expressions for the traveltime are more suited to accurately describe and incorporate the response from highly curved subsurface features and diffracting structures. Consequently, a range of extensions to the zero-offset subset of operator 5 have been introduced in recent years (e.g. Landa et al., 2010;Fomel and Kazinnik, 2013;Schwarz et al., 2014). While these differ in the way the wavefront attributes appear in the expressions, they share the overall double-square-root shape and are equally accurate for the limiting diffraction case (Schwarz and Gajewski, 2017c). In fact, if parametrized consistently they were demonstrated to provide equivalent accuracy and can be considered largely equivalent descriptions (Schwarz and Gajewski, 2017a;. For an effective auxiliary medium they incorporate the conventional Kirchhoff time migration formula where ∆x m i = x m i − x 0 is the midpoint displacement and h i denotes half the source-receiver offset of the ith trace with respect to the time image location ξ = (x 0 , t 0 ). The migration velocity v corresponds to the root-mean-square velocity which can be directly calculated from the wavefront attributes estimated during coherence analysis. Analogous to conventional Kirchhoff depth migration, the diffraction traveltime defined in equation 8 can be divided into two distinct contributions, one corresponding to the source the other to the receiver leg of the full ray-path.
In academia, target depths typically exceed the maximum source-receiver distances by a factor of two or even considerably more. As a result, crucial lateral illumination is not recorded, thereby intrinsically limiting lateral resolution. Also, detailed a-priori information required for successful prestack depth imaging is rarely available, which makes the related computational demands difficult to justify. In general, due to the aforementioned limitations, pragmatic yet stable and computationally efficient imaging approaches with a good Seismic diffraction performance with respect to low signal-to-noise ratios are ideally suited in these scenarios. Despite its unique properties, the diffracted wavefield is still rarely utilized in academic practice. The main benefit of a distinct time-domain separation strategy is that firstly, velocity uncertainties do not immediately affect the extraction, and secondly, the subsequent imaging workflow can either operate in time, which decreases accuracy but increases computational efficiency and stability, or in depth, if the circumstances are favorable. To ensure comparability and to demonstrate that diffraction imaging does not have to be extremely sophisticated to produce highly resolved and useful images, the following examples utilize a time migration kernel and, consequently, make use of expressions 1 and 8.

Diffraction vs. reflection images
As a first simple example of data-driven diffraction imaging, Figure 10 shows a close-up of the Sigbsee 2a synthetic dataset which was made available to the public by the Subsalt Multiples Attenuation and Reduction Joint Venture (SMAART JV). Aside from an extended salt structure mimicking the Sigsbee escarpment below the Gulf of Mexico, the model depicted in Figure 10(a) reveals complicated localized structures in the form of faults, and a horizon of small-scale velocity perturbations which are deemed to cause diffraction patterns in the data. After wavefield separation is performed through the coherent strategy described above, the (conventional) full-wavefield image, comprised of reflected and diffracted contributions (Figure 10(b)), is compared with a reflection-only migration (Figure 10(c)).
The reflection image lacks resolution in essentially all regions that deviate from the picture of a smoothly layered medium. Although this reconstruction is clean in appearance and might add value in tracking horizons, it lacks sufficient detail on where the geodynamically important faults are located in the model. diffraction imaging are supposed to complement each other and that diffraction images alone will not necessarily be of value to an interpreter (e.g. Khaidukov et al., 2004;Moser and Howard, 2008). Confirming this notion of complementarity, the corresponding diffraction image, presented in Figure 11(a), reveals detailed features of the fault system, provides an account of the localized perturbations in the lower half of the model, and even manages to accurately reproduce a detailed pattern of distinguishable diffracting edges (indicated, for example, by the apparent phase reversal near the edge). Although these diffractions were unintentionally introduced through the too coarse definition of the synthetic model, it can be concluded that the diffraction separation was sufficiently successful to enable imaging of features at a resolution approaching the model's spatial discretization. Of course, the diffracted wavefield has not only virtue on its own, but also contributes to the perceived sharpness of the full-wavefield imaging result shown in Figure 10(b). However, the identification and detailed structural characterization of the model coarseness and the highly relevant extent of the fault system is much easier to accomplish when the diffraction image is considered.

Imaging attributes
When only structural imaging is concerned, a noticeably improved signal-to-noise ratio can be achieved when the semblance norm (equation 4) rather than the conventional wavefield sum is estimated during migration (e.g. Landa and Keydar, 1998;Müller, 2000;Heincke et al., 2006). This is generally possible for Kirchhoff-type migrations, as summation over an aperture is directly related to the definition of semblance (Neidell and Taner, 1971). It is important to note the critical difference in treatment that reflected and diffracted contributions experience during Kirchhoff migration. Reflections are merely repositioned and the summation operator is only tangent to the wavefield (meaning that only a subset of the fold within the aperture contributes constructively). Diffractions on the other hand -accurate velocity information provided -fully contribute to the resulting migrated amplitude. This seems reasonable, as diffracted energy, due to its more uniform radiation, was distributed over a larger area and, consequently, must be gathered and refocused using a larger aperture. Semblance, however, represents a normalized quantity and turns out to be largely independent of the absolute amplitude of the signal, which leads to a natural suppression of reflected energy during migration.
Upon closer inspection, the migration of diffraction-only data strongly resembles the procedure of passive-source back-projection which relies on focusing and provides a robust and highly resolved means of detection and localization of natural seismicity in well-instrumented areas in earthquake seismology studies (Ishii et al., 2005). Inspecting the earthquake seismologists' toolbox, the notion of the so-called N-th root stack, which is often used for the detection of weak events in very challenging low signal-to-noise scenarios (Rost and Thomas, 2002), complementary to the varimax norm (Wiggins, 1978;Fomel et al., 2007), can serve to introduce non-linear selectivity and improve robustness: In Figure 11 this N-th root semblance is compared with conventional diffracted wavefield focusing and the linear semblance norm. As expected, both coherence measures reveal clearer delineation of the conjugate faults and the highly localized features in the bottom half of the model (compare Figure 10(a)), while indicating to be less prone to artifacts. In the 10-th root semblance image (Figure 11(c)), an overall higher signal-to-noise ratio and improved reconstruction of the edge diffractions due to the model discretization in the shallower regions can be observed. To prevent undesired destructive interference due to an amplitude reversal at the apex, which typically can be encountered for edge diffraction, both semblance versions are evaluated twice, once with and once without correcting for this phase change. In the final image, for each sample, the maximum value of these two fields is displayed. Aside from normalized quantities like the linear or N-th root semblance norms, tools from complex trace analysis, such as the trace envelope and its first and second order derivatives (Figure 12(a)-(c)) can be employed to additionally sharpen the diffraction imaging results. Due to the perceived robustness in the presented example, in the following, diffraction images were achieved using the newly introduced N-th root semblance measure 9.

Faults, fractures and unconformities
From the previously conducted systematic comparison of full-wavefield and diffraction images it can be concluded that -and this is not entirely unexpected -seismic wave diffraction is responsible for the reconstruction of a majority of the laterally resolved features of the model. In fact, as these comparisons likewise suggest, the diffracted wavefield provides an overall sense of sharpness, which cannot be achieved with the reflections alone. Various other authors have commented on the high-resolution imaging potential of diffractions and essentially all of them conclude in their respective analysis that a successful separation is crucial for these highly resolved features to be distinguishable from artifacts that arise from data imperfections, too coarse spatial wavefield sampling and inaccuracies in the retrieved subsurface velocity model ( and Howard, 2008;Silvestrov et al., 2016). A comparison of full-wavefield and reflection image, as for example presented in Figure 10, suggests that, in a complementary sense, the separated reflected wavefield is overall less artifact prone and appears smoother and somewhat cleaner, which among other things might prove useful for the reliable tracking of horizons in seismic interpretation.
Conversely, diffraction images highlight and emphasize features that are very localized laterally, like abrupt truncations of the aforementioned horizons, erosional surfaces, the carving imprints of river channels, and, most prominently, faults and fracture systems cutting through sedimentary strata or the crystalline basement of the Earth's crust. An integrated (full-wavefield) view may provide a somewhat complete picture when velocity model estimates are accurate and imaging algorithms are sophisticated and well implemented. However, owing to the general faintness of diffractions, their voice tends to be overheard and often does not carry sufficient weight to sustainably have an impact on the image quality and to be accurately honored in the conventional, reflection-biased processing chain. It is only through their separation that targeted reflection and diffraction processing can be fully tailored to the peculiarities and distinct properties and challenges of these insightful wavefields. As a result, a complementary picture of the subsurface emerges, that in some sense, may provide more insights than the mere sum of its parts.

Diffraction in interpretation
From a material science or optics perspective, isolated diffraction patterns are often directly related to material imperfections. In an Earth science context, it is exactly these imperfections that are of interest to an interpreter of seismic images, as they, in some sense, encode the dynamic history of the crust. Conventional and established attributes such as image coherence or image curvature (Bahorich and Farmer, 1995;Marfurt et al., 1998; Seismic diffraction  Chopra and Marfurt, 2007) have proved to bear the potential of being reliable indicators of faults. However, the successful and precise delineation of discontinuous features such as faults and erosional unconformities still remains an area of active research, where recent improvements have been made with the help of conventional and more specialized image processing and analysis techniques. Amongst others these include spectral analysis (Liu and Marfurt, 2007), structural tensor analysis (Wu and Hale, 2016;Wu, 2017), optimal surface voting (Wu and Fomel, 2018) or deep learning (Araya-Polo et al., 2017).
Despite the tremendous value these purely image-driven approaches have delivered, it was recognized early on that they all intrinsically rely in some way on the success of the preceding imaging step (Khaidukov et al., 2004;Moser and Howard, 2008;Berkovitch et al., 2009). Diffraction images are different in that they are naturally coupled to the aforementioned discontinuities. As illustrated in the the second section, diffraction can be viewed as the birth of a new wavefield that focuses -just like a passive source -at the discontinuous structure that caused it. In consequence, diffraction images represent physically derived attribute maps. In addition, the use of robust image functions such as the suggested N-th root semblance norm (equation 9) permit an increase in the signal-to-noise ratio and reduce artifacts directly at the imaging stage, which in turn might improve subsequent image and data-driven processing. Consequently, there exist first attempts to use diffraction imaging as a robust alternative for seismic interpretation (e.g. Tsingas et al., 2011;Sturzu et al., 2014;Schoepp et al., 2014;Decker et al., 2014;Tyiasning et al., 2016;Dell et al., 2018).
In the following, by means of two challenging yet controlled synthetic examples, I want to illustrate how robustly dedicated diffraction imaging can reveal structurally important small-scale subsurface feature whose details normally tend to be drowned out by the strong amplitude reflection foreground. The first example is concerned with another close-up of the Sigsbee 2a synthetic dataset which was also investigated when discussing the mechanism of diffraction focusing. This time, another part of the model near the top of the complicated salt body is investigated. Figure 13 shows the comparison of the actual seismic velocity grid alongside its reconstruction based on the diffracted wavefield that was automatically separated before imaging.
To illustrate the overall robustness of this approach, the image represents a bruteforce constant-velocity (1.6 km/s) focusing result employing the newly introduced N-th root semblance. As before, aside from the clear delineation of a fault at the top of the salt tip, largely every single stair step in the discrete velocity grid can be resolved. In addition, the steep flanks of the rugged top-of-salt in the bottom part of the model are clearly revealed through diffraction imaging. Due to the fact that the N-th root semblance (equation 9), just like its conventional linear counterpart, is a normalized quantity, amplitude-strong features, connected to significant impedance contrasts appear equally well-resolved and bright as, for example, the stair steps of the velocity grid, which are the cause of extremely faint signals that can barely be recognized on individual traces (Meyer, 1934;Heincke et al., 2006;Schwarz, 2019). Although absolute amplitudes can also be used for imaging, destructive interference does not occur for this case, which often results in noticeable migration artifacts. Semblance and its non-linear counterpart are different, in that they uniquely combine the notion of intensity, which is most familiar to what we perceive with our eyes, with phase information leading to destructive interference for non-physical migration trajectories and uncorrelated noise.
The Sigsbee 2a synthetic dataset was simulated based on a moderate macro-gradient Seismic diffraction with small-scale stratigraphic velocity variations. Realistic but controlled variable-density benchmark datasets likewise are publicly available and can be used for the investigation of the potential of diffraction imaging of small-scale density perturbations in the subsurface. One of these challenging datasets is the 2004 BP velocity benchmark (Billette and Brandsberg-Dahl, 2005). In contrast to the Sigsbee example that represents a detailed reproduction of a real-world geologic scenario, the BP model can be viewed as a hybrid in that it synthesizes a variety of different geologies, important in the context of hydrocarbon exploration and likewise challenging to image seismically. The Sigsbee velocity structure can be considered reasonably complicated and fine-structured, whereas the structural complexity of the BP model in large parts is not connected to the relatively simple and well-behaved velocity field, but rather arises from a very detailed and heterogeneous density distribution. In Figure 14 -again for a top-of-salt regime -a close-up of the velocity (Figure 14(a)) and the respective density model (Figure 14(b)), together with the diffraction reconstruction (Figure 14(c)) are displayed. Like before, to emphasize the robustness of coherence-based diffraction imaging, a largely inaccurate constant velocity (1.6 km/s) was chosen for wavefield focusing. Despite the crude assumption of a constant focusing velocity, a stunningly detail-rich reconstruction of the rugged top of the salt body and small-scale discontinuities within the sedimentary overburden, primarily caused by density fluctuations, results, which provides a complementary, highly resolved structural image.

Diffraction imaging in the field
Generally, there is a consensus that seismic wave diffraction routinely occurs in realistic media and that these weak signals are likewise routinely recorded in the field. However, with some impressive exceptions (e.g. Sturzu et al., 2014) it can be argued that diffraction imaging still does not form an integral part of the seismic processor's toolbox (Landa, 2012). One of the reasons is arguably that the overall weakness of their de-focused wavefields, as they appear in the raw recorded data, makes diffractions barely recognizable and their impact on image quality is only perceived indirectly. Additionally, the aforementioned strong interference with the conventionally utilized reflected contributions makes a successful separationdeemed necessary for targeted imaging -very challenging in practice. Also, in particular in the framework of non-commercial surveys with a purely scientific objective, data are generally imperfect and acquisitions only record a limited portion of the back-scattered energy, making focusing strategies prone to artifacts whose amplitudes often approach the strength of the weak diffracted component itself.
Using a non-invasive extraction strategy, Schwarz and Gajewski (2017b) recently demonstrated that a faint background wavefield can be recovered, even if single-channel field data recorded in a highly complex geological setting are considered. Owing to their special illumination capabilities, the detection and characterization of diffractions, even under these unfavorable conditions, can have profound implications for a majority of the reduced measurements that are nowadays conducted in academia. The inherent kinematic symmetries and correlations will be more closely examined in the final section of this chapter. Seismic diffraction

Industry-scale long-offset acquisition in offshore Israel
Here, I want to follow up on previous work (Schwarz and Gajewski, 2017b;Schwarz, 2019) and demonstrate that coherent diffraction separation and focusing, in the simple way it is described in this chapter, can lead to highly resolved images that provide interpretational value that is highly complementary to the full-wavefield migrated images which are normally consulted. To investigate applicability for both extremes of high and low-fold seismic field data, an industry-scale acquisition and the very reduced single-channel acquisition already studied by Schwarz and Gajewski (2017b) are subsequently analyzed. The former constitutes in a sophisticated multi-channel measurement campaign conducted by TGS off the coast of Israel in the Eastern Mediterranean, whereas the latter was performed by academia near Santorini in the Aegean Sea. By means of more conventional processing, both datasets Seismic diffraction were already demonstrated to contain information on very intriguing tectonic features either relating to complicated tectonics that are connected to the formation of massive salt complexes (e.g. Netzeband et al., 2006), or arose from nearby persistent volcanic activity (e.g. Hübscher et al., 2015).
A conventional full-wavefield prestack time migration of the industrial dataset is shown at the top of Figure 15. As can be deduced from this image, underlying structures are very heterogeneous and the subsurface can broadly be divided into at least three distinct regimes. In particular as observed in the left half of the section, the signature of a massive salt body originating from the Messinian salinity crisis (Krijgsman et al., 1999;Gradmann et al., 2005) is overlain by a thick sheet of sediments. From these mostly horizontally stratified sediments upwards to the sea bottom at 1.8 s two-way time a complicated turbulent complex characterized by chaotic but comparably strong reflectivity patterns is followed again by horizontally smooth stratification of sediments. From left to right (corresponding Seismic diffraction to a progression towards the coast of Israel) firstly, a thinning of the salt sheet and secondly, tectonic tilting of different lateral units, separated by large but barely visible faults are revealed. In addition, it may be argued that the overall reflectivity of the subsurface decreases, which may in part be attributed to the fact that structures become more laterally complex and events are generally steeper as the coast is approached.
Displayed at the bottom of Figure 15 is the corresponding diffraction image, resulting from full prestack diffraction separation and focusing based on a simple vertical velocity gradient and the newly introduced N-th root semblance (equation 9). Strikingly, and this becomes even more apparent in the close-up investigated in Figure 16 (indicated by the frame in Figure 15), the diffracted wavefield provides highly resolved structural detail in parts of the section, where the conventional full-wavefield image lacks resolution. This is in full accordance with the previously mentioned overall complementarity of reflection and diffraction imaging. In particular, the chaotic turbulent features and faults, both of which are important indicators of dynamic processes, appear bright and detail-rich in the diffraction image, whereas the undisturbed sedimentary (sub-) units are largely transparent to the diffracted wavefield. For more details on the interaction between the salt body and the overall tectonic regime, I refer the reader to the work of Gradmann et al. (2005) and Reiche et al. (2014).

Academic single-channel acquisition near Santorini
To illustrate the versatility and robustness of coherent diffraction imaging, a very reduced academic marine field dataset, consisting of only a single near-offset channel, is briefly presented and described. The seismic data were recorded above the Anydros Basin, close to the island of Santorini in the Aegean Sea. The region was found to be largely influenced structurally through prevailing volcanism (Hübscher et al., 2015;Nomikou et al., 2016b). One of the dominant features is a large submarine volcano -the Kolumbo seamountwhose activity sustainably shaped the tectonic landscape of the region and is even thought to have been involved in Tsunami generation (Nomikou et al., 2016a). As examples, Figure  17 and Figure 18 show two different sedimentary environments adjacent to the Kolumbo seamount. Owing to the availability of only a single channel, conventional velocity analysis was rendered impossible, so all images in both figures were gained through a brute-force approach with a constant velocity of 1.6 km/s. Like before, the diffraction images show the normalized N-th root semblance, which again proves to be a robust and largely artifact-free diffraction imaging attribute. In Figure 17 the crystalline basement of the crust, rising to the upper few kilometers can be identified as a largely diffuse body. While the strong velocity contrast between the sediments and this basement can clearly be identified in the full-wavefield image, the lower portions are dominated by migration artifacts due to the lack of precise velocities and their severe underestimation during the reconstruction.
The diffraction image on the other hand, owing to the non-linearity of the N-th root semblance norm, is generally less prone to artifacts, so internal unconformities and faults, which are mostly hidden in the conventional image, become distinguishable. In particular in the right part of the section at around 1 s two-way time, a detail-rich system of faults can be delineated in the diffraction image, which remains largely obscured in the conventional migration. Once more supporting the notion of complementarity, Figure 18 shows a comparison of conventional migration and diffraction imaging for a small confined sedimentary Seismic diffraction basin at the flank of the submarine volcano. The conventional full-wavefield image shows its strength in the largely undisturbed sedimentary sequences revealing strong reflectivity, whereas the diffraction image helps to delineate a sequence of erosional surfaces and internal faults and fractures, which were previously hidden.
Both field data examples suggest that, provided the spatial sampling of the emerging seismic wavefield is reasonably dense in at least a single data sub-domain, diffraction imaging becomes feasible and has the potential to provide valuable additional information particularly relevant for the study of dynamic processes related to erosion and tectonics. The lack of laterally resolved velocity information can be considered a fundamental problem in academic studies, where targets are typically deep and lateral illumination is very limited. In addition to permitting high-resolution imaging of faults, fractures and unconformities, the diffracted wavefield radiates uniformly and encodes valuable lateral velocity Seismic diffraction information which can be extracted with tailored inversion strategies (e.g. Fomel et al., 2007;Bauer et al., 2017). Thus, diffractions may help to overcome the traditional shortcomings of these low-fold recordings. In the following, the unique properties of diffractions and their potential impact on seismic data acquisition and imaging will be discussed in more detail.

IMPLICATIONS AND NEW DIRECTIONS
It was previously demonstrated that with the help of non-invasive separation techniques, such as the suggested adaptive subtraction strategy, a rich diffracted wavefield can be utilized for targeted imaging. In many ways, diffractions turn out to be simpler wavefields, in that with conventional Kirchhoff migration, for example, using reasonably accurate velocity information, the utilized operators are not only tangent but coincide with the full response (or a major portion of it). So naturally, Kirchhoff migration utilizes the full redundancy and coherence of the diffracted wavefield, which even for short-offset acquisitions and despite their overall weakness, leads to improved signal-to-noise ratios and a high resolution in the resulting images. In addition to this localized highly resolved nature, diffracted wavefields provide the same illumination that a point source excited at the diffractor location would provide. Accordingly, diffracted wavefields, in the far-field approximation only encode propagation effects from the scatterer location to the receivers, which makes them very suitable for velocity inversion (e.g. Sava et al., 2005;Fomel et al., 2007;Bauer et al., 2017).
These approaches of diffraction-based velocity model building often rely on wavefield focusing -a process that was discussed extensively in the context of imaging in the previous section. Here, I want to draw the attention to additional striking properties of diffractions, which directly follow from their uniform radiation. As diffracted contributions often reveal significant inclination angles when emerging at the registration surface, spatial aliasing constitutes an important and often under-appreciated practical challenge. Accordingly, I Seismic diffraction will briefly start with discussing sampling criteria that need to be met and how unintentional wavefield under-sampling can compromise the reconstruction of subsurface structure, in particular when diffractions are concerned.

Spatial wavefield sampling
Diffraction is a wavefield phenomenon and to optimally utilize its potential this fact has to be honored in how we process the data. For wavefields to be accurately reconstructable, dense spatial recordings are needed. It can be argued that the process of diffraction, aside from surface waves which are commonly of lower frequencies, causes the steepest possible wavefront inclination observable in seismic sections. Accordingly, and because of their overall weakness, diffracted wavefields are the first to suffer from unintentional but sometimes unavoidable wavefield under-sampling. As it was appreciated that diffraction is intimately linked to resolution in the reconstructed image, it was concluded that too coarse sensor spacing also leads to lower image resolution (Grasmueck et al., 2005).
This problem of under-sampling becomes arguably most severe, when the imaging of near-surface structures is concerned, where diffraction, due to medium complexity, becomes particularly likely and constitutes a primary carrier of information. Approaching the limiting case of surface-wave radiation, inclination angles generally become higher the shallower the scattering structure is located. Spatial under-sampling, just like its more commonly discussed temporal counterpart, has a characteristic signature in the frequency-wavenumber (f-k) domain. According to simple geometrical considerations, a Nyquist condition for the spatial reconstructability of a signal containing the smallest wavelength λ min can be formulated where ∆x is the trace spacing, ∆x N is the spatial analog of the Nyquist frequency and α max corresponds to the largest inclination angle of the emerging wavefront measured from the horizintal (e.g. Yilmaz, 2001;Grasmueck et al., 2005). As was convincingly argued by Grasmueck et al. (2005) in the context of ground-penetrating radar imaging, equation 10 finds its direct correspondence in a criterion formulated by Claerbout (1985) describing the maximally achievable lateral resolution in a migrated image. Figure 19 shows thinned out versions of a properly sampled diffracted wavefield and the corresponding f-k spectra. Even if only every fifth trace is considered (j=5) aliasing effects can already be observed.
In particular in earthquake seismology, where typical inter-station distances, even in well-instrumented areas, can reach hundreds of kilometers, spatial aliasing can be considered a major reason why wavefield-based migration-type imaging still rarely leads to convincing results. To illustrate that similar problems are encountered in near-surface prospecting, in Figure 20 a series of vertically aligned diffracted events that could correspond to the signature of a fault are displayed. In particular the response of the shallow diffractors is heavily aliased for insufficient sampling. As a result, when diffracted energy is not isolated but accompanied by other coherent contributions, mis-correlations and unphysical coherent signatures can be identified potentially leading to a flawed interpretation.
In analogy to our eyes, migration is an image forming system, which is prone to miscorrelations when spatial sampling criteria are not met (Gray, 2013). This has to do with Seismic diffraction  Figure 19: Diffraction data (top) and the corresponding frequency-wavenumber (f-k) distributions (bottom) for sufficient spatial sampling (j=1) and for different degrees of undersampling (only every j-th trace is considered). Inclination angles are the highest at the diffraction's flanks, which provide the superior illumination and, accordingly, can be viewed as the most precious part of the event.
the fact that arguably all migration types build on the assumption of local coherence. In particular in Kirchhoff migration (equation 1), where seismic energy is repositioned through constructive and destructive interference by summing along diffraction traveltime trajectories, artifacts due to spatial aliasing represent a common problem and extensive research efforts went into suppressing them (e.g. Abma et al., 1999;Biondi, 2001). Just as illustrated in Figure 20, migration aliasing constitutes a major challenge in near-surface imaging, where the migration trajectories, just like the physical diffracted wavefronts themselves, are most inclined. So, wavefield sampling is a topic not only important in the context of diffraction, but of wavefield imaging as a whole.
To demonstrate the challenges diffraction imaging faces especially in the near-surface zone, Figure 21 shows the result of Kirchhoff migration applied to the set of vertically aligned diffractions introduced in Figure 20. While even for the supposedly well-sampled case (j=1), minor artifacts can be observed, they become increasingly more pronounced if sparser data are considered. In line with what our visual inspection of the datasets already suggested, the severe aliasing resulting from only considering every 10-th or 20-th trace causes strong artifact contamination which is especially pronounced in the near-surface region.
The need of large-aperture recordings and the overall weakness of the wavefield means that diffractions often are almost unrecognizable in the migrated images. In addition to an improved signal-to-noise ratio and a natural enhancement of diffracted energy, the use of coherence measures can likewise enhance image quality, when significantly under-sampled data are considered. To illustrate this, Figure 22 shows the migration result for the extremely unfavorable case of j=20, i.e. when only every 20-th trace is included in the migration (com-Seismic diffraction  Figure 20: The response of vertically aligned diffracting structures (as they could be caused by a sub-vertical fault) for the same degrees of data sparsity as in Figure 19. As the shallowest diffraction caused the steepest inclinations, effects of aliasing are particularly pronounced.  Figure 21: Kirchhoff migration results achieved for dense sampling (j=1), and for reduced acquisitions, where only every j-th trace of the well-sampled example is considered. While for j=5, the deeper diffractors, due to lower dips in the respective response remain wellimaged, artifacts due to spatial aliasing compromise the image quality in the shallower region. In the extreme scenario of j=20, it becomes almost impossible to reliably identify the exact scatterer locations. Seismic diffraction  pare Figure 20 and Figure 21). Even though the migration, in particular in the upper half of the image, results in strong artifact contamination and an overall poor reconstruction of the true diffractor locations, the use of the conventional semblance norm helps to clearly delineate the vertically aligned scatterers. As can especially be observed for the uppermost diffractor, the use of the non-linear (N-th root) semblance norm results in the cleanest image and a reliable reconstruction of all scatterer locations. As both, linear and non-linear semblance represent normalized quantities, i.e. only take values between 0 and 1, they can also be used as a weight to suppress artifacts in conventional migration (see the rightmost panel in Figure 22). Despite these improvements and the availability of powerful wavefield reconstruction algorithms (e.g. Zwartjes and Sacchi, 2006;Hennenfent and Herrmann, 2008;Xie and Gajewski, 2017), it can nevertheless be concluded that sufficiently dense spatial sampling should be ensured already in the field for the potential of diffractions to be fully accessible.

Symmetries and correlations
Seismic imaging using diffracted wavefields has gained some momentum in recent years. However, novel dense acquisition geometries using, for example, fiber-optic strain sensing (e.g. Daley et al., 2013) and improved non-invasive separation and extraction techniques promise the practical realization of even more sophisticated strategies to exploit the unique properties of diffractions, including but not limited to velocity inversion. As was emphasized earlier, diffraction is ultimately linked to resolution and diffraction imaging indeed lets one arrive at a highly resolved reconstruction of the Earth's abrupt property changes and imperfections (see, for example, the illustrative comparison of conventional migration and diffraction focusing in Figure 10 and 11). However, as I will seek to demonstrate, even more profound implications are to be expected when the multi-directional radiation of scattering is Seismic diffraction systematically investigated. Ultimately linked to illumination, provided that the assumption of a point scatterer is reasonably justified at least for a subset of the considered wavefield, the recorded prestack diffraction response turns out to be highly symmetric and redundant (e.g. Bauer et al., 2016;Schwarz and Gajewski, 2017a).
To illustrate this symmetry and the redundancy of information in different data subdomains, the problem of common-source (CS) migration is considered. Due to reciprocity, uniform radiation after scattering also implies that no matter from where the energy arrives before the scattering process, the diffraction response, kinematically encoded in its moveout, will always look the same. Consequently, it can be argued that in contrast to reflection imaging the migration of different CS gathers is not of complementary value, but rather (provided the reconstruction is successful) delivers a multitude of copies of the same image. To illustrate this, information encoded in a word formed by individual diffractors leads to a seemingly different response, when sources located at opposite ends of the acquisition line are considered ( Figure 23). However, the application of CS Kirchhoff migration to the individual records reveals that exactly the same structural image with the same Seismic diffraction information content can be reconstructed. For sufficiently repeatable seismic sources this fundamental redundancy of diffraction in CS gathers, owing to reciprocity, likewise applies to the common-receiver (CR) domain.
As a consequence, if diffractions can be reliably detected, the uniform radiation of point scatterers potentially allows to image and successfully characterize subsurface structures with drastically reduced acquisition effort. As the comparison of the reconstructions based on the left and the right CS gathers suggests, a multi-coverage acquisition does not add more information but rather reproduces the diffraction results gained in individual sub-domains. Although two-way moveouts are concerned in the zero-offset domain, this means that also single-or small-offset acquisitions can potentially lead to diffraction images comparable with industry-scale multi-coverage datasets (compare the academic and industrial data examples discussed in this chapter). To further support this claim of redundancy, the migration results of three additional adjacent CS gathers -this time strongly contaminated with uncorrelated noise -are displayed in Figure 24. Again, equivalent reconstructions could be achieved and noise could effectively be suppressed by employing the N-th root semblance norm. To fully Seismic diffraction appreciate this redundancy of information in the prestack diffraction response, the simple case of a single point diffractor embedded in a smooth background medium is considered. Two CS gathers -with the first source located at the far left and the second at the center of the receiver line -are displayed in Figure 25. The trace-wise cross-correlation of the two gathers not only reveals that the overall shape of the diffraction is independent of the source location, but also that the constant time shift is directly encoded in the event's moveout. So, in principle, a single diffraction record is sufficient to deduce the event's shape and position in time at other locations.
We now consider the full prestack diffraction response shown on the left side of Figure 26. If we follow the same strategy and cross-correlate all CS gathers with the first one (indicated by the red frame), we arrive at perfectly flattened gathers for all considered sources and in the common-receiver (CR) domain the estimated time lags again reproduce the event's moveout (shown on the left side in Figure 27). In analogous fashion, the cross-correlation of all CR gathers with a reference gather has the same behavior as was observed in the CS domain (compare Figure 27). This lets us arrive at a simple recipe for the reconstruction of the full prestack response from a single gather (a similar strategy was formulated for passive events by Diekmann, 2018). For conciseness, only the CS domain is considered in the following. The trace-wise cross-correlation X ij of the i-th with the j-th CS gather can be written as Owing to the invariant shape of diffracted events in different configurations, the crosscorrelation only differs from the auto-correlation by the constant time shift τ ij which can be extracted from the event's moveout in any common-source gather via  Figure 28: Prestack response of a two-diffraction experiment (a) and the corresponding by-products of coherence analysis. Also displayed are (b) the detected coherence and (c) the estimated wavefront emergence angle at the receiver. Very similar to an image segmentation task, the local coherence and continuity arguments can help to automatically detect, distinguish and tag individual diffractions (shown in (d), Bauer et al., 2019).
where x i and x j denote the i-th and the j-th receiver locations, respectively (which coincide with the positions of the i-th and the j-th source). The index k indicates that the traveltime differences do not necessarily have to be estimated in the i-th or the j-th gather, but are encoded in any (k-th) gather available. Following from equations 11 and 12, a combination of cross-correlation (denoted by ) and convolution (denoted by * ) allows us to fully reconstruct the observations of the j-th CS gather solely based on the knowledge of the i-th gather. In practice, this means that kinematically, the full prestack diffraction response can be reconstructed from only a single gather. The applicability of equation 13 is demonstrated in Figure 26, where only the first CS gather indicated by the red frame was used to accurately reconstruct the full prestack diffraction response (shown on the right). While the chosen example suggests the overall practicability of equation 13, complications arise if more than one diffraction is recorded in the CS gather (Figure 28(a)). Seismic diffraction  Bauer et al. (2019). The spurious contributions that arise from inter-event correlations and convolutions become more numerous, the more diffracted events are recorded (compare the cross-correlograms (b) with the actual prestack response (a)).
In order for the reconstruction to be successful, spurious arrivals connected to inter-event correlations need to be prevented. Fortunately -in contrast to the image, where diffractors are represented by small-scale perturbations which can be located very close to each other -each individual diffraction, as it is recorded in the time domain, has a unique shape which makes it distinguishable from other contributions.
In line with the observation that the shape of the diffraction response only contains propagation effects and no information on the internal geometry of the scatterer, it was recently suggested that local coherence and smoothness arguments could be used to efficiently distinguish and tag individual diffractions directly in the data domain (Bauer et al., 2019). Alongside the prestack data cube, Figure 28 shows the detected coherence, the estimated emergence angle, and the aforementioned event tag assigned following the strat-Seismic diffraction Auto-correlation (one diffractor) Figure 30: Recorded passive synthetic wavefield and its auto-correlation for (a) a smooth background medium without scattering structures embedded, and (b) when a single diffractor is present. In both cases, the correlogram is dominated by the auto-correlation of the primary event.
egy described by Bauer et al. (2019). Confirming previous observations, the symmetries contained in the diffracted response clearly express themselves in the estimated angle fields (Figure 28(c)). The automated estimation of the event tags allows for a monogamous application of equation 13, which leads to an accurate reconstruction for the multi-diffractor case ( Figure 29). So depending on the objective and the available data configurations, inherent prestack symmetries connected to the uniform radiation of diffraction either justify imaging based on rapidly reduced acquisitions or data volumes, or help to stabilize and additionally inform multi-coverage imaging and inversion. Strongly related to these findings, the next subsection briefly elaborates on how a similar correlation strategy can help to illuminate processes of multiple scattering.  Figure 31: Displayed are (a) the simulated response of two scatterers, (b) the corresponding auto-correlation, as well as the separated contributions of (c) diffractor 1 and (d) diffractor 2. As contributions with high energy typically dominate the correlogram, the primary event was separated before auto-correlation (compare Figure 30). Red lines indicate the time shifts caused by multiple scattering at the left and the right diffractor, respectively.

Deciphering the coda
The coda of seismic events can be a very lively place in which interesting patterns and systematic behavior can be observed. In particular in passive-source seismology, coda waves, owing to the fact that they encode multiple scattering, arrive later and spend significantly more time in the investigated medium than the primary energy of an earthquake. Based on the concept of stationary phase Snieder (2004) demonstrated that interferometric relations can be employed for Green's function retrieval in the seemingly diffuse coda wavefield, thereby providing important evidence on how ambient noise might (partly) be generated. Despite their overall amplitude weakness these multiply scattered contributions, by means of coda wave interferometry, were demonstrated to provide unique sensitivity with respect to structural changes that occur in the course of an earthquake (Snieder et al., 2002). In the following I will briefly discuss, how the auto-correlation of active or passive trace gathers may Seismic diffraction help in recovering inter-scatterer traveltimes from the coda. As was demonstrated before, it is a special property of point diffraction, that radiation occurs uniformly resulting in a somewhat wasteful redundancy of information in prestack records. A similar redundancy can be observed, when multiple scattering is considered. Figures 30, 31, and 32 show passive seismic data that were modeled using the pseudospectral technique by Tessmer (2011), which is capable of honoring multiple scattering. In Figure 30, for comparison, the two simple cases of a smooth background medium without perturbations and a single diffractor are considered. In the respective auto-correlation, as expected, most of the energy is concentrated in the zero-lag regime, whereas non-zero lag times exist owing to the correlation of the scattered with the primary event. The more interesting case of multiple scattering can only be observed, if more than one scatterer is present. As the simplest possible example, the auto-correlation was performed for two scatterers embedded in a smooth background medium (see Figure 31). Whereas the autocorrelated wavefield is considerably richer, perfectly horizontal features with non-zero time lag appear which correspond to the trace-by-trace correlation of different orders of scattering. Very similar to the multi-coverage case, for multiple scattering the remaining diffractor acts like a secondary source leading to the same, merely time-delayed observable moveout.
To be able to better decipher the noticeably more complicated correlogram the primary event was separated before auto-correlation (Figure 31(a) and (b)). Figure 31(c) and Figure  31 show the separate contributions resulting from correlation with the response of the first and the second scatterer, respectively. Two distinct time lags of τ 1 ≈ 0.5 s and τ 2 ≈ 0.3 s can be observed (indicated by the red lines). As can be deduced from Figure 32, where the case of four diffractors is considered, the auto-correlation's structure becomes increasingly complicated, when more scatterers are involved. To better appreciate the implications of these systematic observations in the autocorrelation, a conceptual sketch of the considered multiple scattering process (for two and four diffractors) and its controlled-source analog is displayed in Figure 33. As was recently suggested by other authors, the fact that the moveout connected to a point diffraction is always the same can be used to systematically investigate the order in which observed scattering occurred (Dylan Mikesell et al., 2012;Meles and Curtis, 2014;Löer et al., 2015). Here, I seek to demonstrate that with the help of wavefield separation and tagging, under certain conditions, inter-scatterer traveltimes can be extracted interferometrically from the coda by means of simple auto-correlation. The conceptual sketch in Figure 33 suggests that the correlation of primary and secondary scattering at diffractors A and B yields the respective lag times τ A and τ B . The redundant leg R r (corresponding to the invariant wavefront curvature originating from the last scatterer) is eliminated in the process, Due to reciprocity we have t BA = t AB for the inter-scatterer propagation time which can be written as the arithmetic mean of the observable individual time lags (for diffractor A and diffractor B) In Figure 33 the distances R i , for a constant velocity or an effective replacement medium, correspond to wavefront radii excited by the actual or by secondary sources. Despite its simplicity, equation 15 constitutes an exact interferometric correspondence of inter-scatterer traveltimes and observables, provided point diffraction is considered and the extraction of the time lags τ A and τ B from the auto-correlation was successful. In contrast to crosscorrelation, physically feasible time lags are always positive, which rightfully implies a positive absolute inter-scatterer time t AB . Like before, the presented findings strongly suggest that in many useful ways, diffracting structures can indeed be viewed as passive sources providing structure and order to the seemingly chaotic coda observed in highly scattering media. Provided that these signals, despite their weakness, can be reliably extracted, the knowledge of inter-scatterer times might help in additionally stabilizing and constraining, for example, wavefront tomographic inversion in depth (Bauer et al., 2017). Again, owing to their uniform radiation, the above considerations likewise apply to passive and controlled-source investigations.

CONCLUSIONS AND OUTLOOK
The aim of this chapter was to provide a conceptual and intuitive introduction to the process of seismic diffraction, the resulting wavefield's unique properties and how they can be utilized for imaging and inversion. Historically, the roots of diffraction and the underlying theoretical framework can be traced back to pioneering works in optics, where the concept (illustrated by the famous slit experiments of Young) can be viewed as the ultimate manifestation of a wave process. Although a fully rigorous treatment might only be possible in a full-waveform picture of seismic wave propagation, powerful far-field approximations exist and can successfully be utilized. This introductory text has also sought to convince the reader that a proper understanding and utilization of the process of diffraction demands the treatment of fields, rather than isolated waveforms as is still often the case in the context of sparse acquisitions.
In that light I briefly reviewed how simple wavefield-based imaging techniques such as Kirchhoff migration intrinsically capture and naturally honor diffraction phenomena, resulting in highly resolved images of small-scale subsurface discontinuities. It was likewise discussed, how insufficient spatial wavefield sampling can compromise diffraction imaginga still prevailing challenge in near-surface geophysical prospecting. A variety of controlled synthetic and field data examples suggest that insightful, highly resolved and complementary images can be gained using different data configurations, shedding new light, for example, on past erosional processes and fault systems, which are normally hidden or hard to distinguish in conventional seismic images. Due to the fact that essentially all conventional migration schemes can be used, the successful separation of diffractions from the more dominant reflected wavefield can be considered a major challenge. Following analogies in optics, I briefly described a simple strategy for wavefield separation that can be easily implemented and was employed for all diffraction imaging examples.
Despite the fact that the first studies concerned with seismic diffraction occurring at crustal faults and rugged unconformities date back at least to the 1950s, its practical utilization for imaging and inversion can still be considered a niche discipline. While this in part can be explained by the overall faintness of the diffracted wavefield and the fact that strong interference with reflections often leaves it obscured and largely invisible in the prestack data, it may also be argued that its overall potential is still often underestimated. Through a simple sequence of cross-correlation and convolution it was demonstrated that the full multi-coverage prestack diffraction response is encoded in and can be reconstructed with a single shot record, suggesting a variety of interesting applications.
In particular when academic objectives are concerned, target depths are typically significant, while measurement campaigns are often limited to very small source-receiver distances and, consequently, suffer from drastically reduced illumination, which renders the successful inversion of laterally resolved velocity structure practically impossible. Following from the revealed and discussed symmetries and redundancies, diffracting structures seem ideally suited to act as highly illuminating secondary sources in the subsurface and in addition to Seismic diffraction recent advances in diffraction imaging and inversion, major advances in the processing of reduced data configurations arising from these poor-man acquisitions are expected in the near future.
Simple controlled simulations likewise suggested that there is order in the coda and that multiple scattering and its systematic study by means of auto-correlation can lead to the interferometric extraction of inter-scatterer traveltimes, which could, for example, be utilized for linking individual measurements and efficiently constraining tomographic inversion in depth. Owing to the multi-directional radiation of diffractions, joint method development and the linking of observations in earthquake seismology and controlled-source exploration might eventually enable an integrated view of the Earth's crust (compare, e.g., Zhu, 2018). While previous studies often suffered from poor data quality, the rise of dense seismic acquisition strategies across the scales, as, for example, provided by fibre-optic strain sensing (Daley et al., 2013;Lindsey et al., 2017;Jousset et al., 2018), promise an exciting future for exploring faults, fractures and unconformities with the weak signatures of seismic diffraction.

APPENDIX: KIRCHHOFF'S DIFFRACTION INTEGRAL
Building on the notion of the Huygens-Fresnel principle of interference, Hagedoorn (1954) was the first to develop a systematic strategy to focus and reposition seismic energy across seismograms to arrive at so-called migrated sections. Despite the suceess this systematic treatment had, the rapid increase of the amount of data recorded, fuelled by the digital revolution, a more quantitative and computational approach to migration was in demand. Building on Kirchhoff's diffraction theory in optics, Schneider (1978) introduced an approximate integral approach to the repositioning of seismic energy, called Kirchhoff migration, which still, up to the present day, is routinely employed in everyday practice (e.g. Etgen et al., 2009). Owing to its persistant importance, its formal beauty, and beacause it engraines the process of diffraction at its heart, in the following, I will briefly re-derive the Kirchhoff integral equation 1 from first principles, i.e. as an approximate solution to the wave equation. To maintain comparability, I will largely follow the derivations provided by Schneider (1978), Wiggins (1984) and Born and Wolf (2013). For a more complete account, the reader is referred to these publications. A time-dependent scalar wavefield D, which may, for example, correspond to the pressure, obeys the scalar wave equation where ∇ 2 is the three-dimensional Laplacian, t is time, v is the wave velocity, and Q encodes the spatial and temporal characteristics of the sources contained in the investigated volume. Following up on the ideas of Huygens and Fresnel, Kirchhoff sought a mathematically sound quantitative theory of wave propagation and image reconstruction, in which the process of diffraction, expressed as secondary, elementary wavefields, forms the central ingredient. He showed that the Huygens-Fresnel principle may be explicitly formulated by utilizing Green's theorem, according to which the wavefield within a closed volume can be expressed in terms of the sum of a volume integral encoding the source terms as well as initial conditions and an integral over the surface enclosing this volume (e.g. Born and Wolf, 2013). For the problem of seismic back-scattering, in good approximation, the homogeneous wave equation (Q = 0) can be investigated. Correspondingly, the only non-vanishing contribution to the spatially and temporally varying wavefield constitutes in an integral over the closed surface S D(r 0 , t 0 ) = 1 4π G(r, t | r 0 , t 0 ) ∂D(r, t) ∂n − D(r, t) ∂G(r, t | r 0 , t 0 ) ∂n dS dt .
In equation 17 the time and space coordinates (r, t) ∈ S indicate measurements on the closed surface, whereas (r 0 , t 0 ) correspond to their counterparts within the enclosed volume. Differentiation ∂/∂n is performed with respect to the surface normal. This representation can be viewed as a form of holography, because measurements of D on a closed surface, i.e. a twodimensional manifold allow for the full reconstruction of the wavefield (in time and space) within the enclosed volume (compare, e.g., Maynard et al., 1985). In a slightly different form, it is commonly also referred to as the Helmholtz-Kirchhoff theorem (Born and Wolf, 2013). A two-dimensional sketch of the configuration can be found in Figure 34(a). For consistency and in contrast to other authors, here, the image or location and time of diffraction are denoted by the subscript 0, whereas the spatial and temporal coordinates corresponding to observations on the surface carry no subscript. This should be kept in mind when other works are consulted. In vertical seismic profiling, measurements are mostly confined to the Seismic diffraction (a) (b) Figure 34: Illustration of (a) the holographic representation theorem and (b) the configuration underlying the derivation of the Kirchhoff integral. In principle, closed-surface data acquisition would permit the full and accurate reconstruction of wavefields in the interior. Kirchhoff migration assumes single-sided illumination and is a quantitative manifestation of the Huygens-Fresnel principle, in which any wavefield can be thought of as a superposition of elementary point diffractions (compare Figures 2 and 3).
Earth's surface, which locally, often in good approximation, can be assumed to be planar ( Figure 34(b)). Consequently, in seismic investigations, the exact correspondence expressed by Green's theorem (equation 17) does not hold and approximations are commonly introduced. One of these approximations is the so-called Sommerfeld radiation condition, which requires that the energy radiated by sources must scatter to infinity and that no energy may be radiated back from infinity. In practice, this means that contributions from the distant hemisphere illustrated in Figure 34(b) can be ignored (Schneider, 1978). If this condition is fulfilled, single-sided illumination as is typical in seismic problems, is deemed sufficient for image reconstruction (Claerbout, 1971;Wiggins, 1984;Born and Wolf, 2013). In the following, we assume that the Sommerfeld radiation condition holds and, for convenience, re-parameterize the location vectors in terms of depth with respect to the registration surface, i.e. we have r = (x, z) and r 0 = (x 0 , z 0 ). The two-component vectors x and x 0 describe the lateral position on the surface and within the investigated subsurface volume, respectively. As a result, in practice, the closed surface S in equation 17 is replaced by a finite surface aperture S x and normal derivatives correspond to derivatives in the z direction (compare Figure 34(b)). As normal wavefield derivatives ∂D/∂n are commonly not available, an appropriate choice of the Green's function can suppress these contributions by vanishing on the measurement surface. According to the original derivation by Schneider (1978), the Green's function of free space (placed at x) subtracted by its reflection in the plane meets this requirement. With R = |(x 0 , z 0 ) − (x, z)| and R r = |(x 0 , z 0 ) − (x, −z)| we have where δ is Dirac's delta function. The radii R and R r correspond to the distance between secondary source and the receiver (located on the registration surface S x ) and its reflected counterpart, respectively. Insertion of 18 into the surface integral representation of the wavefield (equation 17), after simplification, then leads to D(x 0 , z 0 , t 0 ) = 1 2π D(x, t) ∂ ∂n The wavefield D(x, t) = D(x, 0, t) now corresponds to the actual densely spaced seismogram recordings acquired in the field. Despite the assumption of an acquisition plane S x Seismic diffraction (Schneider, 1978), this result is formally equivalent to the generalized solution Wiggins (1984) provided for mildly undulating surfaces and represents one form of the Kirchhoff integral. For practical purposes, the time integration in equation 19 can be carried out explicitly. First, however, the normal differentiation, which in the assumed setup is a differentiation in the vertical (z) direction, can be performed. As a result, one arrives at D(x 0 , z 0 , t 0 ) = 1 2π (I 1 + I 2 ) dS x , with These constituent integrals can be solved in a straight-forward manner through integration by parts and by making use of some of the unique properties of Dirac's delta function δ.
After some minor simplifications, the time integration then leads to where similar to optics t = t 0 + R/v is the advanced time, i.e. the time the field will advance from the present time t 0 . Time advancement captures effects of wave propagation through the second term R/v and, therefore, can be viewed as a process of extrapolation. This reflects once more the formal correspondence of diffraction and the excitation of elementary (secondary) sources.
Equation 22 is in perfect correspondence with the classic results found by Schneider (1978) and Wiggins (1984) in the context of seismic migration. It is interesting to note that exactly the same result can be gained, if the generalized Kirchhoff formula, also valid for non-monochromatic light, is considered instead of Green's theorem (see, e.g., Born and Wolf, 2013). A more generic form of the Kirchhoff integral can be gained, if the time integration is not performed explicitly. Then, equation 20 can be written as For a homogeneous background, R/v denotes the one-way diffraction traveltime t diff , which is the traveltime from the diffractor to the considered measurement location on the surface (compare Figure 1). As, owing to reciprocity, this position can be both, a source (s) or a receiver (g) location, the more practical two-way analog reads In the frame of coherence analysis, it was shown by various authors that this so-called double-square-root equation, in the wavefront picture, lends the basis for powerful traveltime approximations that permit improved data enhancement, attribute determination and timedomain imaging, when reflection at highly curved features or diffractions are considered Seismic diffraction (e.g. Landa et al., 2010;Fomel and Kazinnik, 2013;Schwarz et al., 2014). In homogeneous media, the distances R s and R g correspond to the wavefront radii as they are observed at the source and the receiver location and equation 24 provides exact traveltimes (compare equation 8). For the heterogeneous case (see Figure 1), approximate traveltimes can be gained by considering an effective overburden or optical (straight-ray) projections (Schwarz and Gajewski, 2017a,c).
The first term in equation 23 has a 1/t dependence and therefore decays more rapidly when the scattered far-field is concerned (Wiggins, 1984). Exploiting this observation, common implementations of Kirchhoff migration, therefore, only consider the second contribution containing the temporal derivative of the wavefield ∂D/∂t. Following the image principle introduced by Claerbout (1971), setting t 0 = 0 lets one arrive at the so-called Rayleigh-Sommerfeld approximation of the Kirchhoff integral adapted to the problem of seismic migration which is in exact correspondence with equation 1 provided at the beginning of this chapter. Following modern convention, dS x is now replaced by dx and, for convenience and generality, the factor of the second term in equation 23 is abbreviated with w, as it represents a weight that is applied to the temporal derivative of the wavefield during summation. Like in 1, I now denotes the amplitude of the reconstruction at the image point ξ.
Owing to its generic formulation, equation 25 captures Kirchhoff migration in commonsource gathers, common-receiver gathers (i.e. prestack data) or the poststack (zero-offset) domain. To account for heterogeneity, in practice, the weight function w and the diffraction traveltimes t diff are forward-calculated for a given velocity distribution either using analytical expressions in the case of time migration or by employing ray tracing or eikonal solvers in the depth domain. After its initial introduction to the geophysical community, extensions of equation 25, by means of employing more sophisticated weights, are capable to deliver migrated images with fully preserved amplitudes for complicated media (e.g. Bleistein, 1987;Schleicher et al., 1993a).
As was discussed by Schneider (1978), the Kirchhoff integral represents a three-dimensional convolution of the observed wavefield with a space-time operator related to the point-source solution to the scalar wave equation. In that sense, it can indeed be viewed as a quantitative manifestation of the Huygens-Fresnel principle of interference. It is interesting to see that detectable isolated diffractions, at least from a kinematic perspective, share striking similarities with passive sources and Green's functions. Provided the time of scattering and the scattering location can be determined, for example, by means of source localization techniques, they would represent what Dirk Gajewski once called natural Green's functions that could be used for velocity-independent seismic imaging.