store final PMF in a file

README.md

@@ -4,4 +4,4 @@ TODO
 - Error analysis
 - Better error handling (file read/write)
 - Better File I/O performance (file read in slow)
-- Output to file
+- integration tests

src/cli.yml

@@ -54,4 +54,10 @@ args:
         help: WHAM temperature in Kelvin
         takes_value: true
         required: true
+    - output:
+        short: o
+        long: output
+        help: Free energy output file (defaults to "wham.out")
+        takes_value: true
+        required: false
         

src/io.rs

@@ -8,6 +8,8 @@ use k_B;
 use std::process;
 use std::option::Option;
 use std::path::Path;
+use std::error::Error;
+use std::result::Result;
 
 // Returns the path to path2 relative to path1
 // path1: "path/to/file.dat"
@@ -133,6 +135,16 @@ fn read_window_file(window_file: &str, cfg: &Config) -> Option<Histogram> {
     }
 }
 
+pub fn write_results(out_file: &str, ds: &Dataset, free: &Vec<f32>, prob: &Vec<f32>) -> Result<(), Box<Error>> {
+     let mut output = File::create(out_file)?;
+     writeln!(output, "#{:8}\t{:8}\t{:8}", "x", "Free Energy", "Probability");
+     for bin in 0..free.len() {
+        let x = ds.get_x_for_bin(bin);
+        writeln!(output, "{:8.6}\t{:8.6}\t{:8.6}", x, free[bin], prob[bin])?;
+     }
+     Ok(())
+}
+
 #[cfg(test)]
 mod tests {
     use super::*;
@@ -147,7 +159,8 @@ mod tests {
             tolerance: 0.0,
             max_iterations: 0,
             temperature: 300.0,
-            cyclic: false
+            cyclic: false,
+            output: "qwert".to_string(),
         }
     }
 

src/lib.rs

@@ -27,13 +27,14 @@ pub struct Config {
 	pub max_iterations: usize,
 	pub temperature: f32,
 	pub cyclic: bool,
+	pub output: String,
 }
 
 impl fmt::Display for Config {
 	 fn fmt(&self, f: &mut fmt::Formatter) -> fmt::Result {
-	 	write!(f, "Metadata={}, hist_min={}, hist_max={}, bins={}\nverbose={}, tolerance={}, iterations={}, temperature={}" , self.metadata_file, self.hist_min,
+	 	write!(f, "Metadata={}, hist_min={}, hist_max={}, bins={}\nverbose={}, tolerance={}, iterations={}, temperature={}, cyclic={}" , self.metadata_file, self.hist_min,
 	 		self.hist_max, self.num_bins, self.verbose, self.tolerance,
-	 		self.max_iterations, self.temperature)
+	 		self.max_iterations, self.temperature, self.cyclic)
     }
 }
 
@@ -170,14 +171,15 @@ fn free_energy(ds: &Dataset, P: &mut Vec<f32>, A: &mut Vec<f32>) {
 	for bin in 0..ds.num_bins {
 		A[bin] -= bin_min;
 	}
+
 	// Normalize P
-	// let mut P_sum = 0.0;
-	// for bin in 0..ds.num_bins {
-	// 	P_sum += P[bin];
-	// }
-	// for bin in 0..ds.num_bins {
-	// 	P[bin] /= P_sum;
-	// }
+	let mut P_sum = 0.0;
+	for bin in 0..ds.num_bins {
+		P_sum += P[bin];
+	}
+	for bin in 0..ds.num_bins {
+		P[bin] /= P_sum;
+	}
 
 
 }
@@ -221,6 +223,7 @@ pub fn run(cfg: &Config) -> Result<(), Box<Error>>{
 	}
 	
 	// final free energy calculation and state dump
+	println!("Finished. Dumping final PMF");
 	free_energy(&histograms, &mut P, &mut A);
 	dump_state(&histograms, &F, &F_prev, &P, &A);
 
@@ -228,6 +231,8 @@ pub fn run(cfg: &Config) -> Result<(), Box<Error>>{
 		println!("!!!!! WHAM not converged! (max iterations reached) !!!!!");
 	}
 
+	io::write_results(&cfg.output, &histograms, &A, &P)?;
+	
 	Ok(())
 }
 

src/main.rs

@@ -21,9 +21,10 @@ fn cli() -> Result<Config, Box<Error>> {
 	let tolerance: f32 = matches.value_of("tolerance").unwrap_or("0.000001").parse()?;
 	let max_iterations: usize = matches.value_of("iterations").unwrap_or("100000").parse()?;
 	let cyclic: bool = matches.is_present("cyclic");
+	let output = matches.value_of("output").unwrap_or("wham.out").to_string();
 
 	Ok(wham::Config{metadata_file, hist_min, hist_max, num_bins,
-		verbose, tolerance, max_iterations, temperature, cyclic})
+		verbose, tolerance, max_iterations, temperature, cyclic, output})
 }
 
 fn main() {