imzml/

lib.rs

#![warn(missing_docs)]
#![warn(rustdoc::missing_doc_code_examples)]

//! This library provides a means of accessing mass spectrometry and mass spectrometry imaging data stored in either the
//! .mzML or .imzML data formats.
//!
//! let parser = ImzMLReader::from_path("/path/to/data.imzML").unwrap();
//!
//! for error in parser.errors() {
//!     println!("{:?}", error);
//! }
//!
//! let imzml: ImzML<_> = parser.into();
//!
//! imzml.ion_image(772.573, 100.0);

/// Error handling
pub mod error;
/// ImzML file parsing (relies heavily on `mzml`)
pub mod imzml;
/// MzML file parsing
pub mod mzml;

/// Parsing ontologies stored in the OBO format
pub mod obo;
/// Validating .mzML and .imzML files
pub mod validation;
//mod writer;

use std::io::ErrorKind;
use std::io::Read;
use std::io::Seek;
use std::io::SeekFrom;
use std::sync::Arc;
use std::sync::Mutex;

use byteorder::LittleEndian;
use byteorder::ReadBytesExt;
use error::DataError;
use error::FatalParseError;
use error::ParseError;
use mzml::Tag;

pub use self::imzml::ImzML;
pub use self::imzml::ImzMLReader;
pub use mzml::mzml::MzML;
pub use mzml::MzMLReader;

pub use self::mzml::binarydataarray::*;
pub use self::mzml::mzml::*;
//pub use self::parser::*;
pub use self::mzml::scan::*;
pub use self::mzml::spectrum::*;

#[macro_use]
extern crate lazy_static;

const BUFFER_SIZE: usize = 0x1000;

/// Enables access to data stored in `data_location`. The data could be compressed and/or encoded
/// and stored with various different types. This information is captured in the `data_description`
#[derive(Debug)]
pub struct DataAccess<D: Read + Seek> {
    data_location: Arc<Mutex<D>>,

    array_length: Option<u64>,
    encoded_length: u64,
    offset: u64,

    data_description: DataDescription,
}

/// WritableSpectrum contains information required to be able to write out a spectrum to imzML.
pub struct WritableSpectrum {
    /// Spectrum metadata
    pub spectrum: Spectrum,
    /// m/z values
    pub mz_array: Option<Vec<f64>>,
    /// Intensity values
    pub intensity_array: Option<Vec<f64>>,
}

// impl<'obo> WritableSpectrum<'obo> {
//     fn description(&mut self) -> &mut Spectrum<'obo> {
//         &mut self.spectrum
//     }

//     fn mz_array(&self) -> Option<&Vec<f64>> {
//         self.mz_array.as_ref()
//     }

//     fn intensity_array(&self) -> Option<&Vec<f64>> {
//         self.intensity_array.as_ref()
//     }
// }
// pub trait DataAccess {
//     fn as_f64(&self) -> Result<Vec<f64>, DataError>;
// }

// pub trait SpectrumAccess<'obo> {
//     type DataAccess: DataAccess;

//     fn spectrum(&self) -> &Spectrum<'obo>;
//     fn mz_array(&self) -> Option<&Self::DataAccess>;
//     fn intensity_array(&self) -> Option<&Self::DataAccess>;
// }

impl<D: Read + Seek> DataAccess<D> {
    /// Create a DataAccess instance from a `BinaryDataArray`, using the specified `data_location`
    pub fn from_binary_data_array(
        bda: &BinaryDataArray,
        data_location: Arc<Mutex<D>>,
    ) -> Option<Self> {
        if let (Some(offset), Some(encoded_length)) = (bda.offset(), bda.encoded_length()) {
            Some(DataAccess {
                data_location,
                encoded_length,
                offset,
                array_length: bda.array_length(),
                data_description: DataDescription {
                    is_base64_encoded: !bda.is_data_external(),
                    compression: bda.compression(),
                    binary_type: bda.binary_type(),
                },
            })
        } else {
            None
        }
    }

