{
  "cells": [
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "%matplotlib notebook"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "\n# Calculating a Digital Linear Filter\n\nThis is an example for the add-on ``fdesign``. The example is taken from the\narticle Werthm\u00fcller et al., 2019. Have a look at the article repository on\n[empymod/article-fdesign](https://github.com/emsig/article-fdesign) for\nmany more examples.\n\n**Reference**\n\n- Werthm\u00fcller, D., K. Key, and E. Slob, 2019, **A tool for designing digital\n  filters for the Hankel and Fourier transforms in potential, diffusive, and\n  wavefield modeling**:  Geophysics, 84(2), F47-F56; DOI:\n  [10.1190/geo2018-0069.1](http://doi.org/10.1190/geo2018-0069.1).\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "import empymod\nimport numpy as np\nimport matplotlib.pyplot as plt\nplt.style.use('ggplot')"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "inp = {'r': np.logspace(0, 10, 1000),\n       'r_def': (1, 1, 2),\n       'n': 201,\n       'name': 'test',\n       'full_output': True,\n       'fI': (empymod.fdesign.j0_1(5), empymod.fdesign.j1_1(5))}"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## 1. Rough overview over a wide range\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "filt1, out1 = empymod.fdesign.design(\n        spacing=(0.01, 0.2, 10), shift=(-4, 0, 10), save=False, **inp)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## 2. First focus\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "filt2, out2 = empymod.fdesign.design(\n        spacing=(0.04, 0.1, 10), shift=(-3, -0.5, 10), save=False, **inp)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## 3. Final focus\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "filt, out = empymod.fdesign.design(\n        spacing=(0.047, 0.08, 10), shift=(-2.4, -0.75, 10), finish=False,\n        save=False, **inp)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "### To reproduce exactly the same filter as wer_201_2018:\n\n::\n\n    filt_orig, out_orig = fdesign.load_filter('wer201', True)\n    fdesign.plot_result(filt_orig, out_orig)\n    filt_orig, out_orig = fdesign.design(\n            spacing=out_orig[0][0], shift=out_orig[0][1], **inp)\n\n\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## Plot the result\n\n### Plot function\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "def plotresult(depth, res, zsrc, zrec):\n    x = np.arange(1, 101)*200\n    inp = {\n        'src': [0, 0, depth[1]-zsrc],\n        'rec': [x, x*0, depth[1]-zrec],\n        'depth': depth,\n        'res': res,\n        'ab': 11,\n        'freqtime': 1,\n        'verb': 1,\n    }\n\n    kong241 = empymod.dipole(htarg={'dlf': 'kong_241_2007'}, **inp)\n    key201 = empymod.dipole(htarg={'dlf': 'key_201_2012'}, **inp)\n    and801 = empymod.dipole(htarg={'dlf': 'anderson_801_1982'}, **inp)\n    test = empymod.dipole(htarg={'dlf': filt}, **inp)\n    wer201 = empymod.dipole(htarg={'dlf': 'wer_201_2018'}, **inp)\n    qwe = empymod.dipole(ht='qwe', **inp)\n\n    plt.figure(figsize=(8, 3.5))\n    plt.subplot(121)\n    plt.semilogy(x, qwe.amp(), c='0.5', label='QWE')\n    plt.semilogy(x, kong241.amp(), 'k--', label='Kong241')\n    plt.semilogy(x, key201.amp(), 'k:', label='Key201')\n    plt.semilogy(x, and801.amp(), 'k-.', label='And801')\n    plt.semilogy(x, test.amp(), 'r', label='This filter')\n    plt.semilogy(x, wer201.amp(), 'b', label='Wer201')\n    plt.legend()\n    plt.xticks([0, 5e3, 10e3, 15e3, 20e3])\n    plt.xlim([0, 20e3])\n\n    plt.subplot(122)\n    plt.semilogy(x, np.abs((kong241-qwe)/qwe), 'k--', label='Kong241')\n    plt.semilogy(x, np.abs((key201-qwe)/qwe), 'k:', label='Key201')\n    plt.semilogy(x, np.abs((and801-qwe)/qwe), 'k-.', label='And801')\n    plt.semilogy(x, np.abs((test-qwe)/qwe), 'r', label='This filter')\n    plt.semilogy(x, np.abs((wer201-qwe)/qwe), 'b', label='Wer201')\n    plt.legend()\n    plt.xticks([0, 5e3, 10e3, 15e3, 20e3])\n    plt.xlim([0, 20e3])\n    plt.ylim([1e-14, 1])\n\n    plt.show()"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "### Plot the individual models\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "plotresult([-1e50, 2000], [2e14, 1/3.2, 1], 50, 0)  # KONG model\nplotresult([0, 1000, 2000, 2100], [1/1e-12, 1/3.3, 1, 100, 1], 10, 0)  # KEY m.\nplotresult([0, 1, 1000, 1100], [2e14, 10, 10, 500, 10], 0.5, 0.2)  # LAND model"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "empymod.Report()"
      ]
    }
  ],
  "metadata": {
    "kernelspec": {
      "display_name": "Python 3",
      "language": "python",
      "name": "python3"
    },
    "language_info": {
      "codemirror_mode": {
        "name": "ipython",
        "version": 3
      },
      "file_extension": ".py",
      "mimetype": "text/x-python",
      "name": "python",
      "nbconvert_exporter": "python",
      "pygments_lexer": "ipython3",
      "version": "3.8.6"
    }
  },
  "nbformat": 4,
  "nbformat_minor": 0
}