Adding an interactive way to slice HS over a wavelength range, and view the result spatially!

[0Hughman0]

Describe the functionality you would like to see.

In the context of taking hyperspectral maps, sometimes you may have a sample where in different regions one would expect peaks at different wavelengths to become dominant.

It's useful to be able to quickly scan over a given wavelength (or energy!) range, and see how the intensity of signals within that range vary spatially.

Describe the context

If implemented correctly, this could be used for a number of different hyperspectral signals.

The Dorset software for Chronos CL systems has this feature, but as far as I know this doesn't exist in Lumispy/ Hyperspy.

Proposed Solution

Perhaps this specific sort of analysis is outside of the scope of Lumispy. However, I thought something like this might fit well into the utils section of the library.

I have put together a provisional implementation here:

https://github.com/0Hughman0/lumispy/blob/0fc6fa18e82069af116acdfdc066fb244da05387/lumispy/utils/analysis.py#L4

I've called the feature 'AutoPeakMap'... not sure if this is the best name!

An example would be:

from lumispy.utils.analysis import AutoPeakMap, RangeChannel, plot_auto_peakmap

cl = hs.load(path, signal_type='CL_SEM')

apm = AutoPeakMap(cl,
                  RangeChannel('red', range_=[730, 760]),
                  RangeChannel('green', range_=[300, 400]))
apm.plot()

Renders two plots, the 'navigator', in the hyperspy sense:

and the 'signal':

The ranges in the 'navigator' can be moved around interactively and the signal plot will update accordingly.

For more quick and dirty analysis I implemented a convenience method:

plot_auto_peakmap(cl)

Which chooses sensible (enough) defaults for the RangeChannels.

If this is something you think could slot into LumiSpy, I'll go ahead and write the tests and expand the docstrings etc.

Feedback on the implementation/ naming more than welcome!

Thanks,

Hugh

[ericpre]

This makes sense to have such features; as an alternative of implementing new classes, I would suggest to use already existing functionalities using the following approach:

# Step 1: find peak, use hyperspy find_peaks_ohaver?
peaks =...
width = ...

# Step 2: add roi to figure
# Initialise the roi with the value from the step 1 and the corresponding/given width
roi1_signal = hs.roi.SpanROI(peaks[0]-width/2, peaks[0]+width/2).interactive(s, color="red")
roi2_signal = hs.roi.SpanROI(peaks[1]-width/2, peaks[1]+width/2).interactive(s, color="red")

# Step 3: generate map from signal
map1 = hs.interactive(roi1_signal.sum)
map2 = hs.interactive(roi2_signal.sum)

# Step 4: generate overlay maps figure
hs.plot.plot_images([map1, map2], overlay=True)

In the first two steps need the correct argument to be efficient (avoid redundant computation) and are very similar to: https://github.com/pyxem/pyxem/blob/996651e5c2ac22cbaf0e11fca1eed8cbd5a35616/pyxem/signals/common_diffraction.py#L76-L127

The last step is the overlay of the two (or more) maps and I am not sure if it would be interactive, if not, this should be fixed in hyperspy.

This could be nicely wrapped up in a signal method, which could possibly go in hyperspy as this is a generic functionality and other spectroscopy technique could benefit from it.

[0Hughman0]

Great thanks so much for the feedback.

Indeed the functionality I want to add looks very similar to the one you linked in pyxem.

I do think the ability to have multiple 'channels' would be quite beneficial, and I think interactivity is a must for this to be useful for people to quickly use this tool.

I do think such a thing could be widely useful, but I don't know what the philosophy is hyperspy is in terms of stacking more functionality into the signal classes.

I'm super impressed it looks like this can be accomplished with hs's interactive features. Not sure why, but I sometimes find them a bit bewildering and frequently run into errors, which puts me off using it a bit.

Sadly I have tried to get the snippet you have suggested working, I'm not sure where I'm going wrong:

import hyperspy.api as hs

# generate some appropriately shaped data
s = hs.signals.Signal1D(data=np.random.random((64, 64, 1024)), 
                        axes=[{'name': 'x', 'size': 64}, 
                              {'name': 'y', 'size': 64}, 
                              {'name': 'sig', 'size': 1024}])


peaks = [350, 750]
widths = [100, 50]

