use std::{
fs,
fs::OpenOptions,
io::{self, Write},
path::Path,
process::Command,
};
use crate::{
FrameSlice,
dcd::{DcdFrame, DcdMetadata, DcdTrajectory, read_dcd},
};
fn run_mdconvert(input: &Path, output: &Path) -> io::Result<()> {
let out = Command::new("mdconvert")
.arg(input)
.arg("-o")
.arg(output)
.arg("--force")
.output()
.map_err(|e| {
if e.kind() == io::ErrorKind::NotFound {
io::Error::other("mdconvert not found. Install mdtraj: pip install mdtraj")
} else {
e
}
})?;
if !out.status.success() {
return Err(io::Error::other(format!(
"mdconvert failed: {}",
String::from_utf8_lossy(&out.stderr)
)));
}
Ok(())
}
pub struct XtcMetadata {
pub num_atoms: usize,
pub num_frames: usize,
pub start_time: f32,
pub end_time: f32,
pub dt: f32,
}
impl XtcMetadata {
pub fn read(path: &Path) -> io::Result<Self> {
let pid = std::process::id();
let tmp = std::env::temp_dir();
let dcd_path = tmp.join(format!("bio_xtc_meta_{pid}.dcd"));
run_mdconvert(path, &dcd_path)?;
let md = DcdMetadata::read(&dcd_path);
let _ = fs::remove_file(&dcd_path);
let md = md?;
let start_time = md.start_step * md.dt;
let dt = if md.num_frames > 1 {
md.save_interval_steps as f32 * md.dt
} else {
0.0
};
Ok(Self {
num_atoms: md.num_atoms,
num_frames: md.num_frames,
start_time,
end_time: md.end_time,
dt,
})
}
}
pub fn read_xtc(path: &Path, slice: FrameSlice) -> io::Result<Vec<DcdFrame>> {
let pid = std::process::id();
let tmp = std::env::temp_dir();
let dcd_path = tmp.join(format!("bio_xtc_read_{pid}.dcd"));
run_mdconvert(path, &dcd_path)?;
let frames = read_dcd(&dcd_path, slice);
let _ = fs::remove_file(&dcd_path);
frames
}
pub fn write_xtc(path: &Path, frames: &[DcdFrame]) -> io::Result<()> {
if frames.is_empty() {
return Ok(());
}
let n_atoms = frames[0].atom_posits.len();
if n_atoms == 0 {
return Ok(());
}
for f in frames.iter().skip(1) {
if f.atom_posits.len() != n_atoms {
return Err(io::Error::new(
io::ErrorKind::InvalidInput,
"inconsistent atom counts across frames",
));
}
}
let pid = std::process::id();
let tmp = std::env::temp_dir();
let dcd_path = tmp.join(format!("bio_xtc_write_{pid}.dcd"));
let xtc_tmp = tmp.join(format!("bio_xtc_write_{pid}.xtc"));
let traj = DcdTrajectory {
frames: frames.to_vec(),
};
traj.save(&dcd_path)?;
let result = run_mdconvert(&dcd_path, &xtc_tmp);
let _ = fs::remove_file(&dcd_path);
result?;
let chunk = fs::read(&xtc_tmp);
let _ = fs::remove_file(&xtc_tmp);
let mut out = OpenOptions::new().create(true).append(true).open(path)?;
out.write_all(&chunk?)?;
Ok(())
}