refeff-io 0.2.0

//! FEFF `compton.dat` profile text codec.
//!
//! The COMPTON module writes a comment header describing the integration grid
//! followed by two numeric columns: projected momentum `pq` and the Compton
//! profile `J(pq)`.

use std::fmt::Write as _;
use std::path::Path;

use ndarray::{Array1, Array2};

use crate::error::{IoError, Result};
use crate::format::fortran_list_directed_g15_f64;

const COMPTON_DAT_ROW_WIDTH: usize = 2;
const RHOZZP_DAT_ROW_WIDTH: usize = 2;

#[derive(Default)]
struct ComptonDatHeader {
    ns: Option<usize>,
    nphi: Option<usize>,
    nz: Option<usize>,
    nzp: Option<usize>,
    zpmax: Option<f64>,
    temperature_ev: Option<f64>,
}

/// Parsed FEFF `compton.dat` contents.
#[derive(Debug, Clone, PartialEq)]
pub struct ComptonDatData {
    /// Header and comment lines before the numeric Compton profile table.
    pub header_lines: Vec<String>,
    /// Number of radial integration points in cylindrical radius.
    pub ns: Option<usize>,
    /// Number of azimuthal integration points.
    pub nphi: Option<usize>,
    /// Number of `z` integration points.
    pub nz: Option<usize>,
    /// Number of `z'` integration points.
    pub nzp: Option<usize>,
    /// Maximum `z'` integration coordinate from the header.
    pub zpmax: Option<f64>,
    /// Electronic temperature in eV from the header.
    pub temperature_ev: Option<f64>,
    /// Projected momentum grid `pq`.
    pub momentum: Array1<f64>,
    /// Compton profile values `J(pq)`.
    pub profile: Array1<f64>,
}

impl ComptonDatData {
    /// Number of profile rows.
    #[must_use]
    pub fn point_count(&self) -> usize {
        self.momentum.len()
    }
}

/// Lossless view of a parsed FEFF `compton.dat` profile.
///
/// FEFF emits the numeric rows with compiler-dependent list-directed
/// formatting. This wrapper retains validated source text for exact fixture
/// roundtrips while the semantic codec remains canonical.
#[derive(Debug, Clone, PartialEq)]
pub struct ComptonDatLosslessData {
    pub data: ComptonDatData,
    pub original_text: String,
}

/// Parsed FEFF `rhozzp.dat` diagnostic density slice.
#[derive(Debug, Clone, PartialEq)]
pub struct RhozzpDatData {
    /// Optional comment lines before or around the numeric diagnostic table.
    pub header_lines: Vec<String>,
    /// `z'` coordinate grid.
    pub z_prime: Array1<f64>,
    /// Density matrix slice values `rho(z,z')`.
    pub density: Array1<f64>,
}

impl RhozzpDatData {
    /// Number of diagnostic rows.
    #[must_use]
    pub fn point_count(&self) -> usize {
        self.z_prime.len()
    }
}

/// Lossless view of a parsed FEFF `rhozzp.dat` diagnostic.
///
/// FEFF uses compiler-dependent list-directed formatting for these rows, so
/// exact fixture comparisons retain the validated source representation.
#[derive(Debug, Clone, PartialEq)]
pub struct RhozzpDatLosslessData {
    pub data: RhozzpDatData,
    pub original_text: String,
}

/// Parsed FEFF `jzzp.dat` reusable Compton cache.
#[derive(Debug, Clone, PartialEq)]
pub struct JzzpDatData {
    /// Number of radial integration points in cylindrical radius.
    pub ns: usize,
    /// Number of azimuthal integration points.
    pub nphi: usize,
    /// Number of `z` integration points.
    pub nz: usize,
    /// Number of `z'` integration points.
    pub nzp: usize,
    /// Maximum cylindrical radius used to build the cache.
    pub smax: f64,
    /// Maximum azimuthal angle used to build the cache.
    pub phimax: f64,
    /// Maximum `z` integration coordinate used to build the cache.
    pub zmax: f64,
    /// Maximum `z'` integration coordinate used to build the cache.
    pub zpmax: f64,
    /// Cached `J(z,z')` values with shape `(nz, nzp)`.
    pub values: Array2<f64>,
}

/// Lossless view of a parsed FEFF `jzzp.dat` cache.
///
/// FEFF writes this cache with compiler-dependent list-directed spacing and
/// line wrapping. The semantic codec remains canonical; this wrapper retains
/// the validated source text when byte-exact fixture roundtrips are required.
#[derive(Debug, Clone, PartialEq)]
pub struct JzzpDatLosslessData {
    pub data: JzzpDatData,
    pub original_text: String,
}

impl JzzpDatData {
    /// Number of cached `J(z,z')` values.
    #[must_use]
    pub fn value_count(&self) -> usize {
        self.values.len()
    }
}

/// Render FEFF-compatible `compton.dat` text.
pub fn compton_dat_string(data: &ComptonDatData) -> Result<String> {
    validate_compton_dat(data)?;

    let mut out = String::new();
    for line in &data.header_lines {
        writeln!(out, "{line}")?;
    }
    for (momentum, profile) in data.momentum.iter().zip(data.profile.iter()) {
        writeln!(
            out,
            "{}{}",
            fortran_list_directed_g15_f64(*momentum),
            fortran_list_directed_g15_f64(*profile)
        )?;
    }
    Ok(out)
}

/// Parse FEFF `compton.dat` text.
pub fn parse_compton_dat(text: &str) -> Result<ComptonDatData> {
    let mut header_lines = Vec::new();
    let mut header = ComptonDatHeader::default();
    let mut momentum = Vec::new();
    let mut profile = Vec::new();

    for (index, raw) in text.lines().enumerate() {
        let line_number = index + 1;
        let line = raw.trim_end();
        let tokens = line.split_whitespace().collect::<Vec<_>>();
        if tokens.first().is_some_and(|token| is_numeric_token(token)) {
            if tokens.len() != COMPTON_DAT_ROW_WIDTH {
                return Err(IoError::ComptonDatRowWidth {
                    line: line_number,
                    actual: tokens.len(),
                    expected: COMPTON_DAT_ROW_WIDTH,
                });
            }
            momentum.push(parse_f64(line_number, "momentum", tokens[0])?);
            profile.push(parse_f64(line_number, "profile", tokens[1])?);
        } else {
            parse_header_metadata(line, line_number, &mut header)?;
            header_lines.push(raw.to_string());
        }
    }

    let data = ComptonDatData {
        header_lines,
        ns: header.ns,
        nphi: header.nphi,
        nz: header.nz,
        nzp: header.nzp,
        zpmax: header.zpmax,
        temperature_ev: header.temperature_ev,
        momentum: Array1::from_vec(momentum),
        profile: Array1::from_vec(profile),
    };
    validate_compton_dat(&data)?;
    Ok(data)
}

/// Parse and retain exact validated `compton.dat` source text.
pub fn parse_compton_dat_lossless(text: &str) -> Result<ComptonDatLosslessData> {
    Ok(ComptonDatLosslessData {
        data: parse_compton_dat(text)?,
        original_text: text.to_string(),
    })
}

/// Render a lossless `compton.dat` view.
///
/// Unchanged parsed data reuses the exact validated source. If callers replace
/// the semantic payload, the canonical renderer is used instead.
pub fn compton_dat_lossless_string(data: &ComptonDatLosslessData) -> Result<String> {
    if parse_compton_dat(&data.original_text)? == data.data {
        Ok(data.original_text.clone())
    } else {
        compton_dat_string(&data.data)
    }
}

/// Write FEFF `compton.dat` text to a file.
pub fn write_compton_dat(path: impl AsRef<Path>, data: &ComptonDatData) -> Result<()> {
    let path = path.as_ref();
    std::fs::write(path, compton_dat_string(data)?).map_err(|source| IoError::io(path, source))
}

/// Read FEFF `compton.dat` text from a file.
pub fn read_compton_dat(path: impl AsRef<Path>) -> Result<ComptonDatData> {
    let path = path.as_ref();
    let text = std::fs::read_to_string(path).map_err(|source| IoError::io(path, source))?;
    parse_compton_dat(&text)
}

