"""
Datasets are auxiliary classes used to facilitate the reading and inclusion of
observational data in the IMAGINE pipeline"""
# %% IMPORTS
# Package imports
from astropy import units as u
import numpy as np
# IMAGINE imports
from imagine.tools import BaseClass, distribute_matrix, peye, req_attr
# All declaration
__all__ = ['Dataset', 'TabularDataset', 'HEALPixDataset',
'FaradayDepthHEALPixDataset', 'SynchrotronHEALPixDataset',
'DispersionMeasureHEALPixDataset']
# %% CLASS DEFINITIONS
[docs]class Dataset(BaseClass):
"""
Base class for writing helpers to convert arbitrary
observational datasets into IMAGINE's standardized format
"""
def __init__(self):
# Call super constructor
super().__init__()
self.coords = None
self.Nside = None
self.frequency = 'nan'
self.tag = 'nan'
self._cov = None
self._error = None
self._data = None
@property
@req_attr
def name(self):
return(self.NAME)
@property
def frequency(self):
return self._frequency
@frequency.setter
def frequency(self, frequency):
# Converts the frequency to a string in GHz
if hasattr(frequency, 'unit'):
frequency = frequency.to_value(u.GHz, equivalencies=u.spectral())
if frequency.is_integer():
frequency = int(frequency)
self._frequency = str(frequency)
@property
def data(self):
""" Array in the shape (1, N) """
return self._data[np.newaxis,:]
@property
def key(self):
"""Key used in the Observables_dictionary """
return (self.name,self.frequency, self.Nside, self.tag)
@property
def cov(self):
if self._cov is None:
self._cov = self._error * peye(self._data.size)
return self._cov
[docs]class TabularDataset(Dataset):
"""
Base class for tabular datasets, where the data is input in either
in a Python dictionary-like object
(:py:class:`dict`, :py:class:`astropy.table.Table`,
:py:class:`pandas.DataFrame`, etc).
Parameters
----------
data : dict-like
Can be a :py:class:`dict`, :py:class:`astropy.table.Table`,
:py:class:`pandas.DataFrame`, or similar object
containing the data.
name : str
Standard name of this type of observable. E.g. 'fd', 'sync', 'dm'.
data_column : str
Key used to access the relevant dataset from the provided data
(i.e. data[data_column]).
units : astropy.units.Unit or str
Units used for the data.
coordinates_type : str
Type of coordinates used. Can be 'galactic' or 'cartesian'.
lon_column : str
Key used to access the Galactic longitudes (in deg) from
`data`.
lat_column : str
Key used to access the Galactic latitudes (in deg) from
`data`.
lat_column : str
Key used to access the Galactic latitudes (in deg) from
`data`.
x_column, y_column, z_column : str
Keys used to access the coordinates (in kpc) from
`data`.
frequency : astropy.units.Quantity
Frequency of the measurement (if relevant)
tag : str
Extra information associated with the observable.
* 'I' - total intensity (in unit K-cmb)
* 'Q' - Stokes Q (in unit K-cmb, IAU convention)
* 'U' - Stokes U (in unit K-cmb, IAU convention)
* 'PI' - polarisation intensity (in unit K-cmb)
* 'PA' - polarisation angle (in unit rad, IAU convention)
"""
def __init__(self, data, name, data_column=None, units=None,
coordinates_type='galactic', lon_column='lon', lat_column='lat',
x_column='x', y_column='y', z_column='z',
error_column=None, frequency='nan', tag='nan'):
self.NAME = name
super().__init__()
if data_column is None:
data_column=name
self._data = u.Quantity(data[data_column], units, copy=False)
if coordinates_type == 'galactic':
self.coords = {'type': coordinates_type,
'lon': u.Quantity(data[lon_column], unit=u.deg),
'lat': u.Quantity(data[lat_column], unit=u.deg)}
elif coordinates_type == 'cartesian':
self.coords = {'type': coordinates_type,
'x': u.Quantity(data[x_column], unit=u.kpc),
'y': u.Quantity(data[y_column], unit=u.kpc),
'z': u.Quantity(data[z_column], unit=u.kpc)}
elif coordinates_type == None:
pass
else:
raise ValueError('Unknown coordinates_type!')
self.frequency = frequency
self.tag = tag
if error_column is not None:
self._error = u.Quantity(data[error_column], units, copy=False)
self.Nside = "tab"
[docs]class HEALPixDataset(Dataset):
"""
Base class for HEALPix datasets, which are input as
a simple 1D-array without explicit coordinate information
"""
def __init__(self, data, error=None, cov=None, Nside=None):
super().__init__()
dataNside = np.uint(np.sqrt(data.size/12))
if Nside is None:
Nside = dataNside
try:
assert 12*int(Nside)**2 == data.size
except AssertionError:
print(12*int(Nside)**2, data.size)
raise
self.Nside = str(Nside)
assert len(data.shape)==1
self._data = data
if cov is not None:
assert error is None
self._cov = distribute_matrix(cov)
else:
self._error = error
[docs]class FaradayDepthHEALPixDataset(HEALPixDataset):
r"""
Stores a Faraday depth map into an IMAGINE-compatible
dataset
Parameters
----------
data : numpy.ndarray
1D-array containing the HEALPix map
Nside : int, optional
For extra internal consistency checking. If `Nside` is present,
it will be checked whether :math:`12\times N_{side}^2` matches
Attributes
-------
data
Data in ObservablesDict-compatible shape
key
Standard key associated with this observable
"""
# Class attributes
NAME = 'fd'
[docs]class DispersionMeasureHEALPixDataset(HEALPixDataset):
r"""
Stores a dispersion measure map into an IMAGINE-compatible
dataset
Parameters
----------
data : numpy.ndarray
1D-array containing the HEALPix map
Nside : int, optional
For extra internal consistency checking. If `Nside` is present,
it will be checked whether :math:`12\times N_{side}^2` matches
Attributes
-------
data
Data in ObservablesDict-compatible shape
key
Standard key associated with this observable
"""
# Class attributes
NAME = 'dm'
[docs]class SynchrotronHEALPixDataset(HEALPixDataset):
r"""
Stores a synchrotron emission map into an IMAGINE-compatible
dataset. This can be Stokes parameters, total and polarised
intensity, and polarisation angle.
The parameter `type` and the units of the map in `data` must follow:
* 'I' - total intensity (in unit K-cmb)
* 'Q' - Stokes Q (in unit K-cmb, IAU convention)
* 'U' - Stokes U (in unit K-cmb, IAU convention)
* 'PI' - polarisation intensity (in unit K-cmb)
* 'PA' - polarisation angle (in unit rad, IAU convention)
Parameters
----------
data : numpy.ndarray
1D-array containing the HEALPix map
frequency : astropy.units.Quantity
Frequency of the measurement (if relevant)
Nside : int, optional
For extra internal consistency checking. If `Nside` is present,
it will be checked whether :math:`12\times N_{side}^2` matches data.size
type : str
The type of map being supplied in `data`.
Attributes
-------
data
Data in ObservablesDict-compatible shape
key
Standard key associated with this observable
"""
# Class attributes
NAME = 'sync'
def __init__(self, data, frequency, type,
error=None, cov=None, Nside=None):
super().__init__(data, error=error, cov=cov, Nside=Nside)
self.frequency = frequency
# Checks whether the type is valid
assert type in ['I', 'Q', 'U', 'PI', 'PA']
self.tag = type