version up

Cargo.lock

@@ -362,7 +362,7 @@ checksum = "cccddf32554fecc6acb585f82a32a72e28b48f8c4c1883ddfeeeaa96f7d8e519"
 
 [[package]]
 name = "wham"
-version = "1.1.0"
+version = "1.1.2"
 dependencies = [
  "assert_approx_eq",
  "clap",

Cargo.toml

@@ -1,6 +1,6 @@
 [package]
 name = "wham"
-version = "1.1.0"
+version = "1.1.2"
 authors = ["Daniel Bauer <bauer@cbs.tu-darmstadt.de>"]
 description = "An implementation of the weighted histogram analysis method"
 license = "GPL-3.0"

src/cli.yml

@@ -1,5 +1,5 @@
 name: wham
-version: "1.1.0"
+version: "1.1.2"
 author: D. Bauer <bauer@cbs.tu-darmstadt.de>
 about: |
   wham is a fast implementation of the weighted histogram analysis method (WHAM) written in Rust. It currently supports potential of mean force (PMF) calculations in multiple dimensions at constant temperature.
@@ -112,3 +112,8 @@ args:
         help: "Performs WHAM for slices with the given delta in time and returns an output file for each slice. THis is useful to check the result for convergence. Example: with --convdt 100 and a timeseries ranging from 0-300, free energy surfaces for slices 0-100, 0-200 and 0-300 will be given returned."
         takes_value: true
         required: false
+    - ignore_empty:
+        long: ignore_empty
+        help: If this is set, do not fail if a histogram is empty.
+        takes_value: false
+        required: false

src/io.rs

@@ -103,19 +103,24 @@ pub fn read_data(cfg: &Config) -> Result<Vec<Dataset>> {
         }
     }
 
-    // Histograms are stored as timeseries x convdt right now,
-    // but we need convdt x timeseries to create Datasets
-    // this transposes the data
+    // Datasets are created from histograms.
+    // Empty histograms result in an error when its the final dataset,
+    // and a warning otherwise.
     let num_datasets: usize = histograms.iter().map(|h| h.len()).max().unwrap();
+    let dataset_boundaries: Vec<(f64, f64)> = (0..num_datasets).map(|idx| {
+        (cfg.start, cfg.start+(idx as f64 + 1.0)*cfg.convdt) }
+    ).collect();
+    vprintln(format!("Generating {} datasets from histograms.", num_datasets), cfg.verbose);
     let datasets: Vec<Dataset> = (0..num_datasets).map(|idx| {
         let mut dataset_histograms: Vec<Histogram> = Vec::with_capacity(histograms.len());
         for (hs, path) in histograms.iter().zip(&paths) {
             if hs.len() > idx {
                 dataset_histograms.push(hs[idx].clone())
             } else {
-                let warning = format!("No data points in histogram boundaries: {}", &path);
-                if idx+1 == num_datasets {
-                    bail!(warning);
+                let warning = format!("No data points for interval {}-{} in histogram boundaries: {}.",
+                    dataset_boundaries[idx].0, dataset_boundaries[idx].1 ,&path);
+                if !cfg.ignore_empty && idx+1 == num_datasets {
+                    bail!(warning + " This is the final dataset.");
                 } else {
                     eprintln!("{}", warning);
                 }
@@ -172,12 +177,10 @@ fn get_convdt_boundaries(timeseries: &[f64], cfg: &Config) -> Vec<(f64, f64)> {
     if first_timestep < cfg.start {
         first_timestep = cfg.start;
     }
-    println!("{} to {} with dt={}", first_timestep, last_timestep, cfg.convdt);
     if cfg.convdt == 0.0 {
         vec![(0.0, last_timestep)]
     } else {
         let intervals: usize = ((last_timestep - first_timestep) / cfg.convdt).ceil() as usize;
-        println!("{:?}", intervals);
         (1..intervals+1).map(|i| {
             i as f64 * cfg.convdt + first_timestep
         }).map(|end| { (first_timestep, end) }).collect()

src/lib.rs

@@ -48,6 +48,7 @@ pub struct Config {
     pub end: f64,
     pub uncorr: bool,
     pub convdt: f64,
+    pub ignore_empty: bool
 }
 
 impl fmt::Display for Config {
@@ -55,11 +56,11 @@ impl fmt::Display for Config {
          write!(f, "Metadata={}, hist_min={:?}, hist_max={:?}, bins={:?}, 
             verbose={}, tolerance={}, iterations={}, temperature={},
             cyclic={:?}, uncorr={:?}, bootstrap={:?}, seed={:?},
-            uncorr={:?}, start={:?}, end={:?}, convdt={:?}",
+            uncorr={:?}, start={:?}, end={:?}, convdt={:?}, ignore_empty={:?}",
             self.metadata_file, self.hist_min, self.hist_max, self.num_bins,
             self.verbose, self.tolerance, self.max_iterations, self.temperature,
             self.cyclic, self.uncorr, self.bootstrap, self.bootstrap_seed,
-            self.uncorr, self.start, self.end, self.convdt)
+            self.uncorr, self.start, self.end, self.convdt, self.ignore_empty)
     }
 }
 

src/main.rs

@@ -71,9 +71,11 @@ fn cli() -> Result<Config> {
     let convdt: f64 = matches.value_of("convdt").unwrap_or("0").parse()
         .chain_err(|| "Cannot parse convdt.")?;
 
+    let ignore_empty: bool = matches.is_present("ignore_empty");
+
 	Ok(wham::Config{metadata_file, hist_min, hist_max, num_bins, dimens,
 		verbose, tolerance, max_iterations, temperature, cyclic, output,
-		bootstrap, bootstrap_seed, start, end, uncorr, convdt})
+		bootstrap, bootstrap_seed, start, end, uncorr, convdt, ignore_empty})
 }
 
 fn main() {

tests/examples.rs

@@ -66,4 +66,4 @@ mod integration {
         let output_len = String::from_utf8_lossy(&output.stdout).len();
         assert_eq!(output_len, 0);
     }
-}
+}