nereids-fitting 0.2.1

Optimization engine for resonance fitting (LM, Poisson/KL)
Documentation
//! Quickstart example used by `docs/guide/src/quickstart-rust.md`.
//!
//! Loads ENDF resonance data for U-238, computes a theoretical transmission
//! spectrum, and fits it to recover the areal density. Compile-checked by
//! `cargo check --workspace --examples` so the mdBook page cannot drift out
//! of sync with the live crate APIs.

// ANCHOR: setup
use nereids_core::types::Isotope;
use nereids_endf::parser::parse_endf_file2;
use nereids_endf::retrieval::{EndfLibrary, EndfRetriever, mat_number};
use nereids_fitting::lm::{LmConfig, levenberg_marquardt};
use nereids_fitting::parameters::{FitParameter, ParameterSet};
use nereids_fitting::transmission_model::TransmissionFitModel;
use nereids_physics::transmission::{SampleParams, forward_model};

fn main() -> Result<(), Box<dyn std::error::Error>> {
    // Download and cache ENDF/B-VIII.1 data for U-238 (Z=92, A=238).
    let isotope = Isotope::new(92, 238)?;
    let retriever = EndfRetriever::new();
    let mat = mat_number(&isotope, EndfLibrary::EndfB8_1).expect("U-238 has a known MAT number");
    let (_path, endf_text) = retriever.get_endf_file(&isotope, EndfLibrary::EndfB8_1, mat)?;
    let resonance_data = parse_endf_file2(&endf_text)?;

    println!(
        "U-238: {} resonances, AWR = {:.1}",
        resonance_data.total_resonance_count(),
        resonance_data.awr,
    );
    // ANCHOR_END: setup

    // ANCHOR: forward
    // Energy grid: 1 to 30 eV (covers the 6.67 eV and 20.9 eV resonances).
    let energies: Vec<f64> = (0..2000)
        .map(|i| 1.0 + (i as f64) * 29.0 / 2000.0)
        .collect();

    // Sample: U-238 at 0.001 atoms/barn, room temperature.
    let sample = SampleParams::new(300.0, vec![(resonance_data.clone(), 0.001)])?;

    // No instrument resolution broadening for this example.
    let transmission = forward_model(&energies, &sample, None)?;
    // `transmission[i]` is T(E_i) in [0, 1], with dips at resonance energies.
    // ANCHOR_END: forward

    // ANCHOR: fit
    // Simulate measured data (in practice, load from TIFF/NeXus).
    let measured_t = transmission.clone();
    let sigma: Vec<f64> = vec![0.01; measured_t.len()];

    // Set up the fit model: one density parameter at index 0.
    let model = TransmissionFitModel::new(
        energies.clone(),
        vec![resonance_data],
        300.0,                // temperature_k
        None,                 // no instrument resolution
        (vec![0], vec![1.0]), // density_indices, density_ratios
        None,                 // no temperature fitting
        None,                 // no precomputed cross-sections
    )?;

    // Initial guess: density = 0.0005 atoms/barn (non-negative constraint).
    let mut params = ParameterSet::new(vec![FitParameter::non_negative("U-238 density", 0.0005)]);

    let config = LmConfig::default();
    let result = levenberg_marquardt(&model, &measured_t, &sigma, &mut params, &config)?;

    println!("Fitted density: {:.6} atoms/barn", result.params[0]);
    println!("Reduced chi-squared: {:.3}", result.reduced_chi_squared);
    println!("Converged: {}", result.converged);
    Ok(())
}
// ANCHOR_END: fit