Skip to content

Commit

Permalink
Merge pull request #19 from marshallmcdonnell/safe_divide
Browse files Browse the repository at this point in the history
Fixes for PyStoG command-line program
  • Loading branch information
marshallmcdonnell authored Dec 11, 2018
2 parents 36aeacf + 0bb721a commit 7d22967
Show file tree
Hide file tree
Showing 6 changed files with 82 additions and 88 deletions.
4 changes: 2 additions & 2 deletions bin/pystog_cli
Original file line number Diff line number Diff line change
Expand Up @@ -54,8 +54,8 @@ def main(kwargs=None):
r, gr_out = stog.apply_lorch(q, sq, r)

# Apply final scale number
stog.add_keen_fq(q, sq)
stog.add_keen_gr(r, gr_out)
stog._add_keen_fq(q, sq)
stog._add_keen_gr(r, gr_out)

if "PlotFlag" in kwargs:
if kwargs["PlotFlag"]:
Expand Down
17 changes: 10 additions & 7 deletions pystog/converter.py
Original file line number Diff line number Diff line change
Expand Up @@ -31,6 +31,12 @@ class Converter:
def __init__(self):
pass

def _safe_divide(self, numerator, denominator):
mask = (denominator != 0.0)
out = np.zeros_like(numerator)
out[mask] = numerator[mask] / denominator[mask]
return out

# Reciprocal Space Conversions

def F_to_S(self, q, fq, **kwargs):
Expand All @@ -44,10 +50,7 @@ def F_to_S(self, q, fq, **kwargs):
:return: :math:`S(Q)` vector
:rtype: numpy.array
"""
mask = (q != 0.0)
fq_new = np.zeros_like(fq)
fq_new[mask] = fq[mask] / q[mask]
return fq_new + 1.
return self._safe_divide(fq, q) + 1.

def F_to_FK(self, q, fq, **kwargs):
"""Converts from :math:`Q[S(Q)-1]` to :math:`F(Q)`
Expand All @@ -63,7 +66,7 @@ def F_to_FK(self, q, fq, **kwargs):
mask = (q != 0.0)
fq_new = np.zeros_like(fq)
fq_new[mask] = fq[mask] / q[mask]
return kwargs['<b_coh>^2'] * fq_new
return kwargs['<b_coh>^2'] * self._safe_divide(fq, q)

def F_to_DCS(self, q, fq, **kwargs):
"""Converts from :math:`Q[S(Q)-1]` to :math:`\\frac{d \\sigma}{d \\Omega}(Q)`
Expand Down Expand Up @@ -219,7 +222,7 @@ def G_to_GK(self, r, gr, **kwargs):
:rtype: numpy.array
"""
factor = kwargs['<b_coh>^2'] / (4. * np.pi * kwargs['rho'])
return factor * (gr / r)
return factor * self._safe_divide(gr, r)

def G_to_g(self, r, gr, **kwargs):
"""Convert :math:`G_{PDFFIT}(r)` to :math:`g(r)`
Expand All @@ -233,7 +236,7 @@ def G_to_g(self, r, gr, **kwargs):
:rtype: numpy.array
"""
factor = 4. * np.pi * kwargs['rho']
return gr / (factor * r) + 1.
return self._safe_divide(gr, factor * r) + 1.

# Keen's G(r)
def GK_to_G(self, r, gr, **kwargs):
Expand Down
3 changes: 1 addition & 2 deletions pystog/fourier_filter.py
Original file line number Diff line number Diff line change
Expand Up @@ -75,8 +75,7 @@ def g_using_F(self, r, gr, q, fq, cutoff, **kwargs):
gr_tmp = gr_tmp_initial + 1

# Transform the shifted low-r region to F(Q) to get F(Q)_ft
q_ft = self.transformer._extend_axis_to_low_end(q)
q_ft, fq_ft = self.transformer.g_to_F(r_tmp, gr_tmp, q_ft, **kwargs)
q_ft, fq_ft = self.transformer.g_to_F(r_tmp, gr_tmp, q, **kwargs)
q_ft, fq_ft = self.transformer.apply_cropping(q_ft, fq_ft, qmin, qmax)

# Subtract F(Q)_ft from original F(Q) = delta_F(Q)
Expand Down
16 changes: 0 additions & 16 deletions pystog/transformer.py
Original file line number Diff line number Diff line change
Expand Up @@ -11,7 +11,6 @@
from __future__ import (absolute_import, division, print_function)
import numpy as np

from pystog.utils import create_domain
from pystog.converter import Converter

# -------------------------------------------------------#
Expand Down Expand Up @@ -39,18 +38,6 @@ class Transformer:
def __init__(self):
self.converter = Converter()

def _extend_axis_to_low_end(self, x, decimals=4):
"""Utility to setup axis for the forward transform space.
:param x: vector
:type x: numpy.array or list
:param decimals: max decimals in output
:type decimals: integer
:return: vector
:rtype: numpy.array
"""
return create_domain(min(x), max(x), x[1] - x[0])

