"""Draw the one-pole magnitude response at three corners.

One figure. Equation (5) of the page, evaluated on a logarithmic frequency
axis for the three coefficients equation (6) returns at corners of 250, 1000
and 4000 hertz. The three-decibel point marked on each curve is measured
from the curve by bisection, not taken from the requested corner, because
the two differ by a per cent or two once the corner climbs towards the
sample rate. That difference is the point the caption makes, so drawing the
requested corner would draw the wrong thing.

The colours become custom properties on the way out, so the drawing follows
the page between light and dark. See the note on `save`.

Writes corner.svg beside this file.
"""
import numpy as np
import matplotlib
matplotlib.use("Agg")
matplotlib.rcParams["svg.fonttype"] = "path"
import matplotlib.pyplot as plt

TEAL, INDIGO, ROSE = "#00706B", "#3A49A6", "#A93356"
INK, MUTED, GRID = "#0D1620", "#5A6874", "#D3DBE3"

# Every literal that reaches the file, and the theme property it becomes.
PROPERTY = {TEAL: "--teal", INDIGO: "--indigo", ROSE: "--rose",
            INK: "--ink", MUTED: "--muted", GRID: "--line"}

FS = 48000.0
CORNERS = [(250.0, INDIGO), (1000.0, TEAL), (4000.0, ROSE)]


def coefficient(fc):
    """Equation (6): the coefficient that puts the corner near fc."""
    return 1.0 - np.exp(-2.0 * np.pi * fc / FS)


def magnitude(c, f):
    """Equation (5), as a magnitude rather than a magnitude squared."""
    w = 2.0 * np.pi * f / FS
    return np.sqrt(c ** 2 / (1.0 - 2.0 * (1.0 - c) * np.cos(w) + (1.0 - c) ** 2))


def three_db(c):
    """The frequency where the curve has fallen to 1/sqrt(2), by bisection.

    The response falls monotonically from 1 at zero frequency, so a plain
    bisection on "still above the threshold" converges on the crossing.
    Sixty halvings take the bracket well below a hertz.
    """
    lo, hi = 1.0, FS / 2.0 - 1.0
    for _ in range(60):
        mid = 0.5 * (lo + hi)
        if magnitude(c, mid) > 1.0 / np.sqrt(2.0):
            lo = mid
        else:
            hi = mid
    return 0.5 * (lo + hi)


def save(fig, name):
    """Write the figure, then make it the page's own drawing.

    Matplotlib writes every colour into the SVG as a literal, so a figure
    saved the ordinary way carries fixed values and cannot follow the page
    between light and dark. Each literal becomes var(--prop, literal) here,
    which keeps the file readable on its own and lets the page recolour it.

    Ids are prefixed with the file's stem because several figures may be
    inlined into one page and Matplotlib numbers its own from one in each.
    """
    import re

    fig.tight_layout()
    # transparent, or Matplotlib paints an opaque white rectangle behind
    # everything and the figure reads as a light card on a dark page.
    fig.savefig(name, transparent=True, bbox_inches="tight", pad_inches=0.04)
    text = open(name, encoding="utf-8").read()

    stem = re.sub(r"[^A-Za-z0-9]+", "-", name.rsplit(".", 1)[0]).strip("-")
    text = re.sub(r'(\sid=")([^"]+)"', lambda m: f'{m.group(1)}{stem}-{m.group(2)}"', text)
    text = re.sub(r'((?:xlink:)?href="#)([^"]+)"',
                  lambda m: f'{m.group(1)}{stem}-{m.group(2)}"', text)
    text = re.sub(r'(url\(#)([^)]+)\)', lambda m: f'{m.group(1)}{stem}-{m.group(2)})', text)

    for literal, prop in PROPERTY.items():
        text = re.sub(re.escape(literal), f"var({prop}, {literal.lower()})",
                      text, flags=re.IGNORECASE)

    with open(name, "w", encoding="utf-8") as f:
        f.write(text)
    print(f"{name} written")


f = np.logspace(np.log10(20.0), np.log10(20000.0), 900)

fig, ax = plt.subplots(figsize=(6.4, 3.4))

ax.axhline(-3.0, lw=0.9, color=GRID, zorder=1)
ax.annotate("3 dB down", (21.0, -3.0), ha="left", va="bottom",
            fontsize=9, color=MUTED)

for fc, colour in CORNERS:
    c = coefficient(fc)
    ax.plot(f, 20.0 * np.log10(magnitude(c, f)), lw=1.8, color=colour, zorder=3)

    measured = three_db(c)
    ax.plot([measured], [-3.0], marker="o", ms=4.5, color=colour, zorder=4)
    ax.annotate(f"c = {c:.4f}", (measured * 1.12, -3.6), ha="left", va="top",
                fontsize=9.5, color=colour)
    print(f"corner {fc:.0f} Hz: c = {c:.4f}, measured -3 dB at {measured:.1f} Hz")

ax.set_xscale("log")
ax.set_xlim(20.0, 20000.0)
ax.set_ylim(-42.0, 4.0)
ax.set_xlabel("frequency (Hz)", fontsize=10, color=INK)
ax.set_ylabel("gain (dB)", fontsize=10, color=INK)
ax.set_xticks([20, 100, 1000, 10000])
ax.set_xticklabels(["20", "100", "1k", "10k"])
ax.set_yticks([0, -6, -12, -18, -24, -30, -36])
# which="both", or the minor ticks a logarithmic axis adds keep
# Matplotlib's default black and stay black on a dark page.
ax.tick_params(which="both", colors=MUTED, labelsize=9)
ax.grid(True, which="major", lw=0.6, color=GRID, alpha=0.7)
for side in ("top", "right"):
    ax.spines[side].set_visible(False)
for side in ("left", "bottom"):
    ax.spines[side].set_color(GRID)

save(fig, "corner.svg")