# I guess maybe issue is I'm not plotting the appropriate thing here?
s.plot()

# Step 2: add roi to figure
# Initialise the roi with the value from the step 1 and the corresponding/given width
roi1_signal = hs.roi.SpanROI(peaks[0]-widths[0]/2, peaks[0]+widths[0]/2).interactive(s, color="red")
roi2_signal = hs.roi.SpanROI(peaks[1]-widths[1]/2, peaks[1]+widths[1]/2).interactive(s, color="green")

# Step 3: generate map from signal
map1 = hs.interactive(roi1_signal.sum)
map2 = hs.interactive(roi2_signal.sum)

# Step 4: generate overlay maps figure
hs.plot.plot_images([map1, map2], overlay=True)

From that I get:

---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
~\AppData\Local\Temp/ipykernel_20276/994102482.py in <module>
     12 # Step 2: add roi to figure
     13 # Initialise the roi with the value from the step 1 and the corresponding/given width
---> 14 roi1_signal = hs.roi.SpanROI(peaks[0]-widths[0]/2, peaks[0]+widths[0]/2).interactive(s, color="red")
     15 roi2_signal = hs.roi.SpanROI(peaks[1]-widths[1]/2, peaks[1]+widths[1]/2).interactive(s, color="green")
     16 

C:\ProgramData\Anaconda3\envs\PhD\lib\site-packages\hyperspy\roi.py in interactive(self, signal, navigation_signal, out, color, snap, **kwargs)
    429                 [])
    430         if out is None:
