ReviewsofGeophysics [622552]

ReviewsofGeophysics
Stationary-phase integrals in the cross correlation
of ambient noise
Lapo Boschi1,2andCornelisWeemstra3
1InstitutdesSciencesdelaTerreParis,SorbonneUniversités,PierreandMarieCurieUniversityParis06,UMR7193,Paris,
France,2InstitutdesSciencesdelaTerreParis,CNRS,UMR7193,Paris,France,3DepartmentofGeoscienceand
Engineering,DelftUniversityofTechnology,Delft,Netherlands
Abstract Thecrosscorrelationofambientsignalallowsseismologiststocollectdataevenintheabsence
ofseismicevents.“Seismicinterferometry”showsthatthecrosscorrelationofsimultaneousrecordings
ofarandomwavefieldmadeattwolocationsisformallyrelatedtotheimpulseresponsebetweenthose
locations.Thisideahasfoundmanyapplicationsinseismology,asagrowingnumberofdenseseismic
networks become available: cross-correlating long seismic records, the Green’s function between
instrumentpairsis“reconstructed”andused,justliketheseismicrecordingofanexplosion,intomography,
monitoring,etc.Theseapplicationshavebeenaccompaniedbytheoreticalinvestigationsoftherelationship
betweennoisecrosscorrelationandtheGreen’sfunction;numerousformulationsof“ambientnoise”theory
haveemerged, eachbased ondifferent hypothesesand/or analyticalapproaches. Thepurpose ofthis
studyistopresentmostofthoseapproachestogether,providingacomprehensiveoverviewofthetheory.
Understanding the specific hypotheses behind each Green’s function recipe is critical to its correct
application.Hopingtoguidenonspecialistswhoapproachambientnoisetheoryforthefirsttime,wetreat
thesimplestformulation(thestationary-phaseapproximationappliedtosmoothunboundedmedia)in
detail.Wethenmoveontomoregeneraltreatments ,illustratingthatthe“stationary-phase”and“reciprocity
theorem”approachesleadtothesameformulaewhenappliedtothesamescenario.Weshowthataformal
cross correlation/Green’s function relationship can be found in complex, bounded media and for
nonuniformsourcedistributions.WefinallyprovidethebasesforunderstandinghowtheGreen’sfunctionis
reconstructedinthepresenceofscatteringobstacles.
1. Introduction
In“seismicinterferometry,”theGreen’sfunction,orimpulseresponse,ofamediumcanbedeterminedempir-
ically, based on the background signal recorded by two instruments over some time. The term “Green’s
function,” ubiquitous in ambient noise literature, indicates the response of a medium to an impulsive exci-
tation: a point source [e.g., Morse and Ingard , 1986;Aki and Richards , 2002]. Measuring a Green’s function
amounts to recording the ground oscillations that follow an explosion: exploiting the “ambient noise,” the
same signal can be measured without setting off any explosive. This approach was foreshadowed by sev-
eral early studies in ocean acoustics [e.g., Eckart, 1953;Cox, 1973] and small-scale seismology [ Aki, 1957;
Claerbout ,1968]andlaterappliedsuccessfullytohelioseismology[ Duvalletal. ,1993;Woodard,1997],ultra-
sound[WeaverandLobkis ,2002;Malcolmetal. ,2004],terrestrialseismology[ CampilloandPaul ,2003;Shapiro
etal.,2005],infrasound[ Haney,2009],andengineering[ Sniederand¸ Safak,2006;Kohleretal. ,2007].
Ambient noise seismology on Earth takes advantage of the “ambient,” low-energy signal that seems to
be generated continuously by the coupling between oceans and solid Earth [e.g., Longuet-Higgins , 1950;
Hasselmann ,1963;Stehlyetal. ,2006;Kedaretal. ,2008;Hillersetal. ,2012;Gualtierietal. ,2013;TraerandGerstoft ,
2014].Itsresolutionisnotlimitedbythenonuniformdistributionofearthquakesorbythedifficultiesinherent
insettingupaman-madeseismicsource.Mostambientenergypertainstosurfacewavesofperiodbetween
∼5s and∼30s, which are difficult to observe in the near field of an earthquake, and almost completely
attenuated in the far field, where only longer-period surface waves can be measured. Ambient noise-based
surface-wave observations in this frequency range are thus complementary to traditional observations and
particularlyusefultomapthecrust,lithosphere,andasthenosphere:theshorterthesurface-waveperiod,the
shallower the depth range. Some recent studies [ Ruigroketal. , 2011;ItoandShiomi , 2012;Polietal., 2012a,
2012b;Gorbatovetal. ,2013;Nishida,2013]showthatnotonlysurfacewavesbutalsobodywavestravelingRESEARCHARTICLE
10.1002/2014RG000455
KeyPoints:
• Ensembleaveragecrosscorrelation
notexactly=Green’sfunction
• Relationship between Green’s
function/noisedatadependson
experimentalsetup
• Weprovideapedagogicalreviewof
ambientnoisetheorywithsimple
examples
Correspondenceto:
L.Boschi,
lapo.boschi@upmc.fr
Citation:
Boschi, L., and C. Weemstra
(2015), Stationary-phase integrals
in the cross correlation of ambient
noise,Rev. Geophys. ,53, 411–451,
doi:10.1002/2014RG000455.
Received3APR2014
Accepted30APR2015
Acceptedarticleonline14MAY2015
Publishedonline23JUN2015
Corrected2OCT2015
Thisarticlewascorrectedon2OCT
2015. Seetheendofthefulltextfor
details.
©2015.AmericanGeophysicalUnion.
AllRightsReserved.
BOSCHIANDWEEMSTRA AMBIENT-NOISECROSSCORRELATION 411

ReviewsofGeophysics 10.1002/2014RG000455
overlargedistancesand/orreflectedbyEarthdiscontinuities(theMoho,theupper-to-lowermantleboundary,
etc.)canbeextractedfromambientsignalcrosscorrelations.
Ambientnoisecorrelationallowsonetobuildaseismicdatasetevenintheabsenceofearthquakes.Inprac-
tice,theGreen’sfunctionisreconstructedonlyapproximatelyandwithinalimitedfrequencyband,butthat
isenoughtoestimaterelevantparameterssuchasgroupandphasevelocityofsurfacewaves.Thisisvaluablefor a number of different applications, allowing passive imaging at the reservoir scale and in the context ofhydrocarbonindustry[ DeRidderandDellinger ,2011;Corciuloetal. ,2012;Weemstraetal. ,2013]aswellasmon-
itoringoftime-dependentchangesinmaterialpropertiesaroundanactivefault[ WeglerandSens-Schonfelder ,
2007],intherigidityofalandslide-pronearea[ Mainsantetal. ,2012]orintheshapeandlocationofmagma
[Brenguieretal. ,2011]andhydrocarbon[ DeRidderetal. ,2014]reservoirs.
Thelastdecadehasseenthepublicationofanumberofincreasinglyexhaustivemathematicaldescriptions
of the phenomenon of Green’s function reconstruction from ambient recordings [ LobkisandWeaver , 2001;
Derodeetal. ,2003;Snieder,2004;Wapenaar ,2004;WeaverandLobkis ,2004;Rouxetal. ,2005;Wapenaaretal. ,
2005;Nakahara ,2006;Sanchez-SesmaandCampillo ,2006;Wapenaaretal. ,2006;Tsai,2009;Weaveretal. ,2009].
Ingeneral,amathematicalrelationshipisfoundbetw eentheGreen’sfunctionassociatedwiththelocations
oftworeceivers(i.e.,theresponse,observedatoneofthereceivers,toapointsourcedeployedattheotherreceiver)and thecross correlation[e.g., Smith, 2011], computed over a long time interval, between the ran-
domambientsignalrecordedbythereceivers.(Adetailedaccountofhowambientdataaretreatedisgiven
byBensenetal. [2007].)Differentapproachestoestablishingthisrelationshiphavebeenfollowed,however,
resultingindifferentformulations.Mostauthorsfirstdevelopthetheoryforthesimplecaseofacoustic(scalar)wave propagation in two- or three-dimensional media. Some formulations hold for heterogeneous and/or
boundedmedia,whileothersarelimitedtoinfinitehomogeneousmedia.Definingambientnoise(a“diffuse
field”)mathematicallyisanontrivialproblemperse:someauthorsdescribeitasthesuperpositionofplanewavespropagatinginalldirections;othersprefertosuperposeimpulsiveresponses(Green’sfunctions),some-timesintwodimensions,sometimesinthree;othersyetdefinethenoisefieldasonewhereallnormalmodes
have the same probability of being excited. Finally, the expression Green’s function itself is ambiguous: in
elastodynamics,itmayrefertotheresponseofamediumtoanimpulsiveforce,impulsivestress,orimpulsiveinitial conditions in displacement, or velocity. While specialists are well aware of these differences and theirimplications, interferometry “users” (including the authors of this article) are sometimes confused as to the
theoreticalbasis,andpracticalreliability,ofthemethodstheyapply.
Thegoalofthisstudyistoshowindetailhowseveraldifferentderivationsof“ambientnoisetheory”leadto
apparently different but indeed perfectly coherent results. We first treat some particularly simple scenarios:
ambientnoiseinahomogeneouslossless(i.e.,nonattenuating)2-Dor3-Dmedium,generatedbyazimuthallyuniform,2-Dor3-Ddistributionsofpointsources,orbyplanewavespropagatinginalldirections.Theanal-ogy between acoustic waves in two dimensions and Rayleigh waves is discussed. Analytical expressions for
thecrosscorrelationofambientsignalsinallthesephysicalsettingsaregiveninsection3andimplemented
numericallyinsection4.Inallthesimplecasesweconsider,thecrosscorrelationofambientsignalasderivedinsection3coincideswiththeintegral,overtheareaoccupiedbysources,ofanoverallveryoscillatoryfunc-tion, which becomes slowly varying only around a small set of so-called “stationary” points. We show in
section5howsuchintegralsaresolvedviatheapproximate“stationary-phase”method(AppendixA);differ-
entanalyticalrelationshipsbetweencrosscorrelationandGreen’sfunctionarethusfoundanddiscussed;theyarelatersummarizedinsection9.Insection8weapply,again,thestationary-phasemethodtoamorecom-plicatedmediumincludingonescatteringobstacle,whichallowsustointroducetheconcept,oftenfoundin
ambientnoiseliterature,of“spurious”arrivalincrosscorrelation.Finally,theso-called“reciprocitytheorem”
approach provides an analytical relationship between Green’s function and cross correlation that is valid inarbitrarilycomplex,attenuatingmedia.Wedescribeitindetailinsection6,togetherwithafewother“alterna-tive”approachestoambientnoisetheory;theconsistencybetweenreciprocitytheoremandstationary-phase
resultsisverified.Toavoidambiguityasmuchaspossible,weprovideanoverviewoftheunderlyingtheoryfor
bothacousticsandelastodynamicsinsection2andaderivationofacousticGreen’sfunctionsinAppendixE.Bycollectingallthispreviouslyscatteredmaterialinasinglereview,wehopetoprovideausefultoolforgrad-uatestudentsandnonspecialistsapproachingthetheoryofacousticandseismicambientnoiseforthefirst
time. Most of the results presented here are limited to nonscattering, nondissipative, homogeneous acous-
ticmediaandtosurfacewave(andnotbodywave)propagationinelasticmedia.Thesesimplificationshavetheadvantageofallowingaself-contained,relativelyunclutteredderivation.Weexploremorerealisticsetups
BOSCHIANDWEEMSTRA AMBIENT-NOISECROSSCORRELATION 412

ReviewsofGeophysics 10.1002/2014RG000455
(nonuniform source distributions; scattering) in sections 7 and 8. The thorough understanding of ambient
noisecrosscorrelationinsimplemediawillserveasasolidplatformformoreadvancedinvestigations.Read-ersthataremostlyinterestedinapplicationstorealisticenvironmentsarereferredtothedifferentlyminded
reviewsof Curtisetal. [2006],Laroseetal. [2006],Gouedardetal. [2008],Wapenaaretal. [2010a,2010b], Snieder
andLarose [2013],Ritzwoller [2014],and CampilloandRoux [2014].
2. GoverningEquations
In this study, analytical relationships between Green’s functions and the cross correlation of ambient signal
arederivedinanumberofdifferentscenarios:sphericalacousticwavesintwoandthreedimensions,Rayleigh
waves,andplanewaves.Wefirstsummarizethetheoryunderlyingeachofthesecases.
2.1. AcousticWavesFromaPointSourceinFreeSpace
Pressure pin homogeneous, inviscid fluids occupying a three-dimensional space obeys the linear, lossless
waveequation
∇2p−1
c2𝜕2p
𝜕t2=𝜕q
𝜕t, (1)
where∇is the gradient operator, c2the ratio of adiabatic bulk modulus to density, tdenotes time, and the
forcingterm qisthe(apparent)massproductionperunitvolumeperunittime,representing,e.g.,anexplosion
oraloudspeaker[e.g., Kinsleretal. ,1999,chap.5].
If𝜕q∕𝜕t=𝛿(r−r0)𝛿(t−t0),with𝛿theDiracdistribution, rthepositionvector,and r0,t0thelocationandtimeof
animpulsivesoundsource,respectively,equation(1)isreferredtoas“Green’sproblem”anditssolutionfor p
asGreen’sfunction.OncetheGreen’sproblemissolved,theresponseofthesystemtohowevercomplicatedasourceisfoundbyconvolvingtheGreen’sandsourcefunctions[e.g., MorseandIngard ,1986;AkiandRichards ,
2002].TheGreen’sfunction
G3Dassociatedwith(1)isderivedinAppendixEhere,equations(E21)and(E22),
workinginsphericalcoordinatesandchoosing(withoutlossofgenerality) t0=0andr0=0.
2.2. MembraneWavesFromaPointSource
Equation(1)alsodescribesthemotionofanelasticmembrane.Inthiscase, pcanbeinterpretedasdisplace-
ment in the vertical direction (for a horizontal membrane), ∇is the surface gradient, and cthe ratio of the
membranetensionperunitlengthtoitssurfacedensity[e.g., Kinsleretal. ,1999,section4.2].Theassociated
Green’sfunction G2DisdeterminedanalyticallyinAppendixE,equations(E15)and(E16).
2.3. RayleighWavesFromaPointSource
The “potential” representation of surface wave propagation, adopted, e.g., by Tanimoto [1990], consists of
writingRayleighwavedisplacement
uR=U(z)̂zΨR(x,y;t)+V(z)∇1ΨR(x,y;t), (2)
wherex,y,andzareCartesiancoordinates( zdenotesdepth),the“verticaleigenfunctions” UandVdepend
onzonly,andtheRayleighwavepotential ΨRvariesonlylaterally(butvarieswithtime). ̂x,̂y,and̂zareunit
vectors in the directions of the corresponding coordinates, and the operator ∇1=̂x𝜕
𝜕x+̂y𝜕
𝜕y.L o v ew a v e
displacementisaccordinglywritten
uL=−W(z)̂z×∇1ΨL(x,y;t), (3)
with×denotingavectorproduct.
Substitutingthe Ansätze(2)and(3)intothe3-Dequationofmotion,onefindsthatbothpotentials ΨRandΨL
satisfy
∇2
1ΨR,L−1
c2
R,L𝜕2ΨR,L
𝜕t2=TR,L, (4)
i.e.,equation(1),with cR,Ldenotingphasevelocityand TR,TLscalarforcingterms[e.g., Tanimoto ,1990;Tromp
andDahlen ,1993;Udías,1999].
BOSCHIANDWEEMSTRA AMBIENT-NOISECROSSCORRELATION 413

ReviewsofGeophysics 10.1002/2014RG000455
Figure 1. Sphericalreferenceframeusedthroughoutthis
study.Theoriginisplacedatthelocationofreceiver1( R1),
whilereceiver2( R2)liesonthe 𝜑=0,𝜃=𝜋∕2axis.The
distancebetweentheoriginandapoint Pisdenoted r.The Green’s problem associated with (4) coin-
cides with the 2-D, membrane wave problemof section 2.2 and Appendix E, and the Green’sfunction corresponding to
ΨR,Lis, again, G2Dof
equations (E15) and (E16). The surface waveGreen’s function can be obtained from the scalarGreen’s function
G2Dvia equation (2) [ Trompand
Dahlen,1993,section5].
Importantly,accordingto(2),theverticalcompo-
nent of Rayleigh wave displacement uRis simply
U(z)ΨR(x,y;t);U(z)isafunctionof zthatdoesnot
changewithtimesothatthe phaseofthevertical
componentof uRcoincideswiththatofthepoten-
tialΨR(x,y;t). The analysis of 2-D ambient noise
thatweconducthereholdsthereforeforRayleigh
wavespropagatingina3-Dmedium(intheEarth)andmeasuredontheverticalcomponentofseismograms,asfarasthephaseisconcerned;amplitudeneedstobescaled.
TheRayleigh-wavederivationof Snieder[2004]isalsobasedonthisidea,althoughinthatstudythesepara-
tion between
x,y,tdependence and zdependence follows from a normal-mode formulation [e.g., Snieder,
1986].HallidayandCurtis [2008]exploretheeffectsofsubsurfacesources,orlackthereof,onGreen’sfunction
reconstruction;theirresultsareconsistentwithoursandwiththoseof Snieder[2004]whenevertheirsetupis
thesameasours,i.e.,uniformsourcedistributionovertheEarth’s(laterallyinvariant)surface.
2.4. PlaneWaves
Someauthors[e.g., Aki,1957;Sanchez-SesmaandCampillo ,2006;Tsai,2009;Boschietal. ,2013]havestudied
diffusewavefieldswithcylindricalsymmetry,whichtheydescribedasthesuperpositionofplanewaves(ratherthansphericalorcircularasseensofarinthisarticle)travelingalongmanydifferentazimuths.Ifitsright-handside
𝜕q
𝜕t=0,equation(1)issatisfiedbyamonochromaticplanewave
p(r,t)=S(𝜔)cos(𝜔t+k⋅r), (5)
wherethevector kisconstantanddefinesthedirectionofpropagation,andtheamplitudefunction S(𝜔)is
arbitrary. A plane wave approximates well the response of a medium to a real source, at receivers that aresufficientlyfarfromthesource.
3. CrossCorrelationinaDiffuseWavefield
Weareinterestedinthecrosscorrelationbetweenrecordingsofa diffusewavefieldmadeatapairofreceivers.
Inacoustics,adiffusewavefieldissuchthattheenergyassociatedwithpropagatingwavesisthesameatall
azimuthsofpropagation[e.g., Kinsleretal. ,1999,section12.1].(Theexpression“equipartitionedfield,”often
intended as a synonym of “diffuse,” is used in the literature with slightly different meanings depending onthe context [ Snieder et al. , 2010], and we chose not to employ it here to avoid ambiguity.) While ambient
noiserecordedonEarthisnotstrictlydiffuse[ Mulargia,2012],diffusefieldtheorysuccessfullydescribesmany
seismicobservations,andassuchitisatleastausefulfirstapproximationofrealambientnoise.
Wesimulateanapproximatelydiffusewavefieldbyaveraging(or,inseismologyjargon,“stacking”)crosscor-
relations associated with sources distributed uniformly over a circle or a sphere surrounding the receivers.This is equivalent to the source azimuth being random; i.e., over time, all azimuths are sampled with equalfrequency/probability.
We center our spherical reference frame at the location of “receiver 1” (
R1) and orient it so that “receiver 2”
(R2) lies on the 𝜑=0,𝜃=𝜋∕2axis (Figure 1). We choose sources to be separated in time, that is to say, a
receiverneverrecordssignalfrommorethanonesourceatthesametime:whilediffusewavefieldsintherealworld might result from multiple, simultaneous sources, our simplification is justified by the mathematical
findingthattheso-called“crossterms”,i.e.,thereceiver-receivercrosscorrelationofsignalgeneratedbydif-
ferentsources,arenegligiblewhen“ensembleaveraged.”Aproofisgiven,e.g.,by Weemstraetal. [2014]and
issummarizedhereinAppendixD.
BOSCHIANDWEEMSTRA AMBIENT-NOISECROSSCORRELATION 414

