From a2de08502dc216e2b971415450b48cd598383a60 Mon Sep 17 00:00:00 2001 From: Bradley Dice Date: Wed, 26 May 2021 14:51:13 -0500 Subject: [PATCH 01/29] Rescale k values (validated on cubic box). --- freud/diffraction.pyx | 7 ++++--- 1 file changed, 4 insertions(+), 3 deletions(-) diff --git a/freud/diffraction.pyx b/freud/diffraction.pyx index d76ada390..248f6150d 100644 --- a/freud/diffraction.pyx +++ b/freud/diffraction.pyx @@ -67,6 +67,7 @@ cdef class DiffractionPattern(_Compute): cdef unsigned int _frame_counter cdef double _box_matrix_scale_factor cdef double[:] _view_orientation + cdef double _zoom cdef cbool _k_values_cached cdef cbool _k_vectors_cached @@ -259,8 +260,7 @@ cdef class DiffractionPattern(_Compute): if not self._called_compute: # Create a 1D axis of k-vector magnitudes self._k_values_orig = np.fft.fftshift(np.fft.fftfreq( - n=self.output_size, - d=1/self.output_size)) + n=self.output_size)) # Create a 3D meshgrid of k-vectors with shape # (output_size, output_size, 3) @@ -271,6 +271,7 @@ cdef class DiffractionPattern(_Compute): # lazy evaluation of k-values and k-vectors self._box_matrix_scale_factor = np.max(system.box.to_matrix()) self._view_orientation = view_orientation + self._zoom = zoom self._k_values_cached = False self._k_vectors_cached = False @@ -299,7 +300,7 @@ cdef class DiffractionPattern(_Compute): """(``output_size``, ) :class:`numpy.ndarray`: k-values.""" if not self._k_values_cached: self._k_values = np.asarray( - self._k_values_orig) / self._box_matrix_scale_factor + self._k_values_orig) * (2 * np.pi * self.output_size / (self._box_matrix_scale_factor * self._zoom)) self._k_values_cached = True return np.asarray(self._k_values) From e87fca1a6d746ccfd742e039c2cfd2c3b9f2abdf Mon Sep 17 00:00:00 2001 From: Bradley Dice Date: Wed, 26 May 2021 15:14:47 -0500 Subject: [PATCH 02/29] Use full grid size. --- freud/diffraction.pyx | 6 ++---- 1 file changed, 2 insertions(+), 4 deletions(-) diff --git a/freud/diffraction.pyx b/freud/diffraction.pyx index 248f6150d..e1987f11a 100644 --- a/freud/diffraction.pyx +++ b/freud/diffraction.pyx @@ -131,7 +131,7 @@ cdef class DiffractionPattern(_Compute): """Zoom, shear, and scale diffraction intensities. Args: - img ((``grid_size//zoom, grid_size//zoom``) :class:`numpy.ndarray`): + img ((``grid_size, grid_size``) :class:`numpy.ndarray`): Array of diffraction intensities. box (:class:`~.box.Box`): Simulation box. @@ -219,8 +219,6 @@ cdef class DiffractionPattern(_Compute): view_orientation = freud.util._convert_array( view_orientation, (4,), np.double) - grid_size = int(self.grid_size / zoom) - # Compute the box projection matrix inv_shear = self._calc_proj(view_orientation, system.box) @@ -233,7 +231,7 @@ cdef class DiffractionPattern(_Compute): xy += 0.5 xy %= 1 im, _, _ = np.histogram2d( - xy[:, 0], xy[:, 1], bins=np.linspace(0, 1, grid_size+1)) + xy[:, 0], xy[:, 1], bins=np.linspace(0, 1, self.grid_size+1)) # Compute FFT and convolve with Gaussian cdef double complex[:, :] diffraction_fft From 1ce917a48a238f61ac809c59c9fb715db453b34b Mon Sep 17 00:00:00 2001 From: Bradley Dice Date: Wed, 26 May 2021 15:17:37 -0500 Subject: [PATCH 03/29] Save scale factor, also to k vectors. --- freud/diffraction.pyx | 9 ++++----- 1 file changed, 4 insertions(+), 5 deletions(-) diff --git a/freud/diffraction.pyx b/freud/diffraction.pyx index e1987f11a..2530a878c 100644 --- a/freud/diffraction.pyx +++ b/freud/diffraction.pyx @@ -67,7 +67,7 @@ cdef class DiffractionPattern(_Compute): cdef unsigned int _frame_counter cdef double _box_matrix_scale_factor cdef double[:] _view_orientation - cdef double _zoom + cdef double _k_scale_factor cdef cbool _k_values_cached cdef cbool _k_vectors_cached @@ -269,7 +269,7 @@ cdef class DiffractionPattern(_Compute): # lazy evaluation of k-values and k-vectors self._box_matrix_scale_factor = np.max(system.box.to_matrix()) self._view_orientation = view_orientation - self._zoom = zoom + self._k_scale_factor = 2 * np.pi * self.output_size / (self._box_matrix_scale_factor * zoom) self._k_values_cached = False self._k_vectors_cached = False @@ -297,8 +297,7 @@ cdef class DiffractionPattern(_Compute): def k_values(self): """(``output_size``, ) :class:`numpy.ndarray`: k-values.""" if not self._k_values_cached: - self._k_values = np.asarray( - self._k_values_orig) * (2 * np.pi * self.output_size / (self._box_matrix_scale_factor * self._zoom)) + self._k_values = np.asarray(self._k_values_orig) * self._k_scale_factor self._k_values_cached = True return np.asarray(self._k_values) @@ -311,7 +310,7 @@ cdef class DiffractionPattern(_Compute): if not self._k_vectors_cached: self._k_vectors = rowan.rotate( self._view_orientation, - self._k_vectors_orig) / self._box_matrix_scale_factor + self._k_vectors_orig) * self._k_scale_factor self._k_vectors_cached = True return np.asarray(self._k_vectors) From 49a60ceb16ecbad51718ad74b56fbe29c437d9d7 Mon Sep 17 00:00:00 2001 From: Bradley Dice Date: Wed, 2 Jun 2021 15:04:39 -0500 Subject: [PATCH 04/29] Update docs. --- freud/diffraction.pyx | 9 +++++++-- 1 file changed, 7 insertions(+), 2 deletions(-) diff --git a/freud/diffraction.pyx b/freud/diffraction.pyx index 2530a878c..90306cd8e 100644 --- a/freud/diffraction.pyx +++ b/freud/diffraction.pyx @@ -46,6 +46,10 @@ cdef class DiffractionPattern(_Compute): as a multiplication in Fourier space. The computed diffraction pattern can be accessed as a square array of shape ``(output_size, output_size)``. + The :math:`\vec{k}=0` peak is always located at index + ``(output_size // 2, output_size // 2)`` and is normalized to have a value + of :math:`S(\vec{k}=0) = 1` (not :math:`N`, a common convention). + This method is based on the implementations in the open-source `GIXStapose application `_ and its predecessor, diffractometer :cite:`Jankowski2017`. @@ -148,7 +152,8 @@ cdef class DiffractionPattern(_Compute): # The adjustments to roll and roll_shift ensure that the peak # corresponding to k=0 is located at exactly # (output_size//2, output_size//2), regardless of whether the grid_size - # and output_size are odd or even. + # and output_size are odd or even. This keeps the peak aligned at the + # center of a single pixel, which should always have the maximum value. roll = img.shape[0] / 2 if img.shape[0] % 2 == 1: @@ -289,7 +294,7 @@ cdef class DiffractionPattern(_Compute): def diffraction(self): """ (``output_size``, ``output_size``) :class:`numpy.ndarray`: - diffraction pattern. + Diffraction pattern. """ return np.asarray(self._diffraction) / self._frame_counter From 2d771926bb7ab50f92ef8fe9074f2912498a4981 Mon Sep 17 00:00:00 2001 From: Bradley Dice Date: Wed, 2 Jun 2021 15:06:36 -0500 Subject: [PATCH 05/29] Update tests with comments and instructions. --- tests/test_diffraction_DiffractionPattern.py | 30 ++++++++++++++++++++ 1 file changed, 30 insertions(+) diff --git a/tests/test_diffraction_DiffractionPattern.py b/tests/test_diffraction_DiffractionPattern.py index 140e22030..6b7f2b043 100644 --- a/tests/test_diffraction_DiffractionPattern.py +++ b/tests/test_diffraction_DiffractionPattern.py @@ -160,7 +160,10 @@ def test_k_values_and_k_vectors(self): dp = freud.diffraction.DiffractionPattern() for size in [2, 5, 10]: + # TODO: npoints doesn't actually matter, just use one point. for npoints in [10, 20, 75]: + # TODO: This isn't using the correct signature for + # make_random_system. Flip npoints and size. box, positions = freud.data.make_random_system(npoints, size) dp.compute((box, positions)) @@ -168,3 +171,30 @@ def test_k_values_and_k_vectors(self): npt.assert_allclose(dp.k_values[output_size // 2], 0) center_index = (output_size // 2, output_size // 2) npt.assert_allclose(dp.k_vectors[center_index], [0, 0, 0]) + # Add tests for the left and right bins. Use the math for + # fftfreq and the _k_scale_factor to write an expression for + # the first and last bins' k-values and k-vectors. + + def test_cubic_system(self): + box, positions = freud.data.UnitCell.sc().generate_system(10) + dp = freud.diffraction.DiffractionPattern() + dp.compute((box, positions)) + # Make sure that the peaks are where we expect them. + # Identify the indices of the highest values in dp.diffraction + # and test that k * R == 2*pi*N for some integer N, all peak vectors k, + # and lattice vectors R. This will be inexact because of binning. + # Perhaps try with a larger resolution like 1024. + # Use npt.assert_allclose. + + def test_cubic_system_parameterized(self): + pass + # Second PR: Same as above test but with different grid_size, + # output_size, and zoom values. + for grid_size in (256, 1024): + for output_size in (255, 256, 1023, 1024): + for zoom in (1, 2.5, 4): + # Ensure that peaks are in the correct locations like above. + # Identify the indices of the highest values in + # dp.diffraction and test that k * R == 2*pi*N for some + # integer N, all peak vectors k, and lattice vectors R. + pass From 3fb56a7fad941a973b3dcff098475ab09b33a55a Mon Sep 17 00:00:00 2001 From: Tim Moore Date: Tue, 8 Jun 2021 13:52:12 -0400 Subject: [PATCH 06/29] Refactor diffraction plot function Prior to this commit, the units of the diffraction plot were in pixels...calling something like `ax.axhline(2*np.pi)` would put a horizontal line on the 6th or 7th row of pixels from the top. This commit changes that behavior so that the units in the image space correlate to units in k-space...now calling `ax.axhline(2*np.pi)` will put a horizontal line at y = 2*pi. This enables easy drawing on the plots, which may be useful, e.g., to make sure peaks in s(k) correspond to the actual lattice spacing in a crystalline system. This change does not affect using the diffraction class to find where the bright spots are, as that is more accurately done using the `k_vectors` attribute of the class. --- freud/plot.py | 20 +++++++++++++------- 1 file changed, 13 insertions(+), 7 deletions(-) diff --git a/freud/plot.py b/freud/plot.py index 5f7fa5552..091f64dd3 100644 --- a/freud/plot.py +++ b/freud/plot.py @@ -11,6 +11,7 @@ try: import matplotlib.pyplot as plt from matplotlib.backends.backend_agg import FigureCanvasAgg + from matplotlib.ticker import MaxNLocator, FormatStrFormatter except ImportError: raise ImportError("matplotlib must be installed for freud.plot.") @@ -522,8 +523,14 @@ def diffraction_plot( # Plot the diffraction image and color bar norm = matplotlib.colors.LogNorm(vmin=vmin, vmax=vmax) + extent = (np.min(k_values), np.max(k_values), + np.min(k_values), np.max(k_values)) im = ax.imshow( - np.clip(diffraction, vmin, vmax), interpolation="nearest", cmap=cmap, norm=norm + np.clip(diffraction, vmin, vmax), + interpolation="nearest", + cmap=cmap, + norm=norm, + extent=extent, ) ax_divider = make_axes_locatable(ax) cax = ax_divider.append_axes("right", size="7%", pad="10%") @@ -539,12 +546,11 @@ def diffraction_plot( ticks = np.linspace(0, grid_size, num_ticks) # Set tick locations and labels - tick_labels = np.interp(ticks, range(grid_size), k_values) - tick_labels = [f"{x:.3g}" for x in tick_labels] - ax.xaxis.set_ticks(ticks) - ax.xaxis.set_ticklabels(tick_labels) - ax.yaxis.set_ticks(ticks) - ax.yaxis.set_ticklabels(tick_labels) + ax.xaxis.set_major_locator(MaxNLocator(nbins=6, symmetric=True, min_n_ticks=7)) + ax.yaxis.set_major_locator(MaxNLocator(nbins=6, symmetric=True, min_n_ticks=7)) + formatter = FormatStrFormatter('%.3g') + ax.xaxis.set_major_formatter(formatter) + ax.yaxis.set_major_formatter(formatter) # Set title, limits, aspect ax.set_title("Diffraction Pattern") From 7a0c684b1dc0aa8e40ea8f3a21875302a63f8a4d Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Tue, 8 Jun 2021 17:57:32 +0000 Subject: [PATCH 07/29] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- freud/plot.py | 7 +++---- 1 file changed, 3 insertions(+), 4 deletions(-) diff --git a/freud/plot.py b/freud/plot.py index 091f64dd3..07eb861b8 100644 --- a/freud/plot.py +++ b/freud/plot.py @@ -11,7 +11,7 @@ try: import matplotlib.pyplot as plt from matplotlib.backends.backend_agg import FigureCanvasAgg - from matplotlib.ticker import MaxNLocator, FormatStrFormatter + from matplotlib.ticker import FormatStrFormatter, MaxNLocator except ImportError: raise ImportError("matplotlib must be installed for freud.plot.") @@ -523,8 +523,7 @@ def diffraction_plot( # Plot the diffraction image and color bar norm = matplotlib.colors.LogNorm(vmin=vmin, vmax=vmax) - extent = (np.min(k_values), np.max(k_values), - np.min(k_values), np.max(k_values)) + extent = (np.min(k_values), np.max(k_values), np.min(k_values), np.max(k_values)) im = ax.imshow( np.clip(diffraction, vmin, vmax), interpolation="nearest", @@ -548,7 +547,7 @@ def diffraction_plot( # Set tick locations and labels ax.xaxis.set_major_locator(MaxNLocator(nbins=6, symmetric=True, min_n_ticks=7)) ax.yaxis.set_major_locator(MaxNLocator(nbins=6, symmetric=True, min_n_ticks=7)) - formatter = FormatStrFormatter('%.3g') + formatter = FormatStrFormatter("%.3g") ax.xaxis.set_major_formatter(formatter) ax.yaxis.set_major_formatter(formatter) From e9f2d565031e653fd16e00e768592fadbeff301d Mon Sep 17 00:00:00 2001 From: Andrew Kerr Date: Wed, 9 Jun 2021 08:49:21 -0400 Subject: [PATCH 08/29] Updated k-value and k-vector tests to check for correct values at the ends of each array. Calculations used for comparison are based on np.fft.fftfreq. Any feedback on style and correctness is appreciated, in particular whether it is worth it to add tests for correct k-vectors at indices [0, -1] and [-1, 0] in addition to [0, 0] and [-1, -1]. --- tests/test_diffraction_DiffractionPattern.py | 58 +++++++++++++++----- 1 file changed, 44 insertions(+), 14 deletions(-) diff --git a/tests/test_diffraction_DiffractionPattern.py b/tests/test_diffraction_DiffractionPattern.py index 6b7f2b043..49a99b714 100644 --- a/tests/test_diffraction_DiffractionPattern.py +++ b/tests/test_diffraction_DiffractionPattern.py @@ -160,20 +160,50 @@ def test_k_values_and_k_vectors(self): dp = freud.diffraction.DiffractionPattern() for size in [2, 5, 10]: - # TODO: npoints doesn't actually matter, just use one point. - for npoints in [10, 20, 75]: - # TODO: This isn't using the correct signature for - # make_random_system. Flip npoints and size. - box, positions = freud.data.make_random_system(npoints, size) - dp.compute((box, positions)) - - output_size = dp.output_size - npt.assert_allclose(dp.k_values[output_size // 2], 0) - center_index = (output_size // 2, output_size // 2) - npt.assert_allclose(dp.k_vectors[center_index], [0, 0, 0]) - # Add tests for the left and right bins. Use the math for - # fftfreq and the _k_scale_factor to write an expression for - # the first and last bins' k-values and k-vectors. + box, positions = freud.data.make_random_system(size, 1) + zoom = 4 + view_orientation = np.asarray([1, 0, 0, 0]) + dp.compute((box, positions), view_orientation=view_orientation, zoom=zoom) + + output_size = dp.output_size + npt.assert_allclose(dp.k_values[output_size // 2], 0) + center_index = (output_size // 2, output_size // 2) + npt.assert_allclose(dp.k_vectors[center_index], [0, 0, 0]) + + # Tests for the left and right k_value bins. Uses the math for + # np.fft.n + #fftfreq and the _k_scale_factor to write an expression for + # the first and last bins' k-values and k-vectors. + + scale_factor = 2 * np.pi * output_size / (np.max(box.to_matrix()) * zoom) + + if output_size % 2 == 0: + + first_k_bin = -0.5 * scale_factor + last_k_bin = (0.5 - (1 / output_size)) * scale_factor + + npt.assert_allclose(dp.k_values[0], first_k_bin) + npt.assert_allclose(dp.k_values[-1], last_k_bin) + + first_k_vector = rowan.rotate(view_orientation, [first_k_bin, first_k_bin, 0]) + last_k_vector = rowan.rotate(view_orientation, [last_k_bin, last_k_bin, 0]) + + npt.assert_allclose(dp.k_vectors[0, 0], first_k_vector) + npt.assert_allclose(dp.k_vectors[-1, -1], last_k_vector) + + if output_size % 2 == 1: + + first_k_bin = (-0.5 + 1 / (2 * output_size)) * scale_factor + last_k_bin = (0.5 - 1 / (2 * output_size)) * scale_factor + + npt.assert_allclose(dp.k_values[0], first_k_bin) + npt.assert_allclose(dp.k_values[-1], last_k_bin) + + first_k_vector = rowan.rotate(view_orientation, [first_k_bin, first_k_bin, 0]) + last_k_vector = rowan.rotate(view_orientation, [last_k_bin, last_k_bin, 0]) + + npt.assert_allclose(dp.k_vectors[0, 0], first_k_vector) + npt.assert_allclose(dp.k_vectors[-1, -1], last_k_vector) def test_cubic_system(self): box, positions = freud.data.UnitCell.sc().generate_system(10) From d4a3ee5890f8e52dbcb3e18ca4f34146d8e53da4 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Wed, 9 Jun 2021 12:58:43 +0000 Subject: [PATCH 09/29] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- tests/test_diffraction_DiffractionPattern.py | 18 +++++++++++++----- 1 file changed, 13 insertions(+), 5 deletions(-) diff --git a/tests/test_diffraction_DiffractionPattern.py b/tests/test_diffraction_DiffractionPattern.py index 49a99b714..70655db59 100644 --- a/tests/test_diffraction_DiffractionPattern.py +++ b/tests/test_diffraction_DiffractionPattern.py @@ -172,7 +172,7 @@ def test_k_values_and_k_vectors(self): # Tests for the left and right k_value bins. Uses the math for # np.fft.n - #fftfreq and the _k_scale_factor to write an expression for + # fftfreq and the _k_scale_factor to write an expression for # the first and last bins' k-values and k-vectors. scale_factor = 2 * np.pi * output_size / (np.max(box.to_matrix()) * zoom) @@ -185,8 +185,12 @@ def test_k_values_and_k_vectors(self): npt.assert_allclose(dp.k_values[0], first_k_bin) npt.assert_allclose(dp.k_values[-1], last_k_bin) - first_k_vector = rowan.rotate(view_orientation, [first_k_bin, first_k_bin, 0]) - last_k_vector = rowan.rotate(view_orientation, [last_k_bin, last_k_bin, 0]) + first_k_vector = rowan.rotate( + view_orientation, [first_k_bin, first_k_bin, 0] + ) + last_k_vector = rowan.rotate( + view_orientation, [last_k_bin, last_k_bin, 0] + ) npt.assert_allclose(dp.k_vectors[0, 0], first_k_vector) npt.assert_allclose(dp.k_vectors[-1, -1], last_k_vector) @@ -199,8 +203,12 @@ def test_k_values_and_k_vectors(self): npt.assert_allclose(dp.k_values[0], first_k_bin) npt.assert_allclose(dp.k_values[-1], last_k_bin) - first_k_vector = rowan.rotate(view_orientation, [first_k_bin, first_k_bin, 0]) - last_k_vector = rowan.rotate(view_orientation, [last_k_bin, last_k_bin, 0]) + first_k_vector = rowan.rotate( + view_orientation, [first_k_bin, first_k_bin, 0] + ) + last_k_vector = rowan.rotate( + view_orientation, [last_k_bin, last_k_bin, 0] + ) npt.assert_allclose(dp.k_vectors[0, 0], first_k_vector) npt.assert_allclose(dp.k_vectors[-1, -1], last_k_vector) From 25031b14c71452b7a77fca0939c0a82b7db939a2 Mon Sep 17 00:00:00 2001 From: Andrew Kerr Date: Wed, 9 Jun 2021 22:22:08 -0400 Subject: [PATCH 10/29] Streamlining a few pieces of the k-value and k-vector test per bdice's suggestions. With quite a few tests that are identical except for slight changes to indices, I am wondering if there is someway to put this code in a loop. --- tests/test_diffraction_DiffractionPattern.py | 63 +++++++++++++++----- 1 file changed, 47 insertions(+), 16 deletions(-) diff --git a/tests/test_diffraction_DiffractionPattern.py b/tests/test_diffraction_DiffractionPattern.py index 49a99b714..05f5282a6 100644 --- a/tests/test_diffraction_DiffractionPattern.py +++ b/tests/test_diffraction_DiffractionPattern.py @@ -171,43 +171,74 @@ def test_k_values_and_k_vectors(self): npt.assert_allclose(dp.k_vectors[center_index], [0, 0, 0]) # Tests for the left and right k_value bins. Uses the math for - # np.fft.n - #fftfreq and the _k_scale_factor to write an expression for + # np.fft.fftfreq and the _k_scale_factor to write an expression for # the first and last bins' k-values and k-vectors. scale_factor = 2 * np.pi * output_size / (np.max(box.to_matrix()) * zoom) if output_size % 2 == 0: - first_k_bin = -0.5 * scale_factor - last_k_bin = (0.5 - (1 / output_size)) * scale_factor + first_k_value = -0.5 * scale_factor + last_k_value = (0.5 - (1 / output_size)) * scale_factor - npt.assert_allclose(dp.k_values[0], first_k_bin) - npt.assert_allclose(dp.k_values[-1], last_k_bin) + npt.assert_allclose(dp.k_values[0], first_k_value) + npt.assert_allclose(dp.k_values[-1], last_k_value) - first_k_vector = rowan.rotate(view_orientation, [first_k_bin, first_k_bin, 0]) - last_k_vector = rowan.rotate(view_orientation, [last_k_bin, last_k_bin, 0]) + first_k_vector = rowan.rotate(view_orientation, [first_k_value, first_k_value, 0]) + last_k_vector = rowan.rotate(view_orientation, [last_k_value, last_k_value, 0]) + top_right_k_vector = rowan.rotate(view_orientation, [first_k_value, last_k_value, 0]) + bottom_left_k_vector = rowan.rotate(view_orientation, [last_k_value, first_k_value, 0]) npt.assert_allclose(dp.k_vectors[0, 0], first_k_vector) npt.assert_allclose(dp.k_vectors[-1, -1], last_k_vector) + npt.assert_allclose(dp.k_vectors[0, -1], top_right_k_vector) + npt.assert_allclose(dp.k_vectors[-1, 0], bottom_left_k_vector) - if output_size % 2 == 1: + centre = output_size // 2 + top_centre_k_vector = rowan.rotate(view_orientation, [0, first_k_value, 0]) + bottom_centre_k_vector = rowan.rotate(view_orientation, [0, last_k_value, 0]) + left_centre_k_vector = rowan.rotate(view_orientation, [first_k_value, 0, 0]) + right_centre_k_vector = rowan.rotate(view_orientation, [last_k_value, 0, 0]) - first_k_bin = (-0.5 + 1 / (2 * output_size)) * scale_factor - last_k_bin = (0.5 - 1 / (2 * output_size)) * scale_factor + npt.assert_allclose(dp.k_vectors[centre, 0], top_centre_k_vector) + npt.assert_allclose(dp.k_vectors[centre, -1], bottom_centre_k_vector) + npt.assert_allclose(dp.k_vectors[0, centre], left_centre_k_vector) + npt.assert_allclose(dp.k_vectors[-1, centre], right_centre_k_vector) - npt.assert_allclose(dp.k_values[0], first_k_bin) - npt.assert_allclose(dp.k_values[-1], last_k_bin) + else: - first_k_vector = rowan.rotate(view_orientation, [first_k_bin, first_k_bin, 0]) - last_k_vector = rowan.rotate(view_orientation, [last_k_bin, last_k_bin, 0]) + first_k_value = (-0.5 + 1 / (2 * output_size)) * scale_factor + last_k_value = (0.5 - 1 / (2 * output_size)) * scale_factor + + npt.assert_allclose(dp.k_values[0], first_k_value) + npt.assert_allclose(dp.k_values[-1], last_k_value) + + first_k_vector = rowan.rotate(view_orientation, [first_k_value, first_k_value, 0]) + last_k_vector = rowan.rotate(view_orientation, [last_k_value, last_k_value, 0]) + top_right_k_vector = rowan.rotate(view_orientation, [first_k_value, last_k_value, 0]) + bottom_left_k_vector = rowan.rotate(view_orientation, [last_k_value, first_k_value, 0]) npt.assert_allclose(dp.k_vectors[0, 0], first_k_vector) npt.assert_allclose(dp.k_vectors[-1, -1], last_k_vector) + npt.assert_allclose(dp.k_vectors[0, -1], top_right_k_vector) + npt.assert_allclose(dp.k_vectors[-1, 0], bottom_left_k_vector) + + centre = output_size //2 + top_centre_k_vector = rowan.rotate(view_orientation, [0, first_k_value, 0]) + bottom_centre_k_vector = rowan.rotate(view_orientation, [0, last_k_value, 0]) + left_centre_k_vector = rowan.rotate(view_orientation, [first_k_value, 0, 0]) + right_centre_k_vector = rowan.rotate(view_orientation, [last_k_value, 0, 0]) + + npt.assert_allclose(dp.k_vectors[centre, 0], top_centre_k_vector) + npt.assert_allclose(dp.k_vectors[centre, -1], bottom_centre_k_vector) + npt.assert_allclose(dp.k_vectors[0, centre], left_centre_k_vector) + npt.assert_allclose(dp.k_vectors[-1, centre], right_centre_k_vector) def test_cubic_system(self): + pass + box, positions = freud.data.UnitCell.sc().generate_system(10) - dp = freud.diffraction.DiffractionPattern() + dp = freud.diffraction.DiffractionPattern(grid_size = 1024) dp.compute((box, positions)) # Make sure that the peaks are where we expect them. # Identify the indices of the highest values in dp.diffraction From 2e62ade92232c360cae03073eeaef4a84d4abb9c Mon Sep 17 00:00:00 2001 From: Andrew Kerr Date: Thu, 10 Jun 2021 07:19:20 -0400 Subject: [PATCH 11/29] Cleaning up a few pieces of code to match the current branch. --- tests/test_diffraction_DiffractionPattern.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tests/test_diffraction_DiffractionPattern.py b/tests/test_diffraction_DiffractionPattern.py index 05f5282a6..0bbb64cdb 100644 --- a/tests/test_diffraction_DiffractionPattern.py +++ b/tests/test_diffraction_DiffractionPattern.py @@ -238,7 +238,7 @@ def test_cubic_system(self): pass box, positions = freud.data.UnitCell.sc().generate_system(10) - dp = freud.diffraction.DiffractionPattern(grid_size = 1024) + dp = freud.diffraction.DiffractionPattern() dp.compute((box, positions)) # Make sure that the peaks are where we expect them. # Identify the indices of the highest values in dp.diffraction From c8b839df84cd04a33ed91be8eae4cbdacd1c0311 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Thu, 10 Jun 2021 11:37:41 +0000 Subject: [PATCH 12/29] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- tests/test_diffraction_DiffractionPattern.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tests/test_diffraction_DiffractionPattern.py b/tests/test_diffraction_DiffractionPattern.py index 0435dde19..c912a7ba4 100644 --- a/tests/test_diffraction_DiffractionPattern.py +++ b/tests/test_diffraction_DiffractionPattern.py @@ -247,7 +247,7 @@ def test_k_values_and_k_vectors(self): npt.assert_allclose(dp.k_vectors[0, -1], top_right_k_vector) npt.assert_allclose(dp.k_vectors[-1, 0], bottom_left_k_vector) - centre = output_size //2 + centre = output_size // 2 top_centre_k_vector = rowan.rotate( view_orientation, [0, first_k_value, 0] ) From adbfd3f7cccbde45fb38f58d327a3c2b04c2cec9 Mon Sep 17 00:00:00 2001 From: Andrew Kerr Date: Thu, 10 Jun 2021 08:45:57 -0400 Subject: [PATCH 13/29] Added a very basic test to verify diffraction peak locations in a simple cubic system. It is sensitive to certain parameters right now; the test registers peaks best when sigma_noise = 0.1*lattice_length and vq.kmeans is only passed indices where structure factor is in the greatest 50 values. --- tests/test_diffraction_DiffractionPattern.py | 24 +++++++++++++++----- 1 file changed, 18 insertions(+), 6 deletions(-) diff --git a/tests/test_diffraction_DiffractionPattern.py b/tests/test_diffraction_DiffractionPattern.py index 0435dde19..cdae269a1 100644 --- a/tests/test_diffraction_DiffractionPattern.py +++ b/tests/test_diffraction_DiffractionPattern.py @@ -3,6 +3,7 @@ import numpy.testing as npt import pytest import rowan +from scipy.cluster import vq import freud @@ -267,17 +268,28 @@ def test_k_values_and_k_vectors(self): npt.assert_allclose(dp.k_vectors[-1, centre], right_centre_k_vector) def test_cubic_system(self): - pass - - box, positions = freud.data.UnitCell.sc().generate_system(10) + length = 1 + box, positions = freud.data.UnitCell.sc().generate_system(num_replicas=16, scale=length, sigma_noise=0.1*length) dp = freud.diffraction.DiffractionPattern() dp.compute((box, positions)) + # Make sure that the peaks are where we expect them. # Identify the indices of the highest values in dp.diffraction - # and test that k * R == 2*pi*N for some integer N, all peak vectors k, + # and test that k * R == 2*pi*N for some integer N, for the corresponding peak k vectors # and lattice vectors R. This will be inexact because of binning. - # Perhaps try with a larger resolution like 1024. - # Use npt.assert_allclose. + threshold = np.partition(dp.diffraction, -50)[-50] + xs, ys = np.nonzero(dp.diffraction > threshold) + xy = np.dstack((xs, ys))[0] + + peaks, _ = vq.kmeans(xy.astype(float), 9) + peaks = peaks.astype(int) + + for peak in peaks: + peak_k_vector = dp.k_vectors[peak[0], peak[1]] + lattice_vector = [length, length, length] + dot_prod = np.dot(peak_k_vector, lattice_vector) + + npt.assert_allclose(dot_prod % (2 * np.pi), 0) def test_cubic_system_parameterized(self): pass From 1a1b5da98cdf6022d26898e617dd72bb9e485d39 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Thu, 10 Jun 2021 12:47:14 +0000 Subject: [PATCH 14/29] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- tests/test_diffraction_DiffractionPattern.py | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/tests/test_diffraction_DiffractionPattern.py b/tests/test_diffraction_DiffractionPattern.py index 2c45d0e3d..cd0f2e1c6 100644 --- a/tests/test_diffraction_DiffractionPattern.py +++ b/tests/test_diffraction_DiffractionPattern.py @@ -269,7 +269,9 @@ def test_k_values_and_k_vectors(self): def test_cubic_system(self): length = 1 - box, positions = freud.data.UnitCell.sc().generate_system(num_replicas=16, scale=length, sigma_noise=0.1*length) + box, positions = freud.data.UnitCell.sc().generate_system( + num_replicas=16, scale=length, sigma_noise=0.1 * length + ) dp = freud.diffraction.DiffractionPattern() dp.compute((box, positions)) From 5716ecdc8caebb5b6865e889593aaf49813fe842 Mon Sep 17 00:00:00 2001 From: Andrew Kerr <79798289+AKerr9509@users.noreply.github.com> Date: Thu, 10 Jun 2021 23:24:43 -0400 Subject: [PATCH 15/29] Update tests/test_diffraction_DiffractionPattern.py Co-authored-by: Bradley Dice --- tests/test_diffraction_DiffractionPattern.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tests/test_diffraction_DiffractionPattern.py b/tests/test_diffraction_DiffractionPattern.py index cd0f2e1c6..6d32816a9 100644 --- a/tests/test_diffraction_DiffractionPattern.py +++ b/tests/test_diffraction_DiffractionPattern.py @@ -185,7 +185,7 @@ def test_k_values_and_k_vectors(self): npt.assert_allclose(dp.k_values[0], first_k_value) npt.assert_allclose(dp.k_values[-1], last_k_value) - first_k_vector = rowan.rotate( + top_left_k_vector = rowan.rotate( view_orientation, [first_k_value, first_k_value, 0] ) last_k_vector = rowan.rotate( From d97da8aec411535d27771806ba3a427df834559b Mon Sep 17 00:00:00 2001 From: Andrew Kerr <79798289+AKerr9509@users.noreply.github.com> Date: Thu, 10 Jun 2021 23:24:51 -0400 Subject: [PATCH 16/29] Update tests/test_diffraction_DiffractionPattern.py Co-authored-by: Bradley Dice --- tests/test_diffraction_DiffractionPattern.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tests/test_diffraction_DiffractionPattern.py b/tests/test_diffraction_DiffractionPattern.py index 6d32816a9..5b9a5324e 100644 --- a/tests/test_diffraction_DiffractionPattern.py +++ b/tests/test_diffraction_DiffractionPattern.py @@ -188,7 +188,7 @@ def test_k_values_and_k_vectors(self): top_left_k_vector = rowan.rotate( view_orientation, [first_k_value, first_k_value, 0] ) - last_k_vector = rowan.rotate( + bottom_right_k_vector = rowan.rotate( view_orientation, [last_k_value, last_k_value, 0] ) top_right_k_vector = rowan.rotate( From c9d9c9715ff7045ad6520709e052d46c3bab3f95 Mon Sep 17 00:00:00 2001 From: Andrew Kerr <79798289+AKerr9509@users.noreply.github.com> Date: Thu, 10 Jun 2021 23:24:59 -0400 Subject: [PATCH 17/29] Update tests/test_diffraction_DiffractionPattern.py Co-authored-by: Bradley Dice --- tests/test_diffraction_DiffractionPattern.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tests/test_diffraction_DiffractionPattern.py b/tests/test_diffraction_DiffractionPattern.py index 5b9a5324e..7c447ca00 100644 --- a/tests/test_diffraction_DiffractionPattern.py +++ b/tests/test_diffraction_DiffractionPattern.py @@ -203,7 +203,7 @@ def test_k_values_and_k_vectors(self): npt.assert_allclose(dp.k_vectors[0, -1], top_right_k_vector) npt.assert_allclose(dp.k_vectors[-1, 0], bottom_left_k_vector) - centre = output_size // 2 + center = output_size // 2 top_centre_k_vector = rowan.rotate( view_orientation, [0, first_k_value, 0] ) From 12a55ccbec4d68f2cc9788aaf95a925662ad4541 Mon Sep 17 00:00:00 2001 From: Andrew Kerr Date: Fri, 11 Jun 2021 09:10:01 -0400 Subject: [PATCH 18/29] New test for correct diffraction peak locations in a simple cubic system with various grid size, output size, and zoom values. It isolates the brightest areas of the diffraction pattern (S(q) greater than a certain threshold), and verifies that the ideal peak locations given by k * R = 2 * pi * N are contained in these areas. I opted for this approach over locating peaks with scipy's kmeans algorithm because the random noise in the system and slight inexactness due to binning prevented kmeans from reliably locating the exact pixels which contained diffraction peaks. Parts of this test feel a bit 'hacky' right now (using a dictionary, iterating through it twice), and I am open to any/all style and efficiency suggestions. --- tests/test_diffraction_DiffractionPattern.py | 94 +++++++++++++------- 1 file changed, 62 insertions(+), 32 deletions(-) diff --git a/tests/test_diffraction_DiffractionPattern.py b/tests/test_diffraction_DiffractionPattern.py index 2c45d0e3d..2ea9a6f32 100644 --- a/tests/test_diffraction_DiffractionPattern.py +++ b/tests/test_diffraction_DiffractionPattern.py @@ -203,24 +203,24 @@ def test_k_values_and_k_vectors(self): npt.assert_allclose(dp.k_vectors[0, -1], top_right_k_vector) npt.assert_allclose(dp.k_vectors[-1, 0], bottom_left_k_vector) - centre = output_size // 2 - top_centre_k_vector = rowan.rotate( + center = output_size // 2 + top_center_k_vector = rowan.rotate( view_orientation, [0, first_k_value, 0] ) - bottom_centre_k_vector = rowan.rotate( + bottom_center_k_vector = rowan.rotate( view_orientation, [0, last_k_value, 0] ) - left_centre_k_vector = rowan.rotate( + left_center_k_vector = rowan.rotate( view_orientation, [first_k_value, 0, 0] ) - right_centre_k_vector = rowan.rotate( + right_center_k_vector = rowan.rotate( view_orientation, [last_k_value, 0, 0] ) - npt.assert_allclose(dp.k_vectors[centre, 0], top_centre_k_vector) - npt.assert_allclose(dp.k_vectors[centre, -1], bottom_centre_k_vector) - npt.assert_allclose(dp.k_vectors[0, centre], left_centre_k_vector) - npt.assert_allclose(dp.k_vectors[-1, centre], right_centre_k_vector) + npt.assert_allclose(dp.k_vectors[center, 0], top_center_k_vector) + npt.assert_allclose(dp.k_vectors[center, -1], bottom_center_k_vector) + npt.assert_allclose(dp.k_vectors[0, center], left_center_k_vector) + npt.assert_allclose(dp.k_vectors[-1, center], right_center_k_vector) else: @@ -248,58 +248,88 @@ def test_k_values_and_k_vectors(self): npt.assert_allclose(dp.k_vectors[0, -1], top_right_k_vector) npt.assert_allclose(dp.k_vectors[-1, 0], bottom_left_k_vector) - centre = output_size // 2 - top_centre_k_vector = rowan.rotate( + center = output_size // 2 + top_center_k_vector = rowan.rotate( view_orientation, [0, first_k_value, 0] ) - bottom_centre_k_vector = rowan.rotate( + bottom_center_k_vector = rowan.rotate( view_orientation, [0, last_k_value, 0] ) - left_centre_k_vector = rowan.rotate( + left_center_k_vector = rowan.rotate( view_orientation, [first_k_value, 0, 0] ) - right_centre_k_vector = rowan.rotate( + right_center_k_vector = rowan.rotate( view_orientation, [last_k_value, 0, 0] ) - npt.assert_allclose(dp.k_vectors[centre, 0], top_centre_k_vector) - npt.assert_allclose(dp.k_vectors[centre, -1], bottom_centre_k_vector) - npt.assert_allclose(dp.k_vectors[0, centre], left_centre_k_vector) - npt.assert_allclose(dp.k_vectors[-1, centre], right_centre_k_vector) + npt.assert_allclose(dp.k_vectors[center, 0], top_center_k_vector) + npt.assert_allclose(dp.k_vectors[center, -1], bottom_center_k_vector) + npt.assert_allclose(dp.k_vectors[0, center], left_center_k_vector) + npt.assert_allclose(dp.k_vectors[-1, center], right_center_k_vector) def test_cubic_system(self): length = 1 - box, positions = freud.data.UnitCell.sc().generate_system(num_replicas=16, scale=length, sigma_noise=0.1*length) + box, positions = freud.data.UnitCell.sc().generate_system( + num_replicas=16, scale=length, sigma_noise=0.1*length + ) dp = freud.diffraction.DiffractionPattern() - dp.compute((box, positions)) + dp.compute((box, positions), zoom=2) # Make sure that the peaks are where we expect them. # Identify the indices of the highest values in dp.diffraction # and test that k * R == 2*pi*N for some integer N, for the corresponding peak k vectors # and lattice vectors R. This will be inexact because of binning. - threshold = np.partition(dp.diffraction, -50)[-50] + threshold = np.partition(dp.diffraction.flatten(), -50)[-50] xs, ys = np.nonzero(dp.diffraction > threshold) xy = np.dstack((xs, ys))[0] - peaks, _ = vq.kmeans(xy.astype(float), 9) - peaks = peaks.astype(int) + centers, _ = vq.kmeans(xy.astype(float), 9) + centers = centers.astype(int) - for peak in peaks: - peak_k_vector = dp.k_vectors[peak[0], peak[1]] + for center in centers: + peak_k_vector = dp.k_vectors[center[0], center[1]] lattice_vector = [length, length, length] dot_prod = np.dot(peak_k_vector, lattice_vector) npt.assert_allclose(dot_prod % (2 * np.pi), 0) def test_cubic_system_parameterized(self): - pass - # Second PR: Same as above test but with different grid_size, + length = 1 + box, positions = freud.data.UnitCell.sc().generate_system( + num_replicas=16, scale=length, sigma_noise=0.1*length + ) + #Same as above test but with different grid_size, # output_size, and zoom values. for grid_size in (256, 1024): for output_size in (255, 256, 1023, 1024): for zoom in (1, 2.5, 4): - # Ensure that peaks are in the correct locations like above. - # Identify the indices of the highest values in - # dp.diffraction and test that k * R == 2*pi*N for some - # integer N, all peak vectors k, and lattice vectors R. - pass + dp = freud.diffraction.DiffractionPattern( + grid_size=grid_size, output_size=output_size, + ) + dp.compute((box, positions), zoom=zoom) + + # Locate brightest areas of diffraction pattern (intensity > threshold), + # and check that the ideal diffraction peak locations, given by + # k * R = 2*pi*N for some lattice vector R and integer N, are contained + # within these regions. This test only checks N in range [-2, 2]. + threshold = 0.2 + xs, ys = np.nonzero(dp.diffraction > threshold) + xy = np.dstack((xs, ys))[0] + + ideal_peaks = {-2:'f', -1:'f', 0:'f', 1:'f', 2:'f'} + all_peaks = True + + for peak in ideal_peaks: + for x, y in xy: + k_vector = dp.k_vectors[x, y] + lattice_vector = [1, 1, 1] + dot_prod = np.dot(k_vector, lattice_vector) + + if dot_prod == (peak * 2 * np.pi): + ideal_peaks[peak] = 't' + + for peak in ideal_peaks: + if ideal_peaks[peak] == 'f': + all_peaks = False + + assert all_peaks == True From 7a1e4f3633b4480855bf47adbf9b31d6ab6bf887 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Fri, 11 Jun 2021 13:12:35 +0000 Subject: [PATCH 19/29] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- tests/test_diffraction_DiffractionPattern.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/tests/test_diffraction_DiffractionPattern.py b/tests/test_diffraction_DiffractionPattern.py index 68d47a037..98e64c88a 100644 --- a/tests/test_diffraction_DiffractionPattern.py +++ b/tests/test_diffraction_DiffractionPattern.py @@ -332,10 +332,10 @@ def test_cubic_system_parameterized(self): k_vector = dp.k_vectors[x, y] lattice_vector = [1, 1, 1] dot_prod = np.dot(k_vector, lattice_vector) - + if dot_prod == (peak * 2 * np.pi): ideal_peaks[peak] = 't' - + for peak in ideal_peaks: if ideal_peaks[peak] == 'f': all_peaks = False From 33c85273f95d67bb314fda6e18c376f674abff74 Mon Sep 17 00:00:00 2001 From: Andrew Kerr Date: Fri, 11 Jun 2021 12:16:22 -0400 Subject: [PATCH 20/29] Updated k-vector description in diffraction.pyx docstring. --- freud/diffraction.pyx | 7 ++++++- tests/test_diffraction_DiffractionPattern.py | 8 -------- 2 files changed, 6 insertions(+), 9 deletions(-) diff --git a/freud/diffraction.pyx b/freud/diffraction.pyx index 90306cd8e..97278186d 100644 --- a/freud/diffraction.pyx +++ b/freud/diffraction.pyx @@ -48,7 +48,12 @@ cdef class DiffractionPattern(_Compute): The :math:`\vec{k}=0` peak is always located at index ``(output_size // 2, output_size // 2)`` and is normalized to have a value - of :math:`S(\vec{k}=0) = 1` (not :math:`N`, a common convention). + of :math:`S(\vec{k}=0) = 1` (not :math:`N`, a common convention). The + remaining :math:`\vec{k}`-vectors are computed such that each peak in the + diffraction pattern satisfies the relationship :math:`\vec{k} \cdot + \vec{R} = 2\piN` for some integer :math: `N` and lattice vector of + the system :math: `R` (see + `this article `). This method is based on the implementations in the open-source `GIXStapose application `_ and its diff --git a/tests/test_diffraction_DiffractionPattern.py b/tests/test_diffraction_DiffractionPattern.py index 68d47a037..1f843bef5 100644 --- a/tests/test_diffraction_DiffractionPattern.py +++ b/tests/test_diffraction_DiffractionPattern.py @@ -204,11 +204,7 @@ def test_k_values_and_k_vectors(self): npt.assert_allclose(dp.k_vectors[-1, 0], bottom_left_k_vector) center = output_size // 2 -<<<<<<< HEAD top_center_k_vector = rowan.rotate( -======= - top_centre_k_vector = rowan.rotate( ->>>>>>> c9d9c9715ff7045ad6520709e052d46c3bab3f95 view_orientation, [0, first_k_value, 0] ) bottom_center_k_vector = rowan.rotate( @@ -274,11 +270,7 @@ def test_k_values_and_k_vectors(self): def test_cubic_system(self): length = 1 box, positions = freud.data.UnitCell.sc().generate_system( -<<<<<<< HEAD - num_replicas=16, scale=length, sigma_noise=0.1*length -======= num_replicas=16, scale=length, sigma_noise=0.1 * length ->>>>>>> c9d9c9715ff7045ad6520709e052d46c3bab3f95 ) dp = freud.diffraction.DiffractionPattern() dp.compute((box, positions), zoom=2) From fd61d861dea7133c6ebf15e3c220912c5267a744 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Fri, 11 Jun 2021 16:17:06 +0000 Subject: [PATCH 21/29] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- tests/test_diffraction_DiffractionPattern.py | 15 ++++++++------- 1 file changed, 8 insertions(+), 7 deletions(-) diff --git a/tests/test_diffraction_DiffractionPattern.py b/tests/test_diffraction_DiffractionPattern.py index dd7e4f14c..b483a919b 100644 --- a/tests/test_diffraction_DiffractionPattern.py +++ b/tests/test_diffraction_DiffractionPattern.py @@ -296,16 +296,17 @@ def test_cubic_system(self): def test_cubic_system_parameterized(self): length = 1 box, positions = freud.data.UnitCell.sc().generate_system( - num_replicas=16, scale=length, sigma_noise=0.1*length + num_replicas=16, scale=length, sigma_noise=0.1 * length ) - #Same as above test but with different grid_size, + # Same as above test but with different grid_size, # output_size, and zoom values. for grid_size in (256, 1024): for output_size in (255, 256, 1023, 1024): for zoom in (1, 2.5, 4): dp = freud.diffraction.DiffractionPattern( - grid_size=grid_size, output_size=output_size, - ) + grid_size=grid_size, + output_size=output_size, + ) dp.compute((box, positions), zoom=zoom) # Locate brightest areas of diffraction pattern (intensity > threshold), @@ -316,7 +317,7 @@ def test_cubic_system_parameterized(self): xs, ys = np.nonzero(dp.diffraction > threshold) xy = np.dstack((xs, ys))[0] - ideal_peaks = {-2:'f', -1:'f', 0:'f', 1:'f', 2:'f'} + ideal_peaks = {-2: "f", -1: "f", 0: "f", 1: "f", 2: "f"} all_peaks = True for peak in ideal_peaks: @@ -326,10 +327,10 @@ def test_cubic_system_parameterized(self): dot_prod = np.dot(k_vector, lattice_vector) if dot_prod == (peak * 2 * np.pi): - ideal_peaks[peak] = 't' + ideal_peaks[peak] = "t" for peak in ideal_peaks: - if ideal_peaks[peak] == 'f': + if ideal_peaks[peak] == "f": all_peaks = False assert all_peaks == True From 46328394d5f902d7564399b222cdef6c8ff94fb4 Mon Sep 17 00:00:00 2001 From: Andrew Kerr Date: Mon, 14 Jun 2021 07:58:03 -0400 Subject: [PATCH 22/29] Streamlining and style fixes --- freud/diffraction.pyx | 10 +- tests/test_diffraction_DiffractionPattern.py | 130 ++++++------------- 2 files changed, 47 insertions(+), 93 deletions(-) diff --git a/freud/diffraction.pyx b/freud/diffraction.pyx index 97278186d..2ed0c8246 100644 --- a/freud/diffraction.pyx +++ b/freud/diffraction.pyx @@ -49,11 +49,11 @@ cdef class DiffractionPattern(_Compute): The :math:`\vec{k}=0` peak is always located at index ``(output_size // 2, output_size // 2)`` and is normalized to have a value of :math:`S(\vec{k}=0) = 1` (not :math:`N`, a common convention). The - remaining :math:`\vec{k}`-vectors are computed such that each peak in the + remaining :math:`\vec{k}` vectors are computed such that each peak in the diffraction pattern satisfies the relationship :math:`\vec{k} \cdot - \vec{R} = 2\piN` for some integer :math: `N` and lattice vector of - the system :math: `R` (see - `this article `). + \vec{R} = 2 \pi N` for some integer :math:`N` and lattice vector of + the system :math:`R`. See the `reciprocal lattice Wikipedia page + `__ for more information. This method is based on the implementations in the open-source `GIXStapose application `_ and its @@ -65,7 +65,7 @@ cdef class DiffractionPattern(_Compute): output_size (unsigned int): Resolution of the output diffraction image, uses ``grid_size`` if not provided or ``None`` (Default value = :code:`None`). - """ + """ # noqa: E501 cdef int _grid_size cdef int _output_size cdef double[:] _k_values_orig diff --git a/tests/test_diffraction_DiffractionPattern.py b/tests/test_diffraction_DiffractionPattern.py index dd7e4f14c..0af0e0d80 100644 --- a/tests/test_diffraction_DiffractionPattern.py +++ b/tests/test_diffraction_DiffractionPattern.py @@ -178,94 +178,52 @@ def test_k_values_and_k_vectors(self): scale_factor = 2 * np.pi * output_size / (np.max(box.to_matrix()) * zoom) if output_size % 2 == 0: - first_k_value = -0.5 * scale_factor last_k_value = (0.5 - (1 / output_size)) * scale_factor - npt.assert_allclose(dp.k_values[0], first_k_value) - npt.assert_allclose(dp.k_values[-1], last_k_value) - - top_left_k_vector = rowan.rotate( - view_orientation, [first_k_value, first_k_value, 0] - ) - bottom_right_k_vector = rowan.rotate( - view_orientation, [last_k_value, last_k_value, 0] - ) - top_right_k_vector = rowan.rotate( - view_orientation, [first_k_value, last_k_value, 0] - ) - bottom_left_k_vector = rowan.rotate( - view_orientation, [last_k_value, first_k_value, 0] - ) - - npt.assert_allclose(dp.k_vectors[0, 0], first_k_vector) - npt.assert_allclose(dp.k_vectors[-1, -1], last_k_vector) - npt.assert_allclose(dp.k_vectors[0, -1], top_right_k_vector) - npt.assert_allclose(dp.k_vectors[-1, 0], bottom_left_k_vector) - - center = output_size // 2 - top_center_k_vector = rowan.rotate( - view_orientation, [0, first_k_value, 0] - ) - bottom_center_k_vector = rowan.rotate( - view_orientation, [0, last_k_value, 0] - ) - left_center_k_vector = rowan.rotate( - view_orientation, [first_k_value, 0, 0] - ) - right_center_k_vector = rowan.rotate( - view_orientation, [last_k_value, 0, 0] - ) - - npt.assert_allclose(dp.k_vectors[center, 0], top_center_k_vector) - npt.assert_allclose(dp.k_vectors[center, -1], bottom_center_k_vector) - npt.assert_allclose(dp.k_vectors[0, center], left_center_k_vector) - npt.assert_allclose(dp.k_vectors[-1, center], right_center_k_vector) - else: - first_k_value = (-0.5 + 1 / (2 * output_size)) * scale_factor last_k_value = (0.5 - 1 / (2 * output_size)) * scale_factor - npt.assert_allclose(dp.k_values[0], first_k_value) - npt.assert_allclose(dp.k_values[-1], last_k_value) - - first_k_vector = rowan.rotate( - view_orientation, [first_k_value, first_k_value, 0] - ) - last_k_vector = rowan.rotate( - view_orientation, [last_k_value, last_k_value, 0] - ) - top_right_k_vector = rowan.rotate( - view_orientation, [first_k_value, last_k_value, 0] - ) - bottom_left_k_vector = rowan.rotate( - view_orientation, [last_k_value, first_k_value, 0] - ) - - npt.assert_allclose(dp.k_vectors[0, 0], first_k_vector) - npt.assert_allclose(dp.k_vectors[-1, -1], last_k_vector) - npt.assert_allclose(dp.k_vectors[0, -1], top_right_k_vector) - npt.assert_allclose(dp.k_vectors[-1, 0], bottom_left_k_vector) - - center = output_size // 2 - top_center_k_vector = rowan.rotate( - view_orientation, [0, first_k_value, 0] - ) - bottom_center_k_vector = rowan.rotate( - view_orientation, [0, last_k_value, 0] - ) - left_center_k_vector = rowan.rotate( - view_orientation, [first_k_value, 0, 0] - ) - right_center_k_vector = rowan.rotate( - view_orientation, [last_k_value, 0, 0] - ) - - npt.assert_allclose(dp.k_vectors[center, 0], top_center_k_vector) - npt.assert_allclose(dp.k_vectors[center, -1], bottom_center_k_vector) - npt.assert_allclose(dp.k_vectors[0, center], left_center_k_vector) - npt.assert_allclose(dp.k_vectors[-1, center], right_center_k_vector) + npt.assert_allclose(dp.k_values[0], first_k_value) + npt.assert_allclose(dp.k_values[-1], last_k_value) + + first_k_vector = rowan.rotate( + view_orientation, [first_k_value, first_k_value, 0] + ) + last_k_vector = rowan.rotate( + view_orientation, [last_k_value, last_k_value, 0] + ) + top_right_k_vector = rowan.rotate( + view_orientation, [first_k_value, last_k_value, 0] + ) + bottom_left_k_vector = rowan.rotate( + view_orientation, [last_k_value, first_k_value, 0] + ) + + npt.assert_allclose(dp.k_vectors[0, 0], first_k_vector) + npt.assert_allclose(dp.k_vectors[-1, -1], last_k_vector) + npt.assert_allclose(dp.k_vectors[0, -1], top_right_k_vector) + npt.assert_allclose(dp.k_vectors[-1, 0], bottom_left_k_vector) + + center = output_size // 2 + top_center_k_vector = rowan.rotate( + view_orientation, [0, first_k_value, 0] + ) + bottom_center_k_vector = rowan.rotate( + view_orientation, [0, last_k_value, 0] + ) + left_center_k_vector = rowan.rotate( + view_orientation, [first_k_value, 0, 0] + ) + right_center_k_vector = rowan.rotate( + view_orientation, [last_k_value, 0, 0] + ) + + npt.assert_allclose(dp.k_vectors[center, 0], top_center_k_vector) + npt.assert_allclose(dp.k_vectors[center, -1], bottom_center_k_vector) + npt.assert_allclose(dp.k_vectors[0, center], left_center_k_vector) + npt.assert_allclose(dp.k_vectors[-1, center], right_center_k_vector) def test_cubic_system(self): length = 1 @@ -316,7 +274,7 @@ def test_cubic_system_parameterized(self): xs, ys = np.nonzero(dp.diffraction > threshold) xy = np.dstack((xs, ys))[0] - ideal_peaks = {-2:'f', -1:'f', 0:'f', 1:'f', 2:'f'} + ideal_peaks = {i: False for i in [-2, -1, 0, 1, 2]} all_peaks = True for peak in ideal_peaks: @@ -326,10 +284,6 @@ def test_cubic_system_parameterized(self): dot_prod = np.dot(k_vector, lattice_vector) if dot_prod == (peak * 2 * np.pi): - ideal_peaks[peak] = 't' - - for peak in ideal_peaks: - if ideal_peaks[peak] == 'f': - all_peaks = False + ideal_peaks[peak] = True - assert all_peaks == True + assert all(ideal_peaks.values()) From 2cd6cd391c6a1c1421a711e5016d71e19c41c081 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Mon, 14 Jun 2021 12:01:35 +0000 Subject: [PATCH 23/29] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- tests/test_diffraction_DiffractionPattern.py | 12 +++--------- 1 file changed, 3 insertions(+), 9 deletions(-) diff --git a/tests/test_diffraction_DiffractionPattern.py b/tests/test_diffraction_DiffractionPattern.py index 86b3835c2..d90caa80c 100644 --- a/tests/test_diffraction_DiffractionPattern.py +++ b/tests/test_diffraction_DiffractionPattern.py @@ -207,18 +207,12 @@ def test_k_values_and_k_vectors(self): npt.assert_allclose(dp.k_vectors[-1, 0], bottom_left_k_vector) center = output_size // 2 - top_center_k_vector = rowan.rotate( - view_orientation, [0, first_k_value, 0] - ) + top_center_k_vector = rowan.rotate(view_orientation, [0, first_k_value, 0]) bottom_center_k_vector = rowan.rotate( view_orientation, [0, last_k_value, 0] ) - left_center_k_vector = rowan.rotate( - view_orientation, [first_k_value, 0, 0] - ) - right_center_k_vector = rowan.rotate( - view_orientation, [last_k_value, 0, 0] - ) + left_center_k_vector = rowan.rotate(view_orientation, [first_k_value, 0, 0]) + right_center_k_vector = rowan.rotate(view_orientation, [last_k_value, 0, 0]) npt.assert_allclose(dp.k_vectors[center, 0], top_center_k_vector) npt.assert_allclose(dp.k_vectors[center, -1], bottom_center_k_vector) From 64f3025bc53e644451c394e6e4c1d29b89c7af97 Mon Sep 17 00:00:00 2001 From: Andrew Kerr Date: Mon, 14 Jun 2021 08:28:49 -0400 Subject: [PATCH 24/29] Cleaning up comments / unused variables --- tests/test_diffraction_DiffractionPattern.py | 35 ++++++++++---------- 1 file changed, 17 insertions(+), 18 deletions(-) diff --git a/tests/test_diffraction_DiffractionPattern.py b/tests/test_diffraction_DiffractionPattern.py index 86b3835c2..e2ee1b058 100644 --- a/tests/test_diffraction_DiffractionPattern.py +++ b/tests/test_diffraction_DiffractionPattern.py @@ -3,7 +3,6 @@ import numpy.testing as npt import pytest import rowan -from scipy.cluster import vq import freud @@ -233,23 +232,28 @@ def test_cubic_system(self): dp = freud.diffraction.DiffractionPattern() dp.compute((box, positions), zoom=2) - # Make sure that the peaks are where we expect them. - # Identify the indices of the highest values in dp.diffraction - # and test that k * R == 2*pi*N for some integer N, for the corresponding peak k vectors - # and lattice vectors R. This will be inexact because of binning. - threshold = np.partition(dp.diffraction.flatten(), -50)[-50] + # Locate brightest areas of diffraction pattern + # (intensity > threshold), and check that the ideal + # diffraction peak locations, given by k * R = 2*pi*N + # for some lattice vector R and integer N, are contained + # within these regions. + # This test only checks N in range [-2, 2]. + threshold = 0.2 xs, ys = np.nonzero(dp.diffraction > threshold) xy = np.dstack((xs, ys))[0] - centers, _ = vq.kmeans(xy.astype(float), 9) - centers = centers.astype(int) + ideal_peaks = {i: False for i in [-2, -1, 0, 1, 2]} - for center in centers: - peak_k_vector = dp.k_vectors[center[0], center[1]] - lattice_vector = [length, length, length] - dot_prod = np.dot(peak_k_vector, lattice_vector) + for peak in ideal_peaks: + for x, y in xy: + k_vector = dp.k_vectors[x, y] + lattice_vector = [1, 1, 1] + dot_prod = np.dot(k_vector, lattice_vector) - npt.assert_allclose(dot_prod % (2 * np.pi), 0) + if dot_prod == (peak * 2 * np.pi): + ideal_peaks[peak] = True + + assert all(ideal_peaks.values()) def test_cubic_system_parameterized(self): length = 1 @@ -267,16 +271,11 @@ def test_cubic_system_parameterized(self): ) dp.compute((box, positions), zoom=zoom) - # Locate brightest areas of diffraction pattern (intensity > threshold), - # and check that the ideal diffraction peak locations, given by - # k * R = 2*pi*N for some lattice vector R and integer N, are contained - # within these regions. This test only checks N in range [-2, 2]. threshold = 0.2 xs, ys = np.nonzero(dp.diffraction > threshold) xy = np.dstack((xs, ys))[0] ideal_peaks = {i: False for i in [-2, -1, 0, 1, 2]} - all_peaks = True for peak in ideal_peaks: for x, y in xy: From 87b2947735eb69c1e8d662cd60f39073147e24b5 Mon Sep 17 00:00:00 2001 From: Bradley Dice Date: Mon, 14 Jun 2021 10:52:29 -0500 Subject: [PATCH 25/29] Remove unused code. --- freud/plot.py | 8 -------- 1 file changed, 8 deletions(-) diff --git a/freud/plot.py b/freud/plot.py index 3ff240a52..ef81cfdb0 100644 --- a/freud/plot.py +++ b/freud/plot.py @@ -536,14 +536,6 @@ def diffraction_plot( cb = Colorbar(cax, im) cb.set_label(r"$S(\vec{k})$") - # Determine the number of ticks on the axis - grid_size = diffraction.shape[0] - num_ticks = len([i for i in ax.xaxis.get_ticklocs() if 0 <= i <= grid_size]) - - # Ensure there are an odd number of ticks, so that there is a tick at zero - num_ticks += 1 - num_ticks % 2 - ticks = np.linspace(0, grid_size, num_ticks) - # Set tick locations and labels ax.xaxis.set_major_locator(MaxNLocator(nbins=6, symmetric=True, min_n_ticks=7)) ax.yaxis.set_major_locator(MaxNLocator(nbins=6, symmetric=True, min_n_ticks=7)) From d5fc889396a578cbb4a67d5b608962a360b167e7 Mon Sep 17 00:00:00 2001 From: Bradley Dice Date: Mon, 14 Jun 2021 11:55:37 -0500 Subject: [PATCH 26/29] R vector. --- freud/diffraction.pyx | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/freud/diffraction.pyx b/freud/diffraction.pyx index 2ed0c8246..0f0bd3f17 100644 --- a/freud/diffraction.pyx +++ b/freud/diffraction.pyx @@ -52,7 +52,7 @@ cdef class DiffractionPattern(_Compute): remaining :math:`\vec{k}` vectors are computed such that each peak in the diffraction pattern satisfies the relationship :math:`\vec{k} \cdot \vec{R} = 2 \pi N` for some integer :math:`N` and lattice vector of - the system :math:`R`. See the `reciprocal lattice Wikipedia page + the system :math:`\vec{R}`. See the `reciprocal lattice Wikipedia page `__ for more information. This method is based on the implementations in the open-source From 6871ce7671e32bfe5d79b3873b944d4d0e818568 Mon Sep 17 00:00:00 2001 From: Bradley Dice Date: Mon, 14 Jun 2021 11:56:01 -0500 Subject: [PATCH 27/29] Remove noqa, it's not needed here because the line wrapped within the length limit. --- freud/diffraction.pyx | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/freud/diffraction.pyx b/freud/diffraction.pyx index 0f0bd3f17..d48f354df 100644 --- a/freud/diffraction.pyx +++ b/freud/diffraction.pyx @@ -65,7 +65,7 @@ cdef class DiffractionPattern(_Compute): output_size (unsigned int): Resolution of the output diffraction image, uses ``grid_size`` if not provided or ``None`` (Default value = :code:`None`). - """ # noqa: E501 + """ cdef int _grid_size cdef int _output_size cdef double[:] _k_values_orig From 7b92b126dddecb20ccc5db494b8467980c09f3b2 Mon Sep 17 00:00:00 2001 From: Bradley Dice Date: Mon, 14 Jun 2021 11:56:45 -0500 Subject: [PATCH 28/29] Use more rigorous (imperfect) test case for peaks aligning with k*R. --- tests/test_diffraction_DiffractionPattern.py | 15 ++++++++------- 1 file changed, 8 insertions(+), 7 deletions(-) diff --git a/tests/test_diffraction_DiffractionPattern.py b/tests/test_diffraction_DiffractionPattern.py index b25f34dfa..b442a56a7 100644 --- a/tests/test_diffraction_DiffractionPattern.py +++ b/tests/test_diffraction_DiffractionPattern.py @@ -170,16 +170,15 @@ def test_k_values_and_k_vectors(self): center_index = (output_size // 2, output_size // 2) npt.assert_allclose(dp.k_vectors[center_index], [0, 0, 0]) - # Tests for the left and right k_value bins. Uses the math for + # Tests for the first and last k value bins. Uses the math for # np.fft.fftfreq and the _k_scale_factor to write an expression for - # the first and last bins' k-values and k-vectors. + # the first and last bins' k values and k vectors. scale_factor = 2 * np.pi * output_size / (np.max(box.to_matrix()) * zoom) if output_size % 2 == 0: first_k_value = -0.5 * scale_factor last_k_value = (0.5 - (1 / output_size)) * scale_factor - else: first_k_value = (-0.5 + 1 / (2 * output_size)) * scale_factor last_k_value = (0.5 - 1 / (2 * output_size)) * scale_factor @@ -223,8 +222,10 @@ def test_cubic_system(self): box, positions = freud.data.UnitCell.sc().generate_system( num_replicas=16, scale=length, sigma_noise=0.1 * length ) - dp = freud.diffraction.DiffractionPattern() - dp.compute((box, positions), zoom=2) + # Pick a non-integer value for zoom, to ensure that peaks besides k=0 + # are not perfectly aligned on pixels. + dp = freud.diffraction.DiffractionPattern(grid_size=512) + dp.compute((box, positions), zoom=4.123) # Locate brightest areas of diffraction pattern # (intensity > threshold), and check that the ideal @@ -238,13 +239,13 @@ def test_cubic_system(self): ideal_peaks = {i: False for i in [-2, -1, 0, 1, 2]} + lattice_vector = np.array([length, length, length]) for peak in ideal_peaks: for x, y in xy: k_vector = dp.k_vectors[x, y] - lattice_vector = [1, 1, 1] dot_prod = np.dot(k_vector, lattice_vector) - if dot_prod == (peak * 2 * np.pi): + if np.isclose(dot_prod, peak * 2 * np.pi, atol=1e-2): ideal_peaks[peak] = True assert all(ideal_peaks.values()) From aa80c640276f3d7b38354783e0d21480cf49e3d3 Mon Sep 17 00:00:00 2001 From: Bradley Dice Date: Mon, 14 Jun 2021 11:58:50 -0500 Subject: [PATCH 29/29] Make other DiffractionPattern test use approximate comparisons. --- tests/test_diffraction_DiffractionPattern.py | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/tests/test_diffraction_DiffractionPattern.py b/tests/test_diffraction_DiffractionPattern.py index b442a56a7..28720cd0a 100644 --- a/tests/test_diffraction_DiffractionPattern.py +++ b/tests/test_diffraction_DiffractionPattern.py @@ -259,7 +259,7 @@ def test_cubic_system_parameterized(self): # output_size, and zoom values. for grid_size in (256, 1024): for output_size in (255, 256, 1023, 1024): - for zoom in (1, 2.5, 4): + for zoom in (1, 2.5, 4.123): dp = freud.diffraction.DiffractionPattern( grid_size=grid_size, output_size=output_size, @@ -272,13 +272,13 @@ def test_cubic_system_parameterized(self): ideal_peaks = {i: False for i in [-2, -1, 0, 1, 2]} + lattice_vector = np.array([length, length, length]) for peak in ideal_peaks: for x, y in xy: k_vector = dp.k_vectors[x, y] - lattice_vector = [1, 1, 1] dot_prod = np.dot(k_vector, lattice_vector) - if dot_prod == (peak * 2 * np.pi): + if np.isclose(dot_prod, peak * 2 * np.pi, atol=1e-2): ideal_peaks[peak] = True assert all(ideal_peaks.values())