fixed typo in analyze.py

data/bao_2014.txt

@@ -22,10 +22,10 @@ LOWZ    0.32    8.47    0.17    3
 CMASS   0.57    13.77   0.13    3
 #
 # BOSS LyaQSO (DA/rs) Font-Ribera et al. 1311.1767
-# LyaQSO  2.36   10.8     0.4     5
+ LyaQSO  2.36   10.8     0.4     5
 #
 # BOSS LyaQSO c/(H rs) Font-Ribera et al. 1311.1767
-# LyaQSO	2.36    9.0     0.3     6
+ LyaQSO	2.36    9.0     0.3     6
 #
 # SDSS DR7 MGS, Ross et al. 1409.3242v1
 MGS     0.15   4.47     0.16    3

montepython/analyze.py

@@ -832,6 +832,23 @@ def compute_posterior(information_instances):
                                 "'%s-%s' 2d-plot" % (
                                     info.plotted_parameters[info.native_index],
                                     info.plotted_parameters[info.native_second_index]))
+                        except ValueError as e:
+                            if str(e) == "Contour levels must be increasing":
+                                warnings.warn(
+                                    "The routine could not find the contour of the " +
+                                    "'%s-%s' 2d-plot. \n " % (
+                                        info.plotted_parameters[info.native_index],
+                                        info.plotted_parameters[info.native_second_index]) +
+                                    'The error is: "Contour levels must be increasing"' +
+                                    " but " + str(ctr_level(info.n, info.levels[:2])) +
+                                    " were found. This may happen when most" +
+                                    " points fall in the same bin.")
+                            else:
+                                warnings.warn(
+                                    "The routine could not find the contour of the " +
+                                    "'%s-%s' 2d-plot" % (
+                                        info.plotted_parameters[info.native_index],
+                                        info.plotted_parameters[info.native_second_index]))
 
                         ax2dsub.set_xticks(info.ticks[info.native_second_index])
                         ax2dsub.set_yticks(info.ticks[info.native_index])
@@ -1382,13 +1399,18 @@ def recover_folder_and_files(files):
     # The following list defines the substring that a chain should contain for
     # the code to recognise it as a proper chain.
     substrings = ['.txt', '__']
+    # The following variable defines the substring that identify error_log
+    # files and therefore there must not be taken into account in the analysis.
+    substring_err = 'error_log'
+
     limit = 10
     # If the first element is a folder, grab all chain files inside
     if os.path.isdir(files[0]):
         folder = os.path.normpath(files[0])
         files = [os.path.join(folder, elem) for elem in os.listdir(folder)
                  if not os.path.isdir(os.path.join(folder, elem))
                  and not os.path.getsize(os.path.join(folder, elem)) < limit
+                 and (substring_err not in elem)
                  and all([x in elem for x in substrings])]
     # Otherwise, extract the folder from the chain file-name.
     else:
@@ -1402,6 +1424,7 @@ def recover_folder_and_files(files):
                  if os.path.join(folder, elem) in np.copy(files)
                  and not os.path.isdir(os.path.join(folder, elem))
                  and not os.path.getsize(os.path.join(folder, elem)) < limit
+                 and (substring_err not in elem)
                  and all([x in elem for x in substrings])]
     basename = os.path.basename(folder)
     return folder, files, basename

montepython/cosmo_hammer.py

@@ -22,9 +22,11 @@
 
 import io_mp
 import sampler
-from cosmoHammer.likelihood.chain.LikelihoodComputationChain import (
-    LikelihoodComputationChain)
-from cosmoHammer.sampler.CosmoHammerSampler import CosmoHammerSampler
+from cosmoHammer import MpiCosmoHammerSampler
+from cosmoHammer import LikelihoodComputationChain
+#from cosmoHammer.likelihood.chain.LikelihoodComputationChain import (
+#    LikelihoodComputationChain)
+#from cosmoHammer.sampler.CosmoHammerSampler import CosmoHammerSampler
 from cosmoHammer.util.SampleFileUtil import SampleFileUtil
 
 # Cosmo Hammer subfolder and name separator
