spacetelescope / jdaviz

JWST astronomical data analysis tools in the Jupyter platform
https://jdaviz.readthedocs.io/
BSD 3-Clause "New" or "Revised" License
142 stars 74 forks source link

Model fitting in Specviz2d: AttributeError: 'numpy.ndarray' object has no attribute 'to_value' #2108

Open pllim opened 1 year ago

pllim commented 1 year ago

Workflow to reproduce:

  1. Fire up Specviz2D example notebook.
  2. Select two subsets, preferably one with continuum and one with feature.
  3. Open model fitting plugin.
  4. Add a Const1D (C) model over the continuum.
  5. Add a Gaussian1D (G) model over the feature.
  6. Set equation as "G-C".
  7. Enable "calculate residuals".
  8. Click "Fit Model" and see traceback below.
File ipyvue/VueTemplateWidget.py:60, in Events._handle_event(self, _, content, buffers)
     58     getattr(self, "vue_" + event)(data, buffers)
     59 else:
---> 60     getattr(self, "vue_" + event)(data)

File jdaviz/configs/default/plugins/model_fitting/model_fitting.py:762, in ModelFitting.vue_apply(self, event)
    761 def vue_apply(self, event):
--> 762     self.calculate_fit()

File jdaviz/configs/default/plugins/model_fitting/model_fitting.py:729, in ModelFitting.calculate_fit(self, add_data)
    727     ret = self._fit_model_to_cube(add_data=add_data)
    728 else:
--> 729     ret = self._fit_model_to_spectrum(add_data=add_data)
    731 if ret is None:  # pragma: no cover
    732     # something went wrong in the fitting call and (hopefully) already raised a warning,
    733     # but we don't have anything to add to the table
    734     return ret

File jdaviz/configs/default/plugins/model_fitting/model_fitting.py:795, in ModelFitting._fit_model_to_spectrum(self, add_data)
    793 if add_data:
    794     self.app.fitted_models[self.results_label] = fitted_model
--> 795     self.add_results.add_results_from_plugin(fitted_spectrum)
    797     if self.residuals_calculate:
    798         # NOTE: this will NOT load into the viewer since we have already called
    799         # add_results_from_plugin above.
    800         self.add_results.add_results_from_plugin(masked_spectrum-fitted_spectrum,
    801                                                  label=self.residuals.value,
    802                                                  replace=False)

File jdaviz/core/template_mixin.py:1857, in AddResults.add_results_from_plugin(self, data_item, replace, label)
   1852 self.app.add_data(data_item, label)
   1854 if add_to_viewer_selected != 'None':
   1855     # replace the contents in the selected viewer with the results from this plugin
   1856     # TODO: switch to an instance/classname check?
-> 1857     self.app.add_data_to_viewer(self.viewer.selected_id,
   1858                                 label,
   1859                                 visible=visible, clear_other_data=replace)
   1861 # update overwrite warnings, etc
   1862 self._on_label_changed()

File jdaviz/app.py:1068, in Application.add_data_to_viewer(self, viewer_reference, data_label, visible, clear_other_data)
   1065 if viewer_item is None:
   1066     raise ValueError(f"Could not identify viewer with reference {viewer_reference}")
-> 1068 self.set_data_visibility(viewer_item, data_label, visible=visible, replace=clear_other_data)

File jdaviz/app.py:1458, in Application.set_data_visibility(self, viewer_reference, data_label, visible, replace)
   1455 active_data = self.data_collection[viewer_data_labels[0]]
   1456 if (hasattr(active_data, "_preferred_translation")
   1457         and active_data._preferred_translation is not None):
-> 1458     self._set_plot_axes_labels(viewer_id)
   1460 if self.config == 'imviz':
   1461     viewer.on_limits_change()

File jdaviz/app.py:1081, in Application._set_plot_axes_labels(self, viewer_id)
   1071 """
   1072 Sets the plot axes labels to be the units of the data to be loaded.
   1073 
   (...)
   1077     The UUID associated with the desired viewer item.
   1078 """
   1079 viewer = self._viewer_by_id(viewer_id)
