Computational Problem in Three-Dimensional, Plane-wave Imaging with P to S Converted Waves Recorded by Broadband Arrays

We have developed a set of computer programs that implement the three-dimensional imaging methods described by Poppeliers and Pavlis [2003]. The original work developed the theoretical framework for this technique, but did not implement the fully three-dimensional version due to the requirement that efficient computational procedures are necessary to make the algorithm feasible. We implemented the concepts in two programs linked by a relational database. The first program implements a plane wave decomposition of the wavefield using the pseudostation stacking method to interpolate the wavefield onto a regular spatial grid. We then use a relational database to sort data into plane wave components linked to this interpolated, regular grid. Imaging is completed using sums of plane wave components in an integral equation derived from the inverse generalized Radon transform. This allows the image to accumulate incrementally as a weighted sum of plane wave components with the weights varying in space reducing the memory requirements to practical levels on modern, distributed-memory, parallel systems. This imaging algorithm is readily adaptable on massively parallel machines by distributing back propagation computations for plane wave components among processors and devoting one node to summing the components. The algorithm is built around a new object-oriented, C++ library for manipulating seismograms and the two and three-dimensional data objects that are intermediaries in the algorithm.