Change distribution outputs files

Author: ebachelet

pyLIMA/fits/MCMC_fit.py

@@ -176,12 +176,17 @@ def fit(self, initial_population=[], computational_pool=False):
 
     def reconstruct_chains(self, mcmc_samples, mcmc_prob):
 
+
         rangei, rangej, rangek = mcmc_samples.shape
 
-        #MCMC_chains = np.zeros((rangei, rangej, rangek + 2))
-        #MCMC_chains[:, :, :-2] = mcmc_samples
-        #MCMC_chains[:, :, -2] = mcmc_prob
-        #MCMC_chains[:, :, -1] = np.zeros(mcmc_prob.shape)
+        MCMC_chains = np.zeros((rangei, rangej, rangek + 2))
+        MCMC_chains[:, :, :-2] = mcmc_samples
+        MCMC_chains[:, :, -2] = mcmc_prob
+        MCMC_chains[:, :, -1] = np.zeros(mcmc_prob.shape)
+
+        if self.telescopes_fluxes_method=='fit':
+
+            return MCMC_chains,MCMC_chains
 
         Rangei,Rangej = self.trials_parameters[:,:-2].shape
         #MCMC_chains_with_fluxes = np.zeros((rangei,rangej,Rangej+2))
@@ -204,9 +209,8 @@ def reconstruct_chains(self, mcmc_samples, mcmc_prob):
         MCMC_chains_with_fluxes = np.array(MCMC_FLUXES)[MCMC_rebuild].reshape(rangei,
                                                                               rangej,
                                                                               Rangej + 2)
-        MCMC_chains = mcmc_samples
+        MCMC_chains[:,:,-1] = MCMC_chains_with_fluxes[:,:,-1]
 
-        return MCMC_chains, MCMC_chains_with_fluxes
 
         #MCMC_FLUXES = np.zeros((MCMC_unique.shape[0],trials_unique.shape[1]))
         #match = []
@@ -259,26 +263,24 @@ def reconstruct_chains(self, mcmc_samples, mcmc_prob):
         #        MCMC_chains_with_fluxes[:,j][:,-1] = np.array(unique_priors)[
         #        unique_sample[1].ravel()]
 
-        #columns_to_swap = []
-        #if self.rescale_photometry:
-        #    columns_to_swap += self.rescale_photometry_parameters_index
+        columns_to_swap = []
+        if self.rescale_photometry:
+            columns_to_swap += self.rescale_photometry_parameters_index
 
-        #if self.rescale_astrometry:
-        #    columns_to_swap += self.rescale_photometry_parameters_index
+        if self.rescale_astrometry:
+            columns_to_swap += self.rescale_photometry_parameters_index
 
-        #if (columns_to_swap != []):
+        if (columns_to_swap != []):
 
-        #    old_column = columns_to_swap
-        #    new_column = np.arange(old_column[-1]+1,Rangej-1,1).tolist()
+            old_column = columns_to_swap
+            new_column = np.arange(old_column[-1]+1,Rangej-1,1).tolist()
 
-         #   MCMC_chains_with_fluxes[:, :, old_column + new_column] =
-         #       #   MCMC_chains_with_fluxes[:, :,
-         #                                                        new_column +
-        #                                                        old_column]
+            MCMC_chains_with_fluxes[:, :, old_column + new_column] = MCMC_chains_with_fluxes[:, :,new_column +old_column]
 
         #MCMC_chains[:,:,-1] = np.copy(MCMC_chains_with_fluxes[:,:,-1])
 
 
+        return MCMC_chains, MCMC_chains_with_fluxes
 
     def samples_to_plot(self):
 

pyLIMA/outputs/file_outputs.py

@@ -13,6 +13,8 @@ def json_output(fit_object,json_name='./pyLIMA_fit.json'):
     dictionnary['FitPriors'] = {'Parameters':fit_object.priors_parameters}
     dictionnary['FitResults'] = {}
 
+    distributions_to_output = {}
+
     for key, value in fit_object.fit_results.items():
 
         if isinstance(value, np.ndarray):
@@ -30,17 +32,46 @@ def json_output(fit_object,json_name='./pyLIMA_fit.json'):
 
             try:
 
-                dictionnary['FitResults'][key] = tosave
+                if (key == 'MCMC_chains'):
+
+                    header =  [key for key in fit_object.fit_parameters.keys()]
+                    distributions_to_output[key] = [tosave,header]
+
+                elif (key == 'MCMC_chains_with_fluxes'):
+
+                    header =  [key for key in fit_object.priors_parameters.keys()] +[
+                        'priors',
+                                                                 fit_object.loss_function ]
+                    distributions_to_output[key] = [tosave,header]
+
+                elif (key == 'DE_population'):
+
+                    header =  [key for key in fit_object.priors_parameters.keys()]
+
+                    distributions_to_output[key] = [tosave,header]
+
+                else:
+
+                    dictionnary['FitResults'][key] = tosave
 
             except TypeError:
 
                 pass
 
+
     with open(json_name, 'w') as outfile:
         json.dump(dictionnary,outfile,indent=4)
 
 
+    for key in distributions_to_output.keys():
+
+        dist_name = json_name[:-5]+'_'+key+'.npz'
+        np.savez(dist_name,pyLIMA_distribution=distributions_to_output[key][0],
+                 header=distributions_to_output[key][1]+[
+                     fit_object.loss_function, 'priors'])
 
+   # with open(dist_name, 'w') as outfile:
+   #     json.dump(distributions_to_output,outfile,indent=4)
 
 
 # def json_output(array_parameters, parameters_name, filename='parameters',