Hands-on lab — Day 127: Why We Visualize, and Choosing the Right Chart

"""Chart choice as a function you can call, argue with, and test.

"Which chart should I use?" is normally answered by taste, which is why it
is normally answered badly. Two published results turn most of it into
arithmetic:

* Cleveland and McGill (1984) ranked visual channels by how accurately
  people read magnitudes off them. That ranking is `ENCODING_RANKING`.
* A variable's TYPE decides which channels can carry it honestly at all.
  A nominal variable on a length channel invents an order that is not in
  the data; an ordinal variable on a categorical hue palette destroys the
  order that is.

`best_encoding` combines the two: given a data type and the reader's task,
return the most accurately-judged channel that is still honest. Then
`choose_chart` answers the question one level up -- given the QUESTION,
how many values, and what types they are, name the chart.

These functions are opinionated on purpose. They are not the last word on
visualisation; they are a written-down default you can disagree with
explicitly, which is strictly better than an unwritten one you cannot.
"""

from __future__ import annotations

# --------------------------------------------------------------------------
# The Cleveland-McGill ordering
# --------------------------------------------------------------------------

# Ordered most accurately judged first. From Cleveland and McGill (1984),
# "Graphical Perception: Theory, Experimentation, and Application to the
# Development of Graphical Methods", Journal of the American Statistical
# Association 79(387):531-554. Their elementary perceptual tasks, with the
# names shortened for use as identifiers.
ENCODING_RANKING: tuple[str, ...] = (
    "position_common_scale",  # 1. dots or bar ends on one shared axis
    "position_nonaligned_scales",  # 2. the same, but across separate panels
    "length",  # 3. bar lengths not sharing a baseline
    "angle_slope",  # 4. pie slices, and line steepness
    "area",  # 5. bubble size, treemap tiles
    "volume",  # 6. 3-D bars, spheres
    "color_saturation",  # 7. heatmap intensity, density shading
)

# Channels that carry identity but NOT magnitude. Cleveland and McGill did
# not rank these, because there is no magnitude to read off them: hue is
# how you say "these two lines are different series", never "this one is
# 1.8 times that one".
IDENTITY_CHANNELS: tuple[str, ...] = ("hue", "shape")

# Sequential lightness -- a colormap like viridis, ordered by luminance.
# It carries ORDER faithfully and magnitude only roughly, which places it
# with saturation at the bottom of the accuracy ordering but well above
# categorical hue for anything ordered.
ORDERED_COLOR_CHANNEL = "luminance_sequential"

DATA_TYPES: frozenset[str] = frozenset({"nominal", "ordinal", "quantitative", "temporal"})


def encoding_rank(channel: str) -> int:
    """Position of `channel` in the Cleveland-McGill ordering, 0 = best."""
    try:
        return ENCODING_RANKING.index(channel)
    except ValueError:
        raise ValueError(
            f"{channel!r} is not a ranked magnitude channel; "
            f"ranked channels are {ENCODING_RANKING}"
        ) from None


# `best_encoding`'s task vocabulary. Each task says what the reader is
# trying to DO, because the same variable in the same chart deserves a
# different channel depending on the question being asked of it.
TASKS: frozenset[str] = frozenset(
    {
        "compare",  # read two magnitudes and say which is bigger, by how much
        "compare_across_panels",  # the same, but the values sit in separate small multiples
        "identify_group",  # tell which series a mark belongs to
        "encode_in_color",  # the two spatial axes are taken; colour is all that is left
        "magnitude_on_map",  # position is spent on geography, so magnitude needs another channel
        "trend",  # read direction and rate of change over an ordered axis
    }
)


def best_encoding(data_type: str, task: str) -> str:
    """Most accurately-judged channel that can honestly carry `data_type`.

    The rule in one sentence: take the highest-ranked channel from
    Cleveland and McGill that (a) the task has not already spent on
    something else, and (b) does not claim more structure than the data
    type has.
    """
    if data_type not in DATA_TYPES:
        raise ValueError(f"unknown data type {data_type!r}; expected one of {sorted(DATA_TYPES)}")
    if task not in TASKS:
        raise ValueError(f"unknown task {task!r}; expected one of {sorted(TASKS)}")

    if task == "identify_group":
        # Identity, not magnitude. A nominal variable has no order to
        # preserve and no size to read, so hue is exactly right -- and
        # putting it on a magnitude channel would invent an order.
        return "hue"

    if task == "encode_in_color":
        # Both spatial axes are already spent. Colour is the only channel
        # left, and the data type decides WHICH colour channel.
        if data_type == "nominal":
            return "hue"
        if data_type in ("ordinal", "temporal"):
            # Order must survive. Categorical hue would destroy it; a
            # luminance-ordered ramp preserves it.
            return ORDERED_COLOR_CHANNEL
        return "color_saturation"  # quantitative: magnitude, read roughly

    if task == "magnitude_on_map":
        # Position is spent on latitude and longitude. Of what remains,
        # area is the highest-ranked channel that survives being placed at
        # an arbitrary map location -- length would need a shared baseline
        # the map cannot provide.
        return "area"

    if task == "compare_across_panels":
        # Small multiples: each panel has its own axis, so this is
        # Cleveland and McGill's second task by construction.
        return "position_nonaligned_scales"

    if task == "trend":
        # Direction over an ordered axis. Slope is the thing being read,
        # but it is read off points placed on one common scale, and the
        # accuracy of the reading is the accuracy of those positions.
        return "position_common_scale"

    # task == "compare": one shared axis, the best channel there is. This
    # holds for every data type, including nominal -- a bar chart of
    # nominal categories puts the CATEGORY on the categorical axis and the
    # QUANTITY on the common scale, which invents no order at all.
    return "position_common_scale"


# --------------------------------------------------------------------------
# From the question to the chart
# --------------------------------------------------------------------------

QUESTION_KINDS: frozenset[str] = frozenset(
    {
        "comparison",  # which of these is biggest? by how much?
        "distribution",  # what shape is this variable? where is the mass?
        "relationship",  # does x move with y?
        "composition",  # what are the parts of this whole?
        "change_over_time",  # what happened, and in which direction?
    }
)

# Below this many values, a table is the better instrument. The number is
# a judgement, not a measurement, and the reasoning is written out in
# `choose_chart`'s docstring so you can move it deliberately.
TABLE_MAX_VALUES = 5

# Above this many points, individual marks stop being individually
# readable and the chart is showing density whether you meant it to or
# not. Better to say so and draw density directly.
OVERPLOT_POINT_LIMIT = 2000

# Above this many series, one panel becomes a tangle. Split it.
SMALL_MULTIPLE_LIMIT = 8


def choose_chart(question_kind: str, n_categories: int, data_types: list[str]) -> str:
    """Recommend a chart -- or a table -- for a question.

    `n_categories` is the number of values the reader must take in: bars in
    a bar chart, points in a scatter, lines in a line chart.

    Two recommendations in here are the point of the whole function.

    **It recommends a table below `TABLE_MAX_VALUES` values.** Three
    numbers do not need an axis, a legend and a title in order to be
    compared; they need to be readable. A chart's advantage is that it
    turns comparison into a perceptual judgement instead of an arithmetic
    one, and with three numbers there was never any arithmetic to save.
    The threshold is a judgement call, and 5 is where the trade tips in
    this course's experience -- put it somewhere else if your readers
    differ, but put it somewhere on purpose.

    **It never recommends a pie chart, for anything.** Reading a pie means
    judging angle and area, ranked fourth and fifth by Cleveland and
    McGill, when the identical data on a sorted bar chart would be
    position on a common scale, ranked first. The one case usually offered
    in a pie's defence -- two or three parts of a whole -- falls below
    `TABLE_MAX_VALUES` and gets a table, which answers the question
    exactly rather than approximately.
    """
    if question_kind not in QUESTION_KINDS:
        raise ValueError(
            f"unknown question kind {question_kind!r}; expected one of {sorted(QUESTION_KINDS)}"
        )
    if n_categories < 1:
        raise ValueError(f"n_categories must be at least 1, got {n_categories}")
    unknown = set(data_types) - DATA_TYPES
    if unknown:
        raise ValueError(f"unknown data types {sorted(unknown)}; expected from {sorted(DATA_TYPES)}")

    if question_kind == "change_over_time":
        # A time axis is not optional here: "over time" with no temporal
        # variable is a question about something else.
        if "temporal" not in data_types:
            raise ValueError(
                "change_over_time needs a temporal variable in data_types; "
                f"got {sorted(data_types)}"
            )
        return "small_multiples_line" if n_categories > SMALL_MULTIPLE_LIMIT else "line"

    if question_kind == "relationship":
        if "quantitative" not in data_types:
            raise ValueError(
                "relationship needs at least one quantitative variable; " f"got {sorted(data_types)}"
            )
        # Past the overplot limit the marks are stacked on each other and
        # the reader is looking at ink density, not at points.
        return "hexbin" if n_categories > OVERPLOT_POINT_LIMIT else "scatter"

    if question_kind == "distribution":
        if "quantitative" not in data_types:
            raise ValueError(
                "distribution needs a quantitative variable; " f"got {sorted(data_types)}"
            )
        if n_categories <= 1:
            return "histogram"
        if n_categories <= SMALL_MULTIPLE_LIMIT:
            return "small_multiples_histogram"
        return "boxplot_by_category"

    # comparison and composition share the table rule.
    if n_categories <= TABLE_MAX_VALUES:
        return "table"

    if question_kind == "composition":
        return "stacked_bar"

    return "sorted_horizontal_bar"


# --------------------------------------------------------------------------
# Sorting, measured as reader effort
# --------------------------------------------------------------------------


def comparisons_to_find_max(values: list[float], presented_sorted: bool) -> int:
    """How many comparisons a reader makes to find the largest value.

    This simulates a reader, not a computer. Two behaviours, both real:

    * **Source order.** The reader has no idea where the biggest bar is,
      so they hold a running best and check every remaining bar against
      it: `n - 1` comparisons for `n` bars.
    * **Sorted, descending.** The reader looks at the top row. They still
      make ONE comparison -- top row against the next one -- to confirm
      the chart really is sorted and the leader really is ahead. After
      that they stop, because sorting is a promise about everything below.

    The answer is identical either way. The effort is not, and sorting is
    the cheapest thing you will ever do to a chart.
    """
    n = len(values)
    if n < 2:
        raise ValueError("finding a maximum needs at least two values")
    return 1 if presented_sorted else n - 1


def index_of_max(values: list[float]) -> int:
    """Index of the largest value, for confirming sorting changes only effort."""
    return max(range(len(values)), key=lambda i: values[i])

"""Shared fixtures. pytest finds this file by itself -- nothing imports it.

`points` returns the same seeded cloud on every test, so two tests that
render it are rendering the identical data.

`png_dir` is a temporary directory OUTSIDE the lab, created and removed by
the fixture itself. Every render in this suite writes there and nowhere
else, which is why `tests/run_tests.sh` can assert afterwards that the lab
directory contains no `.png` at all.
"""

from __future__ import annotations

import tempfile
from pathlib import Path

import pytest

import render as R


@pytest.fixture
def points():
    return R.sample_points()


@pytest.fixture
def png_dir():
    with tempfile.TemporaryDirectory(prefix="d127-render-") as d:
        yield Path(d)

"""The measurable half of visual encoding: geometry and colour arithmetic.

Nothing in this module renders anything. Everything here is a pure
function over numbers, which is exactly why it can be asserted on. The
rendering lives in `render.py`.

Two ideas carry the whole file.

1. **Area is quadratic in radius.** If you draw a value as a circle's
   RADIUS, doubling the value quadruples the ink. The reader judges the
   ink, so every ratio in the chart comes out squared. `area_for_radius`
   and the two `radii_*` functions make that concrete enough to assert on.

2. **A colour pair that is far apart for you may be on top of each other
   for a reader with a colour vision deficiency.** `simulate_deuteranopia`
   applies a published transformation matrix; `delta_e_cie76` measures how
   far apart two colours are afterwards. Neither function has an opinion --
   they return numbers, and the tests read them.

Colour maths conventions used throughout:

* An sRGB colour is a 3-tuple of floats in [0, 1] -- the same convention
  matplotlib and seaborn hand you.
* "Linear RGB" is sRGB with its transfer function undone. The colour
  vision deficiency matrix is defined on LINEAR RGB, not on the gamma
  encoded values, so `simulate_deuteranopia` linearises first and
  re-encodes afterwards. Skipping that step is the single most common
  error in hand-rolled deficiency simulators and it visibly changes the
  answer.
* CIELAB is computed against the D65 white point, matching sRGB's own
  reference white.
"""

from __future__ import annotations

import math

# --------------------------------------------------------------------------
# 1. Size encoding -- the square law
# --------------------------------------------------------------------------


def area_for_radius(radius: float) -> float:
    """Area of a circle of this radius. The whole square law in one line."""
    if radius < 0:
        raise ValueError(f"radius must be non-negative, got {radius}")
    return math.pi * radius * radius


def radii_scaled_by_radius(values: list[float], unit_radius: float = 1.0) -> list[float]:
    """Encode each value as a circle whose RADIUS is proportional to it.

    This is what you get from the obvious, wrong implementation --
    `plt.scatter(x, y, s=value)` reads `s` as an area, but a hand-written
    `Circle(xy, radius=value)` or a d3 `.attr("r", value)` does exactly
    this. A value twice as large is drawn with twice the radius and
    therefore FOUR times the area.
    """
    return [unit_radius * float(v) for v in values]


def radii_scaled_by_area(values: list[float], unit_radius: float = 1.0) -> list[float]:
    """Encode each value as a circle whose AREA is proportional to it.

    Take the square root and the distortion disappears: a value twice as
    large is drawn with twice the ink, so the ratio the reader perceives
    is the ratio in the data.
    """
    out = []
    for v in values:
        v = float(v)
        if v < 0:
            raise ValueError(f"cannot encode a negative value as an area: {v}")
        out.append(unit_radius * math.sqrt(v))
    return out


def encoded_area_ratio(values: list[float], mode: str) -> float:
    """Ratio of drawn AREA between the last and first value, under `mode`.

    `mode` is "radius" (the distorting encoding) or "area" (the honest
    one). Returns a plain float so a test can compare it against the
    ratio in the data itself.
    """
    if mode == "radius":
        radii = radii_scaled_by_radius(values)
    elif mode == "area":
        radii = radii_scaled_by_area(values)
    else:
        raise ValueError(f"mode must be 'radius' or 'area', got {mode!r}")
    first, last = area_for_radius(radii[0]), area_for_radius(radii[-1])
    if first == 0:
        raise ValueError("the first value encodes to zero area; ratio is undefined")
    return last / first


# --------------------------------------------------------------------------
# 2. Colour spaces
# --------------------------------------------------------------------------


def _srgb_to_linear_channel(c: float) -> float:
    """Undo sRGB's transfer function for one channel in [0, 1]."""
    return c / 12.92 if c <= 0.04045 else ((c + 0.055) / 1.055) ** 2.4


def _linear_to_srgb_channel(c: float) -> float:
    """Re-apply sRGB's transfer function for one channel, then clamp."""
    c = max(0.0, min(1.0, c))
    v = 12.92 * c if c <= 0.0031308 else 1.055 * (c ** (1 / 2.4)) - 0.055
    return max(0.0, min(1.0, v))


def srgb_to_linear(rgb: tuple[float, float, float]) -> tuple[float, float, float]:
    return tuple(_srgb_to_linear_channel(float(c)) for c in rgb)  # type: ignore[return-value]


def linear_to_srgb(rgb: tuple[float, float, float]) -> tuple[float, float, float]:
    return tuple(_linear_to_srgb_channel(float(c)) for c in rgb)  # type: ignore[return-value]


def relative_luminance(rgb: tuple[float, float, float]) -> float:
    """WCAG relative luminance: the perceived lightness of a colour.

    This is the quantity a greyscale photocopy of your chart keeps and
    everything else throws away. A palette that carries order in hue but
    not in luminance is a palette whose order vanishes in print.
    """
    r, g, b = srgb_to_linear(rgb)
    return 0.2126 * r + 0.7152 * g + 0.0722 * b


# D65 white point in XYZ, scaled so Y = 1. sRGB's own reference white.
_D65 = (0.95047, 1.00000, 1.08883)


def srgb_to_xyz(rgb: tuple[float, float, float]) -> tuple[float, float, float]:
    r, g, b = srgb_to_linear(rgb)
    x = 0.4124564 * r + 0.3575761 * g + 0.1804375 * b
    y = 0.2126729 * r + 0.7151522 * g + 0.0721750 * b
    z = 0.0193339 * r + 0.1191920 * g + 0.9503041 * b
    return (x, y, z)


def _lab_f(t: float) -> float:
    delta = 6 / 29
    return t ** (1 / 3) if t > delta**3 else t / (3 * delta**2) + 4 / 29


