{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## See the signal in the data ##\n",
    "\n",
    "Because GW150914 was relatively loud, it is readily visible in the data. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Requirement already satisfied: pycbc in /home/ahnitz/projects/PyCBC-Tutorials/env/lib/python3.7/site-packages/PyCBC-9rc6331-py3.7-linux-x86_64.egg (9rc6331)\n",
      "Requirement already satisfied: lalsuite in /home/ahnitz/projects/PyCBC-Tutorials/env/lib/python3.7/site-packages (6.70)\n",
      "Requirement already satisfied: ligo-common in /home/ahnitz/projects/PyCBC-Tutorials/env/lib/python3.7/site-packages (1.0.3)\n",
      "Requirement already satisfied: numpy>=1.16.0 in /home/ahnitz/projects/PyCBC-Tutorials/env/lib/python3.7/site-packages (from pycbc) (1.19.0)\n",
      "Requirement already satisfied: Mako>=1.0.1 in /home/ahnitz/projects/PyCBC-Tutorials/env/lib/python3.7/site-packages (from pycbc) (1.1.3)\n",
      "Requirement already satisfied: cython>=0.29 in /home/ahnitz/projects/PyCBC-Tutorials/env/lib/python3.7/site-packages (from pycbc) (0.29.20)\n",
      "Requirement already satisfied: decorator>=3.4.2 in /home/ahnitz/projects/PyCBC-Tutorials/env/lib/python3.7/site-packages (from pycbc) (4.4.2)\n",
      "Requirement already satisfied: matplotlib>=1.5.1 in /home/ahnitz/projects/PyCBC-Tutorials/env/lib/python3.7/site-packages (from pycbc) (3.2.2)\n",
      "Requirement already satisfied: pillow in /home/ahnitz/projects/PyCBC-Tutorials/env/lib/python3.7/site-packages (from pycbc) (7.1.2)\n",
      "Requirement already satisfied: h5py>=2.5 in /home/ahnitz/projects/PyCBC-Tutorials/env/lib/python3.7/site-packages (from pycbc) (2.10.0)\n",
      "Requirement already satisfied: jinja2 in /home/ahnitz/projects/PyCBC-Tutorials/env/lib/python3.7/site-packages (from pycbc) (2.11.2)\n",
      "Requirement already satisfied: mpld3>=0.3 in /home/ahnitz/projects/PyCBC-Tutorials/env/lib/python3.7/site-packages (from pycbc) (0.5.1)\n",
      "Requirement already satisfied: lscsoft-glue>=1.59.3 in /home/ahnitz/projects/PyCBC-Tutorials/env/lib/python3.7/site-packages (from pycbc) (2.0.0)\n",
      "Requirement already satisfied: emcee==2.2.1 in /home/ahnitz/projects/PyCBC-Tutorials/env/lib/python3.7/site-packages (from pycbc) (2.2.1)\n",
      "Requirement already satisfied: requests>=1.2.1 in /home/ahnitz/projects/PyCBC-Tutorials/env/lib/python3.7/site-packages (from pycbc) (2.24.0)\n",
      "Requirement already satisfied: beautifulsoup4>=4.6.0 in /home/ahnitz/projects/PyCBC-Tutorials/env/lib/python3.7/site-packages (from pycbc) (4.9.1)\n",
      "Requirement already satisfied: six>=1.10.0 in /home/ahnitz/projects/PyCBC-Tutorials/env/lib/python3.7/site-packages (from pycbc) (1.15.0)\n"
     ]
    }
   ],
   "source": [
    "import sys\n",
    "!{sys.executable} -m pip install pycbc ligo-common --no-cache-dir"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### See GW150914 in the time series ###"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%matplotlib inline\n",
    "from pycbc.filter import highpass_fir, lowpass_fir\n",
    "from pycbc.psd import welch, interpolate\n",
    "import pylab\n",
    "\n",
    "from pycbc.catalog import Merger\n",
    "\n",
    "for ifo in ['H1', 'L1']:\n",
    "    # Read data: The data is available to the public through losc.ligo.org!\n",
    "    ts = Merger(\"GW150914\").strain(ifo)\n",
    "    \n",
    "    # Estimate the noise spectrum and use it to whiten the data\n",
    "    psd = interpolate(welch(ts), 1.0 / ts.duration)\n",
    "    white_strain = (ts.to_frequencyseries() / psd ** 0.5).to_timeseries()\n",
    "\n",
    "    # remove frequencies below and above where the main portion of the signal lies\n",
    "    smooth = highpass_fir(white_strain, 35, 8)\n",
    "    smooth = lowpass_fir(smooth, 300, 8)\n",
    "\n",
    "    # time shift and flip L1 to match the time of arrival and phase of Hanford\n",
    "    if ifo == 'L1':\n",
    "        smooth *= -1\n",
    "        smooth.roll(int(.007 / smooth.delta_t))\n",
    "\n",
    "    pylab.plot(smooth.sample_times.numpy(), smooth.numpy(), label=ifo)\n",
    "\n",
    "# Plot the region around the signal (time is in gps seconds)\n",
    "pylab.legend()\n",
    "pylab.xlim(1126259462.26, 1126259462.48)\n",
    "pylab.ylabel('Smoothed-Whitened Strain')\n",
    "pylab.grid()\n",
    "pylab.ylim(-100, 100)\n",
    "pylab.xlabel('GPS Time (s)')\n",
    "pylab.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### See the time frequency evolution of GW150914 in Hanford ###"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "ts = Merger(\"GW150914\").strain(\"H1\")\n",
    "\n",
    "ts = ts.whiten(4, 4)\n",
    "zoom = ts.time_slice(1126259462.4 - 1, 1126259462.4 + 1)\n",
    "times, freqs, power = zoom.qtransform(.001, 1, frange=(20, 512), qrange=(4, 64))\n",
    "\n",
    "pylab.figure(figsize=(18, 3))\n",
    "pylab.pcolormesh(times, freqs, power)\n",
    "pylab.ylim(20, 512)\n",
    "pylab.xlabel('Time (s)')\n",
    "pylab.ylabel('Frequency (Hz)')\n",
    "pylab.xlim(times.min(), times.max())\n",
    "pylab.yscale('log')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### See the time frequency evolution of GW150914 in Livingston ###"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from pycbc.frame import read_frame\n",
    "ts = Merger(\"GW150914\").strain(\"L1\")\n",
    "\n",
    "ts = ts.whiten(4, 4)\n",
    "zoom = ts.time_slice(1126259462.4 - 1, 1126259462.4 + 1)\n",
    "times, freqs, power = zoom.qtransform(.001, 1, frange=(20, 512), qrange=(4, 64))\n",
    "\n",
    "pylab.figure(figsize=(18, 3))\n",
    "pylab.pcolormesh(times, freqs, power)\n",
    "pylab.ylim(20, 512)\n",
    "pylab.xlabel('Time (s)')\n",
    "pylab.ylabel('Frequency (Hz)')\n",
    "pylab.xlim(times.min(), times.max())\n",
    "pylab.yscale('log')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "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.7.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}