6.15. סינון וספקטרוגרמות

סינון, החלקה וספקטרום גודל הם המשימות הסמוכות ל-FFT גולמי: החלקה או סינון מעביר-פס של זרם דגימות, חישוב ספקטרום גודל בלולאת זרימה ללא הקצאה, ופירוש מחדש של חוצצי התקנים היקפיים גולמיים כמערכי מספרים בנקודה צפה. הכלים הזמינים:

  • sosfilt() – מסנן דיגיטלי המיושם באמצעות מקטעים מדורגים מסדר שני.

  • spectrogram() – גודל abs(fft(...)) ללא הקצאות ביניים.

  • from_int16_buffer() ושאר פונקציות העזר from_*_buffer של ulab.utils – מושכות מערך מספרים בנקודה צפה מתוך חוצץ שה-dtype שלו אינו מכוסה על ידי frombuffer() המובנית.

6.15.1. סינון באמצעות sosfilt

sosfilt() מיישמת מסנן דיגיטלי בעל תגובה אימפולסיבית אינסופית (IIR) כשרשרת של מקטעים מסדר שני (SOS) – צורה יציבה מבחינה נומרית. sos הוא רצף של מקטעים באורך 6; x הוא הקלט החד-ממדי:

from ulab import numpy as np
from ulab import scipy as sp

x   = np.array([0, 1, 2, 3, 4, 5, 6, 7, 8, 9])
sos = [[1, 2, 3, 1, 5, 6],
       [1, 2, 3, 1, 5, 6]]
y   = sp.signal.sosfilt(sos, x)

כל שורה ב-sos מחזיקה שישה מקדמים [b0, b1, b2, a0, a1, a2] עבור מקטע ביקווד אחד. מערך ה-sos מחושב בדרך כלל מראש על גבי מחשב באמצעות scipy.signal.iirfilter(..., output='sos') ומועתק לתוך סקריפט המצלמה כליטרל של Python.

מילת המפתח האופציונלית zi= נושאת את מצב המסנן בין חוצצים. העבירו מצב התחלתי בצורה (n_sections, 2) והפונקציה מחזירה (y, zf) – הפלט המסונן והמצב הסופי – כך שהמצב הסופי של חוצץ אחד מזין את המצב ההתחלתי של הבא:

y0, zf0 = sp.signal.sosfilt(sos, buffer0, zi=zi)
y1, zf1 = sp.signal.sosfilt(sos, buffer1, zi=zf0)
# ...

זוהי התבנית הסטנדרטית למסנן זרימה על נתונים בחוצץ – קלט מיקרופון הנקרא 1024 דגימות בכל פעם, דגימות ADC הנצברות במנות מונחות-DMA, קריאות IMU הנאספות לאורך חלון.

6.15.2. ספקטרוגרמות

spectrogram() מחשבת את גודל התמרת פורייה. היא שקולה מושגית ל-np.sqrt(real * real + imag * imag) לאחר קריאה ל-fft(), אך מקפלת את העבודה לקריאה אחת – מבלי להחזיק את הביניים real * real, imag * imag, את הסכום או את מערך הגודל המפורש ב-RAM בשום שלב. הדבר הופך אותה לכלי הנכון בכל לולאה שבה מחושב ספקטרום שוב ושוב:

from ulab import numpy as np
from ulab import utils

x        = np.linspace(0, 10, num=1024)
spectrum = utils.spectrogram(x)

צורת הארגומנט משקפת את fft(): מערך ממשי אחד, או זוג (real, imag) כאשר לקלט יש חלק מדומה.

שלוש מילות מפתח מסייעות בהקצאה:

  • scratchpad=None – מערך מספרים בנקודה צפה צפוף חד-ממדי באורך 2 * len(signal) ש-spectrogram() משתמשת בו כמרחב עבודה.

  • out=None – מערך מספרים בנקודה צפה חד-ממדי לכתיבת התוצאה לתוכו.

  • log=False – כאשר הערך True, ביצוע log() של הגודל לפני ההחזרה, מקופל לאותה קריאה.

תבנית הזרימה היא להקצות הכול פעם אחת ולעולם לא להקצות שוב:

from ulab import numpy as np
from ulab import utils

N = 1024
scratch = np.zeros(2 * N)
out     = np.zeros(N)

while True:
    signal = read_samples(N)
    utils.spectrogram(signal, out=out, scratchpad=scratch,
                      log=True)
    # out now holds the log-magnitude spectrum for this window ...

השוו לגרסה המובנת-מאליה-אך-בזבזנית:

while True:
    signal   = read_samples(N)
    spectrum = np.log(utils.spectrogram(signal))   # two allocations

שתיהן מפיקות את אותם מספרים, אך הגרסה הראשונה אינה מקצה דבר בתוך הלולאה – המצלמה משמרת את אותו זיכרון בשימוש בכל איטרציה, והלולאה רצה מהר יותר.

6.15.3. חוצצי התקנים היקפיים רחבים מ-16 סיביות

frombuffer() מטפלת רק ב-dtypes ש-numpy עצמה מגדירה (uint8 / int8, uint16 / int16, float). כאשר התקן היקפי מפיק דגימות שלמות של 32 סיביות – ADC של 24 או 32 סיביות, מיקרופון ברזולוציה גבוהה – ulab.utils חושף פונקציות המרה מפורשות:

כל אחת מקבלת חוצץ דמוי-bytes ומחזירה ndarray של מספרים בנקודה צפה:

from ulab import utils

buf = bytearray([1, 1, 0, 0, 0, 0, 0, 255])
utils.from_uint32_buffer(buf)
# array([257.0, 4278190080.0])

הפונקציות מקבלות את אותם בקרי חיסכון בהקצאה כמו spectrogram():

  • count= ו-offset= לדילוג על כותרת או להגבלת הקריאה.

  • out= לכתיבה לתוך מערך מספרים בנקודה צפה שהוקצה מראש.

  • byteswap=True כאשר ההתקן ההיקפי חולק על ה-MCU בעניין סדר הבתים.

התבנית המשולבת – קריאת from_int32_buffer() אחת ישירות לתוך קריאת spectrogram() אחת, שתיהן עם חוצצי out= מחוץ ללולאה – היא התבנית הנכונה למנתח ספקטרום זרימה הרץ על מיקרופון ברזולוציה גבוהה.