Commit 4edd9dd7 authored by Josué Ortega's avatar Josué Ortega

Import Upstream version 0.2.0

parents
0.2.0 (2016-09-29)
------------------
- Make sure that calling structure_at with an array, list, or tuple all
behave the same. [#98]
- Added support for linked scatter plots and multiple selections.
[#104, #105, #109, #136]
- Added support for custom functions to define what a 'neighbor' is. [#108]
- Fixed a bug that caused the interactive viewer when showing a dendrogram
loaded from a file. [#106, #110]
- Added a 'prune' method to prune dendrograms after computing them. [#111]
- Added support for brightness temperatures in Kelvin. [#112]
- Cache/memoize catalog statistics. [#115]
- Make sure that periodic boundaries (e.g. longitude) are properly
supported. [#121]
- Added progress bar for catalog computation. [#127]
- Better support for image WCS. [#126, #137]
- Improve the performance of dendrogram loading. [#131]
- Include dendrogram parameters in HDF5 files. [#142, #145]
- Give HDUs names in FITS output. [#144]
0.1.0 (2013-11-09)
------------------
Initial release
Computing Astronomical Dendrograms
Copyright (c) 2013 Thomas P. Robitaille, Chris Beaumont, Braden MacDonald,
and Erik Rosolowsky
Permission is hereby granted, free of charge, to any person obtaining a
copy of this software and associated documentation files (the "Software"),
to deal in the Software without restriction, including without limitation
the rights to use, copy, modify, merge, publish, distribute, sublicense,
and/or sell copies of the Software, and to permit persons to whom the
Software is furnished to do so, subject to the following conditions:
The above copyright notice and this permission notice shall be included in
all copies or substantial portions of the Software.
THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
DEALINGS IN THE SOFTWARE.
\ No newline at end of file
include LICENSE
include CHANGES
include runtests.py
recursive-include astrodendro/tests *.npz
recursive-include astrodendro/tests/benchmark_data *.fits
recursive-include docs *
Metadata-Version: 1.1
Name: astrodendro
Version: 0.2.0
Summary: Python package for computation of astronomical dendrograms
Home-page: http://www.dendrograms.org
Author: Thomas Robitaille, Chris Beaumont, Adam Ginsburg, Braden MacDonald, and Erik Rosolowsky
Author-email: thomas.robitaille@gmail.com
License: UNKNOWN
Description: UNKNOWN
Keywords: Scientific/Engineering
Platform: UNKNOWN
Classifier: Development Status :: 4 - Beta
Classifier: Programming Language :: Python
Classifier: License :: OSI Approved :: MIT License
Requires: numpy
Provides: astrodendro
Metadata-Version: 1.1
Name: astrodendro
Version: 0.2.0
Summary: Python package for computation of astronomical dendrograms
Home-page: http://www.dendrograms.org
Author: Thomas Robitaille, Chris Beaumont, Adam Ginsburg, Braden MacDonald, and Erik Rosolowsky
Author-email: thomas.robitaille@gmail.com
License: UNKNOWN
Description: UNKNOWN
Keywords: Scientific/Engineering
Platform: UNKNOWN
Classifier: Development Status :: 4 - Beta
Classifier: Programming Language :: Python
Classifier: License :: OSI Approved :: MIT License
Requires: numpy
Provides: astrodendro
CHANGES
LICENSE
MANIFEST.in
runtests.py
setup.py
astrodendro/__init__.py
astrodendro/analysis.py
astrodendro/dendrogram.py
astrodendro/flux.py
astrodendro/plot.py
astrodendro/progressbar.py
astrodendro/pruning.py
astrodendro/scatter.py
astrodendro/six.py
astrodendro/structure.py
astrodendro/structure_collection.py
astrodendro/viewer.py
astrodendro.egg-info/PKG-INFO
astrodendro.egg-info/SOURCES.txt
astrodendro.egg-info/dependency_links.txt
astrodendro.egg-info/top_level.txt
astrodendro/io/__init__.py
astrodendro/io/fits.py
astrodendro/io/handler.py
astrodendro/io/hdf5.py
astrodendro/io/util.py
astrodendro/tests/__init__.py
astrodendro/tests/_testdata.py
astrodendro/tests/benchmark.py
astrodendro/tests/build_benchmark.py
astrodendro/tests/sample-data-hl.npz
astrodendro/tests/test_analysis.py
astrodendro/tests/test_compute.py
astrodendro/tests/test_flux.py
astrodendro/tests/test_index.py
astrodendro/tests/test_io.py
astrodendro/tests/test_is_independent.py
astrodendro/tests/test_pruning.py
astrodendro/tests/test_recursion.py
astrodendro/tests/test_structure.py
astrodendro/tests/benchmark_data/2d.fits
astrodendro/tests/benchmark_data/2d1.fits
astrodendro/tests/benchmark_data/2d2.fits
astrodendro/tests/benchmark_data/2d3.fits
astrodendro/tests/benchmark_data/3d.fits
astrodendro/tests/benchmark_data/3d1.fits
astrodendro/tests/benchmark_data/3d2.fits
astrodendro/tests/benchmark_data/3d3.fits
docs/L1551_scuba_850mu.fits
docs/Makefile
docs/PerA_Extn2MASS_F_Gal.fits
docs/advanced.rst
docs/algorithm.rst
docs/catalog.rst
docs/conf.py
docs/index.rst
docs/installing.rst
docs/migration.rst
docs/plotting.rst
docs/rtd-pip-requirements
docs/scatter_screenshot.png
docs/scatter_selected_viewer_screenshot.png
docs/schematic_structure_1.png
docs/schematic_structure_2.png
docs/schematic_tree.png
docs/using.rst
docs/viewer_screenshot.png
docs/wcsaxes_docs_screenshot.png
docs/_templates/autosummary/base.rst
docs/_templates/autosummary/class.rst
docs/_templates/autosummary/module.rst
docs/algorithm/large_delta_final.png
docs/algorithm/large_delta_step1.png
docs/algorithm/large_delta_step2.png
docs/algorithm/min_value_final.png
docs/algorithm/simple_final.png
docs/algorithm/simple_step1.png
docs/algorithm/simple_step2.png
docs/algorithm/simple_step3.png
docs/algorithm/simple_step4.png
docs/algorithm/small_delta_final.png
docs/algorithm/small_delta_step1.png
docs/algorithm/small_delta_step2.png
docs/algorithm/small_delta_step3.png
docs/api/astrodendro.analysis.rst
docs/api/astrodendro.dendrogram.Dendrogram.rst
docs/api/astrodendro.dendrogram.periodic_neighbours.rst
docs/api/astrodendro.plot.DendrogramPlotter.rst
docs/api/astrodendro.pruning.rst
docs/api/astrodendro.structure.Structure.rst
docs/diagrams/README.md
docs/diagrams/dendrograms.svg
\ No newline at end of file
# Licensed under an MIT open source license - see LICENSE
from .dendrogram import Dendrogram, periodic_neighbours
from .structure import Structure
from .analysis import ppv_catalog, pp_catalog
from .plot import DendrogramPlotter
from .viewer import BasicDendrogramViewer
__version__ = '0.2.0'
This diff is collapsed.
This diff is collapsed.
import warnings
import numpy as np
from astropy import units as u
from astropy.constants import si
class UnitMetadataWarning(UserWarning):
pass
def quantity_sum(quantities):
"""
In Astropy 0.3, np.sum will do the right thing for quantities, but in the mean time we need a workaround.
"""
return np.sum(quantities.value) * quantities.unit
def compute_flux(input_quantities, output_unit, wavelength=None, spatial_scale=None,
velocity_scale=None, beam_major=None, beam_minor=None):
"""
Given a set of flux values in arbitrary units, find the total flux in a
specific set of units.
Parameters
----------
input_quantities : :class:`~astropy.units.quantity.Quantity` instance
A `astropy.units.quantity.Quantity` instance containing an array of
flux values to be summed.
output_unit : :class:`~astropy.units.core.Unit` instance
The final unit to give the total flux in (should be equivalent to Jy)
wavelength : :class:`~astropy.units.quantity.Quantity` instance
The wavelength of the data (required if converting e.g.
ergs/cm^2/s/micron to Jy)
spatial_scale : :class:`~astropy.units.quantity.Quantity` instance
The pixel scale of the data (should be an angle)
velocity_scale : :class:`~astropy.units.quantity.Quantity` instance
The pixel scale of the data (should be a velocity)
beam_major : :class:`~astropy.units.quantity.Quantity` instance
The beam major full width at half_maximum (FWHM)
beam_minor : :class:`~astropy.units.quantity.Quantity` instance
The beam minor full width at half_maximum (FWHM)
"""
# Start off by finding the total flux in Jy
if input_quantities.unit.is_equivalent(u.Jy): # Fnu
# Simply sum up the values and convert to output unit
total_flux = quantity_sum(input_quantities).to(u.Jy)
elif input_quantities.unit.is_equivalent(u.erg / u.cm ** 2 / u.s / u.m): # Flambda
if wavelength is not None and not wavelength.unit.is_equivalent(u.m):
raise ValueError("wavelength should be a physical length")
# Find the frequency
if wavelength is None:
raise ValueError("wavelength is needed to convert from {0} to Jy".format(input_quantities.unit))
# Find frequency
nu = si.c / wavelength
# Convert input quantity to Fnu in Jy
q = (input_quantities * wavelength / nu).to(u.Jy)
# Find total flux in Jy
total_flux = quantity_sum(q)
elif input_quantities.unit.is_equivalent(u.MJy / u.sr): # surface brightness (Fnu)
if spatial_scale is not None and not spatial_scale.unit.is_equivalent(u.degree):
raise ValueError("spatial_scale should be an angle")
if spatial_scale is None:
raise ValueError("spatial_scale is needed to convert from {0} to Jy".format(input_quantities.unit))
# Find the area of a pixel as a solid angle
pixel_area = (spatial_scale ** 2)
# Convert input quantity to Fnu in Jy
q = (input_quantities * pixel_area).to(u.Jy)
# Find total flux in Jy
total_flux = quantity_sum(q)
elif input_quantities.unit.is_equivalent(u.Jy / u.beam):
if spatial_scale is not None and not spatial_scale.unit.is_equivalent(u.degree):
raise ValueError("spatial_scale should be an angle")
if spatial_scale is None:
raise ValueError("spatial_scale is needed to convert from {0} to Jy".format(input_quantities.unit))
if beam_major is not None and not beam_major.unit.is_equivalent(u.degree):
raise ValueError("beam_major should be an angle")
if beam_major is None:
raise ValueError("beam_major is needed to convert from {0} to Jy".format(input_quantities.unit))
if beam_minor is not None and not beam_minor.unit.is_equivalent(u.degree):
raise ValueError("beam_minor should be an angle")
if beam_minor is None:
raise ValueError("beam_minor is needed to convert from {0} to Jy".format(input_quantities.unit))
# Find the beam area
beams_per_pixel = spatial_scale ** 2 / (beam_minor * beam_major * 1.1331) * u.beam
# Convert input quantity to Fnu in Jy
q = (input_quantities * beams_per_pixel).to(u.Jy)
# Find total flux in Jy
total_flux = quantity_sum(q)
elif input_quantities.unit.is_equivalent(u.K):
if spatial_scale is not None and not spatial_scale.unit.is_equivalent(u.degree):
raise ValueError("spatial_scale should be an angle")
if spatial_scale is None:
raise ValueError("spatial_scale is needed to convert from {0} to Jy".format(input_quantities.unit))
if beam_major is not None and not beam_major.unit.is_equivalent(u.degree):
raise ValueError("beam_major should be an angle")
if beam_major is None:
raise ValueError("beam_major is needed to convert from {0} to Jy".format(input_quantities.unit))
if beam_minor is not None and not beam_minor.unit.is_equivalent(u.degree):
raise ValueError("beam_minor should be an angle")
if beam_minor is None:
raise ValueError("beam_minor is needed to convert from {0} to Jy".format(input_quantities.unit))
if wavelength is not None and not wavelength.unit.is_equivalent(u.m, equivalencies=u.spectral()):
raise ValueError("wavelength should be a physical length")
# Find the frequency
if wavelength is None:
raise ValueError("wavelength is needed to convert from {0} to Jy".format(input_quantities.unit))
warnings.warn("'Kelvin' units interpreted as main beam brightness temperature.",
UnitMetadataWarning)
# Find frequency
nu = wavelength.to(u.Hz, equivalencies=u.spectral())
# Angular area of beam. Conversion between 2D Gaussian FWHM and effective area comes from https://github.com/radio-astro-tools/radio_beam/blob/bc906c38a65e85c6a894ee81519a642665e50f7c/radio_beam/beam.py#L8
omega_beam = np.pi * 2 / (8*np.log(2)) * beam_major * beam_minor
# Find the beam area
beams_per_pixel = spatial_scale ** 2 / omega_beam * u.beam
# Convert input quantity to Fnu in Jy
# Implicitly, this equivalency gives the Janskys in a single beam, so we make this explicit by dividing out a beam
jansky_per_beam = input_quantities.to(u.Jy,
equivalencies=u.brightness_temperature(omega_beam, nu)) / u.beam
q = jansky_per_beam * beams_per_pixel
# Find total flux in Jy
total_flux = quantity_sum(q)
else:
raise ValueError("Flux units {0} not yet supported".format(input_quantities.unit))
if not output_unit.is_equivalent(u.Jy):
raise ValueError("output_unit has to be equivalent to Jy")
else:
return total_flux.to(output_unit)
# Licensed under an MIT open source license - see LICENSE
from .fits import FITSHandler
from .hdf5 import HDF5Handler
IO_FORMATS = {
'fits': FITSHandler,
'hdf5': HDF5Handler
}
def _valid(formats):
return " and ".join(["'{0}'".format(key) for key in formats])
def load_dendrogram(filename, format=None):
if format is not None:
return IO_FORMATS[format].import_dendro(filename)
else:
for handler in IO_FORMATS.values():
if handler.identify(filename, mode='r'):
return handler.import_dendro(filename)
raise IOError("Could not automatically identify file format - use the "
"format= option to specify which format to use (valid "
"options are {0})".format(_valid(IO_FORMATS)))
def save_dendrogram(dendrogram, filename, format=None):
if format is not None:
return IO_FORMATS[format].export_dendro(dendrogram, filename)
else:
for handler in IO_FORMATS.values():
if handler.identify(filename, mode='w'):
return handler.export_dendro(dendrogram, filename)
raise IOError("Could not automatically identify file format - use the "
"format= option to specify which format to use (valid "
"options are {0})".format(_valid(IO_FORMATS)))
# Licensed under an MIT open source license - see LICENSE
import os
import numpy as np
from .util import parse_dendrogram
from .handler import IOHandler
# Import and export
# FITS file signature as per RFC 4047
FITS_SIGNATURE = (b"\x53\x49\x4d\x50\x4c\x45\x20\x20\x3d\x20\x20\x20\x20\x20"
b"\x20\x20\x20\x20\x20\x20\x20\x20\x20\x20\x20\x20\x20\x20"
b"\x20\x54")
def is_fits(filename, mode='r'):
if mode == 'r' and os.path.exists(filename):
fileobj = open(filename, 'rb')
sig = fileobj.read(30)
return sig == FITS_SIGNATURE
elif filename.lower().endswith(('.fits', '.fits.gz', '.fit', '.fit.gz')):
return True
else:
return False
def dendro_export_fits(d, filename):
"""Export the dendrogram 'd' to the FITS file 'filename'"""
from astropy.io import fits
try:
primary_hdu = fits.PrimaryHDU(header=d.wcs.to_header())
except AttributeError:
primary_hdu = fits.PrimaryHDU()
primary_hdu.header["MIN_NPIX"] = (d.params['min_npix'],
"Minimum number of pixels in a leaf.")
primary_hdu.header["MIN_DELT"] = (d.params['min_delta'],
"Minimum branch height.")
primary_hdu.header["MIN_VAL"] = (d.params['min_value'],
"Minimum intensity value.")
hdus = [primary_hdu,
fits.ImageHDU(d.data, name='data'),
fits.ImageHDU(d.index_map, name='index_map'),
fits.ImageHDU(np.array([ord(x) for x in d.to_newick()]), name='newick')]
hdulist = fits.HDUList(hdus)
hdulist.writeto(filename, clobber=True)
def dendro_import_fits(filename):
"""Import 'filename' and construct a dendrogram from it"""
from astropy.io import fits
from astropy.wcs.wcs import WCS
with fits.open(filename) as hdus:
try:
wcs = WCS(hdus[0].header)
except AttributeError:
wcs = None
data = hdus[1].data
index_map = hdus[2].data
newick = ''.join(chr(x) for x in hdus[3].data.flat)
if 'MIN_NPIX' in hdus[0].header:
params = {"min_npix": hdus[0].header['MIN_NPIX'],
"min_value": hdus[0].header['MIN_VAL'],
"min_delta": hdus[0].header['MIN_DELT']}
else:
params = {}
return parse_dendrogram(newick, data, index_map, params, wcs)
FITSHandler = IOHandler(identify=is_fits,
export_dendro=dendro_export_fits,
import_dendro=dendro_import_fits)
from collections import namedtuple
IOHandler = namedtuple("IOHandler", "identify import_dendro export_dendro")
\ No newline at end of file
# Licensed under an MIT open source license - see LICENSE
import os
import numpy as np
from astropy import log
from .util import parse_dendrogram
from .handler import IOHandler
HDF5_SIGNATURE = b'\x89HDF\r\n\x1a\n'
def is_hdf5(filename, mode='r'):
if mode == 'r' and os.path.exists(filename):
fileobj = open(filename, 'rb')
sig = fileobj.read(8)
return sig == HDF5_SIGNATURE
elif filename.lower().endswith(('.hdf5', '.h5')):
return True
else:
return False
def dendro_export_hdf5(d, filename):
"""Export the dendrogram 'd' to the HDF5 file 'filename'"""
import h5py
f = h5py.File(filename, 'w')
f.attrs['n_dim'] = d.n_dim
f.create_dataset('newick', data=d.to_newick())
ds = f.create_dataset('index_map', data=d.index_map, compression=True)
ds.attrs['CLASS'] = 'IMAGE'
ds.attrs['IMAGE_VERSION'] = '1.2'
ds.attrs['IMAGE_MINMAXRANGE'] = [d.index_map.min(), d.index_map.max()]
ds = f.create_dataset('data', data=d.data, compression=True)
ds.attrs['CLASS'] = 'IMAGE'
ds.attrs['IMAGE_VERSION'] = '1.2'
ds.attrs['IMAGE_MINMAXRANGE'] = [d.data.min(), d.data.max()]
for key in d.params.keys():
f.attrs[key] = d.params[key]
try:
f.create_dataset('wcs_header', data=d.wcs.to_header_string())
except AttributeError:
pass
f.close()
def dendro_import_hdf5(filename):
"""Import 'filename' and construct a dendrogram from it"""
import h5py
from astropy.wcs.wcs import WCS
log.debug('Loading HDF5 file from disk...')
with h5py.File(filename, 'r') as h5f:
newick = h5f['newick'].value
data = h5f['data'].value
index_map = h5f['index_map'].value
params = {}
if 'min_value' in h5f.attrs:
params['min_value'] = h5f.attrs['min_value']
params['min_delta'] = h5f.attrs['min_delta']
params['min_npix'] = h5f.attrs['min_npix']
try:
wcs = WCS(h5f['wcs_header'].value)
except KeyError:
wcs = None
log.debug('Parsing dendrogram...')
return parse_dendrogram(newick, data, index_map, params, wcs)
HDF5Handler = IOHandler(identify=is_hdf5,
export_dendro=dendro_export_hdf5,
import_dendro=dendro_import_hdf5)
import numpy as np
from .. import six
from astropy.utils.console import ProgressBar
from astropy import log
def parse_newick(string):
items = {}
# Find maximum level
current_level = 0
max_level = 0
log.debug('String loading...')
for i, c in enumerate(string):
if c == '(':
current_level += 1
if c == ')':
current_level -= 1
max_level = max(max_level, current_level)
# Loop through levels and construct tree
log.debug('Tree loading...')
for level in range(max_level, 0, -1):
pairs = []
current_level = 0
for i, c in enumerate(string):
if c == '(':
current_level += 1
if current_level == level:
start = i
if c == ')':
if current_level == level:
pairs.append((start, i))
current_level -= 1
for pair in pairs[::-1]:
# Extract start and end of branch definition
start, end = pair
# Find the ID of the branch
colon = string.find(":", end)
branch_id = string[end + 1:colon]
if branch_id == '':
branch_id = 'trunk'
else:
branch_id = int(branch_id)
# Add branch definition to overall definition
items[branch_id] = eval("{%s}" % string[start + 1:end])
# Remove branch definition from string
string = string[:start] + string[end + 1:]
new_items = {}
def collect(d):
for item in d:
if item in items:
collect(items[item])
d[item] = (items[item], d[item])
return
collect(items['trunk'])
return items['trunk']
def parse_dendrogram(newick, data, index_map, params, wcs=None):
from ..dendrogram import Dendrogram
from ..structure import Structure
d = Dendrogram()
d.ndim = len(data.shape)
d._structures_dict = {}
d.data = data
d.index_map = index_map
d.params = params
d.wcs = wcs
try:
flux_by_structure, indices_by_structure = _fast_reader(d.index_map, data)
except ImportError:
flux_by_structure, indices_by_structure = _slow_reader(d.index_map, data)
def _construct_tree(repr):
structures = []
for idx in repr:
idx = int(idx)
structure_indices = indices_by_structure[idx]
f = flux_by_structure[idx]
if type(repr[idx]) == tuple:
sub_structures_repr = repr[idx][0] # Parsed representation of sub structures
sub_structures = _construct_tree(sub_structures_repr)
for i in sub_structures:
d._structures_dict[i.idx] = i
b = Structure(structure_indices, f, children=sub_structures, idx=idx, dendrogram=d)
# Correct merge levels - complicated because of the
# order in which we are building the tree.
# What we do is look at the heights of this branch's