/// Render FEFF-compatible `rhozzp.dat` text.
pub fn rhozzp_dat_string(data: &RhozzpDatData) -> Result<String> {
    validate_rhozzp_dat(data)?;

    let mut out = String::new();
    for line in &data.header_lines {
        writeln!(out, "{line}")?;
    }
    for (z_prime, density) in data.z_prime.iter().zip(data.density.iter()) {
        writeln!(
            out,
            "{}{}",
            fortran_list_directed_g15_f64(*z_prime),
            fortran_list_directed_g15_f64(*density)
        )?;
    }
    Ok(out)
}

/// Parse FEFF `rhozzp.dat` diagnostic text.
pub fn parse_rhozzp_dat(text: &str) -> Result<RhozzpDatData> {
    let mut header_lines = Vec::new();
    let mut z_prime = Vec::new();
    let mut density = Vec::new();

    for (index, raw) in text.lines().enumerate() {
        let line_number = index + 1;
        let line = raw.trim_end();
        let tokens = line.split_whitespace().collect::<Vec<_>>();
        if tokens.first().is_some_and(|token| is_numeric_token(token)) {
            if tokens.len() != RHOZZP_DAT_ROW_WIDTH {
                return Err(IoError::RhozzpDatRowWidth {
                    line: line_number,
                    actual: tokens.len(),
                    expected: RHOZZP_DAT_ROW_WIDTH,
                });
            }
            z_prime.push(parse_rhozzp_f64(line_number, "z prime", tokens[0])?);
            density.push(parse_rhozzp_f64(line_number, "density", tokens[1])?);
        } else {
            header_lines.push(raw.to_string());
        }
    }

    let data = RhozzpDatData {
        header_lines,
        z_prime: Array1::from_vec(z_prime),
        density: Array1::from_vec(density),
    };
    validate_rhozzp_dat(&data)?;
    Ok(data)
}

/// Parse and retain exact validated `rhozzp.dat` source text.
pub fn parse_rhozzp_dat_lossless(text: &str) -> Result<RhozzpDatLosslessData> {
    Ok(RhozzpDatLosslessData {
        data: parse_rhozzp_dat(text)?,
        original_text: text.to_string(),
    })
}

/// Render a lossless `rhozzp.dat` view.
pub fn rhozzp_dat_lossless_string(data: &RhozzpDatLosslessData) -> Result<String> {
    if parse_rhozzp_dat(&data.original_text)? == data.data {
        Ok(data.original_text.clone())
    } else {
        rhozzp_dat_string(&data.data)
    }
}

/// Write FEFF `rhozzp.dat` text to a file.
pub fn write_rhozzp_dat(path: impl AsRef<Path>, data: &RhozzpDatData) -> Result<()> {
    let path = path.as_ref();
    std::fs::write(path, rhozzp_dat_string(data)?).map_err(|source| IoError::io(path, source))
}

/// Read FEFF `rhozzp.dat` text from a file.
pub fn read_rhozzp_dat(path: impl AsRef<Path>) -> Result<RhozzpDatData> {
    let path = path.as_ref();
    let text = std::fs::read_to_string(path).map_err(|source| IoError::io(path, source))?;
    parse_rhozzp_dat(&text)
}

/// Render FEFF-compatible `jzzp.dat` cache text.
pub fn jzzp_dat_string(data: &JzzpDatData) -> Result<String> {
    validate_jzzp_dat(data)?;

    let mut out = String::new();
    writeln!(out, "# {} {} {} {}", data.ns, data.nphi, data.nz, data.nzp)?;
    writeln!(
        out,
        "# {smax:24.17E} {phimax:24.17E} {zmax:24.17E} {zpmax:24.17E}",
        smax = data.smax,
        phimax = data.phimax,
        zmax = data.zmax,
        zpmax = data.zpmax
    )?;

    let mut emitted = 0_usize;
    for col in 0..data.nzp {
        for row in 0..data.nz {
            write!(out, " {value:24.17E}", value = data.values[[row, col]])?;
            emitted += 1;
            if emitted.is_multiple_of(3) {
                writeln!(out)?;
            }
        }
    }
    if !emitted.is_multiple_of(3) {
        writeln!(out)?;
    }
    Ok(out)
}

