{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# read src_rec file and check its arrival time and distance from src to rec\n",
    "import numpy as np\n",
    "\n",
    "def read_src_rec_file(fpath):\n",
    "\n",
    "    with open(fpath, \"r\") as f:\n",
    "        lines = f.readlines()\n",
    "\n",
    "        # list to store src lon lat\n",
    "        src_pos = []\n",
    "        # list to store rec lon lat\n",
    "        rec_pos_tmp=[]\n",
    "        rec_pos = []\n",
    "\n",
    "        for line in lines:\n",
    "            # src line if there are 13 elements\n",
    "            if len(line.split()) == 13:\n",
    "                # store source lon lat dep\n",
    "                stlon = float(line.split()[8])\n",
    "                stlat = float(line.split()[7])\n",
    "                src_pos.append([stlon, stlat])\n",
    "\n",
    "                nrec = float(line.split()[11])\n",
    "            # rec line if there are 9 elements\n",
    "            elif len(line.split()) == 9:\n",
    "                # store receiver lon lat dep\n",
    "                rclon = float(line.split()[4])\n",
    "                rclat = float(line.split()[3])\n",
    "                rctime= float(line.split()[8])\n",
    "                # calc epicentral distance from src\n",
    "                dist = np.sqrt((stlon-rclon)**2 + (stlat-rclat)**2)\n",
    "                rec_pos_tmp.append([rclon, rclat, rctime, dist])\n",
    "\n",
    "                nrec-=1\n",
    "            else:\n",
    "                raise ValueError(\"src_rec_test_out.dat file is not correct\")\n",
    "\n",
    "            if nrec==0:\n",
    "                rec_pos.append(rec_pos_tmp)\n",
    "                # remove all rec_pos_tmp\n",
    "                rec_pos_tmp = []\n",
    "\n",
    "\n",
    "\n",
    "    src_pos = np.array(src_pos)\n",
    "    rec_pos = np.array(rec_pos)\n",
    "\n",
    "    return src_pos, rec_pos"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "src_pos_true, rec_pos_true = read_src_rec_file(\"./src_rec_test_out.dat\")\n",
    "src_pos_try, rec_pos_try = read_src_rec_file(\"./src_rec_test_out_out.dat\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# plot rec_pos for each src_pos with color by 3rd element\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "\n",
    "def plot_srcrec(src_pos, rec_pos):\n",
    "\n",
    "    # selected source\n",
    "    id_src = 1\n",
    "    \n",
    "    fig, axs = plt.subplots(2,1,figsize=(10,10))\n",
    "    \n",
    "    # plot arrival time\n",
    "    # src\n",
    "    axs[0].scatter(src_pos[id_src,0], src_pos[id_src,1], c='r', marker='o', s=100)\n",
    "    # rec\n",
    "    axs[0].scatter(rec_pos[id_src][:,0], rec_pos[id_src][:,1], c=rec_pos[id_src][:,2], s=10, marker='o')\n",
    "    \n",
    "    # colorbar\n",
    "    cbar = plt.colorbar(axs[0].scatter(rec_pos[id_src][:,0], rec_pos[id_src][:,1], c=rec_pos[id_src][:,2], s=100, marker='o'))\n",
    "    cbar.set_label('arrival time')\n",
    "    \n",
    "    # plot epicentral distance\n",
    "    # src\n",
    "    axs[1].scatter(src_pos[id_src,0], src_pos[id_src,1], c='r', marker='o', s=100)\n",
    "    # rec\n",
    "    axs[1].scatter(rec_pos[id_src][:,0], rec_pos[id_src][:,1], c=rec_pos[id_src][:,3], s=100, marker='o')\n",
    "    \n",
    "    # colorbar\n",
    "    cbar = plt.colorbar(axs[1].scatter(rec_pos[id_src][:,0], rec_pos[id_src][:,1], c=rec_pos[id_src][:,3], s=10, marker='o'))\n",
    "    cbar.set_label('epicentral distance')\n",
    "    \n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plot_srcrec(src_pos_true, rec_pos_true)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plot_srcrec(src_pos_try, rec_pos_try)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "rec_pos_diff = rec_pos_true[:][:][2]-rec_pos_try[:][:][2]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "rec_pos_diff.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "rec_pos_true.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "interpreter": {
   "hash": "02f83e1f4cd9619657a6845405e2dd67c4de23753ba48bca5dce2ebf57b3813a"
  },
  "kernelspec": {
   "display_name": "Python 3.9.1 64-bit ('3.9.1')",
   "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.9.1"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}