{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Course website**: http://www.leouieda.com/geofisica1\n",
    "\n",
    "**Note**: This notebook is part of the course \"Geofísica 1\" 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 [Jupyter notebook](http://jupyter.org/). É um documento interativo que mistura texto (como esse), código (como abaixo), e o resultado de executar o código (números, texto, figuras, videos, etc)."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Prática 1 - Gravimetria - Interpolação, mapas e a gravidade da Terra"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Objetivos\n",
    "\n",
    "* Entender a influência da interpolação na geração de mapas de dados geofísicos\n",
    "* Visualizar as variações geográficas da gravidade da Terra\n",
    "* Entender como a escala de cores utilizada nos mapas influencia nossa interpretação\n",
    "* Aprender quais são os fatores que devem ser considerados quando visualizamos um dado em mapa"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Instruções\n",
    "\n",
    "O notebook te fornecerá exemplos interativos que trabalham os temas abordados no questionário. Utilize esses exemplos para responder as perguntas.\n",
    "\n",
    "As células com números ao lado, como `In [1]:`, são código [Python](http://python.org/). Algumas dessas células não produzem resultado e servem de preparação para os exemplos interativos. Outras, produzem gráficos interativos. **Você deve executar todas as células, uma de cada vez**, mesmo as que não produzem gráficos.\n",
    "\n",
    "Para executar uma célula, clique em cima dela e aperte `Shift + Enter`. O foco (contorno verde ou cinza em torno da célula) deverá passar para a célula abaixo. Para rodá-la, aperte `Shift + Enter` novamente e assim por diante. Você pode executar células de texto que não acontecerá nada."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Preparação\n",
    "\n",
    "Exectute as células abaixo para carregar as componentes necessárias para nossa prática. Vamos utilizar várias *bibliotecas*, inclusive uma de geofísica chamada [Fatiando a Terra](http://www.fatiando.org)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "%matplotlib inline\n",
    "from __future__ import division\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "from mpl_toolkits.basemap import Basemap\n",
    "import ipywidgets as widgets\n",
    "from IPython.display import display\n",
    "import seaborn\n",
    "from fatiando import utils, gridder\n",
    "import fatiando\n",
    "from icgem import load_icgem_gdf, down_sample"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "print(\"Usando a versão do Fatiando a Terra: {}\".format(fatiando.__version__))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Interpolação"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "O melhor jeito de entendermos o efeito da interpolação é fabricando alguns dados fictícios (sintéticos).\n",
    "Assim, podemos gerar os dados tanto em pontos aleatórios quanto em um grid regular.\n",
    "Isso nos permite comparar os resultados da interpolação com o *verdadeiro*. Nosso verdadeiro será um conjunto de dados medidos em um grid regular. Como se tivéssemos ido ao campo e medido em um grid regular."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Rode a célula abaixo para gerar os dados em pontos aleatórios e em um grid regular."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "area = (-5000., 5000., -5000., 5000.)\n",
    "shape = (100, 100)\n",
    "xp, yp = gridder.scatter(area, 100, seed=6)\n",
    "x, y = [i.reshape(shape) for i in gridder.regular(area, shape)]\n",
    "aletatorio = 50*utils.gaussian2d(xp, yp, 10000, 1000, angle=45)\n",
    "regular = 50*utils.gaussian2d(x, y, 10000, 1000, angle=45).reshape(shape)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Rode as duas células abaixo para gerar um gráfico interativo. Nesse gráfico você poderá controlar:\n",
    "\n",
    "* O número de pontos (em x e y) do grid utilizado na interpolação (`num_pontos`)\n",
    "* O método de interpolação utilizado (`metodo`). Pode ser interpolação cúbica ou linear.\n",
    "* Mostrar ou não os pontos de medição aleatórios no mapa interpolado.\n",
    "\n",
    "**Repare no que acontece com as bordas do mapa e onde não há observações**."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "def interpolacao(num_pontos, metodo, pontos_medidos):\n",
    "    fig, axes = plt.subplots(1, 2, figsize=(14, 6))\n",
    "    ishape = (num_pontos, num_pontos)\n",
    "    tmp = gridder.interp(yp, xp, aletatorio, ishape, area=area, algorithm=metodo, extrapolate=True)\n",
    "    yi, xi, interp = [i.reshape(ishape) for i in tmp]\n",
    "    ranges = np.abs([interp.min(), interp.max()]).max()\n",
    "    kwargs = dict(cmap=\"RdBu_r\", vmin=-ranges, vmax=ranges)\n",
    "    ax = axes[0]\n",
    "    ax.set_title(u'Pontos medidos')\n",
    "    ax.set_aspect('equal')\n",
    "    tmp = ax.scatter(yp*0.001, xp*0.001, s=80, c=aletatorio, **kwargs)\n",
    "    plt.colorbar(tmp, ax=ax, aspect=50, pad=0.01)\n",
    "    ax.set_xlabel('y (km)')\n",
    "    ax.set_ylabel('x (km)')\n",
    "    ax.set_xlim(-5, 5)\n",
    "    ax.set_ylim(-5, 5)\n",
    "    plt.tight_layout(pad=0)\n",
    "    ax = axes[1]\n",
    "    ax.set_title(u'Interpolado')\n",
    "    ax.set_aspect('equal')\n",
    "    tmp = ax.contourf(yi*0.001, xi*0.001, interp, 40, **kwargs)\n",
    "    plt.colorbar(tmp, ax=ax, aspect=50, pad=0.01)\n",
    "    if pontos_medidos:\n",
    "        ax.plot(yp*0.001, xp*0.001, '.k')\n",
    "    ax.set_xlabel('y (km)')\n",
    "    ax.set_ylabel('x (km)')\n",
    "    ax.set_xlim(-5, 5)\n",
    "    ax.set_ylim(-5, 5)\n",
    "    plt.tight_layout(pad=0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "w = widgets.interactive(interpolacao, num_pontos=(5, 100, 5), metodo=['cubic', 'linear'], pontos_medidos=False)\n",
    "display(w)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Vamos verificar se alguma das combinações chegou perto do resultado *verdadeiro*.\n",
    "\n",
    "Rode a célula abaixo para gerar um gráfico dos dados verdadeiros (gerados em um grid regular). Esse deveria ser o resultado observado se a interpolação fosse perfeita."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "fig, ax = plt.subplots(1, 1, figsize=(7, 6))\n",
    "ranges = np.abs([regular.min(), regular.max()]).max()\n",
    "kwargs = dict(cmap=\"RdBu_r\", vmin=-ranges, vmax=ranges)\n",
    "ax.set_title(u'Verdadeiro')\n",
    "ax.set_aspect('equal')\n",
    "tmp = ax.contourf(y*0.001, x*0.001, regular, 40, **kwargs)\n",
    "plt.colorbar(tmp, ax=ax, aspect=50, pad=0.01)\n",
    "ax.plot(yp*0.001, xp*0.001, '.k')\n",
    "ax.set_xlabel('y (km)')\n",
    "ax.set_ylabel('x (km)')\n",
    "plt.tight_layout(pad=0)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Gravidade do mundo"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Vamos visualizar como a gravidade da Terra varia geograficamente. Os dados da gravidade do mundo foram baixados de http://icgem.gfz-potsdam.de/ICGEM/potato/Service.html usando o modelo EIGEN-6c3stat.\n",
    "\n",
    "**As medições foram feitas em cima da superfície da Terra**, ou seja, acompanhando a topografia."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Rode as células abaixo para carregar os dados."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "dados = load_icgem_gdf('../data/eigen-6c3stat-0_5-mundo.gdf')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "lat, lon, grav = dados['latitude'], dados['longitude'], dados['gravity_earth']"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Vamos fazer um mapa da gravidade utilizando a [projeção Mollweid](http://en.wikipedia.org/wiki/Map_projection). Esses dados estão em mGal: 1 mGal = 10⁻⁵ m/s². \n",
    "\n",
    "Rode as duas células abaixo para gerar o gráfico (isso pode demorar um pouco)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "bm = Basemap(projection='moll', lon_0=0, resolution='c')\n",
    "x, y = bm(lon, lat)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "plt.figure(figsize=(18, 10))\n",
    "tmp = bm.contourf(x, y, grav, 100, tri=True, cmap='Reds')\n",
    "plt.colorbar(orientation='horizontal', pad=0.01, aspect=50, shrink=0.5).set_label('mGal')\n",
    "plt.title(\"Gravidade medida na superficie da Terra\", fontsize=16)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Escala de cor\n",
    "\n",
    "A escala de cores que utilizamos para mapear os valores pode ter um impacto grande na nossa interpretação dos resultados. Abaixo, veremos como o nosso dado de gravidade mundial fica quando utilizamos diferentes escalas de cor.\n",
    "\n",
    "As escalas podem ser divididas em 3 categorias:\n",
    "\n",
    "* lineares: as cores variam de um tom claro (geralmente branco) a uma cor (por exemplo, vermelho) de maneira linear\n",
    "* divergente: as cores variam de uma cor escura, passando por um tom claro (geralmente branco), e depois para outra cor escura.\n",
    "* raindow ou qualitativos: as cores variam sem um padrão de intensidade claro. Podem ser as cores do arco-íris ou outra combinação.\n",
    "\n",
    "Nas escalas lineares e divergentes, as cores sempre variam de baixa intensidade para alta intensidade (e vice-versa para escalas divergentes)."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Rode as células abaixo para gerar um mapa interativo da gravidade mundial. Você poderá controlar qual escala de cor você quer usar. Experimente com elas e veja como elas afetam sua percepção.\n",
    "\n",
    "**Para pensar**: Como isso pode afetar alguem que é [daltônico](https://pt.wikipedia.org/wiki/Daltonismo)?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def grav_mundial(escala_de_cor):\n",
    "    plt.figure(figsize=(18, 10))\n",
    "    tmp = bm.contourf(x, y, grav, 100, tri=True, cmap=escala_de_cor)\n",
    "    plt.colorbar(orientation='horizontal', pad=0.01, aspect=50, shrink=0.5).set_label('mGal')\n",
    "    plt.title(\"Escala de cor: {}\".format(escala_de_cor), fontsize=16)  "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "escalas = 'Reds Blues Greys YlOrBr RdBu BrBG PRGn Dark2 jet ocean rainbow gnuplot'.split()\n",
    "w = widgets.interactive(grav_mundial, escala_de_cor=escalas)\n",
    "display(w)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "[]()"
   ]
  }
 ],
 "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.10"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}