/// Parse FEFF `jzzp.dat` cache text.
pub fn parse_jzzp_dat(text: &str) -> Result<JzzpDatData> {
    let mut lines = text
        .lines()
        .enumerate()
        .filter(|(_, line)| !line.trim().is_empty());

    let (grid_line_number, grid_line) = lines
        .next()
        .ok_or(IoError::JzzpDatMissing { field: "grid" })?;
    let grid_line_number = grid_line_number + 1;
    let grid_tokens = header_tokens(grid_line_number, "grid", grid_line, 5)?;
    let ns = parse_jzzp_usize(grid_line_number, "ns", grid_tokens[1])?;
    let nphi = parse_jzzp_usize(grid_line_number, "nphi", grid_tokens[2])?;
    let nz = parse_jzzp_usize(grid_line_number, "nz", grid_tokens[3])?;
    let nzp = parse_jzzp_usize(grid_line_number, "nzp", grid_tokens[4])?;

    let (limits_line_number, limits_line) = lines
        .next()
        .ok_or(IoError::JzzpDatMissing { field: "limits" })?;
    let limits_line_number = limits_line_number + 1;
    let limit_tokens = header_tokens(limits_line_number, "limits", limits_line, 5)?;
    let smax = parse_jzzp_f64(limits_line_number, "smax", limit_tokens[1])?;
    let phimax = parse_jzzp_f64(limits_line_number, "phimax", limit_tokens[2])?;
    let zmax = parse_jzzp_f64(limits_line_number, "zmax", limit_tokens[3])?;
    let zpmax = parse_jzzp_f64(limits_line_number, "zpmax", limit_tokens[4])?;

    let mut payload = Vec::new();
    for (index, raw) in lines {
        let line_number = index + 1;
        for token in raw.split_whitespace() {
            payload.push(parse_jzzp_f64(line_number, "values", token)?);
        }
    }

    let expected = jzzp_expected_len(nz, nzp)?;
    if payload.len() != expected {
        return Err(invalid_jzzp_dat(
            "values",
            format!(
                "expected {expected} value(s) from nz*nzp but found {}",
                payload.len()
            ),
        ));
    }

    let values = Array2::from_shape_fn((nz, nzp), |(row, col)| payload[col * nz + row]);
    let data = JzzpDatData {
        ns,
        nphi,
        nz,
        nzp,
        smax,
        phimax,
        zmax,
        zpmax,
        values,
    };
    validate_jzzp_dat(&data)?;
    Ok(data)
}

/// Parse and retain exact validated `jzzp.dat` source text.
pub fn parse_jzzp_dat_lossless(text: &str) -> Result<JzzpDatLosslessData> {
    Ok(JzzpDatLosslessData {
        data: parse_jzzp_dat(text)?,
        original_text: text.to_string(),
    })
}

/// Render a lossless `jzzp.dat` view.
///
/// Unchanged parsed data reuses the exact validated source. If callers replace
/// the semantic payload, the canonical renderer is used instead.
pub fn jzzp_dat_lossless_string(data: &JzzpDatLosslessData) -> Result<String> {
    if parse_jzzp_dat(&data.original_text)? == data.data {
        Ok(data.original_text.clone())
    } else {
        jzzp_dat_string(&data.data)
    }
}

/// Write FEFF `jzzp.dat` cache text to a file.
pub fn write_jzzp_dat(path: impl AsRef<Path>, data: &JzzpDatData) -> Result<()> {
    let path = path.as_ref();
    std::fs::write(path, jzzp_dat_string(data)?).map_err(|source| IoError::io(path, source))
}

/// Read FEFF `jzzp.dat` cache text from a file.
pub fn read_jzzp_dat(path: impl AsRef<Path>) -> Result<JzzpDatData> {
    let path = path.as_ref();
    let text = std::fs::read_to_string(path).map_err(|source| IoError::io(path, source))?;
    parse_jzzp_dat(&text)
}

