diff --git a/2d_fourie/lattices.py b/2d_fourie/lattices.py index f538676..2d09b71 100644 --- a/2d_fourie/lattices.py +++ b/2d_fourie/lattices.py @@ -1,6 +1,7 @@ import numpy as np from cache import timeit + def deg_2_rad(winkel): return winkel / 180.0 * np.pi @@ -42,7 +43,9 @@ class VO2_Lattice(Lattice): base_c_m = 5.38 base_c_r = 2.856 + base_c_r = 2.8 base_b_r = 4.554 + base_b_r = 4.5 base_a_r = base_b_r alpha_m = 122.64 # degree @@ -64,6 +67,14 @@ class VO2_Lattice(Lattice): offset_a_r, offset_c_r = self._mono_2_rutile(offset_c_m, offset_a_m) + res = 0.05 + offset_a_r = res * int(offset_a_r/res) + offset_c_r = res * int(offset_c_r/res) + + offset_a_r = 0.5 + offset_c_r = 0.5 + print(offset_a_r, offset_c_r) + self.atom_x_mono = offset_a_r + X * \ self.base_c_r + np.mod(Y, 4) * 0.5 * self.base_c_r self.atom_x_mono[np.mod(X, 2) == 0] -= 2 * offset_a_r diff --git a/2d_fourie/main.py b/2d_fourie/main.py index 462980a..5e46d4c 100644 --- a/2d_fourie/main.py +++ b/2d_fourie/main.py @@ -8,6 +8,7 @@ import tqdm from extractors import Rect_Evaluator, Voronoi_Evaluator, Image_Wrapper from cache import timeit + class Plotter: def __init__(self, lat): self.lattice = lat @@ -61,7 +62,7 @@ def rotate(x, y, angle): def test_square(): - LEN = 10 + LEN = 40 lat = SCC_Lattice(LEN, LEN) # lat = VO2_Lattice(LEN, LEN) plot = Plotter(lat) @@ -69,16 +70,26 @@ def test_square(): # pos_x, pos_y = rotate(pos_x, pos_y, 30) si = SpinImage(pos_x, pos_y) fig, axs = plt.subplots(2, 2) - si.plot(axs[0, 0],1) + si.plot(axs[0, 0], 1) # si.gaussian(LEN) # si.blur(3) - si.pad_it_square(int((5./2.)/0.1)) - si.plot(axs[0, 1],1) + si.pad_it_square(size=int((5*LEN)/si.resolution)) + si.plot(axs[0, 1], 1) plt.pause(0.1) + fx, fy, intens = si.fft() plot.plot(fx, fy, intens, axs[1, 0], axs[1, 1]) + plt.tight_layout() plt.savefig("test.pdf") + + plt.figure() + index = np.abs(fy).argmin() + plt.plot(fx, intens[index, :], label="fy=0") + plt.plot(fx, intens[index+1, :], label="fy=df") + plt.plot(fx, intens[index+2, :], label="fy=2df") + plt.legend() + plt.yscale("log") plt.show() @@ -97,87 +108,67 @@ def helper(intens, fx, fy, voro, rect): def test_mixed(): + shift_x = -(5 + 3*28) + shift_y = -43 LEN = 40 lat = VO2_Lattice(LEN, LEN) plot = Plotter(lat) - all_rutile = np.stack(lat.reci()[0]).T - all_mono = np.stack(lat.reci()[1]).T - all_mono2 = np.stack(lat.reci()[2]).T - voro = Voronoi_Evaluator([all_rutile, all_mono, all_mono2]) - rect = Rect_Evaluator([all_rutile, all_mono, all_mono2]) - #voro = Voronoi_Evaluator([all_mono]) - #rect = Rect_Evaluator([all_mono]) - pos_x, pos_y = lat.get_from_mask(np.zeros((40, 40))) - si = SpinImage(pos_x, pos_y) - si.pad_it_square(10, size=2300) - fx, fy, intens_mono = si.fft() + pos_x, pos_y = lat.get_from_mask(np.zeros((LEN, LEN))) + si_mono = SpinImage(pos_x, pos_y) + si_mono.pad_it(int((LEN * 5.712) / si_mono.resolution)+shift_x, + int((LEN * 4.554) / si_mono.resolution)+shift_y) + fx, fy, intens_mono = si_mono.fft() - pos_x, pos_y = lat.get_from_mask(np.ones((40, 40))) - si = SpinImage(pos_x, pos_y) - si.pad_it_square(10, size=2300) - fx, fy, intens_rutile = si.fft() + pos_x, pos_y = lat.get_from_mask(np.ones((LEN, LEN))) + si_rutile = SpinImage(pos_x, pos_y) + print(int((LEN * 5.712) / si_rutile.resolution)+shift_x, + int((LEN * 4.554) / si_rutile.resolution)+shift_y, + si_rutile.img.shape) + si_rutile.pad_it(int((LEN * 5.712) / si_rutile.resolution)+shift_x, + int((LEN * 4.554) / si_rutile.resolution)+shift_y) + fx, fy, intens_rutile = si_rutile.fft() - mask_misk = np.ones((40, 40)) + mask_misk = np.ones((LEN, LEN)) ind = np.arange(mask_misk.size) np.random.shuffle(ind) - mask_misk[np.unravel_index(ind[:800], (40, 40))] = False + mask_misk[np.unravel_index(ind[:800], (LEN, LEN))] = False pos_x, pos_y = lat.get_from_mask(mask_misk) - si = SpinImage(pos_x, pos_y) - si.pad_it_square(10, size=2300) - fx, fy, intens_mixed = si.fft() - img = Image_Wrapper(intens_mono, fx, fy) - voro.generate_mask(img) - rect.generate_mask(Image_Wrapper(intens_mono, fx, fy)) + si_mixed = SpinImage(pos_x, pos_y) + si_mixed.pad_it(int((LEN * 5.712) / si_mixed.resolution)+shift_x, + int((LEN * 4.554) / si_mixed.resolution)+shift_y) + fx, fy, intens_mixed = si_mixed.fft() - fig, axs = plt.subplots(1, 2) - axs[0].imshow(rect.get_mask(), extent=img.ext(), origin="lower") - axs[0].plot(all_rutile[:, 0], all_rutile[:, 1], ".") - axs[1].plot(all_rutile[:, 0], all_rutile[:, 1], ".") - axs[0].plot(all_mono[:, 0], all_mono[:, 1], ".") - axs[1].plot(all_mono[:, 0], all_mono[:, 1], ".") - axs[0].plot(all_mono2[:, 0], all_mono2[:, 1], ".") - axs[1].plot(all_mono2[:, 0], all_mono2[:, 1], ".") - axs[1].imshow(voro.get_mask(), extent=img.ext(), origin="lower") + fig, axs = plt.subplots(3, 3) + si_rutile.plot(axs[0, 0], 1) + si_mono.plot(axs[0, 2], 1) + si_mixed.plot(axs[0, 1], 1) + plot.plot(freqx=fx, freqy=fy, intens=intens_rutile, + ax_log=axs[1, 0], ax_lin=axs[2, 0], vmax=10e7) + plot.plot(freqx=fx, freqy=fy, intens=intens_mono, + ax_log=axs[1, 2], ax_lin=axs[2, 2], vmax=10e7) + plot.plot(freqx=fx, freqy=fy, intens=intens_mixed, + ax_log=axs[1, 1], ax_lin=axs[2, 1], vmax=10e7) - print("mono") - helper(intens=intens_mono, fx=fx, fy=fy, voro=voro, rect=rect) - print("mixed") - helper(intens=intens_mixed, fx=fx, fy=fy, voro=voro, rect=rect) - print("rutile") - helper(intens=intens_rutile, fx=fx, fy=fy, voro=voro, rect=rect) + index = np.abs(fx).argmin() + print(index, "/", fx.shape) + fig, axs = plt.subplots(1, 3) + axs[0].plot(intens_rutile[index, :]) + axs[0].plot(intens_rutile[index-1, :]) + axs[0].plot(intens_rutile[index+1, :]) + axs[1].plot(intens_mixed[index, :]) + axs[2].plot(intens_mono[index, :]) + for a in axs: + a.set_yscale("log") - new_intens_mono = rect.debug(Image_Wrapper(intens_mono, fx, fy)) - new_intens_mixed = rect.debug(Image_Wrapper(intens_mixed, fx, fy)) - new_intens_rutile = rect.debug(Image_Wrapper(intens_rutile, fx, fy)) - - fig, axs = plt.subplots(2, 3) - plot.plot(freqx=fx, freqy=fy, intens=new_intens_rutile, - ax_log=axs[0, 0], ax_lin=axs[1, 0], vmax=10e7) - plot.plot(freqx=fx, freqy=fy, intens=new_intens_mono, - ax_log=axs[0, 2], ax_lin=axs[1, 2], vmax=10e7) - plot.plot(freqx=fx, freqy=fy, intens=new_intens_mixed, - ax_log=axs[0, 1], ax_lin=axs[1, 1], vmax=10e7) - - for ax in axs.flatten(): - ax.set_xlim(-1, 1) - ax.set_ylim(-1, 1) - - new_intens_mono = voro.debug(Image_Wrapper(new_intens_mono, fx, fy)) - new_intens_mixed = voro.debug(Image_Wrapper(intens_mixed, fx, fy)) - new_intens_rutile = voro.debug(Image_Wrapper(intens_rutile, fx, fy)) - - fig, axs = plt.subplots(2, 3) - plot.plot(freqx=fx, freqy=fy, intens=new_intens_rutile, - ax_log=axs[0, 0], ax_lin=axs[1, 0], vmax=10e7) - plot.plot(freqx=fx, freqy=fy, intens=new_intens_mono, - ax_log=axs[0, 2], ax_lin=axs[1, 2], vmax=10e7) - plot.plot(freqx=fx, freqy=fy, intens=new_intens_mixed, - ax_log=axs[0, 1], ax_lin=axs[1, 1], vmax=10e7) - - for ax in axs.flatten(): - ax.set_xlim(-1, 1) - ax.set_ylim(-1, 1) + index = np.abs(fy).argmin() + print(index, "/", fy.shape) + fig, axs = plt.subplots(1, 3) + axs[0].plot(intens_rutile[:, index]) + axs[1].plot(intens_mixed[:, index]) + axs[2].plot(intens_mono[:, index]) + for a in axs: + a.set_yscale("log") def random(seed): @@ -198,7 +189,6 @@ def random(seed): percentage = [] counter = 0 already_inited = False - for i in tqdm.tqdm(ind): maske[np.unravel_index(i, (LEN, LEN))] = True @@ -241,6 +231,7 @@ def sample_index(p): i = np.random.choice(np.arange(p.size), p=p.ravel()) return np.unravel_index(i, p.shape) + @timeit def ising(seed): np.random.seed(seed) @@ -305,8 +296,8 @@ def ising(seed): if __name__ == "__main__": - test_square() - #test_mixed() + # test_square() + test_mixed() plt.show() # random() # np.random.seed(1234) diff --git a/2d_fourie/spin_image.py b/2d_fourie/spin_image.py index 72866c5..d81ca58 100644 --- a/2d_fourie/spin_image.py +++ b/2d_fourie/spin_image.py @@ -3,8 +3,9 @@ import scipy import scipy.fftpack as sfft import scipy.signal + class SpinImage: - resolution = 0.1 + resolution = 0.05 def __init__(self, x_pos, y_pos): x_pos = x_pos - np.min(x_pos) @@ -34,7 +35,7 @@ class SpinImage: def pad_it_square(self, additional_pad=0, size=None): h = self.img.shape[0] w = self.img.shape[1] - xx = np.maximum(h, w) + 2 * additional_pad + xx = np.maximum(h, w) + 2 * additional_pad + 1 if size is not None: xx = np.maximum(xx, size) yy = xx @@ -50,6 +51,23 @@ class SpinImage: self.img = np.pad(self.img, pad_width=( (a, aa), (b, bb)), mode="constant") + def pad_it(self, x_size, y_size): + h = self.img.shape[0] + w = self.img.shape[1] + xx = x_size + yy = y_size + self.length_x = xx * self.resolution + self.length_y = yy * self.resolution + + a = (xx - h) // 2 + aa = xx - a - h + + b = (yy - w) // 2 + bb = yy - b - w + + self.img = np.pad(self.img, pad_width=( + (a, aa), (b, bb)), mode="constant") + def percentage_gaussian(self, mask, sigma): x = np.linspace(-self.length_x / 2, self.length_x / 2, mask.shape[0])