ReviewsofGeophysics 10.1002/2014RG000455
Let us call Sa source location and riSthe distance between Sand theith receiver. If the source at Sis an
impulsive point source, the recorded signals coincide with the Green’s function associated with the source
locationS,evaluatedat R1andR2.In3-D,theGreen’sfunctionisgivenbyequation(E22),andinthefrequency
domainthecrosscorrelationreads
G3D(r1S,𝜔)G∗
3D(r2S,𝜔)=1
(4𝜋c)21
2𝜋ei𝜔
c(r2S−r1S)
r1Sr2S(6)
(here and in the following the superscript ∗denotes complex conjugation). Since the dimension of G3Din
the frequency domain is time over squared distance (Appendix E), that of equation (6) is squared time overdistancetothepowerof4(ortimeoverdistancetothepowerof4inthetimedomain).
In2-D(membranewavesandRayleighwaves),
G2Disgivenby(E16),andthecrosscorrelation
G2D(r1S,𝜔)G∗
2D(r2S,𝜔)=−i
32𝜋c4H(2)
0(𝜔r1S
c)[
−iH(2)
0(𝜔r2S
c)]∗
. (7)
The cumulative effect of multiple sources Sis obtained by rewriting expression (6) or (7) for each S, withr1S
andr2Svaryingasfunctionsofthedistance randazimuth 𝜑(and,in3-D,inclination 𝜃)ofSwithrespecttothe
originandsumming.
Weshallconsiderseveraldifferentscenarios:
(i)In3-Dspace,foracontinuousdistributionofsourcesalongacircleinthe 𝜃=𝜋∕2plane(centered,forthe
sakeofsimplicity,at R1)summingthecrosscorrelationsover Sleadstoanintegralover 𝜑,
IC(𝜔)=1
(4𝜋c)21
4𝜋2∫𝜋
−𝜋d𝜑nC(𝜑)ei𝜔
c(r2S(𝜑)−r1S(𝜑))
r1S(𝜑)r2S(𝜑), (8)
withnCthenumberofsourcesperunitazimuthorsourcedensity.
(ii)Ifsourcesaredistributedonaspherecenteredat R1,thedoubleintegralisoverthesurfaceofthesphere
(4𝜋solidradians),and
IS(𝜔)=1
4𝜋∫𝜋
0d𝜃sin𝜃∫𝜋
−𝜋d𝜑nS(𝜃,𝜑)G3D(r1S,𝜔)G∗
3D(r2S,𝜔)
=1
4𝜋1
(4𝜋c)21
2𝜋∫𝜋
0d𝜃sin𝜃∫𝜋
−𝜋d𝜑nS(𝜃,𝜑)ei𝜔
c(r2S(𝜃,𝜑)−r1S(𝜃,𝜑))
r1S(𝜃,𝜑)r2S(𝜃,𝜑),(9)
wherenSdenotesthenumberofsourcesperunitofsolidangleonthesphere.
(iii)Switchingto2-D,i.e.,toRayleighwavesorelasticwavespropagatingonamembrane,weneedtointegrate
expression(7)overthepositionoccupiedby S.Assumingasourcedistributionanalogoustothatofcase
(i),withsourcesalongacircle,andsourcedensity nMafunctionof 𝜑only,wefindaftersomealgebrathat
crosscorrelationisdescribedby
IMW(𝜔)=1
32𝜋c41
2𝜋∫𝜋
−𝜋d𝜑nM(𝜑)H(2)
0(𝜔r1S(𝜑)
c)
H(2)
0∗(𝜔r2S(𝜑)
c)
≈1
16𝜋2c31
2𝜋∫𝜋
−𝜋d𝜑nM(𝜑)ei𝜔
c(r1S(𝜑)−r2S(𝜑))
𝜔√
r1S(𝜑)r2S(𝜑),(10)
where we have replaced the Hankel function H(2)
0with its asymptotic (high-frequency and/or far-field)
approximation,equation(9.2.4)of AbramowitzandStegun [1964];thisapproximationisnecessarytolater
solve the integral in (10) via the stationary-phase method. Notice that we do not yet require the wave-
field to be isotropic: at this point nC,nS, andnMare arbitrary functions of 𝜃and𝜑. We will show in
section 5, however, that they need to be smooth or constant for the stationary-phase approximation tobeapplicable.
(iv)Intheplanewaveapproachofsection2.4nosourcelocationisspecified,butonecombinesplanewaves
(5)travelingalongallazimuths
𝜑.AtR1,themonochromaticplanewaveoffrequency 𝜔0travelinginthe
direction 𝜑is
p(R1,t)=S(𝜔0)cos(𝜔0t); (11)
atR2thesamesignalisrecordedatadifferenttime,and
p(R2,t)=S(𝜔0)cos[
𝜔0(
t+Δcos𝜑
c)]
. (12)
BOSCHIANDWEEMSTRA AMBIENT-NOISECROSSCORRELATION 415

ReviewsofGeophysics 10.1002/2014RG000455
Based on equations (11) and (12), Boschi et al. [2013] show that the source-averaged, two-station cross
correlationofmonochromaticplanewavesatfrequency 𝜔0canbewritten
IPW(𝜔0,t)=q(𝜔0)
2𝜋∫𝜋
0d𝜑cos[
𝜔0(
t+Δcos𝜑
c)]
, (13)
obtained from equations (18) and (21) of Boschi et al. [2013] through the change of variable
td=Δcos(𝜑)∕cand assuming source density and the amplitude term q(𝜔)(reflecting the frequency
spectrumofthesource)tobeconstantwithrespectto 𝜑.TheFouriertransformof(13)is
IPW(𝜔0,𝜔)=q(𝜔0)
√
8𝜋∫𝜋
0d𝜑[𝛿(𝜔+𝜔0)+𝛿(𝜔−𝜔0)]ei𝜔Δcos𝜑
c. (14)
Equation(14)isvalidforany 𝜔0;usingthepropertiesoftheDirac 𝛿function,wegeneralizeitto
IPW(𝜔)=q(𝜔)√
2𝜋∫𝜋
0d𝜑cos(
𝜔Δcos𝜑
c)
. (15)
Weintegrateequations(8)–(10)and(15)numericallyinsection4;thisexercisewillservetovalidatetheresults,
illustratedinsection5,ofintegratingthe sameequationsanalyticallyviathestationary-phaseapproximation.
Thecaseofaspatially(andnotjustazimuthally)uniformdistributionofsourcesoveraplaneorthe3-Dspace
isnottreatedhereintheinterestofbrevity.Itessentiallyrequiresthatexpression(8),(9),or(10)beadditionally
integrated over source location; in practice, one needs to integrate over the distancebetween source and
origin,alongtheazimuthsofstationarypointsonly(one-dimensionalintegrals)[ Snieder,2004;Sato,2010].
4. NumericalIntegration
Before deriving approximate analytical solutions for the integrals in equations (8)–(10) and (15), we imple-
menttheseequationsnumerically.Inpractice,weevaluateequations(E16)and(E22)at R1,R2,todetermine
G2DorG3Dnumerically for discrete sets of sources and values of 𝜔. For each source, we cross correlate the
Green’s function calculated at R1with that calculated at R2. We stack the resulting cross correlations in the
timedomain,whichisequivalenttocalculatingthe 𝜑and𝜃integralsintheequationsabove.Thestackscan
becomparedwithpredictionsbasedontheanalyticalformulaethatweshallillustratebelow.
4.1. RingofSources
Thesetupassociatedwithequation(8)isimplementednumerically,placing720equallyspacedsourcesalong
the planar ring centered at R1.R1andR2are 20 km away from one another, and wave speed is 2 km/s. We
implement equation (8) directly at frequencies between ±10 Hz, with sampling rate 20 Hz; we taper the
highest and lowest frequencies to avoid ringing artifacts and inverse Fourier transform via numerical fast
Fourier transform [ CooleyandTukey , 1965]. The corresponding cross correlations are shown in the “gather”
plotofFigure2a.Sourcesalignedwiththetworeceivers(azimuth 0∘or180∘)resultinthelongestdelaytime
between arrival of the impulse at R1andR2: cross correlation is nonzero near ±10s, which corresponds to
the propagation time between R1andR2. Sources at azimuth around ∼90∘and∼270∘are approximately
equidistant from R1andR2, resulting in the impulse hitting R1andR2simultaneously, and cross correla-
tion in Figure 2a being nonzero at ∼0s for those source azimuths. The result of averaging over all sources
is shown in Figure 2b. It is clear that the imaginary part of the time domain signal is 0, as it should be. As
for the real part, cross correlations associated with sources at azimuth 0∘or180∘add up constructively;
the cumulative contribution of other sources is smaller but nonnegligible. The source-averaged cross cor-
relation is accordingly dominated by two peaks corr esponding to energy traveling in a straight path from
R2toR1(generated by sources at azimuth 0∘) and from R1toR2(sources at 180∘). The two peaks are
often labeled “causal” and (not quite appropriately) “anticausal” (or “acausal”), respectively. In this particu-
lar case, the causal and anticausal peaks have different amplitude. This asymmetry reflects the asymmetry
in the locations of R1andR2with respect to the source distribution: R2is closer to sources at 𝜑=0∘
thanR1is to sources at 𝜑=180∘, and waves hitting R2have accordingly larger amplitude. The descrip-
tion of seismic/acoustic interferometry in terms of a stationary-phase integral (section 5.1 and Appendix A)explains mathematically why sources at
0∘or180∘are most relevant and provides an analytical relation-
ship between the source-averaged cross correlation in Figure 2b and the Green’s function of the medium.
BOSCHIANDWEEMSTRA AMBIENT-NOISECROSSCORRELATION 416

ReviewsofGeophysics 10.1002/2014RG000455
a
b
Figure 2. Crosscorrelationsassociatedwithacircular,planar
distributionofsourcessurroundingthetworeceivers. R1isat the
centerofacircleofsourceswithradius R=100km.R2liesonthe
planedefinedby R1andthesources, Δ=20km away from R1.P hase
velocityc=2km/s.Sourcedensity nC(𝜑)=2sourcesperdegree,
independent of 𝜑.G3D(r1S,𝜔)andG3D(r2S,𝜔)areevaluated
numericallyandcrosscorrelatedforaringofequallyspacedsources.TheresultisinverseFouriertransformedandaveragedoverallsources.(a)Time-domainsingle-sourcecrosscorrelationsforallazimuths.(b)Real(solidline)andimaginary(dashed) partsofthesource-averaged,timedomaincrosscorrelation,resultingfromstackingthetracesinFigure2a.Thedimension isthat of
G3D(r,𝜔),squaredand
integratedoverfrequency(inverseFouriertransformed).4.2. SourcesonaSphere
We next distribute ∼103approximately
equally spaced sources over the surface
of a sphere of radius R=100km, cen-
tered at R1.R2lies 20 km away from
R1. In analogy with section 4.1, the sig-
nal recorded at R1and that recorded at
R2(equation (E22)) are cross correlated
foreachsource,andcrosscorrelationsare
stacked, implementing equation (9). The
source-averagedfrequency-domaincross
correlationiscomputedbetween ±10Hz,
with sampling rate 20 Hz, and taperedat high frequency to avoid artifacts. The
time-domain result is shown in Figure 3;
in this case, no gather plot was made
because of the difficulty of visualization
when sources are distributed in threedimensions. The source-averaged cross
correlationisreal,asitshouldbe;itisdif-
ferentfromthatofFigure2b,inthatcrosscorrelationsassociatedwithsourcesatall
azimuths add up constructively, result-
ing in a boxcar function. Again, the max-
imum delay time between arrival of the
impulse at
R1andR2(10 s in our setup)
corresponds to sources at 0∘and180∘,
and no energy is recorded at t>10s;
cross correlations are accordingly 0 at
t<−10sandt>10s.Despitetheasym-
metryinreceiverlocationwithrespectto
the sources, the stacked cross correlation is now symmetric. All these results are reproduced analytically in
section5.2.
4.3. MembraneWaves/RayleighWavePotential
Thesetupassociatedwithequation(10)isalsoreproducednumerically,placing720sourcesat 0.5∘azimuth
𝜑intervalsalongtheringcenteredat R1.Frequencyvariesbetween ±10Hzwithsamplingrate=20Hz;high
Figure 3. Real(solidline)andimaginary(dashed)partsoftheaverage
timedomaincrosscorrelation,resultingfromstackingallcrosscorrelationsassociatedwithauniformlydensedistributionofsourcesoverasphere.
R1isatthecenterofasphereofsourceswithradius
R=100km.R2liesΔ=20km away from R1.P hasev elocity c=2km/s.
Thecalculationisconductednumerically(section4.2),byevaluating
G3D(r1S,𝜔)andG3D(r2S,𝜔)foreachsourceviaequation(E22)and
subsequentlyapplyingcrosscorrelation,inverseFouriertransformation,andaveragingoverallsources.frequencies are tapered as above. As
we shall show analytically in section 5.3,
equation(10)leadstoanonconvergent 𝜔
integral when a time domain expression
forIMWis sought; we therefore compute
the time derivative of IMW(multiplication
byi𝜔in the frequency domain) before
inverseFouriertransforming.Intheinter-
est of speed, we implement the asymp-
toticapproximationof H(2)
0(i.e.,thesecond
lineofequation(10))ratherthan H(2)
0itself.
The results are shown in Figure 4. Thegather (Figure 4a) is qualitatively simi-
lar to that of Figure 2a, but Figure 4b
shows that cross correlations associatedwith sources away from azimuth
0∘and
180∘canceloutwhenstacking.Afterstac-
king (Figure4b) we verify that the time
BOSCHIANDWEEMSTRA AMBIENT-NOISECROSSCORRELATION 417

ReviewsofGeophysics 10.1002/2014RG000455
a
b
Figure 4. Crosscorrelations,differentiated withrespecttotime,
associatedwithacircular,planardistributionofsourcessurroundingthetworeceiversintwodimensions, i.e.,onamembrane. ThesetupisthesameasinF igur e2.
G2D(r1S,𝜔)andG2D(r2S,𝜔)areevaluated
numericallyandcrosscorrelatedforaringof720equallyspacedsources.TheresultisinverseFouriertransformedandaveragedoverallsources.(a)Timedomainsingle-sourcecrosscorrelationsforallazimuths.(b)Real(solidline)andimaginary(dashed)partsofthestackedtimederivativeofthecrosscorrelation,resultingfromstackingthetracesinFigure4a.derivative of IMWis real. We also find that
it is antisymmetric with respect to time,
indicatingthat IMWissymmetric,inagree-
ment with the results of section 4.2 butnotwithsection4.1.
5. AnalyticalFormulaeforTwo
ReceiverCrossCorrelationsofaDiffuseWavefield
As pointed out by Snieder[2004], inte-
gration with respect to 𝜃and𝜑in
equations (8)–(10) and (15) can also be
conducted analytically by means of thestationary-phase approximation, thedetailsofwhicharegiveninAppendixA.
The numerical results of section 4 will
serve as a reference to validate theapproximateresultspresentedbelow.
5.1. RingofSourcesinFreeSpace
The integrand at the right-hand side
of equation (8) coincides with that in
(A1), after replacing
𝜆=𝜔,x=𝜑and
defining f(𝜑)=nC(𝜑)∕[4𝜋2r1S(𝜑)r2S(𝜑)]
and𝜓(𝜑)=[r1S(𝜑)−r2S(𝜑)]∕c. We shall
first identify the values of 𝜑such that
𝜓′(𝜑)=0(stationary points), then use
equation(A2)toevaluatethecontributionofeachsta tionarypointtotheintegral,andfinallycombinethem.
Importantly, the stationary-phase approximation is valid at high frequencies 𝜔−→∞(corresponding to
𝜆−→∞insectionA1).Choosingthereferenceframeasdescribedinsection3(Figure1),bydefinition r1Sis
constant(weshallcallit R)and,basedonsomesimplegeometricalconsiderations,
r2S=[(Rcos𝜑−Δ)2+R2sin2𝜑]1
2=(R2+Δ2−2ΔRcos𝜑)1
2, (16)
whereΔisinterstationdistanceand Δ<Rinoursetup.Itfollowsthat
f(𝜑)=1
(4𝜋c)2nC(𝜑)
4𝜋2R(R2+Δ2−2ΔRcos𝜑)1
2, (17)
𝜓(𝜑)=1
c(r1S−r2S)=1
c[(
R2+Δ2−2ΔRcos𝜑)1
2−R]
, (18)
andupondifferentiating 𝜓withrespectto 𝜑,
𝜓′(𝜑)=ΔRsin𝜑
c(R2+Δ2−2ΔRcos𝜑)1
2, (19)
𝜓′′(𝜑)=ΔRcos𝜑
c(R2+Δ2−2ΔRcos𝜑)1
2−Δ2R2sin2𝜑
c(R2+Δ2−2ΔRcos𝜑)3
2. (20)
We infer from (19) that the stationary points of 𝜓(𝜑)within the domain of integration are 𝜑=−𝜋,0,𝜋
(correspondingto sin𝜑=0).Following BenderandOrszag [1978],werewriteequation(8)asasumofintegrals
limitedtothevicinityofstationarypointsandwitha stationarypointasoneoftheintegrationlimits:
IC(𝜔)≈∫−𝜋+𝜀
−𝜋f(𝜑)ei𝜔𝜓(𝜑)d𝜑+∫0
−𝜀f(𝜑)ei𝜔𝜓(𝜑)d𝜑+∫𝜀
0f(𝜑)ei𝜔𝜓(𝜑)d𝜑+∫𝜋
𝜋−𝜀f(𝜑)ei𝜔𝜓(𝜑)d𝜑,(21)
BOSCHIANDWEEMSTRA AMBIENT-NOISECROSSCORRELATION 418

ReviewsofGeophysics 10.1002/2014RG000455
whichisvalidfor 𝜔−→∞andarbitrarilysmall 𝜀.Equation(A2)cannowbeuseddirectlytointegrateeachof
thetermsattheright-handsideof(21),andafternoticingthatboth f(𝜑)and𝜓(𝜑)aresymmetricwithrespect
to𝜑=0,wefind
IC(𝜔)≈2f(𝜋)ei(
𝜔𝜓(𝜋)±𝜋
4)√𝜋
2𝜔|𝜓′′(𝜋)|+2f(0)ei(
𝜔𝜓(0)±𝜋
4)√𝜋
2𝜔|𝜓′′(0)|. (22)
Wenowneedtoevaluate f,𝜓,and𝜓′′at0and𝜋.Itfollowsfromthedefinitionof f,fromequations(18)and
(20)andfromthefactthat R>Δ(andhence|R−Δ|=R−Δ)that
f(0)=1
(4𝜋c)2nC(0)
4𝜋2R(R−Δ), (23)
f(±𝜋)=1
(4𝜋c)2nC(±𝜋)
4𝜋2R(R+Δ), (24)
𝜓(0)=−Δ
c, (25)
𝜓(±𝜋)=Δ
c, (26)
𝜓′′(0)=ΔR
c(R−Δ), (27)
𝜓′′(±𝜋)=−ΔR
c(R+Δ). (28)
To obtain equations (25)–(28), one must also recall that√
R2+Δ2−ΔRcos𝜑=r2Sis, physically, a positive
distance:when 𝜑=0,±𝜋,itfollowsthat√
R2+Δ2−2ΔR=R−Δ.
Substitutingexpressions(23)–(28)intoequation(22),wefindaftersomealgebrathat
IC(𝜔)≈1
(4𝜋c)21
2𝜋√
c
2𝜋ΔR3⎡
⎢
⎢⎣nC(𝜋)ei(𝜔Δ
c−𝜋
4)
√
(R+Δ)𝜔+nC(0)e−i(𝜔Δ
c−𝜋
4)
√
(R−Δ)𝜔⎤
⎥
⎥⎦, (29)
wherethesignof 𝜋∕4intheargumentoftheexponentialfunctionwasselectedbasedonthesignof 𝜓′′as
explainedinsectionA1.Comparingequation(29)with(E17),itisapparentthatbothtermsattheright-handside of (29) are proportional to the high-frequency/far-field form of the Green’s function
G2D(Δ,𝜔)(and,
interestingly,not G3D).
Equation (21) and the subsequent expressions for IC(𝜔)are only valid for large and positive 𝜔so thatIC(𝜔)
remainsundefinedfor 𝜔<0.Weknow,however,thatasumofcrosscorrelationsofreal-valuedfunctionsof
timeshouldberealvaluedinthetimedomain,requiringthat IC(𝜔)=I∗
C(−𝜔),and
F−1[IC(𝜔)] =1√
2𝜋∫+∞
−∞d𝜔IC(𝜔)ei𝜔t=1√
2𝜋(
∫+∞
0d𝜔IC(𝜔)ei𝜔t+∫+∞
0d𝜔I∗
C(𝜔)e−i𝜔t)
(30)
After substituting (29) into (30), the integration over 𝜔can be conducted analytically, making use of the
identity
∫∞
0dxeiax
√
x=⎧
⎪
⎨⎪⎩1√
a√
𝜋
2(1+i)ifa>0,
1√
−a√
𝜋
2(1−i)ifa<0,(31)
BOSCHIANDWEEMSTRA AMBIENT-NOISECROSSCORRELATION 419

