{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Why INCLINATION and SEMIMAJOR_AXIS\n",
    "* Both require external force to change\n",
    "* `SEMIMAJOR_AXIS` ignore variance in `APOAPSIS` and `PERIAPSIS`\n",
    "* Both remains constant in an elliptical orbit\n",
    "* Independent of each other\n",
    "* `INCLINATION` remains fairly constant over long periods of time\n",
    "* Note: `SEMIMAJOR_AXIS` decreases over time due to orbital decay.  In fact, a lot of maneuvers are orbit raising maneuvers to compensate for this.\n",
    "\n",
    "### Problems encountered\n",
    "* Data interval varies a lot - from 1 entry in 3 months to less than 2 hours between 2 TLEs\n",
    "* Duplicates in data\n",
    "* Some error in data (for example, random spikes)\n",
    "* Due to different engine types, some maneuvers are quick (happens between 2 TLEs) and some takes place over a long time (2 weeks for some starlinks).\n",
    "* The closeness of the satellite has a major and non-linear effect on its orbital decay rate\n",
    "* Some satellites constantly maneuvers to stay in orbit (ie: starlink), some let it decay and does big boosts (ISS)\n",
    "* It appears that `INCLINATION` suffers from periodic variance or just inaccuracies in measurement\n",
    "* Solar activity causes increase variance in decay"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pandas as pd\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "import matplotlib\n",
    "from tqdm import tqdm\n",
    "\n",
    "import importlib\n",
    "import detect_maneuver # put it in a .py so it can be used by the SOCRATES detection part also\n",
    "importlib.reload(detect_maneuver);"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "path = f'../../../siads591 data/filtered_raw/payload.pkl.gz' # path to the data file\n",
    "df = pd.read_pickle(path, compression=\"gzip\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "satcat = pd.read_csv(f'../../playground/satcat_all.csv')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "df = df.reset_index().sort_values(by=[\"NORAD_CAT_ID\",\"EPOCH\"]).drop_duplicates(subset=[\"NORAD_CAT_ID\",\"EPOCH\"]).set_index(\"EPOCH\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "maneuver_functions = {\n",
    "    'INCLINATION': [\n",
    "        (\"diff\", lambda x:x - x.shift(), [0.002,0.005,0.008]),\n",
    "#         (\"neighbors_diff\", lambda x:x.shift(-1) - x.shift(), [0.002,0.005,0.008]),\n",
    "        (\"rolling_4_neighbor_diff\", lambda x:x.rolling(4, min_periods=1).mean().shift(-3) - x.rolling(4, min_periods=1).mean(), [0.002,0.005,0.008]),\n",
    "#         (\"rolling_7_neighbor_diff\", lambda x:x.rolling(7, min_periods=1).mean().shift(-6) - x.rolling(7, min_periods=1).mean(), [0.002,0.005,0.008]),\n",
    "        (\"rolling_10_neighbor_diff\", lambda x:x.rolling(10, min_periods=1).mean().shift(-9) - x.rolling(10, min_periods=1).mean(), [0.002,0.005,0.008]),\n",
    "        (\"rolling_20_neighbor_diff\", lambda x:x.rolling(20, min_periods=1).mean().shift(-19) - x.rolling(20, min_periods=1).mean(), [0.002,0.005,0.008]),\n",
    "    ],\n",
    "    'SEMIMAJOR_AXIS': [\n",
    "        (\"diff\", lambda x:x - x.shift(), [0.008, 0.025, 0.06]),\n",
    "        (\"neighbors_diff\", lambda x:x.shift(-1) - x.shift(), [0.008, 0.025, 0.06]),\n",
    "        (\"rolling_3_neighbor_diff\", lambda x:x.rolling(3, min_periods=1).mean().shift(-2) - x.rolling(3, min_periods=1).mean(), [0.008, 0.025, 0.06]),\n",
    "        (\"rolling_5_neighbor_diff\", lambda x:x.rolling(5, min_periods=1).mean().shift(-4) - x.rolling(5, min_periods=1).mean(), [0.008, 0.025, 0.06]),\n",
    "#         (\"rolling_10_neighbor_diff\", lambda x:x.rolling(10, min_periods=1).mean().shift(-9) - x.rolling(10, min_periods=1).mean(), [0.008, 0.015, 0.06]),\n",
    "    ],\n",
    "}\n",
    "\n",
    "combined_maneuver_functions = {\n",
    "    'INCLINATION': [\n",
    "        (\"rolling_10_neighbor_diff\", lambda x:x.rolling(10, min_periods=1).mean().shift(-9) - x.rolling(10, min_periods=1).mean(), [0.008]),\n",
    "    ],\n",
    "    'SEMIMAJOR_AXIS': [\n",
    "        (\"rolling_3_neighbor_diff\", lambda x:x.rolling(3, min_periods=1).mean().shift(-2) - x.rolling(3, min_periods=1).mean(), [0.025]),\n",
    "    ],\n",
    "}\n",
    "\n",
    "satellite_and_range = [\n",
    "    (27424, slice(4300, 5500, None)), # AQUA 2008-10 to 2010-10\n",
    "    (33053, slice(1700, 2150, None)), # FERMI single event April 2012\n",
    "    (25544, slice(-2000, None, None)), # ISS past year\n",
    "    (38337, slice(1500,2500,None)), # GCOM W1 2016-3 to 2017-09\n",
    "    (40059, slice(1800,2200,None)), # OCO-2 2017-12 to 2018-05\n",
    "    (22675, slice(None, None, None)), # COSMOS 2251 - no propulsion\n",
    "    (39210, slice(100, 4000, None)), # PAYLOAD C - 2014 to 2019\n",
    "    (44713, slice(None, None, None)), # STARLINK-1007 - 2019-11 to current\n",
    "]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "100%|██████████| 8/8 [00:55<00:00,  6.88s/it]"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Done\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "\n"
     ]
    }
   ],
   "source": [
    "for norad_id, df_slice in tqdm(satellite_and_range):\n",
    "    raw, fixed, maneuver_results, fig = detect_maneuver.plot_maneuver_results(df, satcat, norad_id, df_slice, maneuver_functions)\n",
    "    plt.savefig(f'../../images/maneuver_analysis/{satcat.loc[(satcat.NORAD_CAT_ID == norad_id),\"SATNAME\"].values[0]} ({norad_id}).png', facecolor='white', transparent=False)\n",
    "    plt.close()\n",
    "    raw, fixed, maneuver_results, fig = detect_maneuver.plot_maneuver_results(df, satcat, norad_id, df_slice, combined_maneuver_functions, True)\n",
    "    plt.savefig(f'../../images/maneuver_analysis/{satcat.loc[(satcat.NORAD_CAT_ID == norad_id),\"SATNAME\"].values[0]} ({norad_id})_combined.png', facecolor='white', transparent=False)\n",
    "    plt.savefig(f'../../images/{norad_id}_maneuver_combined.png', facecolor='white', transparent=False)\n",
    "    plt.close()\n",
    "print(\"Done\")"
   ]
  }
 ],
 "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.8.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}