Initial commit

This commit is contained in:
Sebastian 2022-01-16 16:10:37 +01:00
commit 36c6bd7e46
1 changed files with 409 additions and 0 deletions

409
pressure_wave_data.ipynb Normal file
View File

@ -0,0 +1,409 @@
{
"cells": [
{
"cell_type": "code",
"execution_count": 1,
"id": "4dbd5465",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Defaulting to user installation because normal site-packages is not writeable\n",
"Requirement already satisfied: geopy in /home/sebastian/.local/lib/python3.10/site-packages (2.2.0)\n",
"Requirement already satisfied: geographiclib<2,>=1.49 in /home/sebastian/.local/lib/python3.10/site-packages (from geopy) (1.52)\n",
"Defaulting to user installation because normal site-packages is not writeable\n",
"Requirement already satisfied: ipyleaflet in /home/sebastian/.local/lib/python3.10/site-packages (0.15.0)\n",
"Requirement already satisfied: traittypes<3,>=0.2.1 in /home/sebastian/.local/lib/python3.10/site-packages (from ipyleaflet) (0.2.1)\n",
"Requirement already satisfied: ipywidgets<8,>=7.6.0 in /usr/lib/python3.10/site-packages (from ipyleaflet) (7.6.5)\n",
"Requirement already satisfied: xyzservices>=2021.8.1 in /home/sebastian/.local/lib/python3.10/site-packages (from ipyleaflet) (2021.11.0)\n",
"Requirement already satisfied: traitlets>=4.2.2 in /usr/lib/python3.10/site-packages (from traittypes<3,>=0.2.1->ipyleaflet) (5.1.0)\n",
"Defaulting to user installation because normal site-packages is not writeable\n",
"Collecting tabulate\n",
" Downloading tabulate-0.8.9-py3-none-any.whl (25 kB)\n",
"Installing collected packages: tabulate\n",
"\u001b[33m WARNING: The script tabulate is installed in '/home/sebastian/.local/bin' which is not on PATH.\n",
" Consider adding this directory to PATH or, if you prefer to suppress this warning, use --no-warn-script-location.\u001b[0m\n",
"Successfully installed tabulate-0.8.9\n"
]
}
],
"source": [
"import sys\n",
"!{sys.executable} -m pip install geopy\n",
"!{sys.executable} -m pip install ipyleaflet\n",
"!{sys.executable} -m pip install tabulate"
]
},
{
"cell_type": "code",
"execution_count": 2,
"id": "1ced7adc",
"metadata": {},
"outputs": [],
"source": [
"from datetime import datetime\n",
"from pytz import timezone\n",
"\n",
"data_points = [\n",
" (\"Kaiserslautern\", datetime(2022, 1, 16, 20, 30, 0, tzinfo=timezone(\"Europe/Berlin\")), \"https://chaos.social/@sebastian/107628740679869429\"),\n",
" (\"Aachen\", datetime(2022, 1, 16, 20, 25, 0, tzinfo=timezone(\"Europe/Berlin\")), \"https://chaos.social/@trilader/107628784664752716\"),\n",
" (\"Hamburg\", datetime(2022, 1, 16, 20, 10, 0, tzinfo=timezone(\"Europe/Berlin\")),\"https://chaos.social/@tsia_/107628923762308800\"),\n",
" (\"Edinburgh\", datetime(2022, 1, 16, 18, 45, 0, tzinfo=timezone(\"UTC\")), \"https://toot.cat/@river/107629146900457261\"),\n",
" (\"Copenhagen\", datetime(2022, 1, 16, 19, 55, 0, tzinfo=timezone(\"CET\")), \"https://mcd.dk/@rune/107629163090789143\"),\n",
" (\"Augsburg\", datetime(2022, 1, 16, 20, 38, 0, tzinfo=timezone(\"Europe/Berlin\")), \"https://chaos.social/@phjl/107629214329242643\"),\n",
" (\"Frankfurt\", datetime(2022, 1, 16, 20, 25, 0, tzinfo=timezone(\"Europe/Berlin\")), \"https://mastodon.social/@hko/107629722849924036\"),\n",
" (\"Melbourne\", datetime(2022, 1, 16, 8, 0, 0, tzinfo=timezone(\"UTC\")), \"https://aus.social/@futzle/107631253370573653\"),\n",
"]"
]
},
{
"cell_type": "markdown",
"id": "1f554c69",
"metadata": {},
"source": [
"Since we know roughly where and when the explosion happened we can do a quick sanity check of our data points."
]
},
{
"cell_type": "code",
"execution_count": 26,
"id": "6b525fc1",
"metadata": {},
"outputs": [],
"source": [
"from geopy.distance import distance\n",
"from geopy.geocoders import Nominatim\n",
"\n",
"geolocator = Nominatim(user_agent=\"pressure_wave_triangulation\")\n",
"\n",
"ORIGIN = (-20.5, -175.4)\n",
"T_ZERO = datetime(2022, 1, 16, 4, 27, 0, tzinfo=timezone(\"UTC\"))\n",
"\n",
"SPEED_OF_SOUND = 343.0 # Speed of sound in air at 20°C\n",
"\n",
"enhanced_data_points = []\n",
"\n",
"for location, ts, url in data_points:\n",
" geolocation = geolocator.geocode(location)\n",
" coordinates = (geolocation.latitude, geolocation.longitude)\n",
" dist = distance(ORIGIN, coordinates)\n",
" \n",
" delta_t = ts - T_ZERO\n",
" \n",
" speed = dist.meters / delta_t.total_seconds()\n",
" mach = speed / SPEED_OF_SOUND\n",
"\n",
" enhanced_data_points += [(location, coordinates, ts, dist, delta_t, speed, mach, url)]\n"
]
},
{
"cell_type": "code",
"execution_count": 34,
"id": "842f682a",
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
"<table>\n",
"<tbody>\n",
"<tr><td>Location </td><td>Coordinates </td><td>Timestamp </td><td>Distance </td><td>Time difference</td><td>Speed [m/s] </td><td>Mach Number </td><td>Toot </td></tr>\n",
"<tr><td>Kaiserslautern</td><td>49.443217 7.768995 </td><td>2022-01-16 20:30:00+00:53</td><td>16780.572549711476 km</td><td>15:10:00 </td><td>307.3364935844593 </td><td>0.8960247626369076</td><td>https://chaos.social/@sebastian/107628740679869429</td></tr>\n",
"<tr><td>Aachen </td><td>50.776351 6.083862 </td><td>2022-01-16 20:25:00+00:53</td><td>16641.951165318445 km</td><td>15:05:00 </td><td>306.4816052544834 </td><td>0.8935323768352286</td><td>https://chaos.social/@trilader/107628784664752716 </td></tr>\n",
"<tr><td>Hamburg </td><td>53.550341 10.000654 </td><td>2022-01-16 20:10:00+00:53</td><td>16307.144415681687 km</td><td>14:50:00 </td><td>305.3772362487207 </td><td>0.8903126421245501</td><td>https://chaos.social/@tsia_/107628923762308800 </td></tr>\n",
"<tr><td>Edinburgh </td><td>55.953346 -3.188375 </td><td>2022-01-16 18:45:00+00:00</td><td>16015.614674425336 km</td><td>14:18:00 </td><td>311.1036261543383 </td><td>0.9070076564266423</td><td>https://toot.cat/@river/107629146900457261 </td></tr>\n",
"<tr><td>Copenhagen </td><td>55.686724 12.570072 </td><td>2022-01-16 19:55:00+01:00</td><td>16042.101561186751 km</td><td>14:28:00 </td><td>308.02806377086694</td><td>0.8980410022474254</td><td>https://mcd.dk/@rune/107629163090789143 </td></tr>\n",
"<tr><td>Augsburg </td><td>48.366804 10.898697 </td><td>2022-01-16 20:38:00+00:53</td><td>16861.97644995329 km </td><td>15:18:00 </td><td>306.13610112478744</td><td>0.8925250761655611</td><td>https://chaos.social/@phjl/107629214329242643 </td></tr>\n",
"<tr><td>Frankfurt </td><td>50.110644 8.682092 </td><td>2022-01-16 20:25:00+00:53</td><td>16699.02015177992 km </td><td>15:05:00 </td><td>307.53259948029324</td><td>0.8965964999425459</td><td>https://mastodon.social/@hko/107629722849924036 </td></tr>\n",
"<tr><td>Melbourne </td><td>-37.814218 144.963161</td><td>2022-01-16 08:00:00+00:00</td><td>4264.645211600706 km </td><td>3:33:00 </td><td>333.6968084194606 </td><td>0.9728769924765615</td><td>https://aus.social/@futzle/107631253370573653 </td></tr>\n",
"</tbody>\n",
"</table>"
],
"text/plain": [
"'<table>\\n<tbody>\\n<tr><td>Location </td><td>Coordinates </td><td>Timestamp </td><td>Distance </td><td>Time difference</td><td>Speed [m/s] </td><td>Mach Number </td><td>Toot </td></tr>\\n<tr><td>Kaiserslautern</td><td>49.443217 7.768995 </td><td>2022-01-16 20:30:00+00:53</td><td>16780.572549711476 km</td><td>15:10:00 </td><td>307.3364935844593 </td><td>0.8960247626369076</td><td>https://chaos.social/@sebastian/107628740679869429</td></tr>\\n<tr><td>Aachen </td><td>50.776351 6.083862 </td><td>2022-01-16 20:25:00+00:53</td><td>16641.951165318445 km</td><td>15:05:00 </td><td>306.4816052544834 </td><td>0.8935323768352286</td><td>https://chaos.social/@trilader/107628784664752716 </td></tr>\\n<tr><td>Hamburg </td><td>53.550341 10.000654 </td><td>2022-01-16 20:10:00+00:53</td><td>16307.144415681687 km</td><td>14:50:00 </td><td>305.3772362487207 </td><td>0.8903126421245501</td><td>https://chaos.social/@tsia_/107628923762308800 </td></tr>\\n<tr><td>Edinburgh </td><td>55.953346 -3.188375 </td><td>2022-01-16 18:45:00+00:00</td><td>16015.614674425336 km</td><td>14:18:00 </td><td>311.1036261543383 </td><td>0.9070076564266423</td><td>https://toot.cat/@river/107629146900457261 </td></tr>\\n<tr><td>Copenhagen </td><td>55.686724 12.570072 </td><td>2022-01-16 19:55:00+01:00</td><td>16042.101561186751 km</td><td>14:28:00 </td><td>308.02806377086694</td><td>0.8980410022474254</td><td>https://mcd.dk/@rune/107629163090789143 </td></tr>\\n<tr><td>Augsburg </td><td>48.366804 10.898697 </td><td>2022-01-16 20:38:00+00:53</td><td>16861.97644995329 km </td><td>15:18:00 </td><td>306.13610112478744</td><td>0.8925250761655611</td><td>https://chaos.social/@phjl/107629214329242643 </td></tr>\\n<tr><td>Frankfurt </td><td>50.110644 8.682092 </td><td>2022-01-16 20:25:00+00:53</td><td>16699.02015177992 km </td><td>15:05:00 </td><td>307.53259948029324</td><td>0.8965964999425459</td><td>https://mastodon.social/@hko/107629722849924036 </td></tr>\\n<tr><td>Melbourne </td><td>-37.814218 144.963161</td><td>2022-01-16 08:00:00+00:00</td><td>4264.645211600706 km </td><td>3:33:00 </td><td>333.6968084194606 </td><td>0.9728769924765615</td><td>https://aus.social/@futzle/107631253370573653 </td></tr>\\n</tbody>\\n</table>'"
]
},
"execution_count": 34,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"import tabulate\n",
"\n",
"table_data = [(\"Location\", \"Coordinates\", \"Timestamp\", \"Distance\", \"Time difference\", \"Speed [m/s]\", \"Mach Number\", \"Toot\")]\n",
"for location, coordinates, ts, dist, delta_t, speed, mach, url in enhanced_data_points:\n",
" table_data += [(location, \"%f %f\" % coordinates, ts, dist, delta_t, speed, mach, url)]\n",
"\n",
"table = tabulate.tabulate(table_data, tablefmt='html')\n",
"table\n"
]
},
{
"cell_type": "markdown",
"id": "5341fa32",
"metadata": {},
"source": [
"That looks quite consistent, considering the we pretty much guestimated the timestamps from screenshots."
]
},
{
"cell_type": "code",
"execution_count": 120,
"id": "5ff0ee40",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Kaiserslautern 322.0800766318187\n",
"Aachen 306.3924371618201\n",
"Hamburg 282.48002118519315\n",
"Edinburgh 315.16138946596874\n",
"Copenhagen 283.03327184485505\n",
"Augsburg 314.29962074367353\n",
"Frankfurt 319.0927106279109\n"
]
}
],
"source": [
"relative_speed_data = []\n",
"\n",
"for location, coordinates, ts, _, _, _, _, url in enhanced_data_points:\n",
" count = 0\n",
" speed = 0.0\n",
" for _, other_coordinates, other_ts, _, _, _, _, _ in enhanced_data_points:\n",
" relative_delay = other_ts - ts\n",
" if abs(relative_delay.total_seconds()) < 1:\n",
" continue\n",
" \n",
" relative_dist = distance(coordinates, other_coordinates)\n",
" relative_speed = relative_dist.meters / abs(relative_delay.total_seconds())\n",
" \n",
" if relative_speed > 350:\n",
" continue\n",
" \n",
" speed += relative_speed\n",
" count += 1.0\n",
" \n",
" \n",
" if count < 1.0:\n",
" continue \n",
" \n",
" speed = speed / count\n",
" print(location, speed)\n",
" \n",
" relative_speed_data += [(location, coordinates, ts, speed)]\n",
" \n",
" "
]
},
{
"cell_type": "markdown",
"id": "5d40d782",
"metadata": {},
"source": [
"Some cleanup was required, so I ended up dropping everthing that is much faster then the speed of sound in air.\n",
"That should give use something to work with.\n",
"\n",
"Also the leaftlet map does something weird with the scaling.\n",
"If we cicle each of our locations with the distance we got from geopy, the circles never really intersect."
]
},
{
"cell_type": "code",
"execution_count": 121,
"id": "1a43cbc5",
"metadata": {},
"outputs": [
{
"data": {
"application/vnd.jupyter.widget-view+json": {
"model_id": "175f37ff4cf048d0b886861621cb4598",
"version_major": 2,
"version_minor": 0
},
"text/plain": [
"Map(center=[-37.8142176, 144.9631608], controls=(ZoomControl(options=['position', 'zoom_in_text', 'zoom_in_tit…"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"from ipyleaflet import Map, Marker, Circle\n",
"\n",
"\n",
"MELBOURNE = geolocator.geocode(\"Melbourne\")\n",
"\n",
"\n",
"leaflet_map = Map(center=(MELBOURNE.latitude, MELBOURNE.longitude), zoom=2)\n",
"marker = Marker(location=ORIGIN, draggable=False)\n",
"leaflet_map.add_layer(marker);\n",
"\n",
"for _, coordinates, _, dist, _, _, _, _ in enhanced_data_points:\n",
" marker = Marker(location=coordinates, draggable=False)\n",
" leaflet_map.add_layer(marker);\n",
" \n",
" circle = Circle(location=coordinates, radius=int(dist.meters), fill=False, weight=1)\n",
" leaflet_map.add_layer(circle)\n",
"\n",
"\n",
"display(leaflet_map)\n"
]
},
{
"cell_type": "markdown",
"id": "b62287c9",
"metadata": {},
"source": [
"However we can fiddle with that so that things end up in roughly the right spot. (Highly scientific. I know...)"
]
},
{
"cell_type": "code",
"execution_count": 122,
"id": "1f81d6a2",
"metadata": {},
"outputs": [
{
"data": {
"application/vnd.jupyter.widget-view+json": {
"model_id": "8f4b34318c5b4d9299ad23184f29b1c0",
"version_major": 2,
"version_minor": 0
},
"text/plain": [
"Map(center=[-37.8142176, 144.9631608], controls=(ZoomControl(options=['position', 'zoom_in_text', 'zoom_in_tit…"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"from ipyleaflet import Map, Marker, Circle\n",
"\n",
"\n",
"FUDGE_FACTOR = 0.821\n",
"\n",
"\n",
"leaflet_map = Map(center=(MELBOURNE.latitude, MELBOURNE.longitude), zoom=2)\n",
"marker = Marker(location=ORIGIN, draggable=False)\n",
"leaflet_map.add_layer(marker);\n",
"\n",
"for _, coordinates, _, dist, _, _, _, _ in enhanced_data_points:\n",
" marker = Marker(location=coordinates, draggable=False)\n",
" leaflet_map.add_layer(marker);\n",
" \n",
" circle = Circle(location=coordinates, radius=int(dist.meters * FUDGE_FACTOR), fill=False, weight=1)\n",
" leaflet_map.add_layer(circle)\n",
"\n",
"\n",
"display(leaflet_map)"
]
},
{
"cell_type": "markdown",
"id": "79dd3ef0",
"metadata": {},
"source": [
"Also I don't want to mess with the math required to actually propagate back the circles projected on a sphere until they actually intersec.\n",
"\n",
"So let's just assume we have a rough idea when the event occured and we can use that timestamp +-5min to draw some circles and use our eyes to check for intersections."
]
},
{
"cell_type": "code",
"execution_count": 135,
"id": "0cb443cc",
"metadata": {},
"outputs": [
{
"data": {
"application/vnd.jupyter.widget-view+json": {
"model_id": "19841666b00e498e9cf91b0588a12f7c",
"version_major": 2,
"version_minor": 0
},
"text/plain": [
"Map(center=[-20.5, -175.4], controls=(ZoomControl(options=['position', 'zoom_in_text', 'zoom_in_title', 'zoom_…"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"from ipyleaflet import Map, Marker, Circle\n",
"\n",
"\n",
"leaflet_map = Map(center=ORIGIN, zoom=2)\n",
"marker = Marker(location=ORIGIN, draggable=False)\n",
"leaflet_map.add_layer(marker);\n",
"\n",
"\n",
"T_MINUS = datetime(2022, 1, 16, 4, 22, 0, tzinfo=timezone(\"UTC\"))\n",
"T_PLUS = datetime(2022, 1, 16, 4, 32, 0, tzinfo=timezone(\"UTC\"))\n",
"\n",
"\n",
"for location, coordinates, ts, speed in relative_speed_data:\n",
" marker = Marker(location=coordinates, draggable=False)\n",
" leaflet_map.add_layer(marker);\n",
" \n",
" for start, color in [(T_MINUS, \"green\"), (T_ZERO, \"red\"), (T_PLUS, \"blue\")]:\n",
" dist = speed * (ts - start).total_seconds()\n",
"\n",
" circle = Circle(location=coordinates, radius=int(dist * FUDGE_FACTOR), fill=False, weight=1, color=color)\n",
" leaflet_map.add_layer(circle)\n",
"\n",
"\n",
"display(leaflet_map)"
]
},
{
"cell_type": "markdown",
"id": "457f9ad9",
"metadata": {},
"source": [
"While that looks close, just zoom out until you can see the entire map.\n",
"\n",
"It's a mess.\n",
"Not enough data points, the locations are not spread out wide enough, the timestamps are off, and there is some weird map projection stuff going on...."
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "5ed3a84b",
"metadata": {},
"outputs": [],
"source": []
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3 (ipykernel)",
"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.10.1"
}
},
"nbformat": 4,
"nbformat_minor": 5
}