ReviewsofGeophysics 10.1002/2014RG000455
Figure 5. Source-averagedcrosscorrelationsassociatedwithacircular,planardistributionofsourcessurroundingthe
tworeceivers,obtained(greenline)implementing equation(32),comparedwith(black)therealpartofthatresultingfromthenumericalapproachofsection4.1(shownalreadyinFigure2).Thesource/stationsetupisthesameasforFigure 2.
andintegratingseparatelyforthethreecases t<−Δ
c,−Δ
c<t<Δ
candt>Δ
c.Afteraconsiderableamountof
algebrawefind
IC(t)=1
(4𝜋c)21
2𝜋√
c
𝜋ΔR3×⎧
⎪
⎪
⎨
⎪
⎪⎩nC(0)√
R−Δ1√
Δ
c−tift<−Δ
c,
nC(𝜋)√
R+Δ1√
Δ
c+t+nC(0)√
R−Δ1√
Δ
c−tif−Δ
c<t<Δ
c,
nC(𝜋)√
R+Δ1√
Δ
c+tift>Δ
c,(32)
which,asrequired,haszeroimaginarypart.Weinferfromequation(32)thatinthecurrentsetup,thereisno
explicitrelationshipbetweenthesource-averagedcrosscorrelation IC(t)andtheGreen’sfunctions G2DorG3D
(AppendixE).Figure5showsthat IC(t)asobtainedfromequation(32)isconsistentwiththenumericalresult
of Figure 2. The discussion of section 4.1 remains valid. A slight discrepancy between the “analytical” and
“numerical”resultsinFigure5isexplainedbythestationary-phaseapproximationbeingstrictlyvalidonlyat
highfrequency( 𝜔−→∞).
5.2. SourcesOveraSphericalSurface
Theintegralinequation(9)canbesolvedanalyticallywiththehelpofequation(A12).Equation(9)isindeed
aparticularcaseof(A3),with x=𝜃,y=𝜑,𝜆=𝜔,f(𝜃,𝜑)=nSsin𝜃∕[8𝜋2(r2S(𝜃,𝜑)−r1S(𝜃,𝜑))],and𝜓(𝜃,𝜑)=[r2S(𝜃,𝜑)−r1S(𝜃,𝜑)]∕c. Again, the reference frame is centered on R1so thatr1S(𝜃,𝜑)=Rfor all values of 𝜃
and𝜑.Aftersomealgebra,wefind
r2S=(R2+Δ2−2ΔRsin𝜃cos𝜑)1
2, (33)
andconsequently
f(𝜃,𝜑)=1
(4𝜋c)2nS(𝜃,𝜑)sin𝜃
8𝜋2R(R2+Δ2−2ΔRsin𝜃cos𝜑)1
2, (34)
𝜓(𝜑)=1
c(r2S−r1S)=1
c[(R2+Δ2−2ΔRsin𝜃cos𝜑)1
2−R]
. (35)
Differentiating(35)withrespectto 𝜃and𝜑,wefind
𝜓𝜃=−ΔRcos𝜃cos𝜑
c(R2+Δ2−2ΔRsin𝜃cos𝜑)1
2, (36)
wherethecompact 𝜓𝜃standsfor𝜕𝜓
𝜕𝜃.Likewise,
𝜓𝜑=ΔRsin𝜃sin𝜑
c(R2+Δ2−2ΔRsin𝜃cos𝜑)1
2. (37)
BOSCHIANDWEEMSTRA AMBIENT-NOISECROSSCORRELATION 420

ReviewsofGeophysics 10.1002/2014RG000455
Ifwecontinuedifferentiating,
𝜓𝜃𝜃=−ΔRcos2𝜃cos2𝜑
c(R2+Δ2−2ΔRsin𝜃cos𝜑)3
2+ΔRsin𝜃cos𝜑
c(R2+Δ2−2ΔRsin𝜃cos𝜑)1
2, (38)
𝜓𝜑𝜑=−ΔRsin2𝜃sin2𝜑
c(R2+Δ2−2ΔRsin𝜃cos𝜑)3
2+ΔRsin𝜃cos𝜑
c(R2+Δ2−2ΔRsin𝜃cos𝜑)1
2, (39)
𝜓𝜃𝜑=ΔRsin𝜃sin𝜑cos𝜃cos𝜑
c(R2+Δ2−2ΔRsin𝜃cos𝜑)3
2+ΔRcos𝜃sin𝜑
c(R2+Δ2−2ΔRsin𝜃cos𝜑)1
2. (40)
Equations(36)and(37)allowustoidentifythestationarypoints( 𝜃,𝜑suchthat 𝜓𝜃,𝜓𝜑=0,0)oftheintegrand
at the right-hand side of (9): namely, ( 𝜃=0,𝜑=±𝜋∕2), (𝜋,±𝜋∕2), (𝜋∕2,0), and (𝜋∕2,𝜋). It is sufficient to
evaluate (A12) at these points, and sum the results, to find an analytical expression for the integral (9). Wenoticefirstofallthat
f(𝜃,𝜑)=0if𝜃=0,±𝜋:thecorrespondingstationarypointsgivenocontributiontothe
integralandwillbeneglectedinthefollowing.Weareleftwiththestationarypointsat 𝜃=𝜋∕2and𝜑=0,𝜋.
Letusevaluate f,𝜓,𝜓𝜃𝜃,𝜓𝜑𝜑,and𝜓𝜃𝜑atthosepoints.
f(𝜋
2,0)
=1
(4𝜋c)2nS(
𝜋
2,0)
8𝜋2R(R−Δ), (41)
f(𝜋
2,𝜋)
=1
(4𝜋c)2nS(
𝜋
2,𝜋)
8𝜋2R(R+Δ), (42)
𝜓(𝜋
2,0)
=−Δ
c, (43)
𝜓(𝜋
2,𝜋)
=Δ
c, (44)
𝜓𝜃𝜃(𝜋
2,0)
=ΔR
c(R−Δ), (45)
𝜓𝜃𝜃(𝜋
2,𝜋)
=−ΔR
c(R+Δ), (46)
𝜓𝜑𝜑(𝜋
2,0)
=ΔR
c(R−Δ), (47)
𝜓𝜑𝜑(𝜋
2,𝜋)
=−ΔR
c(R+Δ), (48)
anditfollowsimmediatelyfrom(40)that 𝜓𝜃𝜑=0atallstationarypoints.
Theaboveexpressionscanbesubstitutedinto(A12)tofind
IS(𝜔)≈1
(4𝜋c)2ic
4𝜋R2Δ[
nS(𝜋
2,0)e−i𝜔Δ
c
𝜔−nS(𝜋
2,𝜋)ei𝜔Δ
c
𝜔]
, (49)
whichsatisfies IS(−𝜔)=−I∗
S(𝜔).Noticethatsubstitutingequation(E22)into(49),
IS(𝜔)≈√
2𝜋
(4𝜋R)2i
𝜔[
nS(𝜋
2,0)
G3D(Δ,𝜔)−nS(𝜋
2,𝜋)
G∗
3D(Δ,𝜔)]
, (50)
BOSCHIANDWEEMSTRA AMBIENT-NOISECROSSCORRELATION 421

ReviewsofGeophysics 10.1002/2014RG000455
Figure 6. Stackedcrosscorrelationresultingfromaspherical
distributionofsourcessurroundingthetworeceivers,predicted(greenline)bytheapproximateanalyticalformula(53)and
computed(black)numericallyasdescribedinsection4.2and
alreadyillustratedinFigure3(onlytherealpartisshownhere).Thesource/stationsetupisthesameasforFigure3.and ifnS(
𝜋
2,0)
=nS(
𝜋
2,𝜋)
, equation (50)
canbefurthersimplifiedto
IS(𝜔)≈−√
2𝜋
8𝜋2R2nS(𝜋
2,0)1
𝜔ℑ[
G3D(Δ,𝜔)]
.
(51)
We next inverse Fourier transform
equation (49) to the time domain. In ourconvention,
ℱ−1(1
𝜔)
=i√
𝜋
2sgn(t)(52)
(Dirichlet integral), where the sign function
sgnis+1or−1forpositiveandnegativeval-
uesoftheargument,respectively.Inthetime
domain,ISisthen
IS(t)≈1
(4𝜋c)2c
4√
2𝜋R2Δ[
nS(𝜋
2,𝜋)
sgn(
t+Δ
c)
−nS(𝜋
2,0)
sgn(
t−Δ
c)]
, (53)
illustrated in Figure 6. This result is consistent with equation (15) of Rouxetal. [2005] and equation (27) of
Nakahara [2006],whotreatedthesamephysicalproblemindifferentways,andwithFigure1aof Harmonetal.
[2008].Bycomparisonwithequation(E21),itisapparentthatthetermsattheright-handsideofequation(53)
aretimeintegralsofthetimedomainGreen’sfunction G3D(Δ,t).Equations(49)and(53)meanthattheGreen’s
function G3Dcorrespondingtopropagationfromreceiver R1toR2andviceversacanbefoundby(i)record-
ing at both receivers the signal emitted by a dense distribution of sources covering all azimuths in three
dimensions;(ii)crosscorrelating,foreachsource,thesignalrecordedat R1withthatrecordedat R2;and(iii)
integratingthecrosscorrelationoverthesourcelocation r.
Ourequation(49)coincideswithequation(11)of Snieder[2004],exceptthatherethedifferencebetweenthe
right-handsidetermsistaken,whilein Snieder[2004]theyaresummed.InverseFouriertransformingthefor-
mulaofSnieder[2004]wouldleadtoanonphysicalcrosscorrelation,nonzeroattimeswherecrosscorrelation
is necessarily zero at all 𝜑. We infer that equations (11) and (12) in Snieder[2004] are wrong. However, this
errordidnotaffectthesubsequent,Rayleighwavederivationof Snieder[2004].
5.3. Membrane(orRayleigh)WavesFrom2-DDistributionofPointSources
5.3.1. Source-AveragedCrossCorrelation
Wetreattheintegralinequation(10)withthestationary-phasemethodasillustratedinprevioussections.The
phaseterm 𝜓(𝜑)in(10)istheadditiveinverseofthatin(8),andthetwointegrandssharethesamestationary
points𝜑=0and𝜑=±𝜋.Atthosepoints,
f(0)=nM(0)
32𝜋3c3𝜔√
R(R−Δ), (54)
f(±𝜋)=nM(𝜋)
32𝜋3c3𝜔√
R(R+Δ), (55)
𝜓(0)=Δ
c, (56)
𝜓(±𝜋)=−Δ
c, (57)
𝜓′′(0)=−ΔR
c(R−Δ), (58)
BOSCHIANDWEEMSTRA AMBIENT-NOISECROSSCORRELATION 422

ReviewsofGeophysics 10.1002/2014RG000455
𝜓′′(±𝜋)=ΔR
c(R+Δ). (59)
Substitutingintoequation(A2)(whichismultipliedby2foreachstationarypoint)andsummingoverourset
oftwostationarypoints,
IMW(𝜔)≈1
16𝜋3c3R√
𝜋c
2Δ𝜔−3
2[
nM(0)ei(𝜔Δ
c−𝜋
4)
+nM(𝜋)e−i(𝜔Δ
c−𝜋
4)]
, (60)
coherently with equations (23) and (24) of Snieder[2004], which, however, involve summation over normal
modes.Comparingequation(60)with(E17),wecanalsowrite
IMW(𝜔)≈i
4𝜋2𝜔cR√
𝜋
2[
nM(𝜋)G2D(Δ,𝜔)−nM(0)G∗
2D(Δ,𝜔)]
, (61)
validintheasymptotic(high-frequency/far-field)approximation.
Equations(60)and(61)arevalidinthestationary-phaseapproximation,i.e.,forlargeandpositive 𝜔only.We
know,however,that IMW(t)mustbereal:asinsection5.1, IMWat𝜔<0isthusdefinedby IMW(𝜔)=I∗
MW(−𝜔),
andequation(30)remainsvalidafterreplacing ICwithIMW,i.e.,
ℱ−1[IMW(𝜔)] =1√
2𝜋(
∫+∞
0d𝜔IMW(𝜔)ei𝜔t+∫+∞
0d𝜔I∗
MW(𝜔)e−i𝜔t)
. (62)
After substituting the expression (60) for IMWinto (62), it becomes apparent that finding a time domain
expressionfor IMWrequiresthesolutionofanonconvergentintegral,namely, ∫∞
0dxx−3
2cos(x).
5.3.2. DerivativeoftheSource-AveragedCrossCorrelation
Onecanstillusethepresenttheoreticalformulationtointerpret2-Ddatainthetimedomain,bysimplytak-
ingthetimederivativeofboththeobservedsource-averagedcrosscorrelationanditsanalyticalexpression
IMW(𝜔).BasedonthepropertiesofFouriertransforms,thelatterisachievedbymultiplying IMW(𝜔)byi𝜔,
I′
MW(𝜔)≈i
16𝜋3c3R√
𝜋c
2Δ⎡
⎢
⎢⎣nM(0)ei(𝜔Δ
c−𝜋
4)
√
𝜔+nM(𝜋)e−i(𝜔Δ
c−𝜋
4)
√
𝜔⎤
⎥
⎥⎦. (63)
We infer from equation (63) that for acoustic waves in 2-D, and within the stationary-phase approximation,
the derivative of the source-averaged cross correlation is proportional to the sum of causal and anticausalGreen’sfunctions
G2D(𝜔)givenbyequation(E17)[ Snieder,2004].
Like(60),equation(63)isonlyvalidforlargepositive 𝜔,butthefactthat I′
MW(t)isrealcanbeusedtodefine
I′
MW(𝜔)at𝜔<0via equation (30). The inverse Fourier transform of the resulting expression for I′
MW(𝜔)
involvesaconvergentintegral(equation(31))andcanbefoundanalytically.Theprocedureissimilartothat
ofsection5.1,equations(30)–(32).Theresultisquitedifferent,ascanbeexpectedsincetheimaginaryunitmultipliestheright-handsideofequation(63).Namely,
I′
MW(t)≈1
16𝜋3c3R√
𝜋c
Δ×⎧
⎪
⎪
⎨
⎪
⎪⎩nM(0)√
−t−Δ
cift<−Δ
c,
0if−Δ
c<t<Δ
c,
−nM(𝜋)√
t−Δ
cift>Δ
c.(64)
Figure7showsthatequation(64)isvalidatedb ythenumericalresultsofsection4.3.Noexplicitmathemat-
icalrelationshipbetween I′
MWandG2Dcanbeinferredfromequation(64),althoughFigure7showsthatthe
behaviorofthetwofunctionsisqualitativelysimilar.5.3.3. SymmetricSourceDistribution
InverseFouriertransformationturnsouttobeeasierwhen
nM(𝜋)=nM(0):settingthembothto1forsimplicity,
equation(61)collapsesto
IMW(𝜔)≈i
4𝜋2𝜔cR√
𝜋
2[G2D(Δ,𝜔)−G∗
2D(Δ,𝜔)]
≈−1
2𝜋2𝜔cR√
𝜋
2ℑ[G2D(Δ,𝜔)],(65)
BOSCHIANDWEEMSTRA AMBIENT-NOISECROSSCORRELATION 423

ReviewsofGeophysics 10.1002/2014RG000455
Figure 7. Timederivative I′
MW(t)ofthestackedcrosscorrelationassociatedwithacirculardistributionofsources
surroundingthetworeceiversonamembrane, predicted(greenline)bytheapproximateanalytical formula(64)andcomputed(blackline)numericallyasdescribedinsection4.3(Figure4).Onlytherealpartsareshownheresincewehaveverifiedthatimaginarypartsarenullasrequiredbyphysics.Forcomparison,theGreen’sfunction
−G2Dfrom
equation(E15)with x=Δisalsoshown(redline);itisnormalizedsothatitsamplitudecoincideswiththatofthe
numericalstack.Thephasevelocityandsource/stationsetuparethesameasforFigure4.
and
IMW(t)≈−1
2𝜋2cR√
𝜋
2ℱ−1{
ℑ[
G2D(Δ,𝜔)]
𝜔}
. (66)
Letusdenote Go
2Dtheoddfunction Go
2D(t)=1
2G2D(t)−1
2G2D(−t).WeknowfromAppendixB,equation(B9),
thatℑ[G2D(Δ,𝜔)]=−iGo
2D(Δ,𝜔)sothat
IMW(t)≈−1
2𝜋2cR√
𝜋
2ℱ−1[−iGo
2D(Δ,𝜔)
𝜔]
. (67)
ItalsofollowsfromAppendixB,equation(B4),that
ℱ−1[−iGo
2D(Δ,𝜔)
𝜔]
=∫t
−∞Go
2D(𝜏)d𝜏, (68)
andbydefinitionof Go
2D,
IMW(t)≈−1
4𝜋2cR√
𝜋
2{
∫t
−∞[G2D(𝜏)−G2D(−𝜏)]d𝜏}
. (69)
Differentiatingwithrespecttotime,
I′
MW(t)≈1
4𝜋2cR√
𝜋
2[G2D(−t)−G2D(t)], (70)
consistentwithFigure7.Theapparentdiscrepancybetweenequations(64)and(70)isexplainedbythefact
that the derivation of (70) involved identifying G2Dwith its asymptotic approximation (E17), which was not
thecaseforthederivationofequation(64).
5.4. PlaneWaves
5.4.1. ExactIntegration
An elegant way of reducing equation (15) to a simple and useful identity is to compare it with the integral
formofthezeroth-orderBesselfunctionofthefirstkind,
J0(z)=1
𝜋∫𝜋
0d𝜑cos(zcos𝜑) (71)
[AbramowitzandStegun ,1964,equation(9.1.18)].Substituting(71)into(15),
IPW(𝜔)=√
𝜋
2q(𝜔)J0(𝜔Δ
c)
, (72)
BOSCHIANDWEEMSTRA AMBIENT-NOISECROSSCORRELATION 424

ReviewsofGeophysics 10.1002/2014RG000455
anexactequalitythatdoesnotrequirethestationary-phaseapproximation.Importantly,equation(72)was
originallyobtainedintheearlystudyof Aki[1957],providingthebasisformuchofthelaterworkinambient
noise seismology. Comparing (72) with (E16), we infer that the source-averaged cross correlation IPW(𝜔)is
proportionaltotheimaginarypartofthemembranewaveGreen’sfunction G2D(
𝜔r
c)
ifr=Δ.
5.4.2. ApproximateIntegration
In analogy with previous sections, the integral in equation (15) can also be solved by means of the
stationary-phaseapproximation.Letusrewriteit
IPW(𝜔)=q(𝜔)√
2𝜋ℜ[
∫𝜋
0d𝜑ei𝜔Δcos𝜑
c]
. (73)
Theintegral at theright-hand sideof equation(73) coincides withthat in(A1) after replacing a=0,b=𝜋,
x=𝜑,𝜆=𝜔,and𝜓(𝜑)=Δcos𝜑
c.Takingthe 𝜑derivativeof 𝜓,
𝜓𝜑(𝜑)=−Δ
csin𝜑, (74)
andweseeimmediatelythattherearetwostationarypoin tswithintheintegrationdomain,attheintegration
limits𝜑=0,𝜋.Differentiatingagain,wefind
𝜓𝜑𝜑(𝜑)=−Δ
ccos𝜑. (75)
Atthestationarypoint 𝜑=0wehave𝜓(0)=Δ
cand𝜓𝜑𝜑(0)=−Δ
c.At𝜑=𝜋,𝜓(0)=−Δ
cand𝜓𝜑𝜑(0)=Δ
c.For
eachstationarypoint,wesubstitutethecorrespondingvaluesintoequation(A2),choosing,asusual,thesignof
𝜋∕4in the argument of the exponential based on that of 𝜓𝜑𝜑. We next sum the contributions of both
stationarypoints,finding
∫𝜋
0d𝜑ei𝜔Δcos𝜑
c≈√
2𝜋c
𝜔Δcos(𝜔Δ
c−𝜋
4)
. (76)
Substituting,inturn,(76)into(73),
IPW(𝜔)≈q(𝜔)√
c
𝜔Δcos(𝜔Δ
c−𝜋
4)
≈√
𝜋
2q(𝜔)J0(𝜔Δ
c)
, (77)
wherewehaveusedtheasymptoticapproximationof J0[AbramowitzandStegun ,1964,equation(9.2.1)],valid
inthehigh-frequency(and/orfar-field)limit,i.e.,intherangeofvalidityofthestationary-phaseapproxima-tion. Equation (77) shows that the stationary-phase appr oximation leads to an estimate of source-averaged
crosscorrelation
IPW(𝜔)consistent,atlarge 𝜔,withtheresult(72).
The relationship (72) between observed stacked cross correlations and the Bessel’s function J0has been
applied, e.g., by Ekströmetal. [2009] and Ekström[2014] to analyze ambient noise surface wave data in the
frequencydomainandmeasuretheirvelocity.TheinverseFouriertransformof(72)isobtainedanddiscussedbyNakahara [2006],whoalsoshowsthattheHilberttransformofthestackedcrosscorrelationcoincideswith
the(causalminusanticausal)
G2D[Nakahara ,2006,equation(19)].
5.4.3. MonochromaticPlaneWavesintheTimeDomain
Thetreatmentsof Sanchez-SesmaandCampillo [2006],Tsai[2009],and Boschietal. [2013]areslightlydifferent
fromtheplanewaveformulationpresentedhere,inthattheyworkwithmonochromaticplanewavesinthetime domain. Following Tsai[2009],Boschietal. [2013] make use of the properties of the Bessel and Struve
functions to solve the integral in equation (13) and are eventually able to write the source-averaged cross
correlationofplanewavesoffrequency
𝜔0travelingalongallazimuthsas
ITD(𝜔0,t)=q(𝜔0)√c
8𝜋𝜔0Δ{
cos[
𝜔0(
t+Δ
c)
−𝜋
4]
+cos[
𝜔0(
t−Δ
c)
+𝜋
4]}
(78)
[Boschietal. , 2013, equations (35) and (41)], valid in the far-field (large Δ) and/or high-frequency (large 𝜔0)
approximations.
BOSCHIANDWEEMSTRA AMBIENT-NOISECROSSCORRELATION 425

