ra dec qso

Signed-off-by: Anto Idicherian Lonappan <[email protected]>

Author: antolonappan

marcia/database.py

@@ -193,18 +193,25 @@ def get_SNE(self):
     @load_data_once
     def get_QSO_data(self):
         data = loadtxt(os.path.join(__datapath__, 'Quasars','table3.dat'), usecols=(3,4,5,6,7,9,10, 11, 12))
+        position = loadtxt(os.path.join(__datapath__, 'Quasars','table3.dat'), usecols=(1,2))
+        ra = position[:,0]
+        dec = position[:,1]
         z, lnFUV, lnFUV_err, lnFX, lnFX_err = data[:,0], data[:,1], data[:,2], data[:,3], data[:,4]
         DM, dDM = data[:,7], data[:,8]
-        return z, lnFUV, lnFUV_err, lnFX, lnFX_err, DM, dDM
+        return z, lnFUV, lnFUV_err, lnFX, lnFX_err, DM, dDM, ra, dec
 
     def get_QSO_full(self):
-        z, lnFUV, lnFUV_err, lnFX, lnFX_err, _, _ = self.get_QSO_data()
+        z, lnFUV, lnFUV_err, lnFX, lnFX_err, _, _ ,_,_= self.get_QSO_data()
         return z, lnFUV, lnFUV_err, lnFX, lnFX_err
 
     def get_QSO(self):
-        x,_,_,_,_,y,sigma = self.get_QSO_data()
+        x,_,_,_,_,y,sigma,_,_ = self.get_QSO_data()
         covar = np.diag(sigma**2)
         return x,y,covar
+    
+    def get_QSO_position(self):
+        _,_,_,_,_,_,_,ra,dec = self.get_QSO_data()
+        return ra,dec
 
     def get_CMB_planck_TT(self):
         # This is the Planck 2018 data special case where we want 

marcia/sampler.py

@@ -89,17 +89,23 @@ def get_burnin(self):
             print(f'Burn-in: {burnin} and thin: {thin}')
         except:
             print('Autocorrelation time could not be calculated, increase the number of iterations')
-            burnin = 0
-            thin = 1
+            burnin = 100
+            thin = 10
             print(f'Burn-in: {burnin} and thin: {thin}[DEFAULT VALUES]')
         return burnin, thin
 
-    def get_chain(self, getdist=False):
+    def get_chain(self, getdist=False,burnin=None,thin=None):
         sampler = self.HDFBackend
         if sampler.iteration < self.max_n:
             print(f'Only {sampler.iteration} iterations completed')
             print(f'You should run the sampler to finsih the sampling of {self.max_n} iterations')
-        burnin, thin = self.get_burnin()
+        if (burnin is None) and (thin is None):
+            burnin, thin = self.get_burnin()
+        elif burnin is None:
+            burnin,_ = self.get_burnin()
+        elif thin is None:
+            _,thin = self.get_burning()
+
         samples = sampler.get_chain(discard=burnin, thin=thin, flat=True)
         if getdist:
             lnprob = sampler.get_log_prob(discard=burnin, thin=thin, flat=True)
@@ -109,8 +115,8 @@ def get_chain(self, getdist=False):
             samples = np.concatenate((lnprob[:, None], lnprior[:, None], samples), axis=1)
         return samples
 
-    def corner_plot(self, getdist=False):
-        chains = self.get_chain()
+    def corner_plot(self, getdist=False,burnin=None,thin=None):
+        chains = self.get_chain(getdist,burnin,thin)
         names = self.likelihood.theory.param.parameters
         labels = [p.replace('$', '') for p in self.likelihood.theory.labels]
         if getdist: