//! Minimal VCF allele extraction and left-normalisation.
//!
//! Converts full-sequence path representations (100+ bp ref and alt) to the
//! minimal VCF allele form required for standard variant comparison. For SNVs
//! this trims the flanking shared context. For indels it additionally applies
//! VCF left-normalisation so that deletions or insertions in repetitive regions
//! are always anchored as far left as possible — matching the convention used by
//! bcftools, GATK, and standard truth sets.
//!
//! # Coordinate convention
//!
//! Target IDs in the form `chrN:START-END` use a 1-based closed interval where
//! START is the genomic position of the first base in the stored sequence.
//! A base at sequence index `i` (0-based) is therefore at VCF position
//! `START + i` (1-based).

use std::fmt;

use crate::targeting::parse_target_id;

// ─── Public types ─────────────────────────────────────────────────────────────

/// A minimal, left-normalised variant allele ready for VCF output.
#[derive(Debug, Clone, PartialEq, Eq)]
pub struct MinimalAllele {
    /// Chromosome or sequence name.
    pub chrom: String,
    /// 1-based genomic position of the first base of the REF allele.
    pub pos: u64,
    /// Reference allele string.
    pub ref_allele: String,
    /// Alternate allele string.
    pub alt_allele: String,
}

impl fmt::Display for MinimalAllele {
    fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
        write!(
            f,
            "{}:{} {}->{}",
            self.chrom, self.pos, self.ref_allele, self.alt_allele
        )
    }
}

// ─── Public functions ─────────────────────────────────────────────────────────

/// Extract a minimal left-normalised VCF allele from full target sequences.
///
/// `target_id` must match the format `chrN:START-END` where START is the
/// 1-based genomic position of the first base in `ref_seq`. Returns `None`
/// if `target_id` cannot be parsed or if `ref_seq` and `alt_seq` are identical.
///
/// # Algorithm
///
/// 1. Find the first position where `ref_seq` and `alt_seq` differ.
/// 2. Trim the common suffix of both full sequences to isolate the minimal
///    differing region.
/// 3. For indels: strip any residual inner common prefix and suffix of the
///    trimmed region to isolate the inserted or deleted bases.
/// 4. Left-normalise: shift the anchor left while the preceding base in the
///    reference matches the last base of the inserted/deleted sequence.
/// 5. Prepend the anchor base as required by the VCF specification.
///
/// # Example
///
/// ```
/// use kam_call::allele::extract_minimal_allele;
///
/// // Deletion of one G from a "GG" run.
/// // ref: ...T(50) G(51) G(52) T(53)...
/// // alt: ...T(50) G(51) T(52)...     (G at index 51 deleted)
/// // First diff is at index 52; left-normalisation shifts to anchor at index 50.
/// let ref_seq: Vec<u8> = (0u8..101).map(|i| if i == 51 || i == 52 { b'G' }
///     else if i == 50 || i == 53 { b'T' } else { b'A' }).collect();
/// let alt_seq: Vec<u8> = {
///     let mut v = ref_seq[..51].to_vec();
///     v.extend_from_slice(&ref_seq[52..]);
///     v
/// };
/// let a = extract_minimal_allele("chr2:1000-1100", &ref_seq, &alt_seq).unwrap();
/// assert_eq!(a.pos, 1050);
/// assert_eq!(a.ref_allele, "TG");
/// assert_eq!(a.alt_allele, "T");
/// ```
pub fn extract_minimal_allele(
    target_id: &str,
    ref_seq: &[u8],
    alt_seq: &[u8],
) -> Option<MinimalAllele> {
    let (chrom, target_start_i64) = parse_target_id(target_id)?;
    // Cast i64 to u64; genomic positions are always non-negative.
    let target_start = target_start_i64 as u64;

    // Step 1: find the first differing position.
    let diff_pos = ref_seq
        .iter()
        .zip(alt_seq.iter())
        .position(|(r, a)| r != a)?;

    // Step 2: trim the outer common suffix (must not consume diff_pos).
    let mut outer_cs = 0usize;
    while outer_cs < ref_seq.len().saturating_sub(diff_pos + 1)
        && outer_cs < alt_seq.len().saturating_sub(diff_pos + 1)
        && ref_seq[ref_seq.len() - 1 - outer_cs] == alt_seq[alt_seq.len() - 1 - outer_cs]
    {
        outer_cs += 1;
    }
    let ref_trimmed = &ref_seq[diff_pos..ref_seq.len() - outer_cs];
    let alt_trimmed = &alt_seq[diff_pos..alt_seq.len() - outer_cs];

    if ref_seq.len() == alt_seq.len() {
        // SNV / MNV: no indel normalisation needed.
        Some(MinimalAllele {
            chrom,
            pos: target_start + diff_pos as u64,
            ref_allele: bytes_to_string(ref_trimmed),
            alt_allele: bytes_to_string(alt_trimmed),
        })
    } else {
        let is_deletion = ref_seq.len() > alt_seq.len();
        indel_allele(
            chrom,
            target_start,
            ref_seq,
            ref_trimmed,
            alt_trimmed,
            diff_pos,
            is_deletion,
        )
    }
}

// ─── Private helpers ──────────────────────────────────────────────────────────

/// Compute the left-normalised indel allele from trimmed region slices.
fn indel_allele(
    chrom: String,
    target_start: u64,
    ref_seq: &[u8],
    ref_trimmed: &[u8],
    alt_trimmed: &[u8],
    diff_pos: usize,
    is_deletion: bool,
) -> Option<MinimalAllele> {
    // Step 3: strip residual inner common suffix (allows alt_min to be empty).
    let mut inner_cs = 0usize;
    while inner_cs < ref_trimmed.len()
        && inner_cs < alt_trimmed.len()
        && ref_trimmed[ref_trimmed.len() - 1 - inner_cs]
            == alt_trimmed[alt_trimmed.len() - 1 - inner_cs]
    {
        inner_cs += 1;
    }
    let ref_min = &ref_trimmed[..ref_trimmed.len() - inner_cs];
    let alt_min = &alt_trimmed[..alt_trimmed.len() - inner_cs];

    // Strip inner common prefix.
    let inner_cp = ref_min
        .iter()
        .zip(alt_min.iter())
        .take_while(|(r, a)| r == a)
        .count();
    let event_start = diff_pos + inner_cp;

    if is_deletion {
        // Deleted bases: portion of ref_min beyond the shared prefix.
        let mut del_seq: Vec<u8> = ref_min[inner_cp..].to_vec();
        normalise_and_anchor(
            chrom,
            target_start,
            ref_seq,
            &mut del_seq,
            event_start,
            true,
        )
    } else {
        // Inserted bases: portion of alt_min beyond the shared prefix.
        let mut ins_seq: Vec<u8> = alt_min[inner_cp..].to_vec();
        normalise_and_anchor(
            chrom,
            target_start,
            ref_seq,
            &mut ins_seq,
            event_start,
            false,
        )
    }
}

/// Shift `event_seq` left and build the anchor-prefixed alleles.
///
/// For a deletion: REF = anchor + del_seq, ALT = anchor.
/// For an insertion: REF = anchor, ALT = anchor + ins_seq.
///
/// Left-normalisation: while `anchor_pos > 0` and the preceding reference base
/// equals the last base of `event_seq`, rotate `event_seq` right by one and
/// decrement `anchor_pos`.
fn normalise_and_anchor(
    chrom: String,
    target_start: u64,
    ref_seq: &[u8],
    event_seq: &mut Vec<u8>,
    event_start: usize,
    is_deletion: bool,
) -> Option<MinimalAllele> {
    // The anchor position is the base immediately before the event.
    let mut anchor_pos = event_start.saturating_sub(1);

    // Step 4: left-normalise.
    while anchor_pos > 0 && !event_seq.is_empty() && ref_seq[anchor_pos] == *event_seq.last()? {
        let last = event_seq.remove(event_seq.len() - 1);
        event_seq.insert(0, last);
        anchor_pos -= 1;
    }

    // Step 5: prepend the anchor base.
    let anchor = ref_seq[anchor_pos] as char;
    let (ref_allele, alt_allele) = if is_deletion {
        (
            format!("{anchor}{}", bytes_to_string(event_seq)),
            anchor.to_string(),
        )
    } else {
        (
            anchor.to_string(),
            format!("{anchor}{}", bytes_to_string(event_seq)),
        )
    };

    Some(MinimalAllele {
        chrom,
        pos: target_start + anchor_pos as u64,
        ref_allele,
        alt_allele,
    })
}