ReviewsofGeophysics 10.1002/2014RG000455
Toverifythatourequation(77)isconsistentwith(78),let ustaketheFouriertransformofthelatter.Fromthe
equality ℱ{cos[𝜔0(t+k)]}=√
𝜋
2[𝛿(𝜔+𝜔0)+𝛿(𝜔−𝜔0)]e−i𝜔kitfollowsthat
ITD(𝜔0,𝜔)=q(𝜔0)
2√c
𝜔0Δ[𝛿(𝜔+𝜔0)+𝛿(𝜔−𝜔0)]cos[
𝜔(
𝜋
4𝜔0−Δ
c)]
. (79)
Like(14),equation(79)isvalidforany 𝜔0.Implicitlyrepeatingthemonochromaticwave,timedomainanalysis
atall𝜔0’s,
ITD(𝜔)=q(𝜔)√
c
𝜔Δcos(𝜔Δ
c−𝜋
4)
, (80)
whichcoincides,asexpected,withourexpression(77)for IPW(𝜔).
6. OtherDerivations
6.1. TimeDomainApproach
We loosely follow the treatment of Rouxetal. [2005, section II], with acoustic sources uniformly distributed
throughout an unbounded, infinite 3-D medium. This setup is similar to section 5.2 here. Although in
section5.2sourcesarealllyingonthesurfaceofasphere(azimuthallyuniformdistribution),weshowinthefollowing that the formula we find for the source-averaged cross correlation is proportional to that derivedbyRouxetal. [2005].
Infree3-Dspace,thetimedomaincrosscorrelationbetweenimpulsivesignalsemittedat
Sandrecordedat
R1andR2(i.e.,thetimedomainversionofequation(6))reads
1
T∫+T
2
−T
2d𝜏G3D(r1S,𝜏)G3D(r2S,t+𝜏)=1
(4𝜋c)21
T∫+T
2
−T
2d𝜏𝛿(
𝜏−r1S
c)
𝛿(
𝜏+t−r2S
c)
r1Sr2S
=1
(4𝜋c)21
T𝛿(
t+r1S
c−r2S
c)
r1Sr2S,(81)
where we have used expression (E21) for G3D, and we have limited cross correlation to a finite time interval
(−T∕2,T∕2),whoselength Tisrelated,inpractice,tointerstationdistanceandwavespeed.
Asdiscussedinsection3,ifcrosstermsareneglected(AppendixD),thecrosscorrelationofadiffusewavefield
recorded at R1andR2is estimated by summing expression (81) over many uniformly distributed sources S.
Forsourcesdenselydistributedthroughouttheentire3-Dspace R3,
I3D(t)=1
(4𝜋c)21
T∫R3d3r𝛿(
t+r1S
c−r2S
c)
r1Sr2S, (82)
whereris the source location, and, said r1andr2, respectively, the locations of R1andR2,r1S=|r−r1|,
r2S=|r−r2|.Equation(82)isequivalenttoequation(9)of Rouxetal. [2005].
Asnotedby Rouxetal. [2005],equation(82)showsthatawiggleinthecrosscorrelation I3Dattimetisneces-
sarilyassociatedtooneormoresourceswhoselocation rsatisfiesct=|r−r2|−|r−r1|.Intwodimensions,
the latter is the equation of a hyperbola with foci at R1atR2; in three dimensions, it is the equation of the
single-sheet hyperboloid obtained by rotation of said hyperbola around the vertical axis. It follows that theintegralin(82)canbereducedtoasurfaceintegraloverthehyperboloid.Afterthissimplification, Rouxetal.
[2005] are able to solve the integral in (82) analytically: they find
I3D(t)to be a boxcar function between
t=±Δ ∕c,equivalenttoequation(49)insection5.2here.
6.2. ReciprocityTheoremApproach
All descriptions of seismic and acoustic interferometry that we discussed so far rest on equation (1) and on
theGreen’sfunctionformulaeofAppendixE,whicharestrictlyvalidonlyfor g,infinitemedia.Moregeneral
formulations have been developed by Wapenaar [2004],WeaverandLobkis [2004],vanManenetal. [2005],
BOSCHIANDWEEMSTRA AMBIENT-NOISECROSSCORRELATION 426

ReviewsofGeophysics 10.1002/2014RG000455
WapenaarandFokkema [2006],Snieder[2007], and others, based on the reciprocity or Betti’s theorem [e.g.,
AkiandRichards ,2002,section2.3.2].
6.2.1. AcousticWaves
Thepropagationofacousticwavesinanonhomogeneous,losslessstagnantgasofdensity 𝜌andcompress-
ibility𝜅isdescribedbyequation(6.2.7)of MorseandIngard [1986],plustheforcingterm qcorrespondingto
massinjection.Inthefrequencydomain,
∇p+i𝜔𝜌v=0, (83)
∇⋅v+i𝜔𝜅p−q=0. (84)
Inthemoregeneralcaseofasound-absorbingmedium, 𝜅iscomplexanditsimaginarypartisproportional
to the rate of energy loss (attenuation) [ Kinsleretal. , 1999, chap. 8]. Equations (83) and (84) are equivalent
to equations (2) and (3) of WapenaarandFokkema [2006] or equations (1) and (2) of Snieder[2007] where a
differentFouriertransformconventionapplies.Equation(83)implies v=1
i𝜔𝜌∇p,andthisexpressionfor vcan
besubstitutedinto(84),toobtain
∇⋅(
1
𝜌∇p)
+𝜔2𝜅p=−i𝜔q. (85)
Forhomogeneousorsmoothmedia,wherethespatialderivativesof 𝜌and𝜅arezeroorapproximatelyzero,
(85)furthersimplifiestotheFouriertransformofequation(1),
1
𝜌∇2p+𝜔2𝜅p=−i𝜔q. (86)
Following WapenaarandFokkema [2006],letuscall 𝒢thesolutionof(86)when q(x,t)=𝛿(x)𝛿(t):
1
𝜌∇2𝒢+𝜔2𝜅𝒢=−i𝜔𝛿(x). (87)
𝒢canbeinterpretedastheGreen’sfunctionassociatedwithequation(84)plusthecondition(83).
The relationship between 𝒢andG3Das defined in Appendix E can be determined if one considers that
equation (87) is the usual scalar wave equation (whose Green’s function is G3D) with forcing term −i𝜔𝛿(x);
based on equation (E34), the time domain solution to (87) is then the convolution of G3Dwith the inverse
Fourier transform of the forcing term −i𝜔𝛿(x); in the frequency domain the convolution reduces to the
productofthefunctionsinquestion,and
𝒢=−i𝜔G3D. (88)
WapenaarandFokkema [2006] also introduce a “modified Green’s function” (see their equation (32)) which
coincideswithour G3Dexceptforthesign.
Let us next consider a volume Vbounded by a surface 𝜕V.(𝜕Vis just an arbitrary closed surface within a
medium and generally does not represent a physical boundary.) Let qA(r,𝜔),pA(r,𝜔), andvA(r,𝜔)denote a
possible combination of forcing, pressure, and velocity, respectively, coexisting at rinVand𝜕V. A different
forcingqBwould give rise, through equations (83) and (84), to a different “state” B, defined by pB(r,𝜔)and
vB(r,𝜔).
Ausefulrelationshipbetweenthestates AandB,knownas“reciprocitytheorem,”isobtainedbycombining
equations(83)and(84)asfollows:
∫Vd3r[(83)A⋅v∗
B+(83)∗
B⋅vA+(84)Ap∗
B+(84)∗
BpA]=0 (89)
[e.g.,WapenaarandFokkema , 2006;Snieder, 2007], where (83)Ais short for the expression one obtains after
substituting p=pA(r,𝜔)andv=vA(r,𝜔)intotheleft-handsideofequation(83),etc.Namely,
(83)A⋅v∗
B=∇pA⋅v∗
B+i𝜔𝜌vA⋅v∗
B(90)
BOSCHIANDWEEMSTRA AMBIENT-NOISECROSSCORRELATION 427

ReviewsofGeophysics 10.1002/2014RG000455
(83)∗
B⋅vA=∇p∗
B⋅vA−i𝜔𝜌v∗
B⋅vA (91)
(84)Ap∗
B=∇⋅vAp∗
B+i𝜔𝜅pAp∗
B−qAp∗
B(92)
(84)∗
BpA=∇⋅v∗
BpA−i𝜔𝜅∗p∗
BpA−q∗
BpA. (93)
Since(83)A=(83)B=(84)A=(84)B=0byvirtueofequation(83)and(84),itfollowsthattheexpressionat
theleft-handsideof(89)equals0asanticipated.Aftersubstituting(90)–(93)into(89),
∫Vd3r(∇pA⋅v∗
B+∇⋅v∗
BpA)+∫Vd3r(∇p∗
B⋅vA+∇⋅vAp∗
B)
−∫Vd3ri𝜔(𝜅∗−𝜅)pAp∗
B−∫Vd3r(qAp∗
B+q∗
BpA)=0.(94)
It is convenient to apply the divergence theorem [e.g., Hildebrand , 1976] to the first two integrals at the
left-handsideofequation(94),whichthensimplifiesto
∫𝜕Vd2r(pAv∗
B+p∗
BvA)⋅̂n+2i𝜔∫Vd3rℑ(𝜅)pAp∗
B−∫Vd3r(q∗
BpA+qAp∗
B)=0, (95)
wherênis the unit vector normal to 𝜕Vandℑ(𝜅)denotes the imaginary part of 𝜅. Equation (95) is equiva-
lent to the “reciprocity theorem of the convolution type,” equation (5) of WapenaarandFokkema [2006] or
equation(4)of Snieder[2007].
An equation relating the cross correlation of ambient signal to Green’s functions can be found from (95) by
considering the case qA,B(r)=𝛿(r−rA,B), withrA,Barbitrary locations in V. It follows that pA,B=𝒢(r,rA,B).
Substitutingtheseexpressionsfor pandqinto(95)andusing(83)toeliminatethevelocity,
𝒢(rB,rA)+𝒢∗(rA,rB)=
=1
i𝜔∫𝜕Vd2r1
𝜌[𝒢∗(r,rB)∇𝒢(r,rA)−𝒢(r,rA)∇𝒢∗(r,rB)]⋅̂n
+2i𝜔∫Vd3rℑ(𝜅)𝒢(r,rA)𝒢∗(r,rB)(96)
[e.g.,WapenaarandFokkema ,2006;Snieder,2007;CampilloandRoux ,2014].
Notice that the treatment that leads from equa tion (89) to (96) remains valid for heterogeneous 𝜅and𝜌;
equation(96)holdsforaheterogeneous,attenuatingmediumthatcouldbeboundedorunbounded.Itisthus
more general than similar equations (32), (53), (63), and (72), which are only strictly valid if the propagationmediumishomogeneous.
Provided that the medium be smooth or homogeneous at and near
𝜕V,𝒢still coincides with −i𝜔G3Don
𝜕V, withG3Dthe homogeneous medium Green’s function. ∇𝒢=−i𝜔∇G3Dcan then be computed through
equation(E22),whichimplies
∇G3D(r,𝜔)=r
r[2
rG3D(r,𝜔)−i𝜔
cG3D(r,𝜔)]
. (97)
Ifwemakethefurtherassumptionthatallsourcesarefarfrom 𝜕V(risalwayslarge),itfollowsthat2
rG3Dcan
beneglectedinequation(97),whiler
r⋅̂n≈1sothat
∇G3D(r,𝜔)⋅̂n≈−i𝜔
cG3D(r,𝜔), (98)
and
∇𝒢(r,𝜔)⋅̂n≈−i𝜔
c𝒢(r,𝜔) (99)
BOSCHIANDWEEMSTRA AMBIENT-NOISECROSSCORRELATION 428

ReviewsofGeophysics 10.1002/2014RG000455
[e.g.,Wapenaar and Fokkema , 2006;Snieder, 2007;Campillo and Roux , 2014]. The surface integral in
equation(96)isaccordinglysimplified,
𝒢(rB,rA)+𝒢∗(rA,rB)=
≈−2
𝜌c∫𝜕Vd2r[
𝒢∗(r,rB)𝒢(r,rA)]
+2i𝜔∫Vd3rℑ(𝜅)𝒢(r,rA)𝒢∗(r,rB).(100)
Our treatment in section (5.2) was limited to nonattenuating media, where ℑ(𝜅)=0. In this case
equation(100)reducesto
𝒢(rB,rA)+𝒢∗(rA,rB)≈−2
𝜌c∫𝜕Vd2r[𝒢∗(r,rB)𝒢(r,rA)]. (101)
Applyingequation(88),i.e.,replacing 𝒢=−i𝜔G3D,
i𝜌c
2𝜔[G∗
3D(rA,rB)−G3D(rB,rA)]=−∫𝜕Vd2r[G∗
3D(r,rB)G3D(r,rA)], (102)
andbyvirtueofthereciprocity G3D(rB,rA)=G3D(rA,rB),etc.,
𝜌c
𝜔ℑ[G3D(rA,rB)] = −∫𝜕Vd2r[G∗
3D(r,rB)G3D(r,rA)]. (103)
Theright-handsideofequation(103)issimply IS(𝜔)asdefinedbyequation(9);since nSinequation(51)isarbi-
trary,equations(103)and(51)areequivalent.Inasmooth,losslessmediumilluminatedfromallazimuths,thestationary-phaseandreciprocitytheoremapproachesleadtothesamerelationshipbetweenGreen’sfunctionand cross correlation, establishing, in practice, that the Green’s function between
rAandrBcan be recon-
structedfromobservationsaslongasthemediumisilluminatedbyadensedistributionofsourcescoveringitsboundary
𝜕V(section5.2).
Themoregeneral,reciprocitytheorem-basedresults(96)and/or(100)applytoattenuatingmedia. ℑ(𝜅)≠0
impliesthatthevolumeintegralsattheirright-handsidescannotbeneglected;reconstructionofanattenu-
atingmedium’sGreen’sfunctionfromthedatarequiresthatthemediumbeilluminatedbysources withinV
[e.g.,CampilloandRoux ,2014].
6.2.2. SeismicWaves
In a series of articles, Kees Wapenaar and coworkers have applied the above ideas to the case of an elas-
tic medium, where both compressional and shear deformation exist [e.g., Wapenaar , 2004;Wapenaar and
Fokkema, 2006;Wapenaar et al. , 2006]. Their procedure, based on applying the reciprocity theorem to a
pair of states both excited by impulsive point sources, is qualitatively similar to the acoustic wave formu-
lation of Snieder[2007], illustrated here in section 6.2.1. The most complete description of the reciprocity
theoremapproachisthatof Wapenaaretal. [2006],whoallowformediuminhomogeneityandattenuation.
Wapenaar et al. [2006] show that in analogy with equations (100)–(103) for the acoustic case, the Green’s
functioncanbereconstructedfromnoisecrosscorrelationprovidedthatthemediumisilluminatedbynoisesources densely distributed throughout a volume
V, where receivers are immersed. If the medium is loss-
less,sourceswithin Vareunnecessary,butilluminationfromsourcesdistributedthroughoutthesurfacethat
boundsVisstillneeded.
6.3. Normal-ModeApproach
LobkisandWeaver [2001]useanormal-modeapproachtofindananalyticalexpressionfordiffusefieldcross
correlation in a bounded medium. Following Snieder et al. [2010], we briefly repeat their treatment for the
simplercaseofalosslessboundedmedium.Normalmodesaredefinedastherealfunctions pn(r)suchthat
pn(r)cos(𝜔nt),pn(r)sin(𝜔nt)(n=1,2,3,…,∞), with eigenfrequencies 𝜔n, form a complete set of solutions to
the homogenous version of the scalar wave equation (1). Any wavefield p(r,t)propagating in the medium
underconsiderationcanbewrittenasalinearcombinationofmodeswithcoefficients an,bn
p(r,t)=∑
n[anpn(r)cos(𝜔nt)+bnpn(r)sin(𝜔nt)], (104)
where∑
ndenotessummationovertheintegervaluesof nfrom1toinfinity.
BOSCHIANDWEEMSTRA AMBIENT-NOISECROSSCORRELATION 429

ReviewsofGeophysics 10.1002/2014RG000455
6.3.1. Green’sFunctionasaLinearCombinationofModes
Let us define the Green’s function GMas the solution of the homogeneous version (𝜕F
𝜕t=0) of equation (1),
withinitialconditions p=0and𝜕p
𝜕t=𝛿(r−s)(pointsourceat s).GM(r,t)canbewrittenasalinearcombination
ofmodes,
GM(r,t)=∑
k[
𝛼kcos(𝜔kt)+𝛽ksin(𝜔kt)]
pk(r). (105)
Substitutinginto(105)theinitialconditionon p,wefind
∑
k𝛼kpk(r)=0; (106)
afterdifferentiating(105)withrespecttotime,theinitialconditionon𝜕p
𝜕tgives
∑
k𝜔k𝛽kpk(r)=𝛿(r−s). (107)
We multiply both sides of equations (106) and (107) by the eigenfunction pn(r)and integrate over r.F r o m
equation(107)wefind
∑
k𝜔k𝛽k∫R3d3rpk(r)pn(r)=∫R3d3r𝛿(r−s)pn(r), (108)
which,aftertakingadvantageoftheorthonormalityofthemodes(left-handside)andapplyingtheproperties
ofthe𝛿function(right-handside),collapsesto
𝛽n𝜔n=pn(s). (109)
Equation(106)likewisereducesto 𝛼k=0,andequation(105)becomes
GM(r,t)=∑
npn(s)pn(r)sin(𝜔nt)
𝜔n, (110)
validfort>0,consistentwithequation(4)of LobkisandWeaver [2001]andequation(1)of Sniederetal. [2010].
6.3.2. Ensemble-AveragedCrossCorrelationasaLinearCombinationofModes
Let us write the ambient signals recorded at R1andR2as linear combinations of modes with random
coefficients anandbn.Theircrosscorrelationthenreads
∫+T∕2
−T∕2d𝜏p∗(r1,𝜏)p(r2,t+𝜏)=∑
n,k{
akanpk(r1)pn(r2)∫+T∕2
−T∕2d𝜏cos(𝜔k𝜏)cos[𝜔n(t+𝜏)]
+akbnpk(r1)pn(r2)∫+T∕2
−T∕2d𝜏cos(𝜔k𝜏)sin[𝜔n(t+𝜏)]
+bkanpk(r1)pn(r2)∫+T∕2
−T∕2d𝜏sin(𝜔k𝜏)cos[𝜔n(t+𝜏)]
+bkbnpk(r1)pn(r2)∫+T∕2
−T∕2d𝜏sin(𝜔k𝜏)sin[𝜔n(t+𝜏)]}
.(111)
LobkisandWeaver [2001]maketheassumptionthat“modalamplitudesareuncorrelatedrandomvariables,”
whichisequivalenttonoisesourcesbeingspatiallyandtemporally uncorrelatedsothat“crossterms”canbe
neglected(AppendixD).Thismeansinpracticethatifonerepeatsthecrosscorrelation(111)atmanydifferent
times, the normal modes of the medium stay the same (medium properties do not change), but the coeffi-cients
anandbnchangeinarandomfashionateachrealization.Whentheaverageofallrealizationsistaken,
theproducts akanandbkbnbothaverageto 𝛿knM(𝜔n),withthefunction Mindicatinghowstronglydifferent
BOSCHIANDWEEMSTRA AMBIENT-NOISECROSSCORRELATION 430

