Files
ALF/camerata/main_simple.py
2023-10-23 19:03:37 +02:00

117 lines
2.3 KiB
Python
Executable File

#!/usr/bin/env python
from time import sleep
from collections import deque
import numpy as np
import scipy.signal
from zoetrope import aniplot as plt
from bstrd import BS, bsstream
plt.blit = False
# config
cam_name = "SARES11-XMI125-C4P1.jet_projection"
#xmin = 0
#xmax = 400
#ymin = 400
#ymax = 700
#proj_axis = 0
um_per_px = 2.28
buffer_length = 1000
def get_data(cam):
jet = 0
for i in range(1):
jet += cam.get()
next(bsstream)
proj = jet
return proj
def find_spacing(py, extend=5):
py -= py.mean()
if extend > 1:
py = np.concatenate((py, np.zeros(len(py) * extend)))
spectrum = np.abs(np.fft.fft(py))
freq = np.fft.fftfreq(len(spectrum))
which = (freq > 0.005)
spectrum = spectrum[which]
freq = freq[which]
index_peaks = scipy.signal.find_peaks(spectrum)[0]
peaks_pos = freq[index_peaks]
peak_height = spectrum[index_peaks]
cut_freq = peaks_pos[np.argmax(peak_height)]
return freq, spectrum, cut_freq
def make_xs(arr):
return np.arange(len(arr))
# create channel
cam = BS(cam_name)
proj = get_data(cam)
pixel_proj = make_xs(proj)
plt.fig.set_figheight(12)
plt.fig.set_figwidth(9)
plt.suptitle(cam_name)
plt.subplot(311)
plt.title("projection")
pln_proj = plt.plot(proj, pixel_proj)
plt.grid()
plt.subplot(312)
plt.title("FFT")
plt.xlabel("frequency in 1/pixel")
pln_fft_spec = plt.plot([0])
plt.xlim(0,0.2)
plt.grid()
pln_fft_cut = plt.plot([0])
plt.subplot(313)
plt.title(f"spacing assuming {um_per_px} μm per pixel")
pln_spac = plt.plot([0], 'o')
pln_spac.ax.axhline(y=80, color='r')
plt.ylim(0,250)
plt.grid()
plt.tight_layout()
spacings = deque(maxlen=buffer_length)
for counter in plt.show():
print(counter)
next(bsstream)
proj = get_data(cam)
pln_proj.set(proj, pixel_proj)
freq, spectrum, cut_freq = find_spacing(proj, extend=10)
spac = um_per_px / cut_freq
spacings.append(spac)
pln_spac.set(make_xs(spacings), spacings)
pln_fft_spec.set(freq, spectrum)
pln_fft_cut.set([cut_freq]*2, [spectrum.min(), spectrum.max()])
pln_fft_spec.ax.set_title(f"FFT\ncurrent spacing: {spac:.1f} μm")
# this, I need to move into the library
ps = [pln_proj, pln_spac, pln_fft_spec]
for p in ps:
p.ax.relim()
p.ax.autoscale_view()
bsstream.close()