--> 431             return interactive(self.__call__,
    432                                event=self.events.changed,
    433                                signal=signal,

C:\ProgramData\Anaconda3\envs\PhD\lib\site-packages\hyperspy\interactive.py in interactive(f, event, recompute_out_event, *args, **kwargs)
    128 
    129 def interactive(f, event="auto", recompute_out_event="auto", *args, **kwargs):
--> 130     cls = Interactive(f, event, recompute_out_event, *args, **kwargs)
    131     return cls.out
    132 

C:\ProgramData\Anaconda3\envs\PhD\lib\site-packages\hyperspy\interactive.py in __init__(self, f, event, recompute_out_event, *args, **kwargs)
     79             self.out = self.kwargs.pop('out')
     80         else:
---> 81             self.out = self.f(*self.args, **self.kwargs)
     82         # Reuse the `_plot_kwargs` for the roi if available
     83         if _plot_kwargs and 'signal' in self.kwargs:

C:\ProgramData\Anaconda3\envs\PhD\lib\site-packages\hyperspy\roi.py in __call__(self, signal, out, axes)
    213 
    214         natax = signal.axes_manager._get_axes_in_natural_order()
--> 215         slices = self._make_slices(natax, axes)
    216         nav_axes = [ax.navigate for ax in axes]
    217         nav_dim = signal.axes_manager.navigation_dimension

C:\ProgramData\Anaconda3\envs\PhD\lib\site-packages\hyperspy\roi.py in _make_slices(self, axes_collection, axes, ranges)
    170                 i = axes.index(ax)
    171                 try:
--> 172                     ilow = ax.value2index(ranges[i][0])
    173                 except ValueError:
    174                     if ranges[i][0] < ax.low_value:

C:\ProgramData\Anaconda3\envs\PhD\lib\site-packages\hyperspy\axes.py in value2index(self, value, rounding)
   1231                 return index
   1232             else:
-> 1233                 raise ValueError("The value is out of the axis limits")
   1234 
   1235     def update_axis(self):

ValueError: The value is out of the axis limits

Using the latest '1.7.1' version of hyperspy.

What's worse, I just updated from a dev version of 1.7.0dev0, and now my own snippet is broken 😅 - serves me right I guess!

[0Hughman0]

Ok, I've kept going with the branch here:

https://github.com/0Hughman0/lumispy/blob/auto_peak_map/lumispy/utils/analysis.py

I am trying to make a bit more use of the roi functionality. I did try chaining together events using hs.interactive, but I couldn't get the plots to update, even if the callbacks called sig._plot.navigator_plot.update(), not sure what I was doing wrong.

[ericpre]

Following the approach used in https://github.com/pyxem/pyxem/blob/996651e5c2ac22cbaf0e11fca1eed8cbd5a35616/pyxem/signals/common_diffraction.py#L76-L127, below is an example that works ignoring the part with hs.plot.plot_images that needs to be made interactive:

import hyperspy.api as hs
import numpy as np

# generate some appropriately shaped data
s = hs.signals.Signal1D(data=np.random.random((64, 64, 1024)), 
                        axes=[{'name': 'x', 'size': 64}, 
                              {'name': 'y', 'size': 64}, 
                              {'name': 'sig', 'size': 1024}])


peaks = [350, 750]
widths = [100, 50]

s.plot()

# Initialise the roi with the value from the step 1 and the corresponding/given width
roi1 = hs.roi.SpanROI(peaks[0]-widths[0]/2, peaks[0]+widths[0]/2)
roi2 = hs.roi.SpanROI(peaks[1]-widths[1]/2, peaks[1]+widths[1]/2)

roi1_signal = roi1.interactive(s, color="red", axes=s.axes_manager.signal_axes)
roi2_signal = roi2.interactive(s, color="green", axes=s.axes_manager.signal_axes)

# placehodler for the sum
roi1_sum = roi1_signal.sum(axis=roi1_signal.axes_manager.signal_axes).T
roi2_sum = roi2_signal.sum(axis=roi2_signal.axes_manager.signal_axes).T

hs.interactive(
    roi1_signal.nansum,
    axis=roi1_signal.axes_manager.signal_axes,
    event=roi1.events.changed,
    recompute_out_event=None,
    out=roi1_sum,
)
hs.interactive(
    roi2_signal.nansum,
    axis=roi2_signal.axes_manager.signal_axes,
    event=roi2.events.changed,
    recompute_out_event=None,
    out=roi2_sum,
)


# # Plot the result
roi1_sum.plot()
roi2_sum.plot()

# This part is the one that needs improvement to be interactive
# hs.plot.plot_images([roi1_sum, roi1_sum])

[ericpre]

Regarding implementing it in hs.plot.plot_images, this could be done in a similar way as its Signal1D counterpart: hyperspy/hyperspy#2414

[jordiferrero]

@0Hughman0 thanks for the contribution! Let us know if you managed to make it work using hyperspy interactive method with 2d images. I have also tried in the past and got very confused with the interactivity of 2 rois in parallel.
If it does not work, I believe lumispy could benefit from such functionality.

[jlaehne]

Indeed, where possible to use existing HyperSpy functions that makes sense. As brought up by @jordiferrero in #55 we can consider to still have a set of dedicated plot functions in LumiSpy amending those of HyperSpy - just because it is easier for some types of plots to have a dedicated function and not to have to work with a complicated set of options to a more generic function.

[0Hughman0]

Thanks for the encouragement.

I will have a go at achieving the same thing with built-in hyperspy functionality. Although it sounds like from what you are saying, @ericpre, it will not be possible to overlay the maps and have them automatically update until there's a some tweaks to the hs.plot.plot_images function?

From briefly looking at the source, I think either hs.plot.plot_images needs to provide some sort of events argument, or if it were wrapped in some sort of object (akin to MPL_HyperExplorer), maybe users could hook into the internals to setup their own callbacks.

I don't know how people feel about having a navigation plot (where you select your wavelength ranges), but multiple maps of the integrated counts for each 'channel'? It's not quite a 1 to 1 replication of the Dorset functionality, but probably still useful.

[0Hughman0]

ok please see #148, I haven't gone as far as writing tests and docstrings (excited to try pytest-mpl!)

[ericpre]

I will have a go at achieving the same thing with built-in hyperspy functionality. Although it sounds like from what you are saying, @ericpre, it will not be possible to overlay the maps and have them automatically update until there's a some tweaks to the hs.plot.plot_images function?

From briefly looking at the source, I think either hs.plot.plot_images needs to provide some sort of events argument, or if it were wrapped in some sort of object (akin to MPL_HyperExplorer), maybe users could hook into the internals to setup their own callbacks.

As it has been done for hs.plot.plot_spectra (see hyperspy/hyperspy#2414), it should be implemented in a way that the plot update automatically when the data changed. This changed is needed in the hs.plot.plot_images and the user doesn't need to do anything with events.

[jlaehne]

Upstreamed to hyperspy/hyperspy#3224