ReviewsofGeophysics 10.1002/2014RG000455
eigenfrequenciesareexcitedonaverage,and 𝛿kn=1ifk=n,0otherwise; akbn,ontheotherhand,averages
to0forallvaluesof k,n.Using⟨•⟩todenotetheaveragingprocedure,itfollowsfrom(111)that
⟨∫+T∕2
−T∕2d𝜏p∗(r1,𝜏)p(r2,t+𝜏)⟩=∑
n,k{
⟨akan⟩pk(r1)pn(r2)∫+T∕2
−T∕2d𝜏cos(𝜔k𝜏)cos[𝜔n(t+𝜏)]
+⟨akbn⟩pk(r1)pn(r2)∫+T∕2
−T∕2d𝜏cos(𝜔k𝜏)sin[𝜔n(t+𝜏)]
+⟨bkan⟩pk(r1)pn(r2)∫+T∕2
−T∕2d𝜏sin(𝜔k𝜏)cos[𝜔n(t+𝜏)]
+⟨bkbn⟩pk(r1)pn(r2)∫+T∕2
−T∕2d𝜏sin(𝜔k𝜏)sin[𝜔n(t+𝜏)]}
=∑
nM(𝜔n)pn(r1)pn(r2)∫+T∕2
−T∕2d𝜏cos(𝜔nt)
=∑
nM(𝜔n)pn(r1)pn(r2)Tcos(𝜔nt),(112)
wherethetrigonometricidentity cos(𝜔n𝜏)cos[𝜔n(t+𝜏)]+sin(𝜔n𝜏)sin[𝜔n(t+𝜏)] =cos(𝜔nt)isused[Snieder
etal.,2010].Equation(111)isequivalenttoequations(9)or(10)of Sniederetal. [2010],exceptthatwehave
chosennottonormalizethecoefficients anandbnbythecorrespondingeigenfrequency 𝜔n.
Provided that Mis constant with respect to 𝜔(all modes are equally excited), the time derivative of
the right-hand side of equation (112) is proportional to the right-hand side of (110): in other words, theensemble-averagedcrosscorrelationofadiffusefieldrecordedat
r1andr2isproportionaltothetimederiva-
tiveoftheGreen’sfunction GM(equation(110)insection6.3.1),forasourcelocatedat r1andareceiverlocated
atr2.Thisisequivalenttoequations(50)and(103),validfornonsteadystate3-Dacousticmedia.
6.4. AnalogyBetweenDiffuseFieldandTime-ReversalMirror
Derodeetal. [2003] explain the relationship between diffuse field cross correlation and Green’s function via
the concept of time-reversal mirror [e.g., Fink, 1999, 2006]. A time-reversal mirror can be thought of as an
array of transducers that record sound, reverse it with respect to time, and emit the time-reversed acousticsignal;ifthearrayissufficientlylargeanddense,itwilltimereversetheentirepropagatingwavefield,focusingtime-reversedwavesbacktotheoriginoftheinitialsignal.
Following Stehly[2007],wenextsummarizethereasoningof Derodeetal. [2003]inthreesimplesteps.Firstof
all(section6.4.1),crosscorrelatingtwosignalsisequivalenttotimereversingthefirstsignalthenconvolvingitwiththesecond.Wethenshow(section6.4.2)thattheconvolutionoftwoimpulsivesignalsemittedfrom
r
andrecordedatthelocations r1andr2coincideswiththesignalrecordedat r2afterbeingemittedbyasource
atr1and then time reversed and reemitted by a transducer at r. It follows (section 6.4.3) that if instead of a
singletransducerat ranentirearrayoftransducersforminga time-reversalmirrorarepresent,thementioned
convolution coincides with the Green’s function: through the relationship between convolution and crosscorrelation,anequationconnectingGreen’sfunctionandcrosscorrelationisthusdetermined.6.4.1. Convolution,CrossCorrelation,andTimeReversal
Letusfirstrecallthedefinitionoftheconvolution
f⊗goftworeal-valuedfunctions f(t),g(t):
f⊗g(t)=∫+∞
−∞f(𝜏)g(t−𝜏)d𝜏. (113)
Itcanbeshownthattheconvolutionoperatoriscommutative,
f⊗g=g⊗f, (114)
andassociative,
(f⊗g)⊗h=f⊗(g⊗h). (115)
BOSCHIANDWEEMSTRA AMBIENT-NOISECROSSCORRELATION 431

ReviewsofGeophysics 10.1002/2014RG000455
Figure 8. Setupofsection6.4.1,whereasignal e(t)emittedbyasourceat r(star) is recorded by receivers at r1andr2
(triangles).
Convolution and cross correlation are closely related: if one denotes f−(t)≡f(−t)the function found by
systematicallychangingthesignoftheargumentof f,thatistosay,bytime-reversing f,itfollowsthat
f−⊗g(t)=∫+∞
−∞f(−𝜏)g(t−𝜏)d𝜏
=−∫−∞
+∞f(𝜏)g(t+𝜏)d𝜏
=∫+∞
−∞f(𝜏)g(t+𝜏)d𝜏,(116)
i.e.,theconvolutionof f−withgcoincideswiththecrosscorrelationof fwithg[e.g.,Smith,2011].
Nowconsiderasourceatthelocation remittingthe(real-valued)signal e(t)andtworeceiversatthelocations
r1andr2.BydefinitionofGreen’sfunction G,thereceiverat r1recordsasignal e⊗G(r1,r)(t),whilethatat r2
recordse⊗G(r2,r)(t)(Figure8).Letuscrosscorrelatethesignalrecordedat r1withthatrecordedat r2and
denoteC12thecrosscorrelation.Makinguseofequation(116)andofthefactthatconvolutioniscommutative
andassociative,
C12=[e−⊗G−(r1,r)(t)]⊗[e⊗G(r2,r)(t)]=[G−(r1,r)⊗G(r2,r)(t)]⊗[e−⊗e(t)].(117)
6.4.2. PropagationFromOneSourcetoaTransducerandFromtheTransducertoaReceiver
Letusnowplaceatthelocation r1asoundsourcethatemitsanimpulsivesignal,andat ratransducerthat
Figure 9. Setupofsection6.4.2,whereanimpulseemittedbyasource
atr1(triangle)isrecordedbyatransducerat r(star),whichthentime
reversesitandemitsitback.Areceiverat r2(triangle)wouldthen
recordtheconvolutionofthetime-reversedinitialimpulsewiththeGreen’sfunctioncorrespondingtothelocations
randr2.recordssound,timereversesit,andemits
itback(Figure9).Areceiverisstillplaced
atr2. The signal recorded at rcoincides
with the Green’s function G(r,r1,t).T h e
signalp(r2,t)recorded at r2at a time tis
the convolution of the signal emitted by
r, that is, the time-reversed Green’s func-
tionG(r,r1,−t),withtheGreen’sfunction
G(r2,r,t),
p(r2,t)=G−(r,r1)⊗G(r2,r)(t).(118)
(This is valid in the assumption that the
emitted wavelet is short enough and/or
ris far enough from both r1andr2for
the time-reversed wave packet p(r2,t)
to be easily isolated from the “direct”
arrival at r2.) By the reciprocity of G,
G−(r,r1,t)=G−(r1,r,t).Wecanthussub-
stituteequation(118)into(117),and
C12=p(r2,t)⊗[e−⊗e(t)],(119)
BOSCHIANDWEEMSTRA AMBIENT-NOISECROSSCORRELATION 432

ReviewsofGeophysics 10.1002/2014RG000455
i.e., the cross correlation of the recordings, made at r1andr2, of the signal emitted by a source at rcoin-
cideswiththerecordingmadeat r2ofthesamesignal,emittedat r1andrecordedandtimereversedat r.For
impulsivesignals e(t)=𝛿(t),equation(119)collapsesto
C12=p(r2,t). (120)
6.4.3. MultipleSourcesandTime-ReversalMirror
We next consider a dense, uniform distribution of transducers surrounding the locations r1andr2. Such a
setoftransducersformsatime-reversalmirror[e.g., Fink,1999,2006].Equations(117)and(118)remainvalid,
after replacing G(r1,r,t),G(r2,r,t)with the sum of Green’s functions associated with all locations rwhere
transducersarenowplaced;atthelimitofacontinuoussourcedistribution,
C12=∫Ωdr[G−(r1,r)⊗G(r2,r)(t)]⊗[e−⊗e(t)], (121)
whereΩistheentiresolidangle.
In this setup (Figure 10), an impulsive signal emitted at r1is first recorded at r2asG(r2,r1); it then hits the
transducer array, which, by definition of time-reversal mirror, sends the signal back in such a way that thesame wavefield propagates backward in time: the receiver at
r2records the time-reversed Green’s function
G(r2,r1,−t),andtheback-propagatedsignaleventuallyfocusesbackon r1intheformofanimpulseattime
−t=0.
This means, in practice, that if the individual transducer at r(Figure 8) is replaced by a time-reversal mirror,
thenp(r2,t)inequation(120)canbereplacedbytheGreen’sfunction G(r2,r1,−t),
C12=G(r2,r1,−t). (122)
Now recall from section 6.4.1 that C12is also the cross correlation of the recordings made at r1andr2of
signal generated at r; the individual source at rmust now be replaced by a set of sources occupying the
transducer locations so that cross correlations associated with individual sources are summed. Derodeetal.
[2003]correctlyinferthatifasignalisgeneratedbyadistributionofsourceswiththegeometryofaneffec-
tive time-reversal mirror (i.e., energy propagating along all azimuths and diffuse field), the Green’s function
betweenanytwopoints r1andr2canbefoundbycrosscorrelatingtherecordingsofthesaidsignalmadeat
r1andr2.Thisstatementis,again,equivalenttoequations(50)and(103).
7. UnevenSourceDistributionsand“SpuriousArrivals”
Our derivation so far is based on the hypothesis that the geographic distribution of noise sources be closetouniformwithrespecttosource-receiverazimuth.Thestationary-phaseformulaeofAppendixAonlyholdif
fis asmoothfunction of xin equation (A1) and of both xandyin equation (A3); the source distributions
nC,nM, andnSmust accordingly be smooth with respect to 𝜑and/or𝜃for the treatment of section 5 to be
valid.Theintegralinequation(103)likewiseextendstothewholeboundaryofthevolume Vcontainingthe
receivers:if 𝜕Visnotcovereddenselyanduniformlybysources,noisecrosscorrelationdoesnotcoincidewith
theright-handsideofequation(103),and Gisnotproperlyreconstructed.
Noise sources are generally not uniformly distributed in practical applications, and we know, e.g., from
Mulargia [2012],thatseismicambientnoiseonEarthisnotstrictlydiffuse.Weillustratetheconsequencesof
significant inhomogeneities in source distribution with a simple model. As in sections 4.3 and 5.3, receivers
R1andR2,lying20kmfromoneanotheronamembraneofinfiniteextension,aresurroundedbyacircleof
sourceswhosecenteris R1andwhoseradiusis100km(Figure11).WenumericallyconvolveaRickerwavelet
(centralfrequencyof1Hz)withtheGreen’sfunction G2Dforeachofthesourcesinquestion.Usingawavelet
rather than an impulse allows to better visualize the effects we are interested in. For each location of the
source,wecrosscorrelatethecorrespondingsignalsat R1andR2andplotthecrosscorrelationsinFigure12a.
Theresultofstackingthecrosscorrelations,showninFigure12b,isconsistentwiththeresultsofsection5.3,
after modulating the Green’s function with the Ricker wavelet (we shall speak of “Ricker response” insteadofGreen’s function). We next average only the cross correlations associated with sources denoted in green
BOSCHIANDWEEMSTRA AMBIENT-NOISECROSSCORRELATION 433

ReviewsofGeophysics 10.1002/2014RG000455
Figure10. Setupsofsection6.4.3.(a)Signalsemittedbyacircleofsources(stars)arerecordedbyareceiverpairat r1
andr2(triangles);thisisthesamesetupasinFigure8(section6.4.1),butnowthereismorethanonesource.(b)An
impulseemittedbyasourceat r1isrecordedbyasetoftransducerswhichtimereverseitandemititback.Thisisthe
samesetupatinFigure9(section6.4.2)butwithmorethanonetransducer.Thetransducersoccupythesamelocations(stars)asthesourcesinFigure10aandthusformatime-reversalmirror(intwodimensions). Derodeetal. [2003]and
Stehlyetal. [2006]usethepropertiesofacoustictimereversaltoprovethatinasetupsuchasFigure10a,thecross
correlationbetween thesignalsrecordedat
r1andr2canbeassociatedwiththeGreen’sfunctionbetweenthesame
twolocations.
BOSCHIANDWEEMSTRA AMBIENT-NOISECROSSCORRELATION 434

ReviewsofGeophysics 10.1002/2014RG000455
Figure11. Distributionofsources(stars)andreceivers(triangles)inthesetupofsection7.Only10%ofthesimulated
sourcesareshown.Insomeofthesimulationsdiscussedinthefollowing,theeffectsofunevensourcedistributionsincludingonlythesourcesdenotedbygreenversusyellowstarsaremodeled.
inFigure11and,finally,onlythoseassociatedwiththe“yellow”sourcesofFigure11.Twoinferencescanbe
madefromFigure12c,wherebothaveragesareshown:(i)ifonlytheyellowsourcesare“on,”andenergyonlytravels in the direction
R2−→R1, only the anticausal Ricker’s response between R1andR2emerges from
averaging; likewise, only the causal part shows up if only sources to the left of R1are active. (ii) While both
causal and anticausal arrivals in Figure 12b approximately coincide with those of Figure 12c, the curves in
Figure12ccontaintwoadditionalarrivals,correspondingtothetwoazimuthswherebothsourcedistributions
inFigure11abruptlyend.Thesearrivals,usuallyreferredtoasspurious,havenorelationtotheRickerresponse;theyareartifactscausedbystronginhomogeneitiesinthesourcedistribution.Spuriousarrivalsarelikelyto
affectfielddataandcanbeidentifiedinlaboratory(physicalacoustics)data.
8. Ambient Signal Cross Correlation in the Presence of a Scatterer:
AStationary-PhaseDerivation
The stationary-phase derivations carried out above have established mathematical relationships between
two-station ambient signal cross correlation and a medium’s Green’s function, in the simple case ofhomogeneous , unbounded media. The same approach can also usefully be applied to a homogeneous
mediumincludingalimitednumberofpointscatterers.Following Sniederetal. [2008],weshalltreatinsome
detailthecaseofahomogeneous,3-Dacousticmediumcontainingasinglepointscatterer:extensiontomorescatterers[ Fleuryetal. ,2010]isthenstraightforward,albeitcumbersome.
Itisconvenienttoplacetheoriginofthecoordinatesystematthelocationofthescattererandtochoosethe
xandzaxessothattheplanetheyidentifycontainsthelocations r1,r2ofreceivers R1andR2.Inthissetup,
theGreen’sfunction GS
3Disthesumof G3Dfromequation(E22)plusanadditional,“scattered”terminvolving
propagationfromthesource(locatedatapoint s)tothescattererandfromthescatterertothereceiver,
GS
3D(r1,2,s,𝜔)=G3D(r1,2S,𝜔)+G3D(r1,2,𝜔)G3D(s,𝜔)h(̂r1,2,̂s), (123)
wherer1,2andsdenotethemoduliof r1,2ands,respectively,and ̂r1,2,̂sareunitvectorsparallelto r1,2ands.
Thescatteringfunction(or“matrix”) h(̂r1,2,̂s)accountsfortheazimuthdependenceoftheamountofscattered
energy[e.g., Snieder,1986;GroenenboomandSnieder ,1995;Marston,2001].
BOSCHIANDWEEMSTRA AMBIENT-NOISECROSSCORRELATION 435

ReviewsofGeophysics 10.1002/2014RG000455
Figure12. Crosscorrelationsassociatedwithacircular,planardistributionofsourcessurroundingthetworeceiversas
sketchedinFigure11.EachsourcegeneratesaRickerwavelet(centralfrequency=1Hz);thewaveletisconvolvedwiththeGreen’sfunction(phasevelocity
c=2Km/s)toevaluate thesignalsobservedatthetworeceivers,whicharethen
crosscorrelated.(a)Single-sourcecrosscorrelationsforallsourceazimuths;thedashedlinesmarktheazimuthsseparatingthetwo(yellowversusgreen)subsetsofsourcesasdefinedinFigure11.(b)StackedcrosscorrelationresultingfromFigure12a.(c)StackedcrosscorrelationsthatonewouldobtainifsignalwasgeneratedonlyatthelocationsidentifiedbygreenstarsinFigure11(greenline)versustheyellowstars(yellowline).
Thefrequencydomaincrosscorrelationofthesignalsrecordedat R1andR2reads
GS
3D(r1,s,𝜔)GS∗
3D(r2,s,𝜔)=G3D(r1S,𝜔)G∗
3D(r2S,𝜔)
+G3D(r1S,𝜔)G∗
3D(r2,𝜔)G∗
3D(s,𝜔)h∗(̂r2,̂s)
+G3D(r1,𝜔)G3D(s,𝜔)G∗
3D(r2S,𝜔)h(̂r1,̂s)
+G3D(r1,𝜔)G3D(s,𝜔)G∗
3D(r2,𝜔)G∗
3D(s,𝜔)h(̂r1,̂s)h∗(̂r2,̂s).(124)
Wenextassumesourcestobedistributedoverasphericalsurfaceofradius s=Rsurroundingthereceivers,
asinsection3,scenario(ii);forthesakeofsimplicity,welimitourselvestouniformsourcedensity nS(𝜃,𝜑)=1.
Neglectingcrossterms(AppendixD),thecumulativeeffectofsuchasourcedistributionisobtainedbyinte-
grating equation (124) in sover the whole solid angle. The integral IS(𝜔)of the first term at the right-hand
side of (124) has been treated in detail in section 5.2, and its analytical form is given, e.g., by equation (50);
BOSCHIANDWEEMSTRA AMBIENT-NOISECROSSCORRELATION 436

ReviewsofGeophysics 10.1002/2014RG000455
substituting nS(𝜃,𝜑)=1,
IS(𝜔)≈√
2𝜋
(4𝜋R)2i
𝜔[G3D(Δ,𝜔)−G∗
3D(Δ,𝜔)], (125)
whereΔdenotes, as usual, the distance between R1andR2. Let us call I2(𝜔),I3(𝜔), andI4(𝜔)the source
averagesoftheremainingthreetermsattheright-handsideof(124).Thefirstoftheseintegrals,
I2(𝜔)=1
4𝜋∫𝜋
0d𝜃sin𝜃∫𝜋
−𝜋d𝜑G3D(r1S(𝜃,𝜑),𝜔)G∗
3D(r2,𝜔)G∗
3D(s(𝜃,𝜑),𝜔)h∗(̂r2,̂s(𝜃,𝜑))
=1
(4𝜋)41
(
c√
2𝜋)3∫𝜋
0d𝜃sin𝜃∫𝜋
−𝜋d𝜑ei𝜔
c(R+r2−r1S(𝜃,𝜑))
Rr2r1S(𝜃,𝜑)h∗(̂r2,̂s(𝜃,𝜑)),(126)
can be simplified by making the hypothesis that the source be very far from both receivers and from the
scatterer[ Sniederetal. ,2008],i.e., R≫r1,2,whichimpliesthat r1,2Scanbereplacedwith Ratthedenominator
of(126).Attheexponentofthenumeratorcaremustbetakentoevaluatethedifference R−r1Swhichisof
thesameorderas r2;inthe3-DCartesianreferenceframe
s−r1=(Rsin𝜃cos𝜑−r1sin𝜃1,Rsin𝜃sin𝜑,Rcos𝜃−r1cos𝜃1), (127)
andconsequently
r1S=|s−r1|=[R2+r2
1−2Rr1(sin𝜃cos𝜑sin𝜃1+cos𝜃cos𝜃1)]1
2. (128)
Equation(128)impliesthat
r1S≈R−r1(sin𝜃cos𝜑sin𝜃1+cos𝜃cos𝜃1) (129)
tofirstorderin r1∕R[Sniederetal. ,2008].Substitutingintoequation(126),
I2(𝜔)=1
(4𝜋)41
(
c√
2𝜋)3
R2r2ei𝜔
cr2
∫𝜋
0d𝜃sin𝜃∫𝜋
−𝜋d𝜑ei𝜔
cr1(sin𝜃cos𝜑sin𝜃1+cos𝜃cos𝜃1)h∗(̂r2,̂s(𝜃,𝜑)).(130)
The integral at the right-hand side of equation (130) coincides with that in equation (A3), after replacing
𝜆=𝜔,x=𝜃,y=𝜑,f(𝜃,𝜑)=h∗(̂r2,̂s(𝜃,𝜑))sin𝜃, and𝜓(𝜃,𝜑)=r1(sin𝜃cos𝜑sin𝜃1+cos𝜃cos𝜃1)∕c; the
stationary-phaseformula(A2)canthenbeapplied,providedthat hbeasmoothfunctionof 𝜃and𝜑.
In analogy with the procedure of section 5.2, we differentiate 𝜓(𝜃,𝜑)with respect to 𝜃and𝜑to find the
stationarypointsof(130),
𝜓𝜃=r1
c(cos𝜃cos𝜑sin𝜃1−sin𝜃cos𝜃1), (131)
𝜓𝜑=−r1
csin𝜃sin𝜑sin𝜃1. (132)
Equation (132) establishes that stationary points can be found at either 𝜑=0or𝜑=𝜋, i.e., on the plane
wherethescattererandbothreceiversare.Substituting cos𝜙=±1intoequation(131),wefurtheridentify
thetwostationarypoints( 𝜃=𝜃1,𝜑=0)and(𝜃=𝜋−𝜃1,𝜑=𝜋).Ifwecontinuedifferentiating,
𝜓𝜃𝜃=−r1
c(sin𝜃cos𝜑sin𝜃1+cos𝜃cos𝜃1), (133)
𝜓𝜑𝜑=−r1
csin𝜃cos𝜑sin𝜃1, (134)
𝜓𝜃𝜑=r1
csin𝜃sin𝜑sin𝜃1. (135)
BOSCHIANDWEEMSTRA AMBIENT-NOISECROSSCORRELATION 437

