Skip to content
Merged
Show file tree
Hide file tree
Changes from 1 commit
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
23 changes: 23 additions & 0 deletions news/extrap-warnings.rst
Original file line number Diff line number Diff line change
@@ -0,0 +1,23 @@
**Added:**

* Enable ``diffpy.morph`` to detect extrapolation.

**Changed:**

* <news item>

**Deprecated:**

* <news item>

**Removed:**

* <news item>

**Fixed:**

* <news item>

**Security:**

* <news item>
23 changes: 13 additions & 10 deletions src/diffpy/morph/morph_io.py
Original file line number Diff line number Diff line change
Expand Up @@ -410,27 +410,30 @@ def tabulate_results(multiple_morph_results):

def handle_warnings(squeeze_morph):
if squeeze_morph is not None:
eil = squeeze_morph.extrap_index_low
eih = squeeze_morph.extrap_index_high

if eil is not None or eih is not None:
if eih is None:
extrapolation_info = squeeze_morph.extrapolation_info
is_extrap_low = extrapolation_info["is_extrap_low"]
is_extrap_high = extrapolation_info["is_extrap_high"]
cutoff_low = extrapolation_info["cutoff_low"]
cutoff_high = extrapolation_info["cutoff_high"]

if is_extrap_low or is_extrap_high:
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This is getting better, but I still think it would be more readable with something like:

if low and high:
    low|high flow
elif low:
   low flow
elif high:
   high flow

