Skip to content

Commit

Permalink
Merge branch 'fix_plot_component' into epsic_workshop
Browse files Browse the repository at this point in the history
  • Loading branch information
ericpre committed Mar 22, 2018
2 parents d9bbe17 + f2415b3 commit d8d9039
Show file tree
Hide file tree
Showing 21 changed files with 190 additions and 56 deletions.
8 changes: 6 additions & 2 deletions doc/user_guide/eels.rst
Original file line number Diff line number Diff line change
Expand Up @@ -149,14 +149,18 @@ Define the chemical composition of the sample
>>> s.add_elements(('B', 'N'))
In order to include the effect of plural scattering we provide the low-loss
spectrum to :py:meth:`~._signals.eels.EELSSpectrum_mixin.create_model`:
In order to include the effect of plural scattering, the component of the model is convolved with the loss loss spectrum in which case the low loss spectrum needs to be provided to :py:meth:`~._signals.eels.EELSSpectrum_mixin.create_model`:

.. code-block:: python
>>> m = s.create_model(ll=ll)
.. NOTE::

In case of the convolution with the low loss spectrum, the fitted value of some of the parameter may seem to be unexpectedly small which can be misleading. Those paramater values are actually correct since the values of the low loss spectrum are usually much larger than 1. Indeed, the `convolution <https://en.wikipedia.org/wiki/Convolution>`_ gives the integral of the pointwise multiplication of the component function with the low-loss spectrum shifted over the component and it results that the function values of the convolved component are much larger than in the non-convolved case, which means that the fitted parameter values need to be very small to fit the data.


HyperSpy has created the model and configured it automatically:

.. code-block:: python
Expand Down
2 changes: 1 addition & 1 deletion hyperspy/_components/arctan.py
Original file line number Diff line number Diff line change
Expand Up @@ -59,7 +59,7 @@ def __init__(self, A=1., k=1., x0=1., minimum_at_zero=False):
self._whitelist['minimum_at_zero'] = ('init', minimum_at_zero)

self.isbackground = False
self.isconvolved = False
self.convolved = False
self._position = self.x0

def function(self, x):
Expand Down
1 change: 0 additions & 1 deletion hyperspy/_components/eels_cl_edge.py
Original file line number Diff line number Diff line change
Expand Up @@ -23,7 +23,6 @@
import numpy as np
from scipy.interpolate import splev

from hyperspy.defaults_parser import preferences
from hyperspy.component import Component
from hyperspy.misc.eels.hartree_slater_gos import HartreeSlaterGOS
from hyperspy.misc.eels.hydrogenic_gos import HydrogenicGOS
Expand Down
17 changes: 4 additions & 13 deletions hyperspy/component.py
Original file line number Diff line number Diff line change
Expand Up @@ -17,8 +17,6 @@
# along with HyperSpy. If not, see <http://www.gnu.org/licenses/>.

import os
import functools
import warnings

import numpy as np
from dask.array import Array as dArray
Expand All @@ -27,15 +25,13 @@
import sympy
from sympy.utilities.lambdify import lambdify

from hyperspy.defaults_parser import preferences
from hyperspy.misc.utils import slugify
from hyperspy.misc.io.tools import (incremental_filename,
append2pathname,)
from hyperspy.exceptions import NavigationDimensionError
from hyperspy.misc.export_dictionary import export_to_dictionary, \
load_from_dictionary
from hyperspy.events import Events, Event
from hyperspy.ui_registry import add_gui_method, register_toolkey
from hyperspy.ui_registry import add_gui_method

import logging

