Implementation of numcore and small changes

numcore: deleted c1, c2, c3, name changes in phase_mag_real
analytic: name change beta to phi
phasemapper: changes to include numcore
test_compliance: also checkes .pyx files now
setup: now works properly (hopefully on other systems, too)
compare_discs: added script to compare discs

desktop.ini

 [.ShellClassInfo]
-IconResource=C:\Users\Jan\PhD Thesis\Pyramid\icon.ico,0
+IconResource=C:\Users\Jan\Daten\PhD Thesis\Pyramid\icon.ico,0
 [ViewState]
 Mode=
 Vid=

pyramid/analytic.py

@@ -11,10 +11,12 @@ from numpy import pi
 PHI_0 = -2067.83  # magnetic flux in T*nm²
 
 
-def phase_mag_slab(dim, res, beta, center, width, b_0=1):
+def phase_mag_slab(dim, res, phi, center, width, b_0=1):
     '''Get the analytic solution for a phase map of a slab of specified dimensions
     Arguments:
         dim    - the dimensions of the grid, shape(z, y, x)
+        res    - the resolution of the grid (grid spacing) in nm
+        phi    - the azimuthal angle, describing the direction of the magnetized object
         center - the center of the slab in pixel coordinates, shape(z, y, x)
         width  - the width of the slab in pixel coordinates, shape(z, y, x)
         b_0    - magnetic induction corresponding to a magnetization Mo in T (default: 1)
@@ -28,11 +30,11 @@ def phase_mag_slab(dim, res, beta, center, width, b_0=1):
             a = np.log(x**2 + y**2 + 1E-30)
             b = np.arctan(x / (y+1E-30))
             return x*a - 2*x + 2*y*b
-        return coeff * Lz * (- np.cos(beta) * (F0(x-x0-Lx/2, y-y0-Ly/2)
+        return coeff * Lz * (- np.cos(phi) * (F0(x-x0-Lx/2, y-y0-Ly/2)
                                              - F0(x-x0+Lx/2, y-y0-Ly/2)
                                              - F0(x-x0-Lx/2, y-y0+Ly/2)
                                              + F0(x-x0+Lx/2, y-y0+Ly/2))
-                             + np.sin(beta) * (F0(y-y0-Ly/2, x-x0-Lx/2)
+                             + np.sin(phi) * (F0(y-y0-Ly/2, x-x0-Lx/2)
                                              - F0(y-y0+Ly/2, x-x0-Lx/2)
                                              - F0(y-y0-Ly/2, x-x0+Lx/2)
                                              + F0(y-y0+Ly/2, x-x0+Lx/2)))
@@ -50,10 +52,12 @@ def phase_mag_slab(dim, res, beta, center, width, b_0=1):
     return phiMag(xx, yy)
 
 
-def phase_mag_disc(dim, res, beta, center, radius, height, b_0=1):
+def phase_mag_disc(dim, res, phi, center, radius, height, b_0=1):
     '''Get the analytic solution for a phase map of a disc of specified dimensions
     Arguments:
         dim    - the dimensions of the grid, shape(z, y, x)
+        res    - the resolution of the grid (grid spacing) in nm
+        phi    - the azimuthal angle, describing the direction of the magnetized object
         center - the center of the disc in pixel coordinates, shape(z, y, x)
         radius - the radius of the disc in pixel coordinates (scalar value)
         height - the height of the disc in pixel coordinates (scalar value)
@@ -66,7 +70,7 @@ def phase_mag_disc(dim, res, beta, center, radius, height, b_0=1):
     def phiMag(x, y):
         r = np.hypot(x - x0, y - y0)
         r[center[1], center[2]] = 1E-30
-        result = coeff * Lz * ((y - y0) * np.cos(beta) - (x - x0) * np.sin(beta))
+        result = coeff * Lz * ((y - y0) * np.cos(phi) - (x - x0) * np.sin(phi))
         result *= np.where(r <= R, 1, (R / r) ** 2)
         return result
     # Process input parameters:
@@ -84,10 +88,12 @@ def phase_mag_disc(dim, res, beta, center, radius, height, b_0=1):
     return phiMag(xx, yy)
 
 
-def phase_mag_sphere(dim, res, beta, center, radius, b_0=1):
+def phase_mag_sphere(dim, res, phi, center, radius, b_0=1):
     '''Get the analytic solution for a phase map of a sphere of specified dimensions
     Arguments:
         dim    - the dimensions of the grid, shape(z, y, x)
+        res    - the resolution of the grid (grid spacing) in nm
+        phi    - the azimuthal angle, describing the direction of the magnetized object
         center - the center of the sphere in pixel coordinates, shape(z, y, x)
         radius - the radius of the sphere in pixel coordinates (scalar value)
         b_0    - magnetic induction corresponding to a magnetization Mo in T (default: 1)
@@ -99,7 +105,7 @@ def phase_mag_sphere(dim, res, beta, center, radius, b_0=1):
     def phiMag(x, y):
         r = np.hypot(x - x0, y - y0)
         r[center[1], center[2]] = 1E-30
-        result = coeff * R ** 3 / r ** 2 * ((y - y0) * np.cos(beta) - (x - x0) * np.sin(beta))
+        result = coeff * R ** 3 / r ** 2 * ((y - y0) * np.cos(phi) - (x - x0) * np.sin(phi))
         result *= np.where(r > R, 1, (1 - (1 - (r / R) ** 2) ** (3. / 2.)))
         return result
     # Process input parameters:

pyramid/magcreator.py

@@ -243,4 +243,4 @@ def create_mag_dist_vortex(mag_shape, center=None, magnitude=1):
     z_mag = np.zeros(dim)
     y_mag = -np.ones(dim) * np.sin(phi) * mag_shape * magnitude
     x_mag = -np.ones(dim) * np.cos(phi) * mag_shape * magnitude
-    return z_mag, y_mag, x_mag
\ No newline at end of file
+    return z_mag, y_mag, x_mag

pyramid/numcore/c1.pyx

-import math
-
-def great_circle(float lon1, float lat1, float lon2, float lat2):
-    cdef float radius = 3956.0 
-    cdef float pi = 3.14159265
-    cdef float x = pi / 180.0
-    cdef float a, b, theta,c
-
-    a = (90.0 - lat1) * (x)
-    b = (90.0 - lat2) * (x)
-    theta = (lon2 - lon1) * (x)
-    c = math.acos((math.cos(a)*math.cos(b)) + (math.sin(a)*math.sin(b)*math.cos(theta)))
-    return radius*c

pyramid/numcore/c2.pyx

-cdef extern from "math.h":
-    float cosf(float theta)
-    float sinf(float theta)
-    float acosf(float theta)
-
-def great_circle(float lon1,float lat1,float lon2,float lat2):
-    cdef float radius = 3956.0 
-    cdef float pi = 3.14159265
-    cdef float x = pi/180.0
-    cdef float a,b,theta,c
-
-    a = (90.0-lat1)*(x)
-    b = (90.0-lat2)*(x)
-    theta = (lon2-lon1)*(x)
-    c = acosf((cosf(a)*cosf(b)) + (sinf(a)*sinf(b)*cosf(theta)))
-    return radius*c
\ No newline at end of file

pyramid/numcore/c3.pyx

-cdef extern from "math.h":
-    float cosf(float theta)
-    float sinf(float theta)
-    float acosf(float theta)
-
-cdef float _great_circle(float lon1,float lat1,float lon2,float lat2):
-    cdef float radius = 3956.0 
-    cdef float pi = 3.14159265
-    cdef float x = pi/180.0
-    cdef float a,b,theta,c
-
-    a = (90.0-lat1)*(x)
-    b = (90.0-lat2)*(x)
-    theta = (lon2-lon1)*(x)
-    c = acosf((cosf(a)*cosf(b)) + (sinf(a)*sinf(b)*cosf(theta)))
-    return radius*c
-
-def great_circle(float lon1,float lat1,float lon2,float lat2,int num):
-    cdef int i
-    cdef float x
-    for i from 0 <= i < num:
-        x = _great_circle(lon1,lat1,lon2,lat2)
-    return x
\ No newline at end of file

pyramid/numcore/hello.pyx

-def say_hello_to(name):
-    print 'Hello %s!' % name
\ No newline at end of file

pyramid/numcore/phase_mag_real.pyx

@@ -6,27 +6,24 @@ cimport numpy as np
 
 @cython.boundscheck(False)
 @cython.wraparound(False)
-def phase_mag_real_helper_1(
+def phase_mag_real_helper(
         unsigned int v_dim, unsigned int u_dim,
-        double[:, :] phi_u,
-        double[:, :] phi_v,
-        double[:, :] u_mag,
-        double[:, :] v_mag,
+        double[:, :] phi_u, double[:, :] phi_v,
+        double[:, :] u_mag, double[:, :] v_mag,
         double[:, :] phase, float threshold):
-    cdef unsigned int i, j, ii, jj, iii, jjj
-    cdef double u, v
+    cdef unsigned int i, j, p, q, p_c, q_c
+    cdef double u_m, v_m
     for j in range(v_dim):
         for i in range(u_dim):
-            u = u_mag[j, i]
-            v = v_mag[j, i]
-            iii = u_dim - 1 - i
-            jjj = v_dim - 1 - j
-            if abs(u) > threshold:
-               for jj in range(phase.shape[0]):
-                    for ii in range(phase.shape[1]):
-                        phase[jj, ii] += u * phi_u[jjj + jj, iii + ii]
-            if abs(v) > threshold:
-               for jj in range(phase.shape[0]):
-                    for ii in range(phase.shape[1]):
-                        phase[jj, ii] -= v * phi_v[jjj + jj, iii + ii]
-
+            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 * phi_u[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 * phi_v[q_c + q, p_c + p]

pyramid/phasemapper.py

@@ -103,19 +103,77 @@ def phase_mag_real(res, projection, method, b_0=1, jacobi=None):
                     phase -= v_mag[j, i] * phase_v
         ############################### TODO: NUMERICAL CORE  #####################################
     else:  # Without Jacobi matrix (faster)
-##        phasecopy = phase.copy()
-##        start_time = time.time()
-##        numcore.phase_mag_real_helper_1(v_dim, u_dim, phi_u, phi_v, u_mag, v_mag, phasecopy, threshold)
-##        print time.time() - start_time
-##        start_time = time.time()
+        numcore.phase_mag_real_helper(v_dim, u_dim, phi_u, phi_v, u_mag, v_mag, phase, threshold)
+    # Return the phase:
+    return phase
+
+
+def phase_mag_real2(res, projection, method, b_0=1, jacobi=None):
+    '''Calculate phasemap from magnetization data (real space approach).
+    Arguments:
+        res        - the resolution of the grid (grid spacing) in nm
+        projection - projection of a magnetic distribution (created with pyramid.projector)
+        method     - String, describing the method to use for the pixel field ('slab' or 'disc')
+        b_0        - magnetic induction corresponding to a magnetization Mo in T (default: 1)
+        jacobi     - matrix in which to save the jacobi matrix (default: None, faster computation)
+    Returns:
+        the phasemap as a 2 dimensional array
+
+    '''
+    # Function for creating the lookup-tables:
+    def phi_lookup(method, n, m, res, b_0):
+        if method == 'slab':
+            def F_h(n, m):
+                a = np.log(res**2 * (n**2 + m**2))
+                b = np.arctan(n / m)
+                return n*a - 2*n + 2*m*b
+            return coeff * 0.5 * (F_h(n-0.5, m-0.5) - F_h(n+0.5, m-0.5)
+                                - F_h(n-0.5, m+0.5) + F_h(n+0.5, m+0.5))
+        elif method == 'disc':
+            in_or_out = np.logical_not(np.logical_and(n == 0, m == 0))
+            return coeff * m / (n**2 + m**2 + 1E-30) * in_or_out
+    # Process input parameters:
+    v_dim, u_dim = np.shape(projection[0])
+    v_mag, u_mag = projection
+    coeff = -b_0 * res**2 / (2*PHI_0)
+    # Create lookup-tables for the phase of one pixel:
+    u = np.linspace(-(u_dim-1), u_dim-1, num=2*u_dim-1)
+    v = np.linspace(-(v_dim-1), v_dim-1, num=2*v_dim-1)
+    uu, vv = np.meshgrid(u, v)
+    phi_u = phi_lookup(method, uu, vv, res, b_0)
+    phi_v = phi_lookup(method, vv, uu, res, b_0)
+    # Calculation of the phase:
+    phase = np.zeros((v_dim, u_dim))
+    threshold = 0
+    if jacobi is not None:  # With Jacobian matrix (slower)
+        jacobi[:] = 0  # Jacobi matrix --> zeros
+        ############################### TODO: NUMERICAL CORE  #####################################
+        for j in range(v_dim):
+            for i in range(u_dim):
+                phase_u = phi_u[v_dim-1-j:(2*v_dim-1)-j, u_dim-1-i:(2*u_dim-1)-i]
+                jacobi[:, i+u_dim*j] = phase_u.reshape(-1)
+                if abs(u_mag[j, i]) > threshold:
+                    phase += u_mag[j, i] * phase_u
+                phase_v = phi_v[v_dim-1-j:(2*v_dim-1)-j, u_dim-1-i:(2*u_dim-1)-i]
+                jacobi[:, u_dim*v_dim+i+u_dim*j] = -phase_v.reshape(-1)
+                if abs(v_mag[j, i]) > threshold:
+                    phase -= v_mag[j, i] * phase_v
+        ############################### TODO: NUMERICAL CORE  #####################################
+    else:  # Without Jacobi matrix (faster)
+#        phasecopy = phase.copy()
+#        start_time = time.time()
+#        numcore.phase_mag_real_helper_1(v_dim, u_dim, phi_u, phi_v,
+#                                        u_mag, v_mag, phasecopy, threshold)
+#        print 'with numcore   : ', time.time() - start_time
+#        start_time = time.time()
         for j in range(v_dim):
             for i in range(u_dim):
                 if abs(u_mag[j, i]) > threshold:
                     phase += u_mag[j, i] * phi_u[v_dim-1-j:(2*v_dim-1)-j, u_dim-1-i:(2*u_dim-1)-i]
                 if abs(v_mag[j, i]) > threshold:
                     phase -= v_mag[j, i] * phi_v[v_dim-1-j:(2*v_dim-1)-j, u_dim-1-i:(2*u_dim-1)-i]
-##        print time.time() - start_time
-##        print ((phase - phasecopy) ** 2).sum()
+#        print 'without numcore: ', time.time() - start_time
+#        print 'num. difference: ', ((phase - phasecopy) ** 2).sum()
     # Return the phase:
     return phase
 

pyramid/test/test_compliance.py

@@ -24,7 +24,7 @@ class TestCaseCompliance(unittest.TestCase):
         filepaths = []
         for root, dirs, files in os.walk(rootdir):
             for filename in files:
-                if filename.endswith('.py'):
+                if filename.endswith('.py') or filename.endswith('.pyx'):
                     filepaths.append(os.path.join(root, filename))
         return filepaths
 

scripts/compare_discs.py

+# -*- coding: utf-8 -*-
+"""Create the Pyramid-Logo."""
+
+
+import pdb, traceback, sys
+import numpy as np
+from numpy import pi
+import matplotlib.pyplot as plt
+
+import pyramid.magcreator  as mc
+import pyramid.projector   as pj
+import pyramid.phasemapper as pm
+import pyramid.holoimage   as hi
+import pyramid.analytic    as an
+from pyramid.magdata  import MagData
+from pyramid.phasemap import PhaseMap
+
+
+PHI_0 = -2067.83  # magnetic flux in T*nm²
+
+
+def compare_vortices():
+    '''Calculate and display the Pyramid-Logo.
+    Arguments:
+        None
+    Returns:
+        None
+    
+    '''
+    # Input parameters:  
+    dim_list = [(16, 256, 256), (8, 128, 128), (4, 64, 64), (2, 32, 32), (1, 16, 16)]
+    res_list = [4., 8., 16., 32., 64.]  # in nm
+    density = 1
+    phi = pi/2
+    
+    x = []
+    y = []
+    
+    # Analytic solution:
+    L = 1024.  # in nm
+    Lz = 0.5 * 64.  # in nm
+    R = 0.25 * L
+    x0 = L / 2
+    def F(x):
+        coeff = - pi * Lz / (2*PHI_0)
+        result = coeff * (- (x - x0) * np.sin(phi))
+        result *= np.where(np.abs(x - x0) <= R, 1, (R / (x - x0)) ** 2)
+        return result
+    x_an = np.linspace(0, L, 1000)
+    y_an = F(x_an)
+    
+    # Starting magnetic distribution:
+    dim_start = (2*dim_list[0][0], 2*dim_list[0][1], 2*dim_list[0][2])
+    res_start = res_list[0] / 2
+    center = (dim_start[0]/2 - 0.5, dim_start[1]/2 - 0.5, dim_start[2]/2 - 0.5)
+    radius = 0.25 * dim_start[1]
+    height = 0.5* dim_start[0]
+    mag_shape = mc.Shapes.disc(dim_start, center, radius, height)
+    mag_data = MagData(res_start, mc.create_mag_dist(mag_shape, phi))
+    
+    for i, (dim, res) in enumerate(zip(dim_list, res_list)):
+        # Create coarser grid for the magnetization:
+        print 'dim = ', dim, 'res = ', res
+        z_mag = mag_data.magnitude[0].reshape(dim[0], 2, dim[1], 2, dim[2], 2)
+        y_mag = mag_data.magnitude[1].reshape(dim[0], 2, dim[1], 2, dim[2], 2)
+        x_mag = mag_data.magnitude[2].reshape(dim[0], 2, dim[1], 2, dim[2], 2)
+        magnitude = (z_mag.mean(axis=5).mean(axis=3).mean(axis=1),
+                     y_mag.mean(axis=5).mean(axis=3).mean(axis=1),
+                     x_mag.mean(axis=5).mean(axis=3).mean(axis=1))
+        mag_data = MagData(res, magnitude)
+        projection = pj.simple_axis_projection(mag_data)
+        phase_map = PhaseMap(res, pm.phase_mag_real(res, projection, 'slab'))    
+        hi.display_combined(phase_map, density, 'Vortex State, res = {}'.format(res))
+        x.append(np.linspace(0, dim[1]*res, dim[1]))
+        y.append(phase_map.phase[dim[1]/2, :])
+        fig = plt.figure()
+        fig.add_subplot(1, 1, 1)
+        plt.plot(x[i], y[i])
+
+    # Plot:
+    fig = plt.figure()
+    axis = fig.add_subplot(1, 1, 1)
+    axis.plot(x[0], y[0], 'r',
+              x[1], y[1], 'm',
+              x[2], y[2], 'y',
+              x[3], y[3], 'g',
+              x[4], y[4], 'c',
+              x_an, y_an, 'k')
+    axis.set_xlabel('x [nm]')
+    axis.set_ylabel('phase')
+  
+    
+if __name__ == "__main__":
+    try:
+        compare_vortices()
+    except:
+        type, value, tb = sys.exc_info()
+        traceback.print_exc()
+        pdb.post_mortem(tb)
\ No newline at end of file

scripts/compare_method_errors_res.py

@@ -47,23 +47,28 @@ def compare_method_errors_res():
     data_sl_p_fourier0 = np.zeros((3, len(res_list)))
     data_sl_w_fourier0 = np.zeros((3, len(res_list)))
     data_disc_fourier0 = np.zeros((3, len(res_list)))
-    
+    data_vort_fourier0 = np.zeros((3, len(res_list)))
+
     data_sl_p_fourier1 = np.zeros((3, len(res_list)))
     data_sl_w_fourier1 = np.zeros((3, len(res_list)))
     data_disc_fourier1 = np.zeros((3, len(res_list)))
-    
+    data_vort_fourier1 = np.zeros((3, len(res_list)))
+
     data_sl_p_fourier20 = np.zeros((3, len(res_list)))
     data_sl_w_fourier20 = np.zeros((3, len(res_list)))
     data_disc_fourier20 = np.zeros((3, len(res_list)))
-    
+    data_vort_fourier20 = np.zeros((3, len(res_list)))
+
     data_sl_p_real_s = np.zeros((3, len(res_list)))
     data_sl_w_real_s = np.zeros((3, len(res_list)))
     data_disc_real_s = np.zeros((3, len(res_list)))
-    
+    data_vort_real_s = np.zeros((3, len(res_list)))
+
     data_sl_p_real_d= np.zeros((3, len(res_list)))
     data_sl_w_real_d = np.zeros((3, len(res_list)))
     data_disc_real_d = np.zeros((3, len(res_list)))
-    
+    data_vort_real_d = np.zeros((3, len(res_list)))
+
     
     '''CREATE DATA ARRAYS'''
         

scripts/compare_vortices.py

@@ -15,6 +15,9 @@ from pyramid.magdata  import MagData
 from pyramid.phasemap import PhaseMap
 
 
+PHI_0 = -2067.83  # magnetic flux in T*nm²
+
+
 def compare_vortices():
     '''Calculate and display the Pyramid-Logo.
     Arguments:
@@ -24,59 +27,66 @@ def compare_vortices():
     
     '''
     # Input parameters:  
-    dim_list = [(1, 512, 512), (1, 256, 256), (1, 128, 128), (1, 64, 64), (1, 32, 32), (1, 16, 16)]
-    res_list = [2., 4., 8., 16., 32., 64.]  # in nm
-    density = 10
-    height = 1
+    dim_list = [(16, 256, 256), (8, 128, 128), (4, 64, 64), (2, 32, 32), (1, 16, 16)]
+    res_list = [4., 8., 16., 32., 64.]  # in nm
+    density = 1
     
     x = []
     y = []
     
     # Starting magnetic distribution:
-    dim_start = (1, 2*dim_list[0][1], 2*dim_list[0][2])
+    dim_start = (2*dim_list[0][0], 2*dim_list[0][1], 2*dim_list[0][2])
     res_start = res_list[0] / 2
-    print 'dim_start = ', dim_start
-    print 'res_start = ', res_start
-    center = (0, dim_start[1]/2 - 0.5, dim_start[2]/2 - 0.5)
+    center = (dim_start[0]/2 - 0.5, dim_start[1]/2 - 0.5, dim_start[2]/2 - 0.5)
     radius = 0.25 * dim_start[1]
+    height = 0.5* dim_start[0]
     mag_shape = mc.Shapes.disc(dim_start, center, radius, height)
     mag_data = MagData(res_start, mc.create_mag_dist_vortex(mag_shape))
-    #mag_data.quiver_plot()
+    
+    # Analytic solution:
+    L = 1024.  # in nm
+    Lz = 0.5 * 64.  # in nm
+    R = 0.25 * L
+    x0 = L / 2
+    def F(x):
+        coeff = pi*Lz/PHI_0
+        result = coeff * np.where(np.abs(x - x0) <= R, (np.abs(x-x0)-R), 0)
+        return result
+    x_an = np.linspace(0, L, 1000)
+    y_an = F(x_an)
     
     for i, (dim, res) in enumerate(zip(dim_list, res_list)):
         # Create coarser grid for the magnetization:
-        z_mag = mag_data.magnitude[0].reshape(1, dim[1], 2, dim[1], 2)
-        y_mag = mag_data.magnitude[1].reshape(1, dim[1], 2, dim[1], 2)
-        x_mag = mag_data.magnitude[2].reshape(1, dim[1], 2, dim[1], 2)
-        print 'dim = ', dim
-        print 'res = ', res
-        magnitude = (z_mag.mean(axis=4).mean(axis=2) / 2,
-                     y_mag.mean(axis=4).mean(axis=2) / 2,
-                     x_mag.mean(axis=4).mean(axis=2) / 2)
+        print 'dim = ', dim, 'res = ', res
+        z_mag = mag_data.magnitude[0].reshape(dim[0], 2, dim[1], 2, dim[2], 2)
+        y_mag = mag_data.magnitude[1].reshape(dim[0], 2, dim[1], 2, dim[2], 2)
+        x_mag = mag_data.magnitude[2].reshape(dim[0], 2, dim[1], 2, dim[2], 2)
+        magnitude = (z_mag.mean(axis=5).mean(axis=3).mean(axis=1),
+                     y_mag.mean(axis=5).mean(axis=3).mean(axis=1),
+                     x_mag.mean(axis=5).mean(axis=3).mean(axis=1))
         mag_data = MagData(res, magnitude)
-        print np.shape(mag_data.magnitude[0])
         #mag_data.quiver_plot()
         projection = pj.simple_axis_projection(mag_data)
         phase_map = PhaseMap(res, pm.phase_mag_real(res, projection, 'slab'))    
         hi.display_combined(phase_map, density, 'Vortex State, res = {}'.format(res))
         x.append(np.linspace(0, dim[1]*res, dim[1]))
-        y.append(phase_map.phase[dim[1]/2, :]/pi)
+        y.append(phase_map.phase[dim[1]/2, :])
         fig = plt.figure()
         fig.add_subplot(1, 1, 1)
         plt.plot(x[i], y[i])
-#    
+
     # Plot:
     fig = plt.figure()
     axis = fig.add_subplot(1, 1, 1)
-    axis.plot(x[0], y[0], 'k',
-              x[1], y[1], 'r',
-              x[2], y[2], 'm',
-              x[3], y[3], 'y',
-              x[4], y[4], 'c',
-              x[5], y[5], 'g',
-              x[6], y[6], 'b')
+    axis.plot(x[0], y[0], 'r',
+              x[1], y[1], 'm',
+              x[2], y[2], 'y',
+              x[3], y[3], 'g',
+              x[4], y[4], 'c',#x[5], y[5], 'b',
+              x_an, y_an, 'k')
     axis.set_xlabel('x [nm]')
     axis.set_ylabel('phase')
+   
     
 if __name__ == "__main__":
     try:

setup.py

@@ -2,23 +2,38 @@
 """
 Created on Fri May 03 10:27:04 2013
 
-
 @author: Jan
 """
 
-#call with: python setup.py build_ext --inplace --compiler=mingw32
 
-import os
-import glob
+# Build extensions: 'python setup.py build_ext -i clean'
+# Install package:  'python setup.py install clean'
+
+
+import numpy
 from distutils.core import setup
-from Cython.Build import cythonize
+from distutils.extension import Extension
+from Cython.Distutils import build_ext
+
 
 setup(
-    name = 'Pyramid',
-    version = '0.1',
-    description = 'PYthon based Reconstruction Algorithm for MagnetIc Distributions',
-    author = 'Jan Caron',
-    author_email = 'j.caron@fz-juelich.de',
-    packages = ['pyramid', 'pyramid.numcore'],
-    ext_modules = cythonize(glob.glob(os.path.join('pyramid', 'numcore', '*.pyx')))
+
+      name = 'Pyramid',
+      version = '0.1',
+      description = 'PYthon based Reconstruction Algorithm for MagnetIc Distributions',
+      author = 'Jan Caron',
+      author_email = 'j.caron@fz-juelich.de',
+      
+      packages = ['pyramid', 'pyramid.numcore'],
+      include_dirs = [numpy.get_include()],
+                      
+      cmdclass = {'build_ext': build_ext},              
+      ext_package = 'pyramid/numcore',
+      ext_modules = [
+          Extension('phase_mag_real', ['pyramid/numcore/phase_mag_real.pyx'], 
+                    include_dirs = [numpy.get_include(), numpy.get_numarray_include()],
+                    extra_compile_args=["-march=pentium", "-mtune=pentium"]
+                    )
+          ]
+
 )