Init
This commit is contained in:
commit
caf9571a69
GPG Key ID:
2194FC747048A7FD
202
fft_1d.py
Normal file
202
fft_1d.py
Normal file
@ -0,0 +1,202 @@
|
||||
import numpy as np
|
||||
import matplotlib.pyplot as plt
|
||||
import tqdm
|
||||
|
||||
import multiprocessing as mp
|
||||
|
||||
RESOLUTION = 0.01
|
||||
LENGTH = 500
|
||||
|
||||
|
||||
def generate_image_from_mask(mask: np.array):
|
||||
pos_mono = np.arange(0, mask.size * 2.89 * 2, 2.89)
|
||||
pos_mono[::2][mask] -= 0.27
|
||||
pos_mono[1::2][mask] += 0.27
|
||||
return pos_mono
|
||||
|
||||
|
||||
def pad_zero(img, length):
|
||||
pad = np.zeros(length)
|
||||
img = np.append(img, pad)
|
||||
img = np.append(pad, img)
|
||||
return img
|
||||
|
||||
|
||||
def image_from_pos(pos):
|
||||
length = np.max(pos) + RESOLUTION
|
||||
x = np.arange(0, length, RESOLUTION) # angstrom
|
||||
y = np.zeros_like(x)
|
||||
ind = np.searchsorted(x, pos)
|
||||
y[ind] = 1
|
||||
return y
|
||||
|
||||
|
||||
def beugung(y, resolution):
|
||||
fft = np.fft.fft(y)
|
||||
fft_clean = np.fft.fftshift(fft)
|
||||
fft_freq = np.fft.fftfreq(y.size, resolution)
|
||||
fft_freq_clean = np.fft.fftshift(fft_freq)
|
||||
return fft_freq_clean, np.abs(fft_clean) ** 2
|
||||
|
||||
|
||||
def gaussian_convol(img):
|
||||
sigma = 100 / RESOLUTION
|
||||
mu = img.size/2
|
||||
x = np.arange(0, img.size)
|
||||
gauss = 1/(sigma * np.sqrt(2 * np.pi)) * \
|
||||
np.exp(- (x - mu)**2 / (2 * sigma**2))
|
||||
|
||||
return img*gauss
|
||||
|
||||
|
||||
def analyisis(mask):
|
||||
pos_h = generate_image_from_mask(mask)
|
||||
img = image_from_pos(pos_h)
|
||||
img = gaussian_convol(img)
|
||||
padded = pad_zero(img, int(100 / RESOLUTION))
|
||||
freq, intens = beugung(padded, RESOLUTION)
|
||||
return freq, intens
|
||||
|
||||
|
||||
def get_peaks():
|
||||
orders = np.arange(1, 2, 1)
|
||||
orders = orders / 5.78
|
||||
return np.array(orders)
|
||||
|
||||
|
||||
def eval_peaks(freq, fft):
|
||||
orders = get_peaks()
|
||||
ind = np.searchsorted(freq, orders)
|
||||
return fft[ind]
|
||||
|
||||
|
||||
def basic_test():
|
||||
mask_h = np.zeros(LENGTH).astype(bool)
|
||||
mask_l = np.ones(LENGTH).astype(bool)
|
||||
|
||||
mask_mixed = np.zeros(LENGTH).astype(bool)
|
||||
ind = (np.random.rand(30) * (LENGTH - 1)).astype(int)
|
||||
mask_mixed[ind] = True
|
||||
|
||||
mask_ner = np.zeros(LENGTH).astype(bool)
|
||||
ind = (np.random.rand(1) * (LENGTH - 31)).astype(int)
|
||||
ind = np.arange(ind, ind+30).astype(int)
|
||||
mask_ner[ind] = True
|
||||
|
||||
fig, axs = plt.subplots(4, 1)
|
||||
for mask, ax in zip([mask_h, mask_l, mask_mixed, mask_ner], axs):
|
||||
freq, ffty = analyisis(mask)
|
||||
ax.plot(freq, ffty)
|
||||
for ax in axs:
|
||||
ax.plot([1.0 / 5.78, 1.0 / 2.62, 1.0 / 3.16], [0, 0, 0], "kx")
|
||||
ax.plot([2.0 / 5.78, 2.0 / 2.62, 2.0 / 3.16], [0, 0, 0], "rx")
|
||||
ax.plot([3.0 / 5.78, 3.0 / 2.62, 3.0 / 3.16], [0, 0, 0], "bx")
|
||||
ax.set_xlim(0, 3)
|
||||
plt.show()
|
||||
|
||||
|
||||
def norm(arr):
|
||||
return arr/np.sum(arr)
|
||||
|
||||
|
||||
def norm2(arr):
|
||||
return arr
|
||||
# return arr/np.max(arr)
|
||||
|
||||
|
||||
def next_mask(mask):
|
||||
prob = np.exp((np.roll(mask, 1)*1.0 + np.roll(mask, -1)) / .1)
|
||||
prob[mask] = 0.0
|
||||
prob = norm(prob)
|
||||
|
||||
ind = np.random.choice(LENGTH, p=prob)
|
||||
mask[ind] = True
|
||||
return mask
|
||||
|
||||
|
||||
def random_loop():
|
||||
mask = np.zeros(LENGTH).astype(bool)
|
||||
ind = np.arange(0, LENGTH)
|
||||
np.random.shuffle(ind)
|
||||
percentage = []
|
||||
peaks = []
|
||||
|
||||
masks = []
|
||||
for i in ind:
|
||||
mask[i] = True
|
||||
freq, fft = analyisis(mask)
|
||||
peak = eval_peaks(freq, fft)
|
||||
percentage.append(np.mean(mask))
|
||||
peaks.append(peak)
|
||||
masks.append(mask.copy())
|
||||
masks = np.array(masks)
|
||||
plt.figure()
|
||||
plt.imshow(masks)
|
||||
plt.plot([0, 500], [406, 406])
|
||||
print()
|
||||
percentage = np.array(percentage)
|
||||
peaks = np.array(peaks)
|
||||
return percentage, peaks
|
||||
|
||||
|
||||
def nearest_loop():
|
||||
mask = np.zeros(LENGTH).astype(bool)
|
||||
percentage = []
|
||||
peaks = []
|
||||
for i in range(LENGTH):
|
||||
mask = next_mask(mask)
|
||||
freq, fft = analyisis(mask)
|
||||
peak = eval_peaks(freq, fft)
|
||||
percentage.append(np.mean(mask))
|
||||
peaks.append(peak)
|
||||
percentage = np.array(percentage)
|
||||
peaks = np.array(peaks)
|
||||
return percentage, peaks
|
||||
|
||||
|
||||
def random_helper(seed):
|
||||
np.random.seed(seed)
|
||||
#percentage_near, peaks_near = nearest_loop()
|
||||
percentage_rand, peaks_rand = random_loop()
|
||||
print("done")
|
||||
return percentage_rand, peaks_rand
|
||||
# for i in range(peaks_near.shape[1]):
|
||||
# axs[2].plot(percentage_near, norm2(
|
||||
# peaks_near[:, i]), "-", label="near")
|
||||
|
||||
# for i in range(peaks_rand.shape[1]):
|
||||
# axs[2].plot(percentage_rand, norm2(
|
||||
# peaks_rand[:, i]), ":", label="rand")
|
||||
|
||||
|
||||
def random_increase():
|
||||
fig, axs = plt.subplots(3, 1)
|
||||
|
||||
results = []
|
||||
for i in np.arange(10):
|
||||
results.append(random_helper(i))
|
||||
|
||||
for percentage_rand, peaks_rand in results:
|
||||
for i in range(peaks_rand.shape[1]):
|
||||
axs[2].plot(percentage_rand, norm2(
|
||||
peaks_rand[:, i]), ":", label="rand")
|
||||
|
||||
for ax in [axs[0], axs[1]]:
|
||||
orders = get_peaks()
|
||||
ax.plot(orders, np.zeros_like(orders), "kx")
|
||||
ax.set_xlim(0, 3)
|
||||
|
||||
mask_l = np.ones(LENGTH).astype(bool)
|
||||
mask_h = np.zeros(LENGTH).astype(bool)
|
||||
freq, ffty = analyisis(mask_l)
|
||||
axs[0].plot(freq, ffty)
|
||||
freq, ffty = analyisis(mask_h)
|
||||
axs[1].plot(freq, ffty)
|
||||
plt.xlabel("percentage")
|
||||
plt.ylabel("peak intensity")
|
||||
plt.show()
|
||||
plt.legend()
|
||||
|
||||
|
||||
if __name__ == "__main__":
|
||||
random_increase()
|
||||
137
test_fft.py
Normal file
137
test_fft.py
Normal file
@ -0,0 +1,137 @@
|
||||
import numpy as np
|
||||
import matplotlib.pyplot as plt
|
||||
|
||||
|
||||
def analysis(y, RESOLUTION):
|
||||
fft = np.fft.fft(y)
|
||||
fft_clean = np.fft.fftshift(fft)
|
||||
fft_freq = np.fft.fftfreq(y.size, RESOLUTION)
|
||||
fft_freq_clean = np.fft.fftshift(fft_freq)
|
||||
|
||||
return fft_freq_clean, np.abs(fft_clean) ** 2
|
||||
|
||||
|
||||
def play_1d():
|
||||
RESOLUTION = 0.001
|
||||
LENGTH = 10000
|
||||
|
||||
x = np.arange(0, LENGTH, RESOLUTION) # angstrom
|
||||
y = np.zeros_like(x)
|
||||
|
||||
pos_mono = np.arange(0, x.size, 2890)
|
||||
pos_mono = (
|
||||
pos_mono + np.random.normal(size=pos_mono.shape, loc=0, scale=10)
|
||||
).astype(int)
|
||||
|
||||
pos_rut = np.arange(0, x.size, 5780)
|
||||
|
||||
pos_rut = np.append(pos_rut, pos_rut - 3160)
|
||||
|
||||
# pos_rut = (pos_rut + np.random.normal(size=pos_rut.shape, loc=0, scale=10)).astype(int)
|
||||
|
||||
y[pos_rut] = 1
|
||||
y[pos_rut + 1] = 1
|
||||
y[pos_rut + 2] = 1
|
||||
|
||||
# y = np.sin(x)
|
||||
|
||||
fig, axs = plt.subplots(3, 1)
|
||||
ax = axs[0]
|
||||
ax.plot(x, y)
|
||||
|
||||
ax = axs[1]
|
||||
fft_x, fft_y = analysis(y, RESOLUTION)
|
||||
ax.plot(fft_x, fft_y)
|
||||
ax.plot([1.0 / 5.78, 1.0 / 2.62, 1.0 / 3.16], [0, 0, 0], "kx")
|
||||
ax.plot([2.0 / 5.78, 2.0 / 2.62, 2.0 / 3.16], [0, 0, 0], "rx")
|
||||
ax.plot([3.0 / 5.78, 3.0 / 2.62, 3.0 / 3.16], [0, 0, 0], "bx")
|
||||
ax.set_xlim(0, 3)
|
||||
|
||||
|
||||
def from_mask(mask):
|
||||
pos_mono = np.arange(0, mask.size * 2.89 * 2, 2.89)
|
||||
|
||||
pos_mono[::2][mask] -= 0.27
|
||||
pos_mono[1::2][mask] += 0.27
|
||||
|
||||
return pos_mono
|
||||
|
||||
|
||||
def image_from_pos(pos):
|
||||
RESOLUTION = 0.001
|
||||
LENGTH = 1000000
|
||||
|
||||
x = np.arange(0, LENGTH, RESOLUTION) # angstrom
|
||||
y = np.zeros_like(x)
|
||||
ind = np.searchsorted(x, pos)
|
||||
if np.any(ind > LENGTH):
|
||||
print("overflow")
|
||||
ind = ind[ind < LENGTH]
|
||||
y[ind] = 1
|
||||
|
||||
sigma = 500
|
||||
mu = int(LENGTH / 2)
|
||||
gaussian = (
|
||||
1
|
||||
/ (sigma * np.sqrt(2 * np.pi))
|
||||
* np.exp(-((x - mu) ** 2) / (2 * sigma * sigma))
|
||||
)
|
||||
# y = np.multiply(y, gaussian)
|
||||
|
||||
return x, y
|
||||
|
||||
|
||||
def plot_img(x, y, ax):
|
||||
ax.plot(x, y)
|
||||
|
||||
|
||||
if __name__ == "__main__":
|
||||
RESOLUTION = 0.001
|
||||
print("Done")
|
||||
|
||||
LENGTH = 1000
|
||||
mask_h = np.ones(LENGTH).astype(bool)
|
||||
pos_h = from_mask(mask_h)
|
||||
x, img_h = image_from_pos(pos_h)
|
||||
fftx, ffty_h = analysis(img_h, RESOLUTION)
|
||||
|
||||
mask_l = np.zeros(LENGTH).astype(bool)
|
||||
pos_l = from_mask(mask_l)
|
||||
x, img_l = image_from_pos(pos_l)
|
||||
fftx, ffty_l = analysis(img_l, RESOLUTION)
|
||||
print("Done")
|
||||
|
||||
mask_mixed = np.zeros(LENGTH).astype(bool)
|
||||
ind = (np.random.rand(400) * (LENGTH - 1)).astype(int)
|
||||
mask_mixed[ind] = True
|
||||
pos_mixed = from_mask(mask_mixed)
|
||||
x, img_mixed = image_from_pos(pos_mixed)
|
||||
fftx, ffty_mixed = analysis(img_mixed, RESOLUTION)
|
||||
|
||||
print("Done")
|
||||
mask_near = np.zeros(LENGTH).astype(bool)
|
||||
ind = (np.random.rand(50) * (LENGTH - 1)).astype(int)
|
||||
#for i in range(1, 8):
|
||||
# ind = np.append(ind, ind+i)
|
||||
print("Done")
|
||||
|
||||
|
||||
mask_near[ind] = True
|
||||
pos_near = from_mask(mask_near)
|
||||
x, img_near = image_from_pos(pos_near)
|
||||
fftx, ffty_near = analysis(img_near, RESOLUTION)
|
||||
|
||||
fig, axs = plt.subplots(4, 1)
|
||||
plot_img(fftx, ffty_h, axs[0])
|
||||
plot_img(fftx, ffty_l, axs[1])
|
||||
plot_img(fftx, ffty_mixed, axs[2])
|
||||
plot_img(fftx, ffty_near, axs[3])
|
||||
for ax in axs:
|
||||
ax.plot([1.0 / 5.78, 1.0 / 2.62, 1.0 / 3.16], [0, 0, 0], "kx")
|
||||
ax.plot([2.0 / 5.78, 2.0 / 2.62, 2.0 / 3.16], [0, 0, 0], "rx")
|
||||
ax.plot([3.0 / 5.78, 3.0 / 2.62, 3.0 / 3.16], [0, 0, 0], "bx")
|
||||
ax.set_xlim(0, 3)
|
||||
# play_1d()
|
||||
plt.show()
|
||||
print("Done")
|
||||
pass
|
||||
Loading…
Reference in New Issue
Block a user