fn bytes_to_string(seq: &[u8]) -> String {
    String::from_utf8_lossy(seq).into_owned()
}

// ─── Tests ────────────────────────────────────────────────────────────────────

#[cfg(test)]
mod tests {
    use super::*;

    fn check(target_id: &str, ref_seq: &[u8], alt_seq: &[u8], pos: u64, r: &str, a: &str) {
        let result = extract_minimal_allele(target_id, ref_seq, alt_seq)
            .unwrap_or_else(|| panic!("expected Some for {target_id}"));
        assert_eq!(
            result.pos, pos,
            "pos mismatch: got {}, want {}",
            result.pos, pos
        );
        assert_eq!(
            result.ref_allele, r,
            "REF mismatch: got {}, want {r}",
            result.ref_allele
        );
        assert_eq!(
            result.alt_allele, a,
            "ALT mismatch: got {}, want {a}",
            result.alt_allele
        );
    }

    // Test 1: SNV at index 50 of a 101-bp target.
    #[test]
    fn snv_correct_position() {
        let mut r = vec![b'A'; 101];
        let mut a = r.clone();
        r[50] = b'C';
        a[50] = b'T';
        check("chr1:1000-1100", &r, &a, 1050, "C", "T");
    }

    // Test 2: Identical sequences → None.
    #[test]
    fn identical_returns_none() {
        let seq = vec![b'A'; 101];
        assert!(extract_minimal_allele("chr1:1000-1100", &seq, &seq).is_none());
    }

    // Test 3: Invalid target_id → None.
    #[test]
    fn invalid_id_returns_none() {
        let r = vec![b'A'; 5];
        let mut a = r.clone();
        a[2] = b'T';
        assert!(extract_minimal_allele("notvalid", &r, &a).is_none());
    }

    // Test 4: Deletion in a GG run is left-normalised to the T anchor.
    //
    // ref: ...T(50) G(51) G(52) T(53)...
    // alt: ...T(50) G(51) T(52)...       (G at index 51 deleted)
    // First diff at index 52 (because ref[51]==alt[51] after the shift);
    // left-normalise → anchor at index 50 (T).
    // Expected: pos=target_start+50, REF="TG", ALT="T".
    #[test]
    fn deletion_left_normalised_in_gg_run() {
        let mut r = vec![b'A'; 101];
        r[50] = b'T';
        r[51] = b'G';
        r[52] = b'G';
        r[53] = b'T';
        // Delete r[51] (first G).
        let mut a: Vec<u8> = r[..51].to_vec();
        a.extend_from_slice(&r[52..]);
        check("chr2:1000-1100", &r, &a, 1050, "TG", "T");
    }

    // Test 5: Simple 1-bp insertion with no repeat — no left shift.
    //
    // ref: ...A(50) A(51)... (all A)
    // alt: ...A(50) G A(51)... (G inserted)
    // anchor = A at index 50; A ≠ G, so no shift.
    // Expected: pos=target_start+50, REF="A", ALT="AG".
    #[test]
    fn insertion_no_repeat_anchor_correct() {
        let r: Vec<u8> = vec![b'A'; 100];
        let mut a: Vec<u8> = r[..51].to_vec();
        a.push(b'G');
        a.extend_from_slice(&r[51..]);
        check("chr3:500-599", &r, &a, 550, "A", "AG");
    }

    // Test 6: MNV (two adjacent SNVs) — same length, no normalisation.
    #[test]
    fn mnv_same_length() {
        let mut r = vec![b'A'; 10];
        let mut a = r.clone();
        r[3] = b'C';
        r[4] = b'G';
        a[3] = b'T';
        a[4] = b'T';
        check("chrX:100-109", &r, &a, 103, "CG", "TT");
    }