def srgb_to_lab(rgb: tuple[float, float, float]) -> tuple[float, float, float]:
    """Convert sRGB to CIELAB (L*, a*, b*) against D65.

    CIELAB exists because RGB distance is a poor model of perceived
    difference: two RGB triples the same Euclidean distance apart can look
    obviously different or nearly identical. CIELAB was built so that
    Euclidean distance is at least roughly proportional to perceived
    difference, which is what makes `delta_e_cie76` meaningful at all.
    """
    x, y, z = srgb_to_xyz(rgb)
    fx, fy, fz = _lab_f(x / _D65[0]), _lab_f(y / _D65[1]), _lab_f(z / _D65[2])
    return (116 * fy - 16, 500 * (fx - fy), 200 * (fy - fz))


def delta_e_cie76(a: tuple[float, float, float], b: tuple[float, float, float]) -> float:
    """CIE76 colour difference: Euclidean distance in CIELAB.

    The 1976 formula is the simplest of the family and is used here
    deliberately: it is short enough to read, has no tuning constants to
    get wrong, and the tests only ever ask whether a distance is large or
    small, never for a precise perceptual match. Rough guide, from the
    literature: about 2.3 is the "just noticeable difference" for adjacent
    patches; a difference in the low single digits is two colours a reader
    will struggle to tell apart at all in a legend.
    """
    la, aa, ba = srgb_to_lab(a)
    lb, ab, bb = srgb_to_lab(b)
    return math.sqrt((la - lb) ** 2 + (aa - ab) ** 2 + (ba - bb) ** 2)


# --------------------------------------------------------------------------
# 3. Colour vision deficiency simulation
# --------------------------------------------------------------------------

# Machado, Oliveira and Fernandes (2009), "A Physiologically-based Model
# for Simulation of Color Vision Deficiency", IEEE Transactions on
# Visualization and Computer Graphics 15(6):1291-1298. This is the
# deuteranomaly matrix at severity 1.0 -- that is, deuteranopia, the
# complete absence of a working M cone. The nine coefficients were read
# from the authors' published matrix table, not reconstructed from memory.
#
# The matrix operates on LINEAR RGB. `simulate_deuteranopia` linearises
# before multiplying and re-encodes afterwards.
DEUTERANOPIA_MATRIX_MACHADO_2009 = (
    (0.367322, 0.860646, -0.227968),
    (0.280085, 0.672501, 0.047413),
    (-0.011820, 0.042940, 0.968881),
)


def simulate_deuteranopia(
    rgb: tuple[float, float, float],
    matrix: tuple[tuple[float, float, float], ...] = DEUTERANOPIA_MATRIX_MACHADO_2009,
) -> tuple[float, float, float]:
    """Return an sRGB approximation of how `rgb` appears to a deuteranope.

    A very important limit, stated here rather than buried: this
    APPROXIMATES a deficiency, it does not reproduce an experience. The
    model is a linear transform fitted to a physiological model of cone
    response; it says nothing about how a particular person has learned to
    compensate over a lifetime, it assumes a single severity, and it
    cannot represent the many varieties of anomalous trichromacy at all.
    Treat a small simulated distance as strong evidence that a palette is
    risky, and a large one as weak evidence that it is fine. The reliable
    move is still to carry the distinction in a second channel -- shape,
    position, direct labelling -- so colour is never load-bearing alone.
    """
    r, g, b = srgb_to_linear(rgb)
    out = tuple(m[0] * r + m[1] * g + m[2] * b for m in matrix)
    return linear_to_srgb(out)  # type: ignore[arg-type]


def deuteranopia_collapse(
    a: tuple[float, float, float], b: tuple[float, float, float]
) -> dict[str, float]:
    """Measure how much of a pair's separation survives deuteranopia.

    Returns the normal-vision distance, the simulated distance, and the
    fraction retained. A pair whose separation nearly vanishes is a pair
    that must not be the only thing distinguishing two series.
    """
    normal = delta_e_cie76(a, b)
    simulated = delta_e_cie76(simulate_deuteranopia(a), simulate_deuteranopia(b))
    return {
        "normal_delta_e": normal,
        "simulated_delta_e": simulated,
        "retained_fraction": simulated / normal if normal else 0.0,
    }


# --------------------------------------------------------------------------
# 4. Rank correlation, without scipy
# --------------------------------------------------------------------------


def _ranks(values: list[float]) -> list[float]:
    """Ranks, averaging ties -- the standard midrank treatment."""
    order = sorted(range(len(values)), key=lambda i: values[i])
    ranks = [0.0] * len(values)
    i = 0
    while i < len(order):
        j = i
        while j + 1 < len(order) and values[order[j + 1]] == values[order[i]]:
            j += 1
        midrank = (i + j) / 2 + 1
        for k in range(i, j + 1):
            ranks[order[k]] = midrank
        i = j + 1
    return ranks


def spearman_rho(xs: list[float], ys: list[float]) -> float:
    """Spearman's rank correlation, computed from Pearson on the ranks.

    scipy is not installed in this environment, so this is the honest
    twenty-line version rather than a call to `scipy.stats.spearmanr`.
    Written this way it also handles ties correctly, which the shortcut
    "1 - 6*sum(d^2)/(n^3-n)" formula quietly does not.
    """
    if len(xs) != len(ys):
        raise ValueError("spearman_rho needs two sequences of the same length")
    if len(xs) < 2:
        raise ValueError("spearman_rho needs at least two observations")
    rx, ry = _ranks(list(xs)), _ranks(list(ys))
    mx, my = sum(rx) / len(rx), sum(ry) / len(ry)
    num = sum((a - mx) * (b - my) for a, b in zip(rx, ry))
    dx = math.sqrt(sum((a - mx) ** 2 for a in rx))
    dy = math.sqrt(sum((b - my) ** 2 for b in ry))
    if dx == 0 or dy == 0:
        raise ValueError("a sequence with no variation has no rank correlation")
    return num / (dx * dy)


def luminance_order_correlation(palette: list[tuple[float, float, float]]) -> float:
    """Rank correlation between a palette's position and its luminance.

    A SEQUENTIAL palette is built so that step 1 is darker than step 2 is
    darker than step 3; its correlation is +/-1 exactly. A CATEGORICAL
    palette is built so that neighbouring swatches are as DIFFERENT as
    possible, which says nothing at all about their lightness order --
    so mapping an ordered variable onto one destroys the order.
    """
    positions = [float(i) for i in range(len(palette))]
    luminances = [relative_luminance(c) for c in palette]
    return spearman_rho(positions, luminances)

"""The palettes this lab measures, taken from the libraries themselves.

Nothing here is a hand-picked colour chosen to make a point. Every swatch
is what matplotlib or seaborn hands you by default, which is what makes
the measurements in `test_charts.py` about real tools rather than about a
straw man.

* `PASS_FAIL_RED` and `PASS_FAIL_GREEN` are `tab10`'s red and green --
  matplotlib's default cycle, positions 3 and 2. The reflex "red means
  failed, green means passed" chart is drawn in exactly these two colours
  by anyone who does not stop to think about it.
* `SAFE_BLUE` and `SAFE_ORANGE` are the first two entries of seaborn's
  `colorblind` palette, which is that library's own answer to the same
  problem.
* `viridis_steps` samples matplotlib's default SEQUENTIAL colormap, which
  is built to increase monotonically in luminance.
* `tab10_steps` samples its default CATEGORICAL palette, which is built so
  neighbouring swatches look as different as possible -- and therefore
  carries no luminance order at all.

Both `viridis` and `tab10` ship inside matplotlib. Nothing here touches
the network.
"""

from __future__ import annotations

import matplotlib

matplotlib.use("Agg")

import matplotlib.pyplot as plt  # noqa: E402
import seaborn as sns  # noqa: E402

Color = tuple[float, float, float]


def _rgb(c) -> Color:
    return (float(c[0]), float(c[1]), float(c[2]))


# matplotlib's default categorical cycle, positions 3 (red) and 2 (green).
PASS_FAIL_RED: Color = _rgb(plt.get_cmap("tab10")(3))
PASS_FAIL_GREEN: Color = _rgb(plt.get_cmap("tab10")(2))

# seaborn's own colourblind-safe categorical palette, first two entries.
_CB = sns.color_palette("colorblind")
SAFE_BLUE: Color = _rgb(_CB[0])
SAFE_ORANGE: Color = _rgb(_CB[1])


def viridis_steps(n: int = 5) -> list[Color]:
    """`n` evenly spaced samples of matplotlib's sequential default."""
    cmap = plt.get_cmap("viridis")
    return [_rgb(cmap(i / (n - 1))) for i in range(n)]


def tab10_steps(n: int = 5) -> list[Color]:
    """The first `n` entries of matplotlib's categorical default."""
    cmap = plt.get_cmap("tab10")
    return [_rgb(cmap(i)) for i in range(n)]


# An ordered variable with five levels. Its ORDER is the whole point: a
# palette that does not preserve it has destroyed information the data
# had, which is a different and worse failure than merely looking bad.
SATISFACTION_LEVELS: tuple[str, ...] = (
    "very dissatisfied",
    "dissatisfied",
    "neutral",
    "satisfied",
    "very satisfied",
)

"""Rendering, and measuring what was actually rendered.

"That chart looks better" is not testable. "That chart spends 37% of its
ink on the data where this one spends 93%, and both carry the same eight
numbers" is. Everything in this module exists to turn a claim about a
picture into a number a test can assert on.

Every function here renders through matplotlib's **Agg** backend, which
draws into a memory buffer and needs no display, no window server and no
X11 forwarding. `matplotlib.use("Agg")` is called BEFORE `pyplot` is
imported, because the backend is chosen at import time and switching
afterwards is unreliable. `plt.show()` is never called, and every figure
is closed on the way out so a long test run does not leak them.

Every function takes an explicit output path. Nothing here writes to a
default location, so a test that hands it a `tmp_path` leaves nothing on
disk once pytest cleans up.
"""

from __future__ import annotations

from pathlib import Path

import matplotlib

matplotlib.use("Agg")  # must precede the pyplot import

import matplotlib.pyplot as plt  # noqa: E402
import numpy as np  # noqa: E402
from matplotlib.patches import Circle  # noqa: E402
from PIL import Image  # noqa: E402

WHITE = (255, 255, 255)

# The flat colour the bars are drawn in, as an sRGB 8-bit triple, so a
# test can isolate the data ink from the furniture.
BAR_RGB: tuple[int, int, int] = (0x1D, 0x4E, 0xD8)
BAR_HEX = "#%02x%02x%02x" % BAR_RGB

# The eight-region growth figures the lesson opens with. Deliberately
# close together at the top: 18.9 against 17.4 is a gap a sorted bar chart
# shows instantly and a pie chart cannot resolve at all.
REGION_GROWTH: dict[str, float] = {
    "Nordics": 12.1,
    "Iberia": 18.9,
    "Benelux": 7.4,
    "DACH": 17.4,
    "France": 9.8,
    "Italy": 4.2,
    "Poland": 15.6,
    "Ireland": 11.3,
}


# --------------------------------------------------------------------------
# Pixel measurement
# --------------------------------------------------------------------------


def _load_rgb(path: Path) -> np.ndarray:
    """Load a PNG as an (H, W, 3) uint8 array, discarding any alpha."""
    with Image.open(path) as im:
        return np.asarray(im.convert("RGB"))


def count_non_background_pixels(path: Path, background: tuple[int, int, int] = WHITE) -> int:
    """Count pixels that are not the background colour. This is the ink."""
    arr = _load_rgb(path)
    bg = np.array(background, dtype=np.uint8)
    return int(np.count_nonzero(np.any(arr != bg, axis=-1)))


def count_pixels_of_color(path: Path, rgb: tuple[int, int, int]) -> int:
    """Count pixels exactly matching one colour.

    The bars in `render_region_bar_chart` are drawn in one flat colour and
    nothing else in either figure uses it, so this isolates the DATA ink
    from every other mark on the page -- which is what makes Tufte's
    data-ink ratio a number this lab can measure rather than assert.
    """
    arr = _load_rgb(path)
    target = np.array(rgb, dtype=np.uint8)
    return int(np.count_nonzero(np.all(arr == target, axis=-1)))


def data_ink_ratio(path: Path, data_rgb: tuple[int, int, int]) -> float:
    """Data ink divided by total ink, both counted in pixels.

    Tufte's definition, made literal: of every mark on the page, what
    fraction is the data itself? Erasing non-data ink -- gridlines, a
    tinted panel, a heavy box -- raises this without removing a single
    fact from the chart.
    """
    total = count_non_background_pixels(path)
    if total == 0:
        raise ValueError("the image has no ink at all; a data-ink ratio is undefined")
    return count_pixels_of_color(path, data_rgb) / total


def count_distinct_luminance_levels(path: Path) -> int:
    """Count how many distinct grey levels the image contains.

    This is the measurable trace of density information. Opaque marks that
    overlap produce exactly two levels -- paper and ink -- no matter how
    many marks landed on the same pixel, because the tenth mark on a pixel
    changes nothing. Semi-transparent marks accumulate, so the number of
    levels tells you the image is carrying how MANY marks landed, not just
    whether any did.
    """
    arr = _load_rgb(path).astype(np.float64)
    lum = 0.2126 * arr[..., 0] + 0.7152 * arr[..., 1] + 0.0722 * arr[..., 2]
    return int(np.unique(np.round(lum).astype(np.int64)).size)


# --------------------------------------------------------------------------
# Exercise 1 -- circles drawn at a known radius, then measured
# --------------------------------------------------------------------------


def render_circle(path: Path, radius_px: float, canvas_px: int = 400) -> None:
    """Draw one filled circle of exactly `radius_px` pixels on white.

    The axes fill the whole figure and the data limits are set equal to
    the pixel dimensions, so one data unit is one pixel exactly and the
    circle's radius in data units is its radius on screen. Antialiasing is
    switched off so the edge is a hard boundary and the pixel count is a
    measurement of area rather than of area plus a soft fringe.
    """
    dpi = 100
    fig = plt.figure(figsize=(canvas_px / dpi, canvas_px / dpi), dpi=dpi, facecolor="white")
    ax = fig.add_axes((0, 0, 1, 1))
    ax.set_xlim(0, canvas_px)
    ax.set_ylim(0, canvas_px)
    ax.set_axis_off()
    ax.add_patch(
        Circle(
            (canvas_px / 2, canvas_px / 2),
            radius_px,
            facecolor="black",
            edgecolor="none",
            antialiased=False,
        )
    )
    fig.savefig(path, dpi=dpi, facecolor="white")
    plt.close(fig)


def measure_circle_area_px(path: Path) -> int:
    """The circle's drawn area, in painted pixels."""
    return count_non_background_pixels(path)


# --------------------------------------------------------------------------
# Exercise 7 -- the same chart, decorated and undecorated
# --------------------------------------------------------------------------


def render_region_bar_chart(path: Path, decorated: bool) -> None:
    """Render the eight-region bar chart with or without the furniture.

    Identical data, identical figure size, identical bars. The only
    difference is the non-data ink: a tinted plot background, a full box
    of spines, gridlines on both axes, and a heavy frame. Everything the
    reader needs -- the eight labels, the eight lengths, a common baseline
    -- is present in BOTH.
    """
    dpi = 100
    fig = plt.figure(figsize=(6, 4), dpi=dpi, facecolor="white")
    ax = fig.add_subplot(111)

    names = list(REGION_GROWTH)
    values = [REGION_GROWTH[n] for n in names]
    order = sorted(range(len(values)), key=lambda i: values[i])
    names = [names[i] for i in order]
    values = [values[i] for i in order]

    ax.barh(names, values, color=BAR_HEX)
    ax.set_xlabel("growth %")

    if decorated:
        ax.set_facecolor("#e2e8f0")
        ax.grid(True, which="both", axis="both", color="#64748b", linewidth=1.0)
        ax.set_axisbelow(False)
        for spine in ax.spines.values():
            spine.set_visible(True)
            spine.set_linewidth(3.0)
            spine.set_color("#334155")
    else:
        ax.set_facecolor("white")
        ax.grid(False)
        for name, spine in ax.spines.items():
            spine.set_visible(name == "bottom")
        ax.tick_params(left=False)

    fig.tight_layout()
    fig.savefig(path, dpi=dpi, facecolor="white")
    plt.close(fig)


# --------------------------------------------------------------------------
# Exercise 8 -- overplotting, and two ways out of it
# --------------------------------------------------------------------------


def sample_points(n: int = 10_000, seed: int = 127) -> tuple[np.ndarray, np.ndarray]:
    """A fixed, seeded cloud of `n` correlated points.

    Seeded, so every machine plots the same cloud and every pixel count in
    this lab is reproducible rather than merely typical.
    """
    rng = np.random.default_rng(seed)
    x = rng.normal(0.0, 1.0, n)
    y = 0.6 * x + rng.normal(0.0, 0.8, n)
    return x, y


# The scatter canvas. Both limits are wide enough that every one of the
# 10,000 sampled points falls inside the axes -- nothing is clipped, so
# the painted-pixel count can be compared against the full point count
# without an "and some were off the edge" caveat.
SCATTER_INCHES = 3
SCATTER_LIMIT = 5


