[AstroPy] SpectralCube and ytcube.quick_isocontour with FITS
salome
philippe.salome at obspm.fr
Tue Jun 4 09:38:43 EDT 2019
Hello,
I have troubles with trying to export a FITS cube (118, 480, 480) into Sketchlab with SpectralCube.
I am afraid I don’t understand what’s the array dimension related to. If anyone has an idea, it would
be very helpful. Thanks a lot !
Here is the example
===================================================================
import astropy.units as u
import numpy as np
from spectral_cube import SpectralCube
from yt.mods import ColorTransferFunction, write_bitmap
import astropy.units as u
# Read in spectral cube
filename = '/Users/salome/uid___A001_X88f_X25b.HH212_sci.spw25.cube.I.pbcor.fits'
cube = SpectralCube.read(filename, format='fits')
cube.min()
cube.max()
# Extract the yt object from the SpectralCube instance
ytcube = cube.to_yt(spectral_factor=0.75)
WARNING: StokesWarning: Cube is a Stokes cube, returning spectral cube for I component [spectral_cube.spectral_cube]
yt : [WARNING ] 2019-06-04 15:31:36,176 Cannot find time
yt : [INFO ] 2019-06-04 15:31:36,177 Detected these axes: RA---SIN DEC--SIN FREQ
yt : [WARNING ] 2019-06-04 15:31:36,181 No length conversion provided. Assuming 1 = 1 cm.
yt : [INFO ] 2019-06-04 15:31:36,197 Parameters: current_time = 0.0
yt : [INFO ] 2019-06-04 15:31:36,197 Parameters: domain_dimensions = [480 480 118]
yt : [INFO ] 2019-06-04 15:31:36,198 Parameters: domain_left_edge = [0.5 0.5 0.5]
yt : [INFO ] 2019-06-04 15:31:36,199 Parameters: domain_right_edge = [480.5 480.5 89. ]
yt : [INFO ] 2019-06-04 15:31:36,202 Parameters: cosmological_simulation = 0.0
WARNING: PossiblySlowWarning: This function (<function BaseSpectralCube.min at 0x1207c40d0>) requires loading the entire cube into memory and may therefore be slow. [spectral_cube.utils]
WARNING: PossiblySlowWarning: This function (<function BaseSpectralCube.max at 0x1207afea0>) requires loading the entire cube into memory and may therefore be slow. [spectral_cube.utils]
ytcube.quick_isocontour(export_to='ply', filename='meshes.ply', level=0.02)
===================================================================
—>
WARNING: PossiblySlowWarning: This function (<function BaseSpectralCube.std at 0x1207afb70>) requires loading the entire cube into memory and may therefore be slow. [spectral_cube.utils]
yt : [INFO ] 2019-06-04 15:20:21,909 Adding field flux to the list of fields.
---------------------------------------------------------------------------
IndexError Traceback (most recent call last)
<ipython-input-23-26d45438b450> in <module>
----> 1 ytcube.quick_isocontour()
/anaconda3/lib/python3.6/site-packages/spectral_cube/ytcube.py in quick_isocontour(self, level, title, description, color_map, color_log, export_to, filename, **kwargs)
229 description=description,
230 color_map=color_map,
--> 231 color_log=color_log, **kwargs)
232 elif export_to == 'obj':
233 if filename is None:
/anaconda3/lib/python3.6/site-packages/yt/data_objects/construction_data_containers.py in export_sketchfab(self, title, description, api_key, color_field, color_map, color_log, bounds, no_ghost)
1924 ply_file = TemporaryFile()
1925 self.export_ply(ply_file, bounds, color_field, color_map, color_log,
-> 1926 sample_type = "vertex", no_ghost = no_ghost)
1927 ply_file.seek(0)
1928 # Greater than ten million vertices and we throw an error but dump
/anaconda3/lib/python3.6/site-packages/yt/data_objects/construction_data_containers.py in export_ply(self, filename, bounds, color_field, color_map, color_log, sample_type, no_ghost)
1764 if color_map is None:
1765 color_map = ytcfg.get("yt", "default_colormap")
-> 1766 if self.vertices is None:
1767 self.get_data(color_field, sample_type, no_ghost=no_ghost)
1768 elif color_field is not None:
/anaconda3/lib/python3.6/site-packages/yt/data_objects/construction_data_containers.py in vertices(self)
1300 def vertices(self):
1301 if self._vertices is None:
-> 1302 self.get_data()
1303 return self._vertices
1304
/anaconda3/lib/python3.6/site-packages/yt/data_objects/construction_data_containers.py in get_data(self, fields, sample_type, no_ghost)
1170 my_verts = self._extract_isocontours_from_grid(
1171 block, self.surface_field, self.field_value,
-> 1172 mask, fields, sample_type, no_ghost=no_ghost)
1173 if fields is not None:
1174 my_verts, svals = my_verts
/anaconda3/lib/python3.6/site-packages/yt/data_objects/construction_data_containers.py in _extract_isocontours_from_grid(self, grid, field, value, mask, sample_values, sample_type, no_ghost)
1194 no_ghost = False):
1195 # TODO: check if multiple fields can be passed here
-> 1196 vals = grid.get_vertex_centered_data([field], no_ghost=no_ghost)[field]
1197 if sample_values is not None:
1198 # TODO: is no_ghost=False correct here?
/anaconda3/lib/python3.6/site-packages/yt/data_objects/grid_patch.py in get_vertex_centered_data(self, fields, smoothed, no_ghost)
308 new_fields[field], output_left)
309 else:
--> 310 cg = self.retrieve_ghost_zones(1, fields, smoothed=smoothed)
311 for field in fields:
312 np.add(new_fields[field], cg[field][1: ,1: ,1: ], new_fields[field])
/anaconda3/lib/python3.6/site-packages/yt/data_objects/grid_patch.py in retrieve_ghost_zones(self, n_zones, fields, all_levels, smoothed)
270 level, new_left_edge,
271 field_parameters = field_parameters,
--> 272 **kwargs)
273 else:
274 cube = self.ds.covering_grid(level, new_left_edge,
/anaconda3/lib/python3.6/site-packages/yt/data_objects/construction_data_containers.py in __init__(self, *args, **kwargs)
925 ds.domain_dimensions.astype("float64"))
926 self.global_endindex = None
--> 927 YTCoveringGrid.__init__(self, *args, **kwargs)
928 self._final_start_index = self.global_startindex
929
/anaconda3/lib/python3.6/site-packages/yt/data_objects/construction_data_containers.py in __init__(self, level, left_edge, dims, fields, ds, num_ghost_zones, use_pbar, field_parameters)
553 (self.left_edge-self.ds.domain_left_edge)/self.dds).astype('int64')
554 self._setup_data_source()
--> 555 self.get_data(fields)
556
557 @property
/anaconda3/lib/python3.6/site-packages/yt/data_objects/construction_data_containers.py in get_data(self, fields)
640 raise
641 if len(part) > 0: self._fill_particles(part)
--> 642 if len(fill) > 0: self._fill_fields(fill)
643 for a, f in sorted(alias.items()):
644 if f.particle_type:
/anaconda3/lib/python3.6/site-packages/yt/data_objects/construction_data_containers.py in _fill_fields(self, fields)
1001 domain_dims = domain_dims.astype("int64")
1002 tot = ls.current_dims.prod()
-> 1003 for chunk in ls.data_source.chunks(fields, "io"):
1004 chunk[fields[0]]
1005 input_fields = [chunk[field] for field in fields]
/anaconda3/lib/python3.6/site-packages/yt/data_objects/data_containers.py in chunks(self, fields, chunking_style, **kwargs)
1274 continue
1275 with self._chunked_read(chunk):
-> 1276 self.get_data(fields)
1277 # NOTE: we yield before releasing the context
1278 yield self
/anaconda3/lib/python3.6/site-packages/yt/data_objects/data_containers.py in get_data(self, fields)
1368 # need to be generated.
1369 read_fluids, gen_fluids = self.index._read_fluid_fields(
-> 1370 fluids, self, self._current_chunk)
1371 for f, v in read_fluids.items():
1372 self.field_data[f] = self.ds.arr(v, input_units = finfos[f].units)
/anaconda3/lib/python3.6/site-packages/yt/geometry/geometry_handler.py in _read_fluid_fields(self, fields, dobj, chunk)
243 selector,
244 fields_to_read,
--> 245 chunk_size)
246 return fields_to_return, fields_to_generate
247
/anaconda3/lib/python3.6/site-packages/yt/frontends/fits/io.py in _read_fluid_selection(self, chunks, selector, fields, size)
98 data[np.isnan(data)] = self.ds.nan_mask["all"]
99 data = bzero + bscale*data
--> 100 ind += g.select(selector, data.astype("float64"), rv[field], ind)
101 return rv
/anaconda3/lib/python3.6/site-packages/yt/data_objects/grid_patch.py in select(self, selector, source, dest, offset)
417 slices = get_nodal_slices(source.shape, nodal_flag, dim)
418 for i , sl in enumerate(slices):
--> 419 dest[offset:offset+count, i] = source[sl][np.squeeze(mask)]
420 return count
421
IndexError: boolean index did not match indexed array along dimension 2; dimension is 5 but corresponding boolean dimension is 21
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mail.python.org/pipermail/astropy/attachments/20190604/4e35680b/attachment-0001.html>
More information about the AstroPy
mailing list