Full data access API

HDF5 query functions

pycasso.h5datacube.listRuns(synthesisFile, baseId=None, qbickRunId=None)

Prints STARLIGHT runs contained in synthesis dataset.

Parameters :

synthesisFile : string

Path to an HDF5 file containing the synthesis.

baseId : string

Show results which use base given by baseId.

qbickRunId : string

Show results which use qbick run given by qbickRunId.

pycasso.h5datacube.listBases(synthesisFile)

Prints bases used in STARLIGHT runs contained in synthesis dataset..

Parameters :

synthesisFile : string

Path to an HDF5 file containing the synthesis.

pycasso.h5datacube.listQbickRuns(synthesisFile)

Prints qbick runs contained in synthesis dataset.

Parameters :

synthesisFile : string

Path to an HDF5 file containing the synthesis.

Synthesis data abstraction class

class pycasso.h5datacube.h5Q3DataCube(synthesisFile, runId, califaId, smooth=True)

Interface to CALIFA HDF5 synthesis dataset.

baseId = None

Identifier of the base used in the synthesis.

qbickRunId = None

Identifier of the qbick run used in the synthesis.

pycassoVersion = None

Version of PyCASSO used to import the synthesis.

qVersion = None

Version of qbick.

galaxyName = None

NED name of the galaxy.

califaID = None

CALIFA name of the galaxy.

N_age = None

Number of ages in base.

N_met = None

Number of metallicities in base.

N_base = None

Number of elements in base.

N_x = None

Number of pixels in X direction.

N_y = None

Number of pixels in Y direction.

N_zone = None

Number of voronoi zones.

Nl_obs = None

Number of wavelengths.

x0 = None

X position of the galaxy center.

y0 = None

Y position of the galaxy center.

HLR_pix = None

Half light radius in pixels.

HLR_pc = None

Half light radius in parsecs.

distance_Mpc = None

Distance to galaxy in megaparsecs.

parsecPerPixel = None

Pixel scale in parsecs.

PIXSIZE = None

Pixel scale in arcseconds.

flux_unit = None

Half light radius in pixels.

q_norm = None

\(A_\lambda / A_V\). Used in reddening calculations.

l_ini = None

Initial wavelength of spectra, in angstroms.

l_fin = None

Final wavelength of spectra, in angstroms.

dl = None

Wavelength step, in angstroms.

masterListData = None

Dictionary containing the masterlist data for this galaxy.

keywords = None

Dictionary containing all synthesis, base and qbick keywords. Keys are uppercase.

header = None

Alias for keywords. For backwards compatibility.

integrated_keywords = None

Dictionary containing additional synthesis data for integrated spectra. Keys are uppercase.

pa = None

Position angle of the ellipse used in radial profiles.

ba = None

Ellipticity \(b/a\) of the ellipse used in radial profiles.

pixelDistance__yx = None

Distance of each pixel to the galaxy center.

  • Units: [pixel]
  • Shape: (N_y, N_x)
pixelAngle__yx = None

Angle of each pixel relative to the semimajor axis, corrected for a disk perspective. See pa and ba. Ranges from \(-\pi\) to \(+\pi\).

  • Units: [radians]
  • Shape: (N_y, N_x)
zoneArea_pix = None

Area in pixels of each zone.

  • Units: [pixel]
  • Shape: (N_Zone)
zoneArea_pc2 = None

Area in square parsecs of each zone.

fill_value = None

Fill value used when building images from zones.

loadGalaxy(runId, califaId, smooth=True)

Load a galaxy from a run, and update the metadata.

Parameters :

runId : string

Id of the synthesis run.

califaId : string

CALIFA identifier of the galaxy. Ex.: K0001.

TODO: test LoadGalaxy() :

ageBase

Ages of the base.

  • Units: \([Yr]\)
  • Shape: (N_zone)
metBase

Metalicities of the base.

  • Units: dimensionless
  • Shape: (N_zone)
zonePos

A recarray (x, y) containing the position of the center of each zone.

  • Units: pixel
  • Shape: (N_zone)
qSignal

Signal at wavelength WINDOWSN.

  • Units: \([erg / s / cm^2 / \overset{\circ}{A}]\)
  • Shape: (N_y, N_x)
qNoise

TODO: Add description of cubes.

qSn

S/N at wavelength WINDOWSN.

  • Units: dimensionless
  • Shape: (N_y, N_x)
qSignalUnmasked

Image at wavelength WINDOWSN (quick, not masked).

  • Units: \([erg / s / cm^2 / \overset{\circ}{A}]\)
  • Shape: (N_y, N_x)
qNoiseUnmasked

TODO: Add description of qNoiseUnmasked.

qPipeNoise

TODO: Add description of qPipeNoise.

qZones

Voronoi/segmentation zones (bins).

  • Units: index
  • Shape: (N_y, N_x)
qZonesSn

S/N in Voronoi zones.

  • Units: dimensionless
  • Shape: (N_y, N_x)
qZonesNoise

Noise (RMS) in Voronoi zones.

  • Units: \([erg / s / cm^2 / \overset{\circ}{A}]\)
  • Shape: (N_y, N_x)
qFlagRatio

Ratio of flags in window at wavelength WINDOWSN

  • Units: \([\%]\)
  • Shape: (N_y, N_x)
qZonesSnOrig

S/N in Voronoi zones (unresampled, beta).

  • Units: dimensionless
  • Shape: (N_y, N_x)
qZonesNoiseOrig

Noise (RMS) in Voronoi zones (unresampled, beta).

  • Units: \([erg / s / cm^2 / \overset{\circ}{A}]\)
  • Shape: (N_y, N_x)
qSegmentationSn

Voronoi program output S/N.

  • Units: dimensionless
  • Shape: (N_y, N_x)
qPipeNoiseOrig

Noise image derived from formal errors (unresampled).

  • Units: \([erg / s / cm^2 / \overset{\circ}{A}]\)
  • Shape: (N_y, N_x)
qPipeZonesNoiseOrig

Noise image derived from formal errors in Voronoi zones (unresampled).

  • Units: \([erg / s / cm^2 / \overset{\circ}{A}]\)
  • Shape: (N_y, N_x)
