矩形波をDFTまたは​FFTを利用して周波​数ドメインのスペクト​ラムに変換する方法

矩形波をDFTまたは、FFTを利用して周波数ドメインのスペクトラムに変換してプロットしたいです。 ただし、矩形波の立ち上がり時間と立下り時間を無限小ではなく、有限時間に設定してその波形のスペクトラムの検証を行いたいです。 どのようにすればよろしいでしょうか? 具体的なコード等をご提示いただきたく思います。 宜しくお願い致します。

Answer 1

[Tohru Kikawada]

ご自身でどのようなことを試されたのか、コード例なども交えて記載いただけるとよろしいかと思います。

参考までに矩形波をなまらせて周波数解析した例をご紹介します。

%%周期10Hzの矩形波を生成
Fs = 100; % サンプリング周波数
f = 2;
t = 0:1/Fs:1;
x = double(sin(2*pi*f*t)>=0);
plot(t,x);
xlabel('Time (s)'); ylabel('Amplitude');

%%ローパスフィルタをかけてなまらせる
alpha = 0.4;
b = alpha;
a = [1 alpha-1];
xhat = filter(b,a,x);
hold on;
plot(t,xhat);
legend('Square wave','Filtered square wave');

%%FFTをかけて周波数解析
figure;
periodogram(x,hamming(length(x)),512,Fs);
hold on;
periodogram(xhat,hamming(length(x)),512,Fs);
hAxes = gca;
hAxes.Children(2).Color = 'red';
legend('Square wave','Filtered square wave');

Answer 2

[Atsushi Matsumoto]

立ち上がりの時間で波形生成できたら良いのかもしれませんが、まずは簡単に三角波を飽和させて矩形波を生成してみました。

Fs = 10000;  % サンプリング周波数
dt = 1/Fs;  % サンプル時間
t = 0:dt:1-dt;  % 1sの時間軸
Fsig = 50;  % 矩形波の周波数

for A = 1:5:21
    x = A*sawtooth(2*pi*Fsig*t, 0.5); % 三角波生成してA倍
    x(x>=1) = 1;x(x<=-1)=-1;    % +-1で飽和処理=>矩形波
    figure(1)
    plot(t,x), hold on          % 時間軸応答を追加プロット
    figure(2)
    [p, f] = pspectrum(x);      % スペクトル計算
    plot(f, 20*log10(p)), hold on   % 周波数応答プロット

end

figure(1)
xlim([0 0.2])   % X軸範囲指定
legend('1', '6', '11', '16', '21')
figure(2)
ylim([-100 0])  % Y軸範囲指定