-> 1081 viewer.set_plot_axes()

File jdaviz/configs/specviz/plugins/viewers.py:480, in SpecvizProfileView.set_plot_axes(self)
    478 def set_plot_axes(self):
    479     # Get data to be used for axes labels
--> 480     data = self.data()[0]
    482     # Set axes labels for the spectrum viewer
    483     spectral_axis_unit_type = str(data.spectral_axis.unit.physical_type).title()

File jdaviz/configs/specviz/plugins/viewers.py:94, in SpecvizProfileView.data(self, cls)
     92 handler, _ = data_translator.get_handler_for(_class)
     93 try:
---> 94     layer_data = handler.to_object(layer_data, statistic=statistic)
     95 except IncompatibleAttribute:
     96     continue

File glue_astronomy/translators/spectrum1d.py:284, in Specutils1DHandler.to_object(self, data_or_subset, attribute, statistic)
    279         data_kwargs.update({attribute_label: values,
    280                            'mask': mask})
    282     return data_kwargs
--> 284 data_kwargs = parse_attributes(
    285     [attribute] if not hasattr(attribute, '__len__') else attribute)
    287 return Spectrum1D(**data_kwargs, **kwargs)

File glue_astronomy/translators/spectrum1d.py:251, in Specutils1DHandler.to_object.<locals>.parse_attributes(attributes)
    249     mask = None
    250 else:
--> 251     mask = data.get_mask(subset_state=subset_state)
    252     mask = ~mask
    254 # Collapse values and mask to profile

File glue/core/data.py:1437, in Data.get_mask(self, subset_state, view)
   1435 def get_mask(self, subset_state, view=None):
   1436     try:
-> 1437         return subset_state.to_mask(self, view=view)
   1438     except IncompatibleAttribute:
   1439         return get_mask_with_key_joins(self, self._key_joins, subset_state, view=view)

File glue/core/subset.py:793, in RangeSubsetState.to_mask(self, data, view)
    791 @contract(data='isinstance(Data)', view='array_view')
    792 def to_mask(self, data, view=None):
--> 793     x = data[self.att, view]
    794     result = (x >= self.lo) & (x <= self.hi)
    795     return result

File glue/core/data.py:592, in BaseCartesianData.__getitem__(self, key)
    589     if key is None:
    590         raise IncompatibleAttribute(_k)
--> 592 return self.get_data(key, view=view)

File glue/core/data.py:1418, in Data.get_data(self, cid, view)
   1416     result = comp[view]
   1417 else:
-> 1418     result = comp.data
   1420 return result

File glue/core/component.py:199, in DerivedComponent.data(self)
    196 @property
    197 def data(self):
    198     """Return the numerical data as a numpy array"""
--> 199     return self._link.compute(self._data)

File glue/core/component_link.py:166, in ComponentLink.compute(self, data, view)
    145 """
    146 For a given data set, compute the component comp_to given the data
    147 associated with each comp_from and the ``using`` function
   (...)
    162     The data associated with comp_to component
    163 """
    165 # First we get the values of all the 'from' components.
--> 166 args = [data[join_component_view(f, view)] for f in self._from]
    168 # We keep track of the original shape of the arguments
    169 original_shape = args[0].shape

File glue/core/component_link.py:166, in <listcomp>(.0)
    145 """
    146 For a given data set, compute the component comp_to given the data
    147 associated with each comp_from and the ``using`` function
   (...)
    162     The data associated with comp_to component
    163 """
    165 # First we get the values of all the 'from' components.
--> 166 args = [data[join_component_view(f, view)] for f in self._from]
    168 # We keep track of the original shape of the arguments
    169 original_shape = args[0].shape

File glue/core/data.py:592, in BaseCartesianData.__getitem__(self, key)
    589     if key is None:
    590         raise IncompatibleAttribute(_k)
--> 592 return self.get_data(key, view=view)

File glue/core/data.py:1418, in Data.get_data(self, cid, view)
   1416     result = comp[view]
   1417 else:
-> 1418     result = comp.data
   1420 return result