def render_scatter(path: Path, x: np.ndarray, y: np.ndarray, alpha: float) -> None:
    """Plot the cloud with one-pixel marks at the given opacity.

    The `,` marker is matplotlib's single-pixel marker and antialiasing is
    off, so each point paints exactly one pixel. That makes the painted
    pixel count a direct measurement of how many DISTINCT screen positions
    the data occupies, with no marker-size confound: any shortfall against
    the point count is overplotting and nothing else.
    """
    dpi = 100
    fig = plt.figure(figsize=(SCATTER_INCHES, SCATTER_INCHES), dpi=dpi, facecolor="white")
    ax = fig.add_axes((0, 0, 1, 1))
    ax.set_xlim(-SCATTER_LIMIT, SCATTER_LIMIT)
    ax.set_ylim(-SCATTER_LIMIT, SCATTER_LIMIT)
    ax.set_axis_off()
    ax.plot(x, y, ",", color="black", alpha=alpha, antialiased=False, linestyle="none")
    fig.savefig(path, dpi=dpi, facecolor="white")
    plt.close(fig)


def render_hexbin(path: Path, x: np.ndarray, y: np.ndarray, gridsize: int = 30) -> None:
    """Plot the same cloud as a hexagonal density map.

    Where the scatter throws density away by painting the same pixel over
    and over, hexbin counts the points per cell and encodes the count as
    luminance -- the information the scatter destroyed, recovered by
    aggregating before drawing instead of after.
    """
    dpi = 100
    lim = SCATTER_LIMIT
    fig = plt.figure(figsize=(SCATTER_INCHES, SCATTER_INCHES), dpi=dpi, facecolor="white")
    ax = fig.add_axes((0, 0, 1, 1))
    ax.set_xlim(-lim, lim)
    ax.set_ylim(-lim, lim)
    ax.set_axis_off()
    ax.hexbin(x, y, gridsize=gridsize, cmap="Greys", extent=(-lim, lim, -lim, lim))
    fig.savefig(path, dpi=dpi, facecolor="white")
    plt.close(fig)


def count_painted_pixels(path: Path) -> int:
    """Distinct screen positions carrying at least one mark."""
    return count_non_background_pixels(path)


def points_inside_axes(x: np.ndarray, y: np.ndarray, limit: float = SCATTER_LIMIT) -> int:
    """How many of the sampled points fall inside the scatter's axes.

    Used to confirm that nothing is clipped, so a painted-pixel count
    below the point count means overplotting and never "the rest fell off
    the edge of the picture".
    """
    return int(np.count_nonzero((np.abs(x) <= limit) & (np.abs(y) <= limit)))

"""The worked reference suite for Day 127 -- "Charts That Answer the Question".

Nine exercises. The hard part of a visualisation lab is that "looks
better" is not testable, so this suite never asserts it. It asserts only
things that genuinely are measurable: drawn pixel areas, colour distances
in CIELAB before and after a colour vision deficiency transform, rank
correlation between a palette's order and its luminance, the number of
comparisons a reader performs, the fraction of a chart's ink that is data,
and the number of distinct grey levels an image contains.

Run it:

    pytest examples

`encoding.py` holds the geometry and colour arithmetic, `charts.py` the
two decision functions, `palettes.py` the swatches taken from matplotlib
and seaborn, and `render.py` everything that draws and everything that
measures a drawing. Read `starter/00_brief.md` for the exercise-by-exercise
explanation; this file is the answer key.

Every render writes into `png_dir`, a temporary directory outside the lab
that the fixture removes on the way out. The lab leaves no image behind.
"""

from __future__ import annotations

import math

import pytest

import charts as C
import encoding as E
import palettes as PAL
import render as R

# --------------------------------------------------------------------------
# Exercise 1 -- the square law, measured.
#
# Encode a value as a circle's RADIUS and every ratio in the chart comes
# out squared: a value twice as large is drawn four times as large. Encode
# it as the circle's AREA and the ratio the reader perceives is the ratio
# in the data. First the arithmetic, then the same claim measured off real
# rendered pixels, because arithmetic about a picture is not a picture.
# --------------------------------------------------------------------------

VALUES = [50.0, 100.0]  # the second is exactly twice the first


def test_1_radius_encoding_squares_every_ratio():
    data_ratio = VALUES[-1] / VALUES[0]
    assert data_ratio == 2.0

    by_radius = E.encoded_area_ratio(VALUES, mode="radius")
    by_area = E.encoded_area_ratio(VALUES, mode="area")

    # The distortion is exactly the square of the data ratio.
    assert by_radius == pytest.approx(4.0)
    assert by_radius == pytest.approx(data_ratio**2)

    # Encoding by area removes it entirely.
    assert by_area == pytest.approx(2.0)
    assert by_area == pytest.approx(data_ratio)

    # And the distortion is a factor of two, not a rounding difference.
    assert by_radius / by_area == pytest.approx(2.0)


def test_1_rendered_pixel_areas_confirm_the_square_law(png_dir):
    # Radius encoding: value 50 -> radius 40 px, value 100 -> radius 80 px.
    small = png_dir / "r_small.png"
    big = png_dir / "r_big.png"
    R.render_circle(small, radius_px=40)
    R.render_circle(big, radius_px=80)
    area_small = R.measure_circle_area_px(small)
    area_big = R.measure_circle_area_px(big)

    # Measured on this machine: 5156 px and 20368 px. The rasteriser lands
    # a couple of percent off the ideal pi*r^2 (5026.5 and 20106.2) because
    # a circle's boundary does not fall on pixel edges; the RATIO is what
    # the exercise is about and it survives that intact.
    assert area_small == pytest.approx(math.pi * 40**2, rel=0.03)
    assert area_big == pytest.approx(math.pi * 80**2, rel=0.03)
    assert area_big / area_small == pytest.approx(4.0, rel=0.02)

    # Area encoding: value 100 gets radius 40*sqrt(2), so twice the ink.
    honest = png_dir / "a_big.png"
    R.render_circle(honest, radius_px=40 * math.sqrt(2))
    area_honest = R.measure_circle_area_px(honest)
    assert area_honest / area_small == pytest.approx(2.0, rel=0.02)

    # The two encodings of the SAME value differ by a factor of two on the
    # page. This is the bubble chart that exaggerates without lying about
    # a single number.
    assert area_big / area_honest == pytest.approx(2.0, rel=0.02)


# --------------------------------------------------------------------------
# Exercise 2 -- the perceptual ranking as a decision function.
#
# Cleveland and McGill measured how accurately people read magnitudes off
# each visual channel. That ranking is not taste, and `best_encoding` is
# what it looks like when you actually use it to decide something.
# --------------------------------------------------------------------------


def test_2_ranking_is_in_cleveland_mcgill_order():
    # Position beats length beats angle beats area beats volume beats
    # saturation. Asserting the INDICES rather than the list literal means
    # the test is about the ordering, not about the spelling.
    assert C.encoding_rank("position_common_scale") == 0
    assert C.encoding_rank("position_common_scale") < C.encoding_rank("position_nonaligned_scales")
    assert C.encoding_rank("position_nonaligned_scales") < C.encoding_rank("length")
    assert C.encoding_rank("length") < C.encoding_rank("angle_slope")
    assert C.encoding_rank("angle_slope") < C.encoding_rank("area")
    assert C.encoding_rank("area") < C.encoding_rank("volume")
    assert C.encoding_rank("volume") < C.encoding_rank("color_saturation")

    # Hue is not on the ladder at all, and asking for its rank is an error
    # rather than a large number: there is no magnitude to read off a hue,
    # so "how accurately" is not a question hue can be asked.
    with pytest.raises(ValueError, match="not a ranked magnitude channel"):
        C.encoding_rank("hue")


def test_2_best_encoding_case_table():
    # Each case, and why it is the answer it is.
    cases = [
        # Nothing is competing for the axes, so take the top of the ladder.
        (("quantitative", "compare"), "position_common_scale"),
        # Small multiples: each panel carries its own axis, which IS
        # Cleveland and McGill's second task.
        (("quantitative", "compare_across_panels"), "position_nonaligned_scales"),
        # A map has already spent both spatial axes on geography, so the
        # best REMAINING channel is area -- bubbles, and now you know
        # exactly why they must be scaled by area and not by radius.
        (("quantitative", "magnitude_on_map"), "area"),
        # Both axes taken by other variables: colour is what is left, and
        # for a quantity that means intensity, read roughly.
        (("quantitative", "encode_in_color"), "color_saturation"),
        # An ordinal variable on an axis reads exactly like a quantitative
        # one: the categories are in order and position preserves it.
        (("ordinal", "compare"), "position_common_scale"),
        # But an ordinal variable pushed into colour must keep its order,
        # so it needs a luminance ramp -- never a categorical palette.
        (("ordinal", "encode_in_color"), "luminance_sequential"),
        # Time is ordered, so the same rule applies to it.
        (("temporal", "encode_in_color"), "luminance_sequential"),
        (("temporal", "trend"), "position_common_scale"),
        # Nominal has no order and no magnitude. Hue says "different",
        # which is the entire claim nominal data supports.
        (("nominal", "identify_group"), "hue"),
        (("nominal", "encode_in_color"), "hue"),
        # Comparing counts BY nominal category still puts the count on the
        # common scale; the category goes on the categorical axis and no
        # order is invented.
        (("nominal", "compare"), "position_common_scale"),
    ]
    for (data_type, task), expected in cases:
        assert C.best_encoding(data_type, task) == expected, (data_type, task)

    # The load-bearing negative: ordinal data must NOT land on a
    # categorical hue palette. Exercise 5 measures why.
    assert C.best_encoding("ordinal", "encode_in_color") != "hue"


def test_2_best_encoding_rejects_what_it_does_not_understand():
    with pytest.raises(ValueError, match="unknown data type"):
        C.best_encoding("categorical", "compare")
    with pytest.raises(ValueError, match="unknown task"):
        C.best_encoding("quantitative", "look_nice")


# --------------------------------------------------------------------------
# Exercise 3 -- from the question to the chart.
#
# `choose_chart` takes the reader's question, the number of values, and
# the data types, and names an instrument. Two of its answers are the
# point of the whole function: below a stated number of values it returns
# a TABLE, and it never returns a pie chart for anything.
# --------------------------------------------------------------------------


def test_3_choose_chart_case_table():
    cases = [
        # Comparison. Eight regions is past the table threshold, and
        # sorting is what turns "find the largest" into one glance.
        (("comparison", 8, ["nominal", "quantitative"]), "sorted_horizontal_bar"),
        (("comparison", 30, ["nominal", "quantitative"]), "sorted_horizontal_bar"),
        # Three numbers do not need a chart. Print them.
        (("comparison", 3, ["nominal", "quantitative"]), "table"),
        (("comparison", 5, ["nominal", "quantitative"]), "table"),
        (("comparison", 6, ["nominal", "quantitative"]), "sorted_horizontal_bar"),
        # Distribution. One variable is a histogram; a handful of groups
        # is small multiples; many groups is a box plot grid, where the
        # summary is the only thing that still fits.
        (("distribution", 1, ["quantitative"]), "histogram"),
        (("distribution", 4, ["quantitative", "nominal"]), "small_multiples_histogram"),
        (("distribution", 40, ["quantitative", "nominal"]), "boxplot_by_category"),
        # Relationship. A scatter until the marks stop being individually
        # readable, then a density map that admits what it is showing.
        (("relationship", 200, ["quantitative"]), "scatter"),
        (("relationship", 10_000, ["quantitative"]), "hexbin"),
        # Composition. Small enough is a table; larger is a stacked bar --
        # never a pie, see the next test.
        (("composition", 3, ["nominal", "quantitative"]), "table"),
        (("composition", 12, ["nominal", "quantitative"]), "stacked_bar"),
        # Change over time. One line, or small multiples once the tangle
        # would beat the reader.
        (("change_over_time", 3, ["temporal", "quantitative"]), "line"),
        (("change_over_time", 20, ["temporal", "quantitative"]), "small_multiples_line"),
    ]
    for (kind, n, types), expected in cases:
        assert C.choose_chart(kind, n, types) == expected, (kind, n, types)


def test_3_choose_chart_never_recommends_a_pie():
    recommendations = set()
    for kind in sorted(C.QUESTION_KINDS):
        types = ["temporal", "quantitative"] if kind == "change_over_time" else [
            "nominal",
            "quantitative",
        ]
        for n in (1, 2, 3, 5, 6, 8, 9, 30, 2001, 10_000):
            recommendations.add(C.choose_chart(kind, n, types))

    assert "pie" not in recommendations
    assert "donut" not in recommendations
    # Specifically for ranking, the case where a pie is worst: the answer
    # is always the sorted bar chart, at every size past the table rule.
    for n in (6, 8, 30, 200):
        assert C.choose_chart("comparison", n, ["nominal", "quantitative"]) == (
            "sorted_horizontal_bar"
        )


def test_3_choose_chart_validates_its_inputs():
    with pytest.raises(ValueError, match="unknown question kind"):
        C.choose_chart("pretty", 8, ["quantitative"])
    with pytest.raises(ValueError, match="unknown data types"):
        C.choose_chart("comparison", 8, ["categorical"])
    with pytest.raises(ValueError, match="n_categories must be at least 1"):
        C.choose_chart("comparison", 0, ["quantitative"])
    # "Over time" with no time in the data is a question about something
    # else, and guessing at what would be worse than refusing.
    with pytest.raises(ValueError, match="needs a temporal variable"):
        C.choose_chart("change_over_time", 5, ["nominal", "quantitative"])


# --------------------------------------------------------------------------
# Exercise 4 -- colour deficiency, simulated and measured.
#
# The claim "a red/green pass-fail chart is unreadable for a substantial
# minority of your audience" is usually asserted. Here it is measured: run
# both colours through a published deuteranopia transform and see how much
# of the separation survives.
# --------------------------------------------------------------------------

# CIE76 distance below which two swatches in a legend are, for practical
# purposes, the same colour. The literature puts the just-noticeable
# difference for adjacent patches near 2.3; 10 is a deliberately generous
# threshold, chosen so passing this test is not a close-run thing.
COLLAPSE_THRESHOLD = 10.0

# And the distance a pair must keep to count as genuinely distinguishable.
SURVIVAL_THRESHOLD = 25.0


def test_4_red_green_pair_collapses_under_deuteranopia():
    red, green = PAL.PASS_FAIL_RED, PAL.PASS_FAIL_GREEN
    result = E.deuteranopia_collapse(red, green)

    # To a reader with typical colour vision these are about as far apart
    # as two colours get: measured 119.77 on this machine.
    assert result["normal_delta_e"] > 100.0

    # After the transform they are 7.31 apart -- closer together than many
    # people can reliably separate in a legend, and about 6% of the
    # separation a normal-vision reader gets.
    assert result["simulated_delta_e"] < COLLAPSE_THRESHOLD
    assert result["retained_fraction"] < 0.10

    # Simulation approximates a deficiency; it does not reproduce anyone's
    # experience. What this number licenses is "do not let colour alone
    # carry this distinction", not "this is what they see".
    assert 0.0 < result["simulated_delta_e"] < result["normal_delta_e"]


def test_4_colorblind_safe_pair_survives_the_same_transform():
    result = E.deuteranopia_collapse(PAL.SAFE_BLUE, PAL.SAFE_ORANGE)

    # seaborn's blue and orange start about as far apart as the red/green
    # pair (measured 115.70) and stay there: 116.51 after the transform,
    # so essentially all of the separation survives.
    assert result["normal_delta_e"] > 100.0
    assert result["simulated_delta_e"] > SURVIVAL_THRESHOLD
    assert result["retained_fraction"] > 0.90

    # The two pairs are close to equally distinguishable to a normal
    # -vision reader and nowhere near it afterwards. Same starting point,
    # opposite outcome -- which is the whole argument for choosing the
    # palette on purpose.
    unsafe = E.deuteranopia_collapse(PAL.PASS_FAIL_RED, PAL.PASS_FAIL_GREEN)
    assert result["simulated_delta_e"] > 10 * unsafe["simulated_delta_e"]


# --------------------------------------------------------------------------
# Exercise 5 -- an ordered variable on a categorical palette.
#
# Order is information the data HAS. A palette either carries it or
# destroys it, and the measurement is the rank correlation between the
# variable's order and the palette's luminance order.
# --------------------------------------------------------------------------


def test_5_sequential_palette_preserves_the_order():
    palette = PAL.viridis_steps(len(PAL.SATISFACTION_LEVELS))
    assert len(palette) == 5

    luminances = [E.relative_luminance(c) for c in palette]
    # Measured: 0.0190, 0.0885, 0.2234, 0.4511, 0.7826 -- strictly rising.
    assert luminances == sorted(luminances)

    rho = E.luminance_order_correlation(palette)
    assert rho == pytest.approx(1.0)

    # So a greyscale photocopy of a viridis-coded chart still reads in the
    # right order. That is not a nicety; it is the difference between a
    # legend being needed and being merely helpful.


