{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Demonstration\n",
    "\n",
    "This notebook demostrates possible usages of the code and displays some interesting phenomenology of the memory."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import gwmemory\n",
    "\n",
    "import matplotlib as mpl\n",
    "import matplotlib.pyplot as plt\n",
    "import numpy as np\n",
    "\n",
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "mpl.rcParams[\"font.family\"] = \"serif\"\n",
    "mpl.rcParams[\"font.serif\"] = \"Computer Modern Roman\"\n",
    "mpl.rcParams[\"font.size\"] = 20\n",
    "mpl.rcParams[\"text.usetex\"] = True\n",
    "mpl.rcParams[\"grid.alpha\"] = 0"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Importance of Higher-Order Modes\n",
    "\n",
    "This figure demonstrates the importance of including higher-order oscillatory modes in the evaluation of the memory."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fig = plt.figure(figsize=(12, 6))\n",
    "\n",
    "qs = [1.0, 2.0]\n",
    "S1s = [[0.0, 0.0, 0.0], [0.8, 0.0, 0.0]]\n",
    "S2s = [[0.0, 0.0, 0.0], [0.0, 0.8, 0.0]]\n",
    "\n",
    "colours = [\"r\", \"b\", \"g\", \"k\"]\n",
    "\n",
    "times = np.linspace(-0.08, 0.02, 10001)\n",
    "\n",
    "labels = [\"Equal-mass, non-spinning\"]\n",
    "\n",
    "ax = [fig.add_subplot(2, 1, 1), fig.add_subplot(2, 1, 2)]\n",
    "\n",
    "parameters = dict(\n",
    "    total_mass=60, distance=400, model=\"NRSur7dq2\", inc=np.pi / 2, phase=0, times=times\n",
    ")\n",
    "\n",
    "for ii, q in enumerate(qs):\n",
    "    parameters[\"q\"] = q\n",
    "    for jj, (S1, S2) in enumerate(zip(S1s, S2s)):\n",
    "        colour = colours[ii + jj * 2]\n",
    "\n",
    "        parameters[\"spin_1\"] = S1\n",
    "        parameters[\"spin_2\"] = S2\n",
    "\n",
    "        h_mem, times = gwmemory.time_domain_memory(**parameters)\n",
    "\n",
    "        ax[0].plot(times, h_mem[\"plus\"] * 1e22, linestyle=\"-\", color=colour)\n",
    "        ax[1].plot(times, h_mem[\"cross\"] * 1e22, linestyle=\"-\", color=colour)\n",
    "\n",
    "        h_mem, times = gwmemory.time_domain_memory(\n",
    "            **parameters, Lmax=2, modes=[(2, 2), (2, -2)]\n",
    "        )\n",
    "\n",
    "        ax[0].plot(times, h_mem[\"plus\"] * 1e22, linestyle=\":\", color=colour)\n",
    "        ax[1].plot(times, h_mem[\"cross\"] * 1e22, linestyle=\":\", color=colour)\n",
    "\n",
    "plt.xlabel(\"$t (s)$\")\n",
    "ax[0].set_ylabel(\"$\\delta h_{+} \\, \\left[10^{-22}\\\\right]$\")\n",
    "ax[1].set_ylabel(\"$\\delta h_{\\\\times} \\, \\left[10^{-22}\\\\right]$\")\n",
    "\n",
    "ax[0].set_xlim(-0.02, 0.02)\n",
    "ax[1].set_xlim(-0.02, 0.02)\n",
    "\n",
    "plt.tight_layout()\n",
    "plt.show()\n",
    "plt.close()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Oscillatory Mode Dependance\n",
    "\n",
    "This figure demonstrates the impact of the different oscillatory modes on the memory."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "qs = [1.0, 2.0]\n",
    "S1s = [[0.0, 0.0, 0.0], [0.8, 0.0, 0.0]]\n",
    "S2s = [[0.0, 0.0, 0.0], [0.0, 0.8, 0.0]]\n",
    "\n",
    "colours = [\"r\", \"b\", \"g\", \"k\"]\n",
    "\n",
    "mems_plus = {}\n",
    "mems_cross = {}\n",
    "\n",
    "times = np.linspace(-0.08, 0.02, 10001)\n",
    "\n",
    "parameters = dict(total_mass=60, distance=400, times=times)\n",
    "\n",
    "for ii, q in enumerate(qs):\n",
    "    parameters[\"q\"] = q\n",
    "    for jj, (S1, S2) in enumerate(zip(S1s, S2s)):\n",
    "        parameters[\"spin_1\"] = S1\n",
    "        parameters[\"spin_2\"] = S2\n",
    "\n",
    "        surr = gwmemory.waveforms.Surrogate(**parameters)\n",
    "\n",
    "        hmem = {}\n",
    "\n",
    "        modes = []\n",
    "        for l in range(2, 5):\n",
    "            for m in np.flipud(range(0, l + 1)):\n",
    "                modes += list(set([(l, m), (l, -m)]))\n",
    "                hmem[f\"{l}{m}\"], times = surr.time_domain_memory(\n",
    "                    inc=np.pi / 2, phase=0, modes=modes\n",
    "                )\n",
    "\n",
    "        max_h_mem_plus = []\n",
    "        max_h_mem_cross = []\n",
    "        keys = []\n",
    "        for l in range(2, 5):\n",
    "            for m in np.flipud(range(0, l + 1)):\n",
    "                keys.append(f\"{l}{m}\")\n",
    "                max_h_mem_plus.append(hmem[f\"{l}{m}\"][\"plus\"][-1])\n",
    "                max_h_mem_cross.append(hmem[f\"{l}{m}\"][\"cross\"][-1])\n",
    "\n",
    "        mems_plus[f\"{q}{S1}{S2}\"] = np.array(max_h_mem_plus)\n",
    "        mems_cross[f\"{q}{S1}{S2}\"] = np.array(max_h_mem_cross)\n",
    "\n",
    "colours = [\"r\", \"b\", \"g\", \"k\"]\n",
    "\n",
    "spin_keys = [f\"{S1}{S2}\" for S1, S2 in zip(S1s, S2s)]\n",
    "\n",
    "fig = plt.figure(figsize=(12, 4))\n",
    "\n",
    "for ii, q in enumerate(qs):\n",
    "    for jj, (S1, S2) in enumerate(zip(S1s, S2s)):\n",
    "\n",
    "        key = f\"{q}{S1}{S2}\"\n",
    "        max_h_mems_plus = mems_plus[key]\n",
    "        max_h_mems_cross = mems_cross[key]\n",
    "\n",
    "        plt.plot(\n",
    "            range(len(max_h_mems_plus)),\n",
    "            max_h_mems_plus * 1e22,\n",
    "            color=colours[ii + jj * 2],\n",
    "            linestyle=\"-\",\n",
    "            label=f\"$q={q}$, $S_1={S1}$, $S_2={S2}$\",\n",
    "        )\n",
    "        plt.plot(\n",
    "            range(len(max_h_mems_cross)),\n",
    "            max_h_mems_cross * 1e22,\n",
    "            color=colours[ii + jj * 2],\n",
    "            linestyle=\"--\",\n",
    "        )\n",
    "\n",
    "plt.xticks(range(len(max_h_mems_plus)))\n",
    "keys = [f\"({ell},{delta_m})\" for ell, delta_m in keys]\n",
    "ax = plt.gca()\n",
    "ax.set_xticklabels(keys)\n",
    "plt.xlim(0, len(max_h_mems_plus) - 1)\n",
    "ax.grid(axis=\"both\")\n",
    "\n",
    "plt.xlabel(\"Oscillatory $(\\ell, |m|)_{\\mathrm{last}}$\")\n",
    "plt.ylabel(\"$\\delta h(t=\\infty) \\, \\left[10^{-22}\\\\right]$\")\n",
    "plt.legend(loc=\"lower right\")\n",
    "\n",
    "plt.tight_layout()\n",
    "plt.show()\n",
    "plt.close()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Memory Mode Dependance\n",
    "\n",
    "This figure shows the spin-weighted spherical harmonic decomposition of the memory waveform."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fig = plt.figure(figsize=(12, 4))\n",
    "\n",
    "qs = [1.0, 2.0]\n",
    "S1s = [[0.0, 0.0, 0.0], [0.8, 0.0, 0.0]]\n",
    "S2s = [[0.0, 0.0, 0.0], [0.0, 0.8, 0.0]]\n",
    "\n",
    "colours = [\"r\", \"b\", \"g\", \"k\"]\n",
    "\n",
    "times = np.linspace(-0.08, 0.02, 10001)\n",
    "\n",
    "parameters = dict(total_mass=60, distance=400, times=times, model=\"NRSur7dq2\")\n",
    "\n",
    "for ii, q in enumerate(qs):\n",
    "    parameters[\"q\"] = q\n",
    "    for jj, (S1, S2) in enumerate(zip(S1s, S2s)):\n",
    "        parameters[\"spin_1\"] = S1\n",
    "        parameters[\"spin_2\"] = S2\n",
    "\n",
    "        h_mem, times = gwmemory.time_domain_memory(**parameters)\n",
    "\n",
    "        colour = colours[ii + jj * 2]\n",
    "\n",
    "        for ell, delta_m in h_mem:\n",
    "            if ell <= 4:\n",
    "                if ell == 2 and delta_m == 0:\n",
    "                    plt.scatter(\n",
    "                        ell**2 + delta_m + ell,\n",
    "                        abs(h_mem[(ell, delta_m)][-1]),\n",
    "                        marker=\"x\",\n",
    "                        s=200,\n",
    "                        color=colour,\n",
    "                        label=f\"$q={q}$, $S_1={S1}$, $S_2={S2}$\",\n",
    "                    )\n",
    "                else:\n",
    "                    plt.scatter(\n",
    "                        ell**2 + delta_m + ell,\n",
    "                        abs(h_mem[(ell, delta_m)][-1]),\n",
    "                        marker=\"x\",\n",
    "                        s=200,\n",
    "                        color=colour,\n",
    "                    )\n",
    "                plt.plot(\n",
    "                    [ell**2 + delta_m + ell, ell**2 + delta_m + ell],\n",
    "                    [0, abs(h_mem[(ell, delta_m)][-1])],\n",
    "                    color=colour,\n",
    "                    alpha=0.2,\n",
    "                )\n",
    "\n",
    "plt.xlim(3, 5**2)\n",
    "plt.xticks(range(4, 5**2, 1), rotation=270)\n",
    "keys = [f\"({ell, delta_m})\" for ell, delta_m in gwmemory.harmonics.lmax_modes(4)]\n",
    "ax = plt.gca()\n",
    "ax.set_xticklabels(keys)\n",
    "plt.xlabel(\"Memory $(\\ell, m)$\")\n",
    "\n",
    "plt.ylim(1e-27, 1e-21)\n",
    "plt.yscale(\"log\")\n",
    "plt.ylabel(\"$|\\delta h_{\\ell m}(t=\\infty)|$\")\n",
    "\n",
    "ax = plt.gca()\n",
    "ax.set_yticks(np.logspace(-27, -21, 7))\n",
    "ax.grid(axis=\"y\")\n",
    "\n",
    "plt.tight_layout()\n",
    "plt.show()\n",
    "plt.close()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Orientation Dependance\n",
    "\n",
    "Next we plot the orientation dependance of the + and x polarized memory as a function of source inclination (zenith) and orbital phase (azimuth).\n",
    "\n",
    "*NOTE: this plot requires [`Basemap`](https://github.com/matplotlib/basemap) for the projection.*\n",
    "\n",
    "Install with\n",
    "```console\n",
    "$ conda install -c conda-forge basemap\n",
    "```"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from mpl_toolkits.basemap import Basemap\n",
    "\n",
    "fig = plt.figure(figsize=(24, 10))\n",
    "\n",
    "qs = [1.0]\n",
    "S1s = [[0.0, 0.0, 0.0], [0.8, 0.0, 0.0]]\n",
    "S2s = [[0.0, 0.0, 0.0], [0.0, 0.8, 0.0]]\n",
    "\n",
    "labels = [\"Non-spinning\", \"Precessing\"]\n",
    "\n",
    "parameters = dict(total_mass=60, distance=400, times=times, model=\"NRSur7dq2\")\n",
    "\n",
    "for ii, q in enumerate(qs):\n",
    "    parameters[\"q\"] = q\n",
    "    for jj, (S1, S2) in enumerate(zip(S1s, S2s)):\n",
    "        parameters[\"spin_1\"] = S1\n",
    "        parameters[\"spin_2\"] = S2\n",
    "\n",
    "        h_mem_lm, times = gwmemory.time_domain_memory(**parameters)\n",
    "\n",
    "        inc_array = np.linspace(0, np.pi, 200) - np.pi / 2\n",
    "        pol_array = np.linspace(0, 2 * np.pi, 200) - np.pi\n",
    "        pols, incs = np.meshgrid(pol_array, inc_array)\n",
    "        pols_deg = pols * 180 / np.pi\n",
    "        incs_deg = incs * 180 / np.pi\n",
    "        y_lm = {\n",
    "            (ell, m): gwmemory.harmonics.sYlm(-2, ell, m, incs + np.pi / 2, pols)\n",
    "            for ell, m in h_mem_lm.keys()\n",
    "        }\n",
    "\n",
    "        orientation_map = np.sum(\n",
    "            [y_lm[key] * h_mem_lm[key][-1] for key in y_lm], axis=0\n",
    "        )\n",
    "\n",
    "        ax = plt.subplot(2, 2, ii * 4 + jj * 2 + 1)\n",
    "        m = Basemap(projection=\"moll\", lon_0=-180, resolution=\"c\")\n",
    "        m.contourf(\n",
    "            pols_deg,\n",
    "            incs_deg,\n",
    "            orientation_map.real,\n",
    "            100,\n",
    "            cmap=plt.cm.jet,\n",
    "            latlon=True,\n",
    "            levels=np.linspace(-2.5e-22, 2.5e-22, 201),\n",
    "        )\n",
    "\n",
    "        if jj == 0:\n",
    "            ax.annotate(\n",
    "                text=\"$+$-polarized\",\n",
    "                xy=(0.08, 0.92),\n",
    "                xycoords=\"axes fraction\",\n",
    "                fontsize=30,\n",
    "                horizontalalignment=\"center\",\n",
    "            )\n",
    "\n",
    "        ax = plt.subplot(2, 2, ii * 4 + jj * 2 + 2)\n",
    "        m = Basemap(projection=\"moll\", lon_0=-180, resolution=\"c\")\n",
    "        m.contourf(\n",
    "            pols_deg,\n",
    "            incs_deg,\n",
    "            orientation_map.imag,\n",
    "            100,\n",
    "            cmap=plt.cm.jet,\n",
    "            latlon=True,\n",
    "            levels=np.linspace(-2.5e-22, 2.5e-22, 201),\n",
    "        )\n",
    "        if jj == 0:\n",
    "            ax.annotate(\n",
    "                text=\"$\\\\times$-polarized\",\n",
    "                xy=(0.08, 0.92),\n",
    "                xycoords=\"axes fraction\",\n",
    "                fontsize=30,\n",
    "                horizontalalignment=\"center\",\n",
    "            )\n",
    "        ax.annotate(\n",
    "            text=labels[jj],\n",
    "            xy=(0.9, 0.92),\n",
    "            xycoords=\"axes fraction\",\n",
    "            fontsize=30,\n",
    "            horizontalalignment=\"center\",\n",
    "        )\n",
    "\n",
    "plt.tight_layout()\n",
    "plt.show()\n",
    "plt.close()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Memory of Memory\n",
    "\n",
    "This figure shows the importance of each order of the waveform."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fig = plt.figure(figsize=(12, 4))\n",
    "\n",
    "qs = [1.0, 2.0]\n",
    "S1s = [[0.0, 0.0, 0.0], [0.8, 0.0, 0.0]]\n",
    "S2s = [[0.0, 0.0, 0.0], [0.0, 0.8, 0.0]]\n",
    "\n",
    "colours = [\"r\", \"b\", \"g\", \"k\"]\n",
    "\n",
    "times = np.linspace(-0.08, 0.02, 10001)\n",
    "\n",
    "parameters = dict(total_mass=60, distance=400, times=times)\n",
    "\n",
    "for ii, q in enumerate(qs):\n",
    "    parameters[\"q\"] = q\n",
    "    for jj, (S1, S2) in enumerate(zip(S1s, S2s)):\n",
    "        parameters[\"spin_1\"] = S1\n",
    "        parameters[\"spin_2\"] = S2\n",
    "\n",
    "        surr = gwmemory.waveforms.Surrogate(**parameters)\n",
    "        osc = gwmemory.utils.combine_modes(surr.h_lm, np.pi / 2, 0)\n",
    "        modes = gwmemory.harmonics.lmax_modes(4)\n",
    "        old_h_mem = np.array(0 * (1 + 1j))\n",
    "        delta_h = [max(abs(osc[\"plus\"] - 1j * osc[\"cross\"]))]\n",
    "\n",
    "        for kk in range(1, 5):\n",
    "            h_mem, times = surr.time_domain_memory(np.pi / 2, 0)\n",
    "            delta_h.append(h_mem[\"plus\"][-1] - 1j * h_mem[\"cross\"][-1] - old_h_mem)\n",
    "            old_h_mem = h_mem[\"plus\"][-1] - 1j * h_mem[\"cross\"][-1]\n",
    "            h_mem_lm, times = surr.time_domain_memory()\n",
    "            surr = gwmemory.waveforms.Surrogate(**parameters)\n",
    "            for key in h_mem_lm:\n",
    "                if key[0] <= 4:\n",
    "                    surr.h_lm[key] += h_mem_lm[key]\n",
    "\n",
    "        colour = colours[ii + jj * 2]\n",
    "\n",
    "        plt.semilogy(\n",
    "            range(5),\n",
    "            np.abs(delta_h),\n",
    "            marker=\"x\",\n",
    "            label=f\"$q={q}$, $S_1={S1}$, $S_2={S2}$\",\n",
    "            color=colour,\n",
    "        )\n",
    "        plt.xlabel(\"Memory order, $i$\")\n",
    "        plt.ylabel(\"$\\max(|h^i|)$\")\n",
    "\n",
    "plt.xticks(range(5))\n",
    "plt.legend(loc=\"lower left\")\n",
    "\n",
    "ax = plt.gca()\n",
    "ax.grid(True, axis=\"y\")\n",
    "\n",
    "plt.tight_layout()\n",
    "plt.show()\n",
    "plt.close()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Other Examples\n",
    "\n",
    "Below are some additional examples of using the code."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Compare Waveform models\n",
    "\n",
    "Compare the memory predicted using the available waveforms.\n",
    "\n",
    "*FIXME: add references*"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fig = plt.figure(figsize=(12, 4))\n",
    "\n",
    "q = 1.0\n",
    "S1 = [0.0, 0.0, 0.0]\n",
    "S2 = [0.0, 0.0, 0.0]\n",
    "\n",
    "colours = [\"r\", \"b\", \"g\", \"k\"]\n",
    "\n",
    "parameters = dict(\n",
    "    q=1, spin_1=S1, spin_2=S2, total_mass=60, distance=400, inc=np.pi / 2, phase=0.0\n",
    ")\n",
    "\n",
    "\n",
    "for ii, model in enumerate([\"NRSur7dq2\", \"IMRPhenomD\", \"SEOBNRv4\", \"MWM\"]):\n",
    "    h_mem, times = gwmemory.time_domain_memory(**parameters, model=model)\n",
    "\n",
    "    plt.plot(\n",
    "        times,\n",
    "        h_mem[\"plus\"] - h_mem[\"plus\"][np.argmin(abs(times + 0.25))],\n",
    "        linestyle=\"-\",\n",
    "        color=colours[ii],\n",
    "        label=model,\n",
    "    )\n",
    "    plt.plot(\n",
    "        times,\n",
    "        h_mem[\"cross\"] - h_mem[\"cross\"][np.argmin(abs(times + 0.25))],\n",
    "        linestyle=\"--\",\n",
    "        color=colours[ii],\n",
    "    )\n",
    "\n",
    "plt.xlabel(\"$t (s)$\")\n",
    "plt.ylabel(\"$\\delta h$\")\n",
    "plt.legend(loc=\"upper left\", fontsize=20)\n",
    "\n",
    "plt.xlim(-0.25, 0.02)\n",
    "\n",
    "plt.tight_layout()\n",
    "plt.show()\n",
    "plt.close()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Oscillatory and Memory Waveform\n",
    "\n",
    "In order to generate the memory we necessarily generate the time-domain oscillatory waveform."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "q = 1.0\n",
    "S1 = [0.0, 0.0, 0.0]\n",
    "S2 = [0.0, 0.0, 0.0]\n",
    "\n",
    "times = np.linspace(-0.08, 0.02, 10001)\n",
    "surr = gwmemory.waveforms.Surrogate(\n",
    "    q=q, spin_1=S1, spin_2=S2, total_mass=60, distance=400, times=times\n",
    ")\n",
    "\n",
    "inc = np.pi / 2\n",
    "pol = 0\n",
    "\n",
    "oscillatory, times = surr.time_domain_oscillatory(inc=inc, phase=pol)\n",
    "memory, times = surr.time_domain_memory(inc=inc, phase=pol)\n",
    "\n",
    "fig = plt.figure(figsize=(12, 6))\n",
    "fig.add_subplot(2, 1, 1)\n",
    "plt.plot(times, oscillatory[\"plus\"], linestyle=\"--\", color=\"b\", alpha=0.5)\n",
    "plt.plot(times, oscillatory[\"cross\"], linestyle=\"--\", color=\"r\", alpha=0.5)\n",
    "plt.plot(times, memory[\"plus\"], linestyle=\"-.\", color=\"b\", alpha=0.5)\n",
    "plt.plot(times, memory[\"cross\"], linestyle=\"-.\", color=\"r\", alpha=0.5)\n",
    "plt.plot(times, oscillatory[\"plus\"] + memory[\"plus\"], color=\"b\")\n",
    "plt.plot(times, oscillatory[\"cross\"] + memory[\"cross\"], color=\"r\")\n",
    "plt.xlabel(\"$t [s]$\")\n",
    "plt.ylabel(\"Strain\")\n",
    "plt.axhline(0, linestyle=\":\", color=\"k\")\n",
    "plt.xlim(-0.08, 0.0)\n",
    "\n",
    "fig.add_subplot(2, 1, 2)\n",
    "plt.plot(times, oscillatory[\"plus\"], linestyle=\"--\", color=\"b\", alpha=0.5)\n",
    "plt.plot(times, oscillatory[\"cross\"], linestyle=\"--\", color=\"r\", alpha=0.5)\n",
    "plt.plot(times, memory[\"plus\"], linestyle=\"-.\", color=\"b\", alpha=0.5)\n",
    "plt.plot(times, memory[\"cross\"], linestyle=\"-.\", color=\"r\", alpha=0.5)\n",
    "plt.plot(times, oscillatory[\"plus\"] + memory[\"plus\"], color=\"b\")\n",
    "plt.plot(times, oscillatory[\"cross\"] + memory[\"cross\"], color=\"r\")\n",
    "plt.xlabel(\"$t [s]$\")\n",
    "plt.ylabel(\"Strain\")\n",
    "plt.axhline(0, linestyle=\":\", color=\"k\")\n",
    "plt.xlim(-0.0, 0.02)\n",
    "\n",
    "plt.tight_layout()\n",
    "plt.show()\n",
    "plt.close()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Frequency Domain Memory\n",
    "\n",
    "\n",
    "For graviational-wave parameter estimation we will require frequency domain waveforms. This is currently implemented by performing a fast fourier transform on the time domain waveform."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fig = plt.figure(figsize=(12, 4))\n",
    "\n",
    "q = 1.0\n",
    "S1 = [0, 0.8, 0]\n",
    "S2 = [0.8, 0, 0]\n",
    "\n",
    "times = np.linspace(-0.98, 0.02, 10000)\n",
    "\n",
    "colours = [\"r\", \"b\", \"g\", \"k\"]\n",
    "\n",
    "h_mem, frequencies = gwmemory.frequency_domain_memory(\n",
    "    q=q,\n",
    "    spin_1=S1,\n",
    "    spin_2=S2,\n",
    "    total_mass=60.0,\n",
    "    distance=400.0,\n",
    "    model=\"NRSur7dq2\",\n",
    "    inc=np.pi / 2,\n",
    "    phase=0.0,\n",
    "    times=times,\n",
    ")\n",
    "\n",
    "plt.loglog(frequencies, abs(h_mem[\"plus\"]), linestyle=\"-\", color=\"r\", label=model)\n",
    "plt.loglog(frequencies, abs(h_mem[\"cross\"]), linestyle=\"--\", color=\"r\")\n",
    "\n",
    "plt.xlabel(\"$f$ (Hz)\")\n",
    "plt.ylabel(\"Strain [Hz$^{-1/2}$]\")\n",
    "\n",
    "plt.xlim(10, 2048)\n",
    "plt.ylim(1e-27, 1e-23)\n",
    "\n",
    "ax = plt.gca()\n",
    "ax.grid(True, axis=\"y\")\n",
    "\n",
    "plt.tight_layout()\n",
    "plt.show()\n",
    "plt.close()"
   ]
  }
 ],
 "metadata": {},
 "nbformat": 4,
 "nbformat_minor": 4
}