From 0bfdb8964ac1529c8cb4bbf18f0974b5a0ec3daa Mon Sep 17 00:00:00 2001
From: Jan Caron <j.caron@fz-juelich.de>
Date: Wed, 15 May 2013 09:50:45 +0200
Subject: [PATCH] Middle Stage of converting the program structure (not working
correctly9 For Pushing Purposes
---
pyramid/dataloader.py | 45 ----
pyramid/holoimage.py | 198 +++++++++---------
pyramid/magcreator.py | 94 ++++-----
pyramid/magdata.py | 128 ++++++++++--
pyramid/numcore/setup.py | 30 ---
pyramid/phasemap.py | 288 ++++++++------------------
pyramid/phasemapper.py | 177 ++++++++++++++++
pyramid/projector.py | 17 ++
pyramid/reconstructor.py | 7 +
scripts/create_random_distribution.py | 26 ++-
10 files changed, 564 insertions(+), 446 deletions(-)
delete mode 100644 pyramid/dataloader.py
delete mode 100644 pyramid/numcore/setup.py
create mode 100644 pyramid/phasemapper.py
create mode 100644 pyramid/projector.py
create mode 100644 pyramid/reconstructor.py
diff --git a/pyramid/dataloader.py b/pyramid/dataloader.py
deleted file mode 100644
index 5900c56..0000000
--- a/pyramid/dataloader.py
+++ /dev/null
@@ -1,45 +0,0 @@
-# -*- coding: utf-8 -*-
-"""Load magnetization data from LLG files."""
-
-
-import numpy as np
-
-
-class MagDataLLG:
-
- '''An object storing magnetization data loaded from a LLG-file.'''
-
- def __init__(self, filename):
- '''Load magnetization in LLG-file format.
- Arguments:
- filename - the name of the file where the data are stored
- Returns:
- None
-
- '''
- scale = 1.0E-9 / 1.0E-2 #from cm to nm
- data = np.genfromtxt(filename, skip_header=2)
- x_dim, y_dim, z_dim = np.genfromtxt(filename, dtype=int,
- skip_header=1,
- skip_footer=len(data[:, 0]))
- res = (data[1, 0] - data[0, 0]) / scale
- x_len, y_len, z_len = [data[-1, i]/scale+res/2 for i in range(3)]
- x_mag, y_mag, z_mag = [data[:, i].reshape(z_dim, y_dim, x_dim)
- for i in range(3,6)]
- #Reshape in Python and Igor is different,
- #Python fills rows first, Igor columns!
- self.filename = filename
- self.res = res
- self.dim = (z_dim, y_dim, x_dim)
- self.length = (z_len, y_len, x_len)
- self.magnitude = (z_mag, y_mag, x_mag)
-
- def __str__(self):
- '''Return the filename of the loaded file.
- Arguments:
- None
- Returns:
- the name of the loaded file as a string
-
- '''
- return self.filename
\ No newline at end of file
diff --git a/pyramid/holoimage.py b/pyramid/holoimage.py
index 756bd29..a419b18 100644
--- a/pyramid/holoimage.py
+++ b/pyramid/holoimage.py
@@ -9,111 +9,117 @@ from numpy import pi
from PIL import Image
-CDICT = {'red': [(0.00, 1.0, 0.0),
- (0.25, 1.0, 1.0),
- (0.50, 1.0, 1.0),
- (0.75, 0.0, 0.0),
- (1.00, 0.0, 1.0)],
-
- 'green': [(0.00, 0.0, 0.0),
- (0.25, 0.0, 0.0),
- (0.50, 1.0, 1.0),
- (0.75, 1.0, 1.0),
- (1.00, 0.0, 1.0)],
-
- 'blue': [(0.00, 1.0, 1.0),
- (0.25, 0.0, 0.0),
- (0.50, 0.0, 0.0),
- (0.75, 0.0, 0.0),
- (1.00, 1.0, 1.0)]}
-
-HOLO_CMAP = mpl.colors.LinearSegmentedColormap('my_colormap', CDICT, 256)
-
-
-def holo_image(phase, res, density):
- '''Returns a holography image with color-encoded gradient direction.
- Arguments:
- phase - the phasemap that should be displayed
- res - the resolution of the phasemap
- density - the factor for determining the number of contour lines
- Returns:
- Image
-
- '''
- img_holo = (1 + np.cos(density * phase * pi/2)) /2
+class HoloImage:
- phase_grad_y, phase_grad_x = np.gradient(phase, res, res)
+ '''An object storing magnetization data.'''
- phase_angle = (1 - np.arctan2(phase_grad_y, phase_grad_x)/pi) / 2
-
- # TODO: Delete
-# import pdb; pdb.set_trace()
+
+ CDICT = {'red': [(0.00, 1.0, 0.0),
+ (0.25, 1.0, 1.0),
+ (0.50, 1.0, 1.0),
+ (0.75, 0.0, 0.0),
+ (1.00, 0.0, 1.0)],
- phase_magnitude = np.hypot(phase_grad_x, phase_grad_y)
- phase_magnitude /= phase_magnitude.max()
- phase_magnitude = np.sin(phase_magnitude * pi / 2)
+ 'green': [(0.00, 0.0, 0.0),
+ (0.25, 0.0, 0.0),
+ (0.50, 1.0, 1.0),
+ (0.75, 1.0, 1.0),
+ (1.00, 0.0, 1.0)],
- cmap = HOLO_CMAP
- rgba_img = cmap(phase_angle)
- rgb_img = np.delete(rgba_img, 3, 2)
- red, green, blue = rgb_img[:, :, 0], rgb_img[:, :, 1], rgb_img[:, :, 2]
- red *= 255.999 * img_holo * phase_magnitude
- green *= 255.999 * img_holo * phase_magnitude
- blue *= 255.999 * img_holo * phase_magnitude
- rgb = np.dstack((red, green, blue)).astype(np.uint8)
- # TODO: Which one?
- rgb = (255.999 * img_holo.T * phase_magnitude.T
- * rgba_img[:, :, :3].T).T.astype(np.uint8)
- img = Image.fromarray(rgb)
+ 'blue': [(0.00, 1.0, 1.0),
+ (0.25, 0.0, 0.0),
+ (0.50, 0.0, 0.0),
+ (0.75, 0.0, 0.0),
+ (1.00, 1.0, 1.0)]}
+
+ HOLO_CMAP = mpl.colors.LinearSegmentedColormap('my_colormap', CDICT, 256)
- return img
-
-def make_color_wheel():
- '''Display a color wheel for the gradient direction.
- Arguments:
- None
- Returns:
- None
+ def __init__(self, phase, res, density):
+ '''Returns a holography image with color-encoded gradient direction.
+ Arguments:
+ phase - the phasemap that should be displayed
+ res - the resolution of the phasemap
+ density - the factor for determining the number of contour lines
+ Returns:
+ Image
+
+ '''
+ img_holo = (1 + np.cos(density * phase * pi/2)) /2
- '''
- x = np.linspace(-256, 256, num=512)
- y = np.linspace(-256, 256, num=512)
- xx, yy = np.meshgrid(x, y)
- color_wheel_angle = (1 - np.arctan2(xx, -yy)/pi) / 2
+ phase_grad_y, phase_grad_x = np.gradient(phase, res, res)
+
+ phase_angle = (1 - np.arctan2(phase_grad_y, phase_grad_x)/pi) / 2
+
+ # TODO: Delete
+ # import pdb; pdb.set_trace()
+
+ phase_magnitude = np.hypot(phase_grad_x, phase_grad_y)
+ phase_magnitude /= phase_magnitude.max()
+ phase_magnitude = np.sin(phase_magnitude * pi / 2)
+
+ cmap = HoloImage.HOLO_CMAP
+ rgba_img = cmap(phase_angle)
+ rgb_img = np.delete(rgba_img, 3, 2)
+ red, green, blue = rgb_img[:, :, 0], rgb_img[:, :, 1], rgb_img[:, :, 2]
+ red *= 255.999 * img_holo * phase_magnitude
+ green *= 255.999 * img_holo * phase_magnitude
+ blue *= 255.999 * img_holo * phase_magnitude
+ rgb = np.dstack((red, green, blue)).astype(np.uint8)
+ # TODO: Which one?
+ rgb = (255.999 * img_holo.T * phase_magnitude.T
+ * rgba_img[:, :, :3].T).T.astype(np.uint8)
+ img = Image.fromarray(rgb)
+
+ self.image = img
- r = np.sqrt(xx ** 2 + yy ** 2)
- color_wheel_magnitude = (1 - np.cos(r * pi/360)) / 2
- color_wheel_magnitude *= 0 * (r > 256) + 1 * (r <= 256)
-
- cmap = HOLO_CMAP
- rgba_img = cmap(color_wheel_angle)
- rgb_img = np.delete(rgba_img, 3, 2)
- red, green, blue = rgb_img[:, :, 0], rgb_img[:, :, 1], rgb_img[:, :, 2]
- red *= 255.999 * color_wheel_magnitude
- green *= 255.999 * color_wheel_magnitude
- blue *= 255.999 * color_wheel_magnitude
- #rgb = np.dstack((red, green, blue)).astype(np.uint8)
- # TODO Evtl. einfacher:
- rgb = (255.999 * color_wheel_magnitude.T * rgba_img[:, :, :3].T).T.astype(np.uint8)
- color_wheel = Image.fromarray(rgb)
-
- fig = plt.figure()
- ax = fig.add_subplot(111, aspect='equal')
- ax.imshow(color_wheel)
- ax.set_title('Color Wheel')
- ax.set_xlabel('x-axis')
- ax.set_ylabel('y-axis')
-
-def display_holo(holo, title, axis=None):
- # TODO: Docstring
- if axis == None:
+ @classmethod
+ def make_color_wheel(self):
+ '''Display a color wheel for the gradient direction.
+ Arguments:
+ None
+ Returns:
+ None
+
+ '''
+ x = np.linspace(-256, 256, num=512)
+ y = np.linspace(-256, 256, num=512)
+ xx, yy = np.meshgrid(x, y)
+ color_wheel_angle = (1 - np.arctan2(xx, -yy)/pi) / 2
+
+ r = np.sqrt(xx ** 2 + yy ** 2)
+ color_wheel_magnitude = (1 - np.cos(r * pi/360)) / 2
+ color_wheel_magnitude *= 0 * (r > 256) + 1 * (r <= 256)
+
+ cmap = HoloImage.HOLO_CMAP
+ rgba_img = cmap(color_wheel_angle)
+ rgb_img = np.delete(rgba_img, 3, 2)
+ red, green, blue = rgb_img[:, :, 0], rgb_img[:, :, 1], rgb_img[:, :, 2]
+ red *= 255.999 * color_wheel_magnitude
+ green *= 255.999 * color_wheel_magnitude
+ blue *= 255.999 * color_wheel_magnitude
+ #rgb = np.dstack((red, green, blue)).astype(np.uint8)
+ # TODO Evtl. einfacher:
+ rgb = (255.999 * color_wheel_magnitude.T * rgba_img[:, :, :3].T).T.astype(np.uint8)
+ color_wheel = Image.fromarray(rgb)
+
fig = plt.figure()
- axis = fig.add_subplot(1,1,1, aspect='equal')
+ ax = fig.add_subplot(111, aspect='equal')
+ ax.imshow(color_wheel)
+ ax.set_title('Color Wheel')
+ ax.set_xlabel('x-axis')
+ ax.set_ylabel('y-axis')
- axis.imshow(holo)
- axis.set_title(title)
- axis.set_xlabel('x-axis')
- axis.set_ylabel('y-axis')
\ No newline at end of file
+ def display_holo(holo, title, axis=None):
+ # TODO: Docstring
+ if axis == None:
+ fig = plt.figure()
+ axis = fig.add_subplot(1,1,1, aspect='equal')
+
+ axis.imshow(holo)
+
+ axis.set_title(title)
+ axis.set_xlabel('x-axis')
+ axis.set_ylabel('y-axis')
\ No newline at end of file
diff --git a/pyramid/magcreator.py b/pyramid/magcreator.py
index a497325..7ad98ef 100644
--- a/pyramid/magcreator.py
+++ b/pyramid/magcreator.py
@@ -1,11 +1,12 @@
# -*- coding: utf-8 -*-
"""Create simple LLG Files which describe magnetization in 2D (z-Dim=1)."""
+
import numpy as np
-import matplotlib.pyplot as plt
+from magdata import MagData
-def slab(dim, center, width):
+def shape_slab(dim, center, width):
'''Get the magnetic shape of a slab.
Arguments:
dim - the dimensions of the grid, shape(y, x)
@@ -15,13 +16,14 @@ def slab(dim, center, width):
the magnetic shape as a 2D-boolean-array.
'''
- mag_shape = np.array([[abs(x - center[1]) <= width[1] / 2
- and abs(y - center[0]) <= width[0] / 2
- for x in range(dim[1])] for y in range(dim[0])])
+ mag_shape = np.array([[[abs(x - center[2]) <= width[2] / 2
+ and abs(y - center[1]) <= width[1] / 2
+ and abs(z - center[0]) <= width[0] / 2
+ for x in range(dim[2])] for y in range(dim[1])] for z in range(dim[0])])
return mag_shape
-def disc(dim, center, radius):
+def shape_disc(dim, center, radius, height):
'''Get the magnetic shape of a disc.
Arguments:
dim - the dimensions of the grid, shape(y, x)
@@ -30,13 +32,14 @@ def disc(dim, center, radius):
Returns:
the magnetic shape as a 2D-boolean-array.
- '''
- mag_shape = np.array([[(np.hypot(x-center[1], y-center[0]) <= radius)
- for x in range(dim[1])] for y in range(dim[0])])
+ '''# TODO: up till now only in z-direction
+ mag_shape = np.array([[[np.hypot(x - center[2], y - center[1]) <= radius
+ and abs(z - center[0]) <= height / 2
+ for x in range(dim[2])] for y in range(dim[1])] for z in range(dim[0])])
return mag_shape
-def filament(dim, pos, x_or_y):
+def shape_filament(dim, pos, x_or_y):
'''Get the magnetic shape of a single filament in x or y direction.
Arguments:
dim - the dimensions of the grid, shape(y, x)
@@ -46,16 +49,16 @@ def filament(dim, pos, x_or_y):
Returns:
the magnetic shape as a 2D-boolean-array.
- '''
+ '''# TODO: up till now no z-direction
mag_shape = np.zeros(dim)
if x_or_y == 'y':
- mag_shape[:, pos] = 1
+ mag_shape[:, :, pos] = 1
else:
- mag_shape[pos, :] = 1
+ mag_shape[:, pos, :] = 1
return mag_shape
-def alternating_filaments(dim, spacing, x_or_y):
+def shape_alternating_filaments(dim, spacing, x_or_y):
'''Get the magnetic shape of alternating filaments in x or y direction.
Arguments:
dim - the dimensions of the grid, shape(y, x)
@@ -68,16 +71,15 @@ def alternating_filaments(dim, spacing, x_or_y):
'''
mag_shape = np.zeros(dim)
if x_or_y == 'y':
- # TODO: List comprehension:
for i in range(0, dim[1], spacing):
- mag_shape[:, i] = 1 - 2 * (int(i / spacing) % 2) # 1 or -1
+ mag_shape[:, :, i] = 1 - 2 * (int(i / spacing) % 2) # 1 or -1
else:
for i in range(0, dim[0], spacing):
- mag_shape[i, :] = 1 - 2 * (int(i / spacing) % 2) # 1 or -1
+ mag_shape[:, i, :] = 1 - 2 * (int(i / spacing) % 2) # 1 or -1
return mag_shape
-def single_pixel(dim, pixel):
+def shape_single_pixel(dim, pixel):
'''Get the magnetic shape of a single magnetized pixel.
Arguments:
dim - the dimensions of the grid, shape(y, x)
@@ -87,11 +89,11 @@ def single_pixel(dim, pixel):
'''
mag_shape = np.zeros(dim)
- mag_shape[pixel[0], pixel[1]] = 1
+ mag_shape[pixel] = 1
return mag_shape
-def create_hom_mag(dim, res, beta, mag_shape, filename='output.txt', plot_mag_distr=False):
+def hom_mag(dim, res, mag_shape, beta, magnitude=1):
'''Create homog. magnetization data, saved in a file with LLG-format.
Arguments:
dim - the dimensions of the grid, shape(y, x)
@@ -104,33 +106,27 @@ def create_hom_mag(dim, res, beta, mag_shape, filename='output.txt', plot_mag_di
Returns:
the the magnetic distribution as a 2D-boolean-array.
- '''
- res *= 1.0E-9 / 1.0E-2 # from nm to cm
-
- x = np.linspace(res / 2, dim[1] * res - res / 2, num=dim[1])
- y = np.linspace(res / 2, dim[0] * res - res / 2, num=dim[0])
- xx, yy = np.meshgrid(x, y)
-
- x_mag = np.array(np.ones(dim)) * np.cos(beta) * mag_shape
- y_mag = np.array(np.ones(dim)) * np.sin(beta) * mag_shape
- z_mag = np.array(np.zeros(dim))
-
- if (plot_mag_distr):
- fig = plt.figure()
- fig.add_subplot(111, aspect='equal')
- plt.quiver(x_mag, y_mag, pivot='middle', angles='xy', scale_units='xy',
- scale=1, headwidth=6, headlength=7)
-
- xx = np.reshape(xx,(-1))
- yy = np.reshape(yy,(-1))
- zz = np.array(np.ones(dim[0] * dim[1]) * res / 2)
- x_mag = np.reshape(x_mag,(-1))
- y_mag = np.reshape(y_mag,(-1))
- z_mag = np.array(np.zeros(dim[0] * dim[1]))
+ ''' # TODO: renew Docstring
+ z_mag = np.array(np.zeros(dim))
+ y_mag = np.array(np.ones(dim)) * np.sin(beta) * mag_shape * magnitude
+ x_mag = np.array(np.ones(dim)) * np.cos(beta) * mag_shape * magnitude
+ return z_mag, y_mag, x_mag
- data = np.array([xx, yy, zz, x_mag, y_mag, z_mag]).T
- with open(filename,'w') as mag_file:
- mag_file.write('LLGFileCreator2D: %s\n' % filename.replace('.txt', ''))
- mag_file.write(' %d %d %d\n' % (dim[1], dim[0], 1))
- mag_file.writelines('\n'.join(' '.join('{:7.6e}'.format(cell)
- for cell in row) for row in data) )
\ No newline at end of file
+
+def create_mag_dist(dim, res, mag_shape_list, beta_list, magnitude_list=None):
+ # TODO: Docstring: can take just one or a list of objects! OR JUST A LIST!!!
+ # If no relative magnitude is specified, 1 is assumed for every homog. object:
+ if magnitude_list is None:
+ magnitude_list = np.ones(np.size(beta_list))
+ # For every shape of a homogeneous object a relative magnitude and angle have to be set:
+ assert (np.shape(mag_shape_list)[0],) == np.shape(beta_list) == np.shape(magnitude_list), \
+ 'Lists have not the same length!'
+ x_mag = np.zeros(dim)
+ y_mag = np.zeros(dim)
+ z_mag = np.zeros(dim)
+ for i in range(np.size(beta_list)):
+ pixel_mag = hom_mag(dim, res, mag_shape_list[i], beta_list[i], magnitude_list[i])
+ z_mag += pixel_mag[0]
+ y_mag += pixel_mag[1]
+ x_mag += pixel_mag[2]
+ return MagData(res, z_mag, y_mag, x_mag)
\ No newline at end of file
diff --git a/pyramid/magdata.py b/pyramid/magdata.py
index 3d91b8d..04b883d 100644
--- a/pyramid/magdata.py
+++ b/pyramid/magdata.py
@@ -3,52 +3,136 @@
import numpy as np
+import matplotlib.pyplot as plt
+from mpl_toolkits.mplot3d import Axes3D
+from matplotlib.patches import FancyArrowPatch
+from mpl_toolkits.mplot3d import proj3d
class MagData:
- '''An object storing magnetization data loaded from a LLG-file.'''
+ '''An object storing magnetization data.'''
- def __init__(self, x_mag, y_mag, z_mag, res): # TODO: electrostatic component!
+
+ def __init__(self, res, z_mag, y_mag, x_mag): # TODO: electrostatic component!
'''Load magnetization in LLG-file format.
Arguments:
filename - the name of the file where the data are stored
Returns:
None
- '''
+ '''# TODO: Docstring
assert np.shape(x_mag) == np.shape(y_mag) == np.shape(z_mag), 'Dimensions do not match!'
-
self.res = res
self.dim = np.shape(x_mag)
+ self.magnitude = (z_mag, y_mag, x_mag)
-
-
-
+ @classmethod
+ def load_from_llg(self, filename):
+ # TODO: Docstring
scale = 1.0E-9 / 1.0E-2 #from cm to nm
data = np.genfromtxt(filename, skip_header=2)
x_dim, y_dim, z_dim = np.genfromtxt(filename, dtype=int,
skip_header=1,
skip_footer=len(data[:, 0]))
res = (data[1, 0] - data[0, 0]) / scale
- x_len, y_len, z_len = [data[-1, i]/scale+res/2 for i in range(3)]
x_mag, y_mag, z_mag = [data[:, i].reshape(z_dim, y_dim, x_dim)
for i in range(3,6)]
- #Reshape in Python and Igor is different,
- #Python fills rows first, Igor columns!
- self.filename = filename
- self.res = res
- self.dim = (z_dim, y_dim, x_dim)
- self.length = (z_len, y_len, x_len)
- self.magnitude = (z_mag, y_mag, x_mag)
-
- def __str__(self):
- '''Return the filename of the loaded file.
+ #Reshape in Python and Igor is different, Python fills rows first, Igor columns!
+ return MagData(res, z_mag, y_mag, x_mag)
+
+
+ def save_to_llg(self, filename='output.txt'):
+ '''Create homog. magnetization data, saved in a file with LLG-format.
Arguments:
- None
+ dim - the dimensions of the grid, shape(y, x)
+ res - the resolution of the grid
+ (real space distance between two points)
+ beta - the angle of the magnetization
+ filename - the name of the file in which to save the data
+ mag_shape - an array of shape dim, representing the shape of the magnetic object,
+ a few are supplied in this module
Returns:
- the name of the loaded file as a string
+ the the magnetic distribution as a 2D-boolean-array.
- '''
- return self.filename
\ No newline at end of file
+ ''' # TODO: Renew Docstring
+ dim = self.dim
+ res = self.res * 1.0E-9 / 1.0E-2 # from nm to cm
+
+ zz, yy, xx = np.mgrid[res/2 : (dim[0]*res-res/2) : dim[0]*1j,
+ res/2 : (dim[1]*res-res/2) : dim[1]*1j,
+ res/2 : (dim[2]*res-res/2) : dim[2]*1j]
+
+ xx = np.reshape(xx,(-1))
+ yy = np.reshape(yy,(-1))
+ zz = np.reshape(zz,(-1))
+ x_mag = np.reshape(self.magnitude[2], (-1))
+ y_mag = np.reshape(self.magnitude[1], (-1))
+ z_mag = np.reshape(self.magnitude[0], (-1))
+
+ data = np.array([xx, yy, zz, x_mag, y_mag, z_mag]).T
+ with open(filename,'w') as mag_file:
+ mag_file.write('LLGFileCreator2D: %s\n' % filename.replace('.txt', ''))
+ mag_file.write(' %d %d %d\n' % (dim[2], dim[1], dim[0]))
+ mag_file.writelines('\n'.join(' '.join('{:7.6e}'.format(cell)
+ for cell in row) for row in data) )
+
+
+ def quiver_plot(self, axis='z', ax_slice=0):
+ # TODO: Docstring
+ assert axis == 'z' or axis == 'y' or axis == 'x', 'Axis has to be x, y or z (as string).'
+ if axis == 'z':
+ mag_slice_1 = self.magnitude[2][ax_slice,...]
+ mag_slice_2 = self.magnitude[1][ax_slice,...]
+ elif axis == 'y':
+ mag_slice_1 = self.magnitude[2][:,ax_slice,:]
+ mag_slice_2 = self.magnitude[0][:,ax_slice,:]
+ elif axis == 'x':
+ mag_slice_1 = self.magnitude[1][...,ax_slice]
+ mag_slice_2 = self.magnitude[0][...,ax_slice]
+
+ fig = plt.figure()
+ fig.add_subplot(111, aspect='equal')
+ plt.quiver(mag_slice_1, mag_slice_2, pivot='middle', angles='xy', scale_units='xy',
+ scale=1, headwidth=6, headlength=7)
+
+
+ def quiver_plot_3D(self):
+ # TODO: Docstring
+ res = self.res
+ dim = self.dim
+ fig = plt.figure()
+ ax = fig.add_subplot(111, projection='3d')
+ ax.set_aspect("equal")
+
+ class Arrow3D(FancyArrowPatch):
+ def __init__(self, xs, ys, zs, *args, **kwargs):
+ FancyArrowPatch.__init__(self, (0,0), (0,0), *args, **kwargs)
+ self._verts3d = xs, ys, zs
+ def draw(self, renderer):
+ xs3d, ys3d, zs3d = self._verts3d
+ xs, ys, zs = proj3d.proj_transform(xs3d, ys3d, zs3d, renderer.M)
+ self.set_positions((xs[0],ys[0]),(xs[1],ys[1]))
+ FancyArrowPatch.draw(self, renderer)
+
+ zz, yy, xx = np.mgrid[res/2 : (dim[0]*res-res/2) : dim[0]*1j,
+ res/2 : (dim[1]*res-res/2) : dim[1]*1j,
+ res/2 : (dim[2]*res-res/2) : dim[2]*1j]
+ xx = np.reshape(xx, (-1))
+ yy = np.reshape(yy, (-1))
+ zz = np.reshape(zz, (-1))
+ x_mag = np.reshape(self.magnitude[2], (-1))
+ y_mag = np.reshape(self.magnitude[1], (-1))
+ z_mag = np.reshape(self.magnitude[0], (-1))
+
+ for i in range(np.size(xx)):
+ xs = [xx[i] - x_mag[i]*res/2, xx[i] + x_mag[i]*res/2]
+ ys = [yy[i] - y_mag[i]*res/2, yy[i] + y_mag[i]*res/2]
+ zs = [zz[i] - z_mag[i]*res/2, zz[i] + z_mag[i]*res/2]
+ a = Arrow3D(xs, ys, zs, mutation_scale=10, lw=1, arrowstyle="-|>", color="k")
+ ax.add_artist(a)
+ ax.set_xlim3d(0, xx[-1]+res/2)
+ ax.set_ylim3d(0, yy[-1]+res/2)
+ ax.set_zlim3d(0, zz[-1]+res/2)
+ plt.show()
\ No newline at end of file
diff --git a/pyramid/numcore/setup.py b/pyramid/numcore/setup.py
deleted file mode 100644
index b31924c..0000000
--- a/pyramid/numcore/setup.py
+++ /dev/null
@@ -1,30 +0,0 @@
-"""
-Created on Fri May 03 10:27:04 2013
-
-@author: Jan
-"""
-
-#call with: python setup.py build_ext --inplace --compiler=mingw32
-
-import glob
-import numpy
-from distutils.core import setup
-from distutils.extension import Extension
-from Cython.Build import cythonize
-from Cython.Distutils import build_ext
-
-#setup(
-# cmdclass = {'build_ext': build_ext},
-# ext_modules = [Extension("multiply",
-# sources=["multiply.pyx", "c_multiply.c"],
-# include_dirs=[numpy.get_include()])],
-#)
-
-setup(
- name = 'Pyramex',
- version = '0.1',
- description = 'Pyramid Cython Extensions',
- author = 'Jan Caron',
- author_email = 'j.caron@fz-juelich.de',
- ext_modules = cythonize(glob.glob('*.pyx'))
-)
diff --git a/pyramid/phasemap.py b/pyramid/phasemap.py
index 60c639d..84e4efc 100644
--- a/pyramid/phasemap.py
+++ b/pyramid/phasemap.py
@@ -1,209 +1,105 @@
# -*- coding: utf-8 -*-
-"""Create and display a phase map from magnetization data."""
+"""
+Created on Mon May 13 16:00:57 2013
+@author: Jan
+"""
-import numpy as np
-import matplotlib.pyplot as plt
-from numpy import pi
-
-
-PHI_0 = 2067.83 # magnetic flux in T*nm²
-H_BAR = 6.626E-34 # TODO: unit
-M_E = 9.109E-31 # TODO: unit
-Q_E = 1.602E-19 # TODO: unit
-C = 2.998E8 # TODO: unit
-
-
-def fourier_space(mag_data, b_0=1, padding=0):
- '''Calculate phasemap from magnetization data (fourier approach).
- Arguments:
- mag_data - MagDataLLG object (from magneticimaging.dataloader) storing
- the filename, coordinates and magnetization in x, y and z
- b_0 - magnetic induction corresponding to a magnetization Mo in T
- (default: 1)
- padding - factor for zero padding, the default is 0 (no padding), for
- a factor of n the number of pixels is increase by (1+n)**2,
- padded zeros are cropped at the end
- v_0 - average potential of the sample in V (default: 0)
- v_acc - acceleration voltage of the microscop in V (default: 30000)
- Returns:
- the phasemap as a 2 dimensional array
-
- '''
- res = mag_data.res
- z_dim, y_dim, x_dim = mag_data.dim
- z_mag, y_mag, x_mag = mag_data.magnitude
+class PhaseMap:
- # TODO: interpolation (redefine xDim,yDim,zDim) thickness ramp
+ '''An object storing magnetization data.'''
- x_pad = x_dim/2 * padding
- y_pad = y_dim/2 * padding
- x_mag_big = np.zeros(((1 + padding) * y_dim, (1 + padding) * x_dim))
- y_mag_big = np.zeros(((1 + padding) * y_dim, (1 + padding) * x_dim))
- # TODO: padding so that x_dim and y_dim = 2^n
- x_mag_big[y_pad : y_pad + y_dim, x_pad : x_pad + x_dim] = x_mag
- y_mag_big[y_pad : y_pad + y_dim, x_pad : x_pad + x_dim] = y_mag
- x_mag_fft = np.fft.fftshift(np.fft.rfft2(x_mag_big), axes=0)
- y_mag_fft = np.fft.fftshift(np.fft.rfft2(y_mag_big), axes=0)
- nyq = 1 / res # nyquist frequency
- x = np.linspace( 0, nyq/2, x_mag_fft.shape[1])
- y = np.linspace( -nyq/2, nyq/2, x_mag_fft.shape[0]+1)[:-1]
- xx, yy = np.meshgrid(x, y)
-
- phase_fft = (1j * res * b_0) / (2 * PHI_0) * ((x_mag_fft * yy - y_mag_fft * xx)
- / (xx ** 2 + yy ** 2 + 1e-18))
- phase_big = np.fft.irfft2(np.fft.ifftshift(phase_fft, axes=0))
-
- phase = phase_big[y_pad : y_pad + y_dim, x_pad : x_pad + x_dim]
-
- # TODO: Electrostatic Component
-
- return phase
-
-
-def phi_pixel(method, n, m, res, b_0): # TODO: rename
- 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
- coeff = -b_0 * res**2 / ( 2 * PHI_0 )
- 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':
- coeff = - b_0 * res**2 / ( 2 * PHI_0 )
- 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
-
-
-def real_space(mag_data, method, b_0=1, jacobi=None):
- '''Calculate phasemap from magnetization data (real space approach).
- Arguments:
- mag_data - MagDataLLG object (from magneticimaging.dataloader) storing the filename,
- coordinates and magnetization in x, y and z
- Returns:
- the phasemap as a 2 dimensional array
+ def __init__(self, phase): # TODO: more arguments?
+ '''Load magnetization in LLG-file format.
+ Arguments:
+ filename - the name of the file where the data are stored
+ Returns:
+ None
- '''
- # TODO: Expand docstring!
-
- res = mag_data.res
- z_dim, y_dim, x_dim = mag_data.dim
- z_mag, y_mag, x_mag = mag_data.magnitude
-
- # TODO: proper projection algorithm
- x_mag, y_mag, z_mag = x_mag.mean(0), y_mag.mean(0), z_mag.mean(0)
-
- beta = np.arctan2(y_mag, x_mag)
+ '''# TODO: Docstring
+ self.phase = phase
+
+
+ def load_from_file(self, filename): # TODO: Implement
+ pass
+# # TODO: Docstring
+# scale = 1.0E-9 / 1.0E-2 #from cm to nm
+# data = np.genfromtxt(filename, skip_header=2)
+# x_dim, y_dim, z_dim = np.genfromtxt(filename, dtype=int,
+# skip_header=1,
+# skip_footer=len(data[:, 0]))
+# res = (data[1, 0] - data[0, 0]) / scale
+# x_mag, y_mag, z_mag = [data[:, i].reshape(z_dim, y_dim, x_dim)
+# for i in range(3,6)]
+# #Reshape in Python and Igor is different, Python fills rows first, Igor columns!
+# return MagData(res, z_mag, y_mag, x_mag)
- mag = np.hypot(x_mag, y_mag)
-
- '''CREATE COORDINATE GRIDS'''
- x = np.linspace(0,(x_dim-1),num=x_dim)
- y = np.linspace(0,(y_dim-1),num=y_dim)
- xx, yy = np.meshgrid(x,y)
-
- x_big = np.linspace(-(x_dim-1), x_dim-1, num=2*x_dim-1)
- y_big = np.linspace(-(y_dim-1), y_dim-1, num=2*y_dim-1)
- xx_big, yy_big = np.meshgrid(x_big, y_big)
-
- phi_cos = phi_pixel(method, xx_big, yy_big, res, b_0)
- phi_sin = phi_pixel(method, yy_big, xx_big, res, b_0)
-
- def phi_mag(i, j): # TODO: rename
- return (np.cos(beta[j,i])*phi_cos[y_dim-1-j:(2*y_dim-1)-j, x_dim-1-i:(2*x_dim-1)-i]
- -np.sin(beta[j,i])*phi_sin[y_dim-1-j:(2*y_dim-1)-j, x_dim-1-i:(2*x_dim-1)-i])
-
- def phi_mag_deriv(i, j): # TODO: rename
- return -(np.sin(beta[j,i])*phi_cos[y_dim-1-j:(2*y_dim-1)-j, x_dim-1-i:(2*x_dim-1)-i]
- +np.cos(beta[j,i])*phi_sin[y_dim-1-j:(2*y_dim-1)-j, x_dim-1-i:(2*x_dim-1)-i])
-
- def phi_mag_fd(i, j, h): # TODO: rename
- return ((np.cos(beta[j,i]+h) - np.cos(beta[j,i])) / h
- * phi_cos[y_dim-1-j:(2*y_dim-1)-j, x_dim-1-i:(2*x_dim-1)-i]
- -(np.sin(beta[j,i]+h) - np.sin(beta[j,i])) / h
- * phi_sin[y_dim-1-j:(2*y_dim-1)-j, x_dim-1-i:(2*x_dim-1)-i])
-
- '''CALCULATE THE PHASE'''
- phase = np.zeros((y_dim, x_dim))
-
- # TODO: only iterate over pixels that have a magn. > threshold (first >0)
- if jacobi is not None:
- jacobi_fd = jacobi.copy()
- h = 0.0001
-
- for j in range(y_dim):
- for i in range(x_dim):
- #if (mag[j,i] != 0 ):#or jacobi is not None): # TODO: same result with or without?
- phi_mag_cache = phi_mag(i, j)
- phase += mag[j,i] * phi_mag_cache
- if jacobi is not None:
- jacobi[:,i+x_dim*j] = phi_mag_cache.reshape(-1)
- jacobi[:,x_dim*y_dim+i+x_dim*j] = (mag[j,i]*phi_mag_deriv(i,j)).reshape(-1)
-
- jacobi_fd[:,i+x_dim*j] = phi_mag_cache.reshape(-1)
- jacobi_fd[:,x_dim*y_dim+i+x_dim*j] = (mag[j,i]*phi_mag_fd(i,j,h)).reshape(-1)
-
-
-
- if jacobi is not None:
- jacobi_diff = jacobi_fd - jacobi
- assert (np.abs(jacobi_diff) < 1.0E-8).all(), 'jacobi matrix is not the same'
-
- return phase
-
-def phase_elec(mag_data, b_0=1, v_0=0, v_acc=30000):
- # TODO: Docstring
+ def save_to_file(self, filename='output.txt'): # TODO: Implement
+ pass
+# '''Create homog. magnetization data, saved in a file with LLG-format.
+# Arguments:
+# dim - the dimensions of the grid, shape(y, x)
+# res - the resolution of the grid
+# (real space distance between two points)
+# beta - the angle of the magnetization
+# filename - the name of the file in which to save the data
+# mag_shape - an array of shape dim, representing the shape of the magnetic object,
+# a few are supplied in this module
+# Returns:
+# the the magnetic distribution as a 2D-boolean-array.
+#
+# ''' # TODO: Renew Docstring
+# dim = self.dim
+# res = self.res * 1.0E-9 / 1.0E-2 # from nm to cm
+#
+# zz, yy, xx = np.mgrid[res/2 : (dim[0]*res-res/2) : dim[0]*1j,
+# res/2 : (dim[1]*res-res/2) : dim[1]*1j,
+# res/2 : (dim[2]*res-res/2) : dim[2]*1j]
+#
+# xx = np.reshape(xx,(-1))
+# yy = np.reshape(yy,(-1))
+# zz = np.reshape(zz,(-1))
+# x_mag = np.reshape(self.magnitude[2], (-1))
+# y_mag = np.reshape(self.magnitude[1], (-1))
+# z_mag = np.reshape(self.magnitude[0], (-1))
+#
+# data = np.array([xx, yy, zz, x_mag, y_mag, z_mag]).T
+# with open(filename,'w') as mag_file:
+# mag_file.write('LLGFileCreator2D: %s\n' % filename.replace('.txt', ''))
+# mag_file.write(' %d %d %d\n' % (dim[2], dim[1], dim[0]))
+# mag_file.writelines('\n'.join(' '.join('{:7.6e}'.format(cell)
+# for cell in row) for row in data) )
- # TODO: Delete
-# import pdb; pdb.set_trace()
-
- res = mag_data.res
- z_dim, y_dim, x_dim = mag_data.dim
- z_mag, y_mag, x_mag = mag_data.magnitude
-
- phase = np.logical_or(x_mag, y_mag, z_mag)
-
- lam = H_BAR / np.sqrt(2 * M_E * v_acc * (1 + Q_E*v_acc / (2*M_E*C**2)))
-
- Ce = (2*pi*Q_E/lam * (Q_E*v_acc + M_E*C**2) /
- (Q_E*v_acc * (Q_E*v_acc + 2*M_E*C**2)))
-
- phase *= res * v_0 * Ce
-
- return phase
-
-
-def display_phase(phase, res, title, axis=None):
- '''Display the phasemap as a colormesh.
- Arguments:
- phase - the phasemap that should be displayed
- res - the resolution of the phasemap
- title - the title of the plot
- Returns:
- None
+ def display_phase(self, res, title, axis=None):
+ '''Display the phasemap as a colormesh.
+ Arguments:
+ phase - the phasemap that should be displayed
+ res - the resolution of the phasemap
+ title - the title of the plot
+ Returns:
+ None
+
+ '''
+ if axis == None:
+ fig = plt.figure()
+ axis = fig.add_subplot(1,1,1, aspect='equal')
- '''
- if axis == None:
- fig = plt.figure()
- axis = fig.add_subplot(1,1,1, aspect='equal')
+ im = plt.pcolormesh(phase, cmap='gray')
- im = plt.pcolormesh(phase, cmap='gray')
-
- ticks = axis.get_xticks()*res
- axis.set_xticklabels(ticks)
- ticks = axis.get_yticks()*res
- axis.set_yticklabels(ticks)
-
- axis.set_title(title)
- axis.set_xlabel('x-axis [nm]')
- axis.set_ylabel('y-axis [nm]')
-
- fig = plt.gcf()
- fig.subplots_adjust(right=0.85)
- cbar_ax = fig.add_axes([0.9, 0.15, 0.02, 0.7])
- fig.colorbar(im, cax=cbar_ax)
+ ticks = axis.get_xticks()*res
+ axis.set_xticklabels(ticks)
+ ticks = axis.get_yticks()*res
+ axis.set_yticklabels(ticks)
- plt.show()
\ No newline at end of file
+ axis.set_title(title)
+ axis.set_xlabel('x-axis [nm]')
+ axis.set_ylabel('y-axis [nm]')
+
+ fig = plt.gcf()
+ fig.subplots_adjust(right=0.85)
+ cbar_ax = fig.add_axes([0.9, 0.15, 0.02, 0.7])
+ fig.colorbar(im, cax=cbar_ax)
+
+ plt.show()
\ No newline at end of file
diff --git a/pyramid/phasemapper.py b/pyramid/phasemapper.py
new file mode 100644
index 0000000..53f7888
--- /dev/null
+++ b/pyramid/phasemapper.py
@@ -0,0 +1,177 @@
+# -*- coding: utf-8 -*-
+"""Create and display a phase map from magnetization data."""
+
+
+import numpy as np
+import matplotlib.pyplot as plt
+import pyramid.projector as pj
+from numpy import pi
+from phasemap import PhaseMap
+
+
+PHI_0 = 2067.83 # magnetic flux in T*nm²
+H_BAR = 6.626E-34 # TODO: unit
+M_E = 9.109E-31 # TODO: unit
+Q_E = 1.602E-19 # TODO: unit
+C = 2.998E8 # TODO: unit
+
+
+def fourier_space(mag_data, b_0=1, padding=0):
+ '''Calculate phasemap from magnetization data (fourier approach).
+ Arguments:
+ mag_data - MagDataLLG object (from magneticimaging.dataloader) storing
+ the filename, coordinates and magnetization in x, y and z
+ b_0 - magnetic induction corresponding to a magnetization Mo in T
+ (default: 1)
+ padding - factor for zero padding, the default is 0 (no padding), for
+ a factor of n the number of pixels is increase by (1+n)**2,
+ padded zeros are cropped at the end
+ v_0 - average potential of the sample in V (default: 0)
+ v_acc - acceleration voltage of the microscop in V (default: 30000)
+ Returns:
+ the phasemap as a 2 dimensional array
+
+ '''
+ res = mag_data.res
+ z_dim, y_dim, x_dim = mag_data.dim
+ z_mag, y_mag, x_mag = pj.simple_axis_projection(mag_data)#mag_data.magnitude
+
+ # TODO: interpolation (redefine xDim,yDim,zDim) thickness ramp
+
+ x_pad = x_dim/2 * padding
+ y_pad = y_dim/2 * padding
+ x_mag_big = np.zeros(((1 + padding) * y_dim, (1 + padding) * x_dim))
+ y_mag_big = np.zeros(((1 + padding) * y_dim, (1 + padding) * x_dim))
+ # TODO: padding so that x_dim and y_dim = 2^n
+ x_mag_big[y_pad : y_pad + y_dim, x_pad : x_pad + x_dim] = x_mag
+ y_mag_big[y_pad : y_pad + y_dim, x_pad : x_pad + x_dim] = y_mag
+
+ x_mag_fft = np.fft.fftshift(np.fft.rfft2(x_mag_big), axes=0)
+ y_mag_fft = np.fft.fftshift(np.fft.rfft2(y_mag_big), axes=0)
+ nyq = 1 / res # nyquist frequency
+ x = np.linspace( 0, nyq/2, x_mag_fft.shape[1])
+ y = np.linspace( -nyq/2, nyq/2, x_mag_fft.shape[0]+1)[:-1]
+ xx, yy = np.meshgrid(x, y)
+
+ phase_fft = (1j * res * b_0) / (2 * PHI_0) * ((x_mag_fft * yy - y_mag_fft * xx)
+ / (xx ** 2 + yy ** 2 + 1e-18))
+ phase_big = np.fft.irfft2(np.fft.ifftshift(phase_fft, axes=0))
+
+ phase = phase_big[y_pad : y_pad + y_dim, x_pad : x_pad + x_dim]
+
+ # TODO: Electrostatic Component
+
+ return PhaseMap(phase)
+
+
+def phi_pixel(method, n, m, res, b_0): # TODO: rename
+ 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
+ coeff = -b_0 * res**2 / ( 2 * PHI_0 )
+ 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':
+ coeff = - b_0 * res**2 / ( 2 * PHI_0 )
+ 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
+
+
+def real_space(mag_data, method, b_0=1, jacobi=None):
+ '''Calculate phasemap from magnetization data (real space approach).
+ Arguments:
+ mag_data - MagDataLLG object (from magneticimaging.dataloader) storing the filename,
+ coordinates and magnetization in x, y and z
+ Returns:
+ the phasemap as a 2 dimensional array
+
+ '''
+ # TODO: Expand docstring!
+
+ res = mag_data.res
+ z_dim, y_dim, x_dim = mag_data.dim
+ z_mag, y_mag, x_mag = mag_data.magnitude
+
+ z_mag, y_mag, x_mag = pj.simple_axis_projection(mag_data)
+
+ beta = np.arctan2(y_mag, x_mag)
+
+ mag = np.hypot(x_mag, y_mag)
+
+ '''CREATE COORDINATE GRIDS'''
+ x = np.linspace(0,(x_dim-1),num=x_dim)
+ y = np.linspace(0,(y_dim-1),num=y_dim)
+ xx, yy = np.meshgrid(x,y)
+
+ x_big = np.linspace(-(x_dim-1), x_dim-1, num=2*x_dim-1)
+ y_big = np.linspace(-(y_dim-1), y_dim-1, num=2*y_dim-1)
+ xx_big, yy_big = np.meshgrid(x_big, y_big)
+
+ phi_cos = phi_pixel(method, xx_big, yy_big, res, b_0)
+ phi_sin = phi_pixel(method, yy_big, xx_big, res, b_0)
+
+ def phi_mag(i, j): # TODO: rename
+ return (np.cos(beta[j,i])*phi_cos[y_dim-1-j:(2*y_dim-1)-j, x_dim-1-i:(2*x_dim-1)-i]
+ -np.sin(beta[j,i])*phi_sin[y_dim-1-j:(2*y_dim-1)-j, x_dim-1-i:(2*x_dim-1)-i])
+
+ def phi_mag_deriv(i, j): # TODO: rename
+ return -(np.sin(beta[j,i])*phi_cos[y_dim-1-j:(2*y_dim-1)-j, x_dim-1-i:(2*x_dim-1)-i]
+ +np.cos(beta[j,i])*phi_sin[y_dim-1-j:(2*y_dim-1)-j, x_dim-1-i:(2*x_dim-1)-i])
+
+ def phi_mag_fd(i, j, h): # TODO: rename
+ return ((np.cos(beta[j,i]+h) - np.cos(beta[j,i])) / h
+ * phi_cos[y_dim-1-j:(2*y_dim-1)-j, x_dim-1-i:(2*x_dim-1)-i]
+ -(np.sin(beta[j,i]+h) - np.sin(beta[j,i])) / h
+ * phi_sin[y_dim-1-j:(2*y_dim-1)-j, x_dim-1-i:(2*x_dim-1)-i])
+
+ '''CALCULATE THE PHASE'''
+ phase = np.zeros((y_dim, x_dim))
+
+ # TODO: only iterate over pixels that have a magn. > threshold (first >0)
+ if jacobi is not None:
+ jacobi_fd = jacobi.copy()
+ h = 0.0001
+
+ for j in range(y_dim):
+ for i in range(x_dim):
+ #if (mag[j,i] != 0 ):#or jacobi is not None): # TODO: same result with or without?
+ phi_mag_cache = phi_mag(i, j)
+ phase += mag[j,i] * phi_mag_cache
+ if jacobi is not None:
+ jacobi[:,i+x_dim*j] = phi_mag_cache.reshape(-1)
+ jacobi[:,x_dim*y_dim+i+x_dim*j] = (mag[j,i]*phi_mag_deriv(i,j)).reshape(-1)
+
+ jacobi_fd[:,i+x_dim*j] = phi_mag_cache.reshape(-1)
+ jacobi_fd[:,x_dim*y_dim+i+x_dim*j] = (mag[j,i]*phi_mag_fd(i,j,h)).reshape(-1)
+
+
+
+ if jacobi is not None:
+ jacobi_diff = jacobi_fd - jacobi
+ assert (np.abs(jacobi_diff) < 1.0E-8).all(), 'jacobi matrix is not the same'
+
+ return phase
+
+
+def phase_elec(mag_data, b_0=1, v_0=0, v_acc=30000):
+ # TODO: Docstring
+
+ # TODO: Delete
+# import pdb; pdb.set_trace()
+
+ res = mag_data.res
+ z_dim, y_dim, x_dim = mag_data.dim
+ z_mag, y_mag, x_mag = mag_data.magnitude
+
+ phase = np.logical_or(x_mag, y_mag, z_mag)
+
+ lam = H_BAR / np.sqrt(2 * M_E * v_acc * (1 + Q_E*v_acc / (2*M_E*C**2)))
+
+ Ce = (2*pi*Q_E/lam * (Q_E*v_acc + M_E*C**2) /
+ (Q_E*v_acc * (Q_E*v_acc + 2*M_E*C**2)))
+
+ phase *= res * v_0 * Ce
+
+ return phase
\ No newline at end of file
diff --git a/pyramid/projector.py b/pyramid/projector.py
new file mode 100644
index 0000000..c101733
--- /dev/null
+++ b/pyramid/projector.py
@@ -0,0 +1,17 @@
+# -*- coding: utf-8 -*-
+"""
+Created on Tue May 14 14:44:36 2013
+
+@author: Jan
+"""# TODO: Docstring
+
+from magdata import MagData
+
+def simple_axis_projection(mag_data, axis=0):
+ # TODO: assert isinstance(mag_data, MagData), 'Object is no instance of MagData!'
+ return (mag_data.magnitude[0].mean(axis), # x_mag
+ mag_data.magnitude[1].mean(axis), # y_mag
+ mag_data.magnitude[2].mean(axis)) # z_mag
+
+
+# TODO: proper projection algorithm with two angles and such!
\ No newline at end of file
diff --git a/pyramid/reconstructor.py b/pyramid/reconstructor.py
new file mode 100644
index 0000000..0e3835a
--- /dev/null
+++ b/pyramid/reconstructor.py
@@ -0,0 +1,7 @@
+# -*- coding: utf-8 -*-
+"""
+Created on Tue May 14 15:19:33 2013
+
+@author: Jan
+"""
+
diff --git a/scripts/create_random_distribution.py b/scripts/create_random_distribution.py
index 1ac5c6b..fb7cd80 100644
--- a/scripts/create_random_distribution.py
+++ b/scripts/create_random_distribution.py
@@ -5,30 +5,40 @@ Created on Mon May 13 13:05:40 2013
@author: Jan
"""
-import random
+import random as rnd
import numpy as np
-import matplotlib.pyplot as plt
import pyramid.magcreator as mc
-import time
import pdb, traceback, sys
from numpy import pi
-def phase_from_mag():
+def create_random_distribution():
count = 10
- dim = (128, 128)
+ dim = (1, 128, 128)
+ res = 10 # in nm
- random.seed(42)
+ rnd.seed(42)
+
+ mag_shape_list = np.zeros((count, dim[0], dim[1], dim[2]))
+ beta_list = np.zeros(count)
+ magnitude_list = np.zeros(count)
for i in range(count):
+ pixel = (rnd.randrange(dim[0]), rnd.randrange(dim[1]), rnd.randrange(dim[2]))
+ mag_shape_list[i,...] = mc.shape_single_pixel(dim, pixel)
+ beta_list[i] = 2*pi*rnd.random()
+ magnitude_list[i] = 1#rnd.random()
- x = random.rand
+ mag_data = mc.create_mag_dist(dim, res, mag_shape_list, beta_list, magnitude_list)
+ mag_data.quiver_plot()
+ #mag_data.quiver_plot_3D()
+ mag_data.save_to_llg('../output/mag_dist_random_pixel.txt')
if __name__ == "__main__":
try:
- phase_from_mag()
+ create_random_distribution()
except:
type, value, tb = sys.exc_info()
traceback.print_exc()
--
GitLab