[yt-svn] commit/yt: 11 new changesets
commits-noreply at bitbucket.org
commits-noreply at bitbucket.org
Thu Dec 8 07:29:58 PST 2016
11 new commits in yt:
https://bitbucket.org/yt_analysis/yt/commits/5603bf81bd76/
Changeset: 5603bf81bd76
Branch: yt
User: jzuhone
Date: 2016-12-05 17:50:36+00:00
Summary: Move fits_image to visualization
Affected #: 4 files
diff -r 2668f56ec734f8600cfd4c72b06cc035bdd42d3a -r 5603bf81bd7651f2e976755d3d5ade26efc89792 yt/utilities/fits_image.py
--- a/yt/utilities/fits_image.py
+++ /dev/null
@@ -1,793 +0,0 @@
-"""
-FITSImageData Class
-"""
-
-#-----------------------------------------------------------------------------
-# Copyright (c) 2013, yt Development Team.
-#
-# Distributed under the terms of the Modified BSD License.
-#
-# The full license is in the file COPYING.txt, distributed with this software.
-#-----------------------------------------------------------------------------
-from yt.extern.six import string_types
-import numpy as np
-from yt.funcs import mylog, iterable, fix_axis, ensure_list
-from yt.visualization.fixed_resolution import FixedResolutionBuffer
-from yt.data_objects.construction_data_containers import YTCoveringGrid
-from yt.utilities.on_demand_imports import _astropy
-from yt.units.yt_array import YTQuantity, YTArray
-from yt.units import dimensions
-from yt.utilities.parallel_tools.parallel_analysis_interface import \
- parallel_root_only
-from yt.visualization.volume_rendering.off_axis_projection import off_axis_projection
-import re
-
-class FITSImageData(object):
-
- def __init__(self, data, fields=None, units=None, width=None, wcs=None):
- r""" Initialize a FITSImageData object.
-
- FITSImageData contains a collection of FITS ImageHDU instances and
- WCS information, along with units for each of the images. FITSImageData
- instances can be constructed from ImageArrays, NumPy arrays, dicts
- of such arrays, FixedResolutionBuffers, and YTCoveringGrids. The latter
- two are the most powerful because WCS information can be constructed
- automatically from their coordinates.
-
- Parameters
- ----------
- data : FixedResolutionBuffer or a YTCoveringGrid. Or, an
- ImageArray, an numpy.ndarray, or dict of such arrays
- The data to be made into a FITS image or images.
- fields : single string or list of strings, optional
- The field names for the data. If *fields* is none and *data* has
- keys, it will use these for the fields. If *data* is just a
- single array one field name must be specified.
- units : string
- The units of the WCS coordinates. Defaults to "cm".
- width : float or YTQuantity
- The width of the image. Either a single value or iterable of values.
- If a float, assumed to be in *units*. Only used if this information
- is not already provided by *data*.
- wcs : `astropy.wcs.WCS` instance, optional
- Supply an AstroPy WCS instance. Will override automatic WCS
- creation from FixedResolutionBuffers and YTCoveringGrids.
-
- Examples
- --------
-
- >>> # This example uses a FRB.
- >>> ds = load("sloshing_nomag2_hdf5_plt_cnt_0150")
- >>> prj = ds.proj(2, "kT", weight_field="density")
- >>> frb = prj.to_frb((0.5, "Mpc"), 800)
- >>> # This example just uses the FRB and puts the coords in kpc.
- >>> f_kpc = FITSImageData(frb, fields="kT", units="kpc")
- >>> # This example specifies a specific WCS.
- >>> from astropy.wcs import WCS
- >>> w = WCS(naxis=self.dimensionality)
- >>> w.wcs.crval = [30., 45.] # RA, Dec in degrees
- >>> w.wcs.cunit = ["deg"]*2
- >>> nx, ny = 800, 800
- >>> w.wcs.crpix = [0.5*(nx+1), 0.5*(ny+1)]
- >>> w.wcs.ctype = ["RA---TAN","DEC--TAN"]
- >>> scale = 1./3600. # One arcsec per pixel
- >>> w.wcs.cdelt = [-scale, scale]
- >>> f_deg = FITSImageData(frb, fields="kT", wcs=w)
- >>> f_deg.writeto("temp.fits")
- """
-
- if units is None:
- units = "cm"
- if width is None:
- width = 1.0
-
- exclude_fields = ['x','y','z','px','py','pz',
- 'pdx','pdy','pdz','weight_field']
-
- self.hdulist = _astropy.pyfits.HDUList()
-
- if isinstance(fields, string_types):
- fields = [fields]
-
- if hasattr(data, 'keys'):
- img_data = data
- if fields is None:
- fields = list(img_data.keys())
- elif isinstance(data, np.ndarray):
- if fields is None:
- mylog.warning("No field name given for this array. Calling it 'image_data'.")
- fn = 'image_data'
- fields = [fn]
- else:
- fn = fields[0]
- img_data = {fn: data}
-
- self.fields = []
- for fd in fields:
- if isinstance(fd, tuple):
- self.fields.append(fd[1])
- else:
- self.fields.append(fd)
-
- first = True
- self.field_units = {}
- for key in fields:
- if key not in exclude_fields:
- if hasattr(img_data[key], "units"):
- self.field_units[key] = img_data[key].units
- else:
- self.field_units[key] = "dimensionless"
- mylog.info("Making a FITS image of field %s" % key)
- if first:
- hdu = _astropy.pyfits.PrimaryHDU(np.array(img_data[key]))
- first = False
- else:
- hdu = _astropy.pyfits.ImageHDU(np.array(img_data[key]))
- hdu.name = key
- hdu.header["btype"] = key
- if hasattr(img_data[key], "units"):
- hdu.header["bunit"] = re.sub('()', '', str(img_data[key].units))
- self.hdulist.append(hdu)
-
- self.shape = self.hdulist[0].shape
- self.dimensionality = len(self.shape)
-
- if wcs is None:
- w = _astropy.pywcs.WCS(header=self.hdulist[0].header, naxis=self.dimensionality)
- if isinstance(img_data, FixedResolutionBuffer):
- # FRBs are a special case where we have coordinate
- # information, so we take advantage of this and
- # construct the WCS object
- dx = (img_data.bounds[1]-img_data.bounds[0]).in_units(units).v/self.shape[0]
- dy = (img_data.bounds[3]-img_data.bounds[2]).in_units(units).v/self.shape[1]
- xctr = 0.5*(img_data.bounds[1]+img_data.bounds[0]).in_units(units).v
- yctr = 0.5*(img_data.bounds[3]+img_data.bounds[2]).in_units(units).v
- center = [xctr, yctr]
- cdelt = [dx,dy]
- elif isinstance(img_data, YTCoveringGrid):
- cdelt = img_data.dds.in_units(units).v
- center = 0.5*(img_data.left_edge+img_data.right_edge).in_units(units).v
- else:
- # If img_data is just an array, we assume the center is the origin
- # and use the image width to determine the cell widths
- if not iterable(width):
- width = [width]*self.dimensionality
- if isinstance(width[0], YTQuantity):
- cdelt = [wh.in_units(units).v/n for wh, n in zip(width, self.shape)]
- else:
- cdelt = [float(wh)/n for wh, n in zip(width, self.shape)]
- center = [0.0]*self.dimensionality
- w.wcs.crpix = 0.5*(np.array(self.shape)+1)
- w.wcs.crval = center
- w.wcs.cdelt = cdelt
- w.wcs.ctype = ["linear"]*self.dimensionality
- w.wcs.cunit = [units]*self.dimensionality
- self.set_wcs(w)
- else:
- self.set_wcs(wcs)
-
- def set_wcs(self, wcs, wcsname=None, suffix=None):
- """
- Set the WCS coordinate information for all images
- with a WCS object *wcs*.
- """
- if wcsname is None:
- wcs.wcs.name = 'yt'
- else:
- wcs.wcs.name = wcsname
- h = wcs.to_header()
- if suffix is None:
- self.wcs = wcs
- else:
- setattr(self, "wcs"+suffix, wcs)
- for img in self.hdulist:
- for k, v in h.items():
- kk = k
- if suffix is not None:
- kk += suffix
- img.header[kk] = v
-
- def update_header(self, field, key, value):
- """
- Update the FITS header for *field* with a
- *key*, *value* pair. If *field* == "all", all
- headers will be updated.
- """
- if field == "all":
- for img in self.hdulist:
- img.header[key] = value
- else:
- if field not in self.keys():
- raise KeyError("%s not an image!" % field)
- idx = self.fields.index(field)
- self.hdulist[idx].header[key] = value
-
- def update_all_headers(self, key, value):
- mylog.warning("update_all_headers is deprecated. "+
- "Use update_header('all', key, value) instead.")
- self.update_header("all", key, value)
-
- def keys(self):
- return self.fields
-
- def has_key(self, key):
- return key in self.fields
-
- def values(self):
- return [self.hdulist[k] for k in self.fields]
-
- def items(self):
- return [(k, self.hdulist[k]) for k in self.fields]
-
- def __getitem__(self, item):
- return self.hdulist[item]
-
- def info(self):
- return self.hdulist.info()
-
- @parallel_root_only
- def writeto(self, fileobj, fields=None, clobber=False, **kwargs):
- r"""
- Write all of the fields or a subset of them to a FITS file.
-
- Parameters
- ----------
- fileobj : string
- The name of the file to write to.
- fields : list of strings, optional
- The fields to write to the file. If not specified
- all of the fields in the buffer will be written.
- clobber : boolean, optional
- Whether or not to overwrite a previously existing file.
- Default: False
- All other keyword arguments are passed to the `writeto`
- method of `astropy.io.fits.HDUList`.
- """
- if fields is None:
- hdus = self.hdulist
- else:
- hdus = _astropy.pyfits.HDUList()
- for field in fields:
- hdus.append(self.hdulist[field])
- hdus.writeto(fileobj, clobber=clobber, **kwargs)
-
- def to_glue(self, label="yt", data_collection=None):
- """
- Takes the data in the FITSImageData instance and exports it to
- Glue (http://www.glueviz.org) for interactive analysis. Optionally
- add a *label*. If you are already within the Glue environment, you
- can pass a *data_collection* object, otherwise Glue will be started.
- """
- from glue.core import DataCollection, Data
- from glue.core.coordinates import coordinates_from_header
- from glue.qt.glue_application import GlueApplication
-
- image = Data(label=label)
- image.coords = coordinates_from_header(self.wcs.to_header())
- for k,f in self.hdulist.items():
- image.add_component(f.data, k)
- if data_collection is None:
- dc = DataCollection([image])
- app = GlueApplication(dc)
- app.start()
- else:
- data_collection.append(image)
-
- def to_aplpy(self, **kwargs):
- """
- Use APLpy (http://aplpy.github.io) for plotting. Returns an
- `aplpy.FITSFigure` instance. All keyword arguments are passed to the
- `aplpy.FITSFigure` constructor.
- """
- import aplpy
- return aplpy.FITSFigure(self.hdulist, **kwargs)
-
- def get_data(self, field):
- """
- Return the data array of the image corresponding to *field*
- with units attached.
- """
- return YTArray(self.hdulist[field].data, self.field_units[field])
-
- def set_unit(self, field, units):
- """
- Set the units of *field* to *units*.
- """
- if field not in self.keys():
- raise KeyError("%s not an image!" % field)
- idx = self.fields.index(field)
- new_data = YTArray(self.hdulist[idx].data, self.field_units[field]).in_units(units)
- self.hdulist[idx].data = new_data.v
- self.hdulist[idx].header["bunit"] = units
- self.field_units[field] = units
-
- def pop(self, key):
- """
- Remove a field with name *key*
- and return it as a new FITSImageData
- instance.
- """
- if key not in self.keys():
- raise KeyError("%s not an image!" % key)
- idx = self.fields.index(key)
- im = self.hdulist.pop(idx)
- data = YTArray(im.data, self.field_units[key])
- self.field_units.pop(key)
- self.fields.remove(key)
- return FITSImageData(data, fields=key, wcs=self.wcs)
-
- @classmethod
- def from_file(cls, filename):
- """
- Generate a FITSImageData instance from one previously written to
- disk.
-
- Parameters
- ----------
- filename : string
- The name of the file to open.
- """
- f = _astropy.pyfits.open(filename)
- data = {}
- for hdu in f:
- data[hdu.header["btype"]] = YTArray(hdu.data, hdu.header["bunit"])
- f.close()
- return cls(data, wcs=_astropy.pywcs.WCS(header=hdu.header))
-
- @classmethod
- def from_images(cls, image_list):
- """
- Generate a new FITSImageData instance from a list of FITSImageData
- instances.
-
- Parameters
- ----------
- image_list : list of FITSImageData instances
- The images to be combined.
- """
- w = image_list[0].wcs
- img_shape = image_list[0].shape
- data = {}
- for image in image_list:
- assert_same_wcs(w, image.wcs)
- if img_shape != image.shape:
- raise RuntimeError("Images do not have the same shape!")
- for key in image.keys():
- data[key] = image.get_data(key)
- return cls(data, wcs=w)
-
- def create_sky_wcs(self, sky_center, sky_scale,
- ctype=["RA---TAN","DEC--TAN"],
- crota=None, cd=None, pc=None,
- wcsname="celestial",
- replace_old_wcs=True):
- """
- Takes a Cartesian WCS and converts it to one in a
- celestial coordinate system.
-
- Parameters
- ----------
- sky_center : iterable of floats
- Reference coordinates of the WCS in degrees.
- sky_scale : tuple or YTQuantity
- Conversion between an angle unit and a length unit,
- e.g. (3.0, "arcsec/kpc")
- ctype : list of strings, optional
- The type of the coordinate system to create.
- crota : 2-element ndarray, optional
- Rotation angles between cartesian coordinates and
- the celestial coordinates.
- cd : 2x2-element ndarray, optional
- Dimensioned coordinate transformation matrix.
- pc : 2x2-element ndarray, optional
- Coordinate transformation matrix.
- replace_old_wcs : boolean, optional
- Whether or not to overwrite the default WCS of the
- FITSImageData instance. If false, a second WCS will
- be added to the header. Default: True.
- """
- old_wcs = self.wcs
- naxis = old_wcs.naxis
- crval = [sky_center[0], sky_center[1]]
- if isinstance(sky_scale, YTQuantity):
- scaleq = sky_scale
- else:
- scaleq = YTQuantity(sky_scale[0],sky_scale[1])
- if scaleq.units.dimensions != dimensions.angle/dimensions.length:
- raise RuntimeError("sky_scale %s not in correct dimensions of angle/length!" % sky_scale)
- deltas = old_wcs.wcs.cdelt
- units = [str(unit) for unit in old_wcs.wcs.cunit]
- new_dx = (YTQuantity(-deltas[0], units[0])*scaleq).in_units("deg")
- new_dy = (YTQuantity(deltas[1], units[1])*scaleq).in_units("deg")
- new_wcs = _astropy.pywcs.WCS(naxis=naxis)
- cdelt = [new_dx.v, new_dy.v]
- cunit = ["deg"]*2
- if naxis == 3:
- crval.append(old_wcs.wcs.crval[2])
- cdelt.append(old_wcs.wcs.cdelt[2])
- ctype.append(old_wcs.wcs.ctype[2])
- cunit.append(old_wcs.wcs.cunit[2])
- new_wcs.wcs.crpix = old_wcs.wcs.crpix
- new_wcs.wcs.cdelt = cdelt
- new_wcs.wcs.crval = crval
- new_wcs.wcs.cunit = cunit
- new_wcs.wcs.ctype = ctype
- if crota is not None:
- new_wcs.wcs.crota = crota
- if cd is not None:
- new_wcs.wcs.cd = cd
- if pc is not None:
- new_wcs.wcs.cd = pc
- if replace_old_wcs:
- self.set_wcs(new_wcs, wcsname=wcsname)
- else:
- self.set_wcs(new_wcs, wcsname=wcsname, suffix="a")
-
-class FITSImageBuffer(FITSImageData):
- pass
-
-def sanitize_fits_unit(unit):
- if unit == "Mpc":
- mylog.info("Changing FITS file unit to kpc.")
- unit = "kpc"
- elif unit == "au":
- unit = "AU"
- return unit
-
-axis_wcs = [[1,2],[0,2],[0,1]]
-
-def construct_image(ds, axis, data_source, center, width=None, image_res=None):
- if width is None:
- width = ds.domain_width[axis_wcs[axis]]
- unit = ds.get_smallest_appropriate_unit(width[0])
- mylog.info("Making an image of the entire domain, "+
- "so setting the center to the domain center.")
- else:
- width = ds.coordinates.sanitize_width(axis, width, None)
- unit = str(width[0].units)
- if image_res is None:
- ddims = ds.domain_dimensions*ds.refine_by**ds.index.max_level
- if iterable(axis):
- nx = ddims.max()
- ny = ddims.max()
- else:
- nx, ny = [ddims[idx] for idx in axis_wcs[axis]]
- else:
- if iterable(image_res):
- nx, ny = image_res
- else:
- nx, ny = image_res, image_res
- dx = width[0]/nx
- crpix = [0.5*(nx+1), 0.5*(ny+1)]
- if hasattr(ds, "wcs") and not iterable(axis):
- # This is a FITS dataset, so we use it to construct the WCS
- cunit = [str(ds.wcs.wcs.cunit[idx]) for idx in axis_wcs[axis]]
- ctype = [ds.wcs.wcs.ctype[idx] for idx in axis_wcs[axis]]
- cdelt = [ds.wcs.wcs.cdelt[idx] for idx in axis_wcs[axis]]
- ctr_pix = center.in_units("code_length")[:ds.dimensionality].v
- crval = ds.wcs.wcs_pix2world(ctr_pix.reshape(1, ds.dimensionality))[0]
- crval = [crval[idx] for idx in axis_wcs[axis]]
- else:
- if unit == "unitary":
- unit = ds.get_smallest_appropriate_unit(ds.domain_width.max())
- elif unit == "code_length":
- unit = ds.get_smallest_appropriate_unit(ds.quan(1.0,"code_length"))
- unit = sanitize_fits_unit(unit)
- cunit = [unit]*2
- ctype = ["LINEAR"]*2
- cdelt = [dx.in_units(unit)]*2
- if iterable(axis):
- crval = center.in_units(unit)
- else:
- crval = [center[idx].in_units(unit) for idx in axis_wcs[axis]]
- if hasattr(data_source, 'to_frb'):
- if iterable(axis):
- frb = data_source.to_frb(width[0], (nx, ny), height=width[1])
- else:
- frb = data_source.to_frb(width[0], (nx, ny), center=center, height=width[1])
- else:
- frb = None
- w = _astropy.pywcs.WCS(naxis=2)
- w.wcs.crpix = crpix
- w.wcs.cdelt = cdelt
- w.wcs.crval = crval
- w.wcs.cunit = cunit
- w.wcs.ctype = ctype
- return w, frb
-
-def assert_same_wcs(wcs1, wcs2):
- from numpy.testing import assert_allclose
- assert wcs1.naxis == wcs2.naxis
- for i in range(wcs1.naxis):
- assert wcs1.wcs.cunit[i] == wcs2.wcs.cunit[i]
- assert wcs1.wcs.ctype[i] == wcs2.wcs.ctype[i]
- assert_allclose(wcs1.wcs.crpix, wcs2.wcs.crpix)
- assert_allclose(wcs1.wcs.cdelt, wcs2.wcs.cdelt)
- assert_allclose(wcs1.wcs.crval, wcs2.wcs.crval)
- crota1 = getattr(wcs1.wcs, "crota", None)
- crota2 = getattr(wcs2.wcs, "crota", None)
- if crota1 is None or crota2 is None:
- assert crota1 == crota2
- else:
- assert_allclose(wcs1.wcs.crota, wcs2.wcs.crota)
- cd1 = getattr(wcs1.wcs, "cd", None)
- cd2 = getattr(wcs2.wcs, "cd", None)
- if cd1 is None or cd2 is None:
- assert cd1 == cd2
- else:
- assert_allclose(wcs1.wcs.cd, wcs2.wcs.cd)
- pc1 = getattr(wcs1.wcs, "pc", None)
- pc2 = getattr(wcs2.wcs, "pc", None)
- if pc1 is None or pc2 is None:
- assert pc1 == pc2
- else:
- assert_allclose(wcs1.wcs.pc, wcs2.wcs.pc)
-
-class FITSSlice(FITSImageData):
- r"""
- Generate a FITSImageData of an on-axis slice.
-
- Parameters
- ----------
- ds : :class:`yt.data_objects.api.Dataset`
- The dataset object.
- axis : character or integer
- The axis of the slice. One of "x","y","z", or 0,1,2.
- fields : string or list of strings
- The fields to slice
- center : A sequence of floats, a string, or a tuple.
- The coordinate of the center of the image. If set to 'c', 'center' or
- left blank, the plot is centered on the middle of the domain. If set to
- 'max' or 'm', the center will be located at the maximum of the
- ('gas', 'density') field. Centering on the max or min of a specific
- field is supported by providing a tuple such as ("min","temperature") or
- ("max","dark_matter_density"). Units can be specified by passing in *center*
- as a tuple containing a coordinate and string unit name or by passing
- in a YTArray. If a list or unitless array is supplied, code units are
- assumed.
- width : tuple or a float.
- Width can have four different formats to support windows with variable
- x and y widths. They are:
-
- ================================== =======================
- format example
- ================================== =======================
- (float, string) (10,'kpc')
- ((float, string), (float, string)) ((10,'kpc'),(15,'kpc'))
- float 0.2
- (float, float) (0.2, 0.3)
- ================================== =======================
-
- For example, (10, 'kpc') requests a plot window that is 10 kiloparsecs
- wide in the x and y directions, ((10,'kpc'),(15,'kpc')) requests a
- window that is 10 kiloparsecs wide along the x axis and 15
- kiloparsecs wide along the y axis. In the other two examples, code
- units are assumed, for example (0.2, 0.3) requests a plot that has an
- x width of 0.2 and a y width of 0.3 in code units. If units are
- provided the resulting plot axis labels will use the supplied units.
- image_res : an int or 2-tuple of ints
- Specify the resolution of the resulting image. If not provided, it will be
- determined based on the minimum cell size of the dataset.
- """
- def __init__(self, ds, axis, fields, center="c", width=None, image_res=None, **kwargs):
- fields = ensure_list(fields)
- axis = fix_axis(axis, ds)
- center, dcenter = ds.coordinates.sanitize_center(center, axis)
- slc = ds.slice(axis, center[axis], **kwargs)
- w, frb = construct_image(ds, axis, slc, dcenter, width=width, image_res=image_res)
- super(FITSSlice, self).__init__(frb, fields=fields, wcs=w)
-
-
-class FITSProjection(FITSImageData):
- r"""
- Generate a FITSImageData of an on-axis projection.
-
- Parameters
- ----------
- ds : :class:`yt.data_objects.api.Dataset`
- The dataset object.
- axis : character or integer
- The axis along which to project. One of "x","y","z", or 0,1,2.
- fields : string or list of strings
- The fields to project
- weight_field : string
- The field used to weight the projection.
- center : A sequence of floats, a string, or a tuple.
- The coordinate of the center of the image. If set to 'c', 'center' or
- left blank, the plot is centered on the middle of the domain. If set to
- 'max' or 'm', the center will be located at the maximum of the
- ('gas', 'density') field. Centering on the max or min of a specific
- field is supported by providing a tuple such as ("min","temperature") or
- ("max","dark_matter_density"). Units can be specified by passing in *center*
- as a tuple containing a coordinate and string unit name or by passing
- in a YTArray. If a list or unitless array is supplied, code units are
- assumed.
- width : tuple or a float.
- Width can have four different formats to support windows with variable
- x and y widths. They are:
-
- ================================== =======================
- format example
- ================================== =======================
- (float, string) (10,'kpc')
- ((float, string), (float, string)) ((10,'kpc'),(15,'kpc'))
- float 0.2
- (float, float) (0.2, 0.3)
- ================================== =======================
-
- For example, (10, 'kpc') requests a plot window that is 10 kiloparsecs
- wide in the x and y directions, ((10,'kpc'),(15,'kpc')) requests a
- window that is 10 kiloparsecs wide along the x axis and 15
- kiloparsecs wide along the y axis. In the other two examples, code
- units are assumed, for example (0.2, 0.3) requests a plot that has an
- x width of 0.2 and a y width of 0.3 in code units. If units are
- provided the resulting plot axis labels will use the supplied units.
- image_res : an int or 2-tuple of ints
- Specify the resolution of the resulting image. If not provided, it will be
- determined based on the minimum cell size of the dataset.
- """
- def __init__(self, ds, axis, fields, center="c", width=None,
- weight_field=None, image_res=None, **kwargs):
- fields = ensure_list(fields)
- axis = fix_axis(axis, ds)
- center, dcenter = ds.coordinates.sanitize_center(center, axis)
- prj = ds.proj(fields[0], axis, weight_field=weight_field, **kwargs)
- w, frb = construct_image(ds, axis, prj, dcenter, width=width, image_res=image_res)
- super(FITSProjection, self).__init__(frb, fields=fields, wcs=w)
-
-class FITSOffAxisSlice(FITSImageData):
- r"""
- Generate a FITSImageData of an off-axis slice.
-
- Parameters
- ----------
- ds : :class:`yt.data_objects.api.Dataset`
- The dataset object.
- normal : a sequence of floats
- The vector normal to the projection plane.
- fields : string or list of strings
- The fields to slice
- center : A sequence of floats, a string, or a tuple.
- The coordinate of the center of the image. If set to 'c', 'center' or
- left blank, the plot is centered on the middle of the domain. If set to
- 'max' or 'm', the center will be located at the maximum of the
- ('gas', 'density') field. Centering on the max or min of a specific
- field is supported by providing a tuple such as ("min","temperature") or
- ("max","dark_matter_density"). Units can be specified by passing in *center*
- as a tuple containing a coordinate and string unit name or by passing
- in a YTArray. If a list or unitless array is supplied, code units are
- assumed.
- width : tuple or a float.
- Width can have four different formats to support windows with variable
- x and y widths. They are:
-
- ================================== =======================
- format example
- ================================== =======================
- (float, string) (10,'kpc')
- ((float, string), (float, string)) ((10,'kpc'),(15,'kpc'))
- float 0.2
- (float, float) (0.2, 0.3)
- ================================== =======================
-
- For example, (10, 'kpc') requests a plot window that is 10 kiloparsecs
- wide in the x and y directions, ((10,'kpc'),(15,'kpc')) requests a
- window that is 10 kiloparsecs wide along the x axis and 15
- kiloparsecs wide along the y axis. In the other two examples, code
- units are assumed, for example (0.2, 0.3) requests a plot that has an
- x width of 0.2 and a y width of 0.3 in code units. If units are
- provided the resulting plot axis labels will use the supplied units.
- image_res : an int or 2-tuple of ints
- Specify the resolution of the resulting image.
- north_vector : a sequence of floats
- A vector defining the 'up' direction in the plot. This
- option sets the orientation of the slicing plane. If not
- set, an arbitrary grid-aligned north-vector is chosen.
- """
- def __init__(self, ds, normal, fields, center='c', width=None, image_res=512,
- north_vector=None):
- fields = ensure_list(fields)
- center, dcenter = ds.coordinates.sanitize_center(center, 4)
- cut = ds.cutting(normal, center, north_vector=north_vector)
- center = ds.arr([0.0] * 2, 'code_length')
- w, frb = construct_image(ds, normal, cut, center, width=width, image_res=image_res)
- super(FITSOffAxisSlice, self).__init__(frb, fields=fields, wcs=w)
-
-
-class FITSOffAxisProjection(FITSImageData):
- r"""
- Generate a FITSImageData of an off-axis projection.
-
- Parameters
- ----------
- ds : :class:`yt.data_objects.api.Dataset`
- This is the dataset object corresponding to the
- simulation output to be plotted.
- normal : a sequence of floats
- The vector normal to the projection plane.
- fields : string, list of strings
- The name of the field(s) to be plotted.
- center : A sequence of floats, a string, or a tuple.
- The coordinate of the center of the image. If set to 'c', 'center' or
- left blank, the plot is centered on the middle of the domain. If set to
- 'max' or 'm', the center will be located at the maximum of the
- ('gas', 'density') field. Centering on the max or min of a specific
- field is supported by providing a tuple such as ("min","temperature") or
- ("max","dark_matter_density"). Units can be specified by passing in *center*
- as a tuple containing a coordinate and string unit name or by passing
- in a YTArray. If a list or unitless array is supplied, code units are
- assumed.
- width : tuple or a float.
- Width can have four different formats to support windows with variable
- x and y widths. They are:
-
- ================================== =======================
- format example
- ================================== =======================
- (float, string) (10,'kpc')
- ((float, string), (float, string)) ((10,'kpc'),(15,'kpc'))
- float 0.2
- (float, float) (0.2, 0.3)
- ================================== =======================
-
- For example, (10, 'kpc') requests a plot window that is 10 kiloparsecs
- wide in the x and y directions, ((10,'kpc'),(15,'kpc')) requests a
- window that is 10 kiloparsecs wide along the x axis and 15
- kiloparsecs wide along the y axis. In the other two examples, code
- units are assumed, for example (0.2, 0.3) requests a plot that has an
- x width of 0.2 and a y width of 0.3 in code units. If units are
- provided the resulting plot axis labels will use the supplied units.
- depth : A tuple or a float
- A tuple containing the depth to project through and the string
- key of the unit: (width, 'unit'). If set to a float, code units
- are assumed
- weight_field : string
- The name of the weighting field. Set to None for no weight.
- image_res : an int or 2-tuple of ints
- Specify the resolution of the resulting image.
- depth_res : integer
- Deprecated, this is still in the function signature for API
- compatibility
- north_vector : a sequence of floats
- A vector defining the 'up' direction in the plot. This
- option sets the orientation of the slicing plane. If not
- set, an arbitrary grid-aligned north-vector is chosen.
- method : string
- The method of projection. Valid methods are:
-
- "integrate" with no weight_field specified : integrate the requested
- field along the line of sight.
-
- "integrate" with a weight_field specified : weight the requested
- field by the weighting field and integrate along the line of sight.
-
- "sum" : This method is the same as integrate, except that it does not
- multiply by a path length when performing the integration, and is
- just a straight summation of the field along the given axis. WARNING:
- This should only be used for uniform resolution grid datasets, as other
- datasets may result in unphysical images.
- data_source : yt.data_objects.data_containers.YTSelectionContainer, optional
- If specified, this will be the data source used for selecting regions to project.
- """
- def __init__(self, ds, normal, fields, center='c', width=(1.0, 'unitary'),
- weight_field=None, image_res=512, depth_res=256, data_source=None,
- north_vector=None, depth=(1.0,"unitary"), no_ghost=False, method='integrate'):
- fields = ensure_list(fields)
- center, dcenter = ds.coordinates.sanitize_center(center, 4)
- buf = {}
- width = ds.coordinates.sanitize_width(normal, width, depth)
- wd = tuple(el.in_units('code_length').v for el in width)
- if not iterable(image_res):
- image_res = (image_res, image_res)
- res = (image_res[0], image_res[1])
- if data_source is None:
- source = ds
- else:
- source = data_source
- for field in fields:
- buf[field] = off_axis_projection(source, center, normal, wd, res, field,
- no_ghost=no_ghost, north_vector=north_vector,
- method=method, weight=weight_field).swapaxes(0, 1)
- center = ds.arr([0.0] * 2, 'code_length')
- w, not_an_frb = construct_image(ds, normal, buf, center, width=width, image_res=image_res)
- super(FITSOffAxisProjection, self).__init__(buf, fields=fields, wcs=w)
diff -r 2668f56ec734f8600cfd4c72b06cc035bdd42d3a -r 5603bf81bd7651f2e976755d3d5ade26efc89792 yt/utilities/tests/test_fits_image.py
--- a/yt/utilities/tests/test_fits_image.py
+++ /dev/null
@@ -1,139 +0,0 @@
-"""
-Unit test FITS image creation in yt.
-
-
-
-"""
-
-#-----------------------------------------------------------------------------
-# Copyright (c) 2013, yt Development Team.
-#
-# Distributed under the terms of the Modified BSD License.
-#
-# The full license is in the file COPYING.txt, distributed with this software.
-#-----------------------------------------------------------------------------
-
-import tempfile
-import os
-import shutil
-from yt.testing import fake_random_ds, requires_module
-from yt.convenience import load
-from numpy.testing import \
- assert_equal
-from yt.utilities.fits_image import \
- FITSImageData, FITSProjection, \
- FITSSlice, FITSOffAxisSlice, \
- FITSOffAxisProjection, \
- assert_same_wcs
-from yt.visualization.volume_rendering.off_axis_projection import \
- off_axis_projection
-
-
- at requires_module("astropy")
-def test_fits_image():
- tmpdir = tempfile.mkdtemp()
- curdir = os.getcwd()
- os.chdir(tmpdir)
-
- fields = ("density", "temperature")
- units = ('g/cm**3', 'K',)
- ds = fake_random_ds(64, fields=fields, units=units, nprocs=16,
- length_unit=100.0)
-
- prj = ds.proj("density", 2)
- prj_frb = prj.to_frb((0.5, "unitary"), 128)
-
- fid1 = FITSImageData(prj_frb, fields=["density","temperature"], units="cm")
- fits_prj = FITSProjection(ds, "z", ["density","temperature"], image_res=128,
- width=(0.5,"unitary"))
-
- yield assert_equal, fid1.get_data("density"), fits_prj.get_data("density")
- yield assert_equal, fid1.get_data("temperature"), fits_prj.get_data("temperature")
-
- fid1.writeto("fid1.fits", clobber=True)
- new_fid1 = FITSImageData.from_file("fid1.fits")
-
- yield assert_equal, fid1.get_data("density"), new_fid1.get_data("density")
- yield assert_equal, fid1.get_data("temperature"), new_fid1.get_data("temperature")
-
- ds2 = load("fid1.fits")
- ds2.index
-
- assert ("fits","density") in ds2.field_list
- assert ("fits","temperature") in ds2.field_list
-
- dw_cm = ds2.domain_width.in_units("cm")
-
- assert dw_cm[0].v == 50.
- assert dw_cm[1].v == 50.
-
- slc = ds.slice(2, 0.5)
- slc_frb = slc.to_frb((0.5, "unitary"), 128)
-
- fid2 = FITSImageData(slc_frb, fields=["density","temperature"], units="cm")
- fits_slc = FITSSlice(ds, "z", ["density","temperature"], image_res=128,
- width=(0.5,"unitary"))
-
- yield assert_equal, fid2.get_data("density"), fits_slc.get_data("density")
- yield assert_equal, fid2.get_data("temperature"), fits_slc.get_data("temperature")
-
- dens_img = fid2.pop("density")
- temp_img = fid2.pop("temperature")
-
- # This already has some assertions in it, so we don't need to do anything
- # with it other than just make one
- FITSImageData.from_images([dens_img, temp_img])
-
- cut = ds.cutting([0.1, 0.2, -0.9], [0.5, 0.42, 0.6])
- cut_frb = cut.to_frb((0.5, "unitary"), 128)
-
- fid3 = FITSImageData(cut_frb, fields=["density","temperature"], units="cm")
- fits_cut = FITSOffAxisSlice(ds, [0.1, 0.2, -0.9], ["density","temperature"],
- image_res=128, center=[0.5, 0.42, 0.6],
- width=(0.5,"unitary"))
-
- yield assert_equal, fid3.get_data("density"), fits_cut.get_data("density")
- yield assert_equal, fid3.get_data("temperature"), fits_cut.get_data("temperature")
-
- fid3.create_sky_wcs([30.,45.], (1.0,"arcsec/kpc"))
- fid3.writeto("fid3.fits", clobber=True)
- new_fid3 = FITSImageData.from_file("fid3.fits")
- assert_same_wcs(fid3.wcs, new_fid3.wcs)
- assert new_fid3.wcs.wcs.cunit[0] == "deg"
- assert new_fid3.wcs.wcs.cunit[1] == "deg"
- assert new_fid3.wcs.wcs.ctype[0] == "RA---TAN"
- assert new_fid3.wcs.wcs.ctype[1] == "DEC--TAN"
-
- buf = off_axis_projection(ds, ds.domain_center, [0.1, 0.2, -0.9],
- 0.5, 128, "density").swapaxes(0, 1)
- fid4 = FITSImageData(buf, fields="density", width=100.0)
- fits_oap = FITSOffAxisProjection(ds, [0.1, 0.2, -0.9], "density",
- width=(0.5,"unitary"), image_res=128,
- depth_res=128, depth=(0.5,"unitary"))
-
- yield assert_equal, fid4.get_data("density"), fits_oap.get_data("density")
-
- fid4.create_sky_wcs([30., 45.], (1.0, "arcsec/kpc"), replace_old_wcs=False)
- assert fid4.wcs.wcs.cunit[0] == "cm"
- assert fid4.wcs.wcs.cunit[1] == "cm"
- assert fid4.wcs.wcs.ctype[0] == "linear"
- assert fid4.wcs.wcs.ctype[1] == "linear"
- assert fid4.wcsa.wcs.cunit[0] == "deg"
- assert fid4.wcsa.wcs.cunit[1] == "deg"
- assert fid4.wcsa.wcs.ctype[0] == "RA---TAN"
- assert fid4.wcsa.wcs.ctype[1] == "DEC--TAN"
-
- cvg = ds.covering_grid(ds.index.max_level, [0.25,0.25,0.25],
- [32, 32, 32], fields=["density","temperature"])
- fid5 = FITSImageData(cvg, fields=["density","temperature"])
- assert fid5.dimensionality == 3
-
- fid5.update_header("density", "time", 0.1)
- fid5.update_header("all", "units", "cgs")
-
- assert fid5["density"].header["time"] == 0.1
- assert fid5["temperature"].header["units"] == "cgs"
- assert fid5["density"].header["units"] == "cgs"
-
- os.chdir(curdir)
- shutil.rmtree(tmpdir)
diff -r 2668f56ec734f8600cfd4c72b06cc035bdd42d3a -r 5603bf81bd7651f2e976755d3d5ade26efc89792 yt/visualization/fits_image.py
--- /dev/null
+++ b/yt/visualization/fits_image.py
@@ -0,0 +1,793 @@
+"""
+FITSImageData Class
+"""
+
+#-----------------------------------------------------------------------------
+# Copyright (c) 2013, yt Development Team.
+#
+# Distributed under the terms of the Modified BSD License.
+#
+# The full license is in the file COPYING.txt, distributed with this software.
+#-----------------------------------------------------------------------------
+from yt.extern.six import string_types
+import numpy as np
+from yt.funcs import mylog, iterable, fix_axis, ensure_list
+from yt.visualization.fixed_resolution import FixedResolutionBuffer
+from yt.data_objects.construction_data_containers import YTCoveringGrid
+from yt.utilities.on_demand_imports import _astropy
+from yt.units.yt_array import YTQuantity, YTArray
+from yt.units import dimensions
+from yt.utilities.parallel_tools.parallel_analysis_interface import \
+ parallel_root_only
+from yt.visualization.volume_rendering.off_axis_projection import off_axis_projection
+import re
+
+class FITSImageData(object):
+
+ def __init__(self, data, fields=None, units=None, width=None, wcs=None):
+ r""" Initialize a FITSImageData object.
+
+ FITSImageData contains a collection of FITS ImageHDU instances and
+ WCS information, along with units for each of the images. FITSImageData
+ instances can be constructed from ImageArrays, NumPy arrays, dicts
+ of such arrays, FixedResolutionBuffers, and YTCoveringGrids. The latter
+ two are the most powerful because WCS information can be constructed
+ automatically from their coordinates.
+
+ Parameters
+ ----------
+ data : FixedResolutionBuffer or a YTCoveringGrid. Or, an
+ ImageArray, an numpy.ndarray, or dict of such arrays
+ The data to be made into a FITS image or images.
+ fields : single string or list of strings, optional
+ The field names for the data. If *fields* is none and *data* has
+ keys, it will use these for the fields. If *data* is just a
+ single array one field name must be specified.
+ units : string
+ The units of the WCS coordinates. Defaults to "cm".
+ width : float or YTQuantity
+ The width of the image. Either a single value or iterable of values.
+ If a float, assumed to be in *units*. Only used if this information
+ is not already provided by *data*.
+ wcs : `astropy.wcs.WCS` instance, optional
+ Supply an AstroPy WCS instance. Will override automatic WCS
+ creation from FixedResolutionBuffers and YTCoveringGrids.
+
+ Examples
+ --------
+
+ >>> # This example uses a FRB.
+ >>> ds = load("sloshing_nomag2_hdf5_plt_cnt_0150")
+ >>> prj = ds.proj(2, "kT", weight_field="density")
+ >>> frb = prj.to_frb((0.5, "Mpc"), 800)
+ >>> # This example just uses the FRB and puts the coords in kpc.
+ >>> f_kpc = FITSImageData(frb, fields="kT", units="kpc")
+ >>> # This example specifies a specific WCS.
+ >>> from astropy.wcs import WCS
+ >>> w = WCS(naxis=self.dimensionality)
+ >>> w.wcs.crval = [30., 45.] # RA, Dec in degrees
+ >>> w.wcs.cunit = ["deg"]*2
+ >>> nx, ny = 800, 800
+ >>> w.wcs.crpix = [0.5*(nx+1), 0.5*(ny+1)]
+ >>> w.wcs.ctype = ["RA---TAN","DEC--TAN"]
+ >>> scale = 1./3600. # One arcsec per pixel
+ >>> w.wcs.cdelt = [-scale, scale]
+ >>> f_deg = FITSImageData(frb, fields="kT", wcs=w)
+ >>> f_deg.writeto("temp.fits")
+ """
+
+ if units is None:
+ units = "cm"
+ if width is None:
+ width = 1.0
+
+ exclude_fields = ['x','y','z','px','py','pz',
+ 'pdx','pdy','pdz','weight_field']
+
+ self.hdulist = _astropy.pyfits.HDUList()
+
+ if isinstance(fields, string_types):
+ fields = [fields]
+
+ if hasattr(data, 'keys'):
+ img_data = data
+ if fields is None:
+ fields = list(img_data.keys())
+ elif isinstance(data, np.ndarray):
+ if fields is None:
+ mylog.warning("No field name given for this array. Calling it 'image_data'.")
+ fn = 'image_data'
+ fields = [fn]
+ else:
+ fn = fields[0]
+ img_data = {fn: data}
+
+ self.fields = []
+ for fd in fields:
+ if isinstance(fd, tuple):
+ self.fields.append(fd[1])
+ else:
+ self.fields.append(fd)
+
+ first = True
+ self.field_units = {}
+ for key in fields:
+ if key not in exclude_fields:
+ if hasattr(img_data[key], "units"):
+ self.field_units[key] = img_data[key].units
+ else:
+ self.field_units[key] = "dimensionless"
+ mylog.info("Making a FITS image of field %s" % key)
+ if first:
+ hdu = _astropy.pyfits.PrimaryHDU(np.array(img_data[key]))
+ first = False
+ else:
+ hdu = _astropy.pyfits.ImageHDU(np.array(img_data[key]))
+ hdu.name = key
+ hdu.header["btype"] = key
+ if hasattr(img_data[key], "units"):
+ hdu.header["bunit"] = re.sub('()', '', str(img_data[key].units))
+ self.hdulist.append(hdu)
+
+ self.shape = self.hdulist[0].shape
+ self.dimensionality = len(self.shape)
+
+ if wcs is None:
+ w = _astropy.pywcs.WCS(header=self.hdulist[0].header, naxis=self.dimensionality)
+ if isinstance(img_data, FixedResolutionBuffer):
+ # FRBs are a special case where we have coordinate
+ # information, so we take advantage of this and
+ # construct the WCS object
+ dx = (img_data.bounds[1]-img_data.bounds[0]).in_units(units).v/self.shape[0]
+ dy = (img_data.bounds[3]-img_data.bounds[2]).in_units(units).v/self.shape[1]
+ xctr = 0.5*(img_data.bounds[1]+img_data.bounds[0]).in_units(units).v
+ yctr = 0.5*(img_data.bounds[3]+img_data.bounds[2]).in_units(units).v
+ center = [xctr, yctr]
+ cdelt = [dx,dy]
+ elif isinstance(img_data, YTCoveringGrid):
+ cdelt = img_data.dds.in_units(units).v
+ center = 0.5*(img_data.left_edge+img_data.right_edge).in_units(units).v
+ else:
+ # If img_data is just an array, we assume the center is the origin
+ # and use the image width to determine the cell widths
+ if not iterable(width):
+ width = [width]*self.dimensionality
+ if isinstance(width[0], YTQuantity):
+ cdelt = [wh.in_units(units).v/n for wh, n in zip(width, self.shape)]
+ else:
+ cdelt = [float(wh)/n for wh, n in zip(width, self.shape)]
+ center = [0.0]*self.dimensionality
+ w.wcs.crpix = 0.5*(np.array(self.shape)+1)
+ w.wcs.crval = center
+ w.wcs.cdelt = cdelt
+ w.wcs.ctype = ["linear"]*self.dimensionality
+ w.wcs.cunit = [units]*self.dimensionality
+ self.set_wcs(w)
+ else:
+ self.set_wcs(wcs)
+
+ def set_wcs(self, wcs, wcsname=None, suffix=None):
+ """
+ Set the WCS coordinate information for all images
+ with a WCS object *wcs*.
+ """
+ if wcsname is None:
+ wcs.wcs.name = 'yt'
+ else:
+ wcs.wcs.name = wcsname
+ h = wcs.to_header()
+ if suffix is None:
+ self.wcs = wcs
+ else:
+ setattr(self, "wcs"+suffix, wcs)
+ for img in self.hdulist:
+ for k, v in h.items():
+ kk = k
+ if suffix is not None:
+ kk += suffix
+ img.header[kk] = v
+
+ def update_header(self, field, key, value):
+ """
+ Update the FITS header for *field* with a
+ *key*, *value* pair. If *field* == "all", all
+ headers will be updated.
+ """
+ if field == "all":
+ for img in self.hdulist:
+ img.header[key] = value
+ else:
+ if field not in self.keys():
+ raise KeyError("%s not an image!" % field)
+ idx = self.fields.index(field)
+ self.hdulist[idx].header[key] = value
+
+ def update_all_headers(self, key, value):
+ mylog.warning("update_all_headers is deprecated. "+
+ "Use update_header('all', key, value) instead.")
+ self.update_header("all", key, value)
+
+ def keys(self):
+ return self.fields
+
+ def has_key(self, key):
+ return key in self.fields
+
+ def values(self):
+ return [self.hdulist[k] for k in self.fields]
+
+ def items(self):
+ return [(k, self.hdulist[k]) for k in self.fields]
+
+ def __getitem__(self, item):
+ return self.hdulist[item]
+
+ def info(self):
+ return self.hdulist.info()
+
+ @parallel_root_only
+ def writeto(self, fileobj, fields=None, clobber=False, **kwargs):
+ r"""
+ Write all of the fields or a subset of them to a FITS file.
+
+ Parameters
+ ----------
+ fileobj : string
+ The name of the file to write to.
+ fields : list of strings, optional
+ The fields to write to the file. If not specified
+ all of the fields in the buffer will be written.
+ clobber : boolean, optional
+ Whether or not to overwrite a previously existing file.
+ Default: False
+ All other keyword arguments are passed to the `writeto`
+ method of `astropy.io.fits.HDUList`.
+ """
+ if fields is None:
+ hdus = self.hdulist
+ else:
+ hdus = _astropy.pyfits.HDUList()
+ for field in fields:
+ hdus.append(self.hdulist[field])
+ hdus.writeto(fileobj, clobber=clobber, **kwargs)
+
+ def to_glue(self, label="yt", data_collection=None):
+ """
+ Takes the data in the FITSImageData instance and exports it to
+ Glue (http://www.glueviz.org) for interactive analysis. Optionally
+ add a *label*. If you are already within the Glue environment, you
+ can pass a *data_collection* object, otherwise Glue will be started.
+ """
+ from glue.core import DataCollection, Data
+ from glue.core.coordinates import coordinates_from_header
+ from glue.qt.glue_application import GlueApplication
+
+ image = Data(label=label)
+ image.coords = coordinates_from_header(self.wcs.to_header())
+ for k,f in self.hdulist.items():
+ image.add_component(f.data, k)
+ if data_collection is None:
+ dc = DataCollection([image])
+ app = GlueApplication(dc)
+ app.start()
+ else:
+ data_collection.append(image)
+
+ def to_aplpy(self, **kwargs):
+ """
+ Use APLpy (http://aplpy.github.io) for plotting. Returns an
+ `aplpy.FITSFigure` instance. All keyword arguments are passed to the
+ `aplpy.FITSFigure` constructor.
+ """
+ import aplpy
+ return aplpy.FITSFigure(self.hdulist, **kwargs)
+
+ def get_data(self, field):
+ """
+ Return the data array of the image corresponding to *field*
+ with units attached.
+ """
+ return YTArray(self.hdulist[field].data, self.field_units[field])
+
+ def set_unit(self, field, units):
+ """
+ Set the units of *field* to *units*.
+ """
+ if field not in self.keys():
+ raise KeyError("%s not an image!" % field)
+ idx = self.fields.index(field)
+ new_data = YTArray(self.hdulist[idx].data, self.field_units[field]).in_units(units)
+ self.hdulist[idx].data = new_data.v
+ self.hdulist[idx].header["bunit"] = units
+ self.field_units[field] = units
+
+ def pop(self, key):
+ """
+ Remove a field with name *key*
+ and return it as a new FITSImageData
+ instance.
+ """
+ if key not in self.keys():
+ raise KeyError("%s not an image!" % key)
+ idx = self.fields.index(key)
+ im = self.hdulist.pop(idx)
+ data = YTArray(im.data, self.field_units[key])
+ self.field_units.pop(key)
+ self.fields.remove(key)
+ return FITSImageData(data, fields=key, wcs=self.wcs)
+
+ @classmethod
+ def from_file(cls, filename):
+ """
+ Generate a FITSImageData instance from one previously written to
+ disk.
+
+ Parameters
+ ----------
+ filename : string
+ The name of the file to open.
+ """
+ f = _astropy.pyfits.open(filename)
+ data = {}
+ for hdu in f:
+ data[hdu.header["btype"]] = YTArray(hdu.data, hdu.header["bunit"])
+ f.close()
+ return cls(data, wcs=_astropy.pywcs.WCS(header=hdu.header))
+
+ @classmethod
+ def from_images(cls, image_list):
+ """
+ Generate a new FITSImageData instance from a list of FITSImageData
+ instances.
+
+ Parameters
+ ----------
+ image_list : list of FITSImageData instances
+ The images to be combined.
+ """
+ w = image_list[0].wcs
+ img_shape = image_list[0].shape
+ data = {}
+ for image in image_list:
+ assert_same_wcs(w, image.wcs)
+ if img_shape != image.shape:
+ raise RuntimeError("Images do not have the same shape!")
+ for key in image.keys():
+ data[key] = image.get_data(key)
+ return cls(data, wcs=w)
+
+ def create_sky_wcs(self, sky_center, sky_scale,
+ ctype=["RA---TAN","DEC--TAN"],
+ crota=None, cd=None, pc=None,
+ wcsname="celestial",
+ replace_old_wcs=True):
+ """
+ Takes a Cartesian WCS and converts it to one in a
+ celestial coordinate system.
+
+ Parameters
+ ----------
+ sky_center : iterable of floats
+ Reference coordinates of the WCS in degrees.
+ sky_scale : tuple or YTQuantity
+ Conversion between an angle unit and a length unit,
+ e.g. (3.0, "arcsec/kpc")
+ ctype : list of strings, optional
+ The type of the coordinate system to create.
+ crota : 2-element ndarray, optional
+ Rotation angles between cartesian coordinates and
+ the celestial coordinates.
+ cd : 2x2-element ndarray, optional
+ Dimensioned coordinate transformation matrix.
+ pc : 2x2-element ndarray, optional
+ Coordinate transformation matrix.
+ replace_old_wcs : boolean, optional
+ Whether or not to overwrite the default WCS of the
+ FITSImageData instance. If false, a second WCS will
+ be added to the header. Default: True.
+ """
+ old_wcs = self.wcs
+ naxis = old_wcs.naxis
+ crval = [sky_center[0], sky_center[1]]
+ if isinstance(sky_scale, YTQuantity):
+ scaleq = sky_scale
+ else:
+ scaleq = YTQuantity(sky_scale[0],sky_scale[1])
+ if scaleq.units.dimensions != dimensions.angle/dimensions.length:
+ raise RuntimeError("sky_scale %s not in correct dimensions of angle/length!" % sky_scale)
+ deltas = old_wcs.wcs.cdelt
+ units = [str(unit) for unit in old_wcs.wcs.cunit]
+ new_dx = (YTQuantity(-deltas[0], units[0])*scaleq).in_units("deg")
+ new_dy = (YTQuantity(deltas[1], units[1])*scaleq).in_units("deg")
+ new_wcs = _astropy.pywcs.WCS(naxis=naxis)
+ cdelt = [new_dx.v, new_dy.v]
+ cunit = ["deg"]*2
+ if naxis == 3:
+ crval.append(old_wcs.wcs.crval[2])
+ cdelt.append(old_wcs.wcs.cdelt[2])
+ ctype.append(old_wcs.wcs.ctype[2])
+ cunit.append(old_wcs.wcs.cunit[2])
+ new_wcs.wcs.crpix = old_wcs.wcs.crpix
+ new_wcs.wcs.cdelt = cdelt
+ new_wcs.wcs.crval = crval
+ new_wcs.wcs.cunit = cunit
+ new_wcs.wcs.ctype = ctype
+ if crota is not None:
+ new_wcs.wcs.crota = crota
+ if cd is not None:
+ new_wcs.wcs.cd = cd
+ if pc is not None:
+ new_wcs.wcs.cd = pc
+ if replace_old_wcs:
+ self.set_wcs(new_wcs, wcsname=wcsname)
+ else:
+ self.set_wcs(new_wcs, wcsname=wcsname, suffix="a")
+
+class FITSImageBuffer(FITSImageData):
+ pass
+
+def sanitize_fits_unit(unit):
+ if unit == "Mpc":
+ mylog.info("Changing FITS file unit to kpc.")
+ unit = "kpc"
+ elif unit == "au":
+ unit = "AU"
+ return unit
+
+axis_wcs = [[1,2],[0,2],[0,1]]
+
+def construct_image(ds, axis, data_source, center, width=None, image_res=None):
+ if width is None:
+ width = ds.domain_width[axis_wcs[axis]]
+ unit = ds.get_smallest_appropriate_unit(width[0])
+ mylog.info("Making an image of the entire domain, "+
+ "so setting the center to the domain center.")
+ else:
+ width = ds.coordinates.sanitize_width(axis, width, None)
+ unit = str(width[0].units)
+ if image_res is None:
+ ddims = ds.domain_dimensions*ds.refine_by**ds.index.max_level
+ if iterable(axis):
+ nx = ddims.max()
+ ny = ddims.max()
+ else:
+ nx, ny = [ddims[idx] for idx in axis_wcs[axis]]
+ else:
+ if iterable(image_res):
+ nx, ny = image_res
+ else:
+ nx, ny = image_res, image_res
+ dx = width[0]/nx
+ crpix = [0.5*(nx+1), 0.5*(ny+1)]
+ if hasattr(ds, "wcs") and not iterable(axis):
+ # This is a FITS dataset, so we use it to construct the WCS
+ cunit = [str(ds.wcs.wcs.cunit[idx]) for idx in axis_wcs[axis]]
+ ctype = [ds.wcs.wcs.ctype[idx] for idx in axis_wcs[axis]]
+ cdelt = [ds.wcs.wcs.cdelt[idx] for idx in axis_wcs[axis]]
+ ctr_pix = center.in_units("code_length")[:ds.dimensionality].v
+ crval = ds.wcs.wcs_pix2world(ctr_pix.reshape(1, ds.dimensionality))[0]
+ crval = [crval[idx] for idx in axis_wcs[axis]]
+ else:
+ if unit == "unitary":
+ unit = ds.get_smallest_appropriate_unit(ds.domain_width.max())
+ elif unit == "code_length":
+ unit = ds.get_smallest_appropriate_unit(ds.quan(1.0,"code_length"))
+ unit = sanitize_fits_unit(unit)
+ cunit = [unit]*2
+ ctype = ["LINEAR"]*2
+ cdelt = [dx.in_units(unit)]*2
+ if iterable(axis):
+ crval = center.in_units(unit)
+ else:
+ crval = [center[idx].in_units(unit) for idx in axis_wcs[axis]]
+ if hasattr(data_source, 'to_frb'):
+ if iterable(axis):
+ frb = data_source.to_frb(width[0], (nx, ny), height=width[1])
+ else:
+ frb = data_source.to_frb(width[0], (nx, ny), center=center, height=width[1])
+ else:
+ frb = None
+ w = _astropy.pywcs.WCS(naxis=2)
+ w.wcs.crpix = crpix
+ w.wcs.cdelt = cdelt
+ w.wcs.crval = crval
+ w.wcs.cunit = cunit
+ w.wcs.ctype = ctype
+ return w, frb
+
+def assert_same_wcs(wcs1, wcs2):
+ from numpy.testing import assert_allclose
+ assert wcs1.naxis == wcs2.naxis
+ for i in range(wcs1.naxis):
+ assert wcs1.wcs.cunit[i] == wcs2.wcs.cunit[i]
+ assert wcs1.wcs.ctype[i] == wcs2.wcs.ctype[i]
+ assert_allclose(wcs1.wcs.crpix, wcs2.wcs.crpix)
+ assert_allclose(wcs1.wcs.cdelt, wcs2.wcs.cdelt)
+ assert_allclose(wcs1.wcs.crval, wcs2.wcs.crval)
+ crota1 = getattr(wcs1.wcs, "crota", None)
+ crota2 = getattr(wcs2.wcs, "crota", None)
+ if crota1 is None or crota2 is None:
+ assert crota1 == crota2
+ else:
+ assert_allclose(wcs1.wcs.crota, wcs2.wcs.crota)
+ cd1 = getattr(wcs1.wcs, "cd", None)
+ cd2 = getattr(wcs2.wcs, "cd", None)
+ if cd1 is None or cd2 is None:
+ assert cd1 == cd2
+ else:
+ assert_allclose(wcs1.wcs.cd, wcs2.wcs.cd)
+ pc1 = getattr(wcs1.wcs, "pc", None)
+ pc2 = getattr(wcs2.wcs, "pc", None)
+ if pc1 is None or pc2 is None:
+ assert pc1 == pc2
+ else:
+ assert_allclose(wcs1.wcs.pc, wcs2.wcs.pc)
+
+class FITSSlice(FITSImageData):
+ r"""
+ Generate a FITSImageData of an on-axis slice.
+
+ Parameters
+ ----------
+ ds : :class:`yt.data_objects.api.Dataset`
+ The dataset object.
+ axis : character or integer
+ The axis of the slice. One of "x","y","z", or 0,1,2.
+ fields : string or list of strings
+ The fields to slice
+ center : A sequence of floats, a string, or a tuple.
+ The coordinate of the center of the image. If set to 'c', 'center' or
+ left blank, the plot is centered on the middle of the domain. If set to
+ 'max' or 'm', the center will be located at the maximum of the
+ ('gas', 'density') field. Centering on the max or min of a specific
+ field is supported by providing a tuple such as ("min","temperature") or
+ ("max","dark_matter_density"). Units can be specified by passing in *center*
+ as a tuple containing a coordinate and string unit name or by passing
+ in a YTArray. If a list or unitless array is supplied, code units are
+ assumed.
+ width : tuple or a float.
+ Width can have four different formats to support windows with variable
+ x and y widths. They are:
+
+ ================================== =======================
+ format example
+ ================================== =======================
+ (float, string) (10,'kpc')
+ ((float, string), (float, string)) ((10,'kpc'),(15,'kpc'))
+ float 0.2
+ (float, float) (0.2, 0.3)
+ ================================== =======================
+
+ For example, (10, 'kpc') requests a plot window that is 10 kiloparsecs
+ wide in the x and y directions, ((10,'kpc'),(15,'kpc')) requests a
+ window that is 10 kiloparsecs wide along the x axis and 15
+ kiloparsecs wide along the y axis. In the other two examples, code
+ units are assumed, for example (0.2, 0.3) requests a plot that has an
+ x width of 0.2 and a y width of 0.3 in code units. If units are
+ provided the resulting plot axis labels will use the supplied units.
+ image_res : an int or 2-tuple of ints
+ Specify the resolution of the resulting image. If not provided, it will be
+ determined based on the minimum cell size of the dataset.
+ """
+ def __init__(self, ds, axis, fields, center="c", width=None, image_res=None, **kwargs):
+ fields = ensure_list(fields)
+ axis = fix_axis(axis, ds)
+ center, dcenter = ds.coordinates.sanitize_center(center, axis)
+ slc = ds.slice(axis, center[axis], **kwargs)
+ w, frb = construct_image(ds, axis, slc, dcenter, width=width, image_res=image_res)
+ super(FITSSlice, self).__init__(frb, fields=fields, wcs=w)
+
+
+class FITSProjection(FITSImageData):
+ r"""
+ Generate a FITSImageData of an on-axis projection.
+
+ Parameters
+ ----------
+ ds : :class:`yt.data_objects.api.Dataset`
+ The dataset object.
+ axis : character or integer
+ The axis along which to project. One of "x","y","z", or 0,1,2.
+ fields : string or list of strings
+ The fields to project
+ weight_field : string
+ The field used to weight the projection.
+ center : A sequence of floats, a string, or a tuple.
+ The coordinate of the center of the image. If set to 'c', 'center' or
+ left blank, the plot is centered on the middle of the domain. If set to
+ 'max' or 'm', the center will be located at the maximum of the
+ ('gas', 'density') field. Centering on the max or min of a specific
+ field is supported by providing a tuple such as ("min","temperature") or
+ ("max","dark_matter_density"). Units can be specified by passing in *center*
+ as a tuple containing a coordinate and string unit name or by passing
+ in a YTArray. If a list or unitless array is supplied, code units are
+ assumed.
+ width : tuple or a float.
+ Width can have four different formats to support windows with variable
+ x and y widths. They are:
+
+ ================================== =======================
+ format example
+ ================================== =======================
+ (float, string) (10,'kpc')
+ ((float, string), (float, string)) ((10,'kpc'),(15,'kpc'))
+ float 0.2
+ (float, float) (0.2, 0.3)
+ ================================== =======================
+
+ For example, (10, 'kpc') requests a plot window that is 10 kiloparsecs
+ wide in the x and y directions, ((10,'kpc'),(15,'kpc')) requests a
+ window that is 10 kiloparsecs wide along the x axis and 15
+ kiloparsecs wide along the y axis. In the other two examples, code
+ units are assumed, for example (0.2, 0.3) requests a plot that has an
+ x width of 0.2 and a y width of 0.3 in code units. If units are
+ provided the resulting plot axis labels will use the supplied units.
+ image_res : an int or 2-tuple of ints
+ Specify the resolution of the resulting image. If not provided, it will be
+ determined based on the minimum cell size of the dataset.
+ """
+ def __init__(self, ds, axis, fields, center="c", width=None,
+ weight_field=None, image_res=None, **kwargs):
+ fields = ensure_list(fields)
+ axis = fix_axis(axis, ds)
+ center, dcenter = ds.coordinates.sanitize_center(center, axis)
+ prj = ds.proj(fields[0], axis, weight_field=weight_field, **kwargs)
+ w, frb = construct_image(ds, axis, prj, dcenter, width=width, image_res=image_res)
+ super(FITSProjection, self).__init__(frb, fields=fields, wcs=w)
+
+class FITSOffAxisSlice(FITSImageData):
+ r"""
+ Generate a FITSImageData of an off-axis slice.
+
+ Parameters
+ ----------
+ ds : :class:`yt.data_objects.api.Dataset`
+ The dataset object.
+ normal : a sequence of floats
+ The vector normal to the projection plane.
+ fields : string or list of strings
+ The fields to slice
+ center : A sequence of floats, a string, or a tuple.
+ The coordinate of the center of the image. If set to 'c', 'center' or
+ left blank, the plot is centered on the middle of the domain. If set to
+ 'max' or 'm', the center will be located at the maximum of the
+ ('gas', 'density') field. Centering on the max or min of a specific
+ field is supported by providing a tuple such as ("min","temperature") or
+ ("max","dark_matter_density"). Units can be specified by passing in *center*
+ as a tuple containing a coordinate and string unit name or by passing
+ in a YTArray. If a list or unitless array is supplied, code units are
+ assumed.
+ width : tuple or a float.
+ Width can have four different formats to support windows with variable
+ x and y widths. They are:
+
+ ================================== =======================
+ format example
+ ================================== =======================
+ (float, string) (10,'kpc')
+ ((float, string), (float, string)) ((10,'kpc'),(15,'kpc'))
+ float 0.2
+ (float, float) (0.2, 0.3)
+ ================================== =======================
+
+ For example, (10, 'kpc') requests a plot window that is 10 kiloparsecs
+ wide in the x and y directions, ((10,'kpc'),(15,'kpc')) requests a
+ window that is 10 kiloparsecs wide along the x axis and 15
+ kiloparsecs wide along the y axis. In the other two examples, code
+ units are assumed, for example (0.2, 0.3) requests a plot that has an
+ x width of 0.2 and a y width of 0.3 in code units. If units are
+ provided the resulting plot axis labels will use the supplied units.
+ image_res : an int or 2-tuple of ints
+ Specify the resolution of the resulting image.
+ north_vector : a sequence of floats
+ A vector defining the 'up' direction in the plot. This
+ option sets the orientation of the slicing plane. If not
+ set, an arbitrary grid-aligned north-vector is chosen.
+ """
+ def __init__(self, ds, normal, fields, center='c', width=None, image_res=512,
+ north_vector=None):
+ fields = ensure_list(fields)
+ center, dcenter = ds.coordinates.sanitize_center(center, 4)
+ cut = ds.cutting(normal, center, north_vector=north_vector)
+ center = ds.arr([0.0] * 2, 'code_length')
+ w, frb = construct_image(ds, normal, cut, center, width=width, image_res=image_res)
+ super(FITSOffAxisSlice, self).__init__(frb, fields=fields, wcs=w)
+
+
+class FITSOffAxisProjection(FITSImageData):
+ r"""
+ Generate a FITSImageData of an off-axis projection.
+
+ Parameters
+ ----------
+ ds : :class:`yt.data_objects.api.Dataset`
+ This is the dataset object corresponding to the
+ simulation output to be plotted.
+ normal : a sequence of floats
+ The vector normal to the projection plane.
+ fields : string, list of strings
+ The name of the field(s) to be plotted.
+ center : A sequence of floats, a string, or a tuple.
+ The coordinate of the center of the image. If set to 'c', 'center' or
+ left blank, the plot is centered on the middle of the domain. If set to
+ 'max' or 'm', the center will be located at the maximum of the
+ ('gas', 'density') field. Centering on the max or min of a specific
+ field is supported by providing a tuple such as ("min","temperature") or
+ ("max","dark_matter_density"). Units can be specified by passing in *center*
+ as a tuple containing a coordinate and string unit name or by passing
+ in a YTArray. If a list or unitless array is supplied, code units are
+ assumed.
+ width : tuple or a float.
+ Width can have four different formats to support windows with variable
+ x and y widths. They are:
+
+ ================================== =======================
+ format example
+ ================================== =======================
+ (float, string) (10,'kpc')
+ ((float, string), (float, string)) ((10,'kpc'),(15,'kpc'))
+ float 0.2
+ (float, float) (0.2, 0.3)
+ ================================== =======================
+
+ For example, (10, 'kpc') requests a plot window that is 10 kiloparsecs
+ wide in the x and y directions, ((10,'kpc'),(15,'kpc')) requests a
+ window that is 10 kiloparsecs wide along the x axis and 15
+ kiloparsecs wide along the y axis. In the other two examples, code
+ units are assumed, for example (0.2, 0.3) requests a plot that has an
+ x width of 0.2 and a y width of 0.3 in code units. If units are
+ provided the resulting plot axis labels will use the supplied units.
+ depth : A tuple or a float
+ A tuple containing the depth to project through and the string
+ key of the unit: (width, 'unit'). If set to a float, code units
+ are assumed
+ weight_field : string
+ The name of the weighting field. Set to None for no weight.
+ image_res : an int or 2-tuple of ints
+ Specify the resolution of the resulting image.
+ depth_res : integer
+ Deprecated, this is still in the function signature for API
+ compatibility
+ north_vector : a sequence of floats
+ A vector defining the 'up' direction in the plot. This
+ option sets the orientation of the slicing plane. If not
+ set, an arbitrary grid-aligned north-vector is chosen.
+ method : string
+ The method of projection. Valid methods are:
+
+ "integrate" with no weight_field specified : integrate the requested
+ field along the line of sight.
+
+ "integrate" with a weight_field specified : weight the requested
+ field by the weighting field and integrate along the line of sight.
+
+ "sum" : This method is the same as integrate, except that it does not
+ multiply by a path length when performing the integration, and is
+ just a straight summation of the field along the given axis. WARNING:
+ This should only be used for uniform resolution grid datasets, as other
+ datasets may result in unphysical images.
+ data_source : yt.data_objects.data_containers.YTSelectionContainer, optional
+ If specified, this will be the data source used for selecting regions to project.
+ """
+ def __init__(self, ds, normal, fields, center='c', width=(1.0, 'unitary'),
+ weight_field=None, image_res=512, depth_res=256, data_source=None,
+ north_vector=None, depth=(1.0,"unitary"), no_ghost=False, method='integrate'):
+ fields = ensure_list(fields)
+ center, dcenter = ds.coordinates.sanitize_center(center, 4)
+ buf = {}
+ width = ds.coordinates.sanitize_width(normal, width, depth)
+ wd = tuple(el.in_units('code_length').v for el in width)
+ if not iterable(image_res):
+ image_res = (image_res, image_res)
+ res = (image_res[0], image_res[1])
+ if data_source is None:
+ source = ds
+ else:
+ source = data_source
+ for field in fields:
+ buf[field] = off_axis_projection(source, center, normal, wd, res, field,
+ no_ghost=no_ghost, north_vector=north_vector,
+ method=method, weight=weight_field).swapaxes(0, 1)
+ center = ds.arr([0.0] * 2, 'code_length')
+ w, not_an_frb = construct_image(ds, normal, buf, center, width=width, image_res=image_res)
+ super(FITSOffAxisProjection, self).__init__(buf, fields=fields, wcs=w)
diff -r 2668f56ec734f8600cfd4c72b06cc035bdd42d3a -r 5603bf81bd7651f2e976755d3d5ade26efc89792 yt/visualization/tests/test_fits_image.py
--- /dev/null
+++ b/yt/visualization/tests/test_fits_image.py
@@ -0,0 +1,139 @@
+"""
+Unit test FITS image creation in yt.
+
+
+
+"""
+
+#-----------------------------------------------------------------------------
+# Copyright (c) 2013, yt Development Team.
+#
+# Distributed under the terms of the Modified BSD License.
+#
+# The full license is in the file COPYING.txt, distributed with this software.
+#-----------------------------------------------------------------------------
+
+import tempfile
+import os
+import shutil
+from yt.testing import fake_random_ds, requires_module
+from yt.convenience import load
+from numpy.testing import \
+ assert_equal
+from yt.visualization.fits_image import \
+ FITSImageData, FITSProjection, \
+ FITSSlice, FITSOffAxisSlice, \
+ FITSOffAxisProjection, \
+ assert_same_wcs
+from yt.visualization.volume_rendering.off_axis_projection import \
+ off_axis_projection
+
+
+ at requires_module("astropy")
+def test_fits_image():
+ tmpdir = tempfile.mkdtemp()
+ curdir = os.getcwd()
+ os.chdir(tmpdir)
+
+ fields = ("density", "temperature")
+ units = ('g/cm**3', 'K',)
+ ds = fake_random_ds(64, fields=fields, units=units, nprocs=16,
+ length_unit=100.0)
+
+ prj = ds.proj("density", 2)
+ prj_frb = prj.to_frb((0.5, "unitary"), 128)
+
+ fid1 = FITSImageData(prj_frb, fields=["density","temperature"], units="cm")
+ fits_prj = FITSProjection(ds, "z", ["density","temperature"], image_res=128,
+ width=(0.5,"unitary"))
+
+ yield assert_equal, fid1.get_data("density"), fits_prj.get_data("density")
+ yield assert_equal, fid1.get_data("temperature"), fits_prj.get_data("temperature")
+
+ fid1.writeto("fid1.fits", clobber=True)
+ new_fid1 = FITSImageData.from_file("fid1.fits")
+
+ yield assert_equal, fid1.get_data("density"), new_fid1.get_data("density")
+ yield assert_equal, fid1.get_data("temperature"), new_fid1.get_data("temperature")
+
+ ds2 = load("fid1.fits")
+ ds2.index
+
+ assert ("fits","density") in ds2.field_list
+ assert ("fits","temperature") in ds2.field_list
+
+ dw_cm = ds2.domain_width.in_units("cm")
+
+ assert dw_cm[0].v == 50.
+ assert dw_cm[1].v == 50.
+
+ slc = ds.slice(2, 0.5)
+ slc_frb = slc.to_frb((0.5, "unitary"), 128)
+
+ fid2 = FITSImageData(slc_frb, fields=["density","temperature"], units="cm")
+ fits_slc = FITSSlice(ds, "z", ["density","temperature"], image_res=128,
+ width=(0.5,"unitary"))
+
+ yield assert_equal, fid2.get_data("density"), fits_slc.get_data("density")
+ yield assert_equal, fid2.get_data("temperature"), fits_slc.get_data("temperature")
+
+ dens_img = fid2.pop("density")
+ temp_img = fid2.pop("temperature")
+
+ # This already has some assertions in it, so we don't need to do anything
+ # with it other than just make one
+ FITSImageData.from_images([dens_img, temp_img])
+
+ cut = ds.cutting([0.1, 0.2, -0.9], [0.5, 0.42, 0.6])
+ cut_frb = cut.to_frb((0.5, "unitary"), 128)
+
+ fid3 = FITSImageData(cut_frb, fields=["density","temperature"], units="cm")
+ fits_cut = FITSOffAxisSlice(ds, [0.1, 0.2, -0.9], ["density","temperature"],
+ image_res=128, center=[0.5, 0.42, 0.6],
+ width=(0.5,"unitary"))
+
+ yield assert_equal, fid3.get_data("density"), fits_cut.get_data("density")
+ yield assert_equal, fid3.get_data("temperature"), fits_cut.get_data("temperature")
+
+ fid3.create_sky_wcs([30.,45.], (1.0,"arcsec/kpc"))
+ fid3.writeto("fid3.fits", clobber=True)
+ new_fid3 = FITSImageData.from_file("fid3.fits")
+ assert_same_wcs(fid3.wcs, new_fid3.wcs)
+ assert new_fid3.wcs.wcs.cunit[0] == "deg"
+ assert new_fid3.wcs.wcs.cunit[1] == "deg"
+ assert new_fid3.wcs.wcs.ctype[0] == "RA---TAN"
+ assert new_fid3.wcs.wcs.ctype[1] == "DEC--TAN"
+
+ buf = off_axis_projection(ds, ds.domain_center, [0.1, 0.2, -0.9],
+ 0.5, 128, "density").swapaxes(0, 1)
+ fid4 = FITSImageData(buf, fields="density", width=100.0)
+ fits_oap = FITSOffAxisProjection(ds, [0.1, 0.2, -0.9], "density",
+ width=(0.5,"unitary"), image_res=128,
+ depth_res=128, depth=(0.5,"unitary"))
+
+ yield assert_equal, fid4.get_data("density"), fits_oap.get_data("density")
+
+ fid4.create_sky_wcs([30., 45.], (1.0, "arcsec/kpc"), replace_old_wcs=False)
+ assert fid4.wcs.wcs.cunit[0] == "cm"
+ assert fid4.wcs.wcs.cunit[1] == "cm"
+ assert fid4.wcs.wcs.ctype[0] == "linear"
+ assert fid4.wcs.wcs.ctype[1] == "linear"
+ assert fid4.wcsa.wcs.cunit[0] == "deg"
+ assert fid4.wcsa.wcs.cunit[1] == "deg"
+ assert fid4.wcsa.wcs.ctype[0] == "RA---TAN"
+ assert fid4.wcsa.wcs.ctype[1] == "DEC--TAN"
+
+ cvg = ds.covering_grid(ds.index.max_level, [0.25,0.25,0.25],
+ [32, 32, 32], fields=["density","temperature"])
+ fid5 = FITSImageData(cvg, fields=["density","temperature"])
+ assert fid5.dimensionality == 3
+
+ fid5.update_header("density", "time", 0.1)
+ fid5.update_header("all", "units", "cgs")
+
+ assert fid5["density"].header["time"] == 0.1
+ assert fid5["temperature"].header["units"] == "cgs"
+ assert fid5["density"].header["units"] == "cgs"
+
+ os.chdir(curdir)
+ shutil.rmtree(tmpdir)
https://bitbucket.org/yt_analysis/yt/commits/7b8135f94264/
Changeset: 7b8135f94264
Branch: yt
User: jzuhone
Date: 2016-12-05 17:53:39+00:00
Summary: Move these up to the top
Affected #: 2 files
diff -r 5603bf81bd7651f2e976755d3d5ade26efc89792 -r 7b8135f94264ab1cf66362570b9c683ba6a7b293 yt/__init__.py
--- a/yt/__init__.py
+++ b/yt/__init__.py
@@ -156,7 +156,9 @@
ProjectionPlot, OffAxisProjectionPlot, \
show_colormaps, add_cmap, make_colormap, \
ProfilePlot, PhasePlot, ParticlePhasePlot, \
- ParticleProjectionPlot, ParticleImageBuffer, ParticlePlot
+ ParticleProjectionPlot, ParticleImageBuffer, ParticlePlot, \
+ FITSImageData, FITSSlice, FITSProjection, FITSOffAxisSlice, \
+ FITSOffAxisProjection
from yt.visualization.volume_rendering.api import \
volume_render, create_scene, ColorTransferFunction, TransferFunction, \
diff -r 5603bf81bd7651f2e976755d3d5ade26efc89792 -r 7b8135f94264ab1cf66362570b9c683ba6a7b293 yt/visualization/api.py
--- a/yt/visualization/api.py
+++ b/yt/visualization/api.py
@@ -59,3 +59,10 @@
from .base_plot_types import \
get_multi_plot
+
+from .fits_image import \
+ FITSImageData, \
+ FITSSlice, \
+ FITSOffAxisSlice, \
+ FITSProjection, \
+ FITSOffAxisProjection
https://bitbucket.org/yt_analysis/yt/commits/5eef0ba576ed/
Changeset: 5eef0ba576ed
Branch: yt
User: jzuhone
Date: 2016-12-05 17:59:10+00:00
Summary: Update docs appropriately
Affected #: 1 file
diff -r 7b8135f94264ab1cf66362570b9c683ba6a7b293 -r 5eef0ba576ed6420ca70db05c089b6767aeb7eb6 doc/source/visualizing/FITSImageData.ipynb
--- a/doc/source/visualizing/FITSImageData.ipynb
+++ b/doc/source/visualizing/FITSImageData.ipynb
@@ -15,8 +15,7 @@
},
"outputs": [],
"source": [
- "import yt\n",
- "from yt.utilities.fits_image import FITSImageData, FITSProjection"
+ "import yt"
]
},
{
@@ -73,7 +72,7 @@
},
"outputs": [],
"source": [
- "prj_fits = FITSProjection(ds, \"z\", [\"temperature\"], weight_field=\"density\")"
+ "prj_fits = yt.FITSProjection(ds, \"z\", [\"temperature\"], weight_field=\"density\")"
]
},
{
@@ -236,7 +235,7 @@
"source": [
"slc3 = ds.slice(0, 0.0)\n",
"frb = slc3.to_frb((500.,\"kpc\"), 800)\n",
- "fid_frb = FITSImageData(frb, fields=[\"density\",\"temperature\"], units=\"pc\")"
+ "fid_frb = yt.FITSImageData(frb, fields=[\"density\",\"temperature\"], units=\"pc\")"
]
},
{
@@ -255,7 +254,7 @@
"outputs": [],
"source": [
"cvg = ds.covering_grid(ds.index.max_level, [-0.5,-0.5,-0.5], [64, 64, 64], fields=[\"density\",\"temperature\"])\n",
- "fid_cvg = FITSImageData(cvg, fields=[\"density\",\"temperature\"], units=\"Mpc\")"
+ "fid_cvg = yt.FITSImageData(cvg, fields=[\"density\",\"temperature\"], units=\"Mpc\")"
]
},
{
@@ -280,7 +279,7 @@
},
"outputs": [],
"source": [
- "fid = FITSImageData.from_file(\"sloshing.fits\")\n",
+ "fid = yt.FITSImageData.from_file(\"sloshing.fits\")\n",
"fid.info()"
]
},
@@ -299,8 +298,8 @@
},
"outputs": [],
"source": [
- "prj_fits2 = FITSProjection(ds, \"z\", [\"density\"])\n",
- "prj_fits3 = FITSImageData.from_images([prj_fits, prj_fits2])\n",
+ "prj_fits2 = yt.FITSProjection(ds, \"z\", [\"density\"])\n",
+ "prj_fits3 = yt.FITSImageData.from_images([prj_fits, prj_fits2])\n",
"prj_fits3.info()"
]
},
https://bitbucket.org/yt_analysis/yt/commits/c4f9cefc1435/
Changeset: c4f9cefc1435
Branch: yt
User: jzuhone
Date: 2016-12-05 17:59:29+00:00
Summary: Import this inside the function that uses it
Affected #: 1 file
diff -r 5eef0ba576ed6420ca70db05c089b6767aeb7eb6 -r c4f9cefc1435c3fcf4b56995c8da18ed3a11d7f3 yt/frontends/fits/misc.py
--- a/yt/frontends/fits/misc.py
+++ b/yt/frontends/fits/misc.py
@@ -17,7 +17,6 @@
from yt.funcs import mylog
from yt.utilities.on_demand_imports import _astropy
from yt.units.yt_array import YTQuantity, YTArray
-from yt.utilities.fits_image import FITSImageData
if PY3:
from io import BytesIO as IO
else:
@@ -109,6 +108,7 @@
... nan_mask=0.0)
"""
from spectral_cube import SpectralCube
+ from yt.visualization.fits_image import FITSImageData
from yt.frontends.fits.api import FITSDataset
cube = SpectralCube.read(filename)
if not isinstance(slab_width, YTQuantity):
https://bitbucket.org/yt_analysis/yt/commits/6477f7728cd8/
Changeset: 6477f7728cd8
Branch: yt
User: jzuhone
Date: 2016-12-05 19:26:50+00:00
Summary: Put a fits_image.py in utilities and issue a deprecation warning. Add a helper function for issuing deprecation warnings.
Affected #: 3 files
diff -r c4f9cefc1435c3fcf4b56995c8da18ed3a11d7f3 -r 6477f7728cd8eae4dca45fa8a0f901f28d7abc92 yt/analysis_modules/photon_simulator/api.py
--- a/yt/analysis_modules/photon_simulator/api.py
+++ b/yt/analysis_modules/photon_simulator/api.py
@@ -10,12 +10,10 @@
# The full license is in the file COPYING.txt, distributed with this software.
#-----------------------------------------------------------------------------
-from numpy import VisibleDeprecationWarning
+from yt.funcs import issue_deprecation_warning
-import warnings
-warnings.warn("The photon_simulator module is deprecated. "
- "Please use pyXSIM (http://hea-www.cfa.harvard.edu/~jzuhone/pyxsim) instead.",
- VisibleDeprecationWarning, stacklevel=2)
+issue_deprecation_warning("The photon_simulator module is deprecated. Please use pyXSIM "
+ "(http://hea-www.cfa.harvard.edu/~jzuhone/pyxsim) instead.")
from .photon_models import \
PhotonModel, \
diff -r c4f9cefc1435c3fcf4b56995c8da18ed3a11d7f3 -r 6477f7728cd8eae4dca45fa8a0f901f28d7abc92 yt/funcs.py
--- a/yt/funcs.py
+++ b/yt/funcs.py
@@ -1036,3 +1036,7 @@
return val.decode('utf8')
else:
return val
+
+def issue_deprecation_warning(msg):
+ from numpy import VisibleDeprecationWarning
+ warnings.warn(msg, VisibleDeprecationWarning, stacklevel=3)
diff -r c4f9cefc1435c3fcf4b56995c8da18ed3a11d7f3 -r 6477f7728cd8eae4dca45fa8a0f901f28d7abc92 yt/utilities/fits_image.py
--- /dev/null
+++ b/yt/utilities/fits_image.py
@@ -0,0 +1,9 @@
+from yt.funcs import issue_deprecation_warning
+
+issue_deprecation_warning("The fits_image classes have been moved to yt.visualization "
+ "and are available as top-level imports. yt.utilities.fits_image "
+ "is deprecated.")
+
+from yt.visualization.fits_image import \
+ FITSImageData, FITSSlice, FITSProjection, \
+ FITSOffAxisSlice, FITSOffAxisProjection
\ No newline at end of file
https://bitbucket.org/yt_analysis/yt/commits/e622e127f4c3/
Changeset: e622e127f4c3
Branch: yt
User: jzuhone
Date: 2016-12-05 19:36:03+00:00
Summary: Make this message a bit clearer
Affected #: 1 file
diff -r 6477f7728cd8eae4dca45fa8a0f901f28d7abc92 -r e622e127f4c3dc3fe1ede1114eefcf7697baff19 yt/utilities/fits_image.py
--- a/yt/utilities/fits_image.py
+++ b/yt/utilities/fits_image.py
@@ -1,7 +1,7 @@
from yt.funcs import issue_deprecation_warning
issue_deprecation_warning("The fits_image classes have been moved to yt.visualization "
- "and are available as top-level imports. yt.utilities.fits_image "
+ "and can be imported from yt directly. yt.utilities.fits_image "
"is deprecated.")
from yt.visualization.fits_image import \
https://bitbucket.org/yt_analysis/yt/commits/816b4506d800/
Changeset: 816b4506d800
Branch: yt
User: jzuhone
Date: 2016-12-05 22:53:13+00:00
Summary: These are needed by some analysis modules
Affected #: 1 file
diff -r e622e127f4c3dc3fe1ede1114eefcf7697baff19 -r 816b4506d800adbaf5c3eba9b82e17244d70741e yt/utilities/fits_image.py
--- a/yt/utilities/fits_image.py
+++ b/yt/utilities/fits_image.py
@@ -6,4 +6,5 @@
from yt.visualization.fits_image import \
FITSImageData, FITSSlice, FITSProjection, \
- FITSOffAxisSlice, FITSOffAxisProjection
\ No newline at end of file
+ FITSOffAxisSlice, FITSOffAxisProjection, \
+ assert_same_wcs, sanitize_fits_unit
\ No newline at end of file
https://bitbucket.org/yt_analysis/yt/commits/059ef21284de/
Changeset: 059ef21284de
Branch: yt
User: jzuhone
Date: 2016-12-05 22:59:07+00:00
Summary: Don’t let flake8 check this file
Affected #: 1 file
diff -r 816b4506d800adbaf5c3eba9b82e17244d70741e -r 059ef21284de22c40239d3014271ca4fda4a0fc3 setup.cfg
--- a/setup.cfg
+++ b/setup.cfg
@@ -13,6 +13,6 @@
# unused import errors
# autogenerated __config__.py files
# vendored libraries
-exclude = doc,benchmarks,*/api.py,*/__init__.py,*/__config__.py,yt/visualization/_mpl_imports.py,yt/utilities/lodgeit.py,yt/utilities/lru_cache.py,yt/utilities/poster/*,yt/extern/*,yt/mods.py
+exclude = doc,benchmarks,*/api.py,*/__init__.py,*/__config__.py,yt/visualization/_mpl_imports.py,yt/utilities/lodgeit.py,yt/utilities/lru_cache.py,yt/utilities/poster/*,yt/extern/*,yt/mods.py,yt/utilities/fits_image.py
max-line-length=999
ignore = E111,E121,E122,E123,E124,E125,E126,E127,E128,E129,E131,E201,E202,E211,E221,E222,E227,E228,E241,E301,E203,E225,E226,E231,E251,E261,E262,E265,E266,E302,E303,E402,E502,E701,E703,E731,W291,W292,W293,W391,W503
\ No newline at end of file
https://bitbucket.org/yt_analysis/yt/commits/1b22afa99566/
Changeset: 1b22afa99566
Branch: yt
User: jzuhone
Date: 2016-12-05 23:12:41+00:00
Summary: Make sure we import from the new place. Also removed a line that wasn’t used.
Affected #: 5 files
diff -r 059ef21284de22c40239d3014271ca4fda4a0fc3 -r 1b22afa99566b4d0dcaad8fe6ba898e07c9536fa yt/analysis_modules/photon_simulator/photon_simulator.py
--- a/yt/analysis_modules/photon_simulator/photon_simulator.py
+++ b/yt/analysis_modules/photon_simulator/photon_simulator.py
@@ -31,7 +31,7 @@
from yt.utilities.physical_constants import clight
from yt.utilities.cosmology import Cosmology
from yt.utilities.orientation import Orientation
-from yt.utilities.fits_image import assert_same_wcs
+from yt.visualization.fits_image import assert_same_wcs
from yt.utilities.parallel_tools.parallel_analysis_interface import \
communication_system, parallel_root_only, get_mpi_type, \
parallel_capable
diff -r 059ef21284de22c40239d3014271ca4fda4a0fc3 -r 1b22afa99566b4d0dcaad8fe6ba898e07c9536fa yt/analysis_modules/ppv_cube/ppv_cube.py
--- a/yt/analysis_modules/ppv_cube/ppv_cube.py
+++ b/yt/analysis_modules/ppv_cube/ppv_cube.py
@@ -13,7 +13,7 @@
import numpy as np
from yt.utilities.on_demand_imports import _astropy
from yt.utilities.orientation import Orientation
-from yt.utilities.fits_image import FITSImageData, sanitize_fits_unit
+from yt.visualization.fits_image import FITSImageData, sanitize_fits_unit
from yt.visualization.volume_rendering.off_axis_projection import off_axis_projection
from yt.funcs import get_pbar
from yt.utilities.physical_constants import clight, mh
diff -r 059ef21284de22c40239d3014271ca4fda4a0fc3 -r 1b22afa99566b4d0dcaad8fe6ba898e07c9536fa yt/analysis_modules/sunyaev_zeldovich/projection.py
--- a/yt/analysis_modules/sunyaev_zeldovich/projection.py
+++ b/yt/analysis_modules/sunyaev_zeldovich/projection.py
@@ -375,7 +375,7 @@
>>> sky_center = (30., 45., "deg")
>>> szprj.write_fits("SZbullet.fits", sky_center=sky_center, sky_scale=sky_scale)
"""
- from yt.utilities.fits_image import FITSImageData
+ from yt.visualization.fits_image import FITSImageData
dx = self.dx.in_units("kpc")
dy = dx
diff -r 059ef21284de22c40239d3014271ca4fda4a0fc3 -r 1b22afa99566b4d0dcaad8fe6ba898e07c9536fa yt/visualization/fixed_resolution.py
--- a/yt/visualization/fixed_resolution.py
+++ b/yt/visualization/fixed_resolution.py
@@ -362,7 +362,7 @@
the length units that the coordinates are written in, default 'cm'.
"""
- from yt.utilities.fits_image import FITSImageData
+ from yt.visualization.fits_image import FITSImageData
if fields is None:
fields = list(self.data.keys())
diff -r 059ef21284de22c40239d3014271ca4fda4a0fc3 -r 1b22afa99566b4d0dcaad8fe6ba898e07c9536fa yt/visualization/volume_rendering/image_handling.py
--- a/yt/visualization/volume_rendering/image_handling.py
+++ b/yt/visualization/volume_rendering/image_handling.py
@@ -33,13 +33,12 @@
f.create_dataset("A", data=image[:,:,3])
f.close()
if fits:
- from yt.utilities.fits_image import FITSImageData
+ from yt.visualization.fits_image import FITSImageData
data = {}
data["r"] = image[:,:,0]
data["g"] = image[:,:,1]
data["b"] = image[:,:,2]
data["a"] = image[:,:,3]
- nx, ny = data["r"].shape
fib = FITSImageData(data)
fib.writeto('%s.fits'%fn,clobber=True)
https://bitbucket.org/yt_analysis/yt/commits/0c1251988844/
Changeset: 0c1251988844
Branch: yt
User: jzuhone
Date: 2016-12-05 23:21:58+00:00
Summary: remove these
Affected #: 1 file
diff -r 1b22afa99566b4d0dcaad8fe6ba898e07c9536fa -r 0c1251988844f7138e121a52a16e35037f8a599c yt/utilities/fits_image.py
--- a/yt/utilities/fits_image.py
+++ b/yt/utilities/fits_image.py
@@ -6,5 +6,4 @@
from yt.visualization.fits_image import \
FITSImageData, FITSSlice, FITSProjection, \
- FITSOffAxisSlice, FITSOffAxisProjection, \
- assert_same_wcs, sanitize_fits_unit
\ No newline at end of file
+ FITSOffAxisSlice, FITSOffAxisProjection
\ No newline at end of file
https://bitbucket.org/yt_analysis/yt/commits/959db1b322f8/
Changeset: 959db1b322f8
Branch: yt
User: ngoldbaum
Date: 2016-12-08 15:29:30+00:00
Summary: Merged in jzuhone/yt (pull request #2457)
Move fits_image to visualization
Affected #: 16 files
diff -r 83d2f4c6b6608ea701acf74d326c44125ae93b5e -r 959db1b322f8836bab7a600487819d8cef23dd13 doc/source/visualizing/FITSImageData.ipynb
--- a/doc/source/visualizing/FITSImageData.ipynb
+++ b/doc/source/visualizing/FITSImageData.ipynb
@@ -15,8 +15,7 @@
},
"outputs": [],
"source": [
- "import yt\n",
- "from yt.utilities.fits_image import FITSImageData, FITSProjection"
+ "import yt"
]
},
{
@@ -73,7 +72,7 @@
},
"outputs": [],
"source": [
- "prj_fits = FITSProjection(ds, \"z\", [\"temperature\"], weight_field=\"density\")"
+ "prj_fits = yt.FITSProjection(ds, \"z\", [\"temperature\"], weight_field=\"density\")"
]
},
{
@@ -236,7 +235,7 @@
"source": [
"slc3 = ds.slice(0, 0.0)\n",
"frb = slc3.to_frb((500.,\"kpc\"), 800)\n",
- "fid_frb = FITSImageData(frb, fields=[\"density\",\"temperature\"], units=\"pc\")"
+ "fid_frb = yt.FITSImageData(frb, fields=[\"density\",\"temperature\"], units=\"pc\")"
]
},
{
@@ -255,7 +254,7 @@
"outputs": [],
"source": [
"cvg = ds.covering_grid(ds.index.max_level, [-0.5,-0.5,-0.5], [64, 64, 64], fields=[\"density\",\"temperature\"])\n",
- "fid_cvg = FITSImageData(cvg, fields=[\"density\",\"temperature\"], units=\"Mpc\")"
+ "fid_cvg = yt.FITSImageData(cvg, fields=[\"density\",\"temperature\"], units=\"Mpc\")"
]
},
{
@@ -280,7 +279,7 @@
},
"outputs": [],
"source": [
- "fid = FITSImageData.from_file(\"sloshing.fits\")\n",
+ "fid = yt.FITSImageData.from_file(\"sloshing.fits\")\n",
"fid.info()"
]
},
@@ -299,8 +298,8 @@
},
"outputs": [],
"source": [
- "prj_fits2 = FITSProjection(ds, \"z\", [\"density\"])\n",
- "prj_fits3 = FITSImageData.from_images([prj_fits, prj_fits2])\n",
+ "prj_fits2 = yt.FITSProjection(ds, \"z\", [\"density\"])\n",
+ "prj_fits3 = yt.FITSImageData.from_images([prj_fits, prj_fits2])\n",
"prj_fits3.info()"
]
},
diff -r 83d2f4c6b6608ea701acf74d326c44125ae93b5e -r 959db1b322f8836bab7a600487819d8cef23dd13 setup.cfg
--- a/setup.cfg
+++ b/setup.cfg
@@ -13,6 +13,6 @@
# unused import errors
# autogenerated __config__.py files
# vendored libraries
-exclude = doc,benchmarks,*/api.py,*/__init__.py,*/__config__.py,yt/visualization/_mpl_imports.py,yt/utilities/lodgeit.py,yt/utilities/lru_cache.py,yt/utilities/poster/*,yt/extern/*,yt/mods.py
+exclude = doc,benchmarks,*/api.py,*/__init__.py,*/__config__.py,yt/visualization/_mpl_imports.py,yt/utilities/lodgeit.py,yt/utilities/lru_cache.py,yt/utilities/poster/*,yt/extern/*,yt/mods.py,yt/utilities/fits_image.py
max-line-length=999
ignore = E111,E121,E122,E123,E124,E125,E126,E127,E128,E129,E131,E201,E202,E211,E221,E222,E227,E228,E241,E301,E203,E225,E226,E231,E251,E261,E262,E265,E266,E302,E303,E402,E502,E701,E703,E731,W291,W292,W293,W391,W503
\ No newline at end of file
diff -r 83d2f4c6b6608ea701acf74d326c44125ae93b5e -r 959db1b322f8836bab7a600487819d8cef23dd13 yt/__init__.py
--- a/yt/__init__.py
+++ b/yt/__init__.py
@@ -156,7 +156,9 @@
ProjectionPlot, OffAxisProjectionPlot, \
show_colormaps, add_cmap, make_colormap, \
ProfilePlot, PhasePlot, ParticlePhasePlot, \
- ParticleProjectionPlot, ParticleImageBuffer, ParticlePlot
+ ParticleProjectionPlot, ParticleImageBuffer, ParticlePlot, \
+ FITSImageData, FITSSlice, FITSProjection, FITSOffAxisSlice, \
+ FITSOffAxisProjection
from yt.visualization.volume_rendering.api import \
volume_render, create_scene, ColorTransferFunction, TransferFunction, \
diff -r 83d2f4c6b6608ea701acf74d326c44125ae93b5e -r 959db1b322f8836bab7a600487819d8cef23dd13 yt/analysis_modules/photon_simulator/api.py
--- a/yt/analysis_modules/photon_simulator/api.py
+++ b/yt/analysis_modules/photon_simulator/api.py
@@ -10,12 +10,10 @@
# The full license is in the file COPYING.txt, distributed with this software.
#-----------------------------------------------------------------------------
-from numpy import VisibleDeprecationWarning
+from yt.funcs import issue_deprecation_warning
-import warnings
-warnings.warn("The photon_simulator module is deprecated. "
- "Please use pyXSIM (http://hea-www.cfa.harvard.edu/~jzuhone/pyxsim) instead.",
- VisibleDeprecationWarning, stacklevel=2)
+issue_deprecation_warning("The photon_simulator module is deprecated. Please use pyXSIM "
+ "(http://hea-www.cfa.harvard.edu/~jzuhone/pyxsim) instead.")
from .photon_models import \
PhotonModel, \
diff -r 83d2f4c6b6608ea701acf74d326c44125ae93b5e -r 959db1b322f8836bab7a600487819d8cef23dd13 yt/analysis_modules/photon_simulator/photon_simulator.py
--- a/yt/analysis_modules/photon_simulator/photon_simulator.py
+++ b/yt/analysis_modules/photon_simulator/photon_simulator.py
@@ -31,7 +31,7 @@
from yt.utilities.physical_constants import clight
from yt.utilities.cosmology import Cosmology
from yt.utilities.orientation import Orientation
-from yt.utilities.fits_image import assert_same_wcs
+from yt.visualization.fits_image import assert_same_wcs
from yt.utilities.parallel_tools.parallel_analysis_interface import \
communication_system, parallel_root_only, get_mpi_type, \
parallel_capable
diff -r 83d2f4c6b6608ea701acf74d326c44125ae93b5e -r 959db1b322f8836bab7a600487819d8cef23dd13 yt/analysis_modules/ppv_cube/ppv_cube.py
--- a/yt/analysis_modules/ppv_cube/ppv_cube.py
+++ b/yt/analysis_modules/ppv_cube/ppv_cube.py
@@ -13,7 +13,7 @@
import numpy as np
from yt.utilities.on_demand_imports import _astropy
from yt.utilities.orientation import Orientation
-from yt.utilities.fits_image import FITSImageData, sanitize_fits_unit
+from yt.visualization.fits_image import FITSImageData, sanitize_fits_unit
from yt.visualization.volume_rendering.off_axis_projection import off_axis_projection
from yt.funcs import get_pbar
from yt.utilities.physical_constants import clight, mh
diff -r 83d2f4c6b6608ea701acf74d326c44125ae93b5e -r 959db1b322f8836bab7a600487819d8cef23dd13 yt/analysis_modules/sunyaev_zeldovich/projection.py
--- a/yt/analysis_modules/sunyaev_zeldovich/projection.py
+++ b/yt/analysis_modules/sunyaev_zeldovich/projection.py
@@ -375,7 +375,7 @@
>>> sky_center = (30., 45., "deg")
>>> szprj.write_fits("SZbullet.fits", sky_center=sky_center, sky_scale=sky_scale)
"""
- from yt.utilities.fits_image import FITSImageData
+ from yt.visualization.fits_image import FITSImageData
dx = self.dx.in_units("kpc")
dy = dx
diff -r 83d2f4c6b6608ea701acf74d326c44125ae93b5e -r 959db1b322f8836bab7a600487819d8cef23dd13 yt/frontends/fits/misc.py
--- a/yt/frontends/fits/misc.py
+++ b/yt/frontends/fits/misc.py
@@ -17,7 +17,6 @@
from yt.funcs import mylog
from yt.utilities.on_demand_imports import _astropy
from yt.units.yt_array import YTQuantity, YTArray
-from yt.utilities.fits_image import FITSImageData
if PY3:
from io import BytesIO as IO
else:
@@ -109,6 +108,7 @@
... nan_mask=0.0)
"""
from spectral_cube import SpectralCube
+ from yt.visualization.fits_image import FITSImageData
from yt.frontends.fits.api import FITSDataset
cube = SpectralCube.read(filename)
if not isinstance(slab_width, YTQuantity):
diff -r 83d2f4c6b6608ea701acf74d326c44125ae93b5e -r 959db1b322f8836bab7a600487819d8cef23dd13 yt/funcs.py
--- a/yt/funcs.py
+++ b/yt/funcs.py
@@ -1036,3 +1036,7 @@
return val.decode('utf8')
else:
return val
+
+def issue_deprecation_warning(msg):
+ from numpy import VisibleDeprecationWarning
+ warnings.warn(msg, VisibleDeprecationWarning, stacklevel=3)
diff -r 83d2f4c6b6608ea701acf74d326c44125ae93b5e -r 959db1b322f8836bab7a600487819d8cef23dd13 yt/utilities/fits_image.py
--- a/yt/utilities/fits_image.py
+++ b/yt/utilities/fits_image.py
@@ -1,793 +1,9 @@
-"""
-FITSImageData Class
-"""
-
-#-----------------------------------------------------------------------------
-# Copyright (c) 2013, yt Development Team.
-#
-# Distributed under the terms of the Modified BSD License.
-#
-# The full license is in the file COPYING.txt, distributed with this software.
-#-----------------------------------------------------------------------------
-from yt.extern.six import string_types
-import numpy as np
-from yt.funcs import mylog, iterable, fix_axis, ensure_list
-from yt.visualization.fixed_resolution import FixedResolutionBuffer
-from yt.data_objects.construction_data_containers import YTCoveringGrid
-from yt.utilities.on_demand_imports import _astropy
-from yt.units.yt_array import YTQuantity, YTArray
-from yt.units import dimensions
-from yt.utilities.parallel_tools.parallel_analysis_interface import \
- parallel_root_only
-from yt.visualization.volume_rendering.off_axis_projection import off_axis_projection
-import re
-
-class FITSImageData(object):
-
- def __init__(self, data, fields=None, units=None, width=None, wcs=None):
- r""" Initialize a FITSImageData object.
-
- FITSImageData contains a collection of FITS ImageHDU instances and
- WCS information, along with units for each of the images. FITSImageData
- instances can be constructed from ImageArrays, NumPy arrays, dicts
- of such arrays, FixedResolutionBuffers, and YTCoveringGrids. The latter
- two are the most powerful because WCS information can be constructed
- automatically from their coordinates.
-
- Parameters
- ----------
- data : FixedResolutionBuffer or a YTCoveringGrid. Or, an
- ImageArray, an numpy.ndarray, or dict of such arrays
- The data to be made into a FITS image or images.
- fields : single string or list of strings, optional
- The field names for the data. If *fields* is none and *data* has
- keys, it will use these for the fields. If *data* is just a
- single array one field name must be specified.
- units : string
- The units of the WCS coordinates. Defaults to "cm".
- width : float or YTQuantity
- The width of the image. Either a single value or iterable of values.
- If a float, assumed to be in *units*. Only used if this information
- is not already provided by *data*.
- wcs : `astropy.wcs.WCS` instance, optional
- Supply an AstroPy WCS instance. Will override automatic WCS
- creation from FixedResolutionBuffers and YTCoveringGrids.
-
- Examples
- --------
-
- >>> # This example uses a FRB.
- >>> ds = load("sloshing_nomag2_hdf5_plt_cnt_0150")
- >>> prj = ds.proj(2, "kT", weight_field="density")
- >>> frb = prj.to_frb((0.5, "Mpc"), 800)
- >>> # This example just uses the FRB and puts the coords in kpc.
- >>> f_kpc = FITSImageData(frb, fields="kT", units="kpc")
- >>> # This example specifies a specific WCS.
- >>> from astropy.wcs import WCS
- >>> w = WCS(naxis=self.dimensionality)
- >>> w.wcs.crval = [30., 45.] # RA, Dec in degrees
- >>> w.wcs.cunit = ["deg"]*2
- >>> nx, ny = 800, 800
- >>> w.wcs.crpix = [0.5*(nx+1), 0.5*(ny+1)]
- >>> w.wcs.ctype = ["RA---TAN","DEC--TAN"]
- >>> scale = 1./3600. # One arcsec per pixel
- >>> w.wcs.cdelt = [-scale, scale]
- >>> f_deg = FITSImageData(frb, fields="kT", wcs=w)
- >>> f_deg.writeto("temp.fits")
- """
-
- if units is None:
- units = "cm"
- if width is None:
- width = 1.0
-
- exclude_fields = ['x','y','z','px','py','pz',
- 'pdx','pdy','pdz','weight_field']
-
- self.hdulist = _astropy.pyfits.HDUList()
-
- if isinstance(fields, string_types):
- fields = [fields]
-
- if hasattr(data, 'keys'):
- img_data = data
- if fields is None:
- fields = list(img_data.keys())
- elif isinstance(data, np.ndarray):
- if fields is None:
- mylog.warning("No field name given for this array. Calling it 'image_data'.")
- fn = 'image_data'
- fields = [fn]
- else:
- fn = fields[0]
- img_data = {fn: data}
-
- self.fields = []
- for fd in fields:
- if isinstance(fd, tuple):
- self.fields.append(fd[1])
- else:
- self.fields.append(fd)
-
- first = True
- self.field_units = {}
- for key in fields:
- if key not in exclude_fields:
- if hasattr(img_data[key], "units"):
- self.field_units[key] = img_data[key].units
- else:
- self.field_units[key] = "dimensionless"
- mylog.info("Making a FITS image of field %s" % key)
- if first:
- hdu = _astropy.pyfits.PrimaryHDU(np.array(img_data[key]))
- first = False
- else:
- hdu = _astropy.pyfits.ImageHDU(np.array(img_data[key]))
- hdu.name = key
- hdu.header["btype"] = key
- if hasattr(img_data[key], "units"):
- hdu.header["bunit"] = re.sub('()', '', str(img_data[key].units))
- self.hdulist.append(hdu)
-
- self.shape = self.hdulist[0].shape
- self.dimensionality = len(self.shape)
-
- if wcs is None:
- w = _astropy.pywcs.WCS(header=self.hdulist[0].header, naxis=self.dimensionality)
- if isinstance(img_data, FixedResolutionBuffer):
- # FRBs are a special case where we have coordinate
- # information, so we take advantage of this and
- # construct the WCS object
- dx = (img_data.bounds[1]-img_data.bounds[0]).in_units(units).v/self.shape[0]
- dy = (img_data.bounds[3]-img_data.bounds[2]).in_units(units).v/self.shape[1]
- xctr = 0.5*(img_data.bounds[1]+img_data.bounds[0]).in_units(units).v
- yctr = 0.5*(img_data.bounds[3]+img_data.bounds[2]).in_units(units).v
- center = [xctr, yctr]
- cdelt = [dx,dy]
- elif isinstance(img_data, YTCoveringGrid):
- cdelt = img_data.dds.in_units(units).v
- center = 0.5*(img_data.left_edge+img_data.right_edge).in_units(units).v
- else:
- # If img_data is just an array, we assume the center is the origin
- # and use the image width to determine the cell widths
- if not iterable(width):
- width = [width]*self.dimensionality
- if isinstance(width[0], YTQuantity):
- cdelt = [wh.in_units(units).v/n for wh, n in zip(width, self.shape)]
- else:
- cdelt = [float(wh)/n for wh, n in zip(width, self.shape)]
- center = [0.0]*self.dimensionality
- w.wcs.crpix = 0.5*(np.array(self.shape)+1)
- w.wcs.crval = center
- w.wcs.cdelt = cdelt
- w.wcs.ctype = ["linear"]*self.dimensionality
- w.wcs.cunit = [units]*self.dimensionality
- self.set_wcs(w)
- else:
- self.set_wcs(wcs)
-
- def set_wcs(self, wcs, wcsname=None, suffix=None):
- """
- Set the WCS coordinate information for all images
- with a WCS object *wcs*.
- """
- if wcsname is None:
- wcs.wcs.name = 'yt'
- else:
- wcs.wcs.name = wcsname
- h = wcs.to_header()
- if suffix is None:
- self.wcs = wcs
- else:
- setattr(self, "wcs"+suffix, wcs)
- for img in self.hdulist:
- for k, v in h.items():
- kk = k
- if suffix is not None:
- kk += suffix
- img.header[kk] = v
-
- def update_header(self, field, key, value):
- """
- Update the FITS header for *field* with a
- *key*, *value* pair. If *field* == "all", all
- headers will be updated.
- """
- if field == "all":
- for img in self.hdulist:
- img.header[key] = value
- else:
- if field not in self.keys():
- raise KeyError("%s not an image!" % field)
- idx = self.fields.index(field)
- self.hdulist[idx].header[key] = value
-
- def update_all_headers(self, key, value):
- mylog.warning("update_all_headers is deprecated. "+
- "Use update_header('all', key, value) instead.")
- self.update_header("all", key, value)
-
- def keys(self):
- return self.fields
-
- def has_key(self, key):
- return key in self.fields
-
- def values(self):
- return [self.hdulist[k] for k in self.fields]
-
- def items(self):
- return [(k, self.hdulist[k]) for k in self.fields]
-
- def __getitem__(self, item):
- return self.hdulist[item]
-
- def info(self):
- return self.hdulist.info()
-
- @parallel_root_only
- def writeto(self, fileobj, fields=None, clobber=False, **kwargs):
- r"""
- Write all of the fields or a subset of them to a FITS file.
-
- Parameters
- ----------
- fileobj : string
- The name of the file to write to.
- fields : list of strings, optional
- The fields to write to the file. If not specified
- all of the fields in the buffer will be written.
- clobber : boolean, optional
- Whether or not to overwrite a previously existing file.
- Default: False
- All other keyword arguments are passed to the `writeto`
- method of `astropy.io.fits.HDUList`.
- """
- if fields is None:
- hdus = self.hdulist
- else:
- hdus = _astropy.pyfits.HDUList()
- for field in fields:
- hdus.append(self.hdulist[field])
- hdus.writeto(fileobj, clobber=clobber, **kwargs)
-
- def to_glue(self, label="yt", data_collection=None):
- """
- Takes the data in the FITSImageData instance and exports it to
- Glue (http://www.glueviz.org) for interactive analysis. Optionally
- add a *label*. If you are already within the Glue environment, you
- can pass a *data_collection* object, otherwise Glue will be started.
- """
- from glue.core import DataCollection, Data
- from glue.core.coordinates import coordinates_from_header
- from glue.qt.glue_application import GlueApplication
-
- image = Data(label=label)
- image.coords = coordinates_from_header(self.wcs.to_header())
- for k,f in self.hdulist.items():
- image.add_component(f.data, k)
- if data_collection is None:
- dc = DataCollection([image])
- app = GlueApplication(dc)
- app.start()
- else:
- data_collection.append(image)
-
- def to_aplpy(self, **kwargs):
- """
- Use APLpy (http://aplpy.github.io) for plotting. Returns an
- `aplpy.FITSFigure` instance. All keyword arguments are passed to the
- `aplpy.FITSFigure` constructor.
- """
- import aplpy
- return aplpy.FITSFigure(self.hdulist, **kwargs)
-
- def get_data(self, field):
- """
- Return the data array of the image corresponding to *field*
- with units attached.
- """
- return YTArray(self.hdulist[field].data, self.field_units[field])
-
- def set_unit(self, field, units):
- """
- Set the units of *field* to *units*.
- """
- if field not in self.keys():
- raise KeyError("%s not an image!" % field)
- idx = self.fields.index(field)
- new_data = YTArray(self.hdulist[idx].data, self.field_units[field]).in_units(units)
- self.hdulist[idx].data = new_data.v
- self.hdulist[idx].header["bunit"] = units
- self.field_units[field] = units
-
- def pop(self, key):
- """
- Remove a field with name *key*
- and return it as a new FITSImageData
- instance.
- """
- if key not in self.keys():
- raise KeyError("%s not an image!" % key)
- idx = self.fields.index(key)
- im = self.hdulist.pop(idx)
- data = YTArray(im.data, self.field_units[key])
- self.field_units.pop(key)
- self.fields.remove(key)
- return FITSImageData(data, fields=key, wcs=self.wcs)
-
- @classmethod
- def from_file(cls, filename):
- """
- Generate a FITSImageData instance from one previously written to
- disk.
-
- Parameters
- ----------
- filename : string
- The name of the file to open.
- """
- f = _astropy.pyfits.open(filename)
- data = {}
- for hdu in f:
- data[hdu.header["btype"]] = YTArray(hdu.data, hdu.header["bunit"])
- f.close()
- return cls(data, wcs=_astropy.pywcs.WCS(header=hdu.header))
-
- @classmethod
- def from_images(cls, image_list):
- """
- Generate a new FITSImageData instance from a list of FITSImageData
- instances.
-
- Parameters
- ----------
- image_list : list of FITSImageData instances
- The images to be combined.
- """
- w = image_list[0].wcs
- img_shape = image_list[0].shape
- data = {}
- for image in image_list:
- assert_same_wcs(w, image.wcs)
- if img_shape != image.shape:
- raise RuntimeError("Images do not have the same shape!")
- for key in image.keys():
- data[key] = image.get_data(key)
- return cls(data, wcs=w)
-
- def create_sky_wcs(self, sky_center, sky_scale,
- ctype=["RA---TAN","DEC--TAN"],
- crota=None, cd=None, pc=None,
- wcsname="celestial",
- replace_old_wcs=True):
- """
- Takes a Cartesian WCS and converts it to one in a
- celestial coordinate system.
+from yt.funcs import issue_deprecation_warning
- Parameters
- ----------
- sky_center : iterable of floats
- Reference coordinates of the WCS in degrees.
- sky_scale : tuple or YTQuantity
- Conversion between an angle unit and a length unit,
- e.g. (3.0, "arcsec/kpc")
- ctype : list of strings, optional
- The type of the coordinate system to create.
- crota : 2-element ndarray, optional
- Rotation angles between cartesian coordinates and
- the celestial coordinates.
- cd : 2x2-element ndarray, optional
- Dimensioned coordinate transformation matrix.
- pc : 2x2-element ndarray, optional
- Coordinate transformation matrix.
- replace_old_wcs : boolean, optional
- Whether or not to overwrite the default WCS of the
- FITSImageData instance. If false, a second WCS will
- be added to the header. Default: True.
- """
- old_wcs = self.wcs
- naxis = old_wcs.naxis
- crval = [sky_center[0], sky_center[1]]
- if isinstance(sky_scale, YTQuantity):
- scaleq = sky_scale
- else:
- scaleq = YTQuantity(sky_scale[0],sky_scale[1])
- if scaleq.units.dimensions != dimensions.angle/dimensions.length:
- raise RuntimeError("sky_scale %s not in correct dimensions of angle/length!" % sky_scale)
- deltas = old_wcs.wcs.cdelt
- units = [str(unit) for unit in old_wcs.wcs.cunit]
- new_dx = (YTQuantity(-deltas[0], units[0])*scaleq).in_units("deg")
- new_dy = (YTQuantity(deltas[1], units[1])*scaleq).in_units("deg")
- new_wcs = _astropy.pywcs.WCS(naxis=naxis)
- cdelt = [new_dx.v, new_dy.v]
- cunit = ["deg"]*2
- if naxis == 3:
- crval.append(old_wcs.wcs.crval[2])
- cdelt.append(old_wcs.wcs.cdelt[2])
- ctype.append(old_wcs.wcs.ctype[2])
- cunit.append(old_wcs.wcs.cunit[2])
- new_wcs.wcs.crpix = old_wcs.wcs.crpix
- new_wcs.wcs.cdelt = cdelt
- new_wcs.wcs.crval = crval
- new_wcs.wcs.cunit = cunit
- new_wcs.wcs.ctype = ctype
- if crota is not None:
- new_wcs.wcs.crota = crota
- if cd is not None:
- new_wcs.wcs.cd = cd
- if pc is not None:
- new_wcs.wcs.cd = pc
- if replace_old_wcs:
- self.set_wcs(new_wcs, wcsname=wcsname)
- else:
- self.set_wcs(new_wcs, wcsname=wcsname, suffix="a")
-
-class FITSImageBuffer(FITSImageData):
- pass
-
-def sanitize_fits_unit(unit):
- if unit == "Mpc":
- mylog.info("Changing FITS file unit to kpc.")
- unit = "kpc"
- elif unit == "au":
- unit = "AU"
- return unit
-
-axis_wcs = [[1,2],[0,2],[0,1]]
-
-def construct_image(ds, axis, data_source, center, width=None, image_res=None):
- if width is None:
- width = ds.domain_width[axis_wcs[axis]]
- unit = ds.get_smallest_appropriate_unit(width[0])
- mylog.info("Making an image of the entire domain, "+
- "so setting the center to the domain center.")
- else:
- width = ds.coordinates.sanitize_width(axis, width, None)
- unit = str(width[0].units)
- if image_res is None:
- ddims = ds.domain_dimensions*ds.refine_by**ds.index.max_level
- if iterable(axis):
- nx = ddims.max()
- ny = ddims.max()
- else:
- nx, ny = [ddims[idx] for idx in axis_wcs[axis]]
- else:
- if iterable(image_res):
- nx, ny = image_res
- else:
- nx, ny = image_res, image_res
- dx = width[0]/nx
- crpix = [0.5*(nx+1), 0.5*(ny+1)]
- if hasattr(ds, "wcs") and not iterable(axis):
- # This is a FITS dataset, so we use it to construct the WCS
- cunit = [str(ds.wcs.wcs.cunit[idx]) for idx in axis_wcs[axis]]
- ctype = [ds.wcs.wcs.ctype[idx] for idx in axis_wcs[axis]]
- cdelt = [ds.wcs.wcs.cdelt[idx] for idx in axis_wcs[axis]]
- ctr_pix = center.in_units("code_length")[:ds.dimensionality].v
- crval = ds.wcs.wcs_pix2world(ctr_pix.reshape(1, ds.dimensionality))[0]
- crval = [crval[idx] for idx in axis_wcs[axis]]
- else:
- if unit == "unitary":
- unit = ds.get_smallest_appropriate_unit(ds.domain_width.max())
- elif unit == "code_length":
- unit = ds.get_smallest_appropriate_unit(ds.quan(1.0,"code_length"))
- unit = sanitize_fits_unit(unit)
- cunit = [unit]*2
- ctype = ["LINEAR"]*2
- cdelt = [dx.in_units(unit)]*2
- if iterable(axis):
- crval = center.in_units(unit)
- else:
- crval = [center[idx].in_units(unit) for idx in axis_wcs[axis]]
- if hasattr(data_source, 'to_frb'):
- if iterable(axis):
- frb = data_source.to_frb(width[0], (nx, ny), height=width[1])
- else:
- frb = data_source.to_frb(width[0], (nx, ny), center=center, height=width[1])
- else:
- frb = None
- w = _astropy.pywcs.WCS(naxis=2)
- w.wcs.crpix = crpix
- w.wcs.cdelt = cdelt
- w.wcs.crval = crval
- w.wcs.cunit = cunit
- w.wcs.ctype = ctype
- return w, frb
-
-def assert_same_wcs(wcs1, wcs2):
- from numpy.testing import assert_allclose
- assert wcs1.naxis == wcs2.naxis
- for i in range(wcs1.naxis):
- assert wcs1.wcs.cunit[i] == wcs2.wcs.cunit[i]
- assert wcs1.wcs.ctype[i] == wcs2.wcs.ctype[i]
- assert_allclose(wcs1.wcs.crpix, wcs2.wcs.crpix)
- assert_allclose(wcs1.wcs.cdelt, wcs2.wcs.cdelt)
- assert_allclose(wcs1.wcs.crval, wcs2.wcs.crval)
- crota1 = getattr(wcs1.wcs, "crota", None)
- crota2 = getattr(wcs2.wcs, "crota", None)
- if crota1 is None or crota2 is None:
- assert crota1 == crota2
- else:
- assert_allclose(wcs1.wcs.crota, wcs2.wcs.crota)
- cd1 = getattr(wcs1.wcs, "cd", None)
- cd2 = getattr(wcs2.wcs, "cd", None)
- if cd1 is None or cd2 is None:
- assert cd1 == cd2
- else:
- assert_allclose(wcs1.wcs.cd, wcs2.wcs.cd)
- pc1 = getattr(wcs1.wcs, "pc", None)
- pc2 = getattr(wcs2.wcs, "pc", None)
- if pc1 is None or pc2 is None:
- assert pc1 == pc2
- else:
- assert_allclose(wcs1.wcs.pc, wcs2.wcs.pc)
-
-class FITSSlice(FITSImageData):
- r"""
- Generate a FITSImageData of an on-axis slice.
-
- Parameters
- ----------
- ds : :class:`yt.data_objects.api.Dataset`
- The dataset object.
- axis : character or integer
- The axis of the slice. One of "x","y","z", or 0,1,2.
- fields : string or list of strings
- The fields to slice
- center : A sequence of floats, a string, or a tuple.
- The coordinate of the center of the image. If set to 'c', 'center' or
- left blank, the plot is centered on the middle of the domain. If set to
- 'max' or 'm', the center will be located at the maximum of the
- ('gas', 'density') field. Centering on the max or min of a specific
- field is supported by providing a tuple such as ("min","temperature") or
- ("max","dark_matter_density"). Units can be specified by passing in *center*
- as a tuple containing a coordinate and string unit name or by passing
- in a YTArray. If a list or unitless array is supplied, code units are
- assumed.
- width : tuple or a float.
- Width can have four different formats to support windows with variable
- x and y widths. They are:
-
- ================================== =======================
- format example
- ================================== =======================
- (float, string) (10,'kpc')
- ((float, string), (float, string)) ((10,'kpc'),(15,'kpc'))
- float 0.2
- (float, float) (0.2, 0.3)
- ================================== =======================
-
- For example, (10, 'kpc') requests a plot window that is 10 kiloparsecs
- wide in the x and y directions, ((10,'kpc'),(15,'kpc')) requests a
- window that is 10 kiloparsecs wide along the x axis and 15
- kiloparsecs wide along the y axis. In the other two examples, code
- units are assumed, for example (0.2, 0.3) requests a plot that has an
- x width of 0.2 and a y width of 0.3 in code units. If units are
- provided the resulting plot axis labels will use the supplied units.
- image_res : an int or 2-tuple of ints
- Specify the resolution of the resulting image. If not provided, it will be
- determined based on the minimum cell size of the dataset.
- """
- def __init__(self, ds, axis, fields, center="c", width=None, image_res=None, **kwargs):
- fields = ensure_list(fields)
- axis = fix_axis(axis, ds)
- center, dcenter = ds.coordinates.sanitize_center(center, axis)
- slc = ds.slice(axis, center[axis], **kwargs)
- w, frb = construct_image(ds, axis, slc, dcenter, width=width, image_res=image_res)
- super(FITSSlice, self).__init__(frb, fields=fields, wcs=w)
-
+issue_deprecation_warning("The fits_image classes have been moved to yt.visualization "
+ "and can be imported from yt directly. yt.utilities.fits_image "
+ "is deprecated.")
-class FITSProjection(FITSImageData):
- r"""
- Generate a FITSImageData of an on-axis projection.
-
- Parameters
- ----------
- ds : :class:`yt.data_objects.api.Dataset`
- The dataset object.
- axis : character or integer
- The axis along which to project. One of "x","y","z", or 0,1,2.
- fields : string or list of strings
- The fields to project
- weight_field : string
- The field used to weight the projection.
- center : A sequence of floats, a string, or a tuple.
- The coordinate of the center of the image. If set to 'c', 'center' or
- left blank, the plot is centered on the middle of the domain. If set to
- 'max' or 'm', the center will be located at the maximum of the
- ('gas', 'density') field. Centering on the max or min of a specific
- field is supported by providing a tuple such as ("min","temperature") or
- ("max","dark_matter_density"). Units can be specified by passing in *center*
- as a tuple containing a coordinate and string unit name or by passing
- in a YTArray. If a list or unitless array is supplied, code units are
- assumed.
- width : tuple or a float.
- Width can have four different formats to support windows with variable
- x and y widths. They are:
-
- ================================== =======================
- format example
- ================================== =======================
- (float, string) (10,'kpc')
- ((float, string), (float, string)) ((10,'kpc'),(15,'kpc'))
- float 0.2
- (float, float) (0.2, 0.3)
- ================================== =======================
-
- For example, (10, 'kpc') requests a plot window that is 10 kiloparsecs
- wide in the x and y directions, ((10,'kpc'),(15,'kpc')) requests a
- window that is 10 kiloparsecs wide along the x axis and 15
- kiloparsecs wide along the y axis. In the other two examples, code
- units are assumed, for example (0.2, 0.3) requests a plot that has an
- x width of 0.2 and a y width of 0.3 in code units. If units are
- provided the resulting plot axis labels will use the supplied units.
- image_res : an int or 2-tuple of ints
- Specify the resolution of the resulting image. If not provided, it will be
- determined based on the minimum cell size of the dataset.
- """
- def __init__(self, ds, axis, fields, center="c", width=None,
- weight_field=None, image_res=None, **kwargs):
- fields = ensure_list(fields)
- axis = fix_axis(axis, ds)
- center, dcenter = ds.coordinates.sanitize_center(center, axis)
- prj = ds.proj(fields[0], axis, weight_field=weight_field, **kwargs)
- w, frb = construct_image(ds, axis, prj, dcenter, width=width, image_res=image_res)
- super(FITSProjection, self).__init__(frb, fields=fields, wcs=w)
-
-class FITSOffAxisSlice(FITSImageData):
- r"""
- Generate a FITSImageData of an off-axis slice.
-
- Parameters
- ----------
- ds : :class:`yt.data_objects.api.Dataset`
- The dataset object.
- normal : a sequence of floats
- The vector normal to the projection plane.
- fields : string or list of strings
- The fields to slice
- center : A sequence of floats, a string, or a tuple.
- The coordinate of the center of the image. If set to 'c', 'center' or
- left blank, the plot is centered on the middle of the domain. If set to
- 'max' or 'm', the center will be located at the maximum of the
- ('gas', 'density') field. Centering on the max or min of a specific
- field is supported by providing a tuple such as ("min","temperature") or
- ("max","dark_matter_density"). Units can be specified by passing in *center*
- as a tuple containing a coordinate and string unit name or by passing
- in a YTArray. If a list or unitless array is supplied, code units are
- assumed.
- width : tuple or a float.
- Width can have four different formats to support windows with variable
- x and y widths. They are:
-
- ================================== =======================
- format example
- ================================== =======================
- (float, string) (10,'kpc')
- ((float, string), (float, string)) ((10,'kpc'),(15,'kpc'))
- float 0.2
- (float, float) (0.2, 0.3)
- ================================== =======================
-
- For example, (10, 'kpc') requests a plot window that is 10 kiloparsecs
- wide in the x and y directions, ((10,'kpc'),(15,'kpc')) requests a
- window that is 10 kiloparsecs wide along the x axis and 15
- kiloparsecs wide along the y axis. In the other two examples, code
- units are assumed, for example (0.2, 0.3) requests a plot that has an
- x width of 0.2 and a y width of 0.3 in code units. If units are
- provided the resulting plot axis labels will use the supplied units.
- image_res : an int or 2-tuple of ints
- Specify the resolution of the resulting image.
- north_vector : a sequence of floats
- A vector defining the 'up' direction in the plot. This
- option sets the orientation of the slicing plane. If not
- set, an arbitrary grid-aligned north-vector is chosen.
- """
- def __init__(self, ds, normal, fields, center='c', width=None, image_res=512,
- north_vector=None):
- fields = ensure_list(fields)
- center, dcenter = ds.coordinates.sanitize_center(center, 4)
- cut = ds.cutting(normal, center, north_vector=north_vector)
- center = ds.arr([0.0] * 2, 'code_length')
- w, frb = construct_image(ds, normal, cut, center, width=width, image_res=image_res)
- super(FITSOffAxisSlice, self).__init__(frb, fields=fields, wcs=w)
-
-
-class FITSOffAxisProjection(FITSImageData):
- r"""
- Generate a FITSImageData of an off-axis projection.
-
- Parameters
- ----------
- ds : :class:`yt.data_objects.api.Dataset`
- This is the dataset object corresponding to the
- simulation output to be plotted.
- normal : a sequence of floats
- The vector normal to the projection plane.
- fields : string, list of strings
- The name of the field(s) to be plotted.
- center : A sequence of floats, a string, or a tuple.
- The coordinate of the center of the image. If set to 'c', 'center' or
- left blank, the plot is centered on the middle of the domain. If set to
- 'max' or 'm', the center will be located at the maximum of the
- ('gas', 'density') field. Centering on the max or min of a specific
- field is supported by providing a tuple such as ("min","temperature") or
- ("max","dark_matter_density"). Units can be specified by passing in *center*
- as a tuple containing a coordinate and string unit name or by passing
- in a YTArray. If a list or unitless array is supplied, code units are
- assumed.
- width : tuple or a float.
- Width can have four different formats to support windows with variable
- x and y widths. They are:
-
- ================================== =======================
- format example
- ================================== =======================
- (float, string) (10,'kpc')
- ((float, string), (float, string)) ((10,'kpc'),(15,'kpc'))
- float 0.2
- (float, float) (0.2, 0.3)
- ================================== =======================
-
- For example, (10, 'kpc') requests a plot window that is 10 kiloparsecs
- wide in the x and y directions, ((10,'kpc'),(15,'kpc')) requests a
- window that is 10 kiloparsecs wide along the x axis and 15
- kiloparsecs wide along the y axis. In the other two examples, code
- units are assumed, for example (0.2, 0.3) requests a plot that has an
- x width of 0.2 and a y width of 0.3 in code units. If units are
- provided the resulting plot axis labels will use the supplied units.
- depth : A tuple or a float
- A tuple containing the depth to project through and the string
- key of the unit: (width, 'unit'). If set to a float, code units
- are assumed
- weight_field : string
- The name of the weighting field. Set to None for no weight.
- image_res : an int or 2-tuple of ints
- Specify the resolution of the resulting image.
- depth_res : integer
- Deprecated, this is still in the function signature for API
- compatibility
- north_vector : a sequence of floats
- A vector defining the 'up' direction in the plot. This
- option sets the orientation of the slicing plane. If not
- set, an arbitrary grid-aligned north-vector is chosen.
- method : string
- The method of projection. Valid methods are:
-
- "integrate" with no weight_field specified : integrate the requested
- field along the line of sight.
-
- "integrate" with a weight_field specified : weight the requested
- field by the weighting field and integrate along the line of sight.
-
- "sum" : This method is the same as integrate, except that it does not
- multiply by a path length when performing the integration, and is
- just a straight summation of the field along the given axis. WARNING:
- This should only be used for uniform resolution grid datasets, as other
- datasets may result in unphysical images.
- data_source : yt.data_objects.data_containers.YTSelectionContainer, optional
- If specified, this will be the data source used for selecting regions to project.
- """
- def __init__(self, ds, normal, fields, center='c', width=(1.0, 'unitary'),
- weight_field=None, image_res=512, depth_res=256, data_source=None,
- north_vector=None, depth=(1.0,"unitary"), no_ghost=False, method='integrate'):
- fields = ensure_list(fields)
- center, dcenter = ds.coordinates.sanitize_center(center, 4)
- buf = {}
- width = ds.coordinates.sanitize_width(normal, width, depth)
- wd = tuple(el.in_units('code_length').v for el in width)
- if not iterable(image_res):
- image_res = (image_res, image_res)
- res = (image_res[0], image_res[1])
- if data_source is None:
- source = ds
- else:
- source = data_source
- for field in fields:
- buf[field] = off_axis_projection(source, center, normal, wd, res, field,
- no_ghost=no_ghost, north_vector=north_vector,
- method=method, weight=weight_field).swapaxes(0, 1)
- center = ds.arr([0.0] * 2, 'code_length')
- w, not_an_frb = construct_image(ds, normal, buf, center, width=width, image_res=image_res)
- super(FITSOffAxisProjection, self).__init__(buf, fields=fields, wcs=w)
+from yt.visualization.fits_image import \
+ FITSImageData, FITSSlice, FITSProjection, \
+ FITSOffAxisSlice, FITSOffAxisProjection
\ No newline at end of file
diff -r 83d2f4c6b6608ea701acf74d326c44125ae93b5e -r 959db1b322f8836bab7a600487819d8cef23dd13 yt/utilities/tests/test_fits_image.py
--- a/yt/utilities/tests/test_fits_image.py
+++ /dev/null
@@ -1,139 +0,0 @@
-"""
-Unit test FITS image creation in yt.
-
-
-
-"""
-
-#-----------------------------------------------------------------------------
-# Copyright (c) 2013, yt Development Team.
-#
-# Distributed under the terms of the Modified BSD License.
-#
-# The full license is in the file COPYING.txt, distributed with this software.
-#-----------------------------------------------------------------------------
-
-import tempfile
-import os
-import shutil
-from yt.testing import fake_random_ds, requires_module
-from yt.convenience import load
-from numpy.testing import \
- assert_equal
-from yt.utilities.fits_image import \
- FITSImageData, FITSProjection, \
- FITSSlice, FITSOffAxisSlice, \
- FITSOffAxisProjection, \
- assert_same_wcs
-from yt.visualization.volume_rendering.off_axis_projection import \
- off_axis_projection
-
-
- at requires_module("astropy")
-def test_fits_image():
- tmpdir = tempfile.mkdtemp()
- curdir = os.getcwd()
- os.chdir(tmpdir)
-
- fields = ("density", "temperature")
- units = ('g/cm**3', 'K',)
- ds = fake_random_ds(64, fields=fields, units=units, nprocs=16,
- length_unit=100.0)
-
- prj = ds.proj("density", 2)
- prj_frb = prj.to_frb((0.5, "unitary"), 128)
-
- fid1 = FITSImageData(prj_frb, fields=["density","temperature"], units="cm")
- fits_prj = FITSProjection(ds, "z", ["density","temperature"], image_res=128,
- width=(0.5,"unitary"))
-
- yield assert_equal, fid1.get_data("density"), fits_prj.get_data("density")
- yield assert_equal, fid1.get_data("temperature"), fits_prj.get_data("temperature")
-
- fid1.writeto("fid1.fits", clobber=True)
- new_fid1 = FITSImageData.from_file("fid1.fits")
-
- yield assert_equal, fid1.get_data("density"), new_fid1.get_data("density")
- yield assert_equal, fid1.get_data("temperature"), new_fid1.get_data("temperature")
-
- ds2 = load("fid1.fits")
- ds2.index
-
- assert ("fits","density") in ds2.field_list
- assert ("fits","temperature") in ds2.field_list
-
- dw_cm = ds2.domain_width.in_units("cm")
-
- assert dw_cm[0].v == 50.
- assert dw_cm[1].v == 50.
-
- slc = ds.slice(2, 0.5)
- slc_frb = slc.to_frb((0.5, "unitary"), 128)
-
- fid2 = FITSImageData(slc_frb, fields=["density","temperature"], units="cm")
- fits_slc = FITSSlice(ds, "z", ["density","temperature"], image_res=128,
- width=(0.5,"unitary"))
-
- yield assert_equal, fid2.get_data("density"), fits_slc.get_data("density")
- yield assert_equal, fid2.get_data("temperature"), fits_slc.get_data("temperature")
-
- dens_img = fid2.pop("density")
- temp_img = fid2.pop("temperature")
-
- # This already has some assertions in it, so we don't need to do anything
- # with it other than just make one
- FITSImageData.from_images([dens_img, temp_img])
-
- cut = ds.cutting([0.1, 0.2, -0.9], [0.5, 0.42, 0.6])
- cut_frb = cut.to_frb((0.5, "unitary"), 128)
-
- fid3 = FITSImageData(cut_frb, fields=["density","temperature"], units="cm")
- fits_cut = FITSOffAxisSlice(ds, [0.1, 0.2, -0.9], ["density","temperature"],
- image_res=128, center=[0.5, 0.42, 0.6],
- width=(0.5,"unitary"))
-
- yield assert_equal, fid3.get_data("density"), fits_cut.get_data("density")
- yield assert_equal, fid3.get_data("temperature"), fits_cut.get_data("temperature")
-
- fid3.create_sky_wcs([30.,45.], (1.0,"arcsec/kpc"))
- fid3.writeto("fid3.fits", clobber=True)
- new_fid3 = FITSImageData.from_file("fid3.fits")
- assert_same_wcs(fid3.wcs, new_fid3.wcs)
- assert new_fid3.wcs.wcs.cunit[0] == "deg"
- assert new_fid3.wcs.wcs.cunit[1] == "deg"
- assert new_fid3.wcs.wcs.ctype[0] == "RA---TAN"
- assert new_fid3.wcs.wcs.ctype[1] == "DEC--TAN"
-
- buf = off_axis_projection(ds, ds.domain_center, [0.1, 0.2, -0.9],
- 0.5, 128, "density").swapaxes(0, 1)
- fid4 = FITSImageData(buf, fields="density", width=100.0)
- fits_oap = FITSOffAxisProjection(ds, [0.1, 0.2, -0.9], "density",
- width=(0.5,"unitary"), image_res=128,
- depth_res=128, depth=(0.5,"unitary"))
-
- yield assert_equal, fid4.get_data("density"), fits_oap.get_data("density")
-
- fid4.create_sky_wcs([30., 45.], (1.0, "arcsec/kpc"), replace_old_wcs=False)
- assert fid4.wcs.wcs.cunit[0] == "cm"
- assert fid4.wcs.wcs.cunit[1] == "cm"
- assert fid4.wcs.wcs.ctype[0] == "linear"
- assert fid4.wcs.wcs.ctype[1] == "linear"
- assert fid4.wcsa.wcs.cunit[0] == "deg"
- assert fid4.wcsa.wcs.cunit[1] == "deg"
- assert fid4.wcsa.wcs.ctype[0] == "RA---TAN"
- assert fid4.wcsa.wcs.ctype[1] == "DEC--TAN"
-
- cvg = ds.covering_grid(ds.index.max_level, [0.25,0.25,0.25],
- [32, 32, 32], fields=["density","temperature"])
- fid5 = FITSImageData(cvg, fields=["density","temperature"])
- assert fid5.dimensionality == 3
-
- fid5.update_header("density", "time", 0.1)
- fid5.update_header("all", "units", "cgs")
-
- assert fid5["density"].header["time"] == 0.1
- assert fid5["temperature"].header["units"] == "cgs"
- assert fid5["density"].header["units"] == "cgs"
-
- os.chdir(curdir)
- shutil.rmtree(tmpdir)
diff -r 83d2f4c6b6608ea701acf74d326c44125ae93b5e -r 959db1b322f8836bab7a600487819d8cef23dd13 yt/visualization/api.py
--- a/yt/visualization/api.py
+++ b/yt/visualization/api.py
@@ -59,3 +59,10 @@
from .base_plot_types import \
get_multi_plot
+
+from .fits_image import \
+ FITSImageData, \
+ FITSSlice, \
+ FITSOffAxisSlice, \
+ FITSProjection, \
+ FITSOffAxisProjection
diff -r 83d2f4c6b6608ea701acf74d326c44125ae93b5e -r 959db1b322f8836bab7a600487819d8cef23dd13 yt/visualization/fits_image.py
--- /dev/null
+++ b/yt/visualization/fits_image.py
@@ -0,0 +1,793 @@
+"""
+FITSImageData Class
+"""
+
+#-----------------------------------------------------------------------------
+# Copyright (c) 2013, yt Development Team.
+#
+# Distributed under the terms of the Modified BSD License.
+#
+# The full license is in the file COPYING.txt, distributed with this software.
+#-----------------------------------------------------------------------------
+from yt.extern.six import string_types
+import numpy as np
+from yt.funcs import mylog, iterable, fix_axis, ensure_list
+from yt.visualization.fixed_resolution import FixedResolutionBuffer
+from yt.data_objects.construction_data_containers import YTCoveringGrid
+from yt.utilities.on_demand_imports import _astropy
+from yt.units.yt_array import YTQuantity, YTArray
+from yt.units import dimensions
+from yt.utilities.parallel_tools.parallel_analysis_interface import \
+ parallel_root_only
+from yt.visualization.volume_rendering.off_axis_projection import off_axis_projection
+import re
+
+class FITSImageData(object):
+
+ def __init__(self, data, fields=None, units=None, width=None, wcs=None):
+ r""" Initialize a FITSImageData object.
+
+ FITSImageData contains a collection of FITS ImageHDU instances and
+ WCS information, along with units for each of the images. FITSImageData
+ instances can be constructed from ImageArrays, NumPy arrays, dicts
+ of such arrays, FixedResolutionBuffers, and YTCoveringGrids. The latter
+ two are the most powerful because WCS information can be constructed
+ automatically from their coordinates.
+
+ Parameters
+ ----------
+ data : FixedResolutionBuffer or a YTCoveringGrid. Or, an
+ ImageArray, an numpy.ndarray, or dict of such arrays
+ The data to be made into a FITS image or images.
+ fields : single string or list of strings, optional
+ The field names for the data. If *fields* is none and *data* has
+ keys, it will use these for the fields. If *data* is just a
+ single array one field name must be specified.
+ units : string
+ The units of the WCS coordinates. Defaults to "cm".
+ width : float or YTQuantity
+ The width of the image. Either a single value or iterable of values.
+ If a float, assumed to be in *units*. Only used if this information
+ is not already provided by *data*.
+ wcs : `astropy.wcs.WCS` instance, optional
+ Supply an AstroPy WCS instance. Will override automatic WCS
+ creation from FixedResolutionBuffers and YTCoveringGrids.
+
+ Examples
+ --------
+
+ >>> # This example uses a FRB.
+ >>> ds = load("sloshing_nomag2_hdf5_plt_cnt_0150")
+ >>> prj = ds.proj(2, "kT", weight_field="density")
+ >>> frb = prj.to_frb((0.5, "Mpc"), 800)
+ >>> # This example just uses the FRB and puts the coords in kpc.
+ >>> f_kpc = FITSImageData(frb, fields="kT", units="kpc")
+ >>> # This example specifies a specific WCS.
+ >>> from astropy.wcs import WCS
+ >>> w = WCS(naxis=self.dimensionality)
+ >>> w.wcs.crval = [30., 45.] # RA, Dec in degrees
+ >>> w.wcs.cunit = ["deg"]*2
+ >>> nx, ny = 800, 800
+ >>> w.wcs.crpix = [0.5*(nx+1), 0.5*(ny+1)]
+ >>> w.wcs.ctype = ["RA---TAN","DEC--TAN"]
+ >>> scale = 1./3600. # One arcsec per pixel
+ >>> w.wcs.cdelt = [-scale, scale]
+ >>> f_deg = FITSImageData(frb, fields="kT", wcs=w)
+ >>> f_deg.writeto("temp.fits")
+ """
+
+ if units is None:
+ units = "cm"
+ if width is None:
+ width = 1.0
+
+ exclude_fields = ['x','y','z','px','py','pz',
+ 'pdx','pdy','pdz','weight_field']
+
+ self.hdulist = _astropy.pyfits.HDUList()
+
+ if isinstance(fields, string_types):
+ fields = [fields]
+
+ if hasattr(data, 'keys'):
+ img_data = data
+ if fields is None:
+ fields = list(img_data.keys())
+ elif isinstance(data, np.ndarray):
+ if fields is None:
+ mylog.warning("No field name given for this array. Calling it 'image_data'.")
+ fn = 'image_data'
+ fields = [fn]
+ else:
+ fn = fields[0]
+ img_data = {fn: data}
+
+ self.fields = []
+ for fd in fields:
+ if isinstance(fd, tuple):
+ self.fields.append(fd[1])
+ else:
+ self.fields.append(fd)
+
+ first = True
+ self.field_units = {}
+ for key in fields:
+ if key not in exclude_fields:
+ if hasattr(img_data[key], "units"):
+ self.field_units[key] = img_data[key].units
+ else:
+ self.field_units[key] = "dimensionless"
+ mylog.info("Making a FITS image of field %s" % key)
+ if first:
+ hdu = _astropy.pyfits.PrimaryHDU(np.array(img_data[key]))
+ first = False
+ else:
+ hdu = _astropy.pyfits.ImageHDU(np.array(img_data[key]))
+ hdu.name = key
+ hdu.header["btype"] = key
+ if hasattr(img_data[key], "units"):
+ hdu.header["bunit"] = re.sub('()', '', str(img_data[key].units))
+ self.hdulist.append(hdu)
+
+ self.shape = self.hdulist[0].shape
+ self.dimensionality = len(self.shape)
+
+ if wcs is None:
+ w = _astropy.pywcs.WCS(header=self.hdulist[0].header, naxis=self.dimensionality)
+ if isinstance(img_data, FixedResolutionBuffer):
+ # FRBs are a special case where we have coordinate
+ # information, so we take advantage of this and
+ # construct the WCS object
+ dx = (img_data.bounds[1]-img_data.bounds[0]).in_units(units).v/self.shape[0]
+ dy = (img_data.bounds[3]-img_data.bounds[2]).in_units(units).v/self.shape[1]
+ xctr = 0.5*(img_data.bounds[1]+img_data.bounds[0]).in_units(units).v
+ yctr = 0.5*(img_data.bounds[3]+img_data.bounds[2]).in_units(units).v
+ center = [xctr, yctr]
+ cdelt = [dx,dy]
+ elif isinstance(img_data, YTCoveringGrid):
+ cdelt = img_data.dds.in_units(units).v
+ center = 0.5*(img_data.left_edge+img_data.right_edge).in_units(units).v
+ else:
+ # If img_data is just an array, we assume the center is the origin
+ # and use the image width to determine the cell widths
+ if not iterable(width):
+ width = [width]*self.dimensionality
+ if isinstance(width[0], YTQuantity):
+ cdelt = [wh.in_units(units).v/n for wh, n in zip(width, self.shape)]
+ else:
+ cdelt = [float(wh)/n for wh, n in zip(width, self.shape)]
+ center = [0.0]*self.dimensionality
+ w.wcs.crpix = 0.5*(np.array(self.shape)+1)
+ w.wcs.crval = center
+ w.wcs.cdelt = cdelt
+ w.wcs.ctype = ["linear"]*self.dimensionality
+ w.wcs.cunit = [units]*self.dimensionality
+ self.set_wcs(w)
+ else:
+ self.set_wcs(wcs)
+
+ def set_wcs(self, wcs, wcsname=None, suffix=None):
+ """
+ Set the WCS coordinate information for all images
+ with a WCS object *wcs*.
+ """
+ if wcsname is None:
+ wcs.wcs.name = 'yt'
+ else:
+ wcs.wcs.name = wcsname
+ h = wcs.to_header()
+ if suffix is None:
+ self.wcs = wcs
+ else:
+ setattr(self, "wcs"+suffix, wcs)
+ for img in self.hdulist:
+ for k, v in h.items():
+ kk = k
+ if suffix is not None:
+ kk += suffix
+ img.header[kk] = v
+
+ def update_header(self, field, key, value):
+ """
+ Update the FITS header for *field* with a
+ *key*, *value* pair. If *field* == "all", all
+ headers will be updated.
+ """
+ if field == "all":
+ for img in self.hdulist:
+ img.header[key] = value
+ else:
+ if field not in self.keys():
+ raise KeyError("%s not an image!" % field)
+ idx = self.fields.index(field)
+ self.hdulist[idx].header[key] = value
+
+ def update_all_headers(self, key, value):
+ mylog.warning("update_all_headers is deprecated. "+
+ "Use update_header('all', key, value) instead.")
+ self.update_header("all", key, value)
+
+ def keys(self):
+ return self.fields
+
+ def has_key(self, key):
+ return key in self.fields
+
+ def values(self):
+ return [self.hdulist[k] for k in self.fields]
+
+ def items(self):
+ return [(k, self.hdulist[k]) for k in self.fields]
+
+ def __getitem__(self, item):
+ return self.hdulist[item]
+
+ def info(self):
+ return self.hdulist.info()
+
+ @parallel_root_only
+ def writeto(self, fileobj, fields=None, clobber=False, **kwargs):
+ r"""
+ Write all of the fields or a subset of them to a FITS file.
+
+ Parameters
+ ----------
+ fileobj : string
+ The name of the file to write to.
+ fields : list of strings, optional
+ The fields to write to the file. If not specified
+ all of the fields in the buffer will be written.
+ clobber : boolean, optional
+ Whether or not to overwrite a previously existing file.
+ Default: False
+ All other keyword arguments are passed to the `writeto`
+ method of `astropy.io.fits.HDUList`.
+ """
+ if fields is None:
+ hdus = self.hdulist
+ else:
+ hdus = _astropy.pyfits.HDUList()
+ for field in fields:
+ hdus.append(self.hdulist[field])
+ hdus.writeto(fileobj, clobber=clobber, **kwargs)
+
+ def to_glue(self, label="yt", data_collection=None):
+ """
+ Takes the data in the FITSImageData instance and exports it to
+ Glue (http://www.glueviz.org) for interactive analysis. Optionally
+ add a *label*. If you are already within the Glue environment, you
+ can pass a *data_collection* object, otherwise Glue will be started.
+ """
+ from glue.core import DataCollection, Data
+ from glue.core.coordinates import coordinates_from_header
+ from glue.qt.glue_application import GlueApplication
+
+ image = Data(label=label)
+ image.coords = coordinates_from_header(self.wcs.to_header())
+ for k,f in self.hdulist.items():
+ image.add_component(f.data, k)
+ if data_collection is None:
+ dc = DataCollection([image])
+ app = GlueApplication(dc)
+ app.start()
+ else:
+ data_collection.append(image)
+
+ def to_aplpy(self, **kwargs):
+ """
+ Use APLpy (http://aplpy.github.io) for plotting. Returns an
+ `aplpy.FITSFigure` instance. All keyword arguments are passed to the
+ `aplpy.FITSFigure` constructor.
+ """
+ import aplpy
+ return aplpy.FITSFigure(self.hdulist, **kwargs)
+
+ def get_data(self, field):
+ """
+ Return the data array of the image corresponding to *field*
+ with units attached.
+ """
+ return YTArray(self.hdulist[field].data, self.field_units[field])
+
+ def set_unit(self, field, units):
+ """
+ Set the units of *field* to *units*.
+ """
+ if field not in self.keys():
+ raise KeyError("%s not an image!" % field)
+ idx = self.fields.index(field)
+ new_data = YTArray(self.hdulist[idx].data, self.field_units[field]).in_units(units)
+ self.hdulist[idx].data = new_data.v
+ self.hdulist[idx].header["bunit"] = units
+ self.field_units[field] = units
+
+ def pop(self, key):
+ """
+ Remove a field with name *key*
+ and return it as a new FITSImageData
+ instance.
+ """
+ if key not in self.keys():
+ raise KeyError("%s not an image!" % key)
+ idx = self.fields.index(key)
+ im = self.hdulist.pop(idx)
+ data = YTArray(im.data, self.field_units[key])
+ self.field_units.pop(key)
+ self.fields.remove(key)
+ return FITSImageData(data, fields=key, wcs=self.wcs)
+
+ @classmethod
+ def from_file(cls, filename):
+ """
+ Generate a FITSImageData instance from one previously written to
+ disk.
+
+ Parameters
+ ----------
+ filename : string
+ The name of the file to open.
+ """
+ f = _astropy.pyfits.open(filename)
+ data = {}
+ for hdu in f:
+ data[hdu.header["btype"]] = YTArray(hdu.data, hdu.header["bunit"])
+ f.close()
+ return cls(data, wcs=_astropy.pywcs.WCS(header=hdu.header))
+
+ @classmethod
+ def from_images(cls, image_list):
+ """
+ Generate a new FITSImageData instance from a list of FITSImageData
+ instances.
+
+ Parameters
+ ----------
+ image_list : list of FITSImageData instances
+ The images to be combined.
+ """
+ w = image_list[0].wcs
+ img_shape = image_list[0].shape
+ data = {}
+ for image in image_list:
+ assert_same_wcs(w, image.wcs)
+ if img_shape != image.shape:
+ raise RuntimeError("Images do not have the same shape!")
+ for key in image.keys():
+ data[key] = image.get_data(key)
+ return cls(data, wcs=w)
+
+ def create_sky_wcs(self, sky_center, sky_scale,
+ ctype=["RA---TAN","DEC--TAN"],
+ crota=None, cd=None, pc=None,
+ wcsname="celestial",
+ replace_old_wcs=True):
+ """
+ Takes a Cartesian WCS and converts it to one in a
+ celestial coordinate system.
+
+ Parameters
+ ----------
+ sky_center : iterable of floats
+ Reference coordinates of the WCS in degrees.
+ sky_scale : tuple or YTQuantity
+ Conversion between an angle unit and a length unit,
+ e.g. (3.0, "arcsec/kpc")
+ ctype : list of strings, optional
+ The type of the coordinate system to create.
+ crota : 2-element ndarray, optional
+ Rotation angles between cartesian coordinates and
+ the celestial coordinates.
+ cd : 2x2-element ndarray, optional
+ Dimensioned coordinate transformation matrix.
+ pc : 2x2-element ndarray, optional
+ Coordinate transformation matrix.
+ replace_old_wcs : boolean, optional
+ Whether or not to overwrite the default WCS of the
+ FITSImageData instance. If false, a second WCS will
+ be added to the header. Default: True.
+ """
+ old_wcs = self.wcs
+ naxis = old_wcs.naxis
+ crval = [sky_center[0], sky_center[1]]
+ if isinstance(sky_scale, YTQuantity):
+ scaleq = sky_scale
+ else:
+ scaleq = YTQuantity(sky_scale[0],sky_scale[1])
+ if scaleq.units.dimensions != dimensions.angle/dimensions.length:
+ raise RuntimeError("sky_scale %s not in correct dimensions of angle/length!" % sky_scale)
+ deltas = old_wcs.wcs.cdelt
+ units = [str(unit) for unit in old_wcs.wcs.cunit]
+ new_dx = (YTQuantity(-deltas[0], units[0])*scaleq).in_units("deg")
+ new_dy = (YTQuantity(deltas[1], units[1])*scaleq).in_units("deg")
+ new_wcs = _astropy.pywcs.WCS(naxis=naxis)
+ cdelt = [new_dx.v, new_dy.v]
+ cunit = ["deg"]*2
+ if naxis == 3:
+ crval.append(old_wcs.wcs.crval[2])
+ cdelt.append(old_wcs.wcs.cdelt[2])
+ ctype.append(old_wcs.wcs.ctype[2])
+ cunit.append(old_wcs.wcs.cunit[2])
+ new_wcs.wcs.crpix = old_wcs.wcs.crpix
+ new_wcs.wcs.cdelt = cdelt
+ new_wcs.wcs.crval = crval
+ new_wcs.wcs.cunit = cunit
+ new_wcs.wcs.ctype = ctype
+ if crota is not None:
+ new_wcs.wcs.crota = crota
+ if cd is not None:
+ new_wcs.wcs.cd = cd
+ if pc is not None:
+ new_wcs.wcs.cd = pc
+ if replace_old_wcs:
+ self.set_wcs(new_wcs, wcsname=wcsname)
+ else:
+ self.set_wcs(new_wcs, wcsname=wcsname, suffix="a")
+
+class FITSImageBuffer(FITSImageData):
+ pass
+
+def sanitize_fits_unit(unit):
+ if unit == "Mpc":
+ mylog.info("Changing FITS file unit to kpc.")
+ unit = "kpc"
+ elif unit == "au":
+ unit = "AU"
+ return unit
+
+axis_wcs = [[1,2],[0,2],[0,1]]
+
+def construct_image(ds, axis, data_source, center, width=None, image_res=None):
+ if width is None:
+ width = ds.domain_width[axis_wcs[axis]]
+ unit = ds.get_smallest_appropriate_unit(width[0])
+ mylog.info("Making an image of the entire domain, "+
+ "so setting the center to the domain center.")
+ else:
+ width = ds.coordinates.sanitize_width(axis, width, None)
+ unit = str(width[0].units)
+ if image_res is None:
+ ddims = ds.domain_dimensions*ds.refine_by**ds.index.max_level
+ if iterable(axis):
+ nx = ddims.max()
+ ny = ddims.max()
+ else:
+ nx, ny = [ddims[idx] for idx in axis_wcs[axis]]
+ else:
+ if iterable(image_res):
+ nx, ny = image_res
+ else:
+ nx, ny = image_res, image_res
+ dx = width[0]/nx
+ crpix = [0.5*(nx+1), 0.5*(ny+1)]
+ if hasattr(ds, "wcs") and not iterable(axis):
+ # This is a FITS dataset, so we use it to construct the WCS
+ cunit = [str(ds.wcs.wcs.cunit[idx]) for idx in axis_wcs[axis]]
+ ctype = [ds.wcs.wcs.ctype[idx] for idx in axis_wcs[axis]]
+ cdelt = [ds.wcs.wcs.cdelt[idx] for idx in axis_wcs[axis]]
+ ctr_pix = center.in_units("code_length")[:ds.dimensionality].v
+ crval = ds.wcs.wcs_pix2world(ctr_pix.reshape(1, ds.dimensionality))[0]
+ crval = [crval[idx] for idx in axis_wcs[axis]]
+ else:
+ if unit == "unitary":
+ unit = ds.get_smallest_appropriate_unit(ds.domain_width.max())
+ elif unit == "code_length":
+ unit = ds.get_smallest_appropriate_unit(ds.quan(1.0,"code_length"))
+ unit = sanitize_fits_unit(unit)
+ cunit = [unit]*2
+ ctype = ["LINEAR"]*2
+ cdelt = [dx.in_units(unit)]*2
+ if iterable(axis):
+ crval = center.in_units(unit)
+ else:
+ crval = [center[idx].in_units(unit) for idx in axis_wcs[axis]]
+ if hasattr(data_source, 'to_frb'):
+ if iterable(axis):
+ frb = data_source.to_frb(width[0], (nx, ny), height=width[1])
+ else:
+ frb = data_source.to_frb(width[0], (nx, ny), center=center, height=width[1])
+ else:
+ frb = None
+ w = _astropy.pywcs.WCS(naxis=2)
+ w.wcs.crpix = crpix
+ w.wcs.cdelt = cdelt
+ w.wcs.crval = crval
+ w.wcs.cunit = cunit
+ w.wcs.ctype = ctype
+ return w, frb
+
+def assert_same_wcs(wcs1, wcs2):
+ from numpy.testing import assert_allclose
+ assert wcs1.naxis == wcs2.naxis
+ for i in range(wcs1.naxis):
+ assert wcs1.wcs.cunit[i] == wcs2.wcs.cunit[i]
+ assert wcs1.wcs.ctype[i] == wcs2.wcs.ctype[i]
+ assert_allclose(wcs1.wcs.crpix, wcs2.wcs.crpix)
+ assert_allclose(wcs1.wcs.cdelt, wcs2.wcs.cdelt)
+ assert_allclose(wcs1.wcs.crval, wcs2.wcs.crval)
+ crota1 = getattr(wcs1.wcs, "crota", None)
+ crota2 = getattr(wcs2.wcs, "crota", None)
+ if crota1 is None or crota2 is None:
+ assert crota1 == crota2
+ else:
+ assert_allclose(wcs1.wcs.crota, wcs2.wcs.crota)
+ cd1 = getattr(wcs1.wcs, "cd", None)
+ cd2 = getattr(wcs2.wcs, "cd", None)
+ if cd1 is None or cd2 is None:
+ assert cd1 == cd2
+ else:
+ assert_allclose(wcs1.wcs.cd, wcs2.wcs.cd)
+ pc1 = getattr(wcs1.wcs, "pc", None)
+ pc2 = getattr(wcs2.wcs, "pc", None)
+ if pc1 is None or pc2 is None:
+ assert pc1 == pc2
+ else:
+ assert_allclose(wcs1.wcs.pc, wcs2.wcs.pc)
+
+class FITSSlice(FITSImageData):
+ r"""
+ Generate a FITSImageData of an on-axis slice.
+
+ Parameters
+ ----------
+ ds : :class:`yt.data_objects.api.Dataset`
+ The dataset object.
+ axis : character or integer
+ The axis of the slice. One of "x","y","z", or 0,1,2.
+ fields : string or list of strings
+ The fields to slice
+ center : A sequence of floats, a string, or a tuple.
+ The coordinate of the center of the image. If set to 'c', 'center' or
+ left blank, the plot is centered on the middle of the domain. If set to
+ 'max' or 'm', the center will be located at the maximum of the
+ ('gas', 'density') field. Centering on the max or min of a specific
+ field is supported by providing a tuple such as ("min","temperature") or
+ ("max","dark_matter_density"). Units can be specified by passing in *center*
+ as a tuple containing a coordinate and string unit name or by passing
+ in a YTArray. If a list or unitless array is supplied, code units are
+ assumed.
+ width : tuple or a float.
+ Width can have four different formats to support windows with variable
+ x and y widths. They are:
+
+ ================================== =======================
+ format example
+ ================================== =======================
+ (float, string) (10,'kpc')
+ ((float, string), (float, string)) ((10,'kpc'),(15,'kpc'))
+ float 0.2
+ (float, float) (0.2, 0.3)
+ ================================== =======================
+
+ For example, (10, 'kpc') requests a plot window that is 10 kiloparsecs
+ wide in the x and y directions, ((10,'kpc'),(15,'kpc')) requests a
+ window that is 10 kiloparsecs wide along the x axis and 15
+ kiloparsecs wide along the y axis. In the other two examples, code
+ units are assumed, for example (0.2, 0.3) requests a plot that has an
+ x width of 0.2 and a y width of 0.3 in code units. If units are
+ provided the resulting plot axis labels will use the supplied units.
+ image_res : an int or 2-tuple of ints
+ Specify the resolution of the resulting image. If not provided, it will be
+ determined based on the minimum cell size of the dataset.
+ """
+ def __init__(self, ds, axis, fields, center="c", width=None, image_res=None, **kwargs):
+ fields = ensure_list(fields)
+ axis = fix_axis(axis, ds)
+ center, dcenter = ds.coordinates.sanitize_center(center, axis)
+ slc = ds.slice(axis, center[axis], **kwargs)
+ w, frb = construct_image(ds, axis, slc, dcenter, width=width, image_res=image_res)
+ super(FITSSlice, self).__init__(frb, fields=fields, wcs=w)
+
+
+class FITSProjection(FITSImageData):
+ r"""
+ Generate a FITSImageData of an on-axis projection.
+
+ Parameters
+ ----------
+ ds : :class:`yt.data_objects.api.Dataset`
+ The dataset object.
+ axis : character or integer
+ The axis along which to project. One of "x","y","z", or 0,1,2.
+ fields : string or list of strings
+ The fields to project
+ weight_field : string
+ The field used to weight the projection.
+ center : A sequence of floats, a string, or a tuple.
+ The coordinate of the center of the image. If set to 'c', 'center' or
+ left blank, the plot is centered on the middle of the domain. If set to
+ 'max' or 'm', the center will be located at the maximum of the
+ ('gas', 'density') field. Centering on the max or min of a specific
+ field is supported by providing a tuple such as ("min","temperature") or
+ ("max","dark_matter_density"). Units can be specified by passing in *center*
+ as a tuple containing a coordinate and string unit name or by passing
+ in a YTArray. If a list or unitless array is supplied, code units are
+ assumed.
+ width : tuple or a float.
+ Width can have four different formats to support windows with variable
+ x and y widths. They are:
+
+ ================================== =======================
+ format example
+ ================================== =======================
+ (float, string) (10,'kpc')
+ ((float, string), (float, string)) ((10,'kpc'),(15,'kpc'))
+ float 0.2
+ (float, float) (0.2, 0.3)
+ ================================== =======================
+
+ For example, (10, 'kpc') requests a plot window that is 10 kiloparsecs
+ wide in the x and y directions, ((10,'kpc'),(15,'kpc')) requests a
+ window that is 10 kiloparsecs wide along the x axis and 15
+ kiloparsecs wide along the y axis. In the other two examples, code
+ units are assumed, for example (0.2, 0.3) requests a plot that has an
+ x width of 0.2 and a y width of 0.3 in code units. If units are
+ provided the resulting plot axis labels will use the supplied units.
+ image_res : an int or 2-tuple of ints
+ Specify the resolution of the resulting image. If not provided, it will be
+ determined based on the minimum cell size of the dataset.
+ """
+ def __init__(self, ds, axis, fields, center="c", width=None,
+ weight_field=None, image_res=None, **kwargs):
+ fields = ensure_list(fields)
+ axis = fix_axis(axis, ds)
+ center, dcenter = ds.coordinates.sanitize_center(center, axis)
+ prj = ds.proj(fields[0], axis, weight_field=weight_field, **kwargs)
+ w, frb = construct_image(ds, axis, prj, dcenter, width=width, image_res=image_res)
+ super(FITSProjection, self).__init__(frb, fields=fields, wcs=w)
+
+class FITSOffAxisSlice(FITSImageData):
+ r"""
+ Generate a FITSImageData of an off-axis slice.
+
+ Parameters
+ ----------
+ ds : :class:`yt.data_objects.api.Dataset`
+ The dataset object.
+ normal : a sequence of floats
+ The vector normal to the projection plane.
+ fields : string or list of strings
+ The fields to slice
+ center : A sequence of floats, a string, or a tuple.
+ The coordinate of the center of the image. If set to 'c', 'center' or
+ left blank, the plot is centered on the middle of the domain. If set to
+ 'max' or 'm', the center will be located at the maximum of the
+ ('gas', 'density') field. Centering on the max or min of a specific
+ field is supported by providing a tuple such as ("min","temperature") or
+ ("max","dark_matter_density"). Units can be specified by passing in *center*
+ as a tuple containing a coordinate and string unit name or by passing
+ in a YTArray. If a list or unitless array is supplied, code units are
+ assumed.
+ width : tuple or a float.
+ Width can have four different formats to support windows with variable
+ x and y widths. They are:
+
+ ================================== =======================
+ format example
+ ================================== =======================
+ (float, string) (10,'kpc')
+ ((float, string), (float, string)) ((10,'kpc'),(15,'kpc'))
+ float 0.2
+ (float, float) (0.2, 0.3)
+ ================================== =======================
+
+ For example, (10, 'kpc') requests a plot window that is 10 kiloparsecs
+ wide in the x and y directions, ((10,'kpc'),(15,'kpc')) requests a
+ window that is 10 kiloparsecs wide along the x axis and 15
+ kiloparsecs wide along the y axis. In the other two examples, code
+ units are assumed, for example (0.2, 0.3) requests a plot that has an
+ x width of 0.2 and a y width of 0.3 in code units. If units are
+ provided the resulting plot axis labels will use the supplied units.
+ image_res : an int or 2-tuple of ints
+ Specify the resolution of the resulting image.
+ north_vector : a sequence of floats
+ A vector defining the 'up' direction in the plot. This
+ option sets the orientation of the slicing plane. If not
+ set, an arbitrary grid-aligned north-vector is chosen.
+ """
+ def __init__(self, ds, normal, fields, center='c', width=None, image_res=512,
+ north_vector=None):
+ fields = ensure_list(fields)
+ center, dcenter = ds.coordinates.sanitize_center(center, 4)
+ cut = ds.cutting(normal, center, north_vector=north_vector)
+ center = ds.arr([0.0] * 2, 'code_length')
+ w, frb = construct_image(ds, normal, cut, center, width=width, image_res=image_res)
+ super(FITSOffAxisSlice, self).__init__(frb, fields=fields, wcs=w)
+
+
+class FITSOffAxisProjection(FITSImageData):
+ r"""
+ Generate a FITSImageData of an off-axis projection.
+
+ Parameters
+ ----------
+ ds : :class:`yt.data_objects.api.Dataset`
+ This is the dataset object corresponding to the
+ simulation output to be plotted.
+ normal : a sequence of floats
+ The vector normal to the projection plane.
+ fields : string, list of strings
+ The name of the field(s) to be plotted.
+ center : A sequence of floats, a string, or a tuple.
+ The coordinate of the center of the image. If set to 'c', 'center' or
+ left blank, the plot is centered on the middle of the domain. If set to
+ 'max' or 'm', the center will be located at the maximum of the
+ ('gas', 'density') field. Centering on the max or min of a specific
+ field is supported by providing a tuple such as ("min","temperature") or
+ ("max","dark_matter_density"). Units can be specified by passing in *center*
+ as a tuple containing a coordinate and string unit name or by passing
+ in a YTArray. If a list or unitless array is supplied, code units are
+ assumed.
+ width : tuple or a float.
+ Width can have four different formats to support windows with variable
+ x and y widths. They are:
+
+ ================================== =======================
+ format example
+ ================================== =======================
+ (float, string) (10,'kpc')
+ ((float, string), (float, string)) ((10,'kpc'),(15,'kpc'))
+ float 0.2
+ (float, float) (0.2, 0.3)
+ ================================== =======================
+
+ For example, (10, 'kpc') requests a plot window that is 10 kiloparsecs
+ wide in the x and y directions, ((10,'kpc'),(15,'kpc')) requests a
+ window that is 10 kiloparsecs wide along the x axis and 15
+ kiloparsecs wide along the y axis. In the other two examples, code
+ units are assumed, for example (0.2, 0.3) requests a plot that has an
+ x width of 0.2 and a y width of 0.3 in code units. If units are
+ provided the resulting plot axis labels will use the supplied units.
+ depth : A tuple or a float
+ A tuple containing the depth to project through and the string
+ key of the unit: (width, 'unit'). If set to a float, code units
+ are assumed
+ weight_field : string
+ The name of the weighting field. Set to None for no weight.
+ image_res : an int or 2-tuple of ints
+ Specify the resolution of the resulting image.
+ depth_res : integer
+ Deprecated, this is still in the function signature for API
+ compatibility
+ north_vector : a sequence of floats
+ A vector defining the 'up' direction in the plot. This
+ option sets the orientation of the slicing plane. If not
+ set, an arbitrary grid-aligned north-vector is chosen.
+ method : string
+ The method of projection. Valid methods are:
+
+ "integrate" with no weight_field specified : integrate the requested
+ field along the line of sight.
+
+ "integrate" with a weight_field specified : weight the requested
+ field by the weighting field and integrate along the line of sight.
+
+ "sum" : This method is the same as integrate, except that it does not
+ multiply by a path length when performing the integration, and is
+ just a straight summation of the field along the given axis. WARNING:
+ This should only be used for uniform resolution grid datasets, as other
+ datasets may result in unphysical images.
+ data_source : yt.data_objects.data_containers.YTSelectionContainer, optional
+ If specified, this will be the data source used for selecting regions to project.
+ """
+ def __init__(self, ds, normal, fields, center='c', width=(1.0, 'unitary'),
+ weight_field=None, image_res=512, depth_res=256, data_source=None,
+ north_vector=None, depth=(1.0,"unitary"), no_ghost=False, method='integrate'):
+ fields = ensure_list(fields)
+ center, dcenter = ds.coordinates.sanitize_center(center, 4)
+ buf = {}
+ width = ds.coordinates.sanitize_width(normal, width, depth)
+ wd = tuple(el.in_units('code_length').v for el in width)
+ if not iterable(image_res):
+ image_res = (image_res, image_res)
+ res = (image_res[0], image_res[1])
+ if data_source is None:
+ source = ds
+ else:
+ source = data_source
+ for field in fields:
+ buf[field] = off_axis_projection(source, center, normal, wd, res, field,
+ no_ghost=no_ghost, north_vector=north_vector,
+ method=method, weight=weight_field).swapaxes(0, 1)
+ center = ds.arr([0.0] * 2, 'code_length')
+ w, not_an_frb = construct_image(ds, normal, buf, center, width=width, image_res=image_res)
+ super(FITSOffAxisProjection, self).__init__(buf, fields=fields, wcs=w)
diff -r 83d2f4c6b6608ea701acf74d326c44125ae93b5e -r 959db1b322f8836bab7a600487819d8cef23dd13 yt/visualization/fixed_resolution.py
--- a/yt/visualization/fixed_resolution.py
+++ b/yt/visualization/fixed_resolution.py
@@ -362,7 +362,7 @@
the length units that the coordinates are written in, default 'cm'.
"""
- from yt.utilities.fits_image import FITSImageData
+ from yt.visualization.fits_image import FITSImageData
if fields is None:
fields = list(self.data.keys())
This diff is so big that we needed to truncate the remainder.
Repository URL: https://bitbucket.org/yt_analysis/yt/
--
This is a commit notification from bitbucket.org. You are receiving
this because you have the service enabled, addressing the recipient of
this email.
More information about the yt-svn
mailing list