inc_phot and inc_spec

species/analysis/fit_model.py

@@ -122,7 +122,7 @@ def lnlike(param,
             chisq += (item[0]-flux)**2 / item[1]**2
 
     if spectrum is not None:
-        model = modelspec.get_model(paramdict, ('original', (0.9, 2.5)))
+        model = modelspec.get_model(paramdict, None)
 
         flux_new = spectres.spectres(new_spec_wavs=spectrum[:, 0],
                                      old_spec_wavs=model.wavelength,
@@ -210,7 +210,8 @@ def __init__(self,
                  filters,
                  model,
                  bounds,
-                 data):
+                 inc_phot=True,
+                 inc_spec=True):
         """
         Parameters
         ----------
@@ -224,9 +225,10 @@ def __init__(self,
         bounds : dict
             Parameter boundaries. Full parameter range is used if None or not specified. The
             radius parameter range is set to 0-5 Rjup if not specified.
-        data : tuple(bool, bool)
-            Data to use for the fit ('photometry' and/or 'spectrum'). The first boolean sets
-            the inclusion of photometry and the second boolean the inclusion of a spectrum.
+        inc_phot : bool
+            Include photometry data with the fit.
+        inc_spec : bool
+            Include spectral data with the fit.
 
         Returns
         -------
@@ -239,7 +241,7 @@ def __init__(self,
         self.model = model
         self.bounds = bounds
 
-        if not data[0] and not data[1]:
+        if not inc_phot and not inc_spec:
             raise ValueError('No photometric or spectral data has been selected.')
 
         if self.bounds is not None and 'teff' in self.bounds:
@@ -262,7 +264,7 @@ def __init__(self,
         if 'radius' not in self.bounds:
             self.bounds['radius'] = (0., 5.)
 
-        if data[0]:
+        if inc_phot:
             self.objphot = []
             self.modelphot = []
             self.synphot = []
@@ -288,10 +290,11 @@ def __init__(self,
             self.modelphot = None
             self.synphot = None
 
-        if data[1]:
+        if inc_spec:
             self.spectrum = self.object.get_spectrum()
             self.instrument = self.object.get_instrument()
             self.modelspec = read_model.ReadModel(self.model, (0.9, 2.5), teff_bound)
+
         else:
             self.spectrum = None
             self.instrument = None
@@ -305,8 +308,7 @@ def run_mcmc(self,
                  nsteps,
                  guess,
                  tag,
-                 prior=None,
-                 ncpu=1):
+                 prior=None):
         """
         Function to run the MCMC sampler.
 
@@ -325,9 +327,6 @@ def run_mcmc(self,
             Gaussian prior on one of the parameters. Currently only possible for the mass, e.g.
             ('mass', 13., 3.) for an expected mass of 13 Mjup with an uncertainty of 3 Mjup. Not
             used if set to None.
-        ncpu : int
-            Number of parallel processes. Due to the additional overhead, a value larger than 1
-            will not necessarily result in a decrease in computation time.
 
         Returns
         -------
@@ -367,8 +366,7 @@ def run_mcmc(self,
                                                prior,
                                                self.spectrum,
                                                self.instrument,
-                                               self.modelspec]),
-                                        threads=ncpu)
+                                               self.modelspec]))
 
         progbar = progress.bar.Bar('\rRunning MCMC...',
                                    max=nsteps,

species/data/database.py

@@ -101,7 +101,7 @@ def add_companion(self,
         """
         Parameters
         ----------
-        tuple(str, )
+        name : tuple(str, )
             Companion name. All companions are added if set to None.
 
         Returns
@@ -702,7 +702,9 @@ def get_mcmc_spectra(self,
 
     def get_object(self,
                    object_name,
-                   filter_id=None):
+                   filter_id=None,
+                   inc_phot=True,
+                   inc_spec=True):
         """
         Parameters
         ----------
@@ -711,47 +713,67 @@ def get_object(self,
         filter_id : tuple(str, )
             Filter IDs for which the photometry is selected. All available photometry of the object
             is selected if set to None.
+        inc_phot : bool
+            Include photometry in the box.
+        inc_spec : bool
+            Include spectrum in the box.
 
         Returns
         -------
         species.core.box.ObjectBox
-            Box with the object.
+            Box with the object's data.
         """
 
+        sys.stdout.write('Getting object: '+object_name+'...')
+        sys.stdout.flush()
+
         h5_file = h5py.File(self.database, 'r')
         dset = h5_file['objects/'+object_name]
 
         distance = np.asarray(dset['distance'])
 
-        magnitude = {}
-        flux = {}
+        if inc_phot:
+
+            magnitude = {}
+            flux = {}
+
+            if filter_id:
+                for item in filter_id:
+                    data = dset[item]
 
-        if filter_id:
-            for item in filter_id:
-                data = dset[item]
+                    magnitude[item] = np.asarray(data[0:2])
+                    flux[item] = np.asarray(data[2:4])
 
-                magnitude[item] = np.asarray(data[0:2])
-                flux[item] = np.asarray(data[2:4])
+            else:
+                for key in dset.keys():
+                    if key not in ('distance', 'spectrum'):
+                        for item in dset[key]:
+                            name = key+'/'+item
+
+                            magnitude[name] = np.asarray(dset[name][0:2])
+                            flux[name] = np.asarray(dset[name][2:4])
+
+            filters = tuple(magnitude.keys())
 
         else:
-            for key in dset.keys():
-                if key not in ('distance', 'spectrum'):
-                    for item in dset[key]:
-                        name = key+'/'+item
 
-                        magnitude[name] = np.asarray(dset[name][0:2])
-                        flux[name] = np.asarray(dset[name][2:4])
+            magnitude = None
+            flux = None
+            filters = None
 
-        if 'objects/'+object_name+'/spectrum' in h5_file:
+        if inc_spec and 'objects/'+object_name+'/spectrum' in h5_file:
             spectrum = np.asarray(h5_file['objects/'+object_name+'/spectrum'])
         else:
             spectrum = None
 
         h5_file.close()
 
+        sys.stdout.write(' [DONE]\n')
+        sys.stdout.flush()
+
         return box.create_box('object',
                               name=object_name,
-                              filter=tuple(magnitude.keys()),
+                              filter=filters,
                               magnitude=magnitude,
                               flux=flux,
                               distance=distance,

species/plot/plot_spectrum.py

@@ -308,16 +308,17 @@ def plot_spectrum(boxes,
                          label=boxitem.name, zorder=5)
 
         elif isinstance(boxitem, box.ObjectBox):
-            for item in boxitem.flux:
-                transmission = read_filter.ReadFilter(item)
-                wavelength = transmission.mean_wavelength()
-                fwhm = transmission.filter_fwhm()
+            if boxitem.flux is not None:
+                for item in boxitem.flux:
+                    transmission = read_filter.ReadFilter(item)
+                    wavelength = transmission.mean_wavelength()
+                    fwhm = transmission.filter_fwhm()
 
-                color_obj_phot = colors[j][0]
+                    color_obj_phot = colors[j][0]
 
-                ax1.errorbar(wavelength, boxitem.flux[item][0]/scaling, xerr=fwhm/2.,
-                             yerr=boxitem.flux[item][1]/scaling, marker='s', ms=5, zorder=6,
-                             color=color_obj_phot, markerfacecolor=color_obj_phot)
+                    ax1.errorbar(wavelength, boxitem.flux[item][0]/scaling, xerr=fwhm/2.,
+                                 yerr=boxitem.flux[item][1]/scaling, marker='s', ms=5, zorder=6,
+                                 color=color_obj_phot, markerfacecolor=color_obj_phot)
 
             if boxitem.spectrum is not None:
                 masked = np.ma.array(boxitem.spectrum, mask=np.isnan(boxitem.spectrum))
@@ -368,7 +369,7 @@ def plot_spectrum(boxes,
             if max_tmp > res_max:
                 res_max = max_tmp
 
-        res_lim = math.ceil(max_tmp)
+        res_lim = math.ceil(res_max)
 
         ax3.axhline(0.0, linestyle='--', color='gray', dashes=(2, 4), zorder=1)
         ax3.set_ylim(-res_lim, res_lim)

species/util/phot_util.py

@@ -2,6 +2,8 @@
 Utility functions for photometry.
 """
 
+import sys
+
 import spectres
 import numpy as np
 
@@ -31,6 +33,9 @@ def multi_photometry(datatype,
         Box with synthetic photometry.
     """
 
+    sys.stdout.write('Calculating synthetic photometry...')
+    sys.stdout.flush()
+
     flux = {}
 
     if datatype == 'model':
@@ -43,6 +48,9 @@ def multi_photometry(datatype,
             readcalib = read_calibration.ReadCalibration(spectrum, item)
             flux[item] = readcalib.get_photometry(model_par)
 
+    sys.stdout.write(' [DONE]\n')
+    sys.stdout.flush()
+
     return box.create_box('synphot', name='synphot', flux=flux)
 
 
@@ -67,7 +75,9 @@ def apparent_to_absolute(app_mag,
 def get_residuals(model,
                   model_par,
                   filters,
-                  objectbox):
+                  objectbox,
+                  inc_phot=True,
+                  inc_spec=False):
     """
     Parameters
     ----------
@@ -79,6 +89,10 @@ def get_residuals(model,
         Filter IDs. All available photometry of the object is used if set to None.
     objectbox : species.core.box.ObjectBox
         Box with the photometry and/or spectrum of an object.
+    inc_phot : bool
+        Include photometry.
+    inc_spec : bool
+        Include spectrum.
 
     Returns
     -------
@@ -89,12 +103,12 @@ def get_residuals(model,
     if filters is None:
         filters = objectbox.filter
 
-    model_phot = multi_photometry(datatype='model',
-                                  spectrum=model,
-                                  filters=filters,
-                                  model_par=model_par)
+    if inc_phot:
+        model_phot = multi_photometry(datatype='model',
+                                      spectrum=model,
+                                      filters=filters,
+                                      model_par=model_par)
 
-    if objectbox.flux is not None:
         res_phot = np.zeros((2, len(objectbox.flux)))
 
         for i, item in enumerate(objectbox.flux):
@@ -106,7 +120,10 @@ def get_residuals(model,
     else:
         res_phot = None
 
-    if objectbox.spectrum is not None:
+    sys.stdout.write('Calculating residuals...')
+    sys.stdout.flush()
+
+    if inc_spec:
         wl_range = (0.9*objectbox.spectrum[0, 0], 1.1*objectbox.spectrum[-1, 0])
 
         readmodel = read_model.ReadModel(model, wl_range)
@@ -127,6 +144,9 @@ def get_residuals(model,
     else:
         res_spec = None
 
+    sys.stdout.write(' [DONE]\n')
+    sys.stdout.flush()
+
     return box.create_box(boxtype='residuals',
                           name=objectbox.name,
                           photometry=res_phot,