spitzenfinder 0.4.4

use std::error::Error;
use std::fs::File;
use std::io;
use std::{collections::HashMap, env};

use plotters::prelude::*;
use plotters::{
    element::Circle,
    style::{full_palette::PURPLE, Color, RED, ShapeStyle},
};
use user_query::BinaryQuery;
use unit_conversion::energy::{CONVERT_2_EV_FROM, CONVERT_2_RCM_FROM};

use spitzenfinder::data::FourierData;

fn main() -> Result<(), Box<dyn Error>> {
    // Get command line arguments
    let mut passed_args: env::Args = env::args();
    if passed_args.len() > 3 {
        print_usage();
        panic!(
            "Too many arguments were passed!
               Only name of fourier file and intensity theshold shall be passed."
        )
    }
    let input_file: String = passed_args.nth(1).unwrap();
    let mut threshold: f32 = 0.0_f32;
    if let Some(t) = passed_args.next() {
        threshold = t.parse::<f32>()?;
    }
    println!("Parsing file {}", input_file);
    println!("Setting threshold to {}", threshold);
    // Open the input file
    let file_handle: File = File::open(input_file)?;
    let mut fourier_data: FourierData = FourierData::new(file_handle, threshold).unwrap();
    let mean: f32 = FourierData::mean_intensities(&fourier_data);
    let msd: f32 = FourierData::variance_intensities(fourier_data.clone(), mean);
    let rmsd: f32 = msd.sqrt();
    let mad: f32 = FourierData::mean_abs_error_intensities(fourier_data.clone(), mean);
    // println!("Mean of intensities {}", mean);
    // println!(" MSD of intensities {}", msd);
    // println!("RMSD of intensities {}", rmsd);
    // println!(" MAD of intensities {}", mad);
    let mut local_maxima: HashMap<usize, (f32, f32)> = HashMap::new();
    match FourierData::all_local_maxima(&fourier_data) {
        Ok(maxima) => local_maxima = maxima,
        Err(e) => {
            eprintln!("Error '{}' was raised!", e);
            if let Ok(cons) = BinaryQuery::new(
                "At this point this likely means: \n
            Your data does not end with declining values but with inclining ones. \n
            I will continue by adding (0, (0.0, 0.0)) as the only local maximum is this okay?",
                "yes",
                "no",
                &mut io::stdin().lock(),
            ) {
                if cons.answer.eq("no") {
                    panic!("Ok. Bye then! :-)")
                } else {
                    local_maxima.insert(0, (0.0_f32, 0.0_f32));
                }
            }
        }
    };
    println!("Found {} maxima", local_maxima.len());
    println!("Max intensity {}", fourier_data.max_intensity);
    // Let's count how many values are larger than 4, 3, 2 and 1 times the variance.
    let mut gt_sigma_counter = [0_usize; 4];
    let mut lt_sigma_counter = [0_usize; 4];
    let (mut e_min, mut e_max): (f32, f32) = (0.0_f32, 0.0f32);
    for (e, i) in local_maxima.values() {
        if i > &(4.0_f32 * msd) {
            gt_sigma_counter[0] += 1;
            (e_min, e_max) = (e_min.min(*e), e_max.max(*e))
        } else {
            lt_sigma_counter[0] += 1
        };
        if i > &(3.0_f32 * msd) {
            gt_sigma_counter[1] += 1;
            (e_min, e_max) = (e_min.min(*e), e_max.max(*e))
        } else {
            lt_sigma_counter[0] += 1
        };
        if i > &(2.0_f32 * msd) {
            gt_sigma_counter[2] += 1;
            (e_min, e_max) = (e_min.min(*e), e_max.max(*e))
        } else {
            lt_sigma_counter[0] += 1
        };
        if i > &(1.0_f32 * msd) {
            gt_sigma_counter[3] += 1;
            (e_min, e_max) = (e_min.min(*e), e_max.max(*e))
        } else {
            lt_sigma_counter[0] += 1
        };
    }
    println!(
        "{} values are greater  and {} are smaller than four times the MSD.",
        gt_sigma_counter[0], lt_sigma_counter[0]
    );
    e_min += 0.05 * e_min;
    e_max += 0.05 * e_max;
    let message = format!(
        "Suggesting {} -- {} as x-axis range. D'accord?",
        e_min, e_max
    );
    let consent = BinaryQuery::new(&message, "yes", "no", &mut io::stdin().lock())?;
    if consent.answer.eq(consent.negative_choice) {
        if let Ok(x_lower) = BinaryQuery::new(
            "Enter lower bound:",
            "float to change",
            "k to keep",
            &mut io::stdin().lock(),
        ) {
            if x_lower.answer.ne("k") {
                e_min = x_lower.answer.parse::<f32>().unwrap()
            }
        }
        if let Ok(x_upper) = BinaryQuery::new(
            "Enter upper bound:",
            "float to change",
            "k to keep",
            &mut io::stdin().lock(),
        ) {
            if x_upper.answer.ne("k") {
                e_max = x_upper.answer.parse::<f32>().unwrap()
            }
        }
    }
    let mut energy_shift: f32 = 0.0_f32;
    let shift_query = BinaryQuery::new(
        "Do you wish to enter an energy shift?",
        "yes",
        "no",
        &mut io::stdin().lock(),
    )?;
    if shift_query.answer.eq(shift_query.positive_choice) {
        let user_message = "Please enter shift with unit ( rcm for wavenumbers ) will be converted to rcm automatically.";
        let unser_input = BinaryQuery::new(
            user_message,
            "Float x.x [unit]",
            "k to keep unshifted values.",
            &mut io::stdin().lock(),
        )?;
        if unser_input.answer.ne("k") {
            let mut user_shift = unser_input.answer.split(' ');
            if let Some(energy) = user_shift.next() {
                energy_shift = energy.trim().parse::<f32>()?;
            }
            if let Some(unit) = user_shift.next() {
                match CONVERT_2_RCM_FROM.contains_key(unit.trim()) {
                    true => {
                        energy_shift = CONVERT_2_RCM_FROM[unit.trim()](energy_shift as f64) as f32
                    }
                    // Need to give the user the possibility to re-enter a unit. todo!()
                    false => {
                        panic!("Cannot understand {}", unit)
                    }
                }
            }
        }
    }

    FourierData::shift_energies(&mut fourier_data, energy_shift);
    let max_energy = fourier_data.max_intensity.energy;
    let max_intensity = fourier_data.max_intensity.intensity;
    let mut distances_to_max: Vec<(f32, f32)> = local_maxima
        .values()
        .filter(|(energy, _intensity)| e_min < *energy && *energy < e_max)
        .map(|&(energy, intensity)| (energy - max_energy, intensity))
        .collect();
    distances_to_max.sort_by(|&tuple1, tuple2| tuple2.0.total_cmp(&tuple1.0));

    // let canvas = BitMapBackend::new("fourier_plot.png", (1920, 1080)).into_drawing_area();
    let canvas = SVGBackend::new("fourier_plot.svg", (1920, 1080)).into_drawing_area();
    canvas.fill(&WHITE)?;

    let x_range_rcm = e_min + energy_shift..e_max + energy_shift;
    let energy_shift_ev = CONVERT_2_EV_FROM["rcm"](energy_shift as f64) as f32;
    let lower_bound_ev = CONVERT_2_EV_FROM["rcm"](e_min as f64) as f32;
    let upper_bound_ev = CONVERT_2_EV_FROM["rcm"](e_max as f64) as f32;
    let x_range_ev = (lower_bound_ev+energy_shift_ev)..(upper_bound_ev+energy_shift_ev);

    let canvas = canvas.margin(10, 10, 10, 10);
    let mut chart = ChartBuilder::on(&canvas)
        .x_label_area_size(60)
        .top_x_label_area_size(60)
        .y_label_area_size(60)
        .build_cartesian_2d(
            x_range_rcm,
            0.0_f32..max_intensity,
        )?
        .set_secondary_coord(
            x_range_ev,
            0.0_f32..max_intensity,
        );
    chart
        .configure_mesh()
        .x_desc("Energy in wavenumbers")
        .y_desc("Intensity in arbitrary units")
        .draw()?;

    let line_style = ShapeStyle{ color: Into::<RGBAColor>::into(RED), filled: false, stroke_width: 4 };
    chart.draw_series(LineSeries::new(fourier_data, line_style))?;

    chart.draw_series(
        local_maxima
            .values()
            .filter(|(energy, _intensity)| e_min < *energy && *energy < e_max)
            .map(|(energy, intensity)| {
                Circle::new((*energy + energy_shift, *intensity), 3, PURPLE.filled())
            }),
    )?;
    let difflabel = |energy: f32, intentsity: f32, label: f32| {
        return EmptyElement::at((energy, intentsity))
            + Text::new(
                format!("{:.0}", label),
                (10, 10),
                ("sans-serif", 16)
                    .into_font()
                    .transform(FontTransform::Rotate90),
            );
    };
    distances_to_max.iter().for_each(|&(energy, intensity)| {
        let plot_data = [
            (max_energy, max_intensity),
            (energy + max_energy + energy_shift, max_intensity),
            (energy + max_energy + energy_shift, intensity),
        ];
        let knick: f32 = max_energy + energy + energy_shift;
        let _ = &chart.draw_series(LineSeries::new(plot_data, &BLACK));
        let _ = &chart
            .plotting_area()
            .draw(&difflabel(knick, max_intensity, energy));
    });
    chart.configure_secondary_axes().x_desc("Energy in eV").draw()?;
    canvas.present()?;
    Ok(())
}

fn print_usage() {
    println!("Usage: peak_finder [input_file] [[threshold as float optional]]")
}