stats.percentile

import { type PercentileMethod } from "./stats_percentile_types.ts";

const POW10 = [1, 10, 100, 1e3, 1e4, 1e5, 1e6, 1e7, 1e8, 1e9, 1e10, 1e11, 1e12];

// Half away from zero on the binary64 value, then -0 becomes 0. Only floor,
// * and / are used, so every language rounds the same double the same way.
function roundTo(x: number, decimals: number): number {
  const scale = POW10[decimals];
  const y = Math.abs(x) * scale;
  let r = Math.floor(y);
  if (y - r >= 0.5) r += 1;
  const out = r / scale;
  return (x < 0 ? -out : out) + 0;
}

/**
 * The p-th percentile of `values` by a named method.
 *
 * There is no single "95th percentile": the caller names the definition, so
 * a dashboard and a billing job cannot quietly disagree.
 */
export function percentile(values: readonly number[], p: number, method: PercentileMethod, decimals: number): number {
  if (!Array.isArray(values)) throw new TypeError("values must be a list of numbers");
  if (values.length === 0) throw new RangeError("values must not be empty");
  for (const v of values) {
    if (typeof v !== "number" || !Number.isFinite(v)) {
      throw new TypeError(`values must be finite numbers, received ${v}`);
    }
  }
  if (typeof p !== "number" || !Number.isFinite(p) || p < 0 || p > 100) {
    throw new RangeError(`p must be between 0 and 100, received ${p}`);
  }
  if (!Number.isInteger(decimals) || decimals < 0 || decimals > 12) {
    throw new RangeError(`decimals must be a whole number from 0 to 12, received ${decimals}`);
  }

  const sorted = [...values].sort((a, b) => a - b);
  const n = sorted.length;

  let result: number;
  if (method === "nearest-rank") {
    // The smallest value with at least p% of the sample at or below it.
    let rank = Math.ceil((p * n) / 100);
    if (rank < 1) rank = 1;
    result = sorted[rank - 1];
  } else if (method === "linear") {
    // Hyndman & Fan type 7: zero-based position (n - 1) * p / 100.
    const h = ((n - 1) * p) / 100;
    const lo = Math.floor(h);
    if (lo >= n - 1) {
      result = sorted[n - 1];
    } else {
      result = sorted[lo] + (h - lo) * (sorted[lo + 1] - sorted[lo]);
    }
  } else {
    throw new RangeError(`unknown percentile method "${method}"`);
  }
  return roundTo(result, decimals);
}