@@ -33,6 +35,7 @@
 # Cosmo Hammer file names ending, after the defined 'base_name'
 name_arguments = '.arguments'
 name_chain = 'chain_CH__sampling.txt'
+name_error_chain =  'chain_CH__sampling-error_log.txt'
 
 # Cosmo Hammer option prefix
 CH_prefix = 'CH_'
@@ -48,6 +51,10 @@
     {'help': 'Number of sample iterations',
      'type': int}}
 
+# List of strings that will contain the derived parameters names.
+# Filled in run using the MontePython class data
+# Used in DerivedUtil.persistValues
+CH_derived_names = []
 
 def run(cosmo, data, command_line):
     """
@@ -57,6 +64,9 @@ def run(cosmo, data, command_line):
     # Store the parameters inside the format expected by CosmoHammer
     # TODO: about the derived params?
     parameter_names = data.get_mcmc_parameters(["varying"])
+    # Fill now CH_derived_names to be used in DerivedUtil.persistValues
+    global CH_derived_names
+    CH_derived_names= data.get_mcmc_parameters(["derived"])
 
     # Ensure that their prior is bound and flat
     is_flat, is_bound = sampler.check_flat_bound_priors(
@@ -91,6 +101,7 @@ def run(cosmo, data, command_line):
     # here, since data must be called before cosmo.
     chain.addCoreModule(data)
     chain.addCoreModule(cosmo)
+    chain.addCoreModule(store_cosmo_derived)
 
     # Add each likelihood class as a LikelihoodModule
     for likelihood in data.lkl.itervalues():
@@ -101,7 +112,9 @@ def run(cosmo, data, command_line):
     file_prefix = os.path.join(command_line.folder, CH_subfolder, chain_name)
 
     # Recover the User options
-    data.CH_arguments = {}
+    data.CH_arguments = {'walkersRatio': 50,
+                         'burninIterations': 10,
+                         'sampleIterations': 30}
     for arg in CH_user_arguments:
         value = getattr(command_line, CH_prefix+arg)
         if value != -1:
@@ -127,13 +140,10 @@ def run(cosmo, data, command_line):
         num_threads = 1
 
     # Create the Sampler object
-    sampler_hammer = CosmoHammerSampler(
+    sampler_hammer = MpiCosmoHammerSampler(
         params=params,
         likelihoodComputationChain=chain,
         filePrefix=file_prefix,
-        walkersRatio=50,
-        burninIterations=10,
-        sampleIterations=30,
         storageUtil=derived_util,
         threadCount=num_threads,
         **data.CH_arguments)
@@ -158,6 +168,8 @@ def from_CH_output_to_chains(folder):
 
     This function will be called by the module :mod:`analyze`.
     """
+    chainsClean, lklClean = [], []  # Arrays to store the valid steps.
+    error = []  # Array to store the step did not could not be computed.
 
     chain_name = [a for a in folder.split(os.path.sep) if a][-2]
     base_name = os.path.join(folder, chain_name)
@@ -172,17 +184,25 @@ def from_CH_output_to_chains(folder):
     # multiplicity. This does not mean that the acceptance rate is one, but
     # points where the code stayed are duplicated in the file.
 
-    ## First, reshape the lkl array
-    lkl = np.array([[elem] for elem in lkl])
+    # First, separe errors from accepted points:
+    for x, p in zip(chains, lkl):
+        if p != -np.inf:
+            chainsClean.append(x)
+            lklClean.append([-p])  # MP stores -loglkl in its files.
+        else:
+            error.append([val for val in x if not np.isnan(val)])  # Remove nan's
 
     ## Create the array of ones
