Deprecated numcore module! (No longer needed!)

pyramid/numcore/__init__.py

-# -*- coding: utf-8 -*-
-"""Package for numerical core routines which enable fast computations.
-
-Modules
--------
-phasemapper_core
-    Provides numerical core routines for :class:`~.PhaseMapper` class.
-kernel_core
-    Provides numerical core routines for the :class:`~.Kernel` class.
-
-Notes
------
-Packages are written in `pyx`-format for the use with :mod:`~.cython`.
-
-"""

pyramid/numcore/phasemapper_core.pyx

-# -*- coding: utf-8 -*-
-"""Numerical core routines for the :mod:`~.pyramid.phasemapper` module.
-
-Provides helper functions to speed up phase mapping calculations in the
-:mod:`~pyramid.phasemapper` module, by using C-speed for the for-loops and by omitting boundary
-and wraparound checks.
-
-"""
-
-
-import numpy as np
-import math
-
-cimport cython
-cimport numpy as np
-
-
-@cython.boundscheck(False)
-@cython.wraparound(False)
-def phasemapper_real_convolve(
-        unsigned int v_dim, unsigned int u_dim,
-        float[:, :] v_phi, float[:, :] u_phi,
-        float[:, :] v_mag, float[:, :] u_mag,
-        float[:, :] phase, float threshold):
-    '''Numerical core routine for the phase calculation using the real space approach.
-
-    Parameters
-    ----------
-    v_dim, u_dim : int
-        Dimensions of the projection along the two major axes.
-    v_phi, u_phi : :class:`~numpy.ndarray` (N=2)
-        Lookup tables for the pixel fields oriented in `u`- and `v`-direction.
-    v_mag, u_mag : :class:`~numpy.ndarray` (N=2)
-        Magnetization components in `u`- and `v`-direction.
-    phase : :class:`~numpy.ndarray` (N=2)
-        Matrix in which the resulting magnetic phase map should be stored.
-    threshold : float
-        The `threshold` determines which pixels contribute to the magnetic phase.
-
-    Returns
-    -------
-    None
-
-    '''
-    cdef unsigned int i, j, p, q, p_c, q_c
-    cdef float u_m, v_m
-    for j in range(v_dim):
-        for i in range(u_dim):
-            u_m = u_mag[j, i]
-            v_m = v_mag[j, i]
-            p_c = u_dim - 1 - i
-            q_c = v_dim - 1 - j
-            if abs(u_m) > threshold:
-                for q in range(v_dim):
-                    for p in range(u_dim):
-                        phase[q, p] += u_m * u_phi[q_c + q, p_c + p]
-            if abs(v_m) > threshold:
-                for q in range(v_dim):
-                    for p in range(u_dim):
-                        phase[q, p] -= v_m * v_phi[q_c + q, p_c + p]
-
-
-@cython.boundscheck(False)
-@cython.wraparound(False)
-def jac_dot_real_convolve(
-        unsigned int v_dim, unsigned int u_dim,
-        float[:, :] u_phi, float[:, :] v_phi,
-        float[:] vector,
-        float[:] result):
-    '''Numerical core routine for the Jacobi matrix multiplication for the phase mapper.
-
-    Parameters
-    ----------
-    v_dim, u_dim : int
-        Dimensions of the projection along the two major axes.
-    v_phi, u_phi : :class:`~numpy.ndarray` (N=2)
-        Lookup tables for the pixel fields oriented in `u`- and `v`-direction.
-    vector : :class:`~numpy.ndarray` (N=1)
-        Input vector which should be multiplied by the Jacobi-matrix.
-    result : :class:`~numpy.ndarray` (N=1)
-        Vector in which the result of the multiplication should be stored.
-
-    Returns
-    -------
-    None
-
-    Notes
-    -----
-    The strategy involves iterating over the magnetization first and over the kernel in an inner
-    iteration loop.
-
-    '''
-    cdef unsigned int s, i, j, u_min, u_max, v_min, v_max, r, u, v, size
-    size = u_dim * v_dim  # Number of pixels
-    s = 0  # Current contributing pixel (numbered consecutively)
-    # Iterate over all contributingh pixels:
-    for j in range(v_dim):
-        for i in range(u_dim):
-            u_min = (u_dim-1) - i  # u_max = u_min + u_dim
-            v_min = (v_dim-1) - j  # v_max = v_min + v_dim
-            r = 0  # Current result component / affected pixel (numbered consecutively)
-            # Go over the current kernel cutout [v_min:v_max, u_min:u_max]:
-            for v in range(v_min, v_min + v_dim):
-                for u in range(u_min, u_min + u_dim):
-                    result[r] += vector[s] * u_phi[v, u]
-                    result[r] -= vector[s+size] * v_phi[v, u]
-                    r += 1
-            s += 1
-
-
-@cython.boundscheck(False)
-@cython.wraparound(False)
-def jac_T_dot_real_convolve(
-        unsigned int v_dim, unsigned int u_dim,
-        float[:, :] u_phi, float[:, :] v_phi,
-        float[:] vector,
-        float[:] result):
-    '''Core routine for the transposed Jacobi multiplication for the phase mapper.
-
-    Parameters
-    ----------
-    v_dim, u_dim : int
-        Dimensions of the projection along the two major axes.
-    v_phi, u_phi : :class:`~numpy.ndarray` (N=2)
-        Lookup tables for the pixel fields oriented in `u`- and `v`-direction.
-    vector : :class:`~numpy.ndarray` (N=1)
-        Input vector which should be multiplied by the transposed Jacobi-matrix.
-    result : :class:`~numpy.ndarray` (N=1)
-        Vector in which the result of the multiplication should be stored.
-
-    Returns
-    -------
-    None
-
-    Notes
-    -----
-    The strategy involves iterating over the magnetization first and over the kernel in an inner
-    iteration loop.
-
-    '''
-    cdef unsigned int s, i, j, u_min, u_max, v_min, v_max, r, u, v, size
-    size = u_dim * v_dim  # Number of pixels
-    r = 0  # Current result component / contributing pixel (numbered consecutively)
-    # Iterate over all contributingh pixels:
-    for j in range(v_dim):
-        for i in range(u_dim):
-            u_min = (u_dim-1) - i  # u_max = u_min + u_dim
-            v_min = (v_dim-1) - j  # v_max = v_min + v_dim
-            s = 0  # Current affected pixel (numbered consecutively)
-            # Go over the current kernel cutout [v_min:v_max, u_min:u_max]:
-            for v in range(v_min, v_min + v_dim):
-                for u in range(u_min, u_min + u_dim):
-                    result[r] += vector[s] * u_phi[v, u]
-                    result[r+size] -= vector[s] * v_phi[v, u]
-                    s += 1
-            r += 1

