"""
UDEF 1.0 - reproducible synthetic benchmark and reference algorithm
Python >= 3.9; dependencies: numpy, pandas
IMPORTANT: This is a synthetic proof-of-concept, not industrial validation.
"""
from __future__ import annotations
import numpy as np
import pandas as pd

SEED = 20261105
N_PER_DOMAIN = 5000
DOMAINS = {
    "Bioprocess": {"scale": 1.0, "noise": 0.18},
    "Battery": {"scale": 1.4, "noise": 0.22},
    "Chemical": {"scale": 0.9, "noise": 0.15},
    "Logistics": {"scale": 1.2, "noise": 0.25},
}

def binary_metrics(y_true, y_pred):
    y_true = np.asarray(y_true, dtype=int)
    y_pred = np.asarray(y_pred, dtype=int)
    tp = np.sum((y_true == 1) & (y_pred == 1))
    fp = np.sum((y_true == 0) & (y_pred == 1))
    fn = np.sum((y_true == 1) & (y_pred == 0))
    tn = np.sum((y_true == 0) & (y_pred == 0))
    precision = tp / (tp + fp) if tp + fp else 0.0
    recall = tp / (tp + fn) if tp + fn else 0.0
    f1 = 2 * precision * recall / (precision + recall) if precision + recall else 0.0
    accuracy = (tp + tn) / (tp + tn + fp + fn)
    return {"precision": precision, "recall": recall, "f1": f1, "accuracy": accuracy}

def udef_decision(z_process, qc_invalidity, process_threshold=2.0, qc_threshold=0.5):
    """
    z_process: standardized absolute process residual.
    qc_invalidity: calibrated score in [0,1], higher means less trustworthy measurement.
    Returns independent process-anomaly and measurement-invalidity flags.
    """
    measurement_invalid = np.asarray(qc_invalidity) > qc_threshold
    process_anomaly = (np.asarray(z_process) > process_threshold) & (~measurement_invalid)
    return process_anomaly.astype(int), measurement_invalid.astype(int)

def run(seed=SEED, n_per_domain=N_PER_DOMAIN):
    rng = np.random.default_rng(seed)
    output = []
    for domain, cfg in DOMAINS.items():
        scale, noise = cfg["scale"], cfg["noise"]
        process_bad = rng.binomial(1, 0.25, n_per_domain)
        measurement_bad = rng.binomial(1, 0.20, n_per_domain)
        process_effect = process_bad * rng.normal(3.2*scale, 0.45*scale, n_per_domain)
        random_noise = rng.normal(0, noise*scale, n_per_domain)
        matrix_bias = measurement_bad * rng.normal(2.8*scale, 0.5*scale, n_per_domain)
        residual = process_effect + random_noise + matrix_bias
        z = np.abs(residual) / (noise*scale)
        qc = np.where(measurement_bad,
                      rng.normal(0.82, 0.13, n_per_domain),
                      rng.normal(0.18, 0.12, n_per_domain))
        qc = np.clip(qc, 0, 1)

        naive_process = (z > 2.0).astype(int)
        udef_process, udef_measurement = udef_decision(z, qc)

        for method, p_pred, m_pred in [
            ("Naive single-threshold", naive_process, np.zeros(n_per_domain, dtype=int)),
            ("UDEF dual-engine", udef_process, udef_measurement)
        ]:
            pm = binary_metrics(process_bad, p_pred)
            mm = binary_metrics(measurement_bad, m_pred)
            joint = np.mean((p_pred == process_bad) & (m_pred == measurement_bad))
            output.append({
                "Domain": domain, "Method": method,
                "Process precision": pm["precision"],
                "Process recall": pm["recall"], "Process F1": pm["f1"],
                "Measurement precision": mm["precision"],
                "Measurement recall": mm["recall"], "Measurement F1": mm["f1"],
                "Joint accuracy": joint
            })
    result = pd.DataFrame(output)
    print("Per-domain results")
    print(result.to_string(index=False, float_format=lambda x: f"{x:.3f}"))
    print("\\nMean across domains")
    print(result.groupby("Method").mean(numeric_only=True).to_string(float_format=lambda x: f"{x:.3f}"))
    result.to_csv("UDEF_synthetic_benchmark_results.csv", index=False)
    return result

if __name__ == "__main__":
    run()