//! Substation coordinate parity between the `.pwd` display reader and the
//! aux siblings. The oracle is differential: each decoded substation set is
//! matched 1-1 by (number, name) against the Substation rows of the sibling
//! aux, and the diagram coordinates must correlate with the geography the
//! aux carries (the TAMU and Hawaii layouts are generated by projecting
//! latitude/longitude; hand moved symbols and the June 2016 era transform
//! keep the correlation just under exact). Hawaii40, never hand edited,
//! reproduces the projection bit for bit and pins its constant.
//! Evidence in docs/powerworld.md.

mod common;

use std::collections::BTreeMap;
use std::fs;
use std::path::Path;

use powerio::format::powerworld::{PwdSubstation, parse_aux, parse_pwd};

/// (number, name, latitude, longitude) per aux Substation row. Handles both
/// naming vocabularies (classic SubNum/SubName, 2022 Number/Name).
fn aux_substations(path: &Path) -> Vec<(u32, String, f64, f64)> {
    let aux = parse_aux(&fs::read_to_string(path).unwrap()).unwrap();
    let sub = aux
        .data_of("Substation")
        .next()
        .expect("aux has a Substation section");
    let idx = |names: &[&str]| {
        names
            .iter()
            .find_map(|n| sub.field_index(n))
            .unwrap_or_else(|| panic!("{names:?} not in {:?}", sub.fields))
    };
    let (num, name) = (idx(&["SubNum", "Number"]), idx(&["SubName", "Name"]));
    let (lat, lon) = (idx(&["Latitude"]), idx(&["Longitude"]));
    sub.rows
        .iter()
        .map(|r| {
            (
                r.values[num].parse().unwrap(),
                r.values[name].clone(),
                r.values[lat].parse().unwrap(),
                r.values[lon].parse().unwrap(),
            )
        })
        .collect()
}

/// Join the decoded substations 1-1 with the aux rows by number, assert the
/// names agree, and return (x, y, latitude, longitude) per substation.
fn join(pwd: &[PwdSubstation], aux: &[(u32, String, f64, f64)]) -> Vec<(f64, f64, f64, f64)> {
    assert_eq!(pwd.len(), aux.len(), "substation count");
    let by_num: BTreeMap<u32, &(u32, String, f64, f64)> = aux.iter().map(|r| (r.0, r)).collect();
    pwd.iter()
        .map(|s| {
            let a = by_num
                .get(&s.number)
                .unwrap_or_else(|| panic!("substation {} not in the aux", s.number));
            assert_eq!(s.name, a.1, "substation {} name", s.number);
            (s.x, s.y, a.2, a.3)
        })
        .collect()
}

/// Pearson r² between paired samples.
fn r_squared(pairs: &[(f64, f64)]) -> f64 {
    let n = pairs.len() as f64;
    let (mx, my) = (
        pairs.iter().map(|p| p.0).sum::<f64>() / n,
        pairs.iter().map(|p| p.1).sum::<f64>() / n,
    );
    let (mut sxy, mut sxx, mut syy) = (0.0, 0.0, 0.0);
    for (x, y) in pairs {
        sxy += (x - mx) * (y - my);
        sxx += (x - mx) * (x - mx);
        syy += (y - my) * (y - my);
    }
    sxy * sxy / (sxx * syy)
}

/// The Mercator latitude stretch the auto generated layouts use, in degrees.
fn merc(lat: f64) -> f64 {
    (45.0 + lat / 2.0).to_radians().tan().ln().to_degrees()
}

/// Both correlation checks for one decoded file: x against longitude and y
/// against Mercator latitude.
fn assert_projection_fit(joined: &[(f64, f64, f64, f64)], floor: f64, label: &str) {
    let rx = r_squared(&joined.iter().map(|j| (j.3, j.0)).collect::<Vec<_>>());
    let ry = r_squared(&joined.iter().map(|j| (merc(j.2), j.1)).collect::<Vec<_>>());
    assert!(rx > floor, "{label}: x vs longitude r² {rx}");
    assert!(ry > floor, "{label}: y vs Mercator latitude r² {ry}");
}

/// The vendored 200 bus pair: 111 substations, every one matched by number
/// and name, coordinates correlated with the aux geography.
#[test]
fn activsg200_pwd_matches_its_aux_sibling() {
    let pwd = parse_pwd(&fs::read(common::powerworld_vendored("ACTIVSg200.pwd")).unwrap()).unwrap();
    assert_eq!(pwd.len(), 111);
    let aux = aux_substations(&common::powerworld_vendored("ACTIVSg200.aux"));
    let joined = join(&pwd, &aux);
    assert_projection_fit(&joined, 0.999, "ACTIVSg200");
}

/// The June 2016 era save (no select by criteria block, the oldest probed
/// framing): 1500 substations against the same day aux.
#[test]
fn texas2000_june2016_pwd_matches_its_aux_sibling() {
    let (Some(pwd_path), Some(aux_path)) = (
        common::activsg2000_fetched("Texas2000_June2016.pwd"),
        common::activsg2000_fetched("Texas2000_June2016.AUX"),
    ) else {
        eprintln!("skipped: run benchmarks/fetch_powerworld.sh");
        return;
    };
    let pwd = parse_pwd(&fs::read(pwd_path).unwrap()).unwrap();
    assert_eq!(pwd.len(), 1500);
    let joined = join(&pwd, &aux_substations(&aux_path));
    // The 2016 writer projected with a different transform (linear in
    // latitude, offset axes), so the fit is looser than the Mercator eras.
    assert_projection_fit(&joined, 0.999, "Texas2000 June 2016");
}