fn parse_header_metadata(
    line: &str,
    line_number: usize,
    header: &mut ComptonDatHeader,
) -> Result<()> {
    let lower = line.to_ascii_lowercase();
    if lower.contains("nphi:") {
        header.nphi = Some(parse_usize_header(line_number, "nphi", line)?);
    } else if lower.contains("nzp:") {
        header.nzp = Some(parse_usize_header(line_number, "nzp", line)?);
    } else if lower.contains("ns:") {
        header.ns = Some(parse_usize_header(line_number, "ns", line)?);
    } else if lower.contains("nz:") {
        header.nz = Some(parse_usize_header(line_number, "nz", line)?);
    } else if lower.contains("zpmax:") {
        header.zpmax = Some(parse_f64_header(line_number, "zpmax", line)?);
    } else if lower.contains("temperature") {
        header.temperature_ev = Some(parse_f64_header(line_number, "temperature", line)?);
    }
    Ok(())
}

fn validate_compton_dat(data: &ComptonDatData) -> Result<()> {
    let point_count = data.point_count();
    if point_count == 0 {
        return Err(invalid_compton_dat(
            "rows",
            "at least one Compton profile row is required",
        ));
    }
    validate_len("profile", data.profile.len(), point_count)?;

    validate_positive_header("ns", data.ns)?;
    validate_positive_header("nphi", data.nphi)?;
    validate_positive_header("nz", data.nz)?;
    validate_positive_header("nzp", data.nzp)?;

    if let Some(value) = data.zpmax {
        validate_finite("zpmax", value)?;
    }
    if let Some(value) = data.temperature_ev {
        validate_finite("temperature", value)?;
    }

    for (row, (momentum, profile)) in data.momentum.iter().zip(data.profile.iter()).enumerate() {
        let row = row + 1;
        validate_finite_row("momentum", *momentum, row)?;
        validate_finite_row("profile", *profile, row)?;
    }

    Ok(())
}

fn validate_len(field: &'static str, actual: usize, expected: usize) -> Result<()> {
    if actual == expected {
        Ok(())
    } else {
        Err(IoError::ComptonDatShape {
            field,
            actual,
            expected,
        })
    }
}

fn validate_rhozzp_dat(data: &RhozzpDatData) -> Result<()> {
    let point_count = data.point_count();
    if point_count == 0 {
        return Err(invalid_rhozzp_dat(
            "rows",
            "at least one rhozzp diagnostic row is required",
        ));
    }
    validate_rhozzp_len("density", data.density.len(), point_count)?;

    for (row, (z_prime, density)) in data.z_prime.iter().zip(data.density.iter()).enumerate() {
        let row = row + 1;
        validate_rhozzp_finite_row("z prime", *z_prime, row)?;
        validate_rhozzp_finite_row("density", *density, row)?;
    }

    Ok(())
}

fn validate_jzzp_dat(data: &JzzpDatData) -> Result<()> {
    validate_jzzp_positive("ns", data.ns)?;
    validate_jzzp_positive("nphi", data.nphi)?;
    validate_jzzp_positive("nz", data.nz)?;
    validate_jzzp_positive("nzp", data.nzp)?;
    validate_jzzp_finite("smax", data.smax)?;
    validate_jzzp_finite("phimax", data.phimax)?;
    validate_jzzp_finite("zmax", data.zmax)?;
    validate_jzzp_finite("zpmax", data.zpmax)?;

    let (rows, cols) = data.values.dim();
    if rows != data.nz || cols != data.nzp {
        return Err(IoError::JzzpDatShape {
            field: "values",
            rows,
            cols,
            expected_rows: data.nz,
            expected_cols: data.nzp,
        });
    }

    for (index, value) in data.values.iter().enumerate() {
        let row = index % data.nz + 1;
        validate_jzzp_finite_row("values", *value, row)?;
    }
    Ok(())
}

fn validate_rhozzp_len(field: &'static str, actual: usize, expected: usize) -> Result<()> {
    if actual == expected {
        Ok(())
    } else {
        Err(IoError::RhozzpDatShape {
            field,
            actual,
            expected,
        })
    }
}

