diff --git a/binarycpython/utils/grid.py b/binarycpython/utils/grid.py
index 577b7ba5c24f6c252634cc5bda55390efb44fa68..4c9c1232e7a308b9193b00d6c46b782c2f1b880e 100644
--- a/binarycpython/utils/grid.py
+++ b/binarycpython/utils/grid.py
@@ -67,6 +67,7 @@ from binarycpython.utils.functions import (
     recursive_change_key_to_float,
     custom_sort_dict,
     recursive_change_key_to_string,
+    multiply_values_dict,
 )
 
 
@@ -1017,27 +1018,34 @@ class Population:
                 break
 
         # Extra ensemble result manipulation:
-        # We need to reformat and multiply by a factor
-
         combined_output_dict['ensemble_results']['ensemble'] = self.format_ensemble_results(combined_output_dict['ensemble_results']['ensemble'])
         gc.collect()
 
-        # Multiply the value by a factor:
-        # if self.grid_options['ensemble_factor_in_probability_weighted_mass']:
-            
-
+        # Take into account that we run this on multiple cores
+        combined_output_dict["_total_probability_weighted_mass_run"] = combined_output_dict["_total_probability_weighted_mass_run"]
 
         # Put the values back as object properties
         self.grid_results = combined_output_dict["results"]
+
+        # Put Ensemble results
         self.grid_ensemble_results = combined_output_dict[
             "ensemble_results"
         ]  # Ensemble results are also passed as output from that dictionary
 
         # Add some metadata
-        self.grid_ensemble_results["population_id"] = self.grid_options[
+        self.grid_ensemble_results["metadata"] = {}
+        self.grid_ensemble_results["metadata"]["population_id"] = self.grid_options[
             "_population_id"
         ]
+        self.grid_ensemble_results["metadata"][
+            "total_probability_weighted_mass"
+        ] = combined_output_dict["_total_probability_weighted_mass_run"]
+        self.grid_ensemble_results["metadata"]["factored_in_probability_weighted_mass"] = False
+        if self.grid_options['ensemble_factor_in_probability_weighted_mass']:
+            multiply_values_dict(self.grid_ensemble_results['ensemble'], 1/self.grid_ensemble_results["metadata"]['total_probability_weighted_mass'])
+            self.grid_ensemble_results["metadata"]["factored_in_probability_weighted_mass"] = True
 
+        # 
         self.grid_options["_failed_count"] = combined_output_dict["_failed_count"]
         self.grid_options["_failed_prob"] = combined_output_dict["_failed_prob"]
         self.grid_options["_failed_systems_error_codes"] = list(
@@ -1359,26 +1367,8 @@ class Population:
 
                 # Extra ensemble result manipulation:
                 # We need to reformat and multiply by a factor
-                original_ensemble_results = ensemble_json["ensemble"]
-
-                float_format_ensemble_results = recursive_change_key_to_float(original_ensemble_results)
-                del original_ensemble_results
-                gc.collect()
-
-                # Then sort the dictionary
-                sorted_ensemble_results = custom_sort_dict(float_format_ensemble_results)
-                del float_format_ensemble_results
+                ensemble_json["ensemble"] = self.format_ensemble_results(ensemble_json["ensemble"])
                 gc.collect()
-
-                # Then Change the keys back to a string but with a %g format.
-                reformatted_ensemble_results = recursive_change_key_to_string(sorted_ensemble_results)
-                del sorted_ensemble_results
-                gc.collect()
-
-                # Put back in the dictionary
-                ensemble_json["ensemble"] = reformatted_ensemble_results
-                gc.collect()
-
             else:
                 # If we do not allow this, automatically we will export this to the data_dir, in
                 # some formatted way
@@ -1779,6 +1769,7 @@ class Population:
         # Import packages
         code_string += "import math\n"
         code_string += "import numpy as np\n"
+        code_string += "from collections import OrderedDict\n"
         code_string += "from binarycpython.utils.distribution_functions import *\n"
         code_string += "from binarycpython.utils.spacing_functions import *\n"
         code_string += "from binarycpython.utils.useful_funcs import *\n"
@@ -3514,6 +3505,9 @@ class Population:
         ######
         # Setting up the grid variables
 
+        # Make a normal dict of the options again
+        options = dict(options)
+
         # Multiplicity
         self.add_grid_variable(
             name="multiplicity",
@@ -3558,7 +3552,7 @@ class Population:
             dphasevol="dlnm1",
             precode='M_1 = np.exp(lnm1); options["M1"]=M_1',
             probdist="Moe_de_Stefano_2017_pdf({{{}, {}, {}}}, verbosity=self.grid_options['verbosity'])['total_probdens'] if multiplicity == 1 else 1".format(
-                str(options)[1:-1], "'multiplicity': multiplicity", "'M1': M_1"
+                str(dict(options))[1:-1], "'multiplicity': multiplicity", "'M1': M_1"
             ),
         )
 
@@ -3827,7 +3821,7 @@ eccentricity3=0
         # and add some keys to it:
 
         probdist_addition = "Moe_de_Stefano_2017_pdf({{{}, {}, {}, {}, {}, {}, {}, {}, {}, {}, {}, {}}}, verbosity=self.grid_options['verbosity'])['total_probdens']".format(
-            str(options)[1:-1],
+            str(dict(options))[1:-1],
             '"multiplicity": multiplicity',
             '"M1": M_1',
             '"M2": M_2',
diff --git a/binarycpython/utils/grid_options_defaults.py b/binarycpython/utils/grid_options_defaults.py
index fee8212e03fdb2feb8d25173dce277a36aa848e8..c994fae9dc4b3290103b0503bf9d72c92a02fbf5 100644
--- a/binarycpython/utils/grid_options_defaults.py
+++ b/binarycpython/utils/grid_options_defaults.py
@@ -44,6 +44,7 @@ grid_options_defaults_dict = {
     "_set_ms_grid": False, # Whether the M&S grid has been loaded
     "run_zero_probability_system": True, # Whether to run the zero probability systems
     "_zero_prob_stars_skipped": 0,
+    "ensemble_factor_in_probability_weighted_mass": False, # Whether to multiply the ensemble results by 1/probability_weighted_mass
     ##########################
     # Execution log:
     ##########################
@@ -484,7 +485,8 @@ grid_options_descriptions = {
     "max_queue_size": "Maximum size of the queue that is used to feed the processes. Don't make this too big! Default: 1000. Input: int",
     "_set_ms_grid": "Internal flag whether the M&S grid has been loaded",
     "run_zero_probability_system": "Whether to run the zero probability systems. Default: True. Input: boolean",
-    "_zero_prob_stars_skipped": "Internal counter to track how many systems are skipped because they have 0 probability"
+    "_zero_prob_stars_skipped": "Internal counter to track how many systems are skipped because they have 0 probability",
+    "ensemble_factor_in_probability_weighted_mass": "Flag to multiply all the ensemble results with 1/probability_weighted_mass",
 }
 
 ###