    // Test 7: Multi-base deletion with distinct flanking bases (no repeat).
    // Delete "GCT" at indices 40-42.
    #[test]
    fn multi_base_deletion_no_repeat() {
        let mut r: Vec<u8> = vec![b'A'; 103];
        r[39] = b'C'; // anchor
        r[40] = b'G';
        r[41] = b'C';
        r[42] = b'T';
        r[43] = b'T'; // base after deletion (distinct to prevent shift)
        let mut a: Vec<u8> = r[..40].to_vec();
        a.extend_from_slice(&r[43..]);
        // anchor=r[39]='C', del_seq="GCT", del_seq[-1]='T'; r[39]='C' ≠ 'T' → no shift
        check("chr5:500-602", &r, &a, 539, "CGCT", "C");
    }

    // Test 8: parse_target_id (from targeting) handles standard coordinates.
    #[test]
    fn parse_target_id_standard() {
        use crate::targeting::parse_target_id;
        assert_eq!(
            parse_target_id("chr2:29318295-29318395"),
            Some(("chr2".to_string(), 29_318_295_i64))
        );
    }

    // Test 9: parse_target_id (from targeting) rejects malformed inputs.
    #[test]
    fn parse_target_id_malformed() {
        use crate::targeting::parse_target_id;
        assert!(parse_target_id("no_colon").is_none());
        assert!(parse_target_id("chr1:abc-100").is_none());
        assert!(parse_target_id("chr1:100").is_none()); // no dash
    }

    // Test 10: Insertion with repeat context is left-normalised.
    //
    // Insert T into a "TT" run.
    // ref: ...A(49) T(50) T(51) A(52)...
    // alt: ...A(49) T(50) T(51) T(52) A(53)... (T inserted after index 51)
    // First diff at index 52 (ref[52]=A, alt[52]=T after insert shifts everything).
    // ins_seq = "T"; ref[52]='A'... actually wait, the first diff is where the
    // sequence shifts, which is at index 52. ins_seq[-1]='T', ref[52-1]=ref[51]='T'
    // → shift: anchor_pos=50. ref[50]='T', ins_seq[-1]='T' → shift: anchor_pos=49.
    // ref[49]='A' ≠ 'T' → stop.
    // Expected: pos=target_start+49, REF="A", ALT="AT".
    #[test]
    fn insertion_left_normalised_in_tt_run() {
        let mut r: Vec<u8> = vec![b'A'; 100];
        r[50] = b'T';
        r[51] = b'T';
        // Insert T after index 51.
        let mut a: Vec<u8> = r[..52].to_vec();
        a.push(b'T');
        a.extend_from_slice(&r[52..]);
        check("chr7:1000-1099", &r, &a, 1049, "A", "AT");
    }

    // ── Additional edge-case tests ───────────────────────────────────────────

    // Test 11: SNV at the very first position (index 0) of the target.
    // Verifies that the 1-based output position equals the target start exactly.
    #[test]
    fn snv_at_first_position() {
        let mut r = vec![b'A'; 10];
        let mut a = r.clone();
        r[0] = b'C';
        a[0] = b'T';
        check("chr1:500-509", &r, &a, 500, "C", "T");
    }

    // Test 12: SNV at the very last position of the target.
    // Boundary: common suffix is zero, common prefix is len-1.
    #[test]
    fn snv_at_last_position() {
        let mut r = vec![b'A'; 10];
        let mut a = r.clone();
        r[9] = b'C';
        a[9] = b'T';
        check("chr1:100-109", &r, &a, 109, "C", "T");
    }

    // Test 13: Deletion in a homopolymer context (AAAA).
    // Deleting one A from a run of 4 A's must left-normalise to the leftmost
    // anchor before the run.
    #[test]
    fn deletion_in_homopolymer_left_normalises() {
        // ref: T(0) A(1) A(2) A(3) A(4) G(5)
        // alt: T(0) A(1) A(2) A(3) G(4) — one A deleted from the run.
        // The anchor must shift left to the T at index 0.
        let r = b"TAAAAAG".to_vec();
        let mut a = r[..4].to_vec(); // T A A A
        a.extend_from_slice(&r[5..]); // G
        // After left-normalisation: anchor is at index 0 (T), del_seq = A.
        // REF = "TA", ALT = "T".
        check("chr2:200-206", &r, &a, 200, "TA", "T");
    }