fn validate_jzzp_positive(field: &'static str, value: usize) -> Result<()> {
    if value == 0 {
        Err(invalid_jzzp_dat(field, "value must be positive"))
    } else {
        Ok(())
    }
}

fn validate_positive_header(field: &'static str, value: Option<usize>) -> Result<()> {
    if value.is_some_and(|value| value == 0) {
        Err(invalid_compton_dat(field, "value must be positive"))
    } else {
        Ok(())
    }
}

fn parse_usize_header(line: usize, field: &'static str, text: &str) -> Result<usize> {
    parse_usize(
        line,
        field,
        last_numeric_token(text)
            .ok_or_else(|| invalid_compton_dat(field, "missing numeric header value"))?,
    )
}

fn parse_f64_header(line: usize, field: &'static str, text: &str) -> Result<f64> {
    parse_f64(
        line,
        field,
        last_numeric_token(text)
            .ok_or_else(|| invalid_compton_dat(field, "missing numeric header value"))?,
    )
}

fn parse_f64(line: usize, field: &'static str, token: &str) -> Result<f64> {
    token
        .replace(['D', 'd'], "E")
        .parse::<f64>()
        .map_err(|_| IoError::ComptonDatParse {
            field,
            line,
            token: token.to_string(),
        })
}

fn parse_usize(line: usize, field: &'static str, token: &str) -> Result<usize> {
    token
        .parse::<usize>()
        .map_err(|_| IoError::ComptonDatParse {
            field,
            line,
            token: token.to_string(),
        })
}

fn parse_rhozzp_f64(line: usize, field: &'static str, token: &str) -> Result<f64> {
    token
        .replace(['D', 'd'], "E")
        .parse::<f64>()
        .map_err(|_| IoError::RhozzpDatParse {
            field,
            line,
            token: token.to_string(),
        })
}

fn parse_jzzp_f64(line: usize, field: &'static str, token: &str) -> Result<f64> {
    token
        .replace(['D', 'd'], "E")
        .parse::<f64>()
        .map_err(|_| IoError::JzzpDatParse {
            field,
            line,
            token: token.to_string(),
        })
}

fn parse_jzzp_usize(line: usize, field: &'static str, token: &str) -> Result<usize> {
    token.parse::<usize>().map_err(|_| IoError::JzzpDatParse {
        field,
        line,
        token: token.to_string(),
    })
}

fn validate_finite(field: &'static str, value: f64) -> Result<()> {
    if value.is_finite() {
        Ok(())
    } else {
        Err(invalid_compton_dat(field, "value must be finite"))
    }
}

fn validate_finite_row(field: &'static str, value: f64, row: usize) -> Result<()> {
    if value.is_finite() {
        Ok(())
    } else {
        Err(IoError::InvalidComptonDat {
            field,
            message: format!("row {row} value must be finite"),
        })
    }
}

fn validate_rhozzp_finite_row(field: &'static str, value: f64, row: usize) -> Result<()> {
    if value.is_finite() {
        Ok(())
    } else {
        Err(IoError::InvalidRhozzpDat {
            field,
            message: format!("row {row} value must be finite"),
        })
    }
}

fn validate_jzzp_finite(field: &'static str, value: f64) -> Result<()> {
    if value.is_finite() {
        Ok(())
    } else {
        Err(invalid_jzzp_dat(field, "value must be finite"))
    }
}

fn validate_jzzp_finite_row(field: &'static str, value: f64, row: usize) -> Result<()> {
    if value.is_finite() {
        Ok(())
    } else {
        Err(IoError::InvalidJzzpDat {
            field,
            message: format!("row {row} value must be finite"),
        })
    }
}

fn invalid_compton_dat(field: &'static str, message: impl Into<String>) -> IoError {
    IoError::InvalidComptonDat {
        field,
        message: message.into(),
    }
}

fn invalid_rhozzp_dat(field: &'static str, message: impl Into<String>) -> IoError {
    IoError::InvalidRhozzpDat {
        field,
        message: message.into(),
    }
}

fn invalid_jzzp_dat(field: &'static str, message: impl Into<String>) -> IoError {
    IoError::InvalidJzzpDat {
        field,
        message: message.into(),
    }
}

