Source code for dysmalpy.galaxy

# coding=utf8
# Copyright (c) MPE/IR-Submm Group. See LICENSE.rst for license information. 
#
# Main classes and functions for DYSMALPY for simulating the kinematics of
# a model galaxy and fitting it to observed data.


from __future__ import (absolute_import, division, print_function,
                        unicode_literals)

# Standard library
import logging
import copy
from collections import OrderedDict
import os

# Third party imports
import astropy.cosmology as apy_cosmo

# Package imports
from dysmalpy.observation import Observation
from dysmalpy.models import ModelSet
from dysmalpy.utils_io import write_model_obs_file
from dysmalpy.data_io import load_pickle, dump_pickle


__all__ = ['Galaxy']


# LOGGER SETTINGS
logging.basicConfig(level=logging.INFO)
logger = logging.getLogger('DysmalPy')
logger.setLevel(logging.INFO)


# Default cosmology
_default_cosmo = apy_cosmo.FlatLambdaCDM(H0=70., Om0=0.3)




[docs] class Galaxy: r""" Container for simulating or modelling a galaxy `Galaxy` holds the observed data, model, observing instrument, and general information for a galaxy. This can be a simulated or real galaxy. Parameters ---------- z : float Redshift of the galaxy cosmo : `~astropy.cosmology` object The cosmology to use for modelling. Default is astropy.cosmology.FlatLambdaCDM with H0=70., and Om0=0.3. model : `~dysmalpy.models.ModelSet` object A dysmalpy model to use for simulating and/or fitting data. This generates the intrinsic observables of the galaxy based on the components included in the ModelSet. obs_list : list List of `~dysmalpy.observation.Observation` objects, which hold `~dysmalpy.instrument.Instrument`, `~dysmalpy.observation.ObsModOptions`, and `~dysmalpy.data_classes.Data` instances. For each `obs`, `obs.instrument` and `obs.mod_options` describes how `model` is converted to observed space data. name : str, optional Name of the galaxy. Default is "galaxy." """ def __init__(self, z=0, cosmo=_default_cosmo, obs_list=None, model=None, name='galaxy'): self._z = z self.name = name if model is None: self.model = ModelSet() else: self.model = model self.cosmo = cosmo # v2.0: everything lives inside Observation instances, inside # OrderedDict self.observations self.observations = OrderedDict() if obs_list is not None: for obs in obs_list: self.add_observation(obs) def __setstate__(self, state): # Compatibility hack, to handle the changed galaxy structure # (properties, not attributes for data[*], instrument) self.__dict__ = state # -------------------------------- # v2 Migration: Move everything under observations, # and get rid of separate 1d/2d/3d datas: # Oldest migration: if '_data' in state.keys(): pass else: migrate_keys = ['data', 'data1d', 'data2d', 'data3d', 'instrument', 'dscale'] for mkey in migrate_keys: if (mkey in state.keys()) and ('_{}'.format(mkey) not in state.keys()): self.__dict__['_{}'.format(mkey)] = state[mkey] del self.__dict__[mkey] # Now migrate to the new observation framework: if 'observations' in state.keys(): pass else: inst = self.__dict__['_instrument'] data=self.__dict__['_data'] inst.ndim = data.ndim key_mig_dat2inst = ['apertures', 'slit_width', 'slit_pa', 'smoothing_type', 'smoothing_npix'] for key in key_mig_dat2inst: if key in data.__dict__.keys(): inst.__dict__[key] = data.__dict__[key] if 'apertures' in data.__dict__.keys(): if data.apertures is not None: inst.moment = data.apertures.apertures[0].moment else: inst.apertures = None if 'moment' in data.__dict__.keys(): inst.moment = data.moment else: inst.moment = False else: inst.apertures = None if 'moment' in data.__dict__.keys(): inst.moment = data.moment else: inst.moment = False if 'line_center' not in inst.__dict__.keys(): inst.line_center = None dat_del_keys = ['apertures', 'aper_center_pix_shift', 'slit_width', 'slit_pa', 'smoothing_type', 'smoothing_npix'] for dkey in dat_del_keys: if dkey in data.__dict__.keys(): del data.__dict__[dkey] obs = Observation('OBS', 'LINE', instrument=inst, data=data) obs.model_cube = self.__dict__['model_cube'] obs.model_data = self.__dict__['model_data'] self.observations = OrderedDict() self.add_observation(obs) # Delete old-format keys: delete_keys = ['_data', '_data1d', '_data2d', '_data3d', '_instrument', 'model_cube', 'model_data', 'filename_velocity', 'filename_dispersion'] for dkey in delete_keys: del self.__dict__[dkey] # Astropy migration: need to make new cosmo instance: if 'name' in self._cosmo.__dict__.keys(): kw_cosmo_new = {} keys_params = ['H0', 'Om0', 'Tcmb0', 'Neff', 'm_nu', 'Ob0', 'name', 'meta'] for key in keys_params: if '_{}'.format(key) in self._cosmo.__dict__.keys(): kw_cosmo_new[key] = self._cosmo.__dict__['_{}'.format(key)] cosmo_new = apy_cosmo.FlatLambdaCDM(**kw_cosmo_new) self.cosmo = cosmo_new # --------------------------------
[docs] def copy(self): return copy.deepcopy(self)
@property def z(self): return self._z @z.setter def z(self, value): if value < 0: raise ValueError("Redshift can't be negative!") self._z = value # Reset dscale: self._set_dscale() @property def cosmo(self): return self._cosmo @cosmo.setter def cosmo(self, new_cosmo): if not isinstance(new_cosmo, apy_cosmo.FLRW): raise TypeError("Cosmology must be an astropy.cosmology.FLRW " "instance.") if new_cosmo is None: self._cosmo = _default_cosmo self._cosmo = new_cosmo # Reset dscale: self._set_dscale() @property def dscale(self): return self._dscale def _set_dscale(self): self._dscale = self._cosmo.arcsec_per_kpc_proper(self._z).value
[docs] def add_observation(self, obs): """ Add an observation to the galaxy.observations ordered dict. """ obs_name = obs.name if obs_name in self.observations: logger.warning('Overwriting observation {}!'.format(obs_name)) self.observations[obs_name] = obs
[docs] def get_observation(self, obs_name): """ Retrieve an observation from the galaxy.observations ordered dict. """ try: return self.observations[obs_name] except KeyError: raise KeyError('{} not in self.observations !'.format(obs_name))
[docs] def create_model_data(self, obs_list=None): r""" Function to simulate data for the galaxy The function will initially generate a data cube that will then be optionally reduced to 2D, 1D, or single spectrum data if specified. The generated cube can be accessed via `Galaxy.model_cube` and the generated final data products via `Galaxy.model_data`. Both of these attributes are `data_classes.Data` instances. Calls the per-observation `Observation.create_single_obs_model_data` method Parameters ---------- obs_list : list, optional List of observations to make models for. If omitted, will default to making models for all observations in the galaxy. """ if obs_list is None: obs_list = self.observations.keys() for obs_name in obs_list: # Get the individual observation obs = self.observations[obs_name] obs.create_single_obs_model_data(self.model, self.dscale)
#####################
[docs] def preserve_self(self, filename=None, save_data=True, overwrite=False): # Check for existing file: if (not overwrite) and (filename is not None): if os.path.isfile(filename): logger.warning("overwrite={} & File already exists! Will not save file. \n {}".format(overwrite, filename)) return None if filename is not None: galtmp = copy.deepcopy(self) if not save_data: for obs_name in galtmp.observations: galtmp.observations[obs_name].data = None galtmp.observations[obs_name].model_data = None galtmp.observations[obs_name].model_cube = None # _pickle.dump(galtmp, open(filename, "wb") ) dump_pickle(galtmp, filename=filename, overwrite=overwrite)
[docs] def load_self(self, filename=None): """ Load a saved Galaxy from a pickle file Parameters ---------- filename : str Name of the file with saved Galaxy Returns ------- """ if filename is not None: #galtmp = _pickle.load(open(filename, "rb")) # galtmp = _rename_unpickler(open(filename, "rb")).load() galtmp = load_pickle(filename) return galtmp
[docs] def save_model_data(self, filenames=None, overwrite=False): # Check for existing files: if (not overwrite) and (filenames is not None): for obs_name in self.observations: filename = filenames[obs_name] if os.path.isfile(filename): logger.warning("overwrite={} & File already exists! Will not save file. \n {}".format(overwrite, filename)) return None if filenames is not None: for obs_name in self.observations: obs = self.observations[obs_name] filename = filenames[obs_name] write_model_obs_file(obs=obs, model=self.model, fname=filename, overwrite=overwrite)
def load_galaxy_object(filename=None): """ Load a saved Galaxy from a pickle file Parameters ---------- filename : str Name of the file with saved Galaxy Returns ------- gal: Galaxy object The saved dysmalpy Galaxy object """ gal = Galaxy() gal = gal.load_self(filename=filename) return gal