    fn reader(&self) -> Result<CompressionReader, DataError> {
        let mut reader = self.data_location.lock()?;

        reader.seek(SeekFrom::Start(self.offset))?;
        let mut buf = vec![0; self.encoded_length as usize];
        reader.read_exact(&mut buf)?;

        // We can free the guard now before we move on to decoding
        drop(reader);

        if self.data_description.is_base64_encoded {
            // Base64 decoding
            buf = base64::decode(&buf)?;
        }

        Ok(self.data_description.compression.to_reader(buf))
    }

    /// Return the offset (in the file) of the data
    pub fn offset(&self) -> u64 {
        self.offset
    }

    /// Return the length (in bytes) of the data in the file
    pub fn encoded_length(&self) -> u64 {
        self.encoded_length
    }

    /// Returns the binary type of the data (e.g. 64-bit float)
    pub fn binary_type(&self) -> BinaryDataType {
        self.data_description.binary_type()
    }

    /// Decode, decompress and convert (if necessary) data points to f64
    pub fn as_f64(&self) -> Result<Vec<f64>, DataError> {
        let mut reader = self.reader()?;

        let mut result = match self.array_length {
            Some(array_length) => Vec::with_capacity(array_length as usize),
            None => Vec::new(),
        };

        match self.data_description.binary_type {
            BinaryDataType::Undefined => Err(DataError::UnknownDataType),
            BinaryDataType::Float64 => {
                loop {
                    match reader.read_f64::<LittleEndian>() {
                        Ok(value) => result.push(value),
                        Err(error) => {
                            if error.kind() == ErrorKind::UnexpectedEof {
                                break;
                            }

                            return Err(DataError::IOError(error));
                        }
                    }
                }

                Ok(result)
            }
            BinaryDataType::Float32 => {
                loop {
                    match reader.read_f32::<LittleEndian>() {
                        Ok(value) => result.push(value as f64),
                        Err(error) => {
                            if error.kind() == ErrorKind::UnexpectedEof {
                                break;
                            }

                            return Err(DataError::IOError(error));
                        }
                    }
                }

                Ok(result)
            }
        }
    }
}

// #[derive(Debug)]
// pub enum DataAccess<D: BufRead + Seek> {
//     OnDisk(OnDiskAccess<D>),
//     InMemory(Vec<u8>),
// }

/// Provide ability to access a spectrum (contains `DataAccess` for both m/z and intensity arrays)
pub struct SpectrumAccess<D: Read + Seek> {
    spectrum: Arc<Spectrum>,

    mz_array: Option<DataAccess<D>>,
    intensity_array: Option<DataAccess<D>>,
}

// impl<'obo, D: BufRead + Seek> SpectrumAccess<'obo, D> {
//     type DataAccess = DataAccess<D>;

//     /// Returns the `DataAccess` for the m/z array
//     fn mz_array(&self) -> Option<&DataAccess<D>> {
//         self.mz_array.as_ref()
//     }

//     /// Returns the `DataAccess` for the intensity array
//     fn intensity_array(&self) -> Option<&DataAccess<D>> {
//         self.intensity_array.as_ref()
//     }

//     fn spectrum(&self) -> &Spectrum<'obo> {
//         &self.spectrum
//     }
// }

impl<D: Read + Seek> SpectrumAccess<D> {
    /// Creates a new `SpectrumAccess` from the specified `data_location` and `Spectrum`
    pub fn new(data_location: Arc<Mutex<D>>, spectrum: Arc<Spectrum>) -> Self {
        let mz_array = match spectrum.mz_array() {
            Some(mz_array) => DataAccess::from_binary_data_array(mz_array, data_location.clone()),
            None => None,
        };
        let intensity_array = match spectrum.intensity_array() {
            Some(intensity_array) => {
                DataAccess::from_binary_data_array(intensity_array, data_location)
            }
            None => None,
        };

        Self {
            spectrum,
            mz_array,
            intensity_array,
        }
    }

    /// Returns the `DataAccess` for the m/z array
    pub fn mz_array(&self) -> Option<&DataAccess<D>> {
        self.mz_array.as_ref()
    }

    /// Returns the `DataAccess` for the intensity array
    pub fn intensity_array(&self) -> Option<&DataAccess<D>> {
        self.intensity_array.as_ref()
    }

    /// Returns the spectrum metadata
    pub fn spectrum(&self) -> &Spectrum {
        &self.spectrum
    }