qSpatialMask

TODO: Add description of qSpatialMask.

qSnMask

TODO: Add description of qSnMask.

qFilledMask

TODO: Add description of qFilledMask.

qMask

Boolean image mask used for data.

  • Units: bool
  • Shape: (N_y, N_x)
fobs_norm

Flux in norm window for output spectra, in voronoi zones.

  • Units: \([erg / s / cm^2 / \overset{\circ}{A}]\)
  • Shape: (Nl_obs, N_zone)
f_obs

Observed flux (input spectra for the synthesis), in voronoi zones.

  • Units: \([erg / s / cm^2 / \overset{\circ}{A}]\)
  • Shape: (Nl_obs, N_zone)
f_err

Error in observed spetra, in voronoi zones.

  • Units: \([erg / s / cm^2 / \overset{\circ}{A}]\)
  • Shape: (Nl_obs, N_zone)
f_flag

Flagged spaxels, in voronoi zones.

FIXME: describe flags.

  • Units: dimensionless
  • Shape: (Nl_obs, N_zone)
integrated_f_obs

Observed flux (input spectra for the synthesis), in integrated spectrum.

  • Units: \([erg / s / cm^2 / \overset{\circ}{A}]\)
  • Shape: (Nl_obs)
integrated_f_err

Error in integrated observed spetrum.

  • Units: \([erg / s / cm^2 / \overset{\circ}{A}]\)
  • Shape: (Nl_obs)
integrated_f_flag

Flagged spaxels, in integrated spectrum.

FIXME: describe flags.

  • Units: dimensionless
  • Shape: (Nl_obs)
Mstars

Fraction of the initial stellar mass for a given population that is still trapped inside stars.

  • Units: dimensionless
  • Shape: (N_age, N_met)
fbase_norm

TODO: Add description of cubes.

popx

Light fractions for each population, in voronoi zones.

    • Units: \([\%]\)
  • Shape: (N_age, N_met, N_zone)
popmu_cor

Current mass fractions for each population, in voronoi zones.

  • Units: \([\%]\)
  • Shape: (N_age, N_met, N_zone)
popmu_ini

Initial mass fractions for each population, in voronoi zones.

  • Units: \([\%]\)
  • Shape: (N_age, N_met, N_zone)
popAV_tot

Extinction for each population, in voronoi zones.

  • Units: \([mag]\)
  • Shape: (N_age, N_met, N_zone)
popexAV_flag

TODO: Add description of popexAV_flag.

SSP_chi2r

TODO: Add description of cubes.

SSP_adev

TODO: Add description of cubes.

SSP_AV

TODO: Add description of cubes.

SSP_x

TODO: Add description of cubes.

Lobs_norm

Luminosity density in norm window, in voronoi zones.

  • Units: \([L_\odot/\overset{\circ}{A}]\)
  • Shape: (N_zone)
Mini_tot

Initial mass for each population, in voronoi zones.

  • Units: \([M_\odot]\)
  • Shape: (N_zone)
Mcor_tot

Current mass for each population, in voronoi zones.

  • Units: \([M_\odot]\)
  • Shape: (N_zone)
A_V

Extinction by dust, in voronoi zones.

  • Units: \([mag]\)
  • Shape: (N_zone)
v_0

Velocity displacement, in voronoi zones.

  • Units: \([km/s]\)
  • Shape: (N_zone)
v_d

Velocity dispersion, in voronoi zones.

  • Units: \([km/s]\)
  • Shape: (N_zone)
adev

Mean absolute relative deviation, in percent, only for the Nl_eff points actually used in the synthesis.

  • Units: \([\%]\)
  • Shape: (N_zone)
index_Best_SSP

Best single SSP fit, to use in SSP_*, in voronoi zones.

FIXME: index_Best_SSP is j or (age,met)?

  • Units: index
  • Shape: (N_zone)
NOl_eff

Number of OK wavelengths in input spectrum, in voronoi zones.

  • Units: dimensionless
  • Shape: (N_zone)
Nl_eff

Number of wavelengths actually used in spectral fit, in voronoi zones.

  • Units: dimensionless
  • Shape: (N_zone)
Ntot_clipped

Number of wavelengths clipped, in voronoi zones.

  • Units: dimensionless
  • Shape: (N_zone)
Nglobal_steps

Number of steps in spectral fitting, in voronoi zones.

  • Units: dimensionless
  • Shape: (N_zone)
chi2

\(\chi^2 / Nl_{eff}\) of the fit, in voronoi zones.

  • Units: dimensionless
  • Shape: (N_zone)
chi2_TOT

Total \(\chi^2\) of the fit, in voronoi zones.

  • Units: dimensionless
  • Shape: (N_zone)
chi2_Opt

\(\chi^2\) of the optical fit, in voronoi zones.

  • Units: dimensionless
  • Shape: (N_zone)
chi2_FIR

\(\chi^2\) of the far-IR fit, in voronoi zones.

  • Units: dimensionless
  • Shape: (N_zone)
chi2_QHR

\(\chi^2\) of the QH-related fit in voronoi zones.

  • Units: dimensionless
  • Shape: (N_zone)
chi2_PHO

\(\chi^2\) of the photometric fit, in voronoi zones.

  • Units: dimensionless
  • Shape: (N_zone)
f_syn

Synthetic spectra, in voronoi zones.

  • Units: \([erg / s / cm^2 / \overset{\circ}{A}]\)
  • Shape: (Nl_obs, N_zone)
A_V__yx

Spatially resolved extinction by dust.

  • Units: \([mag]\)
  • Shape: (N_y, N_x)
DeRed_LobnSD__tZyx

Spatially resolved “dereddened” luminosity surface density of each population in normalization window.

  • Units: \([L_\odot / pc^2]\)
  • Shape: (N_age, N_met, N_y, N_x)
DeRed_LobnSD__yx

“Dereddened” luminosity surface density of each population in normalization window.

  • Units: \([L_\odot / pc^2]\)
  • Shape: (N_y, N_x)
DeRed_Lobn__tZz