pyramid/phasemapper.py

@@ -11,7 +11,6 @@ import abc
 import logging
 
 import numpy as np
-import pyramid.numcore.phasemapper_core as nc
 from numpy import pi
 
 from pyramid import fft
@@ -20,7 +19,7 @@ from pyramid.fielddata import VectorData, ScalarData
 from pyramid.phasemap import PhaseMap
 from pyramid.projector import RotTiltProjector, XTiltProjector, YTiltProjector, SimpleProjector
 
-__all__ = ['PhaseMapperRDFC', 'PhaseMapperRDRC', 'PhaseMapperFDFC', 'PhaseMapperMIP', 'pm']
+__all__ = ['PhaseMapperRDFC', 'PhaseMapperFDFC', 'PhaseMapperMIP', 'pm']
 _log = logging.getLogger(__name__)
 
 PHI_0 = 2067.83  # magnetic flux in T*nm²
@@ -193,159 +192,6 @@ class PhaseMapperRDFC(PhaseMapper):
         return result
 
 
-class PhaseMapperRDRC(PhaseMapper):
-    """Class representing a phase mapping strategy using real space discretization.
-
-    The :class:`~.PMReal` class represents a phase mapping strategy involving discretization in
-    real space. It directly takes :class:`~.VectorData` objects and returns :class:`~.PhaseMap`
-    objects.
-
-    Attributes
-    ----------
-    kernel : :class:`~pyramid.Kernel`
-        Convolution kernel, representing the phase contribution of one single magnetized pixel.
-    threshold : float, optional
-        Threshold determining for which voxels the phase contribution will be calculated. The
-        default is 0, which means that all voxels with non-zero magnetization will contribute.
-        Should be above noise level.
-    numcore : boolean, optional
-        Boolean choosing if Cython enhanced routines from the :mod:`~.pyramid.numcore` module
-        should be used. Default is True.
-    m: int
-        Size of the image space.
-    n: int
-        Size of the input space.
-
-    """
-
-    _log = logging.getLogger(__name__ + '.PhaseMapperRDRC')
-
-    def __init__(self, kernel, threshold=0, numcore=True):
-        self._log.debug('Calling __init__')
-        self.kernel = kernel
-        self.threshold = threshold
-        self.numcore = numcore
-        self.m = np.prod(kernel.dim_uv)
-        self.n = 2 * self.m
-        self._log.debug('Created ' + str(self))
-
-    def __repr__(self):
-        self._log.debug('Calling __repr__')
-        return '%s(kernel=%r, threshold=%r, numcore=%r)' % \
-               (self.__class__, self.kernel, self.threshold, self.numcore)
-
-    def __str__(self):
-        self._log.debug('Calling __str__')
-        return 'PhaseMapperRDRC(kernel=%s, threshold=%s, numcore=%s)' % \
-               (self.kernel, self.threshold, self.numcore)
-
-    def __call__(self, mag_data):
-        self._log.debug('Calling __call__')
-        dim_uv = self.kernel.dim_uv
-        assert isinstance(mag_data, VectorData), 'Only VectorData objects can be mapped!'
-        assert mag_data.a == self.kernel.a, 'Grid spacing has to match!'
-        assert mag_data.dim[0] == 1, 'Magnetic distribution must be 2-dimensional!'
-        assert mag_data.dim[1:3] == dim_uv, 'Dimensions do not match!'
-        # Process input parameters:
-        u_mag, v_mag = mag_data.field[0:2, 0, ...]
-        # Get kernel (lookup-tables for the phase of one pixel):
-        u_phi = self.kernel.u
-        v_phi = self.kernel.v
-        # Calculation of the phase:
-        phase = np.zeros(dim_uv, dtype=np.float32)
-        if self.numcore:
-            nc.phasemapper_real_convolve(dim_uv[0], dim_uv[1], u_phi, v_phi,
-                                         v_mag, u_mag, phase, self.threshold)
-        else:
-            for j in range(dim_uv[0]):
-                for i in range(dim_uv[1]):
-                    v_phase = v_phi[dim_uv[0] - 1 - j:(2 * dim_uv[0] - 1) - j,
-                                    dim_uv[1] - 1 - i:(2 * dim_uv[1] - 1) - i]
-                    if abs(v_mag[j, i]) > self.threshold:
-                        phase += u_mag[j, i] * v_phase
-                    u_phase = u_phi[dim_uv[0] - 1 - j:(2 * dim_uv[0] - 1) - j,
-                                    dim_uv[1] - 1 - i:(2 * dim_uv[1] - 1) - i]
-                    if abs(u_mag[j, i]) > self.threshold:
-                        phase -= v_mag[j, i] * u_phase
-        # Return the phase:
-        return PhaseMap(mag_data.a, phase)
-
-    def jac_dot(self, vector):
-        """Calculate the product of the Jacobi matrix with a given `vector`.
-
-        Parameters
-        ----------
-        vector : :class:`~numpy.ndarray` (N=1)
-            Vectorized form of the magnetization in `u`- and `v`-direction of every pixel
-            (row-wise). The first ``N**2`` elements have to correspond to the `u`-, the next
-            ``N**2`` elements to the `v`-component of the magnetization.
-
-        Returns
-        -------
-        result : :class:`~numpy.ndarray` (N=1)
-            Product of the Jacobi matrix (which is not explicitely calculated) with the vector.
-
-        """
-        assert len(vector) == self.n, \
-            'vector size not compatible! vector: {}, size: {}'.format(len(vector), self.n)
-        v_dim, u_dim = self.kernel.dim_uv
-        result = np.zeros(self.m, dtype=np.float32)
-        if self.numcore:
-            if vector.dtype != np.float32:
-                vector = vector.astype(np.float32)
-            nc.jac_dot_real_convolve(v_dim, u_dim, self.kernel.v, self.kernel.u, vector, result)
-        else:
-            # Iterate over all contributing pixels (numbered consecutively)
-            for s in range(self.m):  # column-wise (two columns at a time, u- and v-component)
-                i = s % u_dim  # u-coordinate of current contributing pixel
-                j = int(s / u_dim)  # v-coordinate of current ccontributing pixel
-                u_min = (u_dim - 1) - i  # u_dim-1: center of the kernel
-                u_max = (2 * u_dim - 1) - i  # = u_min + u_dim
-                v_min = (v_dim - 1) - j  # v_dim-1: center of the kernel
-                v_max = (2 * v_dim - 1) - j  # = v_min + v_dim
-                result += vector[s] * self.kernel.v[v_min:v_max, u_min:u_max].reshape(-1)
-                result -= vector[s + self.m] * self.kernel.u[v_min:v_max, u_min:u_max].reshape(-1)
-        return result
-
-    def jac_T_dot(self, vector):
-        """Calculate the product of the transposed Jacobi matrix with a given `vector`.
-
-        Parameters
-        ----------
-        vector : :class:`~numpy.ndarray` (N=1)
-            Vector with ``N**2`` entries which represents a matrix with dimensions like a scalar
-            phasemap.
-
-        Returns
-        -------
-        result : :class:`~numpy.ndarray` (N=1)
-            Product of the transposed Jacobi matrix (which is not explicitely calculated) with
-            the vector, which has ``2*N**2`` entries like a 2D magnetic projection.
-
-        """
-        assert len(vector) == self.m, \
-            'vector size not compatible! vector: {}, size: {}'.format(len(vector), self.m)
-        v_dim, u_dim = self.kernel.dim_uv
-        result = np.zeros(self.n, dtype=np.float32)
-        if self.numcore:
-            if vector.dtype != np.float32:
-                vector = vector.astype(np.float32)
-            nc.jac_T_dot_real_convolve(v_dim, u_dim, self.kernel.v, self.kernel.u, vector, result)
-        else:
-            # Iterate over all contributing pixels (numbered consecutively):
-            for s in range(self.m):  # row-wise (two rows at a time, u- and v-component)
-                i = s % u_dim  # u-coordinate of current contributing pixel
-                j = int(s / u_dim)  # v-coordinate of current contributing pixel
-                u_min = (u_dim - 1) - i  # u_dim-1: center of the kernel
-                u_max = (2 * u_dim - 1) - i  # = u_min + u_dim
-                v_min = (v_dim - 1) - j  # v_dim-1: center of the kernel
-                v_max = (2 * v_dim - 1) - j  # = v_min + v_dim
-                result[s] = np.sum(vector * self.kernel.v[v_min:v_max, u_min:u_max].reshape(-1))
-                result[s + self.m] = np.sum(vector *
-                                            -self.kernel.u[v_min:v_max, u_min:u_max].reshape(-1))
-        return result
-
-
 class PhaseMapperFDFC(PhaseMapper):
     """Class representing a phase mapping strategy using a discretization in Fourier space.
 

pyramid/tests/test_phasemapper.py

@@ -10,8 +10,7 @@ from numpy.testing import assert_allclose
 from pyramid.kernel import Kernel
 from pyramid.fielddata import VectorData, ScalarData
 from pyramid.phasemap import PhaseMap
-from pyramid.phasemapper import PhaseMapperRDFC, PhaseMapperRDRC, PhaseMapperFDFC
-from pyramid.phasemapper import PhaseMapperMIP, pm
+from pyramid.phasemapper import PhaseMapperRDFC, PhaseMapperFDFC, PhaseMapperMIP, pm
 
 
 class TestCasePhaseMapperRDFC(unittest.TestCase):
@@ -52,44 +51,6 @@ class TestCasePhaseMapperRDFC(unittest.TestCase):
                         err_msg='Unexpected behaviour in the the transposed jacobi matrix!')
 
 
-class TestCasePhaseMapperRDRC(unittest.TestCase):
-    def setUp(self):
-        self.path = os.path.join(os.path.dirname(os.path.realpath(__file__)), 'test_phasemapper')
-        self.mag_proj = VectorData.load_from_hdf5(os.path.join(self.path, 'mag_proj.hdf5'))
-        self.mapper = PhaseMapperRDRC(Kernel(self.mag_proj.a, self.mag_proj.dim[1:]))
-
-    def tearDown(self):
-        self.path = None
-        self.mag_proj = None
-        self.mapper = None
-
-    def test_PhaseMapperRDRC_call(self):
-        phase_ref = PhaseMap.load_from_hdf5(os.path.join(self.path, 'phase_map.hdf5'))
-        phase_map = self.mapper(self.mag_proj)
-        assert_allclose(phase_map.phase, phase_ref.phase, atol=1E-7,
-                        err_msg='Unexpected behavior in __call__()!')
-        assert_allclose(phase_map.a, phase_ref.a, err_msg='Unexpected behavior in __call__()!')
-
-    def test_PhaseMapperRDRC_jac_dot(self):
-        phase = self.mapper(self.mag_proj).phase
-        mag_proj_vec = self.mag_proj.field[:2, ...].flatten()
-        phase_jac = self.mapper.jac_dot(mag_proj_vec).reshape(self.mapper.kernel.dim_uv)
-        assert_allclose(phase, phase_jac, atol=1E-7,
-                        err_msg='Inconsistency between __call__() and jac_dot()!')
-        n = self.mapper.n
-        jac = np.array([self.mapper.jac_dot(np.eye(n)[:, i]) for i in range(n)]).T
-        jac_ref = np.load(os.path.join(self.path, 'jac.npy'))
-        assert_allclose(jac, jac_ref, atol=1E-7,
-                        err_msg='Unexpected behaviour in the the jacobi matrix!')
-
-    def test_PhaseMapperRDRC_jac_T_dot(self):
-        m = self.mapper.m
-        jac_T = np.array([self.mapper.jac_T_dot(np.eye(m)[:, i]) for i in range(m)]).T
-        jac_T_ref = np.load(os.path.join(self.path, 'jac.npy')).T
-        assert_allclose(jac_T, jac_T_ref, atol=1E-7,
-                        err_msg='Unexpected behaviour in the the transposed jacobi matrix!')
-
-
 class TestCasePhaseMapperFDFCpad0(unittest.TestCase):
     def setUp(self):
         self.path = os.path.join(os.path.dirname(os.path.realpath(__file__)), 'test_phasemapper')

setup.py

@@ -147,11 +147,8 @@ setup(name=DISTNAME,
       packages=find_packages(exclude=['tests']),
       include_dirs=[numpy.get_include()],
       requires=['numpy', 'scipy', 'matplotlib', 'Pillow',
-                'mayavi', 'pyfftw', 'hyperspy', 'Cython', 'nose'],
+                'mayavi', 'pyfftw', 'hyperspy', 'nose'],
       scripts=get_files('scripts'),
       test_suite='nose.collector',
-      cmdclass={'build_ext': build_ext, 'build': build},
-      ext_package='pyramid/numcore',
-      ext_modules=[Extension('phasemapper_core', ['pyramid/numcore/phasemapper_core.pyx'],
-                             include_dirs=[numpy.get_include()])])
+      cmdclass={'build_ext': build_ext, 'build': build})
 print('-------------------------------------------------------------------------------\n')