    /// Returns the spectral representation (centroid, profile) of this spectrum
    pub fn representation(&self) -> Option<Representation> {
        self.spectrum.representation()
    }

    /// Returns the polarity of this spectrum
    pub fn polarity(&self) -> Option<Polarity> {
        self.spectrum.polarity()
    }

    /// Returns the pixel coordinate of the spectrum
    pub fn coordinate(&self) -> Option<Coordinate> {
        let x = self.spectrum.x_position()?;
        let y = self.spectrum.y_position()?;

        Some(Coordinate { x, y, z: None })
    }
}

/// ScanLocation describes how ImzMLWriter should deal with a particular scan.
/// Either ignore it, or store it with the specified coordinate.
pub enum ScanLocation {
    /// Ignore the scan
    Ignore,
    /// Store the scan with the specified coordinate.
    Location(Coordinate),
}

/// Coordinate describes the pixel location of a spectrum/scan
#[derive(Debug, Hash, PartialEq, Eq)]
pub struct Coordinate {
    /// x-axis coordinate
    pub x: u32,
    /// y-axis coordinate
    pub y: u32,
    /// Optional z-axis coordinate
    pub z: Option<u32>,
}

#[derive(Debug)]
struct DataDescription {
    is_base64_encoded: bool,

    compression: Compression,
    binary_type: BinaryDataType,
}

impl Default for DataDescription {
    fn default() -> Self {
        Self {
            is_base64_encoded: false,
            compression: Compression::None,
            binary_type: BinaryDataType::Float64,
        }
    }
}

impl DataDescription {
    pub fn binary_type(&self) -> BinaryDataType {
        self.binary_type
    }
}

// impl<'obo, 'iter, B: BufRead, D: BufRead + Seek> Iterator for &ImzMLWriter<'obo, 'iter, B, D>  {
//     type Item = Coordinate;

//     fn next(&mut self) -> Option<Self::Item> {
//         self.coord_iterator.lock().unwrap().next()
//     }
// }

// pub trait SpectrumAccessIterator<'obo> {
//     type DataAccess: DataAccess;
//     type SpectrumAccess: SpectrumAccess<'obo, DataAccess = Self::DataAccess>;
//     type AccessIterator: Iterator<Item = Self::SpectrumAccess>;

//     fn iter_access(self) -> Self::AccessIterator;
// }

/// SpectrumAccessIterator provides an iterator over a spectrum
//#[derive(Clone, Copy)]
pub struct SpectrumAccessIterator<D: Read + Seek> {
    data_location: Arc<Mutex<D>>,
    //mzml: MzML<'obo>,
    spectra: Vec<Arc<Spectrum>>,

    next_index: usize,
}

// impl<'obo, D: BufRead + Seek> SpectrumAccessIterator<'obo> for DiskSpectrumAccessIterator<'obo, D> {
//     type DataAccess = DiskDataAccess<D>;
//     type SpectrumAccess = DiskSpectrumAccess<'obo, D>;
//     type AccessIterator = Self;

//     fn iter_access(self) -> <Self as SpectrumAccessIterator<'obo>>::AccessIterator {
//         self
//     }
// }

impl<D: Read + Seek> SpectrumAccessIterator<D> {
    fn new(data_location: Arc<Mutex<D>>, spectrum_list: &SpectrumList) -> Self {
        Self {
            data_location,
            spectra: spectrum_list.spectra().to_vec(),
            next_index: 0,
        }
    }
}

impl<D: Read + Seek> Iterator for SpectrumAccessIterator<D> {
    type Item = SpectrumAccess<D>;

    fn next(&mut self) -> Option<Self::Item> {
        let cur_index = self.next_index;
        self.next_index += 1;

        self.spectra
            .get(cur_index)
            .map(|spectrum| SpectrumAccess::new(self.data_location.clone(), spectrum.clone()))
    }
}

// impl<'a, B: BufRead, D: BufRead + Seek> AsRef<HeaderParser<'a, B>> for Parser<'a, B, D> {
//     fn as_ref(&self) -> &HeaderParser<'a, B> {
//         &self.header_parser
//     }
// }