ReviewsofGeophysics 10.1002/2014RG000455
Atthestationarypoints,
𝜓(𝜃1,0)=r1
c, (136)
𝜓(𝜋−𝜃1,𝜋)=−r1
c, (137)
𝜓𝜃𝜃(𝜃1,0)=−r1
c, (138)
𝜓𝜃𝜃(𝜋−𝜃1,𝜋)=+r1
c, (139)
𝜓𝜑𝜑(𝜃1,0)=−r1
csin2𝜃1, (140)
𝜓𝜑𝜑(𝜋−𝜃1,𝜋)=r1
csin2𝜃1, (141)
𝜓𝜃𝜑(𝜃1,0)=𝜓𝜃𝜑(𝜋−𝜃1,𝜋)=0. (142)
Aftersubstituting(136)–(142)intothestationary-phaseformula(A12),
I2(𝜔)≈1
(4𝜋)41√
2𝜋c2R2ei𝜔
cr2
r2i
𝜔[
h∗(̂r2,̂s(𝜋−𝜃1,𝜋))e−i𝜔
cr1
r1−h∗(̂r2,̂s(𝜃1,0))ei𝜔
cr1
r1]
.(143)
Theintegralin I3(𝜔)issimilartothatin I2(𝜔),andthestationary-phaseapproximationleadsto
I3(𝜔)≈1
(4𝜋)41√
2𝜋c2R2e−i𝜔
cr1
r1i
𝜔[
h(̂r1,̂s(𝜃2,0))e−i𝜔
cr2
r2−h(̂r1,̂s(𝜋−𝜃2,𝜋))ei𝜔
cr2
r2]
. (144)
Theintegralin I4(𝜔),
I4(𝜔)=1
(4𝜋)51
(2𝜋)2c4ei𝜔
c(r2−r1)
R2r2r1∫𝜋
0d𝜃sin𝜃∫𝜋
−𝜋d𝜑h(̂r1,̂s(𝜃,𝜑))h∗(̂r2,̂s(𝜃,𝜑)), (145)
isnotastationary-phaseintegral,ingeneral.
Noticethatbyequation(E22)andsince ̂s(𝜃1,0)=̂r2,theterm
ei𝜔
cr2
r2ei𝜔
cr1
r1h∗(̂r1,̂s(𝜃1,0)) =25𝜋3c2G∗
3D(r2,𝜔)G∗
3D(r1,𝜔)h∗(̂r1,̂r2), (146)
appearing in (143), describes (the complex conjugate of) an impulse propagating from R2toR1via the
scatterer.Likewise,inequation(144),
e−i𝜔
cr1
r1e−i𝜔
cr2
r2h(̂r2,̂s(𝜃2,0)) =25𝜋3c2G3D(r2,𝜔)G3D(r1,𝜔)h(̂r2,̂r1) (147)
BOSCHIANDWEEMSTRA AMBIENT-NOISECROSSCORRELATION 438

ReviewsofGeophysics 10.1002/2014RG000455
isanimpulsepropagatingfrom R1toR2viathescatterer.Theremainingtermsin(143),(144),and(145)donot
haveanimmediatephysicalexplanation.Thesource-averagedcrosscorrelation,coincidingwiththeintegral
of(124),takestheform
IS(𝜔)+I2(𝜔)+I3(𝜔)+I4(𝜔)≈√
2𝜋
(4𝜋R)2i
𝜔[G3D(Δ,𝜔)−G∗
3D(Δ,𝜔)]
−√
2𝜋
(4𝜋R)2i
𝜔G∗
3D(r2,𝜔)G∗
3D(r1,𝜔)h∗(̂r1,̂r2)
+√
2𝜋
(4𝜋R)2i
𝜔G3D(r2,𝜔)G3D(r1,𝜔)h(̂r2,̂r1)
+1
(4𝜋)41√
2𝜋c2i
𝜔ei𝜔
c(r2−r1)
R2r2r1
×{
[h∗(̂r2,−̂r1)−h(̂r1,−̂r2)]−1
4𝜋c2𝜔
(2𝜋)3
2∫𝜋
0d𝜃sin𝜃∫𝜋
−𝜋d𝜑h(̂r1,̂s(𝜃,𝜑))h∗(̂r2,̂s(𝜃,𝜑))}(148)
wherewehaveusedequation(49)with nS=1andtheidentities ̂s(𝜋−𝜃1,𝜋)=−r1and̂s(𝜋−𝜃2,𝜋)=−r2.
Equation (148) was first derived, usin g a slightly different convention/notation, by Sniederetal. [2008], who
observedthatthetermin {…}iszeroasaconsequenceofthegeneralizedopticaltheorem.Weareleftwith
IS(𝜔)+I2(𝜔)+I3(𝜔)+I4(𝜔)≈√
2𝜋
(4𝜋R)2i
𝜔[
GS
3D(Δ,𝜔)−GS∗
3D(Δ,𝜔)]
≈−√
2𝜋
8(𝜋R)21
𝜔ℑ[GS
3D(Δ,𝜔)],(149)
similartoequations(51)and(103).Theprocedureofsection5.3.3couldbeappliedtoshowthatinthetime
domain, the source-averaged cross correlation (149) is proportional to the sum of GS
3D(t)(causal part) and
−GS
3D(−t)(anticausal).
Thisresultconfirmsthat,asfirstpointedoutby WeaverandLobkis [2004],diffuse-fieldcrosscorrelationinhet-
erogeneous(ratherthanjusthomogeneous)mediaallowsinprincipletoreconstructthefullGreen’sfunctionofthemedium,withallreflectionsandscatteringsandpropagationmodes.Thisisimplicitinthereciprocitytheoremformulationofsection6.2andhasbeenverifiedexperimentallyandnumericallyby,e.g., Laroseetal.
[2006],Mikeselletal. [2012],and Colombietal. [2014].
8.1. SpuriousArrivalsandTheirCancelation
Theresult(149)mightbesurprisingifoneconsidersthescattererasa(“secondary”)source;nomatterwhere
theactual(“primary”)sourceis,thescattererisalwaysatthesamelocationrelativeto
R1andR2sothatthe
delaybetweenthearrivalofthescatteredsignalat R1anditsarrivalat R2isalwaysthesame.Thiswouldgive
risetoapeakinthestackedcrosscorrelationthatdoesnotcorrespondwitheitherofthetwoarrivalsin GS
3D
asdefinedbyequation(123).
This is illustrated in Figures 13 and 14. Following Sniederetal. [2008], we convolve a Ricker wavelet (central
frequency of 1 Hz) with the Green’s function associated with a uniformly dense, circular source distribution(Figure13)plusasinglescatterer.Forsimplicityandinanalogywith Sniederetal. [2008],weworkintwodimen-
sions(membranewavesfrompointsources)andassumeisotropicscattering.Weimplementtheexpression
for
hderivedby GroenenboomandSnieder [1995]viatheopticaltheorem[ Newton,2002].Foreachlocationof
thesource,wecrosscorrelatethecorrespondingsignalsat R1andR2,whichwechoosetobeequidistantfrom
thescatterer,andplotthecrosscorrelationsinFigure14a.Asanticipated,whateverthesourceazimuth,apeakinthecrosscorrelationappearsat
t=0,correspondingtothescatteredsignalhittingsimultaneously R1and
R2.Thispeak,or“arrival,”at t=0clearlydoesnotexistin GS
3D,whichwouldseemtocontradictequation(149).
However,whenstackingallsingle-azimuthcrosscorrelations,the t=0peakcancelsoutwiththe“knees”of
someothercrosscorrelationpeaks(Figure14b),andequation(149)isindeedconfirmed.Thefourstationarypointsidentifiedinsection8appearaskneesofthecross-correlationpeaksinFigure14a,correspondinginturntothetwocausalandtwoanticausalpeaksofthesystem’sRickerresponse(Figure14b).
It is critical that source illumination be uniform, for the spurious term associated with the scatterer to
disappear.Toemphasizethispoint,weshowinFigure14cthestacksobtainedbyconsideringonlyhalfofthe
BOSCHIANDWEEMSTRA AMBIENT-NOISECROSSCORRELATION 439

ReviewsofGeophysics 10.1002/2014RG000455
Figure13. Distributionofsources(stars),receivers(triangles),andapointscatterer(redcircle)inthesetupof
section8.1.AsinFigure11,only10%ofthesimulatedsourcesareshown.
availablesources[ Sniederetal. ,2008],oneithersideofthereceiverarray.Dependingonwhichside,onlythe
causaloranticausalpartoftheRickerresponseisreconstructed,asinFigure12.Furthermore,threeartifacts
emerge;theoneat t=0clearlyresultsfromthespuriouswaveletinFigure 14a,whichdoesnotcancelout
since the illumination is (strongly) nonuniform. The other two are simply the spurious arrivals of Figure 12,whichresultfromthesharpinhomogeneitiesinsourcedistribution,andwouldbefoundevenintheabsenceofscatterers.
Scattering thus has a complex effect on Green’s function reconstruction. Each scatterer further complicates
theGreen’sfunction,introducinganadditional,physicaltermtobereconstructed.Italsohindersitsretrieval,generating a spurious term which will only cancel out if the wavefield is sufficiently diffuse. On the otherhand,scatterersthemselvescontributetothewavefield’sdiffusivityandazimuthaluniformityofillumination.Recent and current work [e.g., Fleuryetal. , 2010;Mikeselletal. , 2012;RavasiandCurtis , 2013;Colombietal. ,
2014] aims at disentagling the specific role of scattered ambient signal in reconstructing the main peaks of
theGreen’sfunctionaswellasits“coda.”
9. Summary
Inadiffuseambientwavefield,i.e.,arandomwavefieldwhereenergypropagateswithequalprobabilityinall
directions,crosscorrelatingthesignalrecordedbyapairofreceiversisawaytomeasuretheimpulseresponse
(Green’sfunction)betweenthereceivers.Inpractice,longrecordingsofseismicambientnoiseoftenincludewavestravelingalongapproximatelyallazimuths:thecombinationoftheircrosscorrelationsapproximatesthatofadiffusefield.TherelationshipbetweencrosscorrelationandGreen’sfunctionchangesdependingonsomepropertiesofthemediumandofthewavefield:
(i)Inahomogeneous,lossless,unboundedmembranewherecircularwavesaregeneratedbypointsources,
thecrosscorrelationofdiffusenoiseisgenerallynotexplicitlyrelatedtothecorrespondingGreen’sfunc-tion
G2D(section 5.3). If, however, the source distribution is symmetric along the receiver-receiver axis,
equation(70)stipulatesthat G2Disproportionaltothetimederivativeofthecrosscorrelation.Themem-
branesetupisrelevanttoseismology,becauseitcorr espondsinpracticetothepropagationofRayleigh
waves(section2.3)onEarth.
(ii)In homogeneous, lossless, unbounded 3-D acoustic media, the Green’s function can be exactly recon-
structed from diffuse noise generated at point sources that are smoothly distributed over a sphere
surroundingthereceivers.Moreprecisely,thecrosscorrelationcoincideswiththetimederivativeoftheGreen’sfunction(section5.2).Thisrelationshipholdsforbothcausalandanticausalcontributionstothe
BOSCHIANDWEEMSTRA AMBIENT-NOISECROSSCORRELATION 440

ReviewsofGeophysics 10.1002/2014RG000455
Figure14. Crosscorrelationsassociatedwithacircular,planardistributionofsourcessurroundingthetworeceivers,and
oneisotropicpointscatterer,assketchedinFigure13.EachsourcegeneratesaRickerwavelet(centralfrequency=1Hz);thewaveletisconvolvedwiththeGreen’sfunction(phasevelocity
c=2Km/s),includingthecontributionofscattering,
toevaluatethesignalsobservedatthetworeceivers,whicharethencrosscorrelated.(a)Single-sourcecrosscorrelationsforallsourceazimuths; thecolorscalesaturatesat20%ofthemaximumamplitudetoemphasizethespuriousarrivals;thedashedlinesmarktheazimuthswheretheyellow/“green”sourcedistributionisdiscontinuous(Figure13).(b)Source-averagedcrosscorrelationresultingfromFigure14a.(c)Source-averagecrosscorrelationsthatonewouldobtainifsignalwasgeneratedonlyatthelocationsidentifiedbygreenstarsinFigure13(greenline)versustheyellowstars(yellowline).
Green’sfunction/crosscorrelationandremainsvalidifthesourcedistributionisasymmetric(providedthat
itbesmooth).
(iii)The latter result does not hold if noise in free space is generated by point sources uniformly distributed
along acirclesurrounding the receivers (section 5.1). Still, even in this case, the time of maximum cross
correlationclearlycoincideswiththatof G2D’speak(Figure5).
(iv)Ifnoiseismadeupof planewaves propagatinginalldirections(section5.4),thesameformulationholdsfor
2-Dand3-Dunbounded,losslessmedia.Inthiscase,thesource-averagedcrosscorrelationisproportional,in the frequency domain, to the real part of the 2-D Green’s function
G2D, which in turn is proportional
totheBesselfunction J0(𝜔r∕c),withrtheinterreceiverdistance.Henceequation(72).Thecorresponding
BOSCHIANDWEEMSTRA AMBIENT-NOISECROSSCORRELATION 441

ReviewsofGeophysics 10.1002/2014RG000455
time domain relationship was derived by Nakahara [2006], who showed that the Hilbert transform of
diffusefieldcrosscorrelationcoincideswiththe(causalminusanticausal) G2D.
(v)In bounded, heterogeneous, attenuating media, a general relationship between source-averaged
cross correlation and Green’s function is found via the reciprocity theorem (section 6.2). Essentially,
source-averagedcrosscorrelationcoincideswiththesumofavolume( V)integralandanintegraloverthe
boundary 𝜕Vof such volume. It follows that densely distributed noise sources throughoutavolume are
neededfortheGreen’sfunctiontobereconstructedviacrosscorrelation;sourceswithin Vareunnecessary
ifthemediumislossless,providedthatall azimuths areilluminated.
The results (i) and (iv) are relevant for most seismology applications seen so far, since, provided that atten-
uation can be neglected, they apply to surface waves in the “membrane” approximation [ Tanimoto , 1990]
(section2.1).Forexample Ekströmetal. [2009]haveshownexplicitlythatthephaseofthetheoreticalGreen’s
function is in agreement with cross-correla ted ambient noise data, thus validating the “lossless-medium”
approximation for phase/group-velocity -based seismic imaging. Reconstructing the Green’s function’s
amplitude remainsproblematic[e.g., Weemstraetal. ,2014].
AppendixA:TheStationary-PhaseApproximation
A1. One-DimensionalIntegrals
Thestationary-phaseapproximation[e.g., BenderandOrszag ,1978]appliestointegralsoftheform
I(𝜆)=∫b
af(x)ei𝜆𝜓(x)dx, (A1)
whereaandbarearbitraryconstantsand 𝜆islargeenoughfor ei𝜆𝜓(x)tobearapidlyoscillatingfunctionof x,
withrespecttothesmoothfunctions f(x),𝜓(x).Itisbasedonthefinding[e.g., BenderandOrszag ,1978,section
6.5] that the leading contribution to I(𝜆)comes from a small interval surrounding the “stationary points” of
𝜓(x), i.e., the locations xsuch that 𝜓′(x)=0(where the “prime” denotes differentiation with respect to x).
Elsewhere, the integrand oscillates quickly and its average contribution to the integral is negligible. Bender
andOrszag [1978]demonstratethatforasinglestationarypointat a,
I(𝜆)≈f(a)ei(
𝜆𝜓(a)±𝜋
4)√𝜋
2𝜆|𝜓′′(a)|, (A2)
valid for𝜆−→∞and in the assumption that 𝜓′′(a)≠0. The sign of 𝜋∕4at the exponent of eis positive
if𝜓′′(a)>0andnegativeotherwise.Sinceanyintegraloftheform(A1)canbewrittenasasumofintegrals
withoneoftheintegrationlimitscoincidingwithastationar ypoint,equation(A2)issufficienttosolveall1-D
stationary-phase integrals like (A1), regardless of the number and location of stationary points [ Benderand
Orszag,1978].
A2. Extensionto2-DIntegrals
TheresultofsectionA1canbeusedtofindageneralapproximateformulaforintegralsoftheform
I(𝜆)=∫b
adx∫d
cdyf(x,y)ei𝜆𝜓(x,y). (A3)
If we, again, only consider the limit 𝜆−→∞, the integrand in (A3) turns out to be very strongly oscillatory,
andtheonlynonnegligiblecontributiontotheintegralcomesfromthevicinityofthestationarypoints (xi,yi),
definedby (𝜕𝜓
𝜕x(xi,yi),𝜕𝜓
𝜕y(xi,yi))
=(0,0) (A4)
[e.g.,Wong,1986].Equation(A3)canthusbeapproximatedby
I(𝜆)≈∑
if(xi,yi)∫xi+𝜀
xi−𝜀dx∫yi+𝛿
yi−𝛿dyei𝜆𝜓(x,y), (A5)
wherethesumisextendedtoallstationarypoints i,and𝜀and𝛿aresmall. f(x,y)isapproximatedby f(xi,yi)
sincefvariesmuchmoreslowlythan ei𝜆𝜓when𝜆islarge.WenextconductaTaylorseriesexpansionof 𝜓(x,y)
around(xi,yi),
𝜓(x,y)≈𝜓(xi,yi)+1
2[𝜓xxi(x−xi)2+𝜓yyi(y−yi)2+2𝜓xyi(x−xi)(y−yi)], (A6)
BOSCHIANDWEEMSTRA AMBIENT-NOISECROSSCORRELATION 442

ReviewsofGeophysics 10.1002/2014RG000455
where𝜓xxistands for𝜕2𝜓
𝜕x2(xi,yi), and so on. Substituting (A6) into (A5) and after the changes of variable
ui=x−xiandvi=y−yi,wefind
I(𝜆)≈∑
if(xi,yi)ei𝜆𝜓(xi,yi)
∫+𝜀
−𝜀dui∫+𝛿
−𝛿dviei𝜆
2(𝜓xxu2
i+𝜓yyv2
i+2𝜓xyuivi). (A7)
Theintegralin(A7)canbesolvedanalyticallybyfirstrewriting
𝜓xxiu2
i+𝜓yyiv2
i+2𝜓xyiuivi=𝜓xxi(
ui+𝜓xyi
𝜓xxivi)2
+v2
i(
𝜓yyi−𝜓2
xyi
𝜓xxi)
, (A8)
andsubstitutinginto(A7),
I(𝜆)≈∑
if(xi,yi)ei𝜆𝜓(xi,yi)
∫+𝛿
−𝛿dviei𝜆
2v2
i(𝜓yyi−𝜓2
xyi
𝜓xxi)
∫+𝜀
−𝜀duiei𝜆
2𝜓xxi(
ui+𝜓xyi
𝜓xxivi)2
. (A9)
Since𝜆−→∞,bothintegralsattheright-handsideare1-Dstationary-phaseintegralsasseeninsectionA1.
Theonlystationarypointinthe viintegrationdomainisat vi=0;thephaseterminthe uiintegralislikewise
stationaryat ui=0(providedthat vi≈0).Applicationofequation(A2)tothe uiintegralin(A9)thengives
∫+𝜀
−𝜀duiei𝜆
2𝜓xxi(
ui+𝜓xyi
𝜓xxivi)2
=e±i𝜋
4√
2𝜋
𝜆|𝜓xxi|, (A10)
wherewehaveassumed vi≈0sincethisintegralistobeevaluatednearthestationarypointofthe viintegral,
andwehaveimplicitlymultipliedby2sinceequation(A2)isvalidforasinglestationarypointlocatedatoneoftheintegrationlimits.Similarly,
∫+𝛿
−𝛿dviei𝜆
2v2
i(𝜓yyi−𝜓2
xyi
𝜓xxi)=e±i𝜋
4√√√√√2𝜋
𝜆|𝜓yyi−𝜓2
xyi
𝜓xxi|. (A11)
Substitutingequations(A10)and(A11)into(A9),weareleftwiththefinalformula
I(𝜆)≈2𝜋
𝜆∑
if(xi,yi)ei[
𝜆𝜓(xi,yi)±𝜋
4±𝜋
4]√√√√√1
|𝜓xxi||𝜓yyi−𝜓2
xyi
𝜓xxi|, (A12)
valid for𝜆−→∞and in the assumption that 𝜓xxi≠0,𝜓yyi−𝜓2
xyi
𝜓xxi≠0. The signs of the 𝜋∕4terms in (A12)
dependonthoseof 𝜓xxiand𝜓yyi−𝜓2
xyi∕𝜓xxi(seesectionA1)[ BenderandOrszag ,1978].
AppendixB:FourierTransformConvention
InverseandforwardFouriertransformationsareappliedfrequentlythroughoutthisstudy.Wedenote ℱand
ℱ−1theforwardandinverseFouriertransformoperators,respectively. ℱandℱ−1canbedefinedinvarious
ways.Weadoptherethefollowingconvention:foranyfunction f,
f(𝜔)=ℱ[f(t)] =1√
2𝜋∫+∞
−∞dtf(t)e−i𝜔t, (B1)
andconsequently
f(t)=ℱ−1[f(𝜔)] =1√
2𝜋∫+∞
−∞d𝜔f(𝜔)ei𝜔t, (B2)
wherewehavechosentosimplyspecifytheargument(time torfrequency 𝜔,respectively)todistinguisha
timedomainfunction f(t)fromitsFouriertransform f(𝜔).
BOSCHIANDWEEMSTRA AMBIENT-NOISECROSSCORRELATION 443

