diff --git a/fsl_mrs_mce/mc_sim.py b/fsl_mrs_mce/mc_sim.py
new file mode 100644
index 0000000..f7e745a
--- /dev/null
+++ b/fsl_mrs_mce/mc_sim.py
@@ -0,0 +1,45 @@
+# mc_sim.py - Helpers for Monte-Carlo style simulations
+#
+# Author: Konstantin E Bosbach <konstantin.bosbach@mars.uni-freiburg.de>
+
+
+import numpy as np
+import pandas as pd
+
+
+def get_convergence(
+    output_file_paths, molecule, crit_mean=0.01, crit_std=0.10
+):
+    """Function checks if last two files contain converging datasets."""
+
+    if len(output_file_paths) < 2:
+        print("One iteration, therefore no convergence.")
+        return False
+    else:
+        newer_fit_results = pd.read_csv(output_file_paths[-1])
+        older_fit_results = pd.read_csv(output_file_paths[-2])
+
+    newer_mean = newer_fit_results[molecule].mean()
+    older_mean = older_fit_results[molecule].mean()
+
+    # Normalize with mean of latest dataset, to get deviation
+    norm = newer_mean
+    measure_mean = abs(abs(np.abs(newer_mean) - abs(older_mean))/norm)
+
+    newer_std = newer_fit_results[molecule].std()
+    older_std = older_fit_results[molecule].std()
+    measure_std = abs(abs(np.abs(newer_std) - abs(older_std))/norm)
+
+    if measure_mean <= crit_mean:
+        if measure_std <= crit_std:
+            convergence = True
+    else:
+        convergence = False
+
+    print(
+        "Convergence result for ", output_file_paths[-1], " and ",
+        output_file_paths[-2], "\n\t\t\t",
+        measure_mean, measure_std, "convergence: ", str(convergence)
+    )
+
+    return convergence