-    ones = np.array([[1] for _ in range(len(lkl))])
+    ones = np.array([[1] for _ in range(len(lklClean))])
+
 
     ## Concatenate everything and save to file
-    final = np.concatenate((ones, lkl, chains), axis=1)
+    final = np.concatenate((ones, lklClean, chainsClean), axis=1)
     output_folder = os.path.join(folder, '..')
     output_chain_path = os.path.join(output_folder, name_chain)
+    output_error_chain_path = os.path.join(output_folder, name_error_chain)
     np.savetxt(output_chain_path, final)
+    np.savetxt(output_error_chain_path, error)
 
 
 class DerivedUtil(SampleFileUtil):
@@ -195,8 +215,12 @@ def persistValues(self, posFile, probFile, pos, prob, data):
         """
         # extend the pos array to also contain the value of the derived
         # parameters
+
+        derived_not_computed = [np.NaN] * len(CH_derived_names)
+
         derived = np.array(
-            [[a for a in elem.itervalues()] for elem in data])
+            #[[a for a in elem.itervalues()] if elem else derived_not_computed for elem in data])
+            [[elem[key] for key in sorted(elem.iterkeys())] if elem else derived_not_computed for elem in data])
         final = np.concatenate((pos, derived), axis=1)
 
         posFile.write("\n".join(
@@ -207,3 +231,31 @@ def persistValues(self, posFile, probFile, pos, prob, data):
         probFile.write("\n".join([str(p) for p in prob]))
         probFile.write("\n")
         probFile.flush()
+
+
+def store_cosmo_derived(ctx):
+    """
+    Store derived parameters. Based on sampler.compute_lkl
+    """
+    from classy import CosmoSevereError
+
+    cosmo = ctx.get("cosmo")
+    data = ctx.get("data")
+
+    if data.get_mcmc_parameters(['derived']) != []:
+        try:
+            derived = cosmo.get_current_derived_parameters(
+                data.get_mcmc_parameters(['derived']))
+        except AttributeError:
+            # This happens if the classy wrapper is still using the old
+            # convention, expecting data as the input parameter
+            cosmo.get_current_derived_parameters(data)
+        except CosmoSevereError:
+            raise io_mp.CosmologicalModuleError(
+                "Could not write the current derived parameters")
+
+        ctx.add('key_data', derived)
+
+        # If we want the derived parameters stored rescaled
+        # derived_scaled = {elem: value/data.mcmc_parameters[elem]['scale'] for elem, value in derived.iteritems()}
+        # ctx.add('key_data', derived_scaled)

montepython/data.py

@@ -189,6 +189,11 @@ def __init__(self, command_line, path):
         self.out = None
         self.out_name = ''
 
+        # Create the variable err, and err_name to keep the errors
+        # will be initialised later by the :mod:`io_mp` module
+        self.err = None
+        self.err_name = ''
+
         # If the parameter file is not a log.param, the path will be read
         # before reading the parameter file.
         if self.param.find('log.param') == -1:
@@ -995,8 +1000,13 @@ def __init__(self, array, key):
         dict.__init__(self)
 
         self['initial'] = array[0:4]
-        self['scale'] = array[4]
-        self['role'] = array[-1]
+        if array[4]:
+            self['scale'] = array[4]
+        else:
+            raise io_mp.ParameterError("You chose scale 0 to param {}. Correct\
+                                       it.".format(key))
+
+        self['role'] = array[5]
         self['tex_name'] = io_mp.get_tex_name(key)
         if array[3] == 0:
             self['status'] = 'fixed'

montepython/initialise.py

@@ -84,6 +84,7 @@ def initialise(custom_command=''):
             io_mp.create_output_files(command_line, data)
         # NS: Creating the NS subfolder and the MultiNest arguments
         elif command_line.method == 'NS':
+            io_mp.create_output_files(command_line, data, only_error=True)
             from nested_sampling import initialise as initialise_ns
             initialise_ns(cosmo, data, command_line)