// impl<'a, B: BufRead, D: BufRead + Seek> std::ops::Deref for ImzMLParser<'a, B, D> {
//     type Target = MzMLParser<'a, B>;

//     fn deref(&self) -> &Self::Target {
//         &self.header_parser
//     }
// }

#[cfg(test)]
mod tests {
    use crate::imzml::ImzMLReader;
    use crate::mzml::writer::Writer;
    use crate::mzml::{DataReader, MzMLReader, MzMLTag};

    use std::fs::File;
    use std::io::{BufReader, BufWriter};
    use std::time::Instant;

    #[test]
    fn test_all_in_test_folder() -> std::io::Result<()> {
        let paths = std::fs::read_dir("../test/")?;

        //let ontology = crate::obo::parser::parse("src/obo/imagingMS.obo").unwrap();

        for path in paths {
            let path = path?.path();
            let extension = path.extension().unwrap();
            let extension = extension.to_str().unwrap();

            if !(extension == "imzML" || extension == "mzML") {
                continue;
            }

            println!("Processing {:?}", path);

            let start = Instant::now();

            match extension {
                "imzML" => {
                    let parser = ImzMLReader::from_path(path).unwrap();

                    let duration = start.elapsed();
                    println!("Time elapsed when parsing imzML is: {:?}", duration);

                    for error in parser.errors() {
                        println!("{}", error);
                    }
                }
                "mzML" => {
                    let parser = MzMLReader::from_path(path).unwrap();

                    let duration = start.elapsed();
                    println!("Time elapsed when parsing imzML is: {:?}", duration);

                    for error in parser.errors() {
                        println!("{}", error);
                    }
                }
                _ => {}
            };

            // let parser = MzMLReader::new(header_reader).unwrap();
            //let mut parser = parser::ImzMLParser::from_file(path.path(), Some(&ontology)).unwrap();

            // let start = Instant::now();
            // let mut parser = parser::ImzMLParser::from_file(path, Some(&ontology)).unwrap();

            // let duration = start.elapsed();
            // println!("Time elapsed when parsing (i)mzML is: {:?}", duration);

            // let mut i = 0;
            // while parser.has_errors() {
            //     println!("{:?}", parser.pop_error_front().unwrap());

            //     if i > 10 {
            //         break;
            //     }

            //     i += 1;
            // }
        }

        Ok(())
    }

    #[test]
    fn read_mzml_in_test_folder() -> std::io::Result<()> {
        let paths = std::fs::read_dir("../test/")?;

        //let ontology = crate::obo::parser::parse("src/obo/imagingMS.obo").unwrap();

        for path in paths {
            let path = path?;

            let path = path.path();
            let extension = path.extension().unwrap();
            let extension = extension.to_str().unwrap();

            if extension != "mzML" {
                continue;
            }

            println!("Processing [data] {:?}", path);

            let start = Instant::now();

            let parser = DataReader::from_path(path).unwrap();
            //let mut parser = parser::ImzMLParser::from_file(path.path(), Some(&ontology)).unwrap();

            let duration = start.elapsed();
            println!("Time elapsed when parsing mzML is: {:?}", duration);

            let spectrum = parser.spectrum(0).unwrap();

            let mzs = spectrum.mz_array().unwrap().as_f64().unwrap();

            println!("{:?}", &mzs[..10]);

            // let start = Instant::now();
            // let mut parser = parser::ImzMLParser::from_file(path, Some(&ontology)).unwrap();

            // let duration = start.elapsed();
            // println!("Time elapsed when parsing (i)mzML is: {:?}", duration);

            // let mut i = 0;
            // while parser.has_errors() {
            //     println!("{:?}", parser.pop_error_front().unwrap());

            //     if i > 10 {
            //         break;
            //     }

            //     i += 1;
            // }
        }

        Ok(())
    }