“Dereddened” luminosity of each population in normalization window, in voronoi zones.

  • Units: \([L_\odot]\)
  • Shape: (N_age, N_met, N_zone)
DeRed_Lobn__z

“Dereddened” luminosity in normalization window, in voronoi zones.

  • Units: \([L_\odot]\)
  • Shape: (N_zone)
DeRed_M2L__yx

Spatially resolved “dereddened” mass to light ratio.

  • Units: \([M_\odot / L_\odot]\)
  • Shape: (N_y, N_x)
LobnSD__tZyx

Spatially resolved luminosity surface density of each population in normalization window.

  • Units: \([L_\odot / pc^2]\)
  • Shape: (N_age, N_met, N_y, N_x)
LobnSD__yx

Luminosity surface density of each population in normalization window.

  • Units: \([L_\odot / pc^2]\)
  • Shape: (N_y, N_x)
Lobn__tZz

Luminosity of each population in normalization window, in voronoi zones.

  • Units: \([L_\odot]\)
  • Shape: (N_age, N_met, N_zone)
Lobn__z

Luminosity in normalization window, in voronoi zones.

  • Units: \([L_\odot]\)
  • Shape: (N_zone)
M2L__yx

Spatially resolved mass to light ratio.

  • Units: \([M_\odot / L_\odot]\)
  • Shape: (N_y, N_x)
McorSD__tZyx

Spatially resolved current mass surface density of each population.

  • Units: \([M_\odot / pc^2]\)
  • Shape: (N_age, N_met, N_y, N_x)
McorSD__yx

Spatially resolved current mass surface density.

  • Units: \([M_\odot / pc^2]\)
  • Shape: (N_y, N_x)
Mcor__tZz

Current mass of each population, in voronoi zones.

  • Units: \([M_\odot]\)
  • Shape: (N_age, N_met, N_zone)
Mcor__z

Current mass, in voronoi zones.

  • Units: \([M_\odot]\)
  • Shape: (N_age, N_met, N_zone)
MiniSD__tZyx

Spatially resolved initial mass surface density of each population.

  • Units: \([M_\odot / pc^2]\)
  • Shape: (N_age, N_met, N_y, N_x)
MiniSD__yx

Spatially resolved initial mass surface density.

  • Units: \([M_\odot / pc^2]\)
  • Shape: (N_y, N_x)
Mini__tZz

Initial mass of each population, in voronoi zones.

  • Units: \([M_\odot]\)
  • Shape: (N_age, N_met, N_zone)
Mini__z

Initial mass, in voronoi zones.

  • Units: \([M_\odot]\)
  • Shape: (N_age, N_met, N_zone)
aZ_flux__yx

Spatially resolved, flux-weighted average metallicity.

  • Units: dimensionless
  • Shape: (N_y, N_x)
aZ_flux__z

Flux-weighted average metallicity, in voronoi zones.

  • Units: dimensionless
  • Shape: (N_zone)
aZ_mass__yx

Spatially resolved, mass-weighted average metallicity.

  • Units: dimensionless
  • Shape: (N_y, N_x)
aZ_mass__z

Mass-weighted average metallicity, in voronoi zones.

  • Units: dimensionless
  • Shape: (N_zone)
adevS

From Cid @ 26/05/2012: Here’s my request for pycasso reader: That it defines a new figure of merit analogous to adev, but which uses the synthetic flux in the denominator instead of the observed one. This is the adevS__z thing defined below in awful python. Why? Well, despite all our care there are still some non-flagged pixels with very low fluxes, which screws up adev, and this alternative definition fixes it.

Original code:

>>> adevS__z = np.zeros((self.nZones))
>>> for i_z in np.arange(self.nZones):
>>>    _a = np.abs(self.f_obs[:,i_z] - self.f_syn[:,i_z]) / self.f_syn[:,i_z]
>>>    _b = _a[self.f_wei[:,i_z] > 0]
>>>    adevS__z[i_z] = 100.0 * _b.mean()
>>> return adevS__z
Returns :aDevS : array of length (nZones)
adev__yx

Mean absolute relative deviation, in percent, only for the Nl_eff points actually used in the synthesis.

  • Units: \([\%]\)
  • Shape: (N_y, N_x)
ageBaseMx(n=2)

Return ageBase as a row matrix to simplify multiplications.

Parameters :

n : integer

Dimension of the matrix.

Returns :

ageBaseNd : n-D array

ageBase as a row matrix.

See also

metBaseMx, zoneAreaMx

area_qConvexHull

Area, in pixels, of the convex hull of qMask. This is the area to use when dealing with qConvexHull.

area_qHollowPixels

Aarea, in pixels, of the masked regions. This is the area to use when dealing with qHollowPixels.

area_qMask

Area, in pixels, of the data regions. This is the area to use when dealing with qMask.

at_flux__yx

Spatially resolved, flux-weighted average log. age.

  • Units: \([\log Gyr]\)
  • Shape: (N_y, N_x)
at_flux__z

Flux-weighted average log. age, in voronoi zones.

  • Units: \([\log Gyr]\)
  • Shape: (N_zone)
at_mass__yx

Spatially resolved, mass-weighted average log. age.

  • Units: \([\log Gyr]\)
  • Shape: (N_y, N_x)
at_mass__z

Mass-weighted average log. age, in voronoi zones.

  • Units: \([\log Gyr]\)
  • Shape: (N_zone)
attrs

Attributes attached to this object

azimuthalProfileNd(prop, bin_a, bin_r, rad_scale=None, mask=None, a__yx=None, r__yx=None, mode='mean', return_npts=False)

Calculate the azimuthal profile of a property. The last two dimensions of prop must be of legth N_y and N_x.

Parameters :

prop : array

Image of property to calculate the radial profile.

bin_a : array

Angular bin boundaries in radians.

bin_r : array

Radial bin boundaries in units of rad_scale.

rad_scale : float, optional

Scale of the bins, in pixels. Defaults to HLR_pix.

mask : array, optional

Mask containing the pixels to use in the radial profile. Must have the shape (N_y, N_x). Defaults to qMask.

r__yx : array, optional

