{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Course website**: http://lagex.github.io/geofisica2\n",
    "\n",
    "**Note**: This notebook is part of the course \"Geofísica 2\" of Geology program of the \n",
    "[Universidade do Estado do Rio de Janeiro](http://www.uerj.br/). \n",
    "All content can be freely used and adapted under the terms of the \n",
    "[Creative Commons Attribution 4.0 International License](http://creativecommons.org/licenses/by/4.0/).\n",
    "\n",
    "![Creative Commons License](https://i.creativecommons.org/l/by/4.0/88x31.png)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Esse documento que você está usando é um [IPython notebook](http://ipython.org/notebook.html). É um documento interativo que mistura texto (como esse), código (como abaixo), e o resultado de executar o código (que pode ser números, texto, figuras, videos, etc)."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Prática 1 - Ondas sísmicas\n",
    "\n",
    "Vamos explorar como as ondas sísmicas P e S se comportam, tanto em meios homogêneos como quando há uma interface entre dois meios. Para isso, vamos utilizar as simulações de ondas da biblioteca [Fatiando a Terra](http://www.fatiando.org). Essas simulações utilizam o [método de diferenças finitas](http://en.wikipedia.org/wiki/Finite_difference_method) para calcular soluções da equação da onda.\n",
    "\n",
    "As tarefas abaixo acompanham as questões que vocês devem responder. Usem as simulações abaixo para te ajudar a responder as perguntas."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Preparação\n",
    "\n",
    "Rode a célula abaixo para carregar as componentes necessárias para fazer as simulações. Não se preocupe se aparecer um `:0: FutureWarning: IPython widgets are experimental and may change in the future.` abaixo. Isso é consequência de utilizar tecnologia de ponta."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "%matplotlib inline\n",
    "from __future__ import division\n",
    "import numpy as np\n",
    "from fatiando import utils\n",
    "from fatiando.seismic.wavefd import Ricker, ElasticPSV, ElasticSH\n",
    "from fatiando.vis import mpl\n",
    "import fatiando"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "print(fatiando.__version__)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Frentes de onda\n",
    "\n",
    "Vamos simular uma onda P (esqueçam as ondas S por enquanto). Primeiro vamos ver como ela se comporta em um meio homogêneo, sem camadas.\n",
    "\n",
    "Rode a célula abaixo para criar um modelo da subsuperfície homogêneo."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "shape = (300, 400)\n",
    "spacing = 10\n",
    "extent = [0, shape[1]*spacing, shape[0]*spacing, 0]\n",
    "velocidade = np.ones(shape)*3000\n",
    "densidade = np.ones(shape)*2200"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Agora vamos criar o nosso simulador de ondas."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "sim1 = ElasticSH(velocidade, densidade, spacing=spacing, taper=0.005)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "e colocar uma fonte explosiva na superfície."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "sim1.add_point_source((0, shape[1]//2), Ricker(5, 20, 1/20))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Agora que temos nossa simulação pronta, rode a célcula abaixo para avançar a simulação 400 passos no tempo."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "sim1.run(400)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Quando a barra acima encher, rode a célula abaixo para gerar uma animação da propagação da onda. Isso pode demorar um pouco. A parte de cima da figura é a superfície."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "sim1.animate(every=10, embed=True, dpi=50, cutoff=0.001)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Rode a próxima célula para explorar fotos de cada etapa da simulação, uma de cada vez."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "sim1.explore(every=20)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Nota\n",
    "\n",
    "Se quiser continuar a simulação e avançá-la mais no tempo, rode novamente a célula `sim1.run(400)` para avançar mais 400 iterações. Depois, rode novamente a célula que cria a animação e ela será atualizada."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Onda P incidindo em uma interface (simplificado)\n",
    "\n",
    "Agora que vocês visualizaram a frente de onda, vamos ver o que acontece quando a subsuperfície contem duas camadas. A simulação abaixo considera **somente ondas P e ignora ondas S**. \n",
    "\n",
    "Vamos fazer um modelo com a camada de vaixo mais densa e com velocidade maior. A interface está localizada a 1000m de profundidade."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "densidade = np.ones(shape)*2200\n",
    "velocidade = np.ones(shape)*3000\n",
    "l1 = 100\n",
    "densidade[l1:,:] = 2400\n",
    "velocidade[l1:,:] = 4000"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Novamente, crie a simulação adicione uma fonte explosiva no topo."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "sim2 = ElasticSH(velocidade, densidade, spacing=spacing, taper=0.005)\n",
    "sim2.add_point_source((0, shape[1]//2), Ricker(5, 20, 1/20))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Rode a simulação por 600 iterações."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "sim2.run(600)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Gere a animação abaixo. O que acontece quando a frente de onda chega na interface?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "sim2.animate(every=10, embed=True, dpi=50, cutoff=0.001)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Rode abaixo para ver cada tempo separado, com calma."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "sim2.explore(every=20)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "A célula abaixo **não precisa ser rodada**. Ela somente gera a figura que está no questionário. Para vocês verem que eu não roubei e fiz a figura no Paint."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "fig, axes = mpl.subplots(1, 3)\n",
    "for f, ax in zip([150, 300, 550], axes):\n",
    "    frame = sim2[f,:,:]\n",
    "    ranges = 0.5*np.abs([frame.min(), frame.max()]).max()\n",
    "    ax.imshow(frame, extent=extent, cmap=mpl.cm.Greys, vmin=-ranges, vmax=ranges)\n",
    "    ax.imshow(velocidade, extent=extent, alpha=0.3, cmap=mpl.cm.Greys)\n",
    "    ax.plot(spacing*sim2.sources[0][0][1], 50 + spacing*sim2.sources[0][0][0], '*k', markersize=10)\n",
    "ax.set_xlim(extent[:2])\n",
    "fig.set_size_inches(12, 3.5)\n",
    "fig.tight_layout(pad=0)\n",
    "fig.savefig('../_static/img/ondas-sismicas-tarefa-raios.png')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Ondas P e ondas S\n",
    "\n",
    "Vamos fazer uma simulação mais realista agora que inclui tanto ondas P quanto ondas S.\n",
    "\n",
    "Para isso, precisamos definir qual é a velocidade de cada uma das ondas. Note que a onda P é mais rápida que a S."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "dens = np.ones(shape)*2200\n",
    "vs = np.ones(shape)*3000\n",
    "vp = np.ones(shape)*4000"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Rode a célula abaixo para criar a simulação de ondas P e S. Dessa vez, **não vamos usar uma fonte explosiva**. A fonte será um **puxão com mergulho de 45 graus**. A fonte será localizada em profundidade a 1500m. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "sim3 = ElasticPSV(vp, vs, dens, spacing=spacing)\n",
    "sim3.add_point_source((shape[0]//2, shape[1]//2), dip=45, wavelet=Ricker(5, 20, 1/20))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Agora, avance a simulação 400 iterações no tempo. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "sim3.run(400)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "E gere a animação. Note que essa é diferente da outra. Agora, há vetores na imagem além das frentes de onda coloridas. **Os vetores representam a direção e intensidade das vibrações em todo o espaço**. Os pontos são vetores com tamanho zero. Note que no começo da animação aparece um vetor bem grande no centro. Esse vetor representa a fonte da onda, um puxão na direção de 45 graus.\n",
    "\n",
    "Tente achar as ondas P e as S. Qual frente de onda é mais rápida? Qual é a direção de propagação de cada uma? E a direção das vibrações?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "sim3.animate(every=10, plottype=['vectors', 'wavefield'], scale=150, dpi=70, embed=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Novamente, use a célula abaixo para ver cada etapa da simulação separadamente."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "sim3.explore(every=20, plottype=['vectors', 'wavefield'], scale=150)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "A célula abaixo **não precisa ser rodada**. Ela somente gera a figura que está no questionário."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "fig = mpl.figure()\n",
    "ax = mpl.subplot(111)    \n",
    "ax.plot(spacing*sim3.sources[0][0][1], spacing*sim3.sources[0][0][0], '*k', markersize=20)\n",
    "sim3.snapshot(-1, plottype=['vectors', 'wavefield'], every_particle=5, scale=200, ax=ax)\n",
    "fig.savefig('../_static/img/ondas-sismicas-tarefa-p-s.png')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "## Onda P incidindo em uma interface (completo)\n",
    "\n",
    "Agora que vocês acharam as ondas P e S, vamos fazer a simulação completa do que acontece quando uma onda P encontra uma interface.\n",
    "\n",
    "Rode a célula abaixo para criar um modelo com dois meios com velocidades diferentes. A interface entre os dois meios está a 1000m."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "dens = np.ones(shape)*2000\n",
    "vs = np.ones(shape)*2500\n",
    "vp = np.ones(shape)*4000\n",
    "l = 150\n",
    "dens[l:, :] = 2400\n",
    "vs[l:, :] = 3500\n",
    "vp[l:, :] = 5000"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "A simulação que vamos fazer agora usa uma **fonte explosiva**.  Esse tipo de fonte gera **somente ondas P**. Vamos ver o que acontece quando essa onda encontra a interface. Para variar um pouco, vamos colocar a fonte no fundo do nosso modelo."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "sim4 = ElasticPSV(vp, vs, dens, spacing)\n",
    "sim4.add_blast_source((3*shape[0]//4, shape[1]//2), Ricker(5, 20, 1/20))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Agora, avance a simulação."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "sim4.run(700)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "E gere a animação. Note que no começo há somente ondas P (pergunta, como você sabe que são só ondas P?)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "sim4.animate(every=10, plottype=['vectors', 'wavefield'], cutoff=0.05, every_particle=6, scale=150, dpi=70, embed=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "sim4.explore(every=20, plottype=['vectors', 'wavefield'], scale=150)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Curiosidade extra\n",
    "\n",
    "Avance a simulação no tempo novamente e atualize a animação. O que acontece quando a onda P reflete na superfície?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "A célula abaixo **não precisa ser rodada**. Ela somente gera a figura que está no questionário."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "fig, axes = mpl.subplots(1, 3)\n",
    "for f, ax in zip([250, 400, 650], axes):\n",
    "    sim4.snapshot(f, plottype=['vectors', 'wavefield'], every_particle=5, scale=100, ax=ax, cmap=mpl.cm.Greys)\n",
    "    fig.delaxes(fig.axes[-1])\n",
    "    ax.imshow(vs, extent=extent, alpha=0.2, cmap=mpl.cm.Greys)\n",
    "    ax.plot(spacing*sim4.sources[0][0][1], spacing*sim4.sources[0][0][0], '*k', markersize=10)\n",
    "    ax.set_xlabel('')\n",
    "    ax.set_ylabel('')\n",
    "fig.set_size_inches(12., 3.5)\n",
    "fig.tight_layout(pad=0, w_pad=0)\n",
    "fig.savefig('../_static/img/ondas-sismicas-tarefa-raios-p-s.png')"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 2",
   "language": "python",
   "name": "python2"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.9"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}