    #[ignore]
    #[test]
    fn test_all_in_folder() -> std::io::Result<()> {
        let paths = std::fs::read_dir(
            "/home/alan/Documents/GitProjects/gomsi/metaspace/cmd/metaspacedownloader/output/",
        )
        .unwrap();

        //let ontology = crate::obo::parser::parse("imagingMS.obo").unwrap();

        for path in paths {
            let path = path?;

            if path.path().extension().unwrap() != "imzML" {
                println!("Skipping {:?}.", path.path());
                continue;
            }

            println!("Processing {:?}", path.path());

            let start = Instant::now();

            let file = File::open(path.path()).unwrap();
            let header_reader = BufReader::new(file);

            let parser = MzMLReader::new(header_reader).unwrap();
            //let mut parser = parser::ImzMLParser::from_file(path.path(), Some(&ontology)).unwrap();

            let duration = start.elapsed();
            println!("Time elapsed when parsing imzML is: {:?}", duration);

            for error in parser.errors() {
                println!("{}", error);
            }
        }

        Ok(())
    }

    #[ignore]
    #[test]
    fn test_obo() {
        let start = Instant::now();
        let ontology = crate::obo::parser::parse("imagingMS.obo").unwrap();
        let duration = start.elapsed();
        println!("Time elapsed when parsing ontology is: {:?}", duration);

        let filename = "/home/alan/Documents/GitProjects/gomsi/metaspace/cmd/metaspacedownloader/output/2017-06-30_07h26m26s_Mousebrain_MG08_2017_GruppeF.imzML";
        println!("[Parser] Reading with ontology");
        let parser = ImzMLReader::from_path_with_ontology(filename, ontology.clone()).unwrap();

        let duration = start.elapsed();
        println!("Time elapsed when parsing imzML is: {:?}", duration);

        for error in parser.errors() {
            println!("{}", error);
        }

        let imzml = parser.mzml().unwrap();

        // crate::validation::schema::validate(&imzml);

        let mut errors = Vec::new();

        let mapping = crate::validation::parse("ms-mapping.xml").unwrap();
        mapping.validate(&ontology, imzml, &mut errors);

        let mapping = crate::validation::parse("Ims1.1-mapping.xml").unwrap();
        mapping.validate(&ontology, imzml, &mut errors);

        let file = File::create("tmp.xml").unwrap();

        let mut writer = Writer::new(BufWriter::new(file)).unwrap(); //new(Cursor::new(Vec::new()));

        parser.mzml().unwrap().write_xml(&mut writer).unwrap();
    }

    #[test]
    #[ignore]
    fn test_parser() {
        //let filename = "/home/alan/Documents/GitProjects/gomsi/metaspace/cmd/metaspacedownloader/output/2016-09-22_11h16m17s_Brain01_Bregma-1-46_centroid.imzML";
        let filename = "/home/alan/Documents/GitProjects/gomsi/metaspace/cmd/metaspacedownloader/output/2017-06-30_07h26m26s_Mousebrain_MG08_2017_GruppeF.imzML";
        //let filename = "/home/alan/Documents/GitProjects/gomsi/metaspace/cmd/metaspacedownloader/output/2019-03-28_09h18m57s_brain 25µm raw - root mean square - metaspace.imzML";
        //let filename = "/home/alan/Documents/GitProjects/gomsi/metaspace/cmd/metaspacedownloader/output/2021-01-16_08h47m30s_peptide.imzML";
        //let filename = "/home/alan/Documents/GitProjects/gomsi/metaspace/cmd/metaspacedownloader/output/2016-09-21_16h06m53s_Image5.imzML";

        let start = Instant::now();
        let parser = ImzMLReader::from_path(filename).unwrap();
        let duration = start.elapsed();

        for error in parser.errors() {
            println!("{}", error);
        }

        //println!("width: {:?}", parser.current_mzml.unwrap().get_width());

        let file = File::create("tmp.xml").unwrap();

        let mut writer = Writer::new(BufWriter::new(file)).unwrap(); //new(Cursor::new(Vec::new()));

        parser.mzml().unwrap().write_xml(&mut writer).unwrap();

        println!("Time elapsed when parsing Mousebrain is: {:?}", duration);

        println!("{:?}", parser.mzml().unwrap().spectrum(0));
        println!(
            "{:?}",
            parser.mzml().unwrap().spectrum(0).unwrap().x_position()
        );
        println!(
            "{:?}",
            parser
                .mzml()
                .unwrap()
                .referenceable_param_group_ref("mzArray")
        );
    }
}