def test_5_categorical_palette_destroys_the_order():
    palette = PAL.tab10_steps(len(PAL.SATISFACTION_LEVELS))
    assert len(palette) == 5

    luminances = [E.relative_luminance(c) for c in palette]
    # Measured: 0.1678, 0.3647, 0.2586, 0.1590, 0.1967 -- up, down, down,
    # up. "Satisfied" is DARKER than "dissatisfied", so the picture says
    # the opposite of the data.
    assert luminances != sorted(luminances)

    rho = E.luminance_order_correlation(palette)
    # Measured -0.2 on this machine: not merely weak, but pointing the
    # wrong way. tab10 is not defective -- it is doing its job, which is
    # to make neighbours look DIFFERENT, and difference has no direction.
    assert rho == pytest.approx(-0.2)
    assert abs(rho) < 0.5

    sequential_rho = E.luminance_order_correlation(PAL.viridis_steps(5))
    assert abs(sequential_rho) > abs(rho)


# --------------------------------------------------------------------------
# Exercise 6 -- sorting is an encoding decision, and its cost is the
# reader's, not yours.
# --------------------------------------------------------------------------


def test_6_sorting_changes_the_effort_not_the_answer():
    values = [float(v) for v in (37, 12, 88, 45, 3, 61, 29, 74, 18, 52, 9, 66, 41, 25, 80, 7, 58, 33, 95, 21)]
    assert len(values) == 20

    unsorted_cost = C.comparisons_to_find_max(values, presented_sorted=False)
    sorted_cost = C.comparisons_to_find_max(values, presented_sorted=True)

    # Source order: hold a running best, check all nineteen others.
    assert unsorted_cost == 19
    # Sorted: read the top row, glance at the second to confirm the chart
    # really is sorted, stop.
    assert sorted_cost == 1
    assert unsorted_cost == 19 * sorted_cost

    # And the answer is identical either way. Sorting moved nothing but
    # the reader's effort -- which is exactly why it is free to do and
    # expensive to skip.
    descending = sorted(values, reverse=True)
    assert descending[0] == values[C.index_of_max(values)] == 95.0

    # The same argument at the scale of a real category list.
    assert C.comparisons_to_find_max([0.0] * 200, presented_sorted=False) == 199
    assert C.comparisons_to_find_max([0.0] * 200, presented_sorted=True) == 1


# --------------------------------------------------------------------------
# Exercise 7 -- the data-ink ratio, counted in pixels.
#
# Same eight numbers, same bars, same labels. One chart adds a tinted
# panel, gridlines on both axes and a heavy box; the other does not. Count
# the ink and divide.
# --------------------------------------------------------------------------


def test_7_removing_furniture_raises_the_data_ink_ratio(png_dir):
    decorated = png_dir / "decorated.png"
    plain = png_dir / "plain.png"
    R.render_region_bar_chart(decorated, decorated=True)
    R.render_region_bar_chart(plain, decorated=False)

    total_decorated = R.count_non_background_pixels(decorated)
    total_plain = R.count_non_background_pixels(plain)
    ratio_decorated = R.data_ink_ratio(decorated, R.BAR_RGB)
    ratio_plain = R.data_ink_ratio(plain, R.BAR_RGB)

    # Measured on this machine: 172,351 total ink against 79,107 -- the
    # decorated chart spends more than twice the ink to say the same
    # thing.
    assert total_decorated == pytest.approx(172_351, rel=0.05)
    assert total_plain == pytest.approx(79_107, rel=0.05)
    assert total_decorated > 2 * total_plain

    # And the ratio moves in the expected direction, hard: 0.367 against
    # 0.934. Nearly two thirds of the decorated chart is furniture.
    assert ratio_decorated == pytest.approx(0.3669, abs=0.03)
    assert ratio_plain == pytest.approx(0.9344, abs=0.03)
    assert ratio_plain > ratio_decorated
    assert ratio_plain - ratio_decorated > 0.5

    # The point is not that gridlines are forbidden. It is that every mark
    # is a claim on the reader's attention, and the ones that are not data
    # should have to justify themselves.


# --------------------------------------------------------------------------
# Exercise 8 -- overplotting, and two ways out of it.
#
# Ten thousand points, one pixel each, none of them clipped. Count how
# many distinct pixels end up painted -- the shortfall is data that is on
# the page in principle and invisible in fact.
# --------------------------------------------------------------------------

N_POINTS = 10_000


def test_8_overplotting_hides_a_third_of_the_data(png_dir, points):
    x, y = points
    assert len(x) == N_POINTS
    # Nothing is clipped, so a shortfall below is overplotting and not
    # points falling off the edge of the picture.
    assert R.points_inside_axes(x, y) == N_POINTS

    opaque = png_dir / "opaque.png"
    R.render_scatter(opaque, x, y, alpha=1.0)
    painted = R.count_painted_pixels(opaque)

    # Measured: 6,349 painted pixels for 10,000 points. 3,651 points --
    # 36.5% of the data -- landed on a pixel another point had already
    # blackened and changed nothing about the image.
    assert painted == pytest.approx(6_349, rel=0.05)
    assert painted < N_POINTS
    assert painted / N_POINTS < 0.75

    # The stronger statement: the opaque image contains exactly TWO grey
    # levels, paper and ink. Whether a pixel carries one point or forty,
    # it is the same black -- the density information is not dimmed, it is
    # absent.
    assert R.count_distinct_luminance_levels(opaque) == 2


def test_8_alpha_and_hexbin_recover_the_density(png_dir, points):
    x, y = points

    opaque = png_dir / "opaque.png"
    blended = png_dir / "blended.png"
    hexes = png_dir / "hexbin.png"
    R.render_scatter(opaque, x, y, alpha=1.0)
    R.render_scatter(blended, x, y, alpha=0.05)
    R.render_hexbin(hexes, x, y)

    opaque_levels = R.count_distinct_luminance_levels(opaque)
    blended_levels = R.count_distinct_luminance_levels(blended)
    hexbin_levels = R.count_distinct_luminance_levels(hexes)

    # Measured: 2, 9 and 244. Alpha blending lets a pixel record HOW MANY
    # points landed on it, up to the point where the stack saturates --
    # 9 levels means the busiest pixel in this cloud carries about eight
    # points.
    assert opaque_levels == 2
    assert blended_levels >= 5
    assert blended_levels > opaque_levels

    # Hexbin does better still, because it aggregates BEFORE drawing
    # instead of hoping the compositor will: 244 distinct levels, a
    # genuine density surface rather than a smudge.
    assert hexbin_levels > 50
    assert hexbin_levels > blended_levels

    # And this is why `choose_chart` stops recommending a scatter past
    # `OVERPLOT_POINT_LIMIT` -- the recommendation and the measurement are
    # the same fact seen twice.
    assert C.choose_chart("relationship", N_POINTS, ["quantitative"]) == "hexbin"


# --------------------------------------------------------------------------
# Exercise 9 -- when a table beats a chart.
#
# The reflex to chart everything is itself a failure mode, and the
# threshold where charting starts to pay is a decision you should make on
# purpose rather than by habit.
# --------------------------------------------------------------------------


def test_9_a_table_beats_a_chart_below_the_threshold():
    small = C.choose_chart("comparison", 3, ["nominal", "quantitative"])
    large = C.choose_chart("comparison", 30, ["nominal", "quantitative"])

    assert small == "table"
    assert large == "sorted_horizontal_bar"

    # The boundary, stated: a chart's advantage is that it converts
    # comparison from arithmetic into a perceptual judgement. With three
    # numbers there was no arithmetic to convert -- the reader can just
    # read them, exactly, which a bar length never lets them do. Past
    # TABLE_MAX_VALUES the reader can no longer hold the set in their head
    # at once and the perceptual judgement starts to pay for the precision
    # it costs. Five is where this course puts the line; the number is a
    # judgement and the point is that it is written down.
    assert C.TABLE_MAX_VALUES == 5
    assert C.choose_chart("comparison", C.TABLE_MAX_VALUES, ["quantitative"]) == "table"
    assert C.choose_chart("comparison", C.TABLE_MAX_VALUES + 1, ["quantitative"]) == (
        "sorted_horizontal_bar"
    )

    # Composition obeys the same rule, which is where the last argument
    # for a pie chart goes: "just two or three slices" is precisely the
    # case a table answers exactly.
    assert C.choose_chart("composition", 2, ["nominal", "quantitative"]) == "table"
    assert C.choose_chart("composition", 3, ["nominal", "quantitative"]) == "table"

lesson_id: D127
day: 127
kind: guided-build
languages: [python, bash]
setup_commands:
  - cd labs/sections/math-statistics-and-data/day-127-why-we-visualize-and-choosing-the
  - python3 -m venv .venv
  - .venv/bin/pip install -r requirements/requirements.txt
  - .venv/bin/python3 -c "import matplotlib, seaborn; print(matplotlib.__version__, seaborn.__version__)"
run_commands:
  - .venv/bin/pytest examples
  - .venv/bin/pytest starter
test_commands:
  - bash tests/run_tests.sh
cleanup_commands:
  - "find . -path ./.venv -prune -o -type d -name '__pycache__' -print -exec rm -rf -- {} +"
  - rm -rf .pytest_cache
  - 'rm -rf .venv  # optional: removes the lab virtual environment'
  - 'git checkout -- starter/  # optional: reset your work'