ReviewsofGeophysics 10.1002/2014RG000455
Some properties of the Fourier transform are particularly useful in the context of ambient noise cross
correlation. First of all, it follows from (B1) that the Fourier transform of the derivative of fwith respect
totis
ℱ[df(t)
dt]
=i𝜔f(𝜔). (B3)
TheFouriertransformoftheintegralof fis
ℱ[
∫t
−∞f(𝜏)d𝜏]
=−i
𝜔f(𝜔), (B4)
providedthat f(t)−→0whent−→−∞.
OtherusefulpropertiesoftheFouriertransformconcernevenandoddfunctions.TheGreen’sfunctions G(t)
weworkwiththroughoutthisstudy(AppendixE)arereal(inthetimedomain)andhavetheproperty G(t)=0
ift<0.Letusdefinethereal,evenfunction
Ge(t)=1
2G(t)+1
2G(−t), (B5)
andthereal,oddfunction
Go(t)=1
2G(t)−1
2G(−t). (B6)
Thedefinitions(B5)and(B6)imply
G(t)=Ge(t)+Go(t). (B7)
Itfollowsfrom(B1)thattheFouriertransformofarealevenfunctionisevenandpurelyreal,whiletheFourier
transformofarealoddfunctionisoddandpurelyimaginary.Then,
ℜ[G(𝜔)]=Ge(𝜔) (B8)
and
ℑ[G(𝜔)]=−iGo(𝜔), (B9)
whereℜ[G(𝜔)]andℑ[G(𝜔)]denotetherealandimaginarypartof G(𝜔),respectively.
AppendixC:BesselFunctions
Bessel functions emerge frequently in noise literature, starting with the early works of, e.g., Eckart[1953],
Aki[1957], and Cox[1973]. Their mathematical properties are described in detail by AbramowitzandStegun
[1964]. In seismic interferometry we are in practice only interested in zeroth-order Bessel functions of the
firstandsecondkind,denoted J0(•),Y0(•),respectively,whichtogetherwiththeHankelfunctions H(1)
0(•) =
J0(•)+iY0(•)andH(2)
0(•) =J0(•)−iY0(•)canbedefinedasthesolutionsofthezeroth-orderBesselequation
d2f(x)
dx2+1
xdf(x)
dx+f(x)=0 (C1)
[AbramowitzandStegun ,1964,equation(9.1.1)].
In our implementation we employ the far-field (large r) and/or high-frequency (large 𝜔) approximations for
J0(x),Y0(x),andH(1)
0(x),H(2)
0(x),namely,
J0(x)≈√
2
𝜋xcos(
x−𝜋
4)
, (C2)
Y0(x)≈√
2
𝜋xsin(
x−𝜋
4)
, (C3)
BOSCHIANDWEEMSTRA AMBIENT-NOISECROSSCORRELATION 444

ReviewsofGeophysics 10.1002/2014RG000455
H(1)
0(x)≈√
2
𝜋xei(x−𝜋∕4), (C4)
H(2)
0(x)≈√
2
𝜋xe−i(x−𝜋∕4)(C5)
[AbramowitzandStegun ,1964,equations(9.2.1)–(9.2.4)].Inpractice, J0andY0canroughlybethoughtofas
“damped”sinusoidalfunctions,whoseamplitudedecaysexponentiallyastheirargumenttendstoinfinity;atlarge
x,J0approximatesacosinewitha 𝜋∕4phaseshift,and Y0approximatesasinewitha 𝜋∕4phaseshift.
AppendixD:CancelationofCrossTerms
Throughout this study, simple analytical formulae for the cross correlation of an ambient wavefield are
obtainedneglectingthecontributionoftheso-called“crossterms,”i.e.,thereceiver-receivercrosscorrelationof signal generated by a couple of different sources. We propose in the following a simple proof, similar to
thatofWeemstraetal. [2014],ofthevalidityofthisassumption.
Thepressureduetoanimpulseemittedbysource
jatarandomtime tjisgivenbyequation(E21),afterreplac-
ingt−x∕cwitht−tj−x∕casargumentoftheDiracfunction;inthefrequencydomain,thisisequivalentto
adding a phase anomaly −i𝜔tjto the argument of the exponential in equation (E22). If NSsuch sources are
activeduringthetimeintervaloverwhichacrosscorrelationisconducted,thepressure pirecordedatreceiver
iisalinearcombinationofthesignalsoriginatingfromthesesources,i.e.,in3-D,
pi(𝜔)=NS∑
j=1e−i(𝜔rij
c+𝜔tj)
√
2𝜋4𝜋crij, (D1)
whererijisthedistancebetweenreceiver iandsource j.
Crosscorrelatingthesignalrecordedat R1withthatrecordedat R2,wefind
p1(𝜔)p∗
2(𝜔)=NS∑
j=1NS∑
k=1e−i𝜔
c(r1j−r2k)
2𝜋(4𝜋c)2r1jr2ke−i𝜔(tj−tk), (D2)
andthecross-termcontributioncanbeisolatedbyseparating j≠kterms(crossterms)from j=kones:
p1(𝜔)p∗
2(𝜔)=NS∑
j=1e−i𝜔
c(r1j−r2j)
r1jr2j+NS∑
j=1∑
j≠ke−i𝜔
c(r1j−r2k)
r1jr2ke−i𝜔(tj−tk). (D3)
The first term at the right-hand side of equation (D3) is a sum of single-source cross correlations in 3-D
(equation(6)).(Thealgebraissimilarin2-D(equation(7)), ifthefar-fieldapproximationisapplied.)Itsaverage
is given, e.g., by equations (8)–(10) above, which are obtained from (6) or (7) by neglecting the cross terms
j≠kandreplacingthesumoversources jwithanintegralovertheareaorvolumeoccupiedbythesources.
Whilej≠ktermsin(D3)arenonnegligible,wenextshowthattheircontributiontotheaverageof(D3)over
alargesetofsourcesisnegligible.Wemaketheassumptionsthat(i)ateachrealization,thevalues tjchange
randomlyand(ii)theyareuniformlydistributedbetween0and Ta,i.e.,0<tj<Ta,whereTa=2𝜋∕𝜔.Letus
introducethephase 𝜙j=𝜔tj,whichisrandomlydistributedbetween0and 2𝜋.Theexponent −i𝜔(tj−tk)is
accordinglyreplacedby −i(𝜙j−𝜙k).Intheprocessofaveraging,impulseswillbegeneratedmultipletimesat
eachsourcelocation j,resultinginrandomphases 𝜙j.Thisisequivalenttorequiring,asitisusuallydoneinthe
literature,thatnoisesourcesbespatially andtemporallyuncorrelated[e.g., Snieder,2004;Rouxetal. ,2005].To
takeintoaccounttheeffectsofrandomvariationsinphase,the j≠ktermmustconsequentlyinvolve(besides
theusualsumorintegraloversources)anintegraloverallpossiblevalues(0to 2𝜋)ofeachsourcephase,
Av[NS∑
j=1∑
j≠kei𝜔
c(r1j−r2k)
2𝜋r1jr2kei(𝜙j−𝜙k)]
=1
(2𝜋)NS∫∫…∫2𝜋
0NS∑
j=1∑
j≠kei𝜔
c(r1j−r2k)
2𝜋r1jr2kei(𝜙j−𝜙k)d𝜙1d𝜙2…d𝜙NS.(D4)
AllIntegralsevaluatetozerobecausetheintegrandstraverseacircleinthecomplexplanefrom0to 2𝜋.Con-
sequently,ifcrosscorrelationsareaveragedover sufficientrealizations,thecrosscorrelationinequation(D3)
reducestoasumofsingle-sourcecrosscorrelations(noncrossterms).
BOSCHIANDWEEMSTRA AMBIENT-NOISECROSSCORRELATION 445

ReviewsofGeophysics 10.1002/2014RG000455
The averaging procedure we just described is often referred to as “ensemble-averaging” in ambient noise
theory. This expression is borrowed from statistical mechanics: an “ensemble” is a set of states of a system,
eachdescribedbythesamesetofmicroscopicforcesandsharingsomecommonmacroscopicproperty.The
ensembleconceptthenstatesthatmacroscopicobservablescanbecalculatedbyperformingaveragesoverthe states in the ensemble [ Tuckerman , 1987]. In our case, one state consists of the same acoustic or elastic
medium (its response to a given impulse is always the same) being illuminated by one or more randomly
located sources with random phases: our ensemble aver age is the average over all possible combination
of source locations and phases. Cross-correlating recordings associated with a single state (a unique com-
bination of sources) yield many cross terms; we have just shown, however, that ensemble-averaging cross
correlationsoverdifferentstatesimpliesthatthesecrosstermsstackincoherentlyandhencebecomenegli-
gible.Intherealworld,thelocationandphaseofnoisesources(e.g.,oceanmicroseisms)arenotwellknown;
theassumptionismadethatovertime,asufficientlydiverserangeofsourcesissampledsothataveragingthe
crosscorrelationovertime(i.e.,computingthecrosscorrelationoveraverylongtimewindow)ispractically
equivalenttoensemble-averaging.
AppendixE:Green’sFunctionsoftheScalarWaveEquation(Homogeneous
LosslessMedia)
E1. Green’sProblemasHomogeneousEquation
Following,e.g., Rouxetal. [2005],Sanchez-SesmaandCampillo [2006],and Harmonetal. [2008],wecallGreen’s
functionG=G(x,xS,t)thesolutionof
∇2G−1
c2𝜕2G
𝜕t2=0 (E1)
withinitialconditions
G(x,xS,0)=0, (E2)
𝜕G
𝜕t(x,xS,0)=𝛿(x−xS), (E3)
i.e., an impulsive source at xS. We are only interested in causal Green’s functions, satisfying the radiation
condition,i.e.,vanishingat t<0.
OnceGis known, it can be used to solve rapidly more general initial value problems associated with (E1).
Consider,forexample,
∇2f−1
c2𝜕2f
𝜕t2=0 (E4)
withthemoregeneralinitialconditions
f(x,0)=0, (E5)
𝜕f
𝜕t(x,0)=h(x). (E6)
Itcanbeprovedbydirectsubstitutionthatif Gsolves(E1)–(E3)then
f(x,t)=∫Rdddx′G(x,x′,t)h(x′) (E7)
solves(E4)–(E6),with ddenotingthenumberofdimensions:2or3inourcase.
For the sake of simplicity, we shall set xS=0in the following. G(x,xS,t)can be recovered from G(x,0,t)by
translationofthereferenceframe.
BOSCHIANDWEEMSTRA AMBIENT-NOISECROSSCORRELATION 446

ReviewsofGeophysics 10.1002/2014RG000455
E2. SolutionintheSpatialFourierTransformDomain
Equation(E1)canbesolvedviaaspatialFouriertransform,i.e.,bythe Ansatz
G(x,t)=1
(2𝜋)d∫RdddkG(k,t)eik⋅x, (E8)
which,substitutedinto(E1),gives
1
c2𝜕2G
𝜕t2(k,t)+k2G(k,t)=0, (E9)
withk=|k|.Equation(E9)issolvedby
G(k,t)=A(k)cos(kct)+B(k)sin(kct), (E10)
whereA(k)andB(k)mustbedeterminedbytheinitialconditions.Transforming(E2)and(E3)to kspace,with
𝛿(x)=1
(2𝜋)d∫Rdddkeik⋅x, (E11)
andreplacing Gintheresultingequationswithitsexpression(E10),wefind A(k)=0,B(k)=1
kc,and
G(k,t)=1
kcsin(kct). (E12)
Equation (E8) can now be used to determine G(x,t)from its spatial Fourier transform G(k,t), and the result
differsimportantlydependingon d.
E3. InverseTransformoftheSolutionto2-DSpace
Substituting(E12)into(E8)inthe d=2case,
G2D(x,t)=1
4𝜋2∫+∞
−∞dk1∫+∞
−∞dk2sin(kct)
kceik⋅x. (E13)
Wecall𝜉theanglebetween kandx.Itfollowsthat k⋅x=kxcos𝜉,withx=|x|.Foreach xtheintegration
canbeconductedover kand𝜉,usingdk1dk2=kdkd𝜉,
G2D(x,t)=1
4𝜋2c∫+∞
0dksin(kct)∫2𝜋
0d𝜉eikxcos𝜉
=1
2𝜋2c∫+∞
0dksin(kct)∫𝜋
0d𝜉cos(kxcos𝜉)
=1
2𝜋c∫+∞
0dksin(kct)J0(kx),(E14)
where we have used the symmetry properties of sine and cosine, and the integral form of the zeroth-order
Bessel function of the first kind J0[Abramowitz and Stegun , 1964, equation (9.1.18)] (Appendix C). The
remainingintegralin(E14)issolvedviaequation(11.4.38)of AbramowitzandStegun [1964],resultingin
G2D(x,t)=1
2𝜋c2H(
t−x
c)
√
t2−x2
c2(E15)
[Sanchez-Sesma and Campillo , 2006;Harmon et al. , 2008], where Hdenotes the Heaviside function.
Equation(E15),aswellasothertimedomainformulaefortheGreen’sfunction,isonlyphysicallymeaningful
fort>0.Dimensionalanalysisof(E15)showsthat G2Dinthisformulationhasunitsoftimeoversquareddis-
tance. The Fourier transform of G2D(x,t)is inferred from equation (9.1.24) of AbramowitzandStegun [1964],
afterapplyingourdefinition(B1)toequation(E15):
G2D(x,𝜔)=1√
8𝜋3c2∫∞
x
cdte−i𝜔t
√
t2−x2
c2
=1√
8𝜋3c2∫∞
1dt′e−i𝜔x
ct′
√
t′2−1
=−1
4√
2𝜋c2[
Y0(𝜔x
c)
+iJ0(𝜔x
c)]
=1
4i√
2𝜋c2H(2)
0(𝜔x
c)(E16)
BOSCHIANDWEEMSTRA AMBIENT-NOISECROSSCORRELATION 447

ReviewsofGeophysics 10.1002/2014RG000455
[Sanchez-SesmaandCampillo , 2006;Harmonetal. , 2008], where Y0andH(2)
0denote the zeroth-order Bessel
functionofthesecondkindandthezeroth-orderHankelfunctionofthesecondkind[ AbramowitzandStegun ,
1964](AppendixC).Baseduponthefar-field/high-frequencyasymptoticform(C5)of H(2)
0,wecanalsowrite
G2D(x,𝜔)≈1
4i𝜋c3∕2e−i(𝜔x
c−𝜋
4)
√
𝜔x, (E17)
whichisanalogoustoequation(14)of Snieder[2004].
Equation (9.1.24) of AbramowitzandStegun [1964] is only valid for positive frequency 𝜔>0, and so is, as a
consequence, our relation (E16). We know, however, that G2D(x,t)is a real-valued function: the relationship
G2D(x,𝜔)=G∗
2D(x,−𝜔)thenholdsandallowsustodefine G2Dintheentirefrequencydomain.
E4. InverseTransformoftheSolutionto3-DSpace
Wenowsubstitute(E12)intothe3-Dversionof(E8)andfind
G3D(x,t)=1
8𝜋3∫+∞
−∞dk1∫+∞
−∞dk2∫+∞
−∞dk3sin(kct)
kceik⋅x. (E18)
The integral in (E18) is simplified by switching from k1,k2,k3to spherical coordinates k,𝜉,𝜒, with the 𝜉=0
directioncoincidingwiththatof x.Thendk1dk2dk3=−k2dkd𝜒d(cos𝜉),k⋅x=kxcos𝜉,and
G3D(x,t)=−1
8𝜋3∫+∞
0dkksin(kct)
c ∫2𝜋
0d𝜒∫−1
1d(cos𝜉)eikxcos𝜉
=1
2𝜋2cx∫+∞
0dksin(kct)sin(kx).(E19)
Thekintegralin(E19)issolvedviatheequality sin𝛼sin𝛽=1
2[cos(𝛼−𝛽)−cos(𝛼+𝛽)],whichgives
G3D(x,t)=1
4𝜋2cx∫+∞
0dk{cos[k(ct−x)]−cos[k(ct+x)]}
=1
4𝜋2cxlim
z→∞{[sin[k(ct−x)]
ct−x]k=z
k=0−[sin[k(ct+x)]
ct+x]k=z
k=0}
=1
4𝜋cx[
𝛿(
t−x
c)
−𝛿(
t+x
c)]
,(E20)
where we have used the property of the Dirac 𝛿function that 𝛿(x)=1
𝜋limz→∞sin(zx)
x[e.g.,Weisstein ,
1999–2013].Inourformulationboth tandxarepositivesothat 𝛿(
t+x
c)
=0,andweareleftwith
G3D(x,t)=1
4𝜋c𝛿(
t−x
c)
x(E21)
[e.g.,Aki and Richards , 2002, chap. 4]. Notice that since the dimension of the Dirac 𝛿is the inverse of that
of its argument, that of G3Dis one over squared distance. According to (B1), ℱ[𝛿(t−t0)] =1√
2𝜋e−i𝜔t0, and
consequently,
G3D(x,𝜔)=1√
2𝜋1
4𝜋ce−i𝜔x
c
x. (E22)
E5. Green’sProblemasInhomogeneousEquation
TheGreen’sproblemisalsooftenwrittenasequation(1)plusanimpulsiveforcingterm,i.e.,
∇2g−1
c2𝜕2g
𝜕t2=𝛿(x)𝛿(t), (E23)
withinitialconditions
g(x,0)=0, (E24)
𝜕g
𝜕t(x,0)=0. (E25)
BOSCHIANDWEEMSTRA AMBIENT-NOISECROSSCORRELATION 448

ReviewsofGeophysics 10.1002/2014RG000455
Thesolution gtothisproblemisrelatedtothesolution Gof(E1)–(E3),namely,
g(x,t)=∫t
0dsG(x,t−s). (E26)
Wedemonstratethat(E26)solves(E23)–(E25)viaDuhamel’sprinciple[e.g., Hildebrand ,1976;Strauss,2008]:
letusconsiderthemoregeneralcase
∇2u−1
c2𝜕2u
𝜕t2=h(x,t), (E27)
u(x,0)=0, (E28)
𝜕u
𝜕t(x,0)=0, (E29)
withhan arbitrary forcing term. Suppose that a solution v(x,t)to the following homogeneous problem,
closelyrelatedto(E27)–(E29),canbefound:
∇2v−1
c2𝜕2v
𝜕t2=0, (E30)
v(x,0;s)=0, (E31)
𝜕v
𝜕t(x,0;s)=h(x,s). (E32)
(Inpractice,asolutionofequation(E30)thatsatisfiestheinitialcondition(E32)mustbedeterminedforany
possiblevalueoftheparameter sin(E32);sreplacestintheexpressionof hfirstencounteredin(E27).)Then,
thefollowingrelationbetween uandvholds
u(x,t)=∫t
0dsv(x,t−s;s). (E33)
One canverify that (E33) solves (E30)–(E32)by direct substitution, applying Leibniz’srule for differentiating
undertheintegralsign.Equation(E26)isaparticularcaseof(E33).
Result (E33) can also be combined with (E7) to write the solution of the general inhomogeneous problem
(E27)intermsoftheGreen’sfunction G(e.g.,G2DorG3D).Replacing vinequation(E33)withexpression(E7),
u(x,t)=∫t
0ds∫Rdddx′G(x,x′,t−s;s)h(x′,s)(d=2,3). (E34)
References
Abramowitz,M.,andI.A.Stegun(1964), HandbookofMathematicalFunctions ,Natl.Bur.ofStand.–Appl.Math.Ser.,Washington,D.C.
Aki,K.(1957),Spaceandtimespectraofstationarywaveswithspecialreferencetomicrotremors, Bull.EarthquakeRes.Inst.Univ.Tokyo ,35,
415–456.
Aki,K.,andP.G.Richards(2002), QuantitativeSeismology ,2nded.,Univ.ScienceBooks,Sausalito,Calif.
Bender,C.M.,andS.A.Orszag(1978), AdvancedMathematicalMethodsforScientistsandEngineers ,McGraw-Hill,NewYork.
Bensen, G. D., M. H. Ritzwoller, M. P. Barmin, A. L. Levshin, F. Lin, M. P. Moschetti, N. M. Shapiro, and Y. Yang (2007), Processing
seismic ambient noise data to obtain reliable broad-band surface wave dispersion measurements, Geophys. J. Int. ,169,
1239–1260,doi:10.1111/j.1365-246X.2007.03374.x.
Boschi, L., C. Weemstra, J. Verbeke, G. Ekström, A. Zunino, and D. Giardini (2013), On measuring surface wave phase velocity from
station-stationcross-correlationofambientsignal, Geophys.J.Int. ,192,346–358,doi:10.1093/gji/ggs023.
Brenguier,F.,D.Clarke,Y.Aoki,N.M.Shapiro,M.Campillo,andV.Ferrazzini(2011),Monitoringvolcanoesusingseismicnoisecorrelations,
C.R.Geosci. ,343,633–638.
Campillo,M.,andA.Paul(2003),Long-rangecorrelationsinthediffuseseismiccoda, Science,299,547–549.Acknowledgments
This study greatly benefitted from
our discussions with M. Campillo,L. Campos-Venuti, A. Colombi,
A. Curtis, J. de Rosny, J. P. Montagner,
P. Poli, R. Snieder, P. Roux, J. Verbeke,K. Wapenaar, and four anonymousreviewers. Data and software used
to produce the results of this paper
can be obtained by e-mailing the firstauthor: lapo.boschi@upmc.fr.
The Editor on this paper was Mark
Moldwin. He thanks four anony-mous reviewers for their assistance in
reviewing this manuscript.
BOSCHIANDWEEMSTRA AMBIENT-NOISECROSSCORRELATION 449

