Introduction

Imaging of passive seismic sources provides crucial geomechanical information about the subsurface. In oil and gas exploration, passive sources are useful for monitoring the propagation of fractures when high-pressure fluids are injected into formations, such as hydraulic fracturing (Maxwell, 2014; Montgomery and Smith, 2010). Using the wave equation, Gajewski and Tessmer (2005) and Fink (2006) image passive seismic sources by reverse-time migration (McMechan, 1982), which is also known as time-reversal imaging (TRI). Artman et al. (2010) and Witten and Artman (2011) generalize TRI using P- and S-waves. Artman et al. (2010) also propose to use an auto-correlation time-reversal imaging (ATRI) condition to produce an image of source locations. Zhu (2014) improves the quality of TRI by compensating for attenuation. In order to improve the resolution of TRI, Sun et al. (2015) and Nakata and Beroza (2016) propose a cross-correlation time-reversal imaging (CTRI) condition by performing a zero-lag cross-correlation among all the indenpendently back-propagated receiver wavefields. CTRI is capable of achieving high-resolution images of passive sources. There are two major issues associated with CTRI as well as ATRI. First, the correlation destroys phase information, and leads to exponential growth of relative amplitude of seismic sources. Second, since multiplication is not a linear operation, it does not have an adjoint operator to go back to the data domain so that full waveform inversion (FWI) can be applied for estimating Earth properties.

FWI is a technique commonly used in reflection seismology for obtaining high-resolution Earth properties (Virieux and Operto, 2009; Warner et al., 2013; Plessix, 2009). FWI has also been applied to passive seismic data at the global scale (Zhu et al., 2015; Kim et al., 2011; Tromp et al., 2005) or for exploration purposes (Kamei and Lumley, 2014; Behura, 2015). In both settings, an unknown source function makes the inversion task more difficult. In reflection seismology the source location and excitation time is usually known while one needs to estimate the source signature. Aravkin et al. (2012) and Rickett (2013) have shown that the least-squares inversion for source functions can be incorporated into the nonlinear inversion framework using the method of variable projection (Golub and Pereyra, 1973). For passive seismic data, however, the inversion for source functions is more challenging since the locations and start times of the seismic sources are generally unknown. In fact, the source function spans a four-dimensional space including the three spatial axes and one time axis. Without proper preconditioning, the computational cost of least-squares inversion for the passive source functions alone can be intractable, not to mention joint inversion for both the source function and Earth properties.

In this paper, we formulate a joint inversion framework for both passive seismic sources and seismic velocities. We first review existing passive imaging conditions using wave equation. We then propose a local normalization operator which effectively normalizes all the local maxima that correspond to potential seismic sources. Next, we use the normalized volume as a masking operator to project the time-reversed wavefield into the space of admissible models. This amounts to a simple diagonal weighting of seismic sources, which can be used in iterative inversion of the true seismic sources by the preconditioned conjugate gradient (PCG) method. Since the four-dimensional model space is constrained to a sparse subset using the masking operator, the convergence rate can be vastly improved. Finally, under the framework of variable projection (Golub and Pereyra, 1973; Rickett, 2013), seismic sources and velocity can be jointly inverted using separation of variables. We use several synthetic examples to demonstrate the improved efficiency using the proposed approach.