Expand Down Expand Up @@ -207,7 +203,7 @@ def _load_dictionary(self, dictionary):
load_from_dictionary(self, dictionary)
return dictionary['self']
else:
raise ValueError( "_id_name of parameter and dictionary do not match, \nparameter._id_name = %s\
raise ValueError("_id_name of parameter and dictionary do not match, \nparameter._id_name = %s\
\ndictionary['_id_name'] = %s" % (self._id_name, dictionary['_id_name']))

def __repr__(self):
Expand Down Expand Up @@ -1013,10 +1009,7 @@ def __call__(self):
-------
numpy array
"""

axis = self.model.axis.axis[self.model.channel_switches]
component_array = self.function(axis)
return component_array
return self.model.__call__(component_list=[self])

def _component2plot(self, axes_manager, out_of_range2nans=True):
old_axes_manager = None
Expand All @@ -1027,8 +1020,6 @@ def _component2plot(self, axes_manager, out_of_range2nans=True):
s = self.__call__()
if not self.active:
s.fill(np.nan)
if self.model.signal.metadata.Signal.binned is True:
s *= self.model.signal.axes_manager.signal_axes[0].scale
if old_axes_manager is not None:
self.model.axes_manager = old_axes_manager
self.charge()
Expand Down Expand Up @@ -1183,5 +1174,5 @@ def _load_dictionary(self, dic):
"_id_name of parameters in component and dictionary do not match")
return id_dict
else:
raise ValueError( "_id_name of component and dictionary do not match, \ncomponent._id_name = %s\
raise ValueError("_id_name of component and dictionary do not match, \ncomponent._id_name = %s\
\ndictionary['_id_name'] = %s" % (self._id_name, dic['_id_name']))
5 changes: 0 additions & 5 deletions hyperspy/model.py
Original file line number Diff line number Diff line change
Expand Up @@ -631,11 +631,6 @@ def _close_plot(self):
self._disconnect_parameters2update_plot(components=self)
self._model_line = None

def _update_model_line(self):
if (self._plot_active is True and
self._model_line is not None):
self._model_line.update()

@staticmethod
def _connect_component_line(component):
if hasattr(component, "_model_plot_line"):
Expand Down
2 changes: 0 additions & 2 deletions hyperspy/models/eelsmodel.py
Original file line number Diff line number Diff line change
Expand Up @@ -23,7 +23,6 @@
from hyperspy.models.model1d import Model1D
from hyperspy.components1d import EELSCLEdge
from hyperspy.components1d import PowerLaw
from hyperspy.defaults_parser import preferences
from hyperspy import components1d
from hyperspy._signals.eels import EELSSpectrum

Expand Down Expand Up @@ -494,7 +493,6 @@ def _fit_edge(self, edgenumber, start_energy=None, **kwargs):
ea = self.axis.axis[self.channel_switches]
if start_energy is None:
start_energy = ea[0]
preedge_safe_window_width = self._preedge_safe_window_width
# Declare variables
active_edges = self._active_edges
edge = active_edges[edgenumber]
Expand Down
49 changes: 24 additions & 25 deletions hyperspy/models/model1d.py
Original file line number Diff line number Diff line change
Expand Up @@ -17,18 +17,15 @@
# along with HyperSpy. If not, see <http://www.gnu.org/licenses/>.

import copy
from contextlib import contextmanager

import numpy as np
import traits.api as t

from hyperspy.model import BaseModel, ModelComponents, ModelSpecialSlicers
import hyperspy.drawing.signal1d
from hyperspy.axes import generate_axis
from hyperspy.exceptions import WrongObjectError
from hyperspy.exceptions import WrongObjectError, SignalDimensionError
from hyperspy.decorators import interactive_range_selector
from hyperspy.drawing.widgets import VerticalLineWidget, LabelWidget
from hyperspy.ui_registry import get_gui
from hyperspy.events import EventSuppressor
from hyperspy.signal_tools import SpanSelectorInSignal1D
from hyperspy.ui_registry import add_gui_method, DISPLAY_DT, TOOLKIT_DT
Expand Down Expand Up @@ -377,7 +374,8 @@ def remove(self, things):

remove.__doc__ = BaseModel.remove.__doc__

def __call__(self, non_convolved=False, onlyactive=False):
def __call__(self, non_convolved=False, onlyactive=False,
component_list=None):
"""Returns the corresponding model for the current coordinates
Parameters
Expand All @@ -387,6 +385,9 @@ def __call__(self, non_convolved=False, onlyactive=False):
only_active : bool
If True, only the active components will be used to build the
model.
component_list : list or None
If None, the sum of all the components is returned. If list, only
the provided components are returned
cursor: 1 or 2
Expand All @@ -395,35 +396,33 @@ def __call__(self, non_convolved=False, onlyactive=False):
numpy array
"""

if component_list is None:
component_list = [component for component in self]
if isinstance(component_list, (list, tuple)):
pass
else:
raise ValueError(
"'Component_list' parameter need to be a list or None")

if onlyactive:
component_list = [
component for component in component_list if component.active]

if self.convolved is False or non_convolved is True:
axis = self.axis.axis[self.channel_switches]
sum_ = np.zeros(len(axis))
if onlyactive is True:
for component in self:
if component.active:
sum_ += component.function(axis)
else:
for component in self:
sum_ += component.function(axis)
for component in component_list:
sum_ += component.function(axis)
to_return = sum_

else: # convolved
sum_convolved = np.zeros(len(self.convolution_axis))
sum_ = np.zeros(len(self.axis.axis))
for component in self: # Cut the parameters list
if onlyactive:
if component.active:
if component.convolved:
sum_convolved += component.function(
self.convolution_axis)
else:
sum_ += component.function(self.axis.axis)
for component in component_list:
if component.convolved:
sum_convolved += component.function(self.convolution_axis)
else:
if component.convolved:
sum_convolved += component.function(
self.convolution_axis)
else:
sum_ += component.function(self.axis.axis)
sum_ += component.function(self.axis.axis)

to_return = sum_ + np.convolve(
self.low_loss(self.axes_manager),
Expand Down
Binary file added hyperspy/tests/drawing/data/Cr_L_cl.hspy
Binary file not shown.
Binary file added hyperspy/tests/drawing/data/Cr_L_ll.hspy
Binary file not shown.
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
150 changes: 150 additions & 0 deletions hyperspy/tests/drawing/test_plot_model.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,150 @@
# -*- coding: utf-8 -*-
# Copyright 2007-2018 The HyperSpy developers
#
# This file is part of HyperSpy.
#
# HyperSpy is free software: you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# HyperSpy is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with HyperSpy. If not, see <http://www.gnu.org/licenses/>.

import os
import numpy as np
import pytest
import numpy.testing as nt

import hyperspy.api as hs
from hyperspy.signals import Signal1D, EELSSpectrum
from hyperspy.components1d import Gaussian


my_path = os.path.dirname(__file__)
baseline_dir = 'plot_model'
default_tol = 2.0


def create_ll_signal(signal_shape=1000):
offset = 0
zlp_param = {'A': 10000.0, 'centre': 0.0 + offset, 'sigma': 15.0}
zlp = Gaussian(**zlp_param)
plasmon_param = {'A': 2000.0, 'centre': 200.0 + offset, 'sigma': 75.0}
plasmon = Gaussian(**plasmon_param)
axis = np.arange(signal_shape)
data = zlp.function(axis) + plasmon.function(axis)
ll = EELSSpectrum(data)
ll.axes_manager[-1].offset = -offset
ll.axes_manager[-1].scale = 0.1
return ll


A_value_gaussian = [1000.0, 600.0, 2000.0]
centre_value_gaussian = [50.0, 20.0, 60.0]
sigma_value_gaussian = [5.0, 3.0, 1.0]
scale = 0.1


def create_sum_of_gaussians(convolved=False):
param1 = {'A': A_value_gaussian[0],
'centre': centre_value_gaussian[0]/scale,
'sigma': sigma_value_gaussian[0]/scale}
gs1 = Gaussian(**param1)
param2 = {'A': A_value_gaussian[1],
'centre': centre_value_gaussian[1]/scale,
'sigma': sigma_value_gaussian[1]/scale}
gs2 = Gaussian(**param2)
param3 = {'A': A_value_gaussian[2],
'centre': centre_value_gaussian[2]/scale,
'sigma': sigma_value_gaussian[2]/scale}
gs3 = Gaussian(**param3)

axis = np.arange(1000)
data = gs1.function(axis) + gs2.function(axis) + gs3.function(axis)

if convolved:
to_convolved = create_ll_signal(data.shape[0]).data
data = np.convolve(data, to_convolved) / sum(to_convolved)

s = Signal1D(data[:1000])
s.axes_manager[-1].scale = scale
return s


def _generate_parameters():
parameters = []
for convolved in [True, False]:
for plot_component in [True, False]:
for binned in [True, False]:
parameters.append([convolved, plot_component, binned])
return parameters


@pytest.mark.parametrize(("convolved", "plot_component", "binned"),
_generate_parameters())
@pytest.mark.mpl_image_compare(
baseline_dir=baseline_dir, tolerance=default_tol)
def test_plot_gaussian_eelsmodel(convolved, plot_component, binned):
s = create_sum_of_gaussians(convolved)
s.set_signal_type('EELS')
s.metadata.General.title = 'Convolved: {}, plot_component: {}, binned: {}'.format(
convolved, plot_component, binned)

ll = create_ll_signal(1000) if convolved else None

s.set_microscope_parameters(200, 20, 50)
s.metadata.Signal.binned = binned
m = s.create_model(auto_background=False, ll=ll)

m.extend([Gaussian(), Gaussian(), Gaussian()])

def set_gaussian(gaussian, centre, sigma):
gaussian.centre.value = centre
gaussian.centre.free = False
gaussian.sigma.value = sigma
gaussian.sigma.free = False

for gaussian, centre, sigma in zip(m, centre_value_gaussian,
sigma_value_gaussian):
set_gaussian(gaussian, centre, sigma)

m.fit()
m.plot(plot_components=plot_component)

def A_value(s, component, binned):
if binned:
return component.A.value * scale
else:
return component.A.value

if convolved:
nt.assert_almost_equal(A_value(s, m[0], binned), 0.014034, decimal=5)
nt.assert_almost_equal(A_value(s, m[1], binned), 0.008420, decimal=5)
nt.assert_almost_equal(A_value(s, m[2], binned), 0.028068, decimal=5)
else:
nt.assert_almost_equal(A_value(s, m[0], binned), 100.0)
nt.assert_almost_equal(A_value(s, m[1], binned), 60.0)
nt.assert_almost_equal(A_value(s, m[2], binned), 200.0)

return m._plot.signal_plot.figure


@pytest.mark.parametrize(("convolved"), [False, True])
@pytest.mark.mpl_image_compare(
baseline_dir=baseline_dir, tolerance=default_tol)
def test_fit_EELS_convolved(convolved):
dname = os.path.join(my_path, 'data')
cl = hs.load(os.path.join(dname, 'Cr_L_cl.hspy'))
cl.metadata.Signal.binned = False
cl.metadata.General.title = 'Convolved: {}'.format(convolved)
ll = hs.load(os.path.join(dname, 'Cr_L_ll.hspy')) if convolved else None
m = cl.create_model(auto_background=False, ll=ll, GOS='hydrogenic')
m.fit(kind='smart')
m.plot(plot_components=True)
return m._plot.signal_plot.figure

0 comments on commit d8d9039

Please sign in to comment.