Distance of each pixel to the galaxy center, in pixels. Must have the shape (N_y, N_x). Defaults to pixelDistance__yx.

a__yx : array, optional

Angle associated with each pixel. Must have the shape (N_y, N_x). Defaults to pixelAngle__yx.

mode : string, optional

One of:
  • 'mean': Compute the mean inside the radial bins (default).
  • 'median': Compute the median inside the radial bins.
  • 'sum': Compute the sum inside the radial bins.

return_npts : bool, optional

If set to True, also return the number of points inside each bin. Defaults to False.

Returns :

azProf : array

Array containing the radial and azimuthal profiles as the last dimensions. Note that radProf.shape[-2] == (len(bin_a) - 1) and radProf.shape[-1] == (len(bin_r) - 1).

npts : array, optional

The number of points inside each bin, only if return_npts is set to True.

Examples

Load the dataset:

>>> from pycasso.h5datacube import h5Q3DataCube
>>> K = h5Q3DataCube('qalifa-synthesis.002.h5', 'run001', 'K0277')

Create the radial bins from 0.0 to 3.0 in 0.5 steps, using rad_scale units, and the angular bins as 21 boundaries between -180 and 180 degrees, converted to radians.

>>> import numpy as np
>>> bin_r = np.arange(0.0, 3.0 + 0.5, 0.5)
>>> bin_a = np.linspace(-180.0, 180.0, 21) / 180.0 * np.pi

Calculate the azimuthal profile of the time resolved mass surface density, using a radial scale of 10.5 pixels.

>>> McorSD__tyx = K.McorSD__tZyx.sum(axis=1)
>>> McorSD__taR = K.azimuthalProfileNd(McorSD__tyx, bin_a, bin_r, rad_scale=10.5)

Note that bin_a contains the bin boundaries, it is not fit for plotting along with McorSD__tZaR. We will use the bin centers.

>>> bin_a_center = (bin_a[:-1] + bin_a[1:]) / 2.0

Plot the azimuthal profile for the first and last radial bins, summing in all ages. Note that McorSD__taR.shape is (N_ages, len(bin_a) - 1, len(bin_r) - 1).

