Skip to content
Snippets Groups Projects
Commit 0bfdb896 authored by Jan Caron's avatar Jan Caron
Browse files

Middle Stage of converting the program structure (not working correctly9 

For Pushing Purposes
parent ca664902
No related branches found
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
Showing
with 564 additions and 446 deletions
pyramid/dataloader.py deleted 100644 → 0
+ 0
− 45
View file @ ca664902
# -*- 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
pyramid/holoimage.py
+ 102
− 96
View file @ 0bfdb896
......@@ -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
pyramid/magcreator.py
+ 45
− 49
View file @ 0bfdb896
# -*- 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
pyramid/magdata.py
+ 106
− 22
View file @ 0bfdb896
......@@ -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
pyramid/numcore/setup.py deleted 100644 → 0
+ 0
− 30
View file @ ca664902
"""
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'))
)
pyramid/phasemap.py
+ 92
− 196
View file @ 0bfdb896
# -*- 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
pyramid/phasemapper.py 0 → 100644
+ 177
− 0
View file @ 0bfdb896
# -*- 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
pyramid/projector.py 0 → 100644
+ 17
− 0
View file @ 0bfdb896
# -*- 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
pyramid/reconstructor.py 0 → 100644
+ 7
− 0
View file @ 0bfdb896
# -*- coding: utf-8 -*-
"""
Created on Tue May 14 15:19:33 2013
@author: Jan
"""
scripts/create_random_distribution.py
+ 18
− 8
View file @ 0bfdb896
......@@ -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()
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment