From accc352e3fdc9e264e6055466c4c55a6bd352400 Mon Sep 17 00:00:00 2001 From: Aki Date: Thu, 8 Sep 2022 18:01:35 +0200 Subject: Marked nomogram module as private --- szilagyi/_nomogram.py | 79 +++++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 79 insertions(+) create mode 100644 szilagyi/_nomogram.py (limited to 'szilagyi/_nomogram.py') diff --git a/szilagyi/_nomogram.py b/szilagyi/_nomogram.py new file mode 100644 index 0000000..4dcd141 --- /dev/null +++ b/szilagyi/_nomogram.py @@ -0,0 +1,79 @@ +import math +from collections import deque + +from . import _dataset + + +class NomogramEdgeCase(Exception): + def __init__(self, value): + self.value = value + + +def look_downwards(data, x, start): + for i in range(start, 0, -1): + if data[i - 1].x < 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].x > x: + break + else: + raise IndexError + return i + + +def find_segment(data, x): + width = data[-1].x - data[0].x + relative = x - data[0].x + candidate = math.floor(relative / width * len(data)) + look = look_downwards if data[candidate].x > 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].y > y and data[j].y > y: + segments.append((index, data, i, j)) + break + if data[i].y < y and data[j].y < y: + if segments: + segments.popleft() + segments.append((index, data, i, j)) + if len(segments) == 3: + middle = segments[1][1] + run = middle[j].x - middle[i].x + if run == 0: + raise NomogramEdgeCase(segments[1][0]) # Tidy up the dataset nonetheless + slope = (middle[j].y - middle[i].y) / run + intercept = middle[j].y - slope * middle[j].x + value = slope * x + intercept + if value == y: + raise NomogramEdgeCase(segments[1][0]) + if value < y: + segments.popleft() + else: + segments.pop() + if len(segments) == 1: + raise NomogramEdgeCase(-10) + return segments + + +def calculate_swi(x, y): + try: + low, high = find_boundary_curves(_dataset.INDICES, x, y) + vec = _dataset.Vector(x, y) + dist_to_low = min(abs(vec - p) for p in (low[1][low[2]], low[1][low[2]])) + dist_to_high = min(abs(vec - p) 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] + except NomogramEdgeCase as edge: + return edge.value + except IndexError: + return -10 # Reduce the amount of cases in which it may occur -- cgit v1.1