    // Test 14: Insertion in a homopolymer context (AAAA).
    // Inserting an extra A into a poly-A run must left-normalise.
    #[test]
    fn insertion_in_homopolymer_left_normalises() {
        // ref: G(0) A(1) A(2) A(3) T(4)
        // alt: G(0) A(1) A(2) A(3) A(4) T(5) — one A inserted.
        let r = b"GAAAT".to_vec();
        let mut a = r[..4].to_vec(); // G A A A
        a.push(b'A');
        a.push(b'T');
        // Left-normalise: anchor shifts to index 0 (G).
        // REF = "G", ALT = "GA".
        check("chr3:50-54", &r, &a, 50, "G", "GA");
    }

    // Test 15: Multi-base deletion where the deleted region differs from the
    // flanking context, so the first differing position is detectable.
    // ref: ACGT + CCCC + ACGT = 12bp, alt: ACGT + ACGT = 8bp (CCCC deleted).
    // The diff is at index 4 where ref[4]='C' and alt[4]='A'.
    #[test]
    fn multi_base_deletion_spanning_entire_diff() {
        let r = b"ACGTCCCCACGT".to_vec();
        let a = b"ACGTACGT".to_vec();
        let result =
            extract_minimal_allele("chr1:100-111", &r, &a).expect("should produce an allele");
        // The deletion must have ref_allele longer than alt_allele.
        assert!(
            result.ref_allele.len() > result.alt_allele.len(),
            "deletion must have REF longer than ALT"
        );
        assert_eq!(result.chrom, "chr1", "chrom must be chr1");
    }

    // Test 16: MNV where both bases differ (e.g. CG → TA).
    // Equal-length sequences with >1 differing base, not an inversion.
    #[test]
    fn mnv_two_adjacent_changes() {
        let mut r = vec![b'G'; 20];
        let mut a = r.clone();
        r[10] = b'A';
        r[11] = b'C';
        a[10] = b'T';
        a[11] = b'G';
        check("chr4:300-319", &r, &a, 310, "AC", "TG");
    }

    // Test 17: Very short target (3 bases, shorter than typical k) still works.
    // Ensures the algorithm handles short sequences without panicking.
    #[test]
    fn very_short_target_still_extracts() {
        let r = b"ACG".to_vec();
        let a = b"ATG".to_vec();
        check("chrM:1-3", &r, &a, 2, "C", "T");
    }

    // Test 18: Large deletion at the right boundary of the sequence.
    // The deletion removes several bases near the end, but the sequences
    // differ in the overlapping region so the variant is detectable.
    // ref: ACGT + TTTT + G = 9bp. alt: ACGT + G = 5bp.
    // Diff at index 4: ref[4]='T', alt[4]='G'. Suffix 'G' is shared.
    #[test]
    fn deletion_at_right_boundary() {
        let r = b"ACGTTTTTG".to_vec();
        let a = b"ACGTG".to_vec();
        let result =
            extract_minimal_allele("chr5:10-18", &r, &a).expect("should produce an allele");
        assert!(
            result.ref_allele.len() > result.alt_allele.len(),
            "deletion: REF must be longer than ALT"
        );
    }

    // Test 19: Target ID with large coordinate values (chr1:999999-1000099).
    // Verifies that position arithmetic handles large numbers correctly.
    #[test]
    fn large_coordinate_target_id() {
        let mut r = vec![b'A'; 101];
        let mut a = r.clone();
        r[50] = b'C';
        a[50] = b'G';
        check("chr1:999999-1000099", &r, &a, 1000049, "C", "G");
    }

    // Test 20: MinimalAllele Display trait formats correctly.
    // The Display output is used in logging, so it must be readable.
    #[test]
    fn minimal_allele_display() {
        let allele = MinimalAllele {
            chrom: "chr1".to_string(),
            pos: 12345,
            ref_allele: "A".to_string(),
            alt_allele: "T".to_string(),
        };
        let s = format!("{allele}");
        assert_eq!(s, "chr1:12345 A->T", "Display format mismatch");
    }
}