/// The v19 resave decodes 1250 substations; its substation set differs from
/// every available aux by vintage skew (the known v19 delta family), so the
/// assertion is structural: count and the inline link consistency the
/// parser itself enforces.
#[test]
fn activsg2000_v19_pwd_decodes() {
    let Some(path) = common::activsg2000_fetched("ACTIV_SG_2000_v19.pwd") else {
        eprintln!("skipped: run benchmarks/fetch_powerworld.sh");
        return;
    };
    let pwd = parse_pwd(&fs::read(path).unwrap()).unwrap();
    assert_eq!(pwd.len(), 1250);
}

/// The machine specific corpus siblings (paths via the gitignored
/// manifest): the 2019+ era writers at three scales, including the bit
/// exact Hawaii40 projection.
#[test]
fn local_corpus_pwd_siblings_match_their_auxes() {
    let cases: &[(&str, &str, usize, &str, f64)] = &[
        (
            "ACTIVSg500 export (local only)",
            "ACTIVSg500.pwd",
            208,
            "ACTIVSg500.aux",
            0.9999,
        ),
        (
            "ACTIVSg2000 published set export (local only)",
            "ACTIVSg2000.pwd",
            1250,
            "ACTIVSg2000.aux",
            0.9999,
        ),
        (
            "Hawaii40 2022 export (local only)",
            "Hawaii40.pwd",
            31,
            "Hawaii40_220906.AUX",
            0.9999,
        ),
    ];
    for (label, pwd_name, count, aux_name, floor) in cases {
        let Some(pwb_path) = common::local_corpus_path(label) else {
            eprintln!("skipped {label}: not in the local corpus manifest");
            continue;
        };
        let (pwd_path, aux_path) = (
            pwb_path.with_file_name(pwd_name),
            pwb_path.with_file_name(aux_name),
        );
        if !pwd_path.exists() || !aux_path.exists() {
            eprintln!("skipped {label}: no .pwd/.aux siblings next to the export");
            continue;
        }
        let pwd = parse_pwd(&fs::read(pwd_path).unwrap()).unwrap();
        assert_eq!(pwd.len(), *count, "{label}");
        let joined = join(&pwd, &aux_substations(&aux_path));
        assert_projection_fit(&joined, *floor, label);
        if *label == "Hawaii40 2022 export (local only)" {
            // Never hand edited, so the layout IS the projection: both axes
            // reproduce x = k·longitude, y = k·merc(latitude) to f64
            // rounding and pin the scale constant.
            let k = 535.816_080_362_259_2;
            for (x, y, lat, lon) in &joined {
                assert!((x - k * lon).abs() < 1e-6, "x {x} vs {}", k * lon);
                assert!(
                    (y - k * merc(*lat)).abs() < 1e-6,
                    "y {y} vs {}",
                    k * merc(*lat)
                );
            }
        }
    }
}

/// The malformed input invariant the fuzz harness asserts: every truncation
/// or byte corruption of a real display file parses to Ok or a structured
/// Err, never a panic. Every offset the reader follows derives from file
/// bytes, so the sweep walks the whole file at a stride plus the tail,
/// where a clipped record leaves reads dangling past the end.
#[test]
fn corrupt_display_files_never_panic() {
    let pwd = fs::read(common::powerworld_vendored("ACTIVSg200.pwd")).unwrap();
    let stride = 9973;
    for len in (0..pwd.len())
        .step_by(stride)
        .chain(pwd.len().saturating_sub(48)..pwd.len())
    {
        let _ = parse_pwd(&pwd[..len]);
    }
    let mut copy = pwd.clone();
    for i in (0..pwd.len()).step_by(stride) {
        for byte in [0x00, 0xff] {
            copy[i] = byte;
            let _ = parse_pwd(&copy);
        }
        copy[i] = pwd[i];
    }
}

/// Loud rejections: a binary case is not a display file, garbage is not a
/// display file, and a display header with no substation identity table
/// names what is missing.
#[test]
fn rejects_non_display_inputs() {
    let pwb = fs::read(common::powerworld_vendored("ACTIVSg200.pwb")).unwrap();
    let err = parse_pwd(&pwb).unwrap_err().to_string();
    assert!(
        err.contains("not a recognized PowerWorld display file"),
        "{err}"
    );

    let err = parse_pwd(&[0u8; 16]).unwrap_err().to_string();
    assert!(
        err.contains("not a recognized PowerWorld display file"),
        "{err}"
    );

    let mut headless = Vec::new();
    headless.extend_from_slice(&50u32.to_le_bytes());
    headless.extend_from_slice(&200u16.to_le_bytes());
    headless.extend_from_slice(&200u16.to_le_bytes());
    headless.extend_from_slice(&[0u8; 14]);
    headless.extend_from_slice(&0xa83cu32.to_le_bytes()); // stamp at offset 22
    headless.extend_from_slice(&[0u8; 4096]);
    let err = parse_pwd(&headless).unwrap_err().to_string();
    assert!(err.contains("no substation identity table"), "{err}");
}
