summaryrefslogtreecommitdiff
path: root/plot.py
diff options
context:
space:
mode:
authorAki <please@ignore.pl>2022-08-30 00:15:26 +0200
committerAki <please@ignore.pl>2022-08-30 00:20:01 +0200
commitff0dd2f8c5e9ed5bd63fcc47e09d46a7b2812ff1 (patch)
treee81dfefa6c2f7bf5d56f50816c58055d5b57f29f /plot.py
parent389cbba8e19da1fef16310adcf9e1ca7bfe9dda1 (diff)
downloadszilagyi-ff0dd2f8c5e9ed5bd63fcc47e09d46a7b2812ff1.zip
szilagyi-ff0dd2f8c5e9ed5bd63fcc47e09d46a7b2812ff1.tar.gz
szilagyi-ff0dd2f8c5e9ed5bd63fcc47e09d46a7b2812ff1.tar.bz2
Wrapped script in a setup.py package
Diffstat (limited to 'plot.py')
-rw-r--r--plot.py131
1 files changed, 0 insertions, 131 deletions
diff --git a/plot.py b/plot.py
deleted file mode 100644
index eb2665d..0000000
--- a/plot.py
+++ /dev/null
@@ -1,131 +0,0 @@
-import csv
-import math
-import os
-import re
-from collections import deque
-
-import matplotlib.pyplot as plot
-
-
-def load(directory):
- def _read(iterable):
- for x, y in iterable:
- yield float(x), float(y)
-
- def _load(filename):
- with open(filename) as fd:
- reader = csv.reader(fd)
- return list(_read(reader))
-
- def _files(directory):
- for file in os.listdir(directory):
- match = re.match(r"SWI_(-?\d+)\.csv", file)
- if match:
- yield int(match.group(1)), os.path.join(directory, file)
-
- return [(x, _load(y)) for x, y in sorted(_files(directory), key=lambda x: x[0])]
-
-
-def look_downwards(data, x, start):
- for i in range(start, 0, -1):
- if data[i - 1][0] < x:
- break
- else:
- raise IndexError
- return i - 1
-
-
-def look_upwards(data, x, start):
- for i in range(start, len(data)):
- if data[i + 1][0] > x:
- break
- else:
- raise IndexError
- return i
-
-
-def find_segment(data, x):
- width = data[-1][0] - data[0][0]
- relative = x - data[0][0]
- candidate = math.floor(relative / width * len(data))
- look = look_downwards if data[candidate][0] > x else look_upwards # May raise IndexError
- candidate = look(data, x, candidate)
- return candidate, candidate + 1
-
-
-def find_boundary_curves(swis, x, y):
- segments = deque()
- for index, data in swis:
- i, j = find_segment(data, x)
- if data[i][1] > y and data[j][1] > y:
- segments.append((index, data, i, j))
- break
- if data[i][1] < y and data[j][1] < y:
- if segments:
- segments.popleft()
- segments.append((index, data, i, j))
- if len(segments) == 3:
- middle = segments[1][1]
- run = middle[j][0] - middle[i][0]
- if run == 0:
- raise RuntimeError # tidy up dataset
- slope = (middle[j][1] - middle[i][1]) / run
- intercept = middle[j][1] - slope * middle[j][0]
- value = slope * x + intercept
- if value == y:
- raise RuntimeError # Exactly on point; SWI == index
- if value < y:
- segments.popleft()
- else:
- segments.pop()
- if len(segments) == 1:
- raise RuntimeError # SWI == -10
- return segments
-
-
-def dist(x1, y1, x2, y2):
- return math.sqrt((x2 - x1) ** 2 + (y2 - y1) ** 2)
-
-
-def calculate_swi(segments, x, y):
- low = segments[0]
- high = segments[1]
- dist_to_low = min(dist(p[0], p[1], x, y) for p in (low[1][low[2]], low[1][low[2]]))
- dist_to_high = min(dist(p[0], p[1], x, y) for p in (high[1][high[2]], high[1][high[2]]))
- return dist_to_low / (dist_to_low + dist_to_high) * (high[0] - low[0]) + low[0]
-
-
-swis = load("dataset")
-
-
-def _swi(x, y): # dt, depth
- segments = find_boundary_curves(swis, x, y)
- return calculate_swi(segments, x, y)
-
-
-def grid(start, end, steps):
- step = (end - start) / steps
- i = start
- while i <= end:
- yield i
- i += step
-
-
-C = []
-X = list(grid(0, 40, 1000))
-Y = list(grid(0, 50000, 1000))
-for y_scaled in range(0, 1000):
- y = y_scaled * 50
- row = []
- for x_scaled in range(0, 1000):
- x = x_scaled / 25
- try:
- swi = _swi(x, y)
- except (IndexError, RuntimeError):
- swi = -10
- row.append(swi)
- C.append(row)
-plot.pcolormesh(X, Y, C, cmap='viridis', vmin=-10, vmax=10, rasterized=True)
-for _, data in swis:
- plot.plot([x[0] for x in data], [x[1] for x in data])
-plot.show()