ReviewsofGeophysics 10.1002/2014RG000455
Campillo,M.,andP.Roux(2014),Seismicimagingandmonitoringwithambientnoisecorrelations,in TreatiseofGeophysics ,vol.1,editedby
B.RomanowiczandA.M.Dziewonski,Elsevier,Amsterdam.
Claerbout,J.F.(1968),Synthesisofalayeredmediumfromitsacoustictransmissionresponse, Geophysics ,33,264–269.
Colombi,A.,L.Boschi,P.Roux,andM.Campillo(2014),Green’sfunctionretrievalthroughcross-correlationsinatwo-dimensionalcomplex
reverberatingmedium, J.Acoust.Soc.Am. ,135,1034–1043,doi:10.1121/1.4864485.
Cooley,J.W.,andJ.W.Tukey(1965),AnalgorithmforthemachinecalculationofcomplexFourierseries, Math.Comput. ,19,297–301.
Corciulo,M.,P.Roux,M.Campillo,andD.Dubucq(2012),Instantaneousphasevariationforseismicvelocitymonitoringfromambientnoise
attheexplorationscale, Geophysics ,77,Q37–Q44,doi:10.1190/geo2011–0363.1.
Cox,H.(1973),Spatialcorrelationinarbitrarynoisefieldswithapplicationtoambientseanoise, J.Acoust.Soc.Am. ,54,1289–1301.
Curtis, A., P. Gerstoft, H. Sato, R. Snieder, and K. Wapenaar (2006), Interferometry—Turning noise into signal, The Leading Edge ,25,
1082–1092.
DeRidder,S.,andJ.Dellinger(2011),Ambientseismicnoiseeikonaltomographyfornear-surfaceimagingatValhall, TheLeadingEdge ,30,
506–512,doi:10.1190/1.3589108.
DeRidder,S.A.L.,B.L.Biondi,andR.G.Clapp(2014),Time-lapseseismicnoisecorrelationtomographyatValhall, Geophys.Res.Lett. ,41,
6116–6122,doi:10.1002/2014GL061156.1.
Derode,A.,E.Larose,M.Campillo,andM.Fink(2003),HowtoestimatetheGreen’sfunctionofaheterogeneousmediumbetweentwo
passivesensors?Applicationtoacousticwaves, Appl.Phys.Lett. ,83,3054–3056,doi:10.1063/1.1617373.
Duvall,T.L.,S.M.Jefferies,J.W.Harvey,andM.A.Pomerantz(1993),Time-distancehelioseismology, Nature,362,430–432.
Eckart,C.(1953),Thetheoryofnoiseincontinuousmedia, J.Acoust.Soc.Am. ,25,195–199.
Ekström,G.(2014),LoveandRayleighphase-velocitymaps,5–40s,ofthewesternandcentralUSAfromUSArraydata, EarthPlanet.Sci.
Lett.,402,42–49.
Ekström,G.,G.A.Abers,andS.C.Webb(2009),Determinationofsurface-wavephasevelocitiesacrossUSArrayfromnoiseandAki’sspectral
formulation, Geophys.Res.Lett. ,36,L18301,doi:10.1029/2009GL039131.
Fink,M.(1999),Time-reversedacoustics, Sci.Am.,281,91–97.
Fink,M.(2006),Time-reversalacousticsincomplexenvironments, Geophysics ,71,SI151–SI164,doi:10.1190/1.2215356.
Fleury, C., R. Snieder, and K. Larner (2010), General representation theorem for perturbed media and application to green’s function
retrievalforscatteringproblems, Geophys.J.Int. ,183,1648–1662.
Gorbatov,A.,E.Saygin,andB.L.N.Kennett(2013),Crustalpropertiesfromseismicstationautocorrelograms, Geophys.J.Int. ,192,861–870,
doi:10.1093/gji/ggs064.
Gouedard,P., etal.(2008),Cross-correlationofrandomfields:Mathematicalapproachandapplications, Geophys.Prospect. ,56,375–393.
Groenenboom,J.,andR.Snieder(1995),Attenuation,dispersion,andanisotropybymultiplescatteringoftransmittedwavesthrough
distributionofscattering, J.Acoust.Soc.Am. ,98,3482–3492.
Gualtieri,L.,E.Stutzmann,Y.Capdeville,F.Ardhuin,M.Schimmel,A.Mangeney,andA.Morelli(2013),Modellingsecondarymicroseismic
noisebynormalmodesummation, Geophys.J.Int. ,193,1732–1745,doi:10.1093/gji/ggt090.
Halliday, D., and A. Curtis (2008), Seismic interferometry, surface waves and source distribution, Geophys. J. Int. ,175, 1067–1087,
doi:10.1111/j.1365-246X.2008.03918.x.
Haney, M. M. (2009), Infrasonic ambient noise interferometry from correlations of microbaroms, Geophys. Res. Lett. ,36, L19808,
doi:10.1029/2009GL040179.
Harmon, N., P. Gerstoft, C. A. Rychert, G. A. Abers, M. S. de la Cruz, andK. M. Fischer (2008), Phasevelocities from seismic noise using
beamformingandcrosscorrelationinCostaRicaandNicaragua, Geophys.Res.Lett. ,35,L19303,doi:10.1029/2008GL035387.
Hasselmann,K.(1963),Astatisticalanalysisofthegenerationofmicroseisms, Rev.Geophys. ,1,177–210,doi:10.1029/RG001i002p00177.
Hildebrand,F.B.(1976), AdvancedCalculusforApplications ,2nded.,Prentice-Hall,EnglewoodCliffs,N.J.
Hillers,G.,N.Graham,M.Campillo,S.Kedar,M.Landes,andN.Shapiro(2012),Globaloceanicmicroseismsourcesasseenbyseismicarrays
andpredictedbywaveactionmodels, Geochem.Geophys.Geosyst. ,13,Q01021,doi:10.1029/2011GC003875.
Ito,Y.,andK.Shiomi(2012),Seismicscattererswithinsubductingslabrevealedfromambientnoiseautocorrelation, Geophys.Res.Lett. ,39,
L19303,doi:10.1029/2012GL053321.
Kedar,S.,M.Longuet-Higgins,F.Webb,N.Graham,R.Clayton,andC.Jones(2008),TheoriginofdeepoceanmicroseismsintheNorth
AtlanticOcean, Proc.R.Soc.London,Ser.A ,464,777–793,doi:10.1098/rspa.2007.0277.
Kinsler,L.E.,A.R.Frey,A.B.Coppens,andJ.V.Sanders(1999), FundamentalsofAcoustics ,4thed.,Wiley,Hoboken,N.J.
Kohler,M. D.,T.H. Heaton, and S. C.Bradford (2007),Propagatingwaves in the steel, moment-framefactor buildingrecorded during
earthquakes, Bull.Seismol.Soc.Am. ,97,1334–1345,doi:10.1785/0120060148.
Larose,E.,L.Margerin,A.Derode,B.vanTiggelen,M.Campillo,N.M.Shapiro,A.Paul,L.Stehly,andM.Tanter(2006),Correlationofrandom
wavefields:Aninterdisciplinaryreview, Geophysics ,71,SI11–SI21.
Lobkis,O.I.,andR.L.Weaver(2001),OntheemergenceoftheGreen’sfunctioninthecorrelationsofadiffusefield, J.Acoust.Soc.Am. ,110,
3011–3016,doi:10.1121/1.1417528.
Longuet-Higgins,M.S.(1950),Atheoryoftheoriginofmicroseisms, Philos.Trans.R.Soc.London,Ser.A ,243,1–35.
Mainsant,G.,E.Larose,C.Brönmimann,D.Jongmans,C.Michoud,andM.Jaboyedoff(2012),Ambientseismicnoisemonitoringofaclay
landslide:Towardfailureprediction, J.Geophys.Res. ,117,F01030,doi:10.1029/2011JF002159.
Malcolm,A.E.,J.A.Scales,andB.A.vanTiggelen(2004),Extractingthegreenfunctionfromdiffuse,equipartitionedwaves, Phys.Rev.E ,70,
015601,doi:10.1103/PhysRevE.70.015601.
Marston, P. L. (2001), Generalized optical theorem for scatterers havinginversion symmetry: Applications to acoustic backscattering,
J.Acoust.Soc.Am. ,109,1291–1295.
Mikesell,T.D.,K.vanWijk,T.E.Blum,R.Snieder,andH.Sato(2012),Analyzingthecodafromcorrelatingscatteredsurfacewaves, J.Acoust.
Soc.Am.,131,EL275–EL281,doi:10.1121/1.3687427.
Morse,P.K.,andK.U.Ingard(1986), TheoreticalAcoustics ,PrincetonUniv.Press,Princeton,N.J.
Mulargia,F.(2012),Theseismicnoisewavefieldisnotdiffuse, J.Acoust.Soc.Am. ,131,2853–2858,doi:10.1121/1.3689551.
Nakahara, H. (2006), A systematic study of theoretical relations between spatial correlation and Green’s function in one-, two- and
three-dimensionalrandomscalarwavefields, Geophys.J.Int. ,167,1097–1105,doi:10.1111/j.1365-246X.2006.03170.x.
Newton,R.G.(2002), ScatteringTheoryofWavesandParticles ,Dover,Mineola,NewYork.
Nishida, K. (2013), Global propagation of body waves revealed by cross-correlation analysis of seismic hum, Geophys. Res. Lett. ,40,
1691–1696,doi:10.1002/grl.50269.
Poli,P.,H.A.Pedersen,M.Campillo,andPOLENET/LAPNETWorkingGroup(2012a),Emergenceofbodywavesfromcross-correlationof
shortperiodseismicnoise, Geophys.J.Int. ,188,549–558,doi:10.1111/j.1365-246X.2011.05271.x.
BOSCHIANDWEEMSTRA AMBIENT-NOISECROSSCORRELATION 450

ReviewsofGeophysics 10.1002/2014RG000455
Poli,P.,H.A.Pedersen,M.Campillo,andPOLENET/LAPNETWorkingGroup(2012b),Body-waveimagingofEarth’smantlediscontinuities
fromambientseismicnoise, Science,338,1063–1065,doi:10.1126/science.1228194.
Ravasi,M.,andA.Curtis(2013),Nonlinearscatteringbasedimaginginelasticmedia:Theory,theorems,andimagingconditions, Geophysics ,
78,S137–S155,doi:10.1190/GEO2012-0286.1.
Ritzwoller,M.(2014),Ambientnoiseseismicimaging,in McGrawHillYearBookonScienceandTechnology ,pp.1–13,McGraw-Hill,NewYork.
Roux,P.,K.Sabra,W.Kuperman,andA.Roux(2005),Ambientnoisecrosscorrelationinfreespace:Theoreticalapproach, J.Acoust.Soc.Am. ,
117,79–84,doi:10.1121/1.1830673.
Ruigrok,E.,X.Campman,andK.Wapenaar(2011),ExtractionofP-wavereflectionsfrommicroseisms, C.R.Geosci. ,343,512–525.
Sanchez-Sesma,F.J.,andM.Campillo(2006),RetrievaloftheGreen’sfunctionfromcrosscorrelation:Thecanonicalelasticproblem, Bull.
Seismol.Soc.Am. ,96,1182–1191,doi:10.1785/0120050181.
Sato,H.(2010),RetrievalofGreen’sfunctionhavingcodawavesfromthecross-correlationfunctioninascatteringmediumilluminatedbya
randomlyhomogeneousdistributionofnoisesourcesonthebasisofthefirst-orderBornapproximation, Geophys.J.Int. ,180,759–764,
doi:10.1111/j.1365–246X.2009.04432.x.
Shapiro,N.M.,M.Campillo,L.Stehly,andM.H.Ritzwoller(2005),High-resolutionsurface-wavetomographyfromambientseismicnoise,
Science,307,1615–1618.
Smith,S.W.(2011), TheScientistandEngineer’sGuidetoDigitalSignalProcessing ,Calif.Tech.Publ.[Availableathttp://www.dspguide.com/.]
Snieder,R.(1986),3Dlinearizedscatteringofsurfacewavesandaformalismforsurfacewaveholography, Geophys.J.R.Astron.Soc. ,84,
581–605.
Snieder,R.(2004),ExtractingtheGreen’sfunctionfromthecorrelationofcodawaves:Aderivationbasedonstationaryphase, Phys.Rev.E ,
69,046610,doi:10.1103/PhysRevE.69.046610.
Snieder,R.(2007),ExtractingtheGreen’sfunctionofattenuatingheterogeneousacousticmediafromuncorrelatedwaves, J.Acoust.Soc.
Am.,121,2637–2643,doi:10.1121/1.2713673.
Snieder, R., and E. Larose (2013), Extracting Earth’s elastic wave response from noise measurements, Annu. Rev. Earth Planet. Sci. ,41,
183–206,doi:10.1146/annurev-earth-050212-123936.
Snieder,R.,andE. ¸ Safak(2006),Extractingthebuildingresponseusingseismicinterferometry:TheoryandapplicationtotheMillikan
LibraryinPasadena,California, Bull.Seismol.Soc.Am. ,96,586–598,doi:10.1785/0120050109.
Snieder,R.,K.vanWijk,M.Haney,andR.Calvert(2008),CancellationofspuriousarrivalsinGreen’sfunctionextractionandthegeneralized
opticaltheorem, Phys.Rev.E ,78,036606,doi:10.1103/PhysRevE.78.036606.
Snieder,R.,Y.Fan,E.Slob,andK.Wapenaar(2010),EquipartitioningisnotsufficientforGreen’sfunctionextraction, EarthquakeSci. ,23,
403–415,doi:10.1007/s11589-010-0739-1.
Stehly,L.(2007),Tomographyusingcorrelationsofseismicambientnoise,thesis,Univ.JosephFourier,Grenoble,France.
Stehly,L.,M.Campillo,andN.M.Shapiro(2006),Astudyoftheseismicnoisefromitslong-rangecorrelationproperties, J.Geophys.Res. ,111,
B10306,doi:10.1029/2005JB004237.
Strauss,W.A.(2008), PartialDifferentialEquations ,JohnWiley,Hoboken,N.J.
Tanimoto,T.(1990),Modellingcurvedsurfacewavepaths:Membranesurfacewavesynthetics, Geophys.J.Int. ,102,89–100.
Traer, J., and P. Gerstoft (2014), A unified theory of microseisms and hum, J. Geophys. Res. Solid Earth ,119, 3317–3339,
doi:10.1002/2013JB010504.
Tromp, J., and F. Dahlen (1993), Variational principles for surface wave propagation on a laterally heterogeneous Earth–III. Potential
representation, Geophys.J.Int. ,112,195–209.
Tsai,V.C(2009),Onestablishingtheaccuracyofnoisetomographytravel-timemeasurementsinarealisticmedium, Geophys.J.Int. ,178,
1555–1564,doi:10.1111/j.1365-246X.2009.04239.x.
Tuckerman,M.E.(1987), StatisticalMechanics:TheoryandMolecularSimulation ,OxfordUniv.Press,Oxford,U.K.
Udías,A.(1999), PrinciplesofSeismology ,CambridgeUniv.Press,Cambridge,U.K.
vanManen,D.-J.,J.O.A.Robertsson,andA.Curtis(2005),Modelingofwavepropagationininhomogeneousmedia, Phys.Rev.Lett. ,94,
164301,doi:10.1103/PhysRevLett.94.164301.
Wapenaar,K.(2004),RetrievingtheelastodynamicGreen’sfunctionofanarbitraryinhomogeneousmediumbycrosscorrelation, Phys.Rev.
Lett.,93,254301,doi:10.1103/PhysRevLett.93.254301.
Wapenaar, K., and J. Fokkema (2006), Green’s function representations for seismic interferometry, Geophysics ,71, SI33–SI46,
doi:10.1190/1.2213955.
Wapenaar,K.,J.Fokkema,andR.Snieder(2005),RetrievingtheGreen’sfunctioninanopensystembycrosscorrelation:Acomparisonof
approaches, J.Acoust.Soc.Am. ,118,2783–2786,doi:10.1121/1.2046847.
Wapenaar,K.,E.Slob,andR.Snieder(2006),UnifiedGreen’sfunctionretrievalbyseismicinterferometry, Phys.Rev.Lett. ,97,234301.
Wapenaar,K.,E.Slob,R.Snieder,andA.Curtis(2010a),Tutorialonseismicinterferometry:Part2—Underlyingtheoryandnewadvances,
Geophysics ,75,A211–A227,doi:10.1190/1.3463440.
Wapenaar,K.,D.Draganov,R.Snieder,X.Campman,andA.Verdel(2010b),Tutorialonseismicinterferometry:Part1—Basicprinciplesand
applications, Geophysics ,75,A195–A209,doi:10.1190/1.3457445.
Weaver,R.,B.Froment,andM.Campillo(2009),Onthecorrelationofnon-isotropicallydistributedballisticscalardiffusewaves, J.Acoust.
Soc.Am.,126,1817–1826,doi:10.1121/1.3203359.
Weaver,R.L.,andO.I.Lobkis(2002),OntheemergenceoftheGreen’sfunctioninthecorrelationsofadiffusefield:Pulse-echousing
thermalphonons, Ultrasonics ,40,435–439.
Weaver,R.L.,andO.I.Lobkis(2004),DiffusefieldsinopensystemsandtheemergenceoftheGreen’sfunction, J.Acoust.Soc.Am. ,116,
2731–2734.
Weemstra,C.,L.Boschi,A.Goertz,andB.Artman(2013),Seismicattenuationfromrecordingsofambientnoise, Geophysics ,78,Q1–Q14,
doi:10.1190/geo2012-0132.1.
Weemstra,C.,W.Westra,R.Snieder,andL.Boschi(2014),Onestimatingattenuationfromtheamplitudeofthespectrallywhitenedambient
seismicfield, Geophys.J.Int. ,197,1770–1788,doi:10.1093/gji/ggu088.
Wegler,U.,andC.Sens-Schonfelder(2007),Faultzonemonitoringwithpassiveimageinterferometry, Geophys.J.Int. ,168,1029–1033.
Weisstein, E. W. (1999–2013), Delta Function, From MathWorld—A Wolfram Web Resource. [Available at http://mathworld.
wolfram.com/DeltaFunction.html.]
Wong,R.(1986), AsymptoticApproximationsofIntegrals ,Soc.forInd.andAppl.Math.,Philadephia,Pa.
Woodard,M.F.(1997),Implicationsoflocalized,acousticabsorptionforheliotomographicanalysisofsunspots, Astrophys.J. ,485,890–894.
BOSCHIANDWEEMSTRA AMBIENT-NOISECROSSCORRELATION 451

ReviewsofGeophysics 10.1002/2014RG000455
Erratum
Intheoriginallypublishedversionofthisarticle,EquationsE7andE34areincorrectlyexpressed.Theequations
havesincebeencorrectedandthisversionmaybeconsideredtheauthoritativeversionofrecord.
BOSCHIANDWEEMSTRA AMBIENT-NOISECROSSCORRELATION 452