শর্ট-টাইম ফুরিয়ার ট্রান্সফর্ম (STFT) ও স্পেকট্রোগ্রাম

১ · কেন সাধারণ DFT যথেষ্ট নয়

M২২-এর DFT একটি সম্পূর্ণ সিগন্যালকে একটি একক ফ্রিকোয়েন্সি-স্পেকট্রামে রূপান্তর করে — এটি বলে দেয় সিগন্যালে কোন কোন ফ্রিকোয়েন্সি আছে, কিন্তু কখন তা ঘটেছে সে তথ্য হারিয়ে যায় (গাণিতিকভাবে, DFT সূত্রের সামেশনে সময়-তথ্য সম্পূর্ণ "মিশে" যায়)। বাস্তব সিগন্যাল — কথা বলা, গান, একটি যন্ত্র চালু-বন্ধ হওয়া — প্রায় সবসময় সময়ের সাথে সাথে ফ্রিকোয়েন্সি-বিষয়বস্তু বদলায়।

STFTShort-Time Fourier Transformসিগন্যালকে ছোট ছোট সময়-ফ্রেমে ভাগ করে প্রতিটি ফ্রেমে আলাদাভাবে উইন্ডো-করা DFT প্রয়োগকারী পদ্ধতি, যার ফলাফল সময় ও ফ্রিকোয়েন্সি উভয়ের ফাংশন। এই সমস্যার সমাধান করে সিগন্যালকে ছোট ছোট ফ্রেমে ভাগ করে:

এখানে $L$ ফ্রেমের দৈর্ঘ্য, $H$ হপ-সাইজ (একটি ফ্রেম থেকে পরেরটিতে কত স্যাম্পল এগোনো হয়), $w[n]$ M২৪-এর মতো একটি উইন্ডো ফাংশন (Hamming/Hann), আর $m$ ফ্রেম-ইনডেক্স (সময়)। প্রতিটি $m$-এর জন্য একটি সম্পূর্ণ DFT বের হয় — ফলাফল একটি ২D গ্রিড: সময় ($m$) বনাম ফ্রিকোয়েন্সি ($k$)। এই গ্রিডের ম্যাগনিটিউড ছবি আকারে দেখানো হলে তাকে বলে স্পেকট্রোগ্রাম।

২ · একটি সত্যিকারের ডেমো — ফ্রিকোয়েন্সি-পরিবর্তন শনাক্তকরণ

নিচের কোডে একটি ৬৪-স্যাম্পলের সিগন্যাল তৈরি করা হয়েছে ($f_s=64$ Hz) — প্রথম ৩২ স্যাম্পল একটি $4$ Hz সাইন-তরঙ্গ, শেষ ৩২ স্যাম্পল একটি $12$ Hz সাইন-তরঙ্গ (দুটো অংশ সরাসরি জোড়া দেওয়া)। ফ্রেম-দৈর্ঘ্য $L=16$, হপ $H=16$ (নন-ওভারল্যাপিং, তাই মোট ৪টি ফ্রেম), প্রতিটি ফ্রেমে Hamming উইন্ডো প্রয়োগ করে M৫-এর হাতে-লেখা DFT দিয়ে ম্যাগনিটিউড স্পেকট্রাম বের করা হয়েছে — কোনো ইমেজ লাইব্রেরি ছাড়াই একটি টেক্সট-টেবিল হিসেবে "স্পেকট্রোগ্রাম" প্রিন্ট করা হয়েছে।

import math, cmath

def dft(x):
    N = len(x)
    X = []
    for k in range(N):
        s = 0 + 0j
        for n in range(N):
            s += x[n] * cmath.exp(-2j * math.pi * k * n / N)
        X.append(s)
    return X

def hamming(N):
    return [0.54 - 0.46 * math.cos(2 * math.pi * n / (N - 1)) for n in range(N)]

fs = 64.0
f_low = 4.0
f_high = 12.0
half = 32
low_seg = [math.sin(2 * math.pi * f_low * n / fs) for n in range(half)]
high_seg = [math.sin(2 * math.pi * f_high * n / fs) for n in range(half)]
x = low_seg + high_seg          # প্রথম অর্ধেক নিম্ন-ফ্রিকোয়েন্সি, দ্বিতীয় অর্ধেক উচ্চ-ফ্রিকোয়েন্সি
N_total = len(x)

frame_size = 16
hop = 16
win = hamming(frame_size)

frame_starts = []
start = 0
while start + frame_size <= N_total:
    frame_starts.append(start)
    start += hop

n_bins = frame_size // 2 + 1     # 0 থেকে Nyquist পর্যন্ত বিন
table = []
for fstart in frame_starts:
    seg = x[fstart:fstart + frame_size]
    wseg = [seg[i] * win[i] for i in range(frame_size)]
    X = dft(wseg)
    table.append([abs(X[k]) for k in range(n_bins)])

header = "frame | " + " | ".join(f"{k*fs/frame_size:>5.0f}Hz" for k in range(n_bins))
print(header)
for i, mags in enumerate(table):
    row = f"{i:>5} | " + " | ".join(f"{m:>6.2f}" for m in mags)
    print(row)

print()
for i, mags in enumerate(table):
    peak_bin = max(range(n_bins), key=lambda k: mags[k])
    peak_freq = peak_bin * fs / frame_size
    s0, s1 = frame_starts[i], frame_starts[i] + frame_size - 1
    print(f"frame {i} (samples {s0}-{s1}): পিক bin {peak_bin} (={peak_freq} Hz), ম্যাগনিটিউড {mags[peak_bin]:.2f}")

৩ · সময়-ফ্রিকোয়েন্সি ট্রেড-অফ

STFT-এর একটি মৌলিক সীমাবদ্ধতা আছে — ফ্রেম-দৈর্ঘ্য $L$ একসাথে সময়-রেজোলিউশন ও ফ্রিকোয়েন্সি-রেজোলিউশন উভয়ই নিয়ন্ত্রণ করে, কিন্তু বিপরীত দিকে:

উপরের ডেমোতে $L=16$ ব্যবহার করায় বিন-স্পেসিং $64/16=4$ Hz — যথেষ্ট বড় ফাঁক, কিন্তু যেহেতু ৪ Hz ও ১২ Hz উভয়ই ঠিক বিনের উপর পড়ে (bin ১ ও bin ৩), তাই এই নির্দিষ্ট উদাহরণে সমস্যা হয়নি। বাস্তব অ্যাপ্লিকেশনে (যেমন অডিও স্পেকট্রোগ্রাম, M৫২-এ বিস্তারিত) এই ট্রেড-অফ বিবেচনা করে $L$ বেছে নেওয়া হয়।