This thesis proposes a computationally efficient method for approximating the inverse of the normal operator arising in the linearized inverse problem for reflection seismology. The inversion of the normal operator using direct matrix methods is computationally infeasible. Approximate inverses estimate the solution of the inverse problem or precondition iterative methods. Application of the normal operator requires an expensive solution of large scale PDE problems. However, the normal operator approximately commutes with pseudodifferential operators, hence shares their near diagonality in a frame of localized monochromatic pulses. Estimation of a diagonal representation in this frame encoded in the symbol of the normal operator: (1) follows from its application to a single input vector; (2) suffices to approximate its inverse. I use an efficient algorithm to apply pseudodifferential operators, given their symbol, to construct a rapidly converging optimization algorithm that estimates the symbol of an inverse for the normal operator, thereby approximately solving the inverse problem.