>>> import matplotlib.pyplot as plt
>>> plt.plot(bin_a_center, McorSD__taR[...,0].sum(axis=0)
>>> plt.plot(bin_a_center, McorSD__taR[...,-1].sum(axis=0)
chi2__yx

\(\chi^2 / Nl_{eff}\) of the fit.

  • Units: dimensionless
  • Shape: (N_y, N_x)
close()

Close the file. All open objects become invalid

copy(source, dest, name=None)

Copy an object or group.

The source can be a path, Group, Dataset, or Datatype object. The destination can be either a path or a Group object. The source and destinations need not be in the same file.

If the source is a Group object, all objects contained in that group will be copied recursively.

When the destination is a Group object, by default the target will be created in that group with its current name (basename of obj.name). You can override that by setting “name” to a string.

Example:

>>> f = File('myfile.hdf5')
>>> f.listnames()
['MyGroup']
>>> f.copy('MyGroup', 'MyCopy')
>>> f.listnames()
['MyGroup', 'MyCopy']
create_dataset(name, shape=None, dtype=None, data=None, **kwds)

Create a new HDF5 dataset

name
Name of the dataset (absolute or relative). Provide None to make an anonymous dataset.
shape
Dataset shape. Use “()” for scalar datasets. Required if “data” isn’t provided.
dtype
Numpy dtype or string. If omitted, dtype(‘f’) will be used. Required if “data” isn’t provided; otherwise, overrides data array’s dtype.
data
Provide data to initialize the dataset. If used, you can omit shape and dtype arguments.

Keyword-only arguments:

chunks
(Tuple) Chunk shape, or True to enable auto-chunking.
maxshape
(Tuple) Make the dataset resizable up to this shape. Use None for axes you want to be unlimited.
compression
(String) Compression strategy. Legal values are ‘gzip’, ‘szip’, ‘lzf’. Can also use an integer in range(10) indicating gzip.
compression_opts
Compression settings. This is an integer for gzip, 2-tuple for szip, etc.
shuffle
(T/F) Enable shuffle filter.
fletcher32
(T/F) Enable fletcher32 error detection.
fillvalue
(Scalar) Use this value for uninitialized parts of the dataset.
create_group(name)

Create and return a new subgroup.

Name may be absolute or relative. Fails if the target name already exists.

driver

Low-level HDF5 file driver used to open file

f_err__lyx

Spatially resolved error in observed spetra.

  • Units: \([erg / s / cm^2 / \overset{\circ}{A}]\)
  • Shape: (Nl_obs, N_y, N_x)
f_flag__lyx

Spatially resolved flagged spaxels.

FIXME: describe flags.

  • Units: dimensionless
  • Shape: (Nl_obs, N_y, N_x)
f_obs__lyx

Spatially resolved observed flux (input spectra for the synthesis).

  • Units: \([erg / s / cm^2 / \overset{\circ}{A}]\)
  • Shape: (Nl_obs, N_y, N_x)
f_syn__lyx

Spatially resolved synthetic spectra.

  • Units: \([erg / s / cm^2 / \overset{\circ}{A}]\)
  • Shape: (Nl_obs, N_y, N_x)
f_wei

Weight of the spaxels in the input spectra. This is the weight actually used by the synthesis, after clipping, etc. Values for voronoi zones.

FIXME: describe flags and weights.

  • Units: dimensionless
  • Shape: (Nl_obs, N_zone)
f_wei__lyx

Spatially resolved weight of the spaxels in the input spectra. This is the weight actually used by the synthesis, after clipping, etc.

FIXME: describe flags and weights.

  • Units: dimensionless
  • Shape: (Nl_obs, N_y, N_x)
fid

File ID (backwards compatibility)

file

Return a File instance associated with this object

filename

File name on disk

fillImage(prop, prop__r=None, r=None, r__yx=None, mode='convex')

Fill a 2-D the masked pixels of an image of a property using the values of a radial profile. The pixels to fill are chosen by mode.

Parameters :

prop : array

Property (2-D) to fill the hollow pixels.

prop__r : the radial profile to use to fill the data.

If not set, calculate from prop.

r : array

The radial distances for prop__r values. If not set (or prop__r is not set), use a 1 pixel step.

r__yx : array

A 2-D image containing the geometry of the image. If not set, use pixelDistance__yx.

mode : string

One of:
  • 'convex': fill entire convex hull.
  • 'hollow': fill only hollow pixels.
Returns :

prop_fill : array

A 2-D image of prop, with the missing pixels filled.

mask : array

The effective mask for the filled image.

flush()

Tell the HDF5 library to flush its buffers.

get(name, default=None, getclass=False, getlink=False)

Retrieve an item or other information.

“name” given only:
Return the item, or “default” if it doesn’t exist
“getclass” is True:
Return the class of object (Group, Dataset, etc.), or “default” if nothing with that name exists
“getlink” is True:
Return HardLink, SoftLink or ExternalLink instances. Return “default” if nothing with that name exists.
“getlink” and “getclass” are True:
Return HardLink, SoftLink and ExternalLink classes. Return “default” if nothing with that name exists.

Example:

>>> cls = group.get('foo', getclass=True)
>>> if cls == SoftLink:
...     print '"foo" is a soft link!'
getDezonificationWeight(smooth, prop=None)

Create the weight image for dezonification. If smooth is True, use prop image to weight the pixels in each zone. Otherwise use the zone area. If prop is not set, use qSignal.

Here we use a scheme similar to zoneToYX(), when using smooth dezonification, except that we use numpy.histogram() to calculate the weight of the pixels.

Parameters :

smooth : boolean

Enable or disable smooth dezonification.

prop : array, optional

Image to use as dezonification weights if smooth is True. If set to None, use qSignal.

getEllipseParams()

Estimate ellipticity and orientation of the galaxy using the “Stokes parameters”, as described in: http://adsabs.harvard.edu/abs/2002AJ....123..485S The image used is qSignal.

Returns :

pa : float

Position angle in radians, counter-clockwise relative to the positive X axis.

ba : float

Ellipticity, defined as the ratio between the semiminor axis and the semimajor axis (\(b/a\)).

getHLR_pix(pa=None, ba=None, fill=False)

Find the half light radius using the image at the normalization window (qSignal). Using radial bins of 1 pixel, calculate the cumulative sum of luminosity. The HLR is the radius where the cum. sum reaches 50% of its peak value.

Parameters :

fill : boolean

Whether to fill the hollow areas before calculating the HLR.

Returns :

HLR : float

The half light radius, in pixels.

Notes

This value should be close to $dfrac{HLR_{circular}}{sqrt{b/a}}$ if the isophotes of the galaxy are all ellipses with parameters p.a. and b/a.

getHalfRadius(prop, fill=False)

Find the half radius of the desired property. Using radial bins of 1 pixel, calculate the cumulative sum of prop. The “half prop radius” is the radius where the cumulative sum reaches 50% of its peak value.

Parameters :

prop : array

Image to get the half radius.

fill : boolean

Whether to fill the hollow areas before the calculations. The filling is done using the average value in the respective radial distance of the missing pixels.

Returns :

HXR : float

The half prop radius, in pixels.

Notes

This value should be close to $dfrac{HLR_{circular}}{sqrt{b/a}}$ if the isophotes of the galaxy are all ellipses with parameters p.a. and b/a.

Examples

Load the dataset:

>>> from pycasso.h5datacube import h5Q3DataCube
>>> K = h5Q3DataCube('qalifa-synthesis.002.h5', 'run001', 'K0277')

Calculate the half mass radius, in pixels.

>>> HMR_pix = K.getHalfRadius(K.McorSD)

Calculate the radial profile of extinction by dust, using HMR_pix as radial scale, in bins from 0.0 to 3.0, in 0.1 steps.

>>> import numpy as np
>>> bin_r = np.arange(0.0, 3.0 + 0.1, 0.1)
>>> A_V__r = K.radialProfile(K.A_V__yx, bin_r, rad_scale=HMR_pix)

Note that bin_r is the bin boundaries, it is not fit for plotting along with A_V__r. We need the bin centers.

>>> bin_center = bin_center = (bin_r[:-1] + bin_r[1:]) / 2.0

Plot the Radial profile.

>>> import matplotlib.pyplot as plt
>>> plt.plot(bin_center, A_V__r)
getPixelAngle(units='radians', pa=None, ba=None)

Return an image (numpy.ndarray of same shape as :attr`qSignal`) of the angle in radians (default) of each pixel, relative from the axis of the position angle pa. The projection is fixed assuming the galaxy is a disk, throught the ellipticity parameter ba.

Parameters :

units : string, optional

Angle units to use, one of:
  • 'radians': angles in radians, from \(-\pi\) to \(+\pi\) (default).
  • 'degrees': angles in degrees, from \(-180.0\) to \(+180.0\) (default).

pa : float, optional

Position angle in radians, counter-clockwise relative to the positive X axis.

ba : float, optional

Ellipticity, defined as the ratio between the semiminor axis and the semimajor axis (\(b/a\)).

Returns :

pixelDistance : 2-D array

Image containing the pixel distances.

getPixelDistance(use_HLR_units=True, pixel_scale=None, pa=None, ba=None)

Return an image (numpy.ndarray of same shape as :attr`qSignal`) of the distance from the center of the galaxy (x0, y0) in HLR units (default), assuming a projected disk.

Parameters :

use_HLR_units : boolean, optional

Whether to use units of half light radius or pixels.

pixel_scale : float, optional

Pixel distance scale, used if use_HLR_units is False. If not set, do not scale the distance.

pa : float, optional

Position angle in radians, counter-clockwise relative to the positive X axis.

ba : float, optional

Ellipticity, defined as the ratio between the semiminor axis and the semimajor axis (\(b/a\)).

Returns :

pixelDistance : 2-D array

Image containing the pixel distances.

See also

getPixelAngle

getSFR(logtc_ini=None, logtc_fin=None, logtc_step=0.05, logtc_FWHM=0.5)

Calculate the star formation rate using a smooth-resampled age base, as prescribed by Asari et al. (2007) <http://adsabs.harvard.edu/abs/2007MNRAS.381..263A>.

Parameters :

logtc_ini : float

Logarithm (base 10) of initial age. Defaults to logtb.min().

logtc_fin : float

Logarithm (base 10) of final age. Defaults to logtb.max().

logtc_step : float

Logarithm (base 10) of age step. Defaults to 0.05.

logtc_FWHM : float

Width of the age smoothing kernel. Defaults to 0.5.

Returns :

SFR : array

Star formation rate.

logtc : array

Logarithm (base 10) of smooth-resampled age base.

See also

pystarlight.util.StarlightUtils.calcSFR

growthInAge(prop, mask=None, relative=False)

TODO: move to utils?

Calculate the cumulative growth in age space. In age = 0, we have the sum of prop in age, and as age increases the sum of prop decreases.

Parameters :

prop : array

A property having the age as first dimension.

mask : array

Boolean mask for all dimensions except the first (age). If set, the values will be summed, that is, the return array will only have the age dimension.

relative : bool

If True, normalize by the maximum value of the cumulative sum of prop.

Returns :

cs : array

The growth in age space of prop.

Examples

Load the dataset:

>>> from pycasso.h5datacube import h5Q3DataCube
>>> K = h5Q3DataCube('qalifa-synthesis.002.h5', 'run001', 'K0277')

Get the spatially resolved mass surface density for each age.

>>> MiniSD__tyx = K.MiniSD__tZyx.sum(axis=1)

Calculate the mass build up for the whole galaxy (notice the “all data” mask, K.qMask). The mask could be any selection of data pixels.

>>> MtotGrow__t = K.growthInAge(MiniSD__tyx, mask=K.qMask)

MtotGrow__t is a surface mass density, in \(M_\odot / pc^2\). Converting to absolute mass.

>>> MtotGrow__t *= K.parsecPerPixel ** 2

Plot the mass buildup as a function of log. time.

>>> ages = np.log10(K.ageBase)
>>> import matplotlib.pyplot as plt
>>> plt.plot(ages, MtotGrow__t)
id

Low-level identifier appropriate for this object

items()

Get a list of tuples containing (name, object) pairs

iteritems()

Get an iterator over (name, object) pairs

iterkeys()

Get an iterator over member names

itervalues()

Get an iterator over member objects

keys()

Get a list containing member names

libver

File format version bounds (2-tuple: low, high)

logAgeBaseBins

Age bin edges, in log. Computes the bin edges as the bissection of the bins, expanding the borders accordingly.

  • Units: \([\log age / yr]\)
  • Shape: (N_age + 1)
metBaseMx(n=2)

Return metBase as a row matrix to simplify multiplications.

Parameters :

n : integer

Dimension of the matrix.

Returns :

metBaseNd : n-D array

metBase as a row matrix.

See also

ageBaseMx, zoneAreaMx

mode

Python mode used to open file

name

Return the full name of this object. None if anonymous.

parent

Return the parent group of this object.

This is always equivalent to obj.file[posixpath.dirname(obj.name)]. ValueError if this object is anonymous.

popAV_tot__tZyx

Spatially resolved extinction by dust for each population.

  • Units: \([mag]\)
  • Shape: (N_age, N_met, N_y, N_x)
popmu_cor__tZyx

Spatially resolved corrected mass fractions for each population.

  • Units: \([\%]\)
  • Shape: (N_age, N_met, N_y, N_x)
popmu_ini__tZyx

Spatially resolved initial mass fractions for each population.

  • Units: \([\%]\)
  • Shape: (N_age, N_met, N_y, N_x)
popx__tZyx

Spatially resolved light fractions for each population.

  • Units: \([\%]\)
  • Shape: (N_age, N_met, N_y, N_x)
qConvexHull

Convex hull of data mask.

  • Units: bool
  • Shape: (N-y, N_x)
qHollowPixels

Masked (bad) pixels inside the data mask.

  • Units: bool
  • Shape: (N-y, N_x)
radialProfile(prop, bin_r, rad_scale=None, mask=None, r__yx=None)

Calculate the radial profile of a property image (2-D array).

Parameters :

prop : array

Image of property to calculate the radial profile. Must be a 2-dimensional array.

bin_r : array

Radial bin boundaries in units of rad_scale.

rad_scale : float, optional

Scale of the bins, in pixels. Defaults to HLR_pix.

mask : array, optional

Mask containing the pixels to use in the radial profile. Must have the same dimensions as prop. Defaults to qMask.

r__yx : array, optional

Distance of each pixel to the galaxy center, in pixels. Must have the same dimensions as prop. Defaults to pixelDistance__yx.

Returns :

radProf : array

Array containing the radial profile. Note that len(radProf) == (len(bin_r) - 1)

Examples

Load the dataset:

>>> from pycasso.h5datacube import h5Q3DataCube
>>> K = h5Q3DataCube('qalifa-synthesis.002.h5', 'run001', 'K0277')

Create the bins from 0.0 to 3.0 in 0.1 steps, using rad_scale units.

>>> import numpy as np
>>> bin_r = np.arange(0.0, 3.0 + 0.1, 0.1)

Calculate the radial profile of mass surface density, using a radial scale of 10.5 pixels.

>>> McorSD__r = K.radialProfile(K.McorSD__yx, bin_r, rad_scale=10.5)

Note that bin_r is the bin boundaries, it is not fit for plotting along with McorSD__r. We need the bin centers.

>>> bin_center = bin_center = (bin_r[:-1] + bin_r[1:]) / 2.0

Plot the Radial profile.

>>> import matplotlib.pyplot as plt
>>> plt.plot(bin_center, McorSD__r)
radialProfileNd(prop, bin_r, rad_scale=None, mask=None, r__yx=None, mode='mean', return_npts=False)

Calculate the radial profile of a property. The last two dimensions of prop must be of length N_y and N_x.

Parameters :

prop : array

Image of property to calculate the radial profile.

bin_r : array

Radial bin boundaries in units of rad_scale.

rad_scale : float, optional

Scale of the bins, in pixels. Defaults to HLR_pix.

mask : array, optional

Mask containing the pixels to use in the radial profile. Must have the shape (N_y, N_x). Defaults to qMask.

r__yx : array, optional

Distance of each pixel to the galaxy center, in pixels. Must have the shape (N_y, N_x). Defaults to pixelDistance__yx.

mode : string, optional

One of:
  • 'mean': Compute the mean inside the radial bins (default).
  • 'median': Compute the median inside the radial bins.
  • 'sum': Compute the sum inside the radial bins.

return_npts : bool, optional

If set to True, also return the number of points inside each bin. Defaults to False.

Returns :

radProf : array

Array containing the radial profile as the last dimension. Note that radProf.shape[-1] == (len(bin_r) - 1)

npts : array, optional

The number of points inside each bin, only if return_npts is set to True.

Examples

Load the dataset:

>>> from pycasso.h5datacube import h5Q3DataCube
>>> K = h5Q3DataCube('qalifa-synthesis.002.h5', 'run001', 'K0277')

Create the bins from 0.0 to 3.0 in 0.1 steps, using rad_scale units.

>>> import numpy as np
>>> bin_r = np.arange(0.0, 3.0 + 0.1, 0.1)

Calculate the radial profile of the time resolved mass surface density, using a radial scale of 10.5 pixels.

>>> McorSD__tyx = K.McorSD__tZyx.sum(axis=1)
>>> McorSD__tZr = K.radialProfileNd(McorSD__tyx, bin_r, rad_scale=10.5)

Plot the Radial profile.

>>> import matplotlib.pyplot as plt
>>> plt.pcolormesh(np.log10(K.ageBase), bin_r, McorSD__r)
ref

An (opaque) HDF5 reference to this object

require_dataset(name, shape, dtype, exact=False, **kwds)

Open a dataset, creating it if it doesn’t exist.

If keyword “exact” is False (default), an existing dataset must have the same shape and a conversion-compatible dtype to be returned. If True, the shape and dtype must match exactly.

Other dataset keywords (see create_dataset) may be provided, but are only used if a new dataset is to be created.

Raises TypeError if an incompatible object already exists, or if the shape or dtype don’t match according to the above rules.

require_group(name)

Return a group, creating it if it doesn’t exist.

TypeError is raised if something with that name already exists that isn’t a group.

setGeometry(pa, ba, HLR_pix=None)

Change the geometry of the rings used when calculating radial profiles.

Parameters :

pa : float

Position angle in radians, counter-clockwise relative to the positive X axis.

ba : float

Ellipticity, defined as the ratio between the semiminor axis and the semimajor axis (\(b/a\)).

HLR_pix : float, optional

Effective radius

Examples

Load the dataset:

>>> from pycasso.h5datacube import h5Q3DataCube
>>> K = h5Q3DataCube('qalifa-synthesis.002.h5', 'run001', 'K0277')

Find the ellipse parameters.

>>> pa, ba = K.getEllipseParams()

Set the geometry, using a predefined value for HLR_pix.

>>> K.setGeometry(pa, ba, HLR_pix=10.5)

Get the distance for each pixel from the galaxy center.

>>> dist__yx = K.getPixelDistance()

Plot the distance image. Note that its shape should resemble the ellipticity of the galaxy (if it is well-behaved).

>>> import matplotlib.pyplot as plt
>>> plt.imshow(dist__yx)
setSmoothDezonification(smooth=True)

Enable or disable smooth dezonification. If smooth is True, use qSignal image to weight the pixels in each zone. Otherwise use the zone area.

Parameters :

smooth : boolean, optional

Enable or disable smooth dezonification. Defaults to True.

userblock_size

User block size (in bytes)

v_0__yx

Spatially resolved velocity displacement.

  • Units: \([km/s]\)
  • Shape: (N_y, N_x)
v_d__yx

Spatially resolved velocity dispersion.

  • Units: \([km/s]\)
  • Shape: (N_y, N_x)
values()

Get a list containing member objects

visit(func)

Recursively visit all names in this group and subgroups (HDF5 1.8).

You supply a callable (function, method or callable object); it will be called exactly once for each link in this group and every group below it. Your callable must conform to the signature:

func(<member name>) => <None or return value>

Returning None continues iteration, returning anything else stops and immediately returns that value from the visit method. No particular order of iteration within groups is guranteed.

Example:

>>> # List the entire contents of the file
>>> f = File("foo.hdf5")
>>> list_of_names = []
>>> f.visit(list_of_names.append)
visititems(func)

Recursively visit names and objects in this group (HDF5 1.8).

You supply a callable (function, method or callable object); it will be called exactly once for each link in this group and every group below it. Your callable must conform to the signature:

func(<member name>, <object>) => <None or return value>

Returning None continues iteration, returning anything else stops and immediately returns that value from the visit method. No particular order of iteration within groups is guranteed.

Example:

# Get a list of all datasets in the file >>> mylist = [] >>> def func(name, obj): ... if isinstance(obj, Dataset): ... mylist.append(name) ... >>> f = File(‘foo.hdf5’) >>> f.visititems(func)

zoneAreaMx(n)

Return zoneArea_pc2 as a row matrix to simplify multiplications.

Parameters :

n : integer

Dimension of the matrix.

Returns :

zoneAreaNd : n-D array

zone area as a row matrix.

See also

ageBaseMx, metBaseMx

zoneToRad(prop, bin_r, rad_scale=None, mask=None, r__yx=None, extensive=True)

Calculates the radial profile from a zone array. This method is a wrapper for zoneToYX() and radialProfile().

Parameters :

prop : array

Image of property to calculate the radial profile.

bin_r : array

Radial bin boundaries in units of rad_scale.

rad_scale : float, optional

Scale of the bins, in pixels. Defaults to HLR_pix.

mask : array, optional

Mask containing the pixels to use in the radial profile. Must have the same dimensions as prop. Defaults to qMask.

r__yx : array, optional

Distance of each pixel to the galaxy center, in pixels. Must have the same dimensions as prop. Defaults to pixelDistance__yx.

extensive : boolean, optional

If True, prop is extensive, use dezonification weights. Defaults to True.

Returns :

radProf : array

Array containing the radial profile. Note that len(radProf) == (len(bin_r) - 1)

Examples

Load the dataset:

>>> from pycasso.h5datacube import h5Q3DataCube
>>> K = h5Q3DataCube('qalifa-synthesis.002.h5', 'run001', 'K0277')

Create the bins from 0.0 to 3.0 in 0.1 steps, using rad_scale units.

>>> import numpy as np
>>> bin_r = np.arange(0.0, 3.0 + 0.1, 0.1)

Calculate the radial profile of mass surface density, using a radial scale of 10.5 pixels. The extensive option is the key for “surface density”.

>>> McorSD__r = K.zoneToRad(K.Mcor__z, bin_r, rad_scale=10.5, extensive=True)

Note that bin_r is the bin boundaries, it is not fit for plotting along with McorSD__r. We need the bin centers.

>>> bin_center = (bin_r[:-1] + bin_r[1:]) / 2.0

Plot the Radial profile.

>>> import matplotlib.pyplot as plt
>>> plt.plot(bin_center, McorSD__r)
zoneToYX(prop, extensive=True, surface_density=True, fill_value=None)

Convert a zone array to an image.

This scheme takes advantage of the qZones image, which has, for every pixel (x, y) the index of the corresponding zone. Using this array as a “smart index” for prop, we get to reconstruct the image.

Parameters :

prop : array

Property to be converted to an image. The zone dimension must be the rightmost dimension.

extensive : boolean, optional

If True, prop is extensive, use dezonification weights. Defaults to True.

surface_density : boolean, optional

If True, and extensive is True, divide the return value by the area in parsec^2 of the pixel. Defaults to True.

fill_value : float, optional

Fill value for masked pixels. Defaults to numpy.nan.

Returns :

prop__yx : array

The prop array converted to image. All dimensions are kept the same, except for the rightmost one, which is replaced by y and x.

Examples

Load the dataset:

>>> from pycasso.h5datacube import h5Q3DataCube
>>> K = h5Q3DataCube('qalifa-synthesis.002.h5', 'run001', 'K0277')

The A_V attribute contains the extincion for each zone.

>>> K.A_V.shape
(1638,)

Convert A_V to spatial coordinates. Note that the extinction is not extensive.

>>> A_V__yx = K.zoneToYX(K.A_V, extensive=False)
>>> A_V__yx.shape
(73, 77)

Plot the image.

>>> import matplotlib.pyplot as plt
>>> plt.imshow(A_V__yx)
chains_best_par

TODO: Add description of chains_best_par.

chains_ave_par

TODO: Add description of chains_ave_par.

chains_par

TODO: Add description of chains_par.

chains_best_LAx

TODO: Add description of chains_best_LAx.

chains_ave_LAx

TODO: Add description of chains_ave_LAx.

chains_LAx

TODO: Add description of chains_LAx.

chains_best_mu_cor

TODO: Add description of chains_best_mu_cor.

chains_ave_mu_cor

TODO: Add description of chains_ave_mu_cor.

chains_mu_cor

TODO: Add description of chains_mu_cor.

best_chi2

Best \(\chi^2\) (in chains), in voronoi zones.

  • Units: dimensionless
  • Shape: (N_zone)
ave_chi2

Average \(\chi^2\) (in chains), in voronoi zones.

  • Units: dimensionless
  • Shape: (N_zone)
cha_chi2

TODO: Add description of cubes.

best_Mcor

Best Mcor (in chains), in voronoi zones.

  • Units: \([M_\odot]\)
  • Shape: (N_zone)
ave_Mcor

Average Mcor (in chains), in voronoi zones.

  • Units: \([M_\odot]\)
  • Shape: (N_zone)
cha_Mcor

TODO: Add description of cha_Mcor.

integrated_popx

Light fractions for each population, in integrated spectrum.

  • Units: \([\%]\)
  • Shape: (N_age, N_met)
integrated_popmu_cor

Current mass fractions for each population, in integrated spectrum.

  • Units: \([\%]\)
  • Shape: (N_age, N_met)
integrated_popmu_ini

Current mass fractions for each population, in integrated spectrum.

  • Units: \([\%]\)
  • Shape: (N_age, N_met)
integrated_popAV_tot

Extinction by dust for each population, in integrated spectrum.

  • Units: \([mag]\)
  • Shape: (N_age, N_met)
integrated_popexAV_flag

TODO: Add description of popexAV_flag.

integrated_SSP_chi2r

TODO: Add description of integrated_SSP_chi2r.

integrated_SSP_adev

TODO: Add description of integrated_SSP_adev.

integrated_SSP_AV

TODO: Add description of integrated_SSP_AV.

integrated_SSP_x

TODO: Add description of integrated_SSP_x.

integrated_f_syn

Synthetic integrated spectrum.

  • Units: \([erg / s / cm^2 / \overset{\circ}{A}]\)
  • Shape: (Nl_obs)
integrated_f_wei

Weight of the spaxels in the integrated input spectra. This is the weight actually used by the synthesis, after clipping, etc.

FIXME: describe flags and weights.

  • Units: dimensionless
  • Shape: (Nl_obs, N_y, N_x)
integrated_chains_best_par

TODO: Add description of integrated_chains_best_par.

integrated_chains_ave_par

TODO: Add description of integrated_chains_ave_par.

integrated_chains_par

TODO: Add description of integrated_chains_par.

integrated_chains_best_LAx

TODO: Add description of integrated_chains_best_LAx.

integrated_chains_ave_LAx

TODO: Add description of integrated_chains_ave_LAx.

integrated_chains_LAx

TODO: Add description of integrated_chains_LAx.

integrated_chains_best_mu_cor

TODO: Add description of integrated_chains_best_mu_cor.

integrated_chains_ave_mu_cor

TODO: Add description of integrated_chains_ave_mu_cor.

integrated_chains_mu_cor

TODO: Add description of integrated_chains_mu_cor.

l_obs

Wavelength array for the spectral.

  • Units: \([\overset{\circ}{A}]\)
  • Shape: (Nl_obs)

Table Of Contents

Previous topic

Galaxy property list

Next topic

Utility functions

This Page