if not is_extrap_high:
wmsg = (
"Warning: points with grid value below "
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Change "will be extrapolated" to "are extrapolated"

f"{squeeze_morph.squeeze_cutoff_low} "
f"{cutoff_low} "
f"will be extrapolated."
)
elif eil is None:
elif not is_extrap_low:
wmsg = (
"Warning: points with grid value above "
f"{squeeze_morph.squeeze_cutoff_high} "
f"{cutoff_high} "
f"will be extrapolated."
)
else:
wmsg = (
"Warning: points with grid value below "
f"{squeeze_morph.squeeze_cutoff_low} and above "
f"{squeeze_morph.squeeze_cutoff_high} "
f"{cutoff_low} and above "
f"{cutoff_high} "
f"will be extrapolated."
)
warnings.warn(
Expand Down
5 changes: 4 additions & 1 deletion src/diffpy/morph/morphapp.py
Original file line number Diff line number Diff line change
Expand Up @@ -610,10 +610,12 @@ def single_morph(
config["smear"] = smear_in
# Shift
# Only enable hshift is squeeze is not enabled
shift_morph = None
if (
opts.hshift is not None and squeeze_poly_deg < 0
) or opts.vshift is not None:
chain.append(morphs.MorphShift())
shift_morph = morphs.MorphShift()
chain.append(shift_morph)
if opts.hshift is not None and squeeze_poly_deg < 0:
hshift_in = opts.hshift
config["hshift"] = hshift_in
Expand Down Expand Up @@ -700,6 +702,7 @@ def single_morph(

# THROW ANY WARNINGS HERE
io.handle_warnings(squeeze_morph)
io.handle_warnings(shift_morph)

# Get Rw for the morph range
rw = tools.getRw(chain)
Expand Down
20 changes: 18 additions & 2 deletions src/diffpy/morph/morphs/morph.py
Original file line number Diff line number Diff line change
Expand Up @@ -12,8 +12,7 @@
# See LICENSE.txt for license information.
#
##############################################################################
"""Morph -- base class for defining a morph.
"""
"""Morph -- base class for defining a morph."""


LABEL_RA = "r (A)" # r-grid
Expand Down Expand Up @@ -246,6 +245,23 @@ def plotOutputs(self, xylabels=True, **plotargs):
ylabel(self.youtlabel)
return rv

def checkExtrapolation(self, x_true, x_extrapolate):
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

make this a private function if it is only used within this module, or give it a docstring.

Also, just to chec before we put this in the base class, is this method called in multiple sub-classes?

Finally, the function name seems to be a bit off. It seems as if it is a setter, not checking anything. Something like set_extrapolation_info() or something like that.

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The method will be called in morphsqueeze, morphshift, and morphstretch. I will add its docstring.

Yes, set_extrapolation_info will be better.

import numpy
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

is there a reason we are importing numpy locally? Put it at the top of the file?

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

No. I will put it at the top of the file.


cutoff_low = min(x_true)
cutoff_high = max(x_true)
low_extrap = numpy.where(x_extrapolate < cutoff_low)[0]
high_extrap = numpy.where(x_extrapolate > cutoff_high)[0]
is_extrap_low = False if len(low_extrap) == 0 else True
is_extrap_high = False if len(high_extrap) == 0 else True
extrapolation_info = {
"is_extrap_low": is_extrap_low,
"cutoff_low": cutoff_low,
"is_extrap_high": is_extrap_high,
"cutoff_high": cutoff_high,
}
return extrapolation_info

def __getattr__(self, name):
"""Obtain the value from self.config, when normal lookup fails.

Expand Down
1 change: 1 addition & 0 deletions src/diffpy/morph/morphs/morphshift.py
Original file line number Diff line number Diff line change
Expand Up @@ -57,6 +57,7 @@ def morph(self, x_morph, y_morph, x_target, y_target):
r = self.x_morph_in - hshift
self.y_morph_out = numpy.interp(r, self.x_morph_in, self.y_morph_in)
self.y_morph_out += vshift
self.extrapolation_info = self.checkExtrapolation(self.x_morph_in, r)
return self.xyallout


Expand Down
10 changes: 3 additions & 7 deletions src/diffpy/morph/morphs/morphsqueeze.py
Original file line number Diff line number Diff line change
@@ -1,7 +1,6 @@
"""Class MorphSqueeze -- Apply a polynomial to squeeze the morph
function."""

import numpy as np
from numpy.polynomial import Polynomial
from scipy.interpolate import CubicSpline

Expand Down Expand Up @@ -83,14 +82,11 @@ def morph(self, x_morph, y_morph, x_target, y_target):
coeffs = [self.squeeze[f"a{i}"] for i in range(len(self.squeeze))]
squeeze_polynomial = Polynomial(coeffs)
x_squeezed = self.x_morph_in + squeeze_polynomial(self.x_morph_in)
self.squeeze_cutoff_low = min(x_squeezed)
self.squeeze_cutoff_high = max(x_squeezed)
self.y_morph_out = CubicSpline(x_squeezed, self.y_morph_in)(
self.x_morph_in
)
low_extrap = np.where(self.x_morph_in < self.squeeze_cutoff_low)[0]
high_extrap = np.where(self.x_morph_in > self.squeeze_cutoff_high)[0]
self.extrap_index_low = low_extrap[-1] if low_extrap.size else None
self.extrap_index_high = high_extrap[0] if high_extrap.size else None
self.extrapolation_info = self.checkExtrapolation(
x_squeezed, self.x_morph_in
)

return self.xyallout
26 changes: 14 additions & 12 deletions tests/test_morphsqueeze.py
Original file line number Diff line number Diff line change
Expand Up @@ -46,23 +46,27 @@
@pytest.mark.parametrize("squeeze_coeffs", squeeze_coeffs_dic)
def test_morphsqueeze(x_morph, x_target, squeeze_coeffs):
y_target = np.sin(x_target)
y_morph = np.sin(x_morph)
# expected output
coeffs = [squeeze_coeffs[f"a{i}"] for i in range(len(squeeze_coeffs))]
squeeze_polynomial = Polynomial(coeffs)
x_squeezed = x_morph + squeeze_polynomial(x_morph)
y_morph = np.sin(x_squeezed)
low_extrap = np.where(x_morph < x_squeezed[0])[0]
high_extrap = np.where(x_morph > x_squeezed[-1])[0]
extrap_index_low_expected = low_extrap[-1] if low_extrap.size else None
extrap_index_high_expected = high_extrap[0] if high_extrap.size else None
y_morph_expected = y_morph
x_morph_expected = x_morph
y_morph_expected = np.sin(x_morph)
x_target_expected = x_target
y_target_expected = y_target
# actual output
morph = MorphSqueeze()
y_morph = np.sin(x_squeezed)
morph.squeeze = squeeze_coeffs
x_morph_actual, y_morph_actual, x_target_actual, y_target_actual = morph(
x_morph, y_morph, x_target, y_target
)
extrap_index_low = morph.extrap_index_low
extrap_index_high = morph.extrap_index_high

extrap_low = np.where(x_morph < min(x_squeezed))[0]
Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The code to compute extrap_index_*_expected and extrap_index_*_actual is identical, so we can remove the check and keep only one of them. In this case, extrap_index_*_actual is removed, which makes the implementation in the main program simpler.

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Is there an argument to be made that this work could be done in checkExtrapolation()?

extrap_high = np.where(x_morph > max(x_squeezed))[0]
extrap_index_low = extrap_low[-1] if extrap_low.size else None
extrap_index_high = extrap_high[0] if extrap_high.size else None
if extrap_index_low is None:
extrap_index_low = 0
elif extrap_index_high is None:
Expand All @@ -82,11 +86,9 @@ def test_morphsqueeze(x_morph, x_target, squeeze_coeffs):
y_morph_expected[extrap_index_high:],
atol=1e-3,
)
assert morph.extrap_index_low == extrap_index_low_expected
assert morph.extrap_index_high == extrap_index_high_expected
assert np.allclose(x_morph_actual, x_morph_expected)
assert np.allclose(x_target_actual, x_target)
assert np.allclose(y_target_actual, y_target)
assert np.allclose(x_target_actual, x_target_expected)
assert np.allclose(y_target_actual, y_target_expected)


@pytest.mark.parametrize(
Expand Down
Loading