requires_network: true
requires_api_key: false
estimated_minutes: 50
last_executed: '2026-08-20'
executed_on: >-
  macOS 26.5.2 (Apple Silicon, arm64), Python 3.14.0, matplotlib 3.11.1, seaborn 0.13.2,
  pandas 3.0.5, numpy 2.5.2, Pillow 12.3.0, pytest 9.1.1, bash 3.2.57 -- bash tests/run_tests.sh
  -> 19 checks, 0 failure(s), exit 0. pytest examples -> 17 passed. pytest starter -> 17 skipped
  (untouched checkout). pytest examples starter (one invocation) -> collection aborts with
  "import file mismatch", exit non-zero, confirmed directly rather than assumed, because all six
  of this lab's modules (encoding, charts, palettes, render, conftest, test_charts) are defined
  identically in both directories. Section 6 of the harness solves every exercise in a scratch
  copy (17 passed), deliberately breaks exercise 8's exact "distinct luminance levels == 2"
  assertion, confirms a non-zero exit and a printed failure, restores it, and confirms 17 passed
  again. Separately, the coordinator broke a different assertion directly in
  examples/test_charts.py (exercise 5's measured rank correlation of -0.2 changed to 0.9), not
  in a scratch copy, and re-ran the full harness: pytest reported 1 failed, 16 passed and the
  harness reported 6 of its 19 checks failing with overall exit 1; restoring the line returned
  it to 19 checks, 0 failure(s), exit 0. Everything was run through a real lab-local .venv
  created by the documented setup commands. Honesty notes from this run. FIRST: the harness's
  very first execution genuinely failed, 19 checks 1 failure, exit 1 -- the "nothing calls
  plt.show()" grep was unanchored and matched the sentence in render.py's own docstring
  promising that plt.show is never called. The grep is now anchored to the start of a statement,
  and that accidental failure is itself the first proof this harness can fail. SECOND: every
  pixel count in this lab is a real measurement off a real PNG and lands a couple of percent
  away from the ideal pi*r^2 (5156 px measured against 5026.5 ideal for r=40; 20368 against
  20106.2 for r=80) because a circle's boundary does not fall on pixel edges; the tests assert
  the RATIOS, 4.0 and 2.0 at rel=0.02, which survive rasterisation intact, and never the raw
  counts as if they were exact. THIRD: the deuteranopia matrix is Machado, Oliveira and
  Fernandes (2009) at severity 1.0, and its nine coefficients were read from the authors'
  published matrix table and checked against it rather than recalled; the lesson and the lab
  both state plainly that a simulation approximates a deficiency and does not reproduce anyone's
  experience. FOURTH: the commonly cited prevalence of red-green colour vision deficiency,
  roughly 8% of men of Northern European ancestry, is presented in the lesson as a widely
  reported figure that varies by population, and is not attached to any of this day's five
  verified sources, none of which covers it. FIFTH: scipy is not installed in this environment,
  so Spearman's rank correlation is written out by hand in encoding.py rather than called from
  scipy.stats -- the twenty-line midrank version, which also handles ties correctly where the
  usual shortcut formula does not. SIXTH: Vega-Lite and plotly are described in the lesson's
  Tools section from public documentation only; neither is installed here and no output
  attributed to either is reproduced anywhere in this lab or its lesson.

# What is installed, why, and what it costs

Six packages, all free and open source, installed into a lab-local
virtual environment that `rm -rf .venv` completely undoes.

| Package | Version pinned | Licence | What this lab uses it for |
| --- | --- | --- | --- |
| `matplotlib` | 3.11.1 | matplotlib licence (BSD-style, PSF-derived) | Every render in `render.py`, through the headless `Agg` backend; also the source of the `viridis` and `tab10` palettes exercise 5 measures. |
| `seaborn` | 0.13.2 | BSD 3-Clause | `color_palette("colorblind")` — the colourblind-safe pair exercise 4 measures. Nothing in this lab calls seaborn's plotting functions; it is used as a palette source, offline. |
| `pandas` | 3.0.5 | BSD 3-Clause | Not imported by this lab. It is pinned because seaborn declares it as a hard dependency, and leaving it unpinned would let the resolver pick a different version on a different day. |
| `numpy` | 2.5.2 | BSD 3-Clause | The point cloud in exercise 8, and every pixel count — the images are read as arrays. |
| `pillow` | 12.3.0 | MIT-CMU | Reads the rendered PNGs back off disk so their pixels can be counted. This is what makes "the chart uses more ink" a measurement. |
| `pytest` | 9.1.1 | MIT | The test harness every exercise is written against. |

`math`, `pathlib` and `tempfile` are Python standard library — no
install, no cost. `math` does the square-law arithmetic and the CIELAB
conversion; `tempfile` creates the directory every render is written
into, so nothing lands inside the lab.

Deliberately **not** installed, and described in the lesson from public
documentation only: Vega-Lite and plotly. No output attributed to either
is reproduced anywhere in this lab.

There is no paid tier of anything here, no account, no key and no signup,
personally or commercially.

matplotlib==3.11.1
seaborn==0.13.2
pandas==3.0.5
numpy==2.5.2
pillow==12.3.0
pytest==9.1.1

# Day 127 lab — the brief

Nine exercises, in order. Work top to bottom in `test_charts.py`. The
instruments live in four modules you read but do **not** edit:

- `encoding.py` — the square law (`encoded_area_ratio`), sRGB → CIELAB
  (`delta_e_cie76`), the deuteranopia transform
  (`simulate_deuteranopia`, `deuteranopia_collapse`), WCAG luminance, and
  Spearman's rho written out by hand because scipy is not installed here.
- `charts.py` — `ENCODING_RANKING` and `encoding_rank`, the two decision
  functions `best_encoding` and `choose_chart`, the three thresholds
  (`TABLE_MAX_VALUES`, `OVERPLOT_POINT_LIMIT`, `SMALL_MULTIPLE_LIMIT`),
  and `comparisons_to_find_max`.
- `palettes.py` — the swatches, taken from matplotlib's and seaborn's own
  defaults rather than hand-picked to prove a point.
- `render.py` — everything that draws (`render_circle`,
  `render_region_bar_chart`, `render_scatter`, `render_hexbin`) and
  everything that measures a drawing (`count_non_background_pixels`,
  `count_pixels_of_color`, `data_ink_ratio`,
  `count_distinct_luminance_levels`).

Read all four once before you start. Two fixtures come from
`conftest.py`: `points` (the seeded 10,000-point cloud) and `png_dir` (a
temporary directory outside the lab — every render goes there, and
nowhere else).

Check yourself at any point:

```bash
.venv/bin/pytest starter -v
```

On an untouched checkout that prints `17 skipped`. A **skip** means "not
attempted". Replace a `pytest.skip(...)` line with real assertions and
delete it — when every skip is gone and the suite is green, you are
finished:

```bash
.venv/bin/pytest starter -q
echo $?
```

One rule about tolerances. Pure arithmetic (`encoded_area_ratio`, a rank
correlation, a comparison count) is exact — assert it exactly. Anything
measured off a rendered image is a rasteriser's opinion about where a
circle's edge falls, so give it a tolerance: `pytest.approx(4.0,
rel=0.02)` is honest, and an exact equality that happens to pass is luck
you should not depend on.

Nothing in this lab depends on timing. Never assert a duration.

---

## Exercise 1 — the square law, measured

Two values, 50 and 100. The second is exactly twice the first.

**1a, the arithmetic.** Call `encoding.encoded_area_ratio(VALUES,
mode="radius")` and again with `mode="area"`. The radius encoding gives a
circle twice as wide, which is a circle **four** times as large: assert
`4.0`, and assert it equals the data ratio *squared*. The area encoding
gives `2.0`. Assert the two encodings differ by a factor of exactly 2 —
that factor is the exaggeration a radius-scaled bubble chart applies to
every comparison in it, without misstating a single number.

**1b, the same claim off real pixels.** Render three circles into
`png_dir` with `render.render_circle`: radius 40, radius 80, and radius
`40 * sqrt(2)`. Measure each with `render.measure_circle_area_px`, which
counts painted pixels. Assert the 80-pixel circle covers about four times
the pixels of the 40-pixel one and the `40*sqrt(2)` one about twice, both
with `rel=0.02`. Report both counts in the test's own comments — the
measured numbers on the authoring machine are 5,156 and 20,368 pixels
against an ideal `pi*r^2` of 5,026.5 and 20,106.2, a couple of percent
out because a circle's edge does not land on pixel boundaries. The ratio
survives that intact, which is why the ratio is what the exercise asserts.

## Exercise 2 — the perceptual ranking as a decision function

Cleveland and McGill (1984) ran experiments asking people to judge
magnitudes from each visual channel and ranked the channels by measured
accuracy. `charts.ENCODING_RANKING` is that ranking.

**2a.** Use `charts.encoding_rank` to assert the ordering: position on a
common scale, then position on non-aligned scales, then length, then
angle and slope, then area, then volume, then colour saturation. Then
assert `encoding_rank("hue")` raises `ValueError`. Hue is not a bad
magnitude channel — it is not a magnitude channel at all, and the error
says so.

**2b.** Build a list of `((data_type, task), expected_channel)` cases and
assert `charts.best_encoding` returns each. Cover at least
`quantitative`/`compare`, `quantitative`/`compare_across_panels`,
`quantitative`/`magnitude_on_map`, `ordinal`/`encode_in_color`,
`nominal`/`identify_group` and `nominal`/`compare`. **Justify each
expected answer in a comment** — the justification is the exercise, not
the assertion. Finish with the load-bearing negative: assert
`best_encoding("ordinal", "encode_in_color")` is *not* `"hue"`, which
exercise 5 then measures the reason for.

**2c.** Assert `best_encoding` raises `ValueError` matching
`"unknown data type"` and `"unknown task"` respectively.

## Exercise 3 — from the question to the chart

`choose_chart(question_kind, n_categories, data_types)` takes the reader's
question and names an instrument.

**3a.** Build a case table covering all five question kinds, with at
least one case on each side of `TABLE_MAX_VALUES`,
`OVERPLOT_POINT_LIMIT` and `SMALL_MULTIPLE_LIMIT`, and assert every one.

**3b.** Collect the function's answers over every question kind and a
spread of sizes into a set, and assert that neither `"pie"` nor
`"donut"` is in it — anywhere, for anything. Then assert that for a
ranking question the answer is always `"sorted_horizontal_bar"`.

**3c.** Assert it raises for an unknown question kind, an unknown data
type, `n_categories=0`, and `change_over_time` with no temporal variable
in `data_types`. That last one matters: a function that guesses at a time
axis it was not given is worse than one that refuses.

## Exercise 4 — colour deficiency, simulated and measured

`encoding.simulate_deuteranopia` applies the Machado, Oliveira and
Fernandes (2009) severity-1.0 matrix in linear RGB.
`encoding.deuteranopia_collapse` returns the CIE76 distance before, the
distance after, and the fraction retained.

**4a.** Run it on `PAL.PASS_FAIL_RED` and `PAL.PASS_FAIL_GREEN` —
matplotlib's own `tab10` red and green, the two colours a pass/fail chart
gets by default. Assert the normal-vision distance is above 100, the
simulated distance is below `COLLAPSE_THRESHOLD` (10.0), and the retained
fraction is below 0.10. Report all three.

**4b.** Run it on `PAL.SAFE_BLUE` and `PAL.SAFE_ORANGE`, the first two
entries of seaborn's `colorblind` palette. Assert the simulated distance
stays above `SURVIVAL_THRESHOLD` and the retained fraction above 0.90,
and that the safe pair ends more than ten times further apart than the
red/green pair.

Write into your comments what this does and does not license. A
simulation approximates a deficiency; it does not reproduce anyone's
experience, it assumes one severity, and it cannot represent anomalous
trichromacy at all. A small simulated distance is strong evidence a
palette is risky. A large one is weak evidence it is fine.

## Exercise 5 — an ordered variable on a categorical palette

**5a.** Take `PAL.viridis_steps(5)`, compute
`encoding.relative_luminance` of each swatch, and assert the list comes
out already sorted. Then assert
`encoding.luminance_order_correlation(palette)` is approximately `1.0`.

**5b.** Do the same with `PAL.tab10_steps(5)`. Assert the luminance list
is *not* sorted and the rank correlation is far smaller in absolute
value. Report the number you measure. `tab10` is not defective — it is
doing exactly its job, which is to make neighbouring swatches look as
different as possible, and "different" has no direction. Putting an
ordered variable on it destroys order the data actually had.

## Exercise 6 — sorting is an encoding decision

Build a list of 20 unsorted values.
`charts.comparisons_to_find_max(values, presented_sorted=False)` models a
reader holding a running best and checking every remaining bar: 19
comparisons. With `presented_sorted=True` it is 1 — read the top row,
glance at the second to confirm the chart really is sorted, stop.

Assert both, and then assert the *answer* is identical either way:
`sorted(values, reverse=True)[0] == values[charts.index_of_max(values)]`.
Sorting moved nothing but the reader's effort, which is why skipping it
is expensive and doing it is free.

## Exercise 7 — the data-ink ratio

Render `render.render_region_bar_chart` into `png_dir` twice, once with
`decorated=True` and once with `decorated=False`. Same eight numbers,
same bars, same labels; the decorated version adds a tinted panel,
gridlines on both axes and a heavy box.

Count total ink with `render.count_non_background_pixels` and the data
fraction with `render.data_ink_ratio(path, render.BAR_RGB)` — the bars
are drawn in one flat colour nothing else uses, so that isolates the data
ink exactly. Assert the plain chart's ratio is higher and the gap is more
than 0.5. **Report both counts and both ratios.**

The point is not that gridlines are banned. It is that every mark claims
some of the reader's attention, and the ones that are not data should
have to justify themselves.

## Exercise 8 — overplotting

**8a.** Confirm with `render.points_inside_axes` that none of the 10,000
points is clipped — otherwise a low pixel count could just mean points
fell off the edge. Render with `alpha=1.0` and count painted pixels.
Assert the count is below 75% of `N_POINTS`: every missing pixel is a
point that landed where another point already was and changed nothing.

Then the sharper measurement: assert
`render.count_distinct_luminance_levels` of that image is exactly **2**.
Paper and ink. Whether a pixel carries one point or forty it is the same
black, so the density is not dimmed — it is absent.

**8b.** Render the same cloud three ways: `alpha=1.0`, `alpha=0.05`, and
`render.render_hexbin`. Count distinct luminance levels in each. Assert
the alpha version has more than the opaque one, and the hexbin more than
50. Finish by asserting
`charts.choose_chart("relationship", N_POINTS, ["quantitative"])` is
`"hexbin"` — the recommendation and the measurement are the same fact
seen twice.

## Exercise 9 — when a table beats a chart

Assert `choose_chart("comparison", 3, ...)` is `"table"` and
`choose_chart("comparison", 30, ...)` is `"sorted_horizontal_bar"`.
Assert the boundary sits exactly at `charts.TABLE_MAX_VALUES` by testing
that value and that value plus one.

Then write the comment that is the actual exercise: **why** the boundary
is where it is. A chart's advantage is that it converts comparison from
arithmetic into a perceptual judgement. With three numbers there was no
arithmetic to convert — the reader can simply read them, exactly, which a
bar length never lets them do. Past five values the set stops fitting in
the reader's head at once and the perceptual judgement starts to pay for
the precision it costs. Five is where this course draws the line; the
number is a judgement, and the point is that it is written down where you
can argue with it.

"""Chart choice as a function you can call, argue with, and test.

"Which chart should I use?" is normally answered by taste, which is why it
is normally answered badly. Two published results turn most of it into
arithmetic:

* Cleveland and McGill (1984) ranked visual channels by how accurately
  people read magnitudes off them. That ranking is `ENCODING_RANKING`.
* A variable's TYPE decides which channels can carry it honestly at all.
  A nominal variable on a length channel invents an order that is not in
  the data; an ordinal variable on a categorical hue palette destroys the
  order that is.

`best_encoding` combines the two: given a data type and the reader's task,
return the most accurately-judged channel that is still honest. Then
`choose_chart` answers the question one level up -- given the QUESTION,
how many values, and what types they are, name the chart.

These functions are opinionated on purpose. They are not the last word on
visualisation; they are a written-down default you can disagree with
explicitly, which is strictly better than an unwritten one you cannot.
"""

from __future__ import annotations

# --------------------------------------------------------------------------
# The Cleveland-McGill ordering
# --------------------------------------------------------------------------

# Ordered most accurately judged first. From Cleveland and McGill (1984),
# "Graphical Perception: Theory, Experimentation, and Application to the
# Development of Graphical Methods", Journal of the American Statistical
# Association 79(387):531-554. Their elementary perceptual tasks, with the
# names shortened for use as identifiers.
ENCODING_RANKING: tuple[str, ...] = (
    "position_common_scale",  # 1. dots or bar ends on one shared axis
    "position_nonaligned_scales",  # 2. the same, but across separate panels
    "length",  # 3. bar lengths not sharing a baseline
    "angle_slope",  # 4. pie slices, and line steepness
    "area",  # 5. bubble size, treemap tiles
    "volume",  # 6. 3-D bars, spheres
    "color_saturation",  # 7. heatmap intensity, density shading
)

# Channels that carry identity but NOT magnitude. Cleveland and McGill did
# not rank these, because there is no magnitude to read off them: hue is
# how you say "these two lines are different series", never "this one is
# 1.8 times that one".
IDENTITY_CHANNELS: tuple[str, ...] = ("hue", "shape")

# Sequential lightness -- a colormap like viridis, ordered by luminance.
# It carries ORDER faithfully and magnitude only roughly, which places it
# with saturation at the bottom of the accuracy ordering but well above
# categorical hue for anything ordered.
ORDERED_COLOR_CHANNEL = "luminance_sequential"

DATA_TYPES: frozenset[str] = frozenset({"nominal", "ordinal", "quantitative", "temporal"})


def encoding_rank(channel: str) -> int:
    """Position of `channel` in the Cleveland-McGill ordering, 0 = best."""
    try:
        return ENCODING_RANKING.index(channel)
    except ValueError:
        raise ValueError(
            f"{channel!r} is not a ranked magnitude channel; "
            f"ranked channels are {ENCODING_RANKING}"
        ) from None


# `best_encoding`'s task vocabulary. Each task says what the reader is
# trying to DO, because the same variable in the same chart deserves a
# different channel depending on the question being asked of it.
TASKS: frozenset[str] = frozenset(
    {
        "compare",  # read two magnitudes and say which is bigger, by how much
        "compare_across_panels",  # the same, but the values sit in separate small multiples
        "identify_group",  # tell which series a mark belongs to
        "encode_in_color",  # the two spatial axes are taken; colour is all that is left
        "magnitude_on_map",  # position is spent on geography, so magnitude needs another channel
        "trend",  # read direction and rate of change over an ordered axis
    }
)


def best_encoding(data_type: str, task: str) -> str:
    """Most accurately-judged channel that can honestly carry `data_type`.

    The rule in one sentence: take the highest-ranked channel from
    Cleveland and McGill that (a) the task has not already spent on
    something else, and (b) does not claim more structure than the data
    type has.
    """
    if data_type not in DATA_TYPES:
        raise ValueError(f"unknown data type {data_type!r}; expected one of {sorted(DATA_TYPES)}")
    if task not in TASKS:
        raise ValueError(f"unknown task {task!r}; expected one of {sorted(TASKS)}")

    if task == "identify_group":
        # Identity, not magnitude. A nominal variable has no order to
        # preserve and no size to read, so hue is exactly right -- and
        # putting it on a magnitude channel would invent an order.
        return "hue"

    if task == "encode_in_color":
        # Both spatial axes are already spent. Colour is the only channel
        # left, and the data type decides WHICH colour channel.
        if data_type == "nominal":
            return "hue"
        if data_type in ("ordinal", "temporal"):
            # Order must survive. Categorical hue would destroy it; a
            # luminance-ordered ramp preserves it.
            return ORDERED_COLOR_CHANNEL
        return "color_saturation"  # quantitative: magnitude, read roughly

    if task == "magnitude_on_map":
        # Position is spent on latitude and longitude. Of what remains,
        # area is the highest-ranked channel that survives being placed at
        # an arbitrary map location -- length would need a shared baseline
        # the map cannot provide.
        return "area"

    if task == "compare_across_panels":
        # Small multiples: each panel has its own axis, so this is
        # Cleveland and McGill's second task by construction.
        return "position_nonaligned_scales"

    if task == "trend":
        # Direction over an ordered axis. Slope is the thing being read,
        # but it is read off points placed on one common scale, and the
        # accuracy of the reading is the accuracy of those positions.
        return "position_common_scale"

    # task == "compare": one shared axis, the best channel there is. This
    # holds for every data type, including nominal -- a bar chart of
    # nominal categories puts the CATEGORY on the categorical axis and the
    # QUANTITY on the common scale, which invents no order at all.
    return "position_common_scale"


# --------------------------------------------------------------------------
# From the question to the chart
# --------------------------------------------------------------------------

QUESTION_KINDS: frozenset[str] = frozenset(
    {
        "comparison",  # which of these is biggest? by how much?
        "distribution",  # what shape is this variable? where is the mass?
        "relationship",  # does x move with y?
        "composition",  # what are the parts of this whole?
        "change_over_time",  # what happened, and in which direction?
    }
)

# Below this many values, a table is the better instrument. The number is
# a judgement, not a measurement, and the reasoning is written out in
# `choose_chart`'s docstring so you can move it deliberately.
TABLE_MAX_VALUES = 5

# Above this many points, individual marks stop being individually
# readable and the chart is showing density whether you meant it to or
# not. Better to say so and draw density directly.
OVERPLOT_POINT_LIMIT = 2000

# Above this many series, one panel becomes a tangle. Split it.
SMALL_MULTIPLE_LIMIT = 8


def choose_chart(question_kind: str, n_categories: int, data_types: list[str]) -> str:
    """Recommend a chart -- or a table -- for a question.

    `n_categories` is the number of values the reader must take in: bars in
    a bar chart, points in a scatter, lines in a line chart.

    Two recommendations in here are the point of the whole function.

    **It recommends a table below `TABLE_MAX_VALUES` values.** Three
    numbers do not need an axis, a legend and a title in order to be
    compared; they need to be readable. A chart's advantage is that it
    turns comparison into a perceptual judgement instead of an arithmetic
    one, and with three numbers there was never any arithmetic to save.
    The threshold is a judgement call, and 5 is where the trade tips in
    this course's experience -- put it somewhere else if your readers
    differ, but put it somewhere on purpose.

    **It never recommends a pie chart, for anything.** Reading a pie means
    judging angle and area, ranked fourth and fifth by Cleveland and
    McGill, when the identical data on a sorted bar chart would be
    position on a common scale, ranked first. The one case usually offered
    in a pie's defence -- two or three parts of a whole -- falls below
    `TABLE_MAX_VALUES` and gets a table, which answers the question
    exactly rather than approximately.
    """
    if question_kind not in QUESTION_KINDS:
        raise ValueError(
            f"unknown question kind {question_kind!r}; expected one of {sorted(QUESTION_KINDS)}"
        )
    if n_categories < 1:
        raise ValueError(f"n_categories must be at least 1, got {n_categories}")
    unknown = set(data_types) - DATA_TYPES
    if unknown:
        raise ValueError(f"unknown data types {sorted(unknown)}; expected from {sorted(DATA_TYPES)}")

    if question_kind == "change_over_time":
        # A time axis is not optional here: "over time" with no temporal
        # variable is a question about something else.
        if "temporal" not in data_types:
            raise ValueError(
                "change_over_time needs a temporal variable in data_types; "
                f"got {sorted(data_types)}"
            )
        return "small_multiples_line" if n_categories > SMALL_MULTIPLE_LIMIT else "line"

    if question_kind == "relationship":
        if "quantitative" not in data_types:
            raise ValueError(
                "relationship needs at least one quantitative variable; " f"got {sorted(data_types)}"
            )
        # Past the overplot limit the marks are stacked on each other and
        # the reader is looking at ink density, not at points.
        return "hexbin" if n_categories > OVERPLOT_POINT_LIMIT else "scatter"

    if question_kind == "distribution":
        if "quantitative" not in data_types:
            raise ValueError(
                "distribution needs a quantitative variable; " f"got {sorted(data_types)}"
            )
        if n_categories <= 1:
            return "histogram"
        if n_categories <= SMALL_MULTIPLE_LIMIT:
            return "small_multiples_histogram"
        return "boxplot_by_category"

    # comparison and composition share the table rule.
    if n_categories <= TABLE_MAX_VALUES:
        return "table"

    if question_kind == "composition":
        return "stacked_bar"

    return "sorted_horizontal_bar"


# --------------------------------------------------------------------------
# Sorting, measured as reader effort
# --------------------------------------------------------------------------


def comparisons_to_find_max(values: list[float], presented_sorted: bool) -> int:
    """How many comparisons a reader makes to find the largest value.

    This simulates a reader, not a computer. Two behaviours, both real:

    * **Source order.** The reader has no idea where the biggest bar is,
      so they hold a running best and check every remaining bar against
      it: `n - 1` comparisons for `n` bars.
    * **Sorted, descending.** The reader looks at the top row. They still
      make ONE comparison -- top row against the next one -- to confirm
      the chart really is sorted and the leader really is ahead. After
      that they stop, because sorting is a promise about everything below.

    The answer is identical either way. The effort is not, and sorting is
    the cheapest thing you will ever do to a chart.
    """
    n = len(values)
    if n < 2:
        raise ValueError("finding a maximum needs at least two values")
    return 1 if presented_sorted else n - 1


def index_of_max(values: list[float]) -> int:
    """Index of the largest value, for confirming sorting changes only effort."""
    return max(range(len(values)), key=lambda i: values[i])

"""Shared fixtures. pytest finds this file by itself -- nothing imports it.

`points` returns the same seeded cloud on every test, so two tests that
render it are rendering the identical data.

`png_dir` is a temporary directory OUTSIDE the lab, created and removed by
the fixture itself. Every render in this suite writes there and nowhere
else, which is why `tests/run_tests.sh` can assert afterwards that the lab
directory contains no `.png` at all.
"""

from __future__ import annotations

import tempfile
from pathlib import Path

import pytest

import render as R


@pytest.fixture
def points():
    return R.sample_points()


@pytest.fixture
def png_dir():
    with tempfile.TemporaryDirectory(prefix="d127-render-") as d:
        yield Path(d)

"""The measurable half of visual encoding: geometry and colour arithmetic.

Nothing in this module renders anything. Everything here is a pure
function over numbers, which is exactly why it can be asserted on. The
rendering lives in `render.py`.

Two ideas carry the whole file.

1. **Area is quadratic in radius.** If you draw a value as a circle's
   RADIUS, doubling the value quadruples the ink. The reader judges the
   ink, so every ratio in the chart comes out squared. `area_for_radius`
   and the two `radii_*` functions make that concrete enough to assert on.

2. **A colour pair that is far apart for you may be on top of each other
   for a reader with a colour vision deficiency.** `simulate_deuteranopia`
   applies a published transformation matrix; `delta_e_cie76` measures how
   far apart two colours are afterwards. Neither function has an opinion --
   they return numbers, and the tests read them.

Colour maths conventions used throughout:

* An sRGB colour is a 3-tuple of floats in [0, 1] -- the same convention
  matplotlib and seaborn hand you.
* "Linear RGB" is sRGB with its transfer function undone. The colour
  vision deficiency matrix is defined on LINEAR RGB, not on the gamma
  encoded values, so `simulate_deuteranopia` linearises first and
  re-encodes afterwards. Skipping that step is the single most common
  error in hand-rolled deficiency simulators and it visibly changes the
  answer.
* CIELAB is computed against the D65 white point, matching sRGB's own
  reference white.
"""

from __future__ import annotations

import math

# --------------------------------------------------------------------------
# 1. Size encoding -- the square law
# --------------------------------------------------------------------------


def area_for_radius(radius: float) -> float:
    """Area of a circle of this radius. The whole square law in one line."""
    if radius < 0:
        raise ValueError(f"radius must be non-negative, got {radius}")
    return math.pi * radius * radius


def radii_scaled_by_radius(values: list[float], unit_radius: float = 1.0) -> list[float]:
    """Encode each value as a circle whose RADIUS is proportional to it.

    This is what you get from the obvious, wrong implementation --
    `plt.scatter(x, y, s=value)` reads `s` as an area, but a hand-written
    `Circle(xy, radius=value)` or a d3 `.attr("r", value)` does exactly
    this. A value twice as large is drawn with twice the radius and
    therefore FOUR times the area.
    """
    return [unit_radius * float(v) for v in values]


def radii_scaled_by_area(values: list[float], unit_radius: float = 1.0) -> list[float]:
    """Encode each value as a circle whose AREA is proportional to it.

    Take the square root and the distortion disappears: a value twice as
    large is drawn with twice the ink, so the ratio the reader perceives
    is the ratio in the data.
    """
    out = []
    for v in values:
        v = float(v)
        if v < 0:
            raise ValueError(f"cannot encode a negative value as an area: {v}")
        out.append(unit_radius * math.sqrt(v))
    return out


def encoded_area_ratio(values: list[float], mode: str) -> float:
    """Ratio of drawn AREA between the last and first value, under `mode`.

    `mode` is "radius" (the distorting encoding) or "area" (the honest
    one). Returns a plain float so a test can compare it against the
    ratio in the data itself.
    """
    if mode == "radius":
        radii = radii_scaled_by_radius(values)
    elif mode == "area":
        radii = radii_scaled_by_area(values)
    else:
        raise ValueError(f"mode must be 'radius' or 'area', got {mode!r}")
    first, last = area_for_radius(radii[0]), area_for_radius(radii[-1])
    if first == 0:
        raise ValueError("the first value encodes to zero area; ratio is undefined")
    return last / first


# --------------------------------------------------------------------------
# 2. Colour spaces
# --------------------------------------------------------------------------


def _srgb_to_linear_channel(c: float) -> float:
    """Undo sRGB's transfer function for one channel in [0, 1]."""
    return c / 12.92 if c <= 0.04045 else ((c + 0.055) / 1.055) ** 2.4


def _linear_to_srgb_channel(c: float) -> float:
    """Re-apply sRGB's transfer function for one channel, then clamp."""
    c = max(0.0, min(1.0, c))
    v = 12.92 * c if c <= 0.0031308 else 1.055 * (c ** (1 / 2.4)) - 0.055
    return max(0.0, min(1.0, v))


def srgb_to_linear(rgb: tuple[float, float, float]) -> tuple[float, float, float]:
    return tuple(_srgb_to_linear_channel(float(c)) for c in rgb)  # type: ignore[return-value]


def linear_to_srgb(rgb: tuple[float, float, float]) -> tuple[float, float, float]:
    return tuple(_linear_to_srgb_channel(float(c)) for c in rgb)  # type: ignore[return-value]


def relative_luminance(rgb: tuple[float, float, float]) -> float:
    """WCAG relative luminance: the perceived lightness of a colour.

    This is the quantity a greyscale photocopy of your chart keeps and
    everything else throws away. A palette that carries order in hue but
    not in luminance is a palette whose order vanishes in print.
    """
    r, g, b = srgb_to_linear(rgb)
    return 0.2126 * r + 0.7152 * g + 0.0722 * b


# D65 white point in XYZ, scaled so Y = 1. sRGB's own reference white.
_D65 = (0.95047, 1.00000, 1.08883)


def srgb_to_xyz(rgb: tuple[float, float, float]) -> tuple[float, float, float]:
    r, g, b = srgb_to_linear(rgb)
    x = 0.4124564 * r + 0.3575761 * g + 0.1804375 * b
    y = 0.2126729 * r + 0.7151522 * g + 0.0721750 * b
    z = 0.0193339 * r + 0.1191920 * g + 0.9503041 * b
    return (x, y, z)


def _lab_f(t: float) -> float:
    delta = 6 / 29
    return t ** (1 / 3) if t > delta**3 else t / (3 * delta**2) + 4 / 29


def srgb_to_lab(rgb: tuple[float, float, float]) -> tuple[float, float, float]:
    """Convert sRGB to CIELAB (L*, a*, b*) against D65.

    CIELAB exists because RGB distance is a poor model of perceived
    difference: two RGB triples the same Euclidean distance apart can look
    obviously different or nearly identical. CIELAB was built so that
    Euclidean distance is at least roughly proportional to perceived
    difference, which is what makes `delta_e_cie76` meaningful at all.
    """
    x, y, z = srgb_to_xyz(rgb)
    fx, fy, fz = _lab_f(x / _D65[0]), _lab_f(y / _D65[1]), _lab_f(z / _D65[2])
    return (116 * fy - 16, 500 * (fx - fy), 200 * (fy - fz))


def delta_e_cie76(a: tuple[float, float, float], b: tuple[float, float, float]) -> float:
    """CIE76 colour difference: Euclidean distance in CIELAB.

    The 1976 formula is the simplest of the family and is used here
    deliberately: it is short enough to read, has no tuning constants to
    get wrong, and the tests only ever ask whether a distance is large or
    small, never for a precise perceptual match. Rough guide, from the
    literature: about 2.3 is the "just noticeable difference" for adjacent
    patches; a difference in the low single digits is two colours a reader
    will struggle to tell apart at all in a legend.
    """
    la, aa, ba = srgb_to_lab(a)
    lb, ab, bb = srgb_to_lab(b)
    return math.sqrt((la - lb) ** 2 + (aa - ab) ** 2 + (ba - bb) ** 2)


# --------------------------------------------------------------------------
# 3. Colour vision deficiency simulation
# --------------------------------------------------------------------------

# Machado, Oliveira and Fernandes (2009), "A Physiologically-based Model
# for Simulation of Color Vision Deficiency", IEEE Transactions on
# Visualization and Computer Graphics 15(6):1291-1298. This is the
# deuteranomaly matrix at severity 1.0 -- that is, deuteranopia, the
# complete absence of a working M cone. The nine coefficients were read
# from the authors' published matrix table, not reconstructed from memory.
#
# The matrix operates on LINEAR RGB. `simulate_deuteranopia` linearises
# before multiplying and re-encodes afterwards.
DEUTERANOPIA_MATRIX_MACHADO_2009 = (
    (0.367322, 0.860646, -0.227968),
    (0.280085, 0.672501, 0.047413),
    (-0.011820, 0.042940, 0.968881),
)


def simulate_deuteranopia(
    rgb: tuple[float, float, float],
    matrix: tuple[tuple[float, float, float], ...] = DEUTERANOPIA_MATRIX_MACHADO_2009,
) -> tuple[float, float, float]:
    """Return an sRGB approximation of how `rgb` appears to a deuteranope.

    A very important limit, stated here rather than buried: this
    APPROXIMATES a deficiency, it does not reproduce an experience. The
    model is a linear transform fitted to a physiological model of cone
    response; it says nothing about how a particular person has learned to
    compensate over a lifetime, it assumes a single severity, and it
    cannot represent the many varieties of anomalous trichromacy at all.
    Treat a small simulated distance as strong evidence that a palette is
    risky, and a large one as weak evidence that it is fine. The reliable
    move is still to carry the distinction in a second channel -- shape,
    position, direct labelling -- so colour is never load-bearing alone.
    """
    r, g, b = srgb_to_linear(rgb)
    out = tuple(m[0] * r + m[1] * g + m[2] * b for m in matrix)
    return linear_to_srgb(out)  # type: ignore[arg-type]


def deuteranopia_collapse(
    a: tuple[float, float, float], b: tuple[float, float, float]
) -> dict[str, float]:
    """Measure how much of a pair's separation survives deuteranopia.

    Returns the normal-vision distance, the simulated distance, and the
    fraction retained. A pair whose separation nearly vanishes is a pair
    that must not be the only thing distinguishing two series.
    """
    normal = delta_e_cie76(a, b)
    simulated = delta_e_cie76(simulate_deuteranopia(a), simulate_deuteranopia(b))
    return {
        "normal_delta_e": normal,
        "simulated_delta_e": simulated,
        "retained_fraction": simulated / normal if normal else 0.0,
    }


# --------------------------------------------------------------------------
# 4. Rank correlation, without scipy
# --------------------------------------------------------------------------


def _ranks(values: list[float]) -> list[float]:
    """Ranks, averaging ties -- the standard midrank treatment."""
    order = sorted(range(len(values)), key=lambda i: values[i])
    ranks = [0.0] * len(values)
    i = 0
    while i < len(order):
        j = i
        while j + 1 < len(order) and values[order[j + 1]] == values[order[i]]:
            j += 1
        midrank = (i + j) / 2 + 1
        for k in range(i, j + 1):
            ranks[order[k]] = midrank
        i = j + 1
    return ranks


def spearman_rho(xs: list[float], ys: list[float]) -> float:
    """Spearman's rank correlation, computed from Pearson on the ranks.

    scipy is not installed in this environment, so this is the honest
    twenty-line version rather than a call to `scipy.stats.spearmanr`.
    Written this way it also handles ties correctly, which the shortcut
    "1 - 6*sum(d^2)/(n^3-n)" formula quietly does not.
    """
    if len(xs) != len(ys):
        raise ValueError("spearman_rho needs two sequences of the same length")
    if len(xs) < 2:
        raise ValueError("spearman_rho needs at least two observations")
    rx, ry = _ranks(list(xs)), _ranks(list(ys))
    mx, my = sum(rx) / len(rx), sum(ry) / len(ry)
    num = sum((a - mx) * (b - my) for a, b in zip(rx, ry))
    dx = math.sqrt(sum((a - mx) ** 2 for a in rx))
    dy = math.sqrt(sum((b - my) ** 2 for b in ry))
    if dx == 0 or dy == 0:
        raise ValueError("a sequence with no variation has no rank correlation")
    return num / (dx * dy)


def luminance_order_correlation(palette: list[tuple[float, float, float]]) -> float:
    """Rank correlation between a palette's position and its luminance.

    A SEQUENTIAL palette is built so that step 1 is darker than step 2 is
    darker than step 3; its correlation is +/-1 exactly. A CATEGORICAL
    palette is built so that neighbouring swatches are as DIFFERENT as
    possible, which says nothing at all about their lightness order --
    so mapping an ordered variable onto one destroys the order.
    """
    positions = [float(i) for i in range(len(palette))]
    luminances = [relative_luminance(c) for c in palette]
    return spearman_rho(positions, luminances)

"""The palettes this lab measures, taken from the libraries themselves.

Nothing here is a hand-picked colour chosen to make a point. Every swatch
is what matplotlib or seaborn hands you by default, which is what makes
the measurements in `test_charts.py` about real tools rather than about a
straw man.

* `PASS_FAIL_RED` and `PASS_FAIL_GREEN` are `tab10`'s red and green --
  matplotlib's default cycle, positions 3 and 2. The reflex "red means
  failed, green means passed" chart is drawn in exactly these two colours
  by anyone who does not stop to think about it.
* `SAFE_BLUE` and `SAFE_ORANGE` are the first two entries of seaborn's
  `colorblind` palette, which is that library's own answer to the same
  problem.
* `viridis_steps` samples matplotlib's default SEQUENTIAL colormap, which
  is built to increase monotonically in luminance.
* `tab10_steps` samples its default CATEGORICAL palette, which is built so
  neighbouring swatches look as different as possible -- and therefore
  carries no luminance order at all.

Both `viridis` and `tab10` ship inside matplotlib. Nothing here touches
the network.
"""

from __future__ import annotations

import matplotlib

matplotlib.use("Agg")

import matplotlib.pyplot as plt  # noqa: E402
import seaborn as sns  # noqa: E402

Color = tuple[float, float, float]


def _rgb(c) -> Color:
    return (float(c[0]), float(c[1]), float(c[2]))


# matplotlib's default categorical cycle, positions 3 (red) and 2 (green).
PASS_FAIL_RED: Color = _rgb(plt.get_cmap("tab10")(3))
PASS_FAIL_GREEN: Color = _rgb(plt.get_cmap("tab10")(2))

# seaborn's own colourblind-safe categorical palette, first two entries.
_CB = sns.color_palette("colorblind")
SAFE_BLUE: Color = _rgb(_CB[0])
SAFE_ORANGE: Color = _rgb(_CB[1])


def viridis_steps(n: int = 5) -> list[Color]:
    """`n` evenly spaced samples of matplotlib's sequential default."""
    cmap = plt.get_cmap("viridis")
    return [_rgb(cmap(i / (n - 1))) for i in range(n)]


def tab10_steps(n: int = 5) -> list[Color]:
    """The first `n` entries of matplotlib's categorical default."""
    cmap = plt.get_cmap("tab10")
    return [_rgb(cmap(i)) for i in range(n)]


# An ordered variable with five levels. Its ORDER is the whole point: a
# palette that does not preserve it has destroyed information the data
# had, which is a different and worse failure than merely looking bad.
SATISFACTION_LEVELS: tuple[str, ...] = (
    "very dissatisfied",
    "dissatisfied",
    "neutral",
    "satisfied",
    "very satisfied",
)

"""Rendering, and measuring what was actually rendered.

"That chart looks better" is not testable. "That chart spends 37% of its
ink on the data where this one spends 93%, and both carry the same eight
numbers" is. Everything in this module exists to turn a claim about a
picture into a number a test can assert on.

Every function here renders through matplotlib's **Agg** backend, which
draws into a memory buffer and needs no display, no window server and no
X11 forwarding. `matplotlib.use("Agg")` is called BEFORE `pyplot` is
imported, because the backend is chosen at import time and switching
afterwards is unreliable. `plt.show()` is never called, and every figure
is closed on the way out so a long test run does not leak them.

Every function takes an explicit output path. Nothing here writes to a
default location, so a test that hands it a `tmp_path` leaves nothing on
disk once pytest cleans up.
"""

from __future__ import annotations

from pathlib import Path

import matplotlib

matplotlib.use("Agg")  # must precede the pyplot import

import matplotlib.pyplot as plt  # noqa: E402
import numpy as np  # noqa: E402
from matplotlib.patches import Circle  # noqa: E402
from PIL import Image  # noqa: E402

WHITE = (255, 255, 255)

# The flat colour the bars are drawn in, as an sRGB 8-bit triple, so a
# test can isolate the data ink from the furniture.
BAR_RGB: tuple[int, int, int] = (0x1D, 0x4E, 0xD8)
BAR_HEX = "#%02x%02x%02x" % BAR_RGB

# The eight-region growth figures the lesson opens with. Deliberately
# close together at the top: 18.9 against 17.4 is a gap a sorted bar chart
# shows instantly and a pie chart cannot resolve at all.
REGION_GROWTH: dict[str, float] = {
    "Nordics": 12.1,
    "Iberia": 18.9,
    "Benelux": 7.4,
    "DACH": 17.4,
    "France": 9.8,
    "Italy": 4.2,
    "Poland": 15.6,
    "Ireland": 11.3,
}


# --------------------------------------------------------------------------
# Pixel measurement
# --------------------------------------------------------------------------


def _load_rgb(path: Path) -> np.ndarray:
    """Load a PNG as an (H, W, 3) uint8 array, discarding any alpha."""
    with Image.open(path) as im:
        return np.asarray(im.convert("RGB"))


def count_non_background_pixels(path: Path, background: tuple[int, int, int] = WHITE) -> int:
    """Count pixels that are not the background colour. This is the ink."""
    arr = _load_rgb(path)
    bg = np.array(background, dtype=np.uint8)
    return int(np.count_nonzero(np.any(arr != bg, axis=-1)))


def count_pixels_of_color(path: Path, rgb: tuple[int, int, int]) -> int:
    """Count pixels exactly matching one colour.

    The bars in `render_region_bar_chart` are drawn in one flat colour and
    nothing else in either figure uses it, so this isolates the DATA ink
    from every other mark on the page -- which is what makes Tufte's
    data-ink ratio a number this lab can measure rather than assert.
    """
    arr = _load_rgb(path)
    target = np.array(rgb, dtype=np.uint8)
    return int(np.count_nonzero(np.all(arr == target, axis=-1)))


def data_ink_ratio(path: Path, data_rgb: tuple[int, int, int]) -> float:
    """Data ink divided by total ink, both counted in pixels.

    Tufte's definition, made literal: of every mark on the page, what
    fraction is the data itself? Erasing non-data ink -- gridlines, a
    tinted panel, a heavy box -- raises this without removing a single
    fact from the chart.
    """
    total = count_non_background_pixels(path)
    if total == 0:
        raise ValueError("the image has no ink at all; a data-ink ratio is undefined")
    return count_pixels_of_color(path, data_rgb) / total


def count_distinct_luminance_levels(path: Path) -> int:
    """Count how many distinct grey levels the image contains.

    This is the measurable trace of density information. Opaque marks that
    overlap produce exactly two levels -- paper and ink -- no matter how
    many marks landed on the same pixel, because the tenth mark on a pixel
    changes nothing. Semi-transparent marks accumulate, so the number of
    levels tells you the image is carrying how MANY marks landed, not just
    whether any did.
    """
    arr = _load_rgb(path).astype(np.float64)
    lum = 0.2126 * arr[..., 0] + 0.7152 * arr[..., 1] + 0.0722 * arr[..., 2]
    return int(np.unique(np.round(lum).astype(np.int64)).size)


# --------------------------------------------------------------------------
# Exercise 1 -- circles drawn at a known radius, then measured
# --------------------------------------------------------------------------


def render_circle(path: Path, radius_px: float, canvas_px: int = 400) -> None:
    """Draw one filled circle of exactly `radius_px` pixels on white.

    The axes fill the whole figure and the data limits are set equal to
    the pixel dimensions, so one data unit is one pixel exactly and the
    circle's radius in data units is its radius on screen. Antialiasing is
    switched off so the edge is a hard boundary and the pixel count is a
    measurement of area rather than of area plus a soft fringe.
    """
    dpi = 100
    fig = plt.figure(figsize=(canvas_px / dpi, canvas_px / dpi), dpi=dpi, facecolor="white")
    ax = fig.add_axes((0, 0, 1, 1))
    ax.set_xlim(0, canvas_px)
    ax.set_ylim(0, canvas_px)
    ax.set_axis_off()
    ax.add_patch(
        Circle(
            (canvas_px / 2, canvas_px / 2),
            radius_px,
            facecolor="black",
            edgecolor="none",
            antialiased=False,
        )
    )
    fig.savefig(path, dpi=dpi, facecolor="white")
    plt.close(fig)


def measure_circle_area_px(path: Path) -> int:
    """The circle's drawn area, in painted pixels."""
    return count_non_background_pixels(path)


# --------------------------------------------------------------------------
# Exercise 7 -- the same chart, decorated and undecorated
# --------------------------------------------------------------------------


def render_region_bar_chart(path: Path, decorated: bool) -> None:
    """Render the eight-region bar chart with or without the furniture.

    Identical data, identical figure size, identical bars. The only
    difference is the non-data ink: a tinted plot background, a full box
    of spines, gridlines on both axes, and a heavy frame. Everything the
    reader needs -- the eight labels, the eight lengths, a common baseline
    -- is present in BOTH.
    """
    dpi = 100
    fig = plt.figure(figsize=(6, 4), dpi=dpi, facecolor="white")
    ax = fig.add_subplot(111)

    names = list(REGION_GROWTH)
    values = [REGION_GROWTH[n] for n in names]
    order = sorted(range(len(values)), key=lambda i: values[i])
    names = [names[i] for i in order]
    values = [values[i] for i in order]

    ax.barh(names, values, color=BAR_HEX)
    ax.set_xlabel("growth %")

    if decorated:
        ax.set_facecolor("#e2e8f0")
        ax.grid(True, which="both", axis="both", color="#64748b", linewidth=1.0)
        ax.set_axisbelow(False)
        for spine in ax.spines.values():
            spine.set_visible(True)
            spine.set_linewidth(3.0)
            spine.set_color("#334155")
    else:
        ax.set_facecolor("white")
        ax.grid(False)
        for name, spine in ax.spines.items():
            spine.set_visible(name == "bottom")
        ax.tick_params(left=False)

    fig.tight_layout()
    fig.savefig(path, dpi=dpi, facecolor="white")
    plt.close(fig)


# --------------------------------------------------------------------------
# Exercise 8 -- overplotting, and two ways out of it
# --------------------------------------------------------------------------


def sample_points(n: int = 10_000, seed: int = 127) -> tuple[np.ndarray, np.ndarray]:
    """A fixed, seeded cloud of `n` correlated points.

    Seeded, so every machine plots the same cloud and every pixel count in
    this lab is reproducible rather than merely typical.
    """
    rng = np.random.default_rng(seed)
    x = rng.normal(0.0, 1.0, n)
    y = 0.6 * x + rng.normal(0.0, 0.8, n)
    return x, y


# The scatter canvas. Both limits are wide enough that every one of the
# 10,000 sampled points falls inside the axes -- nothing is clipped, so
# the painted-pixel count can be compared against the full point count
# without an "and some were off the edge" caveat.
SCATTER_INCHES = 3
SCATTER_LIMIT = 5


def render_scatter(path: Path, x: np.ndarray, y: np.ndarray, alpha: float) -> None:
    """Plot the cloud with one-pixel marks at the given opacity.

    The `,` marker is matplotlib's single-pixel marker and antialiasing is
    off, so each point paints exactly one pixel. That makes the painted
    pixel count a direct measurement of how many DISTINCT screen positions
    the data occupies, with no marker-size confound: any shortfall against
    the point count is overplotting and nothing else.
    """
    dpi = 100
    fig = plt.figure(figsize=(SCATTER_INCHES, SCATTER_INCHES), dpi=dpi, facecolor="white")
    ax = fig.add_axes((0, 0, 1, 1))
    ax.set_xlim(-SCATTER_LIMIT, SCATTER_LIMIT)
    ax.set_ylim(-SCATTER_LIMIT, SCATTER_LIMIT)
    ax.set_axis_off()
    ax.plot(x, y, ",", color="black", alpha=alpha, antialiased=False, linestyle="none")
    fig.savefig(path, dpi=dpi, facecolor="white")
    plt.close(fig)


def render_hexbin(path: Path, x: np.ndarray, y: np.ndarray, gridsize: int = 30) -> None:
    """Plot the same cloud as a hexagonal density map.

    Where the scatter throws density away by painting the same pixel over
    and over, hexbin counts the points per cell and encodes the count as
    luminance -- the information the scatter destroyed, recovered by
    aggregating before drawing instead of after.
    """
    dpi = 100
    lim = SCATTER_LIMIT
    fig = plt.figure(figsize=(SCATTER_INCHES, SCATTER_INCHES), dpi=dpi, facecolor="white")
    ax = fig.add_axes((0, 0, 1, 1))
    ax.set_xlim(-lim, lim)
    ax.set_ylim(-lim, lim)
    ax.set_axis_off()
    ax.hexbin(x, y, gridsize=gridsize, cmap="Greys", extent=(-lim, lim, -lim, lim))
    fig.savefig(path, dpi=dpi, facecolor="white")
    plt.close(fig)


def count_painted_pixels(path: Path) -> int:
    """Distinct screen positions carrying at least one mark."""
    return count_non_background_pixels(path)


def points_inside_axes(x: np.ndarray, y: np.ndarray, limit: float = SCATTER_LIMIT) -> int:
    """How many of the sampled points fall inside the scatter's axes.

    Used to confirm that nothing is clipped, so a painted-pixel count
    below the point count means overplotting and never "the rest fell off
    the edge of the picture".
    """
    return int(np.count_nonzero((np.abs(x) <= limit) & (np.abs(y) <= limit)))

"""YOUR test suite for Day 127 -- "Charts That Answer the Question".

Nine exercises. Run it from the lab directory, not from here:

    pytest starter -v

Every exercise below ends in a `pytest.skip(...)` line. pytest reports a
skip as `s` and moves on, so an unfinished suite still exits 0. Replace
each skip with real assertions -- deleting the skip line is part of the
exercise. `starter/00_brief.md` explains each exercise in full.

The four modules beside this file are the lab's instruments, not
exercises. Read them before you start; do not edit them:

  encoding.py  geometry and colour arithmetic -- the square law, sRGB to
               CIELAB, the deuteranopia transform, Spearman's rho
  charts.py    the two decision functions, `best_encoding` and
               `choose_chart`, plus the Cleveland-McGill ordering
  palettes.py  the swatches, taken from matplotlib and seaborn themselves
  render.py    everything that draws, and everything that measures a
               drawing in pixels

Every render must go into the `png_dir` fixture -- a temporary directory
outside the lab. Write a PNG anywhere else and `tests/run_tests.sh` will
catch it.

Assert measured values, and give a rendered measurement a tolerance: a
rasteriser puts a circle's edge a pixel or two either way, so
`pytest.approx(..., rel=0.02)` is honest where an exact equality would be
luck. Nothing in this lab depends on timing.
"""

from __future__ import annotations

import math

import pytest

import charts as C
import encoding as E
import palettes as PAL
import render as R

# --------------------------------------------------------------------------
# EXERCISE 1 -- the square law, measured. See 00_brief.md exercise 1.
#
# Check with:   pytest starter -v -k test_1
# --------------------------------------------------------------------------

VALUES = [50.0, 100.0]  # the second is exactly twice the first


def test_1_radius_encoding_squares_every_ratio():
    pytest.skip(
        "exercise 1a: call encoding.encoded_area_ratio(VALUES, mode='radius') and again with "
        "mode='area'. Assert the radius encoding gives 4.0 -- the square of the data ratio 2.0 -- "
        "and the area encoding gives 2.0, and that the two differ by a factor of exactly 2."
    )


def test_1_rendered_pixel_areas_confirm_the_square_law(png_dir):
    pytest.skip(
        "exercise 1b: render.render_circle three circles into png_dir -- radius 40, radius 80, "
        "and radius 40*sqrt(2) -- then measure each with render.measure_circle_area_px. Assert "
        "the 80px circle covers about 4x the pixels of the 40px one (rel=0.02) and the "
        "40*sqrt(2) one about 2x, and report both counts."
    )


# --------------------------------------------------------------------------
# EXERCISE 2 -- the Cleveland-McGill ranking as a decision function.
#
# Check with:   pytest starter -v -k test_2
# --------------------------------------------------------------------------


def test_2_ranking_is_in_cleveland_mcgill_order():
    pytest.skip(
        "exercise 2a: use charts.encoding_rank to assert the ordering position_common_scale < "
        "position_nonaligned_scales < length < angle_slope < area < volume < color_saturation, "
        "and that asking for the rank of 'hue' raises ValueError."
    )


def test_2_best_encoding_case_table():
    pytest.skip(
        "exercise 2b: build a list of ((data_type, task), expected_channel) cases and assert "
        "charts.best_encoding returns each one. Cover at least: quantitative/compare, "
        "quantitative/compare_across_panels, quantitative/magnitude_on_map, "
        "ordinal/encode_in_color, nominal/identify_group and nominal/compare. Justify each "
        "expected answer in a comment, and assert explicitly that ordinal/encode_in_color is NOT "
        "'hue'."
    )


def test_2_best_encoding_rejects_what_it_does_not_understand():
    pytest.skip(
        "exercise 2c: assert charts.best_encoding raises ValueError matching 'unknown data type' "
        "for a bad data type and 'unknown task' for a bad task."
    )


# --------------------------------------------------------------------------
# EXERCISE 3 -- from the question to the chart.
#
# Check with:   pytest starter -v -k test_3
# --------------------------------------------------------------------------


def test_3_choose_chart_case_table():
    pytest.skip(
        "exercise 3a: build a ((question_kind, n_categories, data_types), expected) case table "
        "and assert charts.choose_chart returns each one. Cover all five question kinds, and "
        "include at least one case on each side of TABLE_MAX_VALUES, OVERPLOT_POINT_LIMIT and "
        "SMALL_MULTIPLE_LIMIT."
    )


def test_3_choose_chart_never_recommends_a_pie():
    pytest.skip(
        "exercise 3b: collect choose_chart's answer over every question kind and a spread of "
        "n_categories into a set, and assert neither 'pie' nor 'donut' is in it. Then assert "
        "that for ranking specifically the answer is always 'sorted_horizontal_bar'."
    )


def test_3_choose_chart_validates_its_inputs():
    pytest.skip(
        "exercise 3c: assert choose_chart raises ValueError for an unknown question kind, an "
        "unknown data type, n_categories of 0, and change_over_time with no temporal variable."
    )


# --------------------------------------------------------------------------
# EXERCISE 4 -- colour deficiency, simulated and measured.
#
# Check with:   pytest starter -v -k test_4
# --------------------------------------------------------------------------

COLLAPSE_THRESHOLD = 10.0
SURVIVAL_THRESHOLD = 25.0


def test_4_red_green_pair_collapses_under_deuteranopia():
    pytest.skip(
        "exercise 4a: call encoding.deuteranopia_collapse(PAL.PASS_FAIL_RED, PAL.PASS_FAIL_GREEN). "
        "Assert the normal-vision CIE76 distance is above 100, the simulated distance is below "
        "COLLAPSE_THRESHOLD, and the retained fraction is below 0.10. Report all three numbers."
    )


def test_4_colorblind_safe_pair_survives_the_same_transform():
    pytest.skip(
        "exercise 4b: run the same measurement on PAL.SAFE_BLUE and PAL.SAFE_ORANGE. Assert the "
        "simulated distance stays above SURVIVAL_THRESHOLD and the retained fraction above 0.90, "
        "and that the safe pair's simulated distance is more than 10x the red/green pair's."
    )


# --------------------------------------------------------------------------
# EXERCISE 5 -- an ordered variable on a categorical palette.
#
# Check with:   pytest starter -v -k test_5
# --------------------------------------------------------------------------


def test_5_sequential_palette_preserves_the_order():
    pytest.skip(
        "exercise 5a: take PAL.viridis_steps(5), compute encoding.relative_luminance of each, "
        "and assert the list is already sorted. Then assert "
        "encoding.luminance_order_correlation(palette) is approximately 1.0."
    )


def test_5_categorical_palette_destroys_the_order():
    pytest.skip(
        "exercise 5b: do the same with PAL.tab10_steps(5). Assert the luminance list is NOT "
        "sorted, that the rank correlation is well below the sequential palette's in absolute "
        "value, and report the number you measure."
    )


# --------------------------------------------------------------------------
# EXERCISE 6 -- sorting is an encoding decision.
#
# Check with:   pytest starter -v -k test_6
# --------------------------------------------------------------------------


def test_6_sorting_changes_the_effort_not_the_answer():
    pytest.skip(
        "exercise 6: build a list of 20 unsorted values. Assert "
        "charts.comparisons_to_find_max(values, presented_sorted=False) is 19 and with "
        "presented_sorted=True is 1, and that sorted(values, reverse=True)[0] equals "
        "values[charts.index_of_max(values)] -- sorting moved the effort, not the answer."
    )


# --------------------------------------------------------------------------
# EXERCISE 7 -- the data-ink ratio, counted in pixels.
#
# Check with:   pytest starter -v -k test_7
# --------------------------------------------------------------------------


def test_7_removing_furniture_raises_the_data_ink_ratio(png_dir):
    pytest.skip(
        "exercise 7: render.render_region_bar_chart into png_dir twice, decorated=True and "
        "decorated=False. Count total ink with render.count_non_background_pixels and the "
        "data-ink fraction with render.data_ink_ratio(path, render.BAR_RGB). Assert the plain "
        "chart's ratio is higher, and that the gap is more than 0.5. Report both counts and "
        "both ratios."
    )


# --------------------------------------------------------------------------
# EXERCISE 8 -- overplotting, and two ways out of it.
#
# Check with:   pytest starter -v -k test_8
# --------------------------------------------------------------------------

N_POINTS = 10_000


def test_8_overplotting_hides_a_third_of_the_data(png_dir, points):
    pytest.skip(
        "exercise 8a: confirm render.points_inside_axes says nothing is clipped, render the "
        "cloud with alpha=1.0, and assert the painted-pixel count is below 75% of N_POINTS. "
        "Then assert render.count_distinct_luminance_levels of that image is exactly 2 -- the "
        "density information is not dimmed, it is absent."
    )


def test_8_alpha_and_hexbin_recover_the_density(png_dir, points):
    pytest.skip(
        "exercise 8b: render the same cloud three ways -- alpha=1.0, alpha=0.05, and "
        "render.render_hexbin -- and count distinct luminance levels in each. Assert the alpha "
        "version has more than the opaque one and the hexbin more than 50. Finish by asserting "
        "charts.choose_chart('relationship', N_POINTS, ['quantitative']) == 'hexbin'."
    )


# --------------------------------------------------------------------------
# EXERCISE 9 -- when a table beats a chart.
#
# Check with:   pytest starter -v -k test_9
# --------------------------------------------------------------------------


def test_9_a_table_beats_a_chart_below_the_threshold():
    pytest.skip(
        "exercise 9: assert choose_chart('comparison', 3, ...) is 'table' and "
        "choose_chart('comparison', 30, ...) is 'sorted_horizontal_bar'. Assert the boundary "
        "sits exactly at charts.TABLE_MAX_VALUES by testing that value and that value plus one. "
        "Then write a comment explaining WHY the boundary is where it is -- what a chart buys "
        "you, and why three numbers do not need it."
    )

#!/usr/bin/env bash
# Tests for the Day 127 lab. Run from the lab directory:
#   bash tests/run_tests.sh
#
# A visualisation lab cannot assert that a chart "looks better", so this
# one asserts only what is genuinely measurable, by running code and
# reading real values:
#
#   * encoding a value as a circle's RADIUS squares every ratio in the
#     chart -- confirmed analytically (4.0 against a data ratio of 2.0)
#     and again by counting the pixels of two rendered circles;
#   * the Cleveland-McGill ordering, used as a decision function, and a
#     chart-choice function that recommends a TABLE below a stated number
#     of values and never recommends a pie chart for anything;
#   * matplotlib's default red and green -- the pass/fail reflex -- start
#     119.8 apart in CIELAB and end 7.3 apart under a published
#     deuteranopia transform, while seaborn's colourblind-safe blue and
#     orange keep essentially all of their separation;
#   * a sequential palette's luminance order matches an ordinal
#     variable's order exactly (rank correlation 1.0) while a categorical
#     palette's does not;
#   * sorting turns 19 reader comparisons into 1 without changing the
#     answer;
#   * the same eight numbers drawn with and without furniture: 37% of the
#     decorated chart's ink is data against 93% of the plain one's;
#   * 10,000 one-pixel points paint only 6,349 distinct pixels, and the
#     opaque image contains exactly TWO grey levels -- density is not
#     dimmed, it is absent -- which alpha blending and hexbin recover;
#   * the reference suite (`examples/`) passes in full;
#   * the exercise suite (`starter/`) is all-skip on an untouched
#     checkout, and the harness proves it can genuinely FAIL by solving
#     every exercise in a scratch copy, breaking one assertion on
#     purpose, confirming a non-zero exit, then restoring it;
#   * matplotlib really is headless (Agg) and nothing calls plt.show();
#   * nothing -- no .png, no __pycache__ -- is left behind by this run.
#
# Everything after the one-time install runs offline. Nothing binds a
# port, nothing writes outside the lab or a temporary directory, nothing
# needs a key. Deterministic, non-interactive, exits 0 only if every
# check passes.
set -u

export PYTHONDONTWRITEBYTECODE=1

# MPLBACKEND is deliberately NOT set here. render.py calls
# matplotlib.use("Agg") itself, and section 2 below checks that it really
# did -- pre-setting the environment variable would make that check pass
# for the wrong reason.
unset MPLBACKEND

lab_dir="$(cd "$(dirname "${BASH_SOURCE[0]}")/.." && pwd)"

# Bytecode left by an EARLIER command is not this run's litter. The
# README documents `pytest starter` and `pytest examples` as separate
# commands, and running either writes .pyc files that would then fail the
# cleanliness check at the end -- failing the reader for following the
# instructions. Clearing them here makes that final check measure what it
# claims to. `.venv` is untouched: the packages' bytecode is theirs.
find "${lab_dir}" -name '.venv' -prune -o -type d -name '__pycache__' -exec rm -rf {} + 2>/dev/null || true
find "${lab_dir}" -name '.venv' -prune -o -type d -name '.pytest_cache' -exec rm -rf {} + 2>/dev/null || true

failures=0
checks=0

check() {
  local label="$1" ok="$2"
  checks=$((checks + 1))
  if [ "${ok}" = "yes" ]; then
    echo "  ok: ${label}"
  else
    echo "  FAIL: ${label}"
    failures=$((failures + 1))
  fi
}

check_eq() {
  # check_eq <label> <expected> <actual>
  if [ "$2" = "$3" ]; then
    check "$1" "yes"
  else
    check "$1 (expected [$2], got [$3])" "no"
  fi
}

# Resolve pytest: an explicit override, then this lab's .venv, then PATH.
# Fails loudly with instructions rather than silently skipping checks.
resolve_tool() {
  local tool="$1" override="$2"
  if [ -n "${override}" ] && [ -x "${override}" ]; then echo "${override}"; return 0; fi
  if [ -x "${lab_dir}/.venv/bin/${tool}" ]; then echo "${lab_dir}/.venv/bin/${tool}"; return 0; fi
  if command -v "${tool}" >/dev/null 2>&1; then command -v "${tool}"; return 0; fi
  return 1
}

pytest_bin="$(resolve_tool pytest "${PYTEST:-}")" || {
  echo "FAIL: pytest not found." >&2
  echo "  Install the lab's dependencies with:" >&2
  echo "    python3 -m venv .venv" >&2
  echo "    .venv/bin/pip install -r requirements/requirements.txt" >&2
  echo "  Or point this suite at an existing pytest:" >&2
  echo "    PYTEST=/path/to/pytest bash tests/run_tests.sh" >&2
  exit 1
}

python_bin="$(dirname "${pytest_bin}")/python3"
if [ ! -x "${python_bin}" ]; then
  python_bin="$(command -v python3 || true)"
fi
if [ -z "${python_bin}" ]; then
  echo "FAIL: python3 not found on PATH." >&2
  exit 1
fi

if ! "${python_bin}" -c "import matplotlib, seaborn, PIL" >/dev/null 2>&1; then
  echo "FAIL: matplotlib, seaborn or Pillow is not importable from ${python_bin}." >&2
  echo "  Install the lab's dependencies with:" >&2
  echo "    python3 -m venv .venv" >&2
  echo "    .venv/bin/pip install -r requirements/requirements.txt" >&2
  exit 1
fi

echo "Day 127 — Charts That Answer the Question"
echo

# --------------------------------------------------------------------------
echo "1. The tools and the versions this lab was written against"
# --------------------------------------------------------------------------

versions="$("${python_bin}" - <<'PY'
import platform
from importlib.metadata import version

print(f"python     {platform.python_version()}")
for name in ("matplotlib", "seaborn", "pandas", "numpy", "pillow", "pytest"):
    try:
        print(f"{name:<10} {version(name)}")
    except Exception as exc:  # pragma: no cover
        print(f"{name:<10} NOT INSTALLED ({exc})")
PY
)"
echo "${versions}"
echo

mpl_version="$("${python_bin}" -c "import matplotlib; print(matplotlib.__version__)" 2>/dev/null || echo "")"
pinned_mpl="$(grep -m1 '^matplotlib==' "${lab_dir}/requirements/requirements.txt" | cut -d= -f3)"
check_eq "installed matplotlib matches requirements.txt exactly" "${pinned_mpl}" "${mpl_version}"

sns_version="$("${python_bin}" -c "import seaborn; print(seaborn.__version__)" 2>/dev/null || echo "")"
pinned_sns="$(grep -m1 '^seaborn==' "${lab_dir}/requirements/requirements.txt" | cut -d= -f3)"
check_eq "installed seaborn matches requirements.txt exactly" "${pinned_sns}" "${sns_version}"
echo

# --------------------------------------------------------------------------
echo "2. Rendering is headless -- Agg, no display, no window server"
# --------------------------------------------------------------------------

backend="$(cd "${lab_dir}/examples" && "${python_bin}" -c "import render, matplotlib; print(matplotlib.get_backend().lower())" 2>/dev/null || echo "")"
check_eq "importing render.py selects the Agg backend" "agg" "${backend}"

# Anchored to the start of a statement on purpose. An unanchored search
# also matches the prose in render.py's docstring that PROMISES plt.show
# is never called -- and a check that fails because the code documents
# itself is a check measuring the wrong thing. (Observed here: the first
# run of this harness reported exactly that failure, which is also the
# first proof this harness can fail.)
show_hits="$(grep -rnE '^[[:space:]]*plt\.show\(' "${lab_dir}/examples" "${lab_dir}/starter" 2>/dev/null || true)"
check "nothing calls plt.show() -- it would hang a headless run" "$( [ -z "${show_hits}" ] && echo yes || echo no )"
echo

# --------------------------------------------------------------------------
echo "3. Reference suite -- examples/ must pass in full"
# --------------------------------------------------------------------------

examples_output="$(cd "${lab_dir}" && "${pytest_bin}" examples -q 2>&1)"
examples_status=$?
echo "${examples_output}" | tail -5
check "examples/ exits 0" "$( [ ${examples_status} -eq 0 ] && echo yes || echo no )"

check "examples/ reports 17 passed, 0 failed" "$( echo "${examples_output}" | grep -qE '^17 passed' && echo yes || echo no )"
echo

# --------------------------------------------------------------------------
echo "4. Exercise suite -- starter/ is all-skip on an untouched checkout"
# --------------------------------------------------------------------------

starter_output="$(cd "${lab_dir}" && "${pytest_bin}" starter -q 2>&1)"
starter_status=$?
echo "${starter_output}" | tail -5
check "starter/ (untouched) exits 0" "$( [ ${starter_status} -eq 0 ] && echo yes || echo no )"
check "starter/ (untouched) reports 17 skipped, 0 failed" "$( echo "${starter_output}" | grep -qE '^17 skipped' && echo yes || echo no )"
echo

# --------------------------------------------------------------------------
echo "5. Never run 'pytest examples starter' in one invocation -- every"
echo "   module name (encoding, charts, palettes, render, conftest,"
echo "   test_charts) is defined identically in both directories, so the"
echo "   second collected collides with the first. Checked below."
# --------------------------------------------------------------------------

both_output="$(cd "${lab_dir}" && "${pytest_bin}" examples starter -q 2>&1)"
both_status=$?
check "pytest examples starter (one invocation) does NOT exit 0" "$( [ ${both_status} -ne 0 ] && echo yes || echo no )"
check "pytest examples starter reports an import file mismatch, not a quiet partial run" "$( echo "${both_output}" | grep -qi 'import file mismatch' && echo yes || echo no )"
echo

# --------------------------------------------------------------------------
echo "6. Prove the suite can genuinely FAIL: solve every exercise in a"
echo "   scratch copy, confirm green, break one assertion on purpose,"
echo "   confirm a non-zero exit and a printed failure, then restore."
# --------------------------------------------------------------------------

scratch_dir="$(mktemp -d "${TMPDIR:-/tmp}/d127-scratch.XXXXXX")"
cleanup_scratch() { rm -rf "${scratch_dir}"; }
trap cleanup_scratch EXIT

for module in test_charts encoding charts palettes render conftest; do
  cp "${lab_dir}/examples/${module}.py" "${scratch_dir}/${module}.py"
done

solved_output="$("${pytest_bin}" "${scratch_dir}" -q 2>&1)"
solved_status=$?
check "scratch copy of the solved suite exits 0" "$( [ ${solved_status} -eq 0 ] && echo yes || echo no )"
check "scratch copy reports 17 passed" "$( echo "${solved_output}" | grep -qE '^17 passed' && echo yes || echo no )"

# Break exercise 8's exact luminance-level assertion on purpose. Two grey
# levels is the measured truth; three is not.
sed -i.bak 's/assert R\.count_distinct_luminance_levels(opaque) == 2/assert R.count_distinct_luminance_levels(opaque) == 3/' "${scratch_dir}/test_charts.py"

broken_output="$("${pytest_bin}" "${scratch_dir}" -q 2>&1)"
broken_status=$?
check "broken scratch copy exits non-zero" "$( [ ${broken_status} -ne 0 ] && echo yes || echo no )"
check "broken scratch copy prints a failure line" "$( echo "${broken_output}" | grep -qiE 'failed|assert' && echo yes || echo no )"

mv "${scratch_dir}/test_charts.py.bak" "${scratch_dir}/test_charts.py"
restored_output="$("${pytest_bin}" "${scratch_dir}" -q 2>&1)"
restored_status=$?
check "restored scratch copy exits 0 again" "$( [ ${restored_status} -eq 0 ] && echo yes || echo no )"
check "restored scratch copy reports 17 passed again" "$( echo "${restored_output}" | grep -qE '^17 passed' && echo yes || echo no )"

cleanup_scratch
trap - EXIT
echo

# --------------------------------------------------------------------------
echo "7. Nothing in examples/ or starter/ opens a network connection"
# --------------------------------------------------------------------------

url_hits="$(grep -rEl 'https?://|ftp://' "${lab_dir}/examples" "${lab_dir}/starter" 2>/dev/null || true)"
check "no URLs inside examples/ or starter/" "$( [ -z "${url_hits}" ] && echo yes || echo no )"
echo

# --------------------------------------------------------------------------
echo "8. A chart-rendering day that litters images would be embarrassing"
echo "   -- confirm no .png, .jpg, .svg or .pdf is left anywhere in the lab"
# --------------------------------------------------------------------------

image_hits="$(find "${lab_dir}" -name '.venv' -prune -o \
  \( -type f \( -name '*.png' -o -name '*.jpg' -o -name '*.jpeg' -o -name '*.svg' -o -name '*.pdf' \) -print \) 2>/dev/null || true)"
check "no image files left under the lab" "$( [ -z "${image_hits}" ] && echo yes || echo no )"
echo

# --------------------------------------------------------------------------
echo "9. Cleanliness -- nothing left behind by THIS run"
# --------------------------------------------------------------------------

find "${lab_dir}" -name '.venv' -prune -o -type d -name '__pycache__' -exec rm -rf {} + 2>/dev/null || true
find "${lab_dir}" -name '.venv' -prune -o -type d -name '.pytest_cache' -exec rm -rf {} + 2>/dev/null || true

stray="$(find "${lab_dir}" -name '.venv' -prune -o \( -type d -name '__pycache__' -print -o -type d -name '.pytest_cache' -print \) 2>/dev/null || true)"
check "no __pycache__ or .pytest_cache left behind" "$( [ -z "${stray}" ] && echo yes || echo no )"
echo

echo "-------------------------------------------------------------"
echo "${checks} checks, ${failures} failure(s)"
if [ "${failures}" -gt 0 ]; then
  exit 1
fi
exit 0