File glue/core/component.py:199, in DerivedComponent.data(self)
    196 @property
    197 def data(self):
    198     """Return the numerical data as a numpy array"""
--> 199     return self._link.compute(self._data)

File glue/core/component_link.py:187, in ComponentLink.compute(self, data, view)
    184 args = np.broadcast_arrays(*args)
    186 # We call the actual linking function
--> 187 result = self._using(*args)
    189 # We call asarray since link functions may return Python scalars in some cases
    190 result = np.asarray(result)

File glue/core/component_link.py:395, in CoordinateComponentLink.using(self, *args)
    393     return pixel2world_single_axis(self.coords, *args2[::-1], world_axis=self.ndim - 1 - self.index)
    394 else:
--> 395     return world2pixel_single_axis(self.coords, *args2[::-1], pixel_axis=self.ndim - 1 - self.index)

File glue/core/coordinate_helpers.py:122, in world2pixel_single_axis(wcs, pixel_axis, *world)
    120     result = result.reshape(world_shape)
    121 else:
--> 122     result = wcs.world_to_pixel_values(*world)
    123     if len(world) > 1:
    124         result = result[pixel_axis]

File glue_astronomy/translators/spectrum1d.py:80, in PaddedSpectrumWCS.world_to_pixel_values(self, *world_arrays)
     76 def world_to_pixel_values(self, *world_arrays):
     77     # The ravel and reshape are needed because of
     78     # https://github.com/astropy/astropy/issues/12154
     79     wx = np.array(world_arrays[0])
---> 80     pixel_arrays = [self.spectral_wcs.world_to_pixel_values(wx.ravel()).reshape(wx.shape),
     81                     *world_arrays[1:]]
     82     return tuple(pixel_arrays)

File gwcs/api.py:142, in GWCSAPIMixin.world_to_pixel_values(self, *world_arrays)
    138 world_arrays = self._add_units_input(world_arrays, self.backward_transform, self.output_frame)
    140 result = self.invert(*world_arrays, with_units=False)
--> 142 return self._remove_quantity_output(result, self.input_frame)

File gwcs/api.py:81, in GWCSAPIMixin._remove_quantity_output(self, result, frame)
     78     if self.output_frame.naxes == 1:
     79         result = [result]
---> 81     result = tuple(r.to_value(unit) for r, unit in zip(result, frame.unit))
     83 # If we only have one output axes, we shouldn't return a tuple.
     84 if self.output_frame.naxes == 1 and isinstance(result, tuple):

File gwcs/api.py:81, in <genexpr>(.0)
     78     if self.output_frame.naxes == 1:
     79         result = [result]
---> 81     result = tuple(r.to_value(unit) for r, unit in zip(result, frame.unit))
     83 # If we only have one output axes, we shouldn't return a tuple.
     84 if self.output_frame.naxes == 1 and isinstance(result, tuple):

AttributeError: 'numpy.ndarray' object has no attribute 'to_value'

🐱

rosteen commented 1 year ago

I confirmed that I see the same error trace when trying to fit a linear model in Specviz2D.

Also, I wanted to point out that the workflow you wrote out isn't actually something we support, as far as I know. You can't select one subset to fit one model component to and set another component to a different subset as part of the same model fitting call. Every component gets fit to the subset that's selected when you press the Fit Model button.

pllim commented 1 year ago

Oh... Hmm... But isn't it natural to want to fit the continuum and the signal separately from different areas, and then calculate the signal parameters with continuum subtracted? Am I thinking about this wrong?

rosteen commented 1 year ago

I agree that's a common thing to want to do, but you currently have to do it as two separate model fits. This is facilitated by the Calculate Residuals toggle, which will add the model-subtracted data as a data entry. This is tangential to the bug here of course.

rosteen commented 1 year ago

I just checked and (somewhat surprisingly) model fitting does work in Mosviz.

kecnry commented 1 year ago

I just checked and (somewhat surprisingly) model fitting does work in Mosviz.

The spectrum in specviz2d is from the spectral extraction plugin (unless a 1d spectrum is loaded manually). This is not the case (yet) in mosviz.