fn header_tokens<'a>(
    line: usize,
    field: &'static str,
    text: &'a str,
    expected: usize,
) -> Result<Vec<&'a str>> {
    let tokens = text.split_whitespace().collect::<Vec<_>>();
    if tokens.len() == expected && tokens.first() == Some(&"#") {
        Ok(tokens)
    } else {
        Err(invalid_jzzp_dat(
            field,
            format!("line {line} must start with # and contain {expected} token(s)"),
        ))
    }
}

fn jzzp_expected_len(nz: usize, nzp: usize) -> Result<usize> {
    nz.checked_mul(nzp)
        .ok_or_else(|| invalid_jzzp_dat("values", "nz*nzp overflows usize"))
}

fn is_numeric_token(token: &str) -> bool {
    token.replace(['D', 'd'], "E").parse::<f64>().is_ok()
}

fn last_numeric_token(line: &str) -> Option<&str> {
    line.split_whitespace()
        .rev()
        .find(|token| is_numeric_token(token))
}

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

    #[test]
    fn parses_feff_compton_reference_shape_and_metadata() -> Result<()> {
        let data = parse_compton_dat(COMPTON_DAT)?;
        assert_eq!(data.point_count(), 3);
        assert_eq!(data.ns, Some(32));
        assert_eq!(data.nphi, Some(32));
        assert_eq!(data.nz, Some(32));
        assert_eq!(data.nzp, Some(120));
        assert_eq!(data.zpmax, Some(10.0));
        assert_eq!(data.temperature_ev, Some(0.0));
        assert_eq!(data.momentum[1], 5.005004815757275e-3);
        assert_eq!(data.profile[2], 2.74462341659279);
        Ok(())
    }

    #[test]
    fn roundtrips_compton_text() -> Result<()> {
        let data = parse_compton_dat(COMPTON_DAT)?;
        let rendered = compton_dat_string(&data)?;
        assert_eq!(rendered, COMPTON_DAT);
        assert_eq!(parse_compton_dat(&rendered)?, data);
        Ok(())
    }

    #[test]
    fn lossless_compton_preserves_list_directed_rows() -> Result<()> {
        let text = "# Compton profile\n   0.0000000000000000        2.7447719580900216     \n";
        let data = parse_compton_dat_lossless(text)?;
        assert_eq!(compton_dat_lossless_string(&data)?, text);
        Ok(())
    }

    #[test]
    fn rejects_bad_compton_inputs() {
        assert!(parse_compton_dat("# no data\n").is_err());
        assert!(parse_compton_dat("1 2 3\n").is_err());
        assert!(parse_compton_dat("1 NaN\n").is_err());
        assert!(parse_compton_dat("# ns: 0\n1 2\n").is_err());

        let bad_shape = ComptonDatData {
            header_lines: Vec::new(),
            ns: None,
            nphi: None,
            nz: None,
            nzp: None,
            zpmax: None,
            temperature_ev: None,
            momentum: Array1::from_vec(vec![1.0, 2.0]),
            profile: Array1::from_vec(vec![1.0]),
        };
        assert!(compton_dat_string(&bad_shape).is_err());
    }

    #[test]
    fn parses_feff_rhozzp_reference_shape() -> Result<()> {
        let data = parse_rhozzp_dat(RHOZZP_DAT)?;
        assert_eq!(data.point_count(), 3);
        assert_eq!(data.z_prime[0], 9.999999776482582e-3);
        assert_eq!(data.density[1], 2.66921255682004);
        assert_eq!(data.density[2], 1.84694165446344);
        Ok(())
    }

    #[test]
    fn roundtrips_rhozzp_text() -> Result<()> {
        let data = parse_rhozzp_dat(RHOZZP_DAT)?;
        let rendered = rhozzp_dat_string(&data)?;
        assert_eq!(rendered, RHOZZP_DAT);
        assert_eq!(parse_rhozzp_dat(&rendered)?, data);
        Ok(())
    }

    #[test]
    fn lossless_rhozzp_preserves_list_directed_rows() -> Result<()> {
        let text = "   9.9999997764825821E-003   3.7109598151487062     \n";
        let data = parse_rhozzp_dat_lossless(text)?;
        assert_eq!(rhozzp_dat_lossless_string(&data)?, text);
        Ok(())
    }

    #[test]
    fn rejects_bad_rhozzp_inputs() {
        assert!(parse_rhozzp_dat("# no data\n").is_err());
        assert!(parse_rhozzp_dat("1 2 3\n").is_err());
        assert!(parse_rhozzp_dat("1 NaN\n").is_err());

        let bad_shape = RhozzpDatData {
            header_lines: Vec::new(),
            z_prime: Array1::from_vec(vec![1.0, 2.0]),
            density: Array1::from_vec(vec![1.0]),
        };
        assert!(rhozzp_dat_string(&bad_shape).is_err());
    }

    #[test]
    fn parses_jzzp_fortran_order_cache() -> Result<()> {
        let data = parse_jzzp_dat(JZZP_DAT)?;
        assert_eq!(data.ns, 2);
        assert_eq!(data.nphi, 3);
        assert_eq!(data.nz, 2);
        assert_eq!(data.nzp, 3);
        assert_eq!(data.value_count(), 6);
        assert_eq!(data.smax, 1.0);
        assert_eq!(data.phimax, 3.125);
        assert_eq!(data.values[[0, 0]], 1.0);
        assert_eq!(data.values[[1, 0]], 2.0);
        assert_eq!(data.values[[0, 1]], 3.0);
        assert_eq!(data.values[[1, 2]], 6.0);
        Ok(())
    }

    #[test]
    fn roundtrips_jzzp_text() -> Result<()> {
        let data = parse_jzzp_dat(JZZP_DAT)?;
        let rendered = jzzp_dat_string(&data)?;
        assert_eq!(parse_jzzp_dat(&rendered)?, data);
        Ok(())
    }

    #[test]
    fn lossless_jzzp_preserves_list_directed_layout() -> Result<()> {
        let text = "#          2           3           2           3\n\
#     1.0000000000000     3.1250000000000     2.0000000000000     4.0000000000000\n\
  1.0  2.0  3.0  4.0  5.0  6.0\n";
        let data = parse_jzzp_dat_lossless(text)?;
        assert_eq!(jzzp_dat_lossless_string(&data)?, text);
        Ok(())
    }

    #[test]
    fn rejects_bad_jzzp_inputs() {
        assert!(parse_jzzp_dat("").is_err());
        assert!(parse_jzzp_dat("# 2 3 2\n# 1 2 3 4\n1 2 3 4\n").is_err());
        assert!(parse_jzzp_dat("# 2 3 2 3\n# 1 2 3 4\n1 2 3\n").is_err());
        assert!(parse_jzzp_dat("# 2 3 2 3\n# 1 2 3 4\n1 2 3 4 5 NaN\n").is_err());

        let bad_shape = JzzpDatData {
            ns: 2,
            nphi: 3,
            nz: 2,
            nzp: 3,
            smax: 1.0,
            phimax: 3.125,
            zmax: 2.0,
            zpmax: 4.0,
            values: Array2::zeros((2, 2)),
        };
        assert!(jzzp_dat_string(&bad_shape).is_err());
    }

    const COMPTON_DAT: &str = r#" # Compton profile, J(pq)
 # ns:            32
 # nphi:          32
 # nz:            32
 # nzp:          120
 # zpmax:   10.0000000000000     
 # temperature (eV):  0.0000000E+00
 #----------------------------
 # pq               J
  0.000000000000000E+000   2.74476734850343     
  5.005004815757275E-003   2.74473136578831     
  1.001000963151455E-002   2.74462341659279     
"#;

    const RHOZZP_DAT: &str = concat!(
        "  9.999999776482582E-003   3.71096344005271     \n",
        "  2.001000978649259E-002   2.66921255682004     \n",
        "  3.002001979650260E-002   1.84694165446344     \n",
    );

    const JZZP_DAT: &str = r#"# 2 3 2 3
# 1.0 3.125 2.0 4.0
 1.0 2.0 3.0
 4.0 5.0 6.0
"#;
}