Example of how to do soft and hard iron calibration from a 3DOF magnetometer. http://psas.pdx.edu/
You can not select more than 25 topics Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.

980 lines
2.1 MiB

4 years ago
4 years ago
4 years ago
  1. {
  2. "cells": [
  3. {
  4. "cell_type": "raw",
  5. "metadata": {},
  6. "source": [
  7. "---\n",
  8. "title: PSAS Magnetometer Calibration\n",
  9. "layout: post.liquid\n",
  10. "published_date: 2019-11-05 00:00:00 +0000\n",
  11. "slug: magnetometer_calibration\n",
  12. "is_draft: false\n",
  13. "description: |\n",
  14. " Solving for the soft and hard iron IMU magnetometer calbiration on a\n",
  15. " amateur student rocket launched by Portland State Aerospace Society.\n",
  16. "data:\n",
  17. " mathjax: True\n",
  18. "---\n",
  19. "\n"
  20. ]
  21. },
  22. {
  23. "cell_type": "markdown",
  24. "metadata": {},
  25. "source": [
  26. "# PSAS Magnetometer Calibration\n",
  27. "\n",
  28. "_Published in November 2019_\n",
  29. "\n",
  30. "For a number of years I was involved with a university rocketry club called PSAS{% capture sidenote %}[Portland State Aerospace Society](http://psas.pdx.edu), a student aerospace engineering project at Portland State University. They build ultra-low-cost, open source rockets that feature very sophisticated amateur rocket avionics systems.{% endcapture %}{% include \"sidenote.liquid\" %}. One of the things I really liked to do was play with the data from the launches and learn how rockets and flight electronics work.\n",
  31. "\n",
  32. "\n",
  33. "Our rockets carry an instrument on them called an **IMU** (Inertial Measument Unit). An IMU typically measures both acceleration and rotation-rate of an object in all directions so with some clever math you can recreate the exact position, velocity, and orientation of the rocket over time. This is the only way to know where something is in space, and very important for rockets. IMUs have a problem though: they're not very precise.\n",
  34. "\n",
  35. "\n",
  36. "Since our IMU is fixed to the rocket, {% capture marginnote %}![diagram of the rocket on it's side showing the layout of the internal components](img/L-12_overview.png) Overview of the rocket \"LV2.3\". The IMU is near the primary flight computer.{% endcapture %}{% include \"marginnote.liquid\" %} which direction is \"up\" or \"left\", etc. relative to the Earth changes constantly as the rocket flies about. In order for the data to be useful we need to know which way we are pointed, which is why IMUs always have some kind of gryoscope to account for rotation. Our particular IMU has rate-gyroscopes that can sense rotation rate, and so we integrate that once to get orientation. Since any integration will give an estimate that drifts from the true value over time, our IMU also includes a 3-axis _magnetometer_ as well.\n",
  37. "\n",
  38. "## 9DOF IMU\n",
  39. "\n",
  40. "This makes what is often what is refered to as a \"9DOF\" IMU, because it has \"nine degrees of freedom\". That would be _x, y, z_ accleration, _x, y, z_ rotation-rate, and _x, y, z_ magnetic field. The reason to have a magnetometer is so you can use Earth's own magnetic field as a kind of guide to the orientation of the rocket. This doesn't instantly solve all problems in life, sadly. But it provides a good reference for the rough orientation of the rocket that can be used to produce a real-time estimate of rate-gyroscape drift, or 'bias', as we fly."
  41. ]
  42. },
  43. {
  44. "cell_type": "markdown",
  45. "metadata": {},
  46. "source": [
  47. "The magnetic field sensor in the rocket is sensitive, but because the Earth's field is so weak it's easily overwhelmed by local effects (metal screws, magnetic fields from nearby wires, etc.). In order to get good orientation data we need to undo{% capture marginnote %}![photo of two men awkwardly holding a large rocket body and an angle](img/L-12_ground_calibration.jpg) Members of the PSAS ground crew lifting and aranging the rocket around as many different orientations as possible before the flight.{% endcapture %}{% include \"marginnote.liquid\" %} these local effects.\n",
  48. "\n",
  49. "\n",
  50. "So a little before the flight we took the nearly complete rocket, powered the electronics up, picked it up and tried to move it around in every direction.\n",
  51. "\n",
  52. "## Magnetometer Calibration\n",
  53. "\n",
  54. "What do we expect good magnetometer data to look like? The Earth's magnetic field shouldn't change much, so it should look like a single vector going through the IMU. If we rotate the rocket one way or another, the angle that the vector goes through will change, but it should stay the same strength. That means that the magnitude of the local magnetic field should be constant, and it should measure it to be exactly the same as Earth's magnetic field.\n",
  55. "\n"
  56. ]
  57. },
  58. {
  59. "cell_type": "markdown",
  60. "metadata": {},
  61. "source": [
  62. "## Earth's Field Strength\n",
  63. "\n",
  64. "\n",
  65. "But what is the strength of Earth's magnetic field? It varies over time and over the surface of the Earth. We know where we launched from{% capture sidenote %}\n",
  66. "Latitude: `43.79613280°` N<br>\n",
  67. "Longitude: `120.65175340°` W<br>\n",
  68. "Elevation: `1390.0` m Mean Sea Level<br>\n",
  69. "{% endcapture %}{% include \"sidenote.liquid\" %} and the date, so we can look up{% capture sidenote %}[NOAA's magnetic field calculator](https://www.ngdc.noaa.gov/geomag/magfield.shtml)\n",
  70. "Model Used: `WMM2015`\n",
  71. " {% endcapture %}{% include \"sidenote.liquid\" %} what the expected magnetic field should be:\n",
  72. "\n",
  73. "Its direction:\n",
  74. "\n",
  75. "| Declination (+E/-W) | Inclination (+D/-U) | Horizontal Intensity |\n",
  76. "| ------------------: | ------------------: | -------------------: | \n",
  77. "| 14.7990° ±0.36° | 66.5386° ±0.22 | 20,754.1 nT ±133 nT| \n",
  78. "\n",
  79. "\n",
  80. "\n",
  81. "And as a vector:\n",
  82. "\n",
  83. "| North Comp (+N/-S) | East Comp (+E/-W) | Vertical Comp (+D/-U) |\n",
  84. "| -------------------: | ----------------: | --------------------: |\n",
  85. "| 20,065.7 nT ±138 nT | 5,301.2 nT ±89 nT | 47,819.4 nT ±165 nT |\n",
  86. "\n",
  87. "\n",
  88. "And finally, the total strength:\n",
  89. "\n",
  90. "| Total Field |\n",
  91. "| ------------------: |\n",
  92. "| 52,129.0 nT ±152 nT | \n"
  93. ]
  94. },
  95. {
  96. "cell_type": "markdown",
  97. "metadata": {},
  98. "source": [
  99. "## Calibration Data\n",
  100. "\n",
  101. "That's what we _expect_ to see. What do we actually get?"
  102. ]
  103. },
  104. {
  105. "cell_type": "code",
  106. "execution_count": 1,
  107. "metadata": {},
  108. "outputs": [],
  109. "source": [
  110. "import h5py\n",
  111. "import numpy as np\n",
  112. "from numpy import mgrid, pi, sin, cos\n",
  113. "from scipy import linalg\n",
  114. "import matplotlib.pyplot as plt\n",
  115. "from matplotlib import gridspec\n",
  116. "from matplotlib.ticker import FuncFormatter\n",
  117. "import matplotlib.patches as patches\n",
  118. "from mpl_toolkits.mplot3d import Axes3D\n",
  119. "from IPython.core.display import Markdown \n",
  120. "%matplotlib inline\n",
  121. "\n",
  122. "def minsec(x, pos):\n",
  123. " m = int(x/60)\n",
  124. " s = int(x - (m*60))\n",
  125. " return f\"{m:02d}:{s:02d}\"\n",
  126. "\n",
  127. "true_magnitude = 52.129\n",
  128. " \n",
  129. "# show a sphere\n",
  130. "u, v = mgrid[0:2*pi:40j, 0:pi:20j]\n",
  131. "radius = true_magnitude\n",
  132. "x=cos(u)*sin(v)*radius\n",
  133. "y=sin(u)*sin(v)*radius\n",
  134. "z=cos(v)*radius"
  135. ]
  136. },
  137. {
  138. "cell_type": "code",
  139. "execution_count": 2,
  140. "metadata": {},
  141. "outputs": [
  142. {
  143. "data": {
  144. "text/markdown": [
  145. "In the 22.1 minutes that we had the flight computer collecting data in our calibration run we recoreded 1,079,342 data points from the IMU."
  146. ],
  147. "text/plain": [
  148. "<IPython.core.display.Markdown object>"
  149. ]
  150. },
  151. "execution_count": 2,
  152. "metadata": {},
  153. "output_type": "execute_result"
  154. }
  155. ],
  156. "source": [
  157. "cal_data = h5py.File(\"data/L-12_calibration.hdf5\", \"r\")\n",
  158. "adis = cal_data[\"ADIS\"]\n",
  159. "\n",
  160. "time = np.array(adis[\"time\"])\n",
  161. "\n",
  162. "m_raw_x = np.array(adis[\"mag_x\"]) * 1e6\n",
  163. "m_raw_y = np.array(adis[\"mag_y\"]) * 1e6\n",
  164. "m_raw_z = np.array(adis[\"mag_z\"]) * 1e6\n",
  165. "\n",
  166. "time_elapsed = time[-1] - time[0]\n",
  167. "Markdown(f\"In the {time_elapsed/60:0.1f} minutes that we had the flight computer collecting data in our calibration run we recoreded {len(time):,} data points from the IMU.\")"
  168. ]
  169. },
  170. {
  171. "cell_type": "code",
  172. "execution_count": 3,
  173. "metadata": {},
  174. "outputs": [],
  175. "source": [
  176. "raw_mag = []\n",
  177. "for i, t in enumerate(time):\n",
  178. " raw_mag.append(np.sqrt((m_raw_x[i]*m_raw_x[i]) + (m_raw_y[i]*m_raw_y[i]) + (m_raw_z[i]*m_raw_z[i])))"
  179. ]
  180. },
  181. {
  182. "cell_type": "markdown",
  183. "metadata": {},
  184. "source": [
  185. "Looking over time at the _x, y, z_ values of the magnetometer and the mangitue compared to the NOAA predicted field we see it vary a lot."
  186. ]
  187. },
  188. {
  189. "cell_type": "code",
  190. "execution_count": 4,
  191. "metadata": {
  192. "class": "fullwidth"
  193. },
  194. "outputs": [
  195. {
  196. "data": {
  197. "image/png": "iVBORw0KGgoAAAANSUhEUgAABYAAAAH3CAYAAAACOftEAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8rg+JYAAAACXBIWXMAAAsTAAALEwEAmpwYAAEAAElEQVR4nOzdd3xb9bn48c/3HE3L247tJI6TkEkSCCGEFTZlJ5TelrZ07/aW0gGdtL+W7nG7C9370tve292wKWHPQAghhJAdxyt2vIf2+f7+0JEiy7ItT0nO8+aVF9bR0dFjWTo65znP9/kqrTVCCCGEEEIIIYQQQgghZh4j2wEIIYQQQgghhBBCCCGEmBqSABZCCCGEEEIIIYQQQogZShLAQgghhBBCCCGEEEIIMUNJAlgIIYQQQgghhBBCCCFmKEkACyGEEEIIIYQQQgghxAwlCWAhhBBCCCGEEEIIIYSYoSQBLIQQQgghhABAKfUOpZRWSi1IWvZbpdTBlPUOKqUem+740lFKLbBjfke2YxFCCCGEyEWSABZCCCGEAJRSF9hJpPckLYsnlrRS6sfDPO61Seu8I2n5LfayxcM87iupibZh1ht3DPlEKVVuv2YXZDuWTOVCzEqp/1BK3aWUalVKhZVSbUqpu5VSb1VKmdmKayoopW7Mxfd3ymc0/q9HKfWYUuqN2Y5PCCGEEEISwEIIIYQQowsAb1BKudLc9zb7/uMhhqlUDnwBuCDLcYxF1mJWSjmVUn8C/gpUAT8A3gd8C3AAvwc+NUlP915g2SRtayJuBN6RZvkhwAv897RGM9Q/gbcS+zx+DagB/qiUemdWoxJCCCHEcU8SwEIIIYQQo9tELNl3ZfJCpVQFcAXwr+MkBjFNlFK+UVb5KvAG4Ata69O01l/VWv9Ga/1fWutLgLOB5smIRWsd1loHJ2NbcRn8fhnTMQGtdXSytjlOO7TWt2ut/1tr/Q1gPdAHfDzLcQkhhBDiOCcJYCGEEEKI0e0GniFW3ZfsOsAC/pyLMSil5iulfqSUekkp1Wf/e1QpdUW6J1BKvU0ptUspFbD//454K4uU9R5SSjXYQ983KaV6lVIdSqmfKqU8abZ7slLqb0qpdnvbLyql3p10/wXAHvvmF5KG0f82aZ25di/aI0qpoFJqp1LqY0opNUxsJyil7rBjO6KU+rKKqVBK/bcdb49S6jdKKe8UxexTSn1NKbVPKRVSSjXbr1H5MDEvs9s59AB3pvsb2evXAB8BHtFafyndOlrrJ7XWv0l6zI1KqYftVhFBpdRe+zVJV1Ge+nxDegAn3XeGUupxpdSAUqrR3qYj099PKXWSUupXSqk99ja6lFL3KqXOSNmGBuYC5ye91gft+9L2AB7Heyaj93OmtNZHgF1Aog3McLEmxfFQ6u+tlLpdKXWRUmqL/V48pJS6YbxxCSGEEOL44xh9FSGEEEIIQWx4+beVUmVa60572duIVeZ25WgM64BXAX8HDgKlwFuAO5VSl2itH4ivqJR6G/A7YBvwGaAI+CbQOEwsXuDfwEPAJ4AzgfcDbcD/S9ru2cD9wH5i7Ql6gY3AL5VSVVrrrwMvE6uS/LYd69/sh++zt1EBPEFsSP1t9rY2AN8FFgEfShPb/cB9wCeB1wCfI1aN+Qb7+T4LnE+spcAR4NOTHLMbeABYCfzSXn8pcD1wllLqDK11ctuOAvv1vM/erpXmNY+7CnABvxlhnVQfB+4C/kGsXcg5wM3AfGLvofGYA9wN/I/97wpir3Ml8J8p6w73+10GnGQ/vgGoBt4NPKSUWqu13mmv91bgh8T+Vl+1l/UNF9g43zOjvp/Hwk6E1wLt43l8ktXAn4BfEPubvwH4oVJqZ/JnWAghhBBiWFpr+Sf/5J/8k3/yT/7Jv+P+H7E+rhp4T9KyBfayrxBLaoWAD9j3LbPv20gsyaqBdyQ99hZ72eJhnu8r9v0LRolrIjEUpNmeh1gy8t6kZU6gBdib/Bh7++HYIeOgbTxkP9dHU5b/A2hNuq2Al4CnAWfKun8BBoAy+/Zie5u3pIn5W/Z9r03Z9t/s5Selie0jKb9fE7Gk43dStr0F6JiCmD9p/61OS1m+wX7MB9LEfHOG79Xv2OuvGcP7O9174QtAFJibtOwdqe9L4LfAwZTHHrTX+880r5EGVmTy+w0TVyXQCvwsZXkD8NAIn5Hk9/543jMjvp8z+Ix+x459FnAKsaS2Tn7PpYs1JY6HUpZp+290atIyN7FE+P9l+veXf/JP/sk/+Sf/5N/x/U9aQAghhBBCZEBrfZRYtWO8BcPbgKPAPbkag9Z6IP6zUspjV0UWEks0rUtadS2xystfJD9Ga/3KcNsmlkz9Wcqyh4FZSqki+/bJwArgD0CJUqoy/o9YCwAvsV61o7ka2Ku1/mtSbBr4L/vmxpFi01qHibXPUMBPU9Z9HChTSpVNcszXAc8CB1O28RTQD1yc5jE/zmC7AMX2/3syXD/xXlBKmUqpUjuWB4m1hFub6XZS9AK/Sln2Xfv/G9KsP+T3S3mPFtjvUYj9vdalrj8GE3rP2FLfz6O5kVjFcCvwPPB6YtXHnxlb6EM8rbXeGr+hY/2YnyJWySyEEEIIMSppASGEEEIIkbn/Bv6slFpMrJXCn7TW4ZSWomOlR19lfDHY/V0/RyxRPH+E511g/38PQ6VbBrHKSH/KsnhbinJiycFl9u0f2P/SqRpmebIFxFoypIq3B1iYJrZAyrIu+//1wywvJxb/ZMW8jFiyuC3DbXRorbvSrZhGPPGbaWISpdTlwOeB04hVRCcrzXQ7KQ5qrUMpy16x/5/6N0n7+ymliolVt19LrF1DsgPjjAvG954Z7f08mj8Qq5Z2EkuqfwYoIVZFPxGH0izrJHaxQgghhBBiVJIAFkIIIYTIXLzX7s+AOmLJ2OHEE5BDJhizFdj/T006TWYM3wc+APwEeAzoIDac/J3Am8b4vKmiI9ynUv7/Rfv509k5zPKJGCm24e6b7JgV8CSxpGs6XSm3x/I+iD//amI9m0cORKmziFUvPw3cQKyVQpDYpGq/ZXomhh7u9/sTcBHwPWAr0E2sGvczTG+Faybv59Ec1Fr/2/75bqXUAeB2YtW6t9nLR7rgYw4Tx2jvWSGEEEKIEUkCWAghhBAiQ1rroFLq/4D3Aa9orZ8ZYfWD9v9PBF5Mc/+JxFoBjGmCqDHGcB3we6319ckLlVLvHibWJWm2sXQs8aXYa//fn5QYG85IibEDwPI0y09Mun+yTFbMe4n1Ch5tG+NxF7H+wu8gNnHfaF5vr39xcpWrUuqyCcaxQCnlSqkCjldQj/o3UUqVEps47ota61tS7vtymoeMpVp+Ot8zaWmt/6CU+ihwi1Lqd1rrPo5VFZelechCjr3/hBBCCCEmjfQAFkIIIYQYm+8Tqw792Cjr/ZtYFfAH7VYMCUqpk4hVPd6ltR6p8nCiMURJqRJUSi0DrklZ7zlifUvfo5QqSFl3IknC54m1BPioUmpW6p1KqeQ2CH32/9MlxjYBi5VSr0l6rAI+bt/81wRiTDVZMf8RWK6UemvqHXYf3vLxBqi1bgZuBS5QSn023TpKqTOVUu+wb0aJJU+NpPtN4BPjjcFWBKReTLjR/v8dGTw+Hlfqe/QC4Iw06/eR/rVOZzrfMyP5OrGJ4f4TQGvdQ+yzdlHySkqp1xOryBZCCCGEmHRSASyEEEIIMQZa65eBWzJY76hS6mZik2I9a1ftthGrSnwfsWTWuCaHyjQG4B/AO5VS/cQSmycQS0S9DKxJ2l5YKfVp4NfA40qp3xObLO5DwHbg1HHGaSml3kmsF+tOpdSvgH3ALOAU4NWA2173iFKqHnijUmo3scroA1rrp4FvAm8A/qiUug3YD1xFrHr0Nq31jvHEN8Uxf9eO8XdKqSuBJ4glOhcBryXWm/m3Ewj108T6On9FKfVqYn/rZqCC2ARzlwHx5PA/iSVmH1BK/TextiRvYOLFIAeAryqlVhJrS3EFscnffqG1HrVNhta6Vym1GfikUspLrN/0ScRalLzE0B7HW4C3KKW+AOwG+rTWm4bZ/LS9Z0bxD2AXcJNS6kd
  198. "text/plain": [
  199. "<Figure size 1728x576 with 1 Axes>"
  200. ]
  201. },
  202. "metadata": {
  203. "needs_background": "light"
  204. },
  205. "output_type": "display_data"
  206. }
  207. ],
  208. "source": [
  209. "fig, ax1 = plt.subplots(figsize=(24,8))\n",
  210. "plt.title(r\"IMU Magnetometer Calibration Run\")\n",
  211. "plt.ylabel(r\"Magnetic Field [$ \\mu$T]\")\n",
  212. "plt.xlabel(r\"Time [mm:ss]\")\n",
  213. "\n",
  214. "plt.plot(time, raw_mag, 'k-', alpha=0.3, label=\"Magnitude\")\n",
  215. "plt.plot([-100,5000], [52.129, 52.129], 'k-.', alpha=0.3, label=\"True Magnitude (NOAA estimate)\")\n",
  216. "plt.plot(time, m_raw_x, lw=0.5, label=\"X ('Up')\")\n",
  217. "plt.plot(time, m_raw_y, lw=0.5, label=\"Y\")\n",
  218. "plt.plot(time, m_raw_z, lw=0.5, label=\"Z\")\n",
  219. "\n",
  220. "plt.xlim([0,1350])\n",
  221. "plt.ylim([-100,100])\n",
  222. "ax1.xaxis.set_major_formatter(FuncFormatter(minsec))\n",
  223. "ax1.legend(loc=4)\n",
  224. "plt.show()"
  225. ]
  226. },
  227. {
  228. "cell_type": "markdown",
  229. "metadata": {},
  230. "source": [
  231. "This is because we have a couple of problems. One is that the effective _center_ of our magnetometer values are pushed off to one side. And the other is that the values are skewed (or \"stretched\") off to one side as well. This is somewhat easier to see in 3D:"
  232. ]
  233. },
  234. {
  235. "cell_type": "code",
  236. "execution_count": 5,
  237. "metadata": {},
  238. "outputs": [
  239. {
  240. "data": {
  241. "image/png": "iVBORw0KGgoAAAANSUhEUgAABYUAAAVhCAYAAADbYOFNAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8rg+JYAAAACXBIWXMAAAsTAAALEwEAmpwYAAEAAElEQVR4nOzddXxk1fnH8c8ZyUwm7rrJStZxd3cWKQ7FWrRIabEWihUoxf1XpNBCcShWCqW4Fy26u6xns3G3Scbv74+ZTJM1NkuyE/m+X6+8krn3zr3PvbFzn3vOc4xlWYiIiIiIiIiIiIjI+GBLdAAiIiIiIiIiIiIisuEoKSwiIiIiIiIiIiIyjigpLCIiIiIiIiIiIjKOKCksIiIiIiIiIiIiMo4oKSwiIiIiIiIiIiIyjigpLCIiIiIiIiIiIjKOKCksIqOGMeYdY8w7Ky2zjDFX9Xu9W2zZqRs6vtUxxlxljLESHYeIiIiIyNqorS0iMr4oKSwiw8IYk2mMudIY86UxptMY4zPGLDbG/NkYs3mi4xtKxpjJsQbpZomOZWV9DeV+H73GmFpjzGvGmPONMZk/cv8/N8b8amiiHbDfScYYrzHmzTWsf8QYEzLGbGWMOTl2buevYdvDVr6hERERERnN1NYeGUZjW9sYs7sxJmKMuXMN67eOtbMfMsZMXOn81vTxzlDGKCIbhiPRAYjI2GOMmQ38CygCngEeBHzAVOBI4BRjTJllWdVDcLhkIDQE+/kxJgNXApXAVyutuxa4fgPHszoXAg2AEygAdgJuAi4yxhxuWdZH67nfnwOlwO1DEWQfy7KWGWMuA241xpxqWdYDfeuMMQcAxwM3W5b1OfC5MeZY4BpjzPOWZS3rt20GcDfwHXDdUMYoIiIikghqaw+gtvYgWZb1tjHmz8DZxpgnLMv6T986Y4wDeABoAn4NBIAT1rK704GdgfU9PxFJICWFRWRIGWNSgX8AKcC2lmX9d6X1vwMuAsxQHM+yLN9Q7KePMcYAbsuyeodif5ZlhUh8QxrgRcuyFvdfYIzZmugNxT+MMbMsy2pMTGhrdAdwDHCTMeZly7LqjDFpwL3AYuCKftueTjTx+2dgr37LbwLygUMsywpsmLBFREREhofa2gOprb3eLgIOAB40xmzWr518MbAJcLhlWW2xZY+ubgfGmD2BHYEPGNguF5FRQuUjRGSonU70af6FKzdSIdpwsyzrj5ZlrQAwxpQbY+4yxsw1xnTHPt43xuy/LgdbS1kAY4y5yBhTGRtO94UxZq+VNugbDnVtrATBXMAPHB1bf7Ix5t+xIWABY0xVLNb0fvs4GXg99vKv/YZQXRVbv9o6Z8aY/Y0xH8ZKJHTGjrPtWuI7yhjzXexcFhpjjlqX67M2lmV9RrQHQA5wTr/jrtP3xBhTSbQhWN5/+Fi/9ecbY941xjQaY/wmOqTxGmNM0jrGFwFOJXrT86fY4huJ9pY4rf/NhGVZy4HfAHsaY34eO/6usfffFjtXERERkdFObW21tfvWr3db27KsTuAMYCZwWWx/04DLgb9blvXc2t5vjCkBngCagaNjyXkRGWXUU1hEhtpPiDb2Hl/H7bcm2rPzeaJDwjKJlgZ42Rizt2VZq60puw7OBjKI9io1RBs9rxhj9rAs64OVtj0EyCWaeGwCFsSWnwMsBG4DOoAtiDbENwF2jW3zHtEha78F7gfejy3/Zk2BxRqZT8aOcxWQBJwJvGuM2dOyrA9Xesv+wM9i59JONNH5hDHmK8uyFv7wpVirp2Jx78v/nvCv6/fkV8ANQDbRBu/KLgReAV4gOqRxJ+BSoBw4cV2CsyzrW2PM9cDlxpjbiH4f/2xZ1jur2fweoj2LbzHGvB07ryWo54KIiIiMHWprR6mt/SPb2pZlvWKMeRT4rTHm78CdQA/9EtirY4xxAk8TTXbvY1lW7Q8dS0RGKMuy9KEPfehjyD6AFuDrQWzvWc0yNzAf+PdKy98B3llpmQVc1e/1brFlHUBhv+VFQBfwSb9lE2Pb+oFJ6xjbSbH3bN9v2V6xZSevZvuron9q468dQC2wAsjst7w0Ft/nq4mvGyjpt7wgFvON63B9r4rto2It23wNtK7n9+QDoHIQ39srgXD/81mHc3AB82LnUQ1krGXbqUQbs81ABNg1Eb8H+tCHPvShD33oQx/D8aG29irbq609cNmg2tpEE7v1sbazBZy4Du+5Lbbt5Yn6PdCHPvQxNB8qHyEiQy0d6FzXjS3L6un72hjjNsbkAKlEG6Vb/4g4nrQsq77fceqIDnHaxhiTv9K2/7L6TU62cmzGGJsxJsMYk0u0twI/IratiDaa77Msq73fsaqJ9vjY0hhTvNJ7/mFZVk2/bRuA74Ep6xnDyrqAtH77H5LvSb/rZzfRGbJzgbeJli7achDxBYneeAB8ZFlWx5o2tCxrEfB7og3cByzLencQxxEREREZ6dTWXju1tQfR1rYsqwU4j2jb+XXLsv62tu2NMUcQ7cH8KtFJ/kRkFFP5CBEZap30a/T8kFjNq8uIDnEqX2n1KvXBBmHBWpZNAvpP9LBkDbFtA/yB6FAs90qrM9czromxz/NXs25ev/j6D8Navppt24gOJRsKafS7uRiq74kxZj+iw+S2IjoTc3+Zg4jvHGA7osMEjzTG7GVZ1htr2f6T2OdPB3EMERERkdFAbe21mxj7rLb2uluntnOs5vBfiPbCPt6yrB/z8yMiI4CSwiIy1OYBWxtjXJZl+ddh+9uJ1vi6h+jwqFaiQ55+Bhw3XEGuZJXZj40xE4k+aa8kOgtvZWw7O9En4xtypEV4Dct/9KzSsUbpdKLD2vrczo/8nhhjtgdeJtrIPJdo2Qc/UAI8xDpeP2NMOdGbhTeJ1tCbC9xnjNnIGqJZq0VERERGEbW1h964bWuvK2OMB3iWaPL+qFgPYxEZ5ZQUFpGh9gLRp/3HEm2Q/JBjgb9ZlnV2/4XGmFN+ZBzT17JsleFrq3EI4AEOtCyrsl9cq9vvYJ6S9x175mrWzVxpmw3hGKI1e1/tt2ww35M1nftRQADYs3/y1hiz7yDju49oo/Z0y7K6jDG/AP4JXA1cNMh9iYiIiIx2L6C29tqorT087gU2As6zLOvjYTqGiGxgqiksIkPtfqJP+m8yxmy68kpjjMMY8xtjTGlsUZiVnsLHGoOH/sg4jjHGFPbbZxHRBthnlmU1rvltcX09BlbuIfDb1WzbHfuctQ77/YLocLXTjTHp/eIrBn4KfGFtoBl8jTFbAbcSnbDk7n6rBvM96QYyjTErX6cw0UZs/P+MMcbOIBK5xpiTiM3UbFnWUgDLsl4mOpv0r40xW6zrvkRERETGCLW1105t7SFmjDkDOAF4xrKsO4d6/yKSOOopLCJDKtab82DgX8Bnxpingf8APqACOAKYDDwae8sLwM+MMV7gy9i6XxCtA7b5jwhlGfAfY8x9sddnEh3udOE6vv/VWMyvxPYRAQ5i9bXF5gI9wC+MMd1EJ5P4zrKs71be0LKskDHm10QTmx8bY/4CJMXicxKd6GE4HGKMaSD6dz8f2Bk4gGi9t8Mty2rqt+0LrPv35DNgP+B2Y8wnQMSyrCeBF4HzgTeNMY8AycDRrHvZiHyijejPiQ6x6+88YB/gAWPM1pZlrWnIn4iIiMiYora22tpD0dZeV8aYjYE7iPZKftsYc/waNu22LOuFoTy2iAw/JYVFZMhZlvVtrAHxK6JDww4l2ghbQbQ27OH9Zvj9FdH6YYcRraP1PXAG0eFdP6ah+n9EJ1g4m+gMxHOBMyzLem9tb+p3DouNMXOA64jWtO0hWrfrOAZOnIFlWd3GmBOB38eO64x9vUpDNbb907EG7e+IlkEIE23MHzWMw7Fujn32E5044zuiPQn+allW20rb/op1/57cTPQG5ASi9cwM0dmo3zfGHEv0HG8mWivtGeDPrOG6rOQuorNrn7Jy0teyrMZYY/9h4ALgxnXYn4iIiMiYoLa22tpD0NZeV1sSLX8B8Ke1bLecaLJbREYRowkjRURERERERERERMYP1RQWERERERERERERGUdGbFLYGPNrY8xcY8x
  242. "text/plain": [
  243. "<Figure size 1728x1728 with 4 Axes>"
  244. ]
  245. },
  246. "metadata": {
  247. "needs_background": "light"
  248. },
  249. "output_type": "display_data"
  250. }
  251. ],
  252. "source": [
  253. "\n",
  254. "fig = plt.figure(figsize=(24,24))\n",
  255. "gs = gridspec.GridSpec(2, 2, width_ratios=[1, 1])\n",
  256. "\n",
  257. "ax1 = plt.subplot(gs[0])\n",
  258. "plt.title(r\"Calibration Data XY\")\n",
  259. "plt.xlabel(r\"Magnetic Field X [$ \\mu$T]\")\n",
  260. "plt.ylabel(r\"Magnetic Field Y [$ \\mu$T]\")\n",
  261. "ax1.plot(m_raw_x, m_raw_y, lw=0.8, label=\"\")\n",
  262. "ax1.add_patch(patches.Circle((0, 0), 53, edgecolor=\"#cccccc\", linewidth=1.0, linestyle='--', fill=False))\n",
  263. "plt.xlim([-80,80])\n",
  264. "plt.ylim([-80,80])\n",
  265. "\n",
  266. "ax2 = plt.subplot(gs[1])\n",
  267. "plt.title(r\"Calibration Data YZ\")\n",
  268. "plt.xlabel(r\"Magnetic Field Y [$ \\mu$T]\")\n",
  269. "plt.ylabel(r\"Magnetic Field Z [$ \\mu$T]\")\n",
  270. "ax2.plot(m_raw_y, m_raw_z, lw=0.8, label=\"\")\n",
  271. "ax2.add_patch(patches.Circle((0, 0), 53, edgecolor=\"#cccccc\", linewidth=1.0, linestyle='--', fill=False))\n",
  272. "plt.xlim([-80,80])\n",
  273. "plt.ylim([-80,80])\n",
  274. "\n",
  275. "ax3 = plt.subplot(gs[2])\n",
  276. "plt.title(r\"Calibration Data XZ\")\n",
  277. "plt.xlabel(r\"Magnetic Field X [$ \\mu$T]\")\n",
  278. "plt.ylabel(r\"Magnetic Field Z [$ \\mu$T]\")\n",
  279. "ax3.plot(m_raw_x, m_raw_z, lw=0.8, label=\"\")\n",
  280. "ax3.add_patch(patches.Circle((0, 0), 53, edgecolor=\"#cccccc\", linewidth=1.0, linestyle='--', fill=False))\n",
  281. "plt.xlim([-80,80])\n",
  282. "plt.ylim([-80,80])\n",
  283. "\n",
  284. "ax4 = plt.subplot(gs[3], projection='3d')\n",
  285. "plt.title(r\"Calibration Data XYZ\")\n",
  286. "plt.xlabel(r\"Field Strength Y [$\\mu$T]\")\n",
  287. "plt.ylabel(r\"Field Strength Z [$\\mu$T]\")\n",
  288. "ax4.set_zlabel('Field Strength X [$\\mu$T]')\n",
  289. "\n",
  290. "ax4.plot_wireframe(x, y, z, color=\"k\", alpha=0.1, lw=0.2)\n",
  291. "\n",
  292. "ax4.plot(m_raw_y, m_raw_z, m_raw_x, '-', lw=0.5)\n",
  293. "ax4.plot([0],[0],[0], 'g.')\n",
  294. "\n",
  295. "ax4.set_xlim(-60, 60)\n",
  296. "ax4.set_ylim(-60, 60)\n",
  297. "ax4.set_zlim(-60, 60)\n",
  298. "\n",
  299. "plt.show()"
  300. ]
  301. },
  302. {
  303. "cell_type": "markdown",
  304. "metadata": {},
  305. "source": [
  306. "## Correction\n",
  307. "\n",
  308. "The two parts of the correction are called \"Hard Iron\" (fixed center offset) and \"Soft Iron\" (streched sphere) corrections.\n",
  309. "\n",
  310. "### Hard Iron\n",
  311. "\n",
  312. "This is the simpler of the two, one can essentially find the midrange value of across the entire calibration dataset and subtract that offset to move the '0' point back to center."
  313. ]
  314. },
  315. {
  316. "cell_type": "markdown",
  317. "metadata": {},
  318. "source": [
  319. "### Soft Iron\n",
  320. "\n",
  321. "Finding the soft iron correction is a bit trickier because we want to fit an matching elongated ellipsoid to the data, and then once we have an approximation for that ellipsoid apply stretch to the data to undo the elongation and get it back to a sphere. Luckily an algorithm for this has been worked out. For a detailed breakdown see\n",
  322. "\n",
  323. "<https://teslabs.com/articles/magnetometer-calibration/>\n",
  324. "\n",
  325. "After doing fitting we end up with both a correction matrix and an offset vector. This is both the soft iron and hard iron correction.\n",
  326. "\n",
  327. "To invert the stretch we multiply a vector representing each magnetometer reading (a 'sample', $\\vec s$) by the correction matrix (after subtracting the center offset).\n",
  328. "\n",
  329. "$$\n",
  330. "\\vec s_\\textrm{corrected} = \\mathbf{A} \\cdot (\\vec s - \\vec b)\n",
  331. "$$\n",
  332. "\n",
  333. "Where\n",
  334. "\n",
  335. " - $\\vec s_\\textrm{corrected}$ is a fully corrected sample at time t\n",
  336. " - $\\mathbf{A}$ is our soft iron correction matrix\n",
  337. " - $\\vec s$ is a raw sample from the IMU at time t\n",
  338. " - $\\vec b$ our hard iron offset vector.\n",
  339. "\n",
  340. "When we solve for $\\mathbf{A}$ on the calibration data we get the matrix:"
  341. ]
  342. },
  343. {
  344. "cell_type": "code",
  345. "execution_count": 6,
  346. "metadata": {},
  347. "outputs": [],
  348. "source": [
  349. "def ellipsoid_fit(s):\n",
  350. " ''' Estimate ellipsoid parameters from a set of points.\n",
  351. "\n",
  352. " Parameters\n",
  353. " ----------\n",
  354. " s : array_like\n",
  355. " The samples (M,N) where M=3 (x,y,z) and N=number of samples.\n",
  356. "\n",
  357. " Returns\n",
  358. " -------\n",
  359. " M, n, d : array_like, array_like, float\n",
  360. " The ellipsoid parameters M, n, d.\n",
  361. "\n",
  362. " References\n",
  363. " ----------\n",
  364. " .. [1] Qingde Li; Griffiths, J.G., \"Least squares ellipsoid specific\n",
  365. " fitting,\" in Geometric Modeling and Processing, 2004.\n",
  366. " Proceedings, vol., no., pp.335-340, 2004\n",
  367. " '''\n",
  368. "\n",
  369. " # D (samples)\n",
  370. " D = np.array([s[0]**2., s[1]**2., s[2]**2.,\n",
  371. " 2.*s[1]*s[2], 2.*s[0]*s[2], 2.*s[0]*s[1],\n",
  372. " 2.*s[0], 2.*s[1], 2.*s[2], np.ones_like(s[0])])\n",
  373. "\n",
  374. " # S, S_11, S_12, S_21, S_22 (eq. 11)\n",
  375. " S = np.dot(D, D.T)\n",
  376. " S_11 = S[:6,:6]\n",
  377. " S_12 = S[:6,6:]\n",
  378. " S_21 = S[6:,:6]\n",
  379. " S_22 = S[6:,6:]\n",
  380. "\n",
  381. " # C (Eq. 8, k=4)\n",
  382. " C = np.array([[-1, 1, 1, 0, 0, 0],\n",
  383. " [ 1, -1, 1, 0, 0, 0],\n",
  384. " [ 1, 1, -1, 0, 0, 0],\n",
  385. " [ 0, 0, 0, -4, 0, 0],\n",
  386. " [ 0, 0, 0, 0, -4, 0],\n",
  387. " [ 0, 0, 0, 0, 0, -4]])\n",
  388. "\n",
  389. " # v_1 (eq. 15, solution)\n",
  390. " E = np.dot(linalg.inv(C), S_11 - np.dot(S_12, np.dot(linalg.inv(S_22), S_21)))\n",
  391. "\n",
  392. " E_w, E_v = linalg.eig(E)\n",
  393. "\n",
  394. " v_1 = E_v[:, np.argmax(E_w)]\n",
  395. " if v_1[0] < 0: v_1 = -v_1\n",
  396. "\n",
  397. " # v_2 (eq. 13, solution)\n",
  398. " v_2 = np.dot(np.dot(-linalg.inv(S_22), S_21), v_1)\n",
  399. "\n",
  400. " # quadric-form parameters\n",
  401. " M = np.array([[v_1[0], v_1[3], v_1[4]],\n",
  402. " [v_1[3], v_1[1], v_1[5]],\n",
  403. " [v_1[4], v_1[5], v_1[2]]])\n",
  404. " n = np.array([[v_2[0]],\n",
  405. " [v_2[1]],\n",
  406. " [v_2[2]]])\n",
  407. " d = v_2[3]\n",
  408. "\n",
  409. " return M, n, d"
  410. ]
  411. },
  412. {
  413. "cell_type": "code",
  414. "execution_count": 7,
  415. "metadata": {},
  416. "outputs": [],
  417. "source": [
  418. "# samples\n",
  419. "s = []\n",
  420. "for i, t in enumerate(time):\n",
  421. " s.append((m_raw_x[i], m_raw_y[i], m_raw_z[i]))\n",
  422. "\n",
  423. "M, n, d = ellipsoid_fit(np.array(s).T)\n",
  424. "\n",
  425. "\n",
  426. "\n",
  427. "F = true_magnitude\n",
  428. "\n",
  429. "M_1 = linalg.inv(M)\n",
  430. "b = np.dot(M_1, n)\n",
  431. "A_1 = np.real(F / np.sqrt(np.dot(n.T, np.dot(M_1, n)) - d) * linalg.sqrtm(M))\n",
  432. "\n",
  433. "#print(A_1)\n",
  434. "#print(b)"
  435. ]
  436. },
  437. {
  438. "cell_type": "code",
  439. "execution_count": 8,
  440. "metadata": {},
  441. "outputs": [
  442. {
  443. "data": {
  444. "text/markdown": [
  445. "\n",
  446. "$$\n",
  447. "\\textbf{A} = \\left[\\begin{array}{ccc}\n",
  448. " 0.870368 & -0.128543 & -0.283684 \\\\\\\\\n",
  449. " -0.128543 & 1.510386 & -0.046543 \\\\\\\\\n",
  450. " -0.283684 & -0.046543 & 1.440805\n",
  451. "\\end{array}\\right]\n",
  452. "$$\n"
  453. ],
  454. "text/plain": [
  455. "<IPython.core.display.Markdown object>"
  456. ]
  457. },
  458. "execution_count": 8,
  459. "metadata": {},
  460. "output_type": "execute_result"
  461. }
  462. ],
  463. "source": [
  464. "Markdown(f\"\"\"\n",
  465. "$$\n",
  466. "\\\\textbf{{A}} = \\\\left[\\\\begin{{array}}{{ccc}}\n",
  467. " {A_1[0][0]:0.6f} & {A_1[0][1]:0.6f} & {A_1[0][2]:0.6f} \\\\\\\\\\\\\\\\\n",
  468. " {A_1[1][0]:0.6f} & {A_1[1][1]:0.6f} & {A_1[1][2]:0.6f} \\\\\\\\\\\\\\\\\n",
  469. " {A_1[2][0]:0.6f} & {A_1[2][1]:0.6f} & {A_1[2][2]:0.6f}\n",
  470. "\\\\end{{array}}\\\\right]\n",
  471. "$$\n",
  472. "\"\"\")"
  473. ]
  474. },
  475. {
  476. "cell_type": "markdown",
  477. "metadata": {},
  478. "source": [
  479. "And a hard iron offset vector:"
  480. ]
  481. },
  482. {
  483. "cell_type": "code",
  484. "execution_count": 9,
  485. "metadata": {},
  486. "outputs": [
  487. {
  488. "data": {
  489. "text/markdown": [
  490. "\n",
  491. "$$\n",
  492. "\\vec b = \\left[ \\begin{array}{ccc}\n",
  493. "-12.019415 & -3.209783 & -1.939041\n",
  494. "\\end{array}\\right]\n",
  495. "$$\n"
  496. ],
  497. "text/plain": [
  498. "<IPython.core.display.Markdown object>"
  499. ]
  500. },
  501. "execution_count": 9,
  502. "metadata": {},
  503. "output_type": "execute_result"
  504. }
  505. ],
  506. "source": [
  507. "Markdown(f\"\"\"\n",
  508. "$$\n",
  509. "\\\\vec b = \\\\left[ \\\\begin{{array}}{{ccc}}\n",
  510. "{b[0][0]:0.6f} & {b[1][0]:0.6f} & {b[2][0]:0.6f}\n",
  511. "\\\\end{{array}}\\\\right]\n",
  512. "$$\n",
  513. "\"\"\")"
  514. ]
  515. },
  516. {
  517. "cell_type": "markdown",
  518. "metadata": {},
  519. "source": [
  520. "## Python\n",
  521. "\n",
  522. "If we want to apply this correction we can make a convienient function to call on all our samples:"
  523. ]
  524. },
  525. {
  526. "cell_type": "code",
  527. "execution_count": 10,
  528. "metadata": {},
  529. "outputs": [],
  530. "source": [
  531. "def apply_mag_correction(sample):\n",
  532. " return np.dot(A_1, sample + b)"
  533. ]
  534. },
  535. {
  536. "cell_type": "code",
  537. "execution_count": 11,
  538. "metadata": {},
  539. "outputs": [
  540. {
  541. "data": {
  542. "text/markdown": [
  543. "\n",
  544. "``` python\n",
  545. "import numpy as np\n",
  546. "\n",
  547. "MAG_CORRECTION_A = np.array((\n",
  548. " ( 0.870367858077, -0.128543320363, -0.283683583608),\n",
  549. " (-0.128543320363, 1.510386103995, -0.046543028701),\n",
  550. " (-0.283683583608, -0.046543028701, 1.440804950101),\n",
  551. "))\n",
  552. "\n",
  553. "MAG_CORRECTION_b = np.array((-12.019414737824, -3.209782771540, -1.939040882716))\n",
  554. "\n",
  555. "def apply_mag_correction(sample):\n",
  556. " \"\"\"Take a sample and correct for hard and soft iron effects\n",
  557. " based on calibration run before Launch 12.\n",
  558. " \n",
  559. " :param sample: Vector (np.array) (x,y,z) instantaintous magnetometer reading\n",
  560. " :returns: vector (x,y,z) corrected magnetometer reading\n",
  561. " \"\"\"\n",
  562. " \n",
  563. " return np.dot(MAG_CORRECTION_A, sample + MAG_CORRECTION_b)\n",
  564. "```\n"
  565. ],
  566. "text/plain": [
  567. "<IPython.core.display.Markdown object>"
  568. ]
  569. },
  570. "execution_count": 11,
  571. "metadata": {},
  572. "output_type": "execute_result"
  573. }
  574. ],
  575. "source": [
  576. "Markdown(f\"\"\"\n",
  577. "``` python\n",
  578. "import numpy as np\n",
  579. "\n",
  580. "MAG_CORRECTION_A = np.array((\n",
  581. " ({A_1[0][0]:15.12f}, {A_1[0][1]:15.12f}, {A_1[0][2]:15.12f}),\n",
  582. " ({A_1[1][0]:15.12f}, {A_1[1][1]:15.12f}, {A_1[1][2]:15.12f}),\n",
  583. " ({A_1[2][0]:15.12f}, {A_1[2][1]:15.12f}, {A_1[2][2]:15.12f}),\n",
  584. "))\n",
  585. "\n",
  586. "MAG_CORRECTION_b = np.array(({b[0][0]:15.12f}, {b[1][0]:15.12f}, {b[2][0]:15.12f}))\n",
  587. "\n",
  588. "def apply_mag_correction(sample):\n",
  589. " \\\"\"\"Take a sample and correct for hard and soft iron effects\n",
  590. " based on calibration run before Launch 12.\n",
  591. " \n",
  592. " :param sample: Vector (np.array) (x,y,z) instantaintous magnetometer reading\n",
  593. " :returns: vector (x,y,z) corrected magnetometer reading\n",
  594. " \\\"\"\"\n",
  595. " \n",
  596. " return np.dot(MAG_CORRECTION_A, sample + MAG_CORRECTION_b)\n",
  597. "```\n",
  598. "\"\"\")"
  599. ]
  600. },
  601. {
  602. "cell_type": "markdown",
  603. "metadata": {},
  604. "source": [
  605. "## Apply Calibration\n",
  606. "\n",
  607. "After we apply the calibration fix above, do we get a better result?"
  608. ]
  609. },
  610. {
  611. "cell_type": "code",
  612. "execution_count": 12,
  613. "metadata": {},
  614. "outputs": [],
  615. "source": [
  616. "mag = []\n",
  617. "m_x = []\n",
  618. "m_y = []\n",
  619. "m_z = []\n",
  620. "for i, t in enumerate(time):\n",
  621. " s = np.array((m_raw_x[i], m_raw_y[i], m_raw_z[i])).reshape(3, 1)\n",
  622. " s = apply_mag_correction(s)\n",
  623. " mx, my, mz = s[0][0], s[1][0], s[2][0]\n",
  624. " m_x.append(mx)\n",
  625. " m_y.append(my)\n",
  626. " m_z.append(mz)\n",
  627. " mag.append(np.sqrt((mx*mx) + (my*my) + (mz*mz)))"
  628. ]
  629. },
  630. {
  631. "cell_type": "code",
  632. "execution_count": 13,
  633. "metadata": {
  634. "class": "fullwidth"
  635. },
  636. "outputs": [
  637. {
  638. "data": {
  639. "image/png": "iVBORw0KGgoAAAANSUhEUgAABYAAAAH3CAYAAAACOftEAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8rg+JYAAAACXBIWXMAAAsTAAALEwEAmpwYAAEAAElEQVR4nOzdd3wc1bn/8c/Zpt4syZZlueGGbUwvoVfTAumdFJJASC75JSHJDaQTcnNJctM7N6SSQBqXJHRM75hubAzuVbZk9b5tzu+P2V2tpJW0kiXvav19v156aXd2yqPVzuzMM+c8x1hrEREREREREREREZHc48l0ACIiIiIiIiIiIiIyOZQAFhEREREREREREclRSgCLiIiIiIiIiIiI5CglgEVERERERERERERylBLAIiIiIiIiIiIiIjlKCWARERERERERERGRHKUEsIiIiIjIAWSMudQYY40x8zK0fWuMuTYT245tf14shkuTpg15T4wxv49N82UizsGMMduMMb/PdBwiIiIiY6UEsIiIiGQ9Y8wZsUTQZUnT4kkka4z5xTDLvT1pnkuTpl8bm7ZwmOX+K50E3f7EMJUYY6bF3rMzMh1LujIZcyxRaIf5Of9Ax7O/jDFHG2P+aIzZbowJGmM6jDFPG2OuMcaUZTq+iWSMeVsmk+MjSfG56jXGrDPGfNkYE8h0fCIiIpK9suJuuoiIiMh+6APebYz5jLU2NOi1D8Zezz8IYphM04Cvxx4/nME4xiLTMW8Avpli+svAKuAvQPCARjQOxpjPAN8H6oGbgY24n+U34L6/ZwHnTsCmbiI73pO3AZcA16Z4bQngHNBohkr+XE0D3gL8F7AY+FCGYhIREZEspwSwiIiITHW3A+8ELgT+GZ9ojKkELgBuA951EMQgB4gxpsha2z3KbPustX8a4fXoRMY0GYwxFwI/BO4C3mmt7Ul6+WfGmKuBy1IuPEbW2igT/J6k+X9Km7U208lpGPS5Msb8FHga+IAx5j+ttY2ZC01ERESylUpAiIiIyFS3AVgNfGDQ9Pfittb7ezbGYIyZa4z5aawLd1fs5zFjzAWpNmCM+aAx5jVjTF/s96XxUhaD5nvYGLMrVp7idmNMpzGmxRjzK2PMkFbIxpjDjTH/Z4xpjq37FWPMR5NePwO31SfA15O6n/8+aZ5ZsXqtDbESAa8aY64yxphhYjvEGHNHLLYGY8w3javSGHNTLN4OY8zvjDEFkxRzkTHmv40xm40xIWPMnth7NG2YmJcYY+4yxnQAd6b6H6XLDKp3a4xZEHsvViW/Z8aYiti21yX/74wxpxtj7jPGtMfKAKw2xrw5xXZmGGNujs3Xboz5uzFm5hhCvR5oB94/KPkLgLV2t7X2G0nbe5Mx5p/GmJ2xz8FeY8wfjDG1Y31PBqmJxR7/O/5sjJk+aPl4WZejY//HRqAr9to0Y8x3jDEvJr1nzxtj3j9oHQ/jtv6N10m2ZuD/aUgNYGNMXuzzuzn2N+8yxvzEDCqNkRTf4bHXG40xPcaYu40xc0d7f4ZjrbXAo4ABFiRtL2W9YjMBxwwRERGZetQCWERERHLBTcD3jDEV1trW2LQP4rbMbcvSGI4DzsFtHbwNKAfeD9xpjFlprX0gPqMx5oPAH4CXgC8CJcB3gN3DxFIA3I9b+uA/cbvrXwHsA76atN6TcMsRbAG+C3QCFwM3GmOmW2uvB9YDnwe+F4v1/2KLb46toxJ4EqgBfh5b10XAD3ATUp9MEdsq4D7gC8Bbga/gJuveHdvel4HTgUuBBuCaCY45D3gAWA7cGJt/MXAlcKIx5gRrbV9SzIWx9/O+2HrTKQPgM8ZUDZoWtta2D57RWrvZuKUWbgQ+g9vqFuBXQDVwUTweY8zbgb/itvr8BhDBvdHwT2PMJdbam5P+xvuBpbH1rAfOx23NOypjzALgcOAPSZ/n0XwENxH5C6AJOBS4HHiDMeaIQe/pWNwB7AS+BCwDPg4sN8Ycn6Lkyh+BvcB1uCUSAA4B3gP8A/hf3BIWbwNuMsYErLW/jc33LcAPnMTAmzn7RojtH7if97/jlso4Avczf7Ix5qQUrYZ/C7TE4qsBPgv8CTh1lPdgJPNjv5v3Yx1pHTNERERkirLW6kc/+tGPfvSjH/1k9Q9wBmCBy5KmzYtN+y+gCggBH4+9tiT22sW4SVYLXJq07LWxaQuH2d5/xV6fN0pc+xNDYYr15eMm6u5NmubHTWhtSl4mtv4wsUaASdMfjm3rM4Om/xNoTHpugHXAM4B/0Lz/AHqAitjzhbF1Xpsi5u/GXnv7oHX/X2z6ihSxfXrQ31ePm1T9/qB1Pwu0TELMX4j9r44dNP2i2DIfTxHzl8bwed0WW2bwz8Ox1y9N9fnCTVb3AYfhJiAtcHXyZwY3sfp/g5bzxt6T3YAnNu3K2PKfGDTvX4d7XwbNd3FsvqvG8Hen+kyfHlvPe1PsN8n7w5D3BPh9bNpfB63z/w3+2+jfp+8BzKD58wDvoGkG9ybA64Om/4lB+9Sg/+vvk55fGNvmTwfN9+nY9E+miO+2QfN+JjZ9WZqfq6dxjzVVwCLcmyMO8PxIsQ6OY9C0+Gd8xGOGfvSjH/3oRz/6mbo/KgEhIiIiU561tgm4m/5Wex/ETZTdk60x2KQu9caY/FhL2mLcZMxxSbMeA8wAfp28jLX29eHWjZsQumHQtEeAamNMSez54bitKf8MlBljquI/uCUOCnBbQo7mTcAma+2tSbFZ4H9iTy8eKTZrbRi3fIbBbama7AmgwhhTMcExvxd4Dtg2aB1PA93A2SmW+UUa6022Flg56OdzoyxzOdCKm6T9GW7X/v9Jen0lUAn8cVDcFbh/fy1ui19w3/dO4DeDtvH9NOMvjf3uSHP+xGfauEpjsa3DbQF/3EjLjuJHg57/L26L8cGfLYBfxT5/yXEFrVtjGGNMwLhlPipxW7wuNsaUplhPOt4U+/2dwTHgvm9vYqhfDnr+SOz3gsEzDuME3Fa5+3BLz1yP+3e8Jc3lh5POMUNERESmKJWAEBERkVxxE/B3Y8xC3FIKf7HWhs3AMrRjZUefZXwxGGMCuKUPPggMrgGavN15sd8bGSrVNHBb7fUOmhbvxj8NNzG4JPb8x7GfVKYPMz3ZPNySDIO9Gvs9f9D0Rju0FEBb7PeOYaZPw41/omJegpssHq5r/+B1tFhr21LNOIJ2a+39Y1nAWttkjPkP3NbT3cAHrbXJ5Sbif/9tI6xmOm7SdR6wzQ4tkfB6muHEE79pJ/+MMYuBb+MmqosHvVye7npSGBCztTZojNnG0M8WxMp8DIrL4LbK/ThuqY/BO2Q5Y0h0J5kHdFprd6WIb8sw8W0f9Dx5v0zHWuAq3LFc5uO2AJ6B2/p9f6RzzBAREZEpSglgERERyRXxWrs3AHNwk7HDiScghwwwFlMY+z04ITKRMfwINyH1S+Bx3LqgUeDDwPvGuN3BoiO8Zgb9/kZs+6m8Osz0/TFSbMO9NtExG+Ap4GvDvN426PlYPwf748LY70LcMhbJCcP43/8J3JIgqbw8QXHE38cj0pk51or2UdzSGt/AvTnRg3sz4y8cuMGnU/2v/hO3le6fcev87sOtnXwh/cnUA2W0z/hoBtxYMMbcj5sU/hXwzqT5hrt55R1jXGOJTURERLKUEsAiIiKSE2Kt7v4GfAy3rufqEWbfFvu9FHglxetLcVtgjmlQpTHG8F7gj9baK5MnGmM+Okysi1KsY/FY4hsknkDsTaOl6kgtobfiDvY12NKk1yfKRMW8CbdW8Jha6E42Y8ybgMtwWzefB/zeGLMiqfVx/O9vTSP2bcBJsUHOklsBLxlm/gGsOzDdWuCtxpjP2BSD1w1yJm5L1DOttQ8n/U0FuCUq9scS3IR9fJ15uK1vH0tz+fcCj1hr35880RiTqtTHWFr9bwXOM8bUJbcCjrXuPwS3LvOkiv2ffg58zhjzBmvt07GXWkn9vh8y2TGJiIhI9lENYBEREcklP8JtfXjVKPPdj9sK+D9iyZoEY8wK4Czgrnjd0EmKIcqglnXGmCUMreX5PNAIXGaMKRw073njiC/uRdyu9Z8
  640. "text/plain": [
  641. "<Figure size 1728x576 with 1 Axes>"
  642. ]
  643. },
  644. "metadata": {
  645. "needs_background": "light"
  646. },
  647. "output_type": "display_data"
  648. }
  649. ],
  650. "source": [
  651. "fig, ax1 = plt.subplots(figsize=(24,8))\n",
  652. "plt.title(r\"IMU Magnetometer Fixed Calibration Run\")\n",
  653. "plt.ylabel(r\"Magnetic Field [$ \\mu$T]\")\n",
  654. "plt.xlabel(r\"Time [mm:ss]\")\n",
  655. "\n",
  656. "plt.plot(time, mag, 'k-', alpha=0.3, label=\"Magnitude\")\n",
  657. "plt.plot([-100,5000], [52.129, 52.129], 'k-.', alpha=0.3, label=\"True Magnitude (NOAA estimate)\")\n",
  658. "plt.plot(time, m_x, lw=0.5, label=\"X ('Up')\")\n",
  659. "plt.plot(time, m_y, lw=0.5, label=\"Y\")\n",
  660. "plt.plot(time, m_z, lw=0.5, label=\"Z\")\n",
  661. "\n",
  662. "plt.xlim([0,1350])\n",
  663. "plt.ylim([-100,100])\n",
  664. "ax1.xaxis.set_major_formatter(FuncFormatter(minsec))\n",
  665. "ax1.legend(loc=4)\n",
  666. "plt.show()"
  667. ]
  668. },
  669. {
  670. "cell_type": "markdown",
  671. "metadata": {},
  672. "source": [
  673. "Yes! Quite a bit better. Notice how the magnitue of the vector now stays very close to constant and is very close to the NOAA estimate!\n",
  674. "\n",
  675. "Again in 3D we can see a much closer to spherical data:"
  676. ]
  677. },
  678. {
  679. "cell_type": "code",
  680. "execution_count": 14,
  681. "metadata": {
  682. "scrolled": true
  683. },
  684. "outputs": [
  685. {
  686. "data": {
  687. "image/png": "iVBORw0KGgoAAAANSUhEUgAABYUAAAVhCAYAAADbYOFNAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8rg+JYAAAACXBIWXMAAAsTAAALEwEAmpwYAAEAAElEQVR4nOzdd3ib1dnH8e/RsOS9t2NnJ2Swyy4b+rJXKLvMsilQZoFSdhkFShll75FSNoVSdhkFSqGsEJKQYcd7b0vWOO8fklU7cRI7sWMn/n2uS5elZ57nsZwc3Trnvo21FhEREREREREREREZGxwj3QARERERERERERERWX8UFBYREREREREREREZQxQUFhERERERERERERlDFBQWERERERERERERGUMUFBYREREREREREREZQxQUFhERERERERERERlDFBQWkT6MMe8bY94foXOfYIyxxpjxI3H+aBuuMsbYFZatdE+i7XxyvTZuFYwxu0bbs+tIt0VEREREVk19bfW1RURGCwWFRcaQXh2a/h7VI92+wTLGOI0xJxlj3jXGNBhjAsaYKmPMi8aYg0a6fUPJGOOIdqIPHum2rKif91XAGFNvjPnEGHOTMWbSOh5/9+i1pw1Rk3uO6zLGfG2MqTXGZPaz/uTo9ZxvjJlujPEZY15ZxbGKjTFtxpj3jDFmKNspIiIiGwb1tTdc6msPbV/bGJNtjKkzxnxpjHH1sz7FGLPcGPODMcZrjFm2mr+d2GOo2iciESv9cYrImPAQ8P4Ky7qiP/dev01ZO8aYFOBlYFci13ITUA/kA/sDLxljjrHWPj0EpxsN98QB/A54DHhphXUfAPFA93pu04p63lcOIB3YAjgTOM8Yc6619t61PO7uwOXAo0DzOrcyylobNMacDHwK3AEc27POGFMA3Ap8BtxhrQ0bY64GbjDGHGGt/csKh/sz4AROsdaqwyoiIjK2qa89OKPhnqivPYR9bWttnTHmV8DTwAVE3j+93QgUAj+11vqMMecBSas43E+BU4FPhqJtIvI/CgqLjE2fWmv7nY5lrR3pzs5APQjsApxorX10hXXXG2MOGKoTDcc9McYkWms7huJY1tow4BuKY62jld5XxphLgL8B9xhjfrTWvj0yTeuftfY/xpg/AhcYY5621r4eXXUPkc7/ydH7C3ALMAf4kzHmbWttA4Ax5mhgX+ACa+3i9XsFIiIiMgqprz0I6msP2AbV17bWPmOMOQq4yhjzgrV2EYAxZkfgdOAua+3H0W1f6u8Yxpg8In3weuDn66XhImOI0keISB9mhZxexpiLo9N1jlthu8Ojy6/stcwV3f776FT7emPM08aY4n7Os7cx5ovodsuMMRcBA5p2b4zZEjgceKKfTioA1tpXrbWvRrePi06J+iw69c0XbeN5A5nqv+I9WWHdfsaY/0aPucQYc24/2ywzxnxkjNnOGPOBMaYTuDu67qfRe7S01z173hgzrdf+44FA9OXxvaZQvR9d32+eMxNJefC8MabRGNMVnb7V5/e4Qvs2i15rpzGm2hhznTFmnf6fsNbWEOnAhYHe75UB/U6MMY8SGbkAsLTXte8aXX+gMeal6PQzf7Tdj0VH+g7UlcAS4F5jTLIx5gjgIOAGa+28XtcSBE4CMoDbo+fPBP5IZETxHwdxThERERmD1Nde8z1ZYZ362quxAfS1Twf8wIMmIg54ACgFfrO6HY0xTuAvQA5wjLW2fIDnFJEB0khhkbEpyRiTtcKyNmutv59t/0BkFORdxpgPrLWlxphC4D4i0+6vB4h2Lv4K7Edk6tEdRKYEnQ3sbIzZwlpbF912F+A1oBy4KnqeMxn4dKVDoj8fGeD2KdHjPwv0fLu+N5HAXga9OlCDtDlwAHBvtC2HA380xnittStOkSoics1PRB8t0eU/B3KBh4EqYDyR6VEfGWNmWmtrgTrgeCLT2T4E7o/uW7OqhhljJhOZYmWAu6LHOBJ43BiTa639wwq75ANvAHOjj32IdBCXEpmqttastUuNMf8EdjHGJFtr2xj47+Q+II1IkPZ8IqMEAOZHf54UvcZ7ouumA78EtjPGbGatXeOoDmttpzHmVODt6HF+BnwH3NDPtl8bY24ErjDGPA0cDaTSd0SxiIiIjG3qa6uvrb52pG2VxpgLiIw8Pw0oADYB9h7ASO7fAzsDV1tr31zDtiKyNqy1euihxxh5EMkJZlfxOCG6zfvA+yvsV0ykE/lPInlT3wLagEm9tvl59DgHrrDv5kAQuLHXss+JdNTyei3Liy6zwPg1XMfz0e3SB3jdTsDTz/JHotcR12vZVZF/Gvts19896blv+/Ra5ibSee8EMnotXxbd9uh+2pDQz7KpRL5R/02vZa7oMR5dze91117LniUyYmDrXsviiIxo7QKy+mnfwSsc9yvg34N4X52ymm3uiG6z6Vr8Tq5b1ftiFfdvl+j2Rw3y7+Oh6H4hYNvVbBcHzAMaottfuTZ/j3rooYceeuihx8b1QH3tFZerr923fQevcNyvGCN97eh7ujV63x8ewPYHR8/xJuAY7N+iHnroMbCH0keIjE23AXut8PjHqja21pYRHYVA5NvzPYHzbN/8qUcBFcC/jDFZPQ8iIxQWAXtALC/U1sBca211r3NUA08NsP0p0Z+tA9nYWhuy0ZEZ0Wl36dG2vUOkoMH0AZ53RQustX/vdZ4A8CciuWj3XGHbRiKjAlZsW2fPc2NMkomkI2gEFgA/WZtGRada7Qu8a639T69zdRNJceBl5YIeVXblXF7/BNapmnEvbdGfydG2DMnvpOf+RaejpUSPMY/IB6vB3r+ekRGVRDrpqzpnN/9LI/E9kVEMIiIiIj3U11ZfW33tvn5JJKjfDPx6dRsaYyYRGQ1fTiRthGbjiQwTpY8QGZvm20EWIbDWPmmMORw4EPibtXbFaU7TiExhq1vFIcqiP8dHfy7oZ5v+lvWnp4OazACnwRljjiVS+XY2kW/Oe0sb4HlXtLCfZT3XMGGF5cv669AYY/KJVN89gEgV4d7qV9x+gLKBRP437au371fRvtJ+tm0iEvgcCsnRn7EPF0PxOzHGTCVy//Zi5YrFAzpG9DhbRtvyDbAp8FvgilVtb639LJqO7b/RDygiIiIiPdTX7ittgOddkfraAzeq+9rW2mXGmBqg3FrbvJrzeYmMVE8E9rXRlCgiMjwUFBaRAYl+K7xN9OUUY0y8tbar9yZEcmKduopDDGXF3u+BQ4HNiHzDvlrGmJ8TyS32d+BOoBroBrYEbmL9FN3sWnFBtLDEm0Q6+LcT+da9nchUtD+up3b1CA3z8WdHz7EUhuZ3YoxJAT6I7nc1kVEynUSmms0dyDGix3ERSR3RCOxOJOfZxcaYudba7wZ8hSIiIiJrSX3tdaa+9ijtaw/S3UTedxdYa/81DMcXkV4UFBaRgbqfyDfZ5xOZEncLkWluPX4kkl/qn2sYObks+nNaP+v6W9afl4iM4jyRAXRUiUy3Wwrs33sEQXRq0rqY2s+ynmtYOoD9ZwOzgBPtCpWdjTEZ9B29YAfRrjqgg0gRhxX1LBtI+4aEMWYCkemQH1tr26OLB/M7WdW170akcMhu1tr3ex0jnpVHgqzORUTy8R1lrW0wxpwVPfaDxpgdNGVNRERE1gP1tVemvvYAbAB97QExxpxEJEXbC9ba24b6+CKyMuUUFpE1MsacSKQK8eXW2j8S+Wb9TGNM71xZzxCZtnRpP/ubngrM0XxmXwBHRnOe9WyTBxwzkPZYa78AXgB+YYw5bhVt3s8Yc0D0Zc8386bX+njgVwM532pMM8bs0+uY7ugxfcBApgyu1K7ocY4nUqE4xlobih53jR2w6LavA7tH0yL0bt+5RAo8rJcKvsaYHOAvRP6/ubbXqsH8Tno6tytee7/3D7iYgY8Snkqk+vLfrLVzIVIlGbgE2Ja+H8ZEREREhpz62qukvvYajPa+9kAZYzYnMkr4RyJfRojIeqCRwiKyWtFvnu8gUhW45xvb3xApnvCIMWa2tbaRSEf1EOAaY8x2RIoY+Ijk0zoIeJpItWGIdCTeBD4xxtwXXXYasITIiM2
  688. "text/plain": [
  689. "<Figure size 1728x1728 with 4 Axes>"
  690. ]
  691. },
  692. "metadata": {
  693. "needs_background": "light"
  694. },
  695. "output_type": "display_data"
  696. }
  697. ],
  698. "source": [
  699. "\n",
  700. "fig = plt.figure(figsize=(24,24))\n",
  701. "gs = gridspec.GridSpec(2, 2, width_ratios=[1, 1])\n",
  702. "\n",
  703. "ax1 = plt.subplot(gs[0])\n",
  704. "plt.title(r\"Fixed Calibration Data XY\")\n",
  705. "plt.xlabel(r\"Magnetic Field X [$ \\mu$T]\")\n",
  706. "plt.ylabel(r\"Magnetic Field Y [$ \\mu$T]\")\n",
  707. "ax1.plot(m_x, m_y, lw=0.8, label=\"\")\n",
  708. "ax1.add_patch(patches.Circle((0, 0), 53, edgecolor=\"#cccccc\", linewidth=1.0, linestyle='--', fill=False))\n",
  709. "plt.xlim([-80,80])\n",
  710. "plt.ylim([-80,80])\n",
  711. "\n",
  712. "ax2 = plt.subplot(gs[1])\n",
  713. "plt.title(r\"Fixed Calibration Data YZ\")\n",
  714. "plt.xlabel(r\"Magnetic Field Y [$ \\mu$T]\")\n",
  715. "plt.ylabel(r\"Magnetic Field Z [$ \\mu$T]\")\n",
  716. "ax2.plot(m_y, m_z, lw=0.8, label=\"\")\n",
  717. "ax2.add_patch(patches.Circle((0, 0), 53, edgecolor=\"#cccccc\", linewidth=1.0, linestyle='--', fill=False))\n",
  718. "plt.xlim([-80,80])\n",
  719. "plt.ylim([-80,80])\n",
  720. "\n",
  721. "ax3 = plt.subplot(gs[2])\n",
  722. "plt.title(r\"Fixed Calibration Data XZ\")\n",
  723. "plt.xlabel(r\"Magnetic Field X [$ \\mu$T]\")\n",
  724. "plt.ylabel(r\"Magnetic Field Z [$ \\mu$T]\")\n",
  725. "ax3.plot(m_x, m_z, lw=0.8, label=\"\")\n",
  726. "ax3.add_patch(patches.Circle((0, 0), 53, edgecolor=\"#cccccc\", linewidth=1.0, linestyle='--', fill=False))\n",
  727. "plt.xlim([-80,80])\n",
  728. "plt.ylim([-80,80])\n",
  729. "\n",
  730. "ax4 = plt.subplot(gs[3], projection='3d')\n",
  731. "plt.title(r\"Fixed Calibration Data XYZ\")\n",
  732. "plt.xlabel(r\"Field Strength Y [$\\mu$T]\")\n",
  733. "plt.ylabel(r\"Field Strength Z [$\\mu$T]\")\n",
  734. "ax4.set_zlabel('Field Strength X [$\\mu$T]')\n",
  735. "\n",
  736. "ax4.plot_wireframe(x, y, z, color=\"k\", alpha=0.1, lw=0.2)\n",
  737. "\n",
  738. "ax4.plot(m_y, m_z, m_x, '-', lw=0.5)\n",
  739. "ax4.plot([0],[0],[0], 'g.')\n",
  740. "\n",
  741. "ax4.set_xlim(-60, 60)\n",
  742. "ax4.set_ylim(-60, 60)\n",
  743. "ax4.set_zlim(-60, 60)\n",
  744. "\n",
  745. "plt.show()"
  746. ]
  747. },
  748. {
  749. "cell_type": "markdown",
  750. "metadata": {},
  751. "source": [
  752. "## Applying to Flight Data\n",
  753. "\n",
  754. "We can now take our calibration matrix and apply it to real flight data! Here is a 3D look at the [Launch-12](https://github.com/psas/Launch-12) **raw** (uncalibrated) magnetometer data from liftoff to apogee:"
  755. ]
  756. },
  757. {
  758. "cell_type": "code",
  759. "execution_count": 15,
  760. "metadata": {},
  761. "outputs": [],
  762. "source": [
  763. "l12_data = h5py.File(\"data/L-12_flight.hdf5\", \"r\")\n",
  764. "l12_adis = l12_data[\"ADIS\"]\n",
  765. "\n",
  766. "l12_time = np.array(l12_adis[\"time\"])\n",
  767. "\n",
  768. "l12_m_raw_x = np.array(l12_adis[\"mag_x\"]) * 1e6\n",
  769. "l12_m_raw_y = np.array(l12_adis[\"mag_y\"]) * 1e6\n",
  770. "l12_m_raw_z = np.array(l12_adis[\"mag_z\"]) * 1e6"
  771. ]
  772. },
  773. {
  774. "cell_type": "code",
  775. "execution_count": 16,
  776. "metadata": {},
  777. "outputs": [
  778. {
  779. "data": {
  780. "image/png": "iVBORw0KGgoAAAANSUhEUgAABYUAAAVhCAYAAADbYOFNAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8rg+JYAAAACXBIWXMAAAsTAAALEwEAmpwYAAEAAElEQVR4nOzdd5icVdnH8e+ZmZ2Zne29ZtMrMQHkBQVBQZDekaI06UgXkF4F6Yh0AqFIFQQFUZogTaSIgALpbZNs72V2+nn/2M2aQBKyyW6e2Z3f57pyZefZp9xzdnb2zP2ccx9jrUVEREREREREREREUoPL6QBEREREREREREREZPNRUlhEREREREREREQkhSgpLCIiIiIiIiIiIpJClBQWERERERERERERSSFKCouIiIiIiIiIiIikECWFRURERERERERERFKIksIiIutgjHnYGGONMR6nYxlMxpg3jTFvfmWbNcZcuZHn+0Hf8ScMRnwiIiIiMrKpn73B51M/W0SGjJLCIjJo1Gn5OmPMKcaY3xtjlvS1zbvr2K/SGHOpMeY9Y0yTMabdGPMvY8xJxhj3Bl7ryr5rrO3fU4P7zDaOMWZcX5xbbuD+3zHGxI0xD67le8YY85YxpsMYU7Xa8z9oHef6Rd/3j920ZyEiIiKyeamf/XXqZ69pIP1sY8yxfbH/Yh3fP2hVMnu11943/Xt4sJ+TiAytEXVXTkQkCV0I5AAfAYXr2e8A4FLgReD3QATYG7gP+AHwkwFc8zyg/ivbln7DMelAbADX2FjjgCvojefTb9rZWvu+MeZO4ExjzBPW2r+t9u1TgZ2A06211caYXwMHA3caY96w1rat2tEYMwa4GnjVWvvw4DwVEREREXGQ+tlr2uB+trX2YWPMEcCvjDF/tNYuWfU9Y0wOcCfwOfBrIA84aj2nuxiYCry3KcGLyOanpLCIyND6PlBtrbXGmKXr2e9NYLS1dvVO5j3GmN8BRxljbrLWfrKB13zeWrtwIEFaa0MD2X8zuwTYH5hljPmWtbbbGFMFXA+8C9wNYK2NGGOOA/4J3AIcv9o57gMscNJmjVxEREREhor62ZvmJHoTv/cDu662/SagGNjfWhuhNwn+2NpO0DcDbyrwe2vtrCGNVkQGncpHiMhmZ4zZzxjzJ2PMcmNM2BhTZ4x5xBhT/pX9Vk1V+sFazrF09SlKxpgxffteY4w51BjzuTEmZIyZb4w5dB1xHGeM+cgY022MaTPG/LPvjvlXFRhjHu3bp7Nvmlr+hjxXa+0ya63dgP0+/0pHdZVn+v7fYkOut7FWTQ/7yrYyY8xTfeUZ2o0xzxhjSte272rHnG6MWdT3c/3UGLPzat87Fnit7+FDq001W+u5VrHWdgEnA2OBa/o23wukASes3r7W2o+AW4HjjDE/7Lvu0cCPgAuttcs2rEVEREREhh/1s9e6n/rZa9HXL74A+GHfwAqMMd8HTgB+09evXt/zmkHv4Iy5fceIyDCjkcIi4oTjAENvJ6IJmAKcCHzHGDNzE++m7wn8jN6kYRu9HZQnjTGfWmvnr9rJ9JYkOI3ekaZXAj3AVsBewJNfOedfgMXARcBk4HR6p52tbxrVYFnVgW8cwDF5xpivTqFrsdYmNvQExhg/8Dd6n+89wBxgN+Cv6znsFCATmEVv+5wNPG+MGW2tbQXepnd074V9+7zTd9x/vikea+0rxphH6S0jYej9OV9srZ23lt2voHea4CxjzC70Jon7RxSLiIiIjGDqZ2849bN7r384cIsx5u99xy4CLv+G55ADPEvvTLxD+gZxiMgwo6SwiDjhJ9ba4OobjDEv0Du160C+3lkciMnAZGvtyr7zPg1U09tp/WXftu/R21F9Ejhy9U5cX8Lxq9611p692j4ApxtjTrPWdmxCrOtljEkHfgHU0Ns2G+rDtWwbyzfXO1vdicA04OTVpoLd3ZeY3Wodx5QAU1d1Cvs6lp8ARwB3W2sXG2Nep7ez+k9r7Vqnoa3HOcAewFn01km7aW07WWt7jDHHA28BHwMZwPEbMpJEREREZJhTP3sDqJ/dq6/0xvHAZ/TWZs4HdrbW9nzDoQ8BE4CjrbVfbMi1RCT5qHyEiGx2qzqqpld23932L+gdcfB/m3j6F1Z1VPuuVU/vlKbxq+2zaprbpV+9q7+OxOFXR5i+BbiB0ZsY6zd5AJgEnGqtDQ/guGPoHW2w+r+6AV57b6CD3g7f6m5bzzGPrj5KwFr7ad85xq/ziIHpAVad/zVr7ToX7LDWvkNvHeEC4Ferj14RERERGanUz95g6mf/71wLgKvo7Tc/YK19a337G2POo/cGwyxr7aOben0RcY5GCovIZmeMmUTv9Kbd6J0GtbrcTTz92mrGttJ713uViUC3tXbxRp6zte//Dap3tjGMMTfTuxLyRdbaFwZ4+HsDXQBjLcYAy6y10a9sX7CeYzak7TfFNfSOxPgvcJYx5mFr7Zfr2f8DeqfarW1Eh4iIiMiIo372N1M/e60+6Pt/vf1mY8yOwHXAv4EzB+naIuIQJYVFZLMyxmTTW/MqQu8d6QVAkN56VE+x5gyG9U33d69je3xdlx5YpEN+znUyxlwBnAvcZK29fiiuMUSGrJ2MMdvSWzbiQXprBn8J3G+M+Z7KQoiIiIion70h1M/eeMaYEnpfR13Ajwc4wlpEkpCSwiKyue1Mb02sna21b67a2FfXK+8r+64aKbDG9r7FGco2IYYFwB7GmHEDGMWwWRhjzqd3QY57rbW/dDCUpcB3jTFpXxnFMGkTzzvgBK4xJg2YDTQA51lrW40xFwF3Aj8H7trEmERERERGAvWz10P97I1njHHTmxAuAw5Mtp+tiGwc1RQWkc1t1V3ur97V/iVff09aCsSAXb6y/UzWPYJhQzzd9/81xpg1rrmOBTA2C2PMacCNwKP0Jjud9Fcgm94Vpld31iaed1UttK9+MFmfi4HpwOl9qytDb/2594DrjDGVmxiTiIiIyEigfvY6qJ+9ya4BfgDcbK19fgjOLyIO0EhhERkKextjStey/R3gH0Aj8DtjzJ1AJ72d0f8Dmlff2VrbYYx5HPh5Xyfyc+C7wI5A08YGZ6191xhzL731ZkcZY56nd2rdloAfOHpjz/1Vxph9gZl9D3N6N5lL+x5/Zq39c99++wN3ALXA68BPv9Jv/o+19j+DFdcGuB84ld6VkKfTu4jIbkBV3/c3diTCF/S29anGmC56f/6fW2s/X9vOxphp9CaF/2itfXbV9r6Vkk8APqV3pPD+GxmPiIiIyHCifnYf9bO/ZkD97A1ljNkduABoAb40xhy5jl3rrbWvbcq1RGTzUlJYRIbCAX3/vupaa+1bxpg9gJuBS+gd0fB3eu88/30tx5xN73vVUfSOcHiD3qlx610VdwP8HPgPcDLwKyAEzAF+u4nn/aqD6V2leJXcvusBPAL8ue/rregd1VEGPLyW81zVF+9mYa3tMcb8kN5VkH8GJICX6F1ReiG97bUx5+0yxhxN7/O5C0jr+/prndW+0SWzgR7gtLWca44x5lrgKmPModbap7+6j4iIiMgIcwDqZ6+ifvaa593gfvYAfZfe9ssHHlrPfm8BSgqLDCNG6/OIiMiGMsZsDXwM/NRa+4TT8YiIiIiIjATqZ4vI5qaawiIislZ9i5J81Xn0jmbY1BEkIiIiIiIpSf1sEUkGSVs+whhzDnACvfV0/kvvtIoyele8LKD3DtpR1tqIY0GKiIxsTxlj2oCPAC+wN7116e601q50MjAREdk06muLiDhK/WwRcVxSlo8wxlQA7wLT+urtPE3vCp17Ac9Za5/qK17/mbX2HidjFREZqYwxp9JbC24svQuDLKa3DttN1tqEg6GJiMgmUF9bRMRZ6meLSDJI5qTw+/SuJNoB/Ine1UIfB0qttTFjzHeBK621uzsWqIiIiIjIMKO+toiIiIgkZU3hvukSNwPVQC3QTu8UtjZrbaxvtxVAhTMRioiIiIgMT+pri4iIiEhS1hQ2xuQB+9M7laINeAbYYwDHnwScBHDyySd/u7S0dAiiHJni8Thut9vpMIYNtdfAqL0GTm02MGq
  781. "text/plain": [
  782. "<Figure size 1728x1728 with 4 Axes>"
  783. ]
  784. },
  785. "metadata": {
  786. "needs_background": "light"
  787. },
  788. "output_type": "display_data"
  789. }
  790. ],
  791. "source": [
  792. "fig = plt.figure(figsize=(24,24))\n",
  793. "gs = gridspec.GridSpec(2, 2, width_ratios=[1, 1])\n",
  794. "\n",
  795. "ax1 = plt.subplot(gs[0])\n",
  796. "plt.title(r\"Launch 12 Flight XY\")\n",
  797. "plt.xlabel(r\"Magnetic Field X [$ \\mu$T]\")\n",
  798. "plt.ylabel(r\"Magnetic Field Y [$ \\mu$T]\")\n",
  799. "ax1.plot(l12_m_raw_x, l12_m_raw_y, lw=0.8, label=\"\")\n",
  800. "ax1.add_patch(patches.Circle((0, 0), 53, edgecolor=\"#cccccc\", linewidth=1.0, linestyle='--', fill=False))\n",
  801. "plt.xlim([-80,80])\n",
  802. "plt.ylim([-80,80])\n",
  803. "\n",
  804. "ax2 = plt.subplot(gs[1])\n",
  805. "plt.title(r\"Launch 12 Flight YZ\")\n",
  806. "plt.xlabel(r\"Magnetic Field Y [$ \\mu$T]\")\n",
  807. "plt.ylabel(r\"Magnetic Field Z [$ \\mu$T]\")\n",
  808. "ax2.plot(l12_m_raw_y, l12_m_raw_z, lw=0.8, label=\"\")\n",
  809. "ax2.add_patch(patches.Circle((0, 0), 53, edgecolor=\"#cccccc\", linewidth=1.0, linestyle='--', fill=False))\n",
  810. "plt.xlim([-80,80])\n",
  811. "plt.ylim([-80,80])\n",
  812. "\n",
  813. "ax3 = plt.subplot(gs[2])\n",
  814. "plt.title(r\"Launch 12 Flight XZ\")\n",
  815. "plt.xlabel(r\"Magnetic Field X [$ \\mu$T]\")\n",
  816. "plt.ylabel(r\"Magnetic Field Z [$ \\mu$T]\")\n",
  817. "ax3.plot(l12_m_raw_x, l12_m_raw_z, lw=0.8, label=\"\")\n",
  818. "ax3.add_patch(patches.Circle((0, 0), 53, edgecolor=\"#cccccc\", linewidth=1.0, linestyle='--', fill=False))\n",
  819. "plt.xlim([-80,80])\n",
  820. "plt.ylim([-80,80])\n",
  821. "\n",
  822. "ax4 = plt.subplot(gs[3], projection='3d')\n",
  823. "plt.title(r\"Launch 12 Flight XYZ\")\n",
  824. "plt.xlabel(r\"Field Strength Y [$\\mu$T]\")\n",
  825. "plt.ylabel(r\"Field Strength Z [$\\mu$T]\")\n",
  826. "ax4.set_zlabel('Field Strength X [$\\mu$T]')\n",
  827. "\n",
  828. "ax4.plot_wireframe(x, y, z, color=\"k\", alpha=0.1, lw=0.2)\n",
  829. "\n",
  830. "ax4.plot(l12_m_raw_y, l12_m_raw_z, l12_m_raw_x, '-', lw=0.5)\n",
  831. "ax4.plot([0],[0],[0], 'g.')\n",
  832. "\n",
  833. "ax4.set_xlim(-60, 60)\n",
  834. "ax4.set_ylim(-60, 60)\n",
  835. "ax4.set_zlim(-60, 60)\n",
  836. "\n",
  837. "plt.show()"
  838. ]
  839. },
  840. {
  841. "cell_type": "markdown",
  842. "metadata": {},
  843. "source": [
  844. "The rocket spins during flight, and we see the magnetic field measurement spiral around the plots. We also see the familiar stretch and offset that we saw in the calibration data.\n",
  845. "\n",
  846. "\n",
  847. "### Calibrated Flight Data\n",
  848. "\n",
  849. "So now lets calibrate the flight data!"
  850. ]
  851. },
  852. {
  853. "cell_type": "code",
  854. "execution_count": 17,
  855. "metadata": {},
  856. "outputs": [],
  857. "source": [
  858. "l12_m_x = []\n",
  859. "l12_m_y = []\n",
  860. "l12_m_z = []\n",
  861. "for i, t in enumerate(l12_time):\n",
  862. " s = np.array((l12_m_raw_x[i], l12_m_raw_y[i], l12_m_raw_z[i])).reshape(3, 1)\n",
  863. " s = apply_mag_correction(s)\n",
  864. " mx, my, mz = s[0][0], s[1][0], s[2][0]\n",
  865. " l12_m_x.append(mx)\n",
  866. " l12_m_y.append(my)\n",
  867. " l12_m_z.append(mz)"
  868. ]
  869. },
  870. {
  871. "cell_type": "code",
  872. "execution_count": 18,
  873. "metadata": {
  874. "class": "fullwidth"
  875. },
  876. "outputs": [
  877. {
  878. "data": {
  879. "image/png": "iVBORw0KGgoAAAANSUhEUgAABYUAAAVhCAYAAADbYOFNAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8rg+JYAAAACXBIWXMAAAsTAAALEwEAmpwYAAEAAElEQVR4nOzdeXxcVf3/8deZPZPJZN+adEu6F2RT+ILKpmyyySYIIruCIiAqoIKCgiyCCiKLgqCAiiiLCqKyiSiiKCBQ2qZb0jRp9n32mfP7I2l+bSnQJcmdZN7PxyOPNJM7974nk2bOfO65n2OstYiIiIiIiIiIiIhIbnA5HUBEREREREREREREJo6KwiIiIiIiIiIiIiI5REVhERERERERERERkRyiorCIiIiIiIiIiIhIDlFRWERERERERERERCSHqCgsIiIiIiIiIiIikkNUFBaRSckYc68xxhpjPE5nGSsjj+dKh459rzFmzRju70pjjN3stueMMc/twD6tMeb+HQ4nIiIiIu9I4+wxP7bG2SKSlVQUFskRxpj9R17sz3Y6S7YwxpxrjHnQGLN65GfzwjtsV2uMudwY8w9jTKcxps8Y87Ix5jPGGPdWHuvKkWNs6eNXY/vIxs9Gv0db+ljvdD4AY4xr5Of98a3cvtAY02KMecsY49/C97898viOMcYcYIzJGGNueYd9fcAYkzLG3Ltjj0JEREQmC42z307j7G031cbZxpgFxpiYMeZ37/D9GcaYAWPMs2bYOz32jT/WjPVjEsllU+bMn4jIdrgMKAT+DZS9y3YfBy4H/gA8CCSAw4E7gf2Bk7fhmF8G2ja7bc3I5zwgtQ37ctLdwHOb3RZ9j/scPD5R3sYFfBP4GfDoe21sre0zxnweeBi4guHnGgBjzC7ApcBvrLWPjNz2E+DzxphfWmtf3GhbD3AX0AF8ccwejYiIiMjko3H29psS42xr7VJjzFXAd4wxJ1prH9xsk9sBN3C2tdYaY059l90dAxwL/GO7k4vI26goLCK5bD+gaWQQsuZdtnsOmGmt3XiQebsx5ufAqcaY71prX9nKYz5mrV2xpW9Ya2NbuY9s8E9r7TZdYmatTYxXmB1lrX3EGPNb4BJjzK+ttf8bmZ1yNzAAnL/R5l8BPgbcbYzZdaPHdQnwPuA4a23PROYXERERyTIaZ2+/qTTO/i5wPHCLMeYpa20XgDHmZIbH01+y1q4EeKfHbIxZzPBJggbg3AlJLZIj1D5CRDZhjDnKGPOoMWatMSZujFlvjPmZMWbaZtttuLxp/y3sY83Gl88bY2aNbHu1MeYTxpg3Ri4lWm6M+cQ75DjTGPNvY8yQMabXGPOiMeaTW9i01Bhz38g2AyOXqZVszWO11jZaa+1WbPfGZgPVDR4a+bx4a473Xsxmvc6MMY+PPP75m213mzEmvfHP3hhTboz50cjzljDGNBpjrjfGBLZwnC8aY1aNPAevGGOOGIv878VsodeZMcZrjLnWDLdviBhj/m6M+b8tbbvRfQ4c+d2IjTzOL2z0vVlAcuTL0za61GyL+9rM+cAQw8VeN8OzTfYALt74+bfW9gOfBRYyMqvYGDOP4VnGv7HWPrwVxxIREZEco3H2FrfTOHsMZOs421qbAs4ESoDvj+ynFPgB8NLI53d7XAXAbwEDHD8yDheRMaKZwiKyuTMZftG9DegEFgDnAP9njNllB8+yHwacAdwB9AJnA780xrxqrV2+YSNjzK3A54EXgCsZvlxqN4bPJv9ys30+DqwCvgrMZ7iwlwDe7fKjsbJhAN+xDfcpNsZsfgldt7U2s4VtzwJeB+43xuxjrU0aYz4GnAd811r7HIwOrP4JhIEfA40MFzO/BOxijDlsw6DcGPN14GrgbwwPwmqAX4zcZ1uEtvA4Bqy18W3cz50M/078DvgTw79vTwDdQPMWtt8F+BXwE+Ae4ESGZx4ssdY+zfBzcRrDl7T9jeGfB7z9UsK3sdauN8Z8ieHZwT8cyfVna+3PtrDtE2Z4MY7LjDG/AW4BImw6o1hERERkYxpnbz2Nszc1acfZ1trXjDHXAZcbY37BcEuQQuCsd3huNnYXw797Z1hr//ce24rItrLW6kMf+siBD4Z7clmGeza923bBLdy238h9P7mF/e2/he3XAPdu9PWskW0HgZqNbq8E4sANG932oZFtfwG4Ntuv2ejf945s94PNtvkBw/3Cwtv481kDvLAN2+cBy4B1gH8rtr9yJO+WPmaNbGOBKze739Ejt18NlAPrgVcB30bb3Ab0MHzp3cb3PX/kvoeOfF0KxIC/A56NtjtsZLs12/B7tKWP0zd/vJvd9znguY2+3nnDc73ZdmeN3P7cZrdbIA3svtFtfoYHor/e6DbPyLb3vtfjeYfH+NRGv68z32W70pHno3Nk+09vz/H0oQ996EMf+tDH5P5A4+z3+vmsQePsNdvwezTlxtmAD3gT6Bq5/ze24j4Xjmx797YcSx/60MfWf2imsIhswlobATDGGKCA//8C3gt8gLfPINgWv7PWrtvoWG3GmKVA/UbbbLjM7XK72Zlja63dwj5v2+zrvzI8gJjJ8Nn/8XIXMA842m7bWfvTgJbNbnvH1YSttY8ZY+5meLGOgxk+q/4RO9I3bOR5OhH4MzC02ayCP498/gjwJHAQw4O7W+3wpVwbjvFHY8xbQHAbHsf3gD9udtub23B/GF5EBODmzW7/GXDjO9znJWvtfzd8Ya2NG2P+yaa/Qzuqa+Tzm9bad5zZYa3tMsZcyPCMir9Ya38+hhlERERkitE4e6tpnD3FxtnW2oQx5kyGZ10vAa59t+2NMXsz3I/4NYZntovIOFBRWEQ2YYZ7o17H8MAmtNm3i3Zw91sqsPUw3GNqg7nAkLV21Xbus2fk81b1O9sexpgbGb7s6avW2t9t493/Yd9hAYx3cRHDsww+AHzZWrvxoLCc4cf6Cf7/QH9zFSOfZ418XraFbZYxfOng1nrLWvvUNmy/JRvyNGx8o7U2ZYxZ/Q73eaffofftYBZguNcfwz/H/wF7GmPOttbe9S53eWnk87/G4vgiIiIydWmc/d40zgam6DjbWvvScJ2dV6y1yXfabqT4/muGW5scbyfXIoEik4qKwiIyyhgTBp5nuFfYVQwPIiIMX7bzKzZdnHJLswk2cL/D7el3OvS2JR33fb4jY8w3Ge4h9l1r7XXjcYwt2A2oGvn3LptHGvn8KPCjd7h/6zhkcsq4Pd/GmEKGZ8S8BezN8GyY7xpj/mCtfcdZJiIiIiLvRePs96ZxtuMm9Pne4oGMcTHc3qQWOG47ivwisg1UFBaRjR3AcP+xA+zI4goAxpg8oHizbTfMFNjk9pFVeKt3IEMDcKgxpm4bZjFMCGPMVxju4XWHtfaSCTpmGPg5w4t8PAp82RjzmLX2tyObdAD9QGArZhSsGfk8H/jvZt+bz8RbM/J5Lv9/xi3GGA8wm+HLxbbHu72Rejc3MPy7+wlrbcQYczbwMsOLzp2wnfsUERERAY2z35XG2WNuzcjnbBlnb60rGZ5J/31r7cPjfCyRnOd6701EJIdsODu8+dngS3j734s1DC80ceBmt1/AO89g2Bq/Hvl89ciZ4lEjfb0cYYz5PMNFw/uAz03goX8ITGd4leevAv8G7jTGVAOM9IN7EDjEGLPf5nc2xgSMMQUjX/6F4QVHzh8ZEG7Y5jBg4bg+ii17YuTzhZvdfho7cAmltTbN8EIfm7/BekfGmP0ZXv37NmvtP0b28xrDPdeOH2krISIiIrK9NM5+Bxpnj4usGWdvrZGf1eXAPxj+fyEi40wzhUVyz+HGmKot3P43hlfL7QB+boy5FRhgeDD6Af7/4lsAWGv7jTEPAJ8bGUS+wfAl9x8GOrc3nLX2BWPMHcC5wHRjzGMMX1q3KxAAPr29+96cMeZI/v9lYoXDN5nLR75+zVr7+5HtjmZ40NgKPA2cstm4+X/W2v+NVa6N8h3P8OP9lrX2nyO3fQp4Bfgpw/3PYHgQux/wF2PMzxienZDH8AIdJwDHMbzCcJcx5jsMX7L4jDHmIaCG4cH3GwwveDJhrLX/M8b8HPi
  880. "text/plain": [
  881. "<Figure size 1728x1728 with 4 Axes>"
  882. ]
  883. },
  884. "metadata": {
  885. "needs_background": "light"
  886. },
  887. "output_type": "display_data"
  888. }
  889. ],
  890. "source": [
  891. "fig = plt.figure(figsize=(24,24))\n",
  892. "gs = gridspec.GridSpec(2, 2, width_ratios=[1, 1])\n",
  893. "\n",
  894. "ax1 = plt.subplot(gs[0])\n",
  895. "plt.title(r\"Launch 12 Fixed Flight XY\")\n",
  896. "plt.xlabel(r\"Magnetic Field X [$ \\mu$T]\")\n",
  897. "plt.ylabel(r\"Magnetic Field Y [$ \\mu$T]\")\n",
  898. "ax1.plot(l12_m_x, l12_m_y, lw=0.8, label=\"\")\n",
  899. "ax1.add_patch(patches.Circle((0, 0), 53, edgecolor=\"#cccccc\", linewidth=1.0, linestyle='--', fill=False))\n",
  900. "plt.xlim([-60,60])\n",
  901. "plt.ylim([-60,60])\n",
  902. "\n",
  903. "ax2 = plt.subplot(gs[1])\n",
  904. "plt.title(r\"Launch 12 Fixed Flight YZ\")\n",
  905. "plt.xlabel(r\"Magnetic Field Y [$ \\mu$T]\")\n",
  906. "plt.ylabel(r\"Magnetic Field Z [$ \\mu$T]\")\n",
  907. "ax2.plot(l12_m_y, l12_m_z, lw=0.8, label=\"\")\n",
  908. "ax2.add_patch(patches.Circle((0, 0), 53, edgecolor=\"#cccccc\", linewidth=1.0, linestyle='--', fill=False))\n",
  909. "plt.xlim([-60,60])\n",
  910. "plt.ylim([-60,60])\n",
  911. "\n",
  912. "ax3 = plt.subplot(gs[2])\n",
  913. "plt.title(r\"Launch 12 Fixed Flight XZ\")\n",
  914. "plt.xlabel(r\"Magnetic Field X [$ \\mu$T]\")\n",
  915. "plt.ylabel(r\"Magnetic Field Z [$ \\mu$T]\")\n",
  916. "ax3.plot(l12_m_x, l12_m_z, lw=0.8, label=\"\")\n",
  917. "ax3.add_patch(patches.Circle((0, 0), 53, edgecolor=\"#cccccc\", linewidth=1.0, linestyle='--', fill=False))\n",
  918. "plt.xlim([-60,60])\n",
  919. "plt.ylim([-60,60])\n",
  920. "\n",
  921. "ax4 = plt.subplot(gs[3], projection='3d')\n",
  922. "plt.title(r\"Launch 12 Fixed Flight XYZ\")\n",
  923. "plt.xlabel(r\"Field Strength Y [$\\mu$T]\")\n",
  924. "plt.ylabel(r\"Field Strength Z [$\\mu$T]\")\n",
  925. "ax4.set_zlabel('Field Strength X [$\\mu$T]')\n",
  926. "\n",
  927. "ax4.plot_wireframe(x, y, z, color=\"k\", alpha=0.1, lw=0.2)\n",
  928. "\n",
  929. "ax4.plot(l12_m_y, l12_m_z, l12_m_x, '-', lw=0.5)\n",
  930. "ax4.plot([0],[0],[0], 'g.')\n",
  931. "\n",
  932. "ax4.set_xlim(-60, 60)\n",
  933. "ax4.set_ylim(-60, 60)\n",
  934. "ax4.set_zlim(-60, 60)\n",
  935. "\n",
  936. "plt.show()"
  937. ]
  938. },
  939. {
  940. "cell_type": "markdown",
  941. "metadata": {},
  942. "source": [
  943. "And it looks like the calibration did a reasonable job. The values now come very close to landing on the nominal Earth field sphere. The XY view is still off a little bit but it might just be that we had some bias in the calibration run. It's still a huge improvement to the original dataset and it now usable in IMU reconstructions of the flight of the rocket."
  944. ]
  945. },
  946. {
  947. "cell_type": "markdown",
  948. "metadata": {},
  949. "source": [
  950. "----------------------------------------------\n",
  951. "\n",
  952. "This post is written as a [jupyter notebook](https://jupyter.org/) and all its code and data can be viewed as a stand-alone own project here:\n",
  953. "\n",
  954. "<https://git.natronics.org/natronics/psas-magnetometer-calibration>"
  955. ]
  956. }
  957. ],
  958. "metadata": {
  959. "celltoolbar": "Edit Metadata",
  960. "kernelspec": {
  961. "display_name": "Python 3 (ipykernel)",
  962. "language": "python",
  963. "name": "python3"
  964. },
  965. "language_info": {
  966. "codemirror_mode": {
  967. "name": "ipython",
  968. "version": 3
  969. },
  970. "file_extension": ".py",
  971. "mimetype": "text/x-python",
  972. "name": "python",
  973. "nbconvert_exporter": "python",
  974. "pygments_lexer": "ipython3",
  975. "version": "3.8.6"
  976. }
  977. },
  978. "nbformat": 4,
  979. "nbformat_minor": 2
  980. }