def _low_x_correction(self, xin, yin, xout, yout, **kwargs):
"""Omitted low-x range correction performed in the
space you have transformed to. Does so by assumming a
Expand Down Expand Up @@ -170,8 +157,6 @@ def fourier_transform(self, xin, yin, xout,

xin, yin = self.apply_cropping(xin, yin, xmin, xmax)

xout = self._extend_axis_to_low_end(xout)

factor = np.full_like(yin, 1.0)
if 'lorch' in kwargs:
if kwargs['lorch']:
Expand Down Expand Up @@ -434,7 +419,6 @@ def G_to_F(self, r, gr, q, **kwargs):
:return: :math:`Q` and :math:`Q[S(Q)-1]` vector pair
:rtype: numpy.array pair
"""
q = self._extend_axis_to_low_end(q)
q, fq = self.fourier_transform(r, gr, q, **kwargs)
return q, fq

Expand Down
124 changes: 69 additions & 55 deletions tests/test_converter.py
Original file line number Diff line number Diff line change
@@ -1,12 +1,26 @@
import unittest
import numpy
import numpy as np
from tests.utils import \
load_data, get_index_of_function
from tests.materials import Nickel, Argon
from pystog.utils import \
RealSpaceHeaders, ReciprocalSpaceHeaders
from pystog.converter import Converter


class TestConverterUtilities(unittest.TestCase):
def setUp(self):
unittest.TestCase.setUp(self)
self.converter = Converter()

def tearDown(self):
unittest.TestCase.tearDown(self)

def test_safe_divide(self):
self.assertTrue(np.array_equal(self.converter._safe_divide(np.arange(10), np.arange(10)),
[0, 1, 1, 1, 1, 1, 1, 1, 1, 1]))


# Real Space Function


Expand Down Expand Up @@ -43,41 +57,41 @@ def tearDown(self):
# g(r) tests
def g_to_G(self):
GofR = self.converter.g_to_G(self.r, self.gofr, **self.kwargs)
self.assertTrue(numpy.allclose(GofR[self.first:self.last],
self.GofR_target,
rtol=self.rtol, atol=self.atol))
self.assertTrue(np.allclose(GofR[self.first:self.last],
self.GofR_target,
rtol=self.rtol, atol=self.atol))

def g_to_GK(self):
GKofR = self.converter.g_to_GK(self.r, self.gofr, **self.kwargs)
self.assertTrue(numpy.allclose(GKofR[self.first:self.last],
self.GKofR_target,
rtol=self.rtol, atol=self.atol))
self.assertTrue(np.allclose(GKofR[self.first:self.last],
self.GKofR_target,
rtol=self.rtol, atol=self.atol))

# G(r) tests
def G_to_g(self):
gofr = self.converter.G_to_g(self.r, self.GofR, **self.kwargs)
self.assertTrue(numpy.allclose(gofr[self.first:self.last],
self.gofr_target,
rtol=self.rtol, atol=self.atol))
self.assertTrue(np.allclose(gofr[self.first:self.last],
self.gofr_target,
rtol=self.rtol, atol=self.atol))

def G_to_GK(self):
GKofR = self.converter.G_to_GK(self.r, self.GofR, **self.kwargs)
self.assertTrue(numpy.allclose(GKofR[self.first:self.last],
self.GKofR_target,
rtol=self.rtol, atol=self.atol))
self.assertTrue(np.allclose(GKofR[self.first:self.last],
self.GKofR_target,
rtol=self.rtol, atol=self.atol))

# GK(r) tests
def GK_to_g(self):
gofr = self.converter.GK_to_g(self.r, self.GKofR, **self.kwargs)
self.assertTrue(numpy.allclose(gofr[self.first:self.last],
self.gofr_target,
rtol=self.rtol, atol=self.atol))
self.assertTrue(np.allclose(gofr[self.first:self.last],
self.gofr_target,
rtol=self.rtol, atol=self.atol))

def GK_to_G(self):
GofR = self.converter.GK_to_G(self.r, self.GKofR, **self.kwargs)
self.assertTrue(numpy.allclose(GofR[self.first:self.last],
self.GofR_target,
rtol=self.rtol, atol=self.atol))
self.assertTrue(np.allclose(GofR[self.first:self.last],
self.GofR_target,
rtol=self.rtol, atol=self.atol))


class TestConverterRealSpaceNickel(TestConverterRealSpaceBase):
Expand Down Expand Up @@ -171,78 +185,78 @@ def tearDown(self):
# S(Q) tests
def S_to_F(self):
fq = self.converter.S_to_F(self.q, self.sq, **self.kwargs)
self.assertTrue(numpy.allclose(fq[self.first:self.last],
self.fq_target,
rtol=self.rtol, atol=self.atol))
self.assertTrue(np.allclose(fq[self.first:self.last],
self.fq_target,
rtol=self.rtol, atol=self.atol))

def S_to_FK(self):
fq_keen = self.converter.S_to_FK(self.q, self.sq, **self.kwargs)
self.assertTrue(numpy.allclose(fq_keen[self.first:self.last],
self.fq_keen_target,
rtol=self.rtol, atol=self.atol))
self.assertTrue(np.allclose(fq_keen[self.first:self.last],
self.fq_keen_target,
rtol=self.rtol, atol=self.atol))

