SyntheticSeismogram.ipynb

{
 "metadata": {
  "name": "",
  "signature": "sha256:61b755300c238685b7228857dace815c244576bc2bf23349ef94f7ff39dcda70"
 },
 "nbformat": 3,
 "nbformat_minor": 0,
 "worksheets": [
  {
   "cells": [
    {
     "cell_type": "heading",
     "level": 1,
     "metadata": {},
     "source": [
      "Using the Notebook"
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "This is the <a href=\"http://ipython.org/notebook.html\">iPython Notebook</a>, an interactive coding and computation environment. For this lab, you do not have to write any code, you will only be running it. \n",
      "\n",
      "To use the notebook:\n",
      "- \"Shift + Enter\" runs the code within the cell (so does the forward arrow button near the top of the document)\n",
      "- You can alter variables and re-run cells\n",
      "- If you want to start with a clean slate, restart the Kernel either by going to the top, clicking on Kernel: Restart, or by \"esc + 00\" (if you do this, you will need to re-run Step 0 before running any other cells in the notebook) \n",
      "\n",
      "Instructions as to how to set up Python and the iPython notebook on your personal computer are attached in the appendix of the lab"
     ]
    },
    {
     "cell_type": "heading",
     "level": 2,
     "metadata": {},
     "source": [
      "Step 0: Import the Necessary Packages"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "# Import the necessary packages\n",
      "%pylab inline                            \n",
      "from syntheticSeismogramImport import *  "
     ],
     "language": "python",
     "metadata": {},
     "outputs": []
    },
    {
     "cell_type": "heading",
     "level": 1,
     "metadata": {},
     "source": [
      "Constructing a Synthetic Seismogram"
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "We will consider a model consisting of two horizontal layers over a half-space. Each layer has density $\\rho_i$ and seismic velocity $v_i$. In this notebook, we will start from the physical property model and walk through how to construct a synthetic, normal incidence seismogram. \n",
      "\n",
      "<img src=\"./Figures/geoModel.png\" style=\"width: 50%; height: 50%\"></img>"
     ]
    },
    {
     "cell_type": "heading",
     "level": 2,
     "metadata": {},
     "source": [
      "Step 1: From Physical Properties to a Reflectivity Series"
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "We begin by constructing a reflectivity series from the physical property model. The physical properties characterizing the model are:\n",
      "- density ($\\rho$) and \n",
      "- seismic velocity ($v$)\n",
      "\n",
      "The acoustic impedance of each layer is:\n",
      "$$\n",
      "Z_i = \\rho_i v_i\n",
      "$$\n",
      "\n",
      "From acoustic impedance, the down-going reflection coefficient for each interface is given by\n",
      "$$\n",
      "r_{i,i+1} = \\frac{Z_{i+1}-Z_i}{Z_{i+1}+Z_i}\n",
      "$$\n",
      "\n",
      "and the transmission coefficient is\n",
      "$$\n",
      "t_{i,i+1} = \\frac{2Z_{i}}{Z_{i+1}+Z_i}\n",
      "$$\n",
      "\n",
      "<img src=\"./Figures/rt1interface.png\" style=\"width: 45%; height: 45%\"></img>\n",
      "\n",
      "The true reflectivity accounts for both reflection and transmission. For the reflection off of the first layer, these are equivalent. For the reflection off of the second interface, the reflection coefficient $r_{2,3}$ in multiplied by $t_{1,2}t_{2,1}$ to get the true reflectivity. In the below plot, this effect of the transmission coefficients can be included or not using the toggle \"usingT.\"\n",
      "\n",
      "<img src=\"./Figures/rt2interface.png\" style=\"width: 50%; height: 50%\"></img>\n",
      "\n"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "# Create Interactive Plot for Logs\n",
      "logs = interact(plotLogsInteract,d2=(0.,100.,5),d3=(100.,200.,5),rho1=(2000.,5000.,50.),rho2=(2000.,5000.,50.),rho3=(2000.,5000.,50.),v1=(300.,4000.,50.),v2=(300.,4000.,50.),v3=(300.,4000.,50.))"
     ],
     "language": "python",
     "metadata": {},
     "outputs": []
    },
    {
     "cell_type": "heading",
     "level": 2,
     "metadata": {},
     "source": [
      "Step 2: Converting from Depth to Time"
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Now we have the reflectivity series as a function of depth. With seismic, we measure a signal as a function of time. So we must apply a conversion from depth to time. We do this by computing the time it takes for a signal to reach a given depth and return to the surface. "
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "# Create depth-time interactive plot\n",
      "interact(plotTimeDepthInteract,d2=(0.,100.,5),d3=(100.,200.,5),v1=(300.,4000.,50.),v2=(300.,4000.,50.),v3=(300.,4000.,50.))"
     ],
     "language": "python",
     "metadata": {},
     "outputs": []
    },
    {
     "cell_type": "heading",
     "level": 2,
     "metadata": {},
     "source": [
      "Step 3: Wavelet * Reflectivity = Seismogram"
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Now that we have the reflectivity series in time, the next step is to choose the input pulse and construct our seismogram. For the following examples, we will use a <a href=\"http://subsurfwiki.org/wiki/Ricker_wavelet\">Ricker Wavelet</a> with peak frequency $f$. \n",
      "\n",
      "A seismogram is the convolution of the wavelet and reflectivity series. Here, you can adjust the peak frequency (wavF) of the wavelet and its amplitude (wavA). \n",
      "\n",
      "The geologic model used is:\n",
      "<img src=\"./Figures/Step3FixedMod.png\" style=\"width: 50%; height: 50%\"></img>"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "# Interactive seismogram plot for a fixed geologic model\n",
      "interact(plotSeismogramInteractFixMod,wavf=(5.,100.,5.),wavA=(-2.,2.,0.25))"
     ],
     "language": "python",
     "metadata": {},
     "outputs": []
    },
    {
     "cell_type": "heading",
     "level": 2,
     "metadata": {},
     "source": [
      "Step 4: Seismic Resolution"
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Here, we have the following geologic model: \n",
      "\n",
      "<img src=\"./Figures/resolnModel.png\" style=\"width: 50%; height: 50%\"></img>\n",
      "\n",
      "When referring to vertical resolution, the question to ask is: \"Can the two arrivals (one from the top, and one from the bottom of the layer) be distinguished?\" \n",
      "\n",
      "Adjust the layer thickness (h2) for the middle layer and the frequency of the input pulse to investigate vertical resolution. You can also add noise to the trace. "
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "interact(plotSeismogramInteractRes,h2=(0.,200,1.),wavf=(5.,100.,2.5),AddNoise=False)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": []
    },
    {
     "cell_type": "heading",
     "level": 1,
     "metadata": {},
     "source": [
      "Step 5: Forward Modelling"
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "You can adjust all of the parameters. Have fun!"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "interact(plotSeismogramInteract,d2=(0.,200.,5),d3=(50.,200.,5),rho1=(2000.,5000.,50.),rho2=(2000.,5000.,50.),rho3=(2000.,5000.,50.),v1=(300.,4000.,50.),v2=(300.,4000.,50.),v3=(300.,4000.,50.),wavf=(5.,100.,2.5),wavA=(-1.,2.,0.5),addNoise=False,usingT=True) "
     ],
     "language": "python",
     "metadata": {},
     "outputs": []
    }
   ],
   "metadata": {}
  }
 ]
}