def S_to_DCS(self):
dcs = self.converter.S_to_DCS(self.q, self.sq, **self.kwargs)
self.assertTrue(numpy.allclose(dcs[self.first:self.last],
self.dcs_target,
rtol=self.rtol, atol=self.atol))
self.assertTrue(np.allclose(dcs[self.first:self.last],
self.dcs_target,
rtol=self.rtol, atol=self.atol))
# Q[S(Q)-1] tests

def F_to_S(self):
sq = self.converter.F_to_S(self.q, self.fq, **self.kwargs)
self.assertTrue(numpy.allclose(sq[self.first:self.last],
self.sq_target,
rtol=self.rtol, atol=self.atol))
self.assertTrue(np.allclose(sq[self.first:self.last],
self.sq_target,
rtol=self.rtol, atol=self.atol))

def F_to_FK(self):
fq_keen = self.converter.F_to_FK(self.q, self.fq, **self.kwargs)
self.assertTrue(numpy.allclose(fq_keen[self.first:self.last],
self.fq_keen_target,
rtol=self.rtol, atol=self.atol))
self.assertTrue(np.allclose(fq_keen[self.first:self.last],
self.fq_keen_target,
rtol=self.rtol, atol=self.atol))

def F_to_DCS(self):
dcs = self.converter.F_to_DCS(self.q, self.fq, **self.kwargs)
self.assertTrue(numpy.allclose(dcs[self.first:self.last],
self.dcs_target,
rtol=self.rtol, atol=self.atol))
self.assertTrue(np.allclose(dcs[self.first:self.last],
self.dcs_target,
rtol=self.rtol, atol=self.atol))
# FK(Q) tests

def FK_to_S(self):
sq = self.converter.FK_to_S(self.q, self.fq_keen, **self.kwargs)
self.assertTrue(numpy.allclose(sq[self.first:self.last],
self.sq_target,
rtol=self.rtol, atol=self.atol))
self.assertTrue(np.allclose(sq[self.first:self.last],
self.sq_target,
rtol=self.rtol, atol=self.atol))

def FK_to_F(self):
fq = self.converter.FK_to_F(self.q, self.fq_keen, **self.kwargs)
self.assertTrue(numpy.allclose(fq[self.first:self.last],
self.fq_target,
rtol=self.rtol, atol=self.atol))
self.assertTrue(np.allclose(fq[self.first:self.last],
self.fq_target,
rtol=self.rtol, atol=self.atol))

def FK_to_DCS(self):
dcs = self.converter.FK_to_DCS(self.q, self.fq_keen, **self.kwargs)
self.assertTrue(numpy.allclose(dcs[self.first:self.last],
self.dcs_target,
rtol=self.rtol, atol=self.atol))
self.assertTrue(np.allclose(dcs[self.first:self.last],
self.dcs_target,
rtol=self.rtol, atol=self.atol))
# DCS(Q) tests

def DCS_to_S(self):
sq = self.converter.DCS_to_S(self.q, self.dcs, **self.kwargs)
self.assertTrue(numpy.allclose(sq[self.first:self.last],
self.sq_target,
rtol=self.rtol, atol=self.atol))
self.assertTrue(np.allclose(sq[self.first:self.last],
self.sq_target,
rtol=self.rtol, atol=self.atol))

def DCS_to_F(self):
fq = self.converter.DCS_to_F(self.q, self.dcs, **self.kwargs)
self.assertTrue(numpy.allclose(fq[self.first:self.last],
self.fq_target,
rtol=self.rtol, atol=self.atol))
self.assertTrue(np.allclose(fq[self.first:self.last],
self.fq_target,
rtol=self.rtol, atol=self.atol))

def DCS_to_FK(self):
fq_keen = self.converter.DCS_to_FK(self.q, self.dcs, **self.kwargs)
self.assertTrue(numpy.allclose(fq_keen[self.first:self.last],
self.fq_keen_target,
rtol=self.rtol, atol=self.atol))
self.assertTrue(np.allclose(fq_keen[self.first:self.last],
self.fq_keen_target,
rtol=self.rtol, atol=self.atol))


class TestConverterReciprocalSpaceNickel(TestConverterReciprocalSpaceBase):
Expand Down
6 changes: 0 additions & 6 deletions tests/test_transformer.py
Original file line number Diff line number Diff line change
Expand Up @@ -73,12 +73,6 @@ def tearDown(self):

# Utilities

def test_extend_axis_to_low_end(self):
xin = numpy.linspace(0.5, 1.0, 11)
xout = self.transformer._extend_axis_to_low_end(xin)
self.assertAlmostEqual(xout[0], 0.5)
self.assertAlmostEqual(xout[-1], 1.0)

def test_apply_cropping(self):
xin = numpy.linspace(0.5, 1.0, 11)
yin = numpy.linspace(4.5, 5.0, 11)
Expand Down

0 comments on commit 7d22967

Please sign in to comment.