goo blog サービス終了のお知らせ 

Take3's blog

日々の記録などをこつこつ

GNU Octave 3.8.2フィルタの周波数特性の確認 Vol. 2

2015-02-06 06:25:08 | アプリケーション

フィルタの周波数特性についてはまだ不安なので、
以下のようなコードでもう一度確認した.

おそらくこれで確認できていると思う。

% Sampling freq. 10 Hz

t = 0:0.1:10; % Sampling freq. 10 Hz
x = sin(t * 3.0 * 2 * pi) + sin(t * 1.0 * 2 * pi);


figure(1);
subplot(2, 1, 1);
[b, a] = butter(2, 0.4);
[h, w] = freqz(b,a, 201);
plot(w/pi * 10/2,abs(h));

filt_x = filter(b, a, x);

hold on;
[b, a] = butter(4, 0.4);
[h, w] = freqz(b,a, 201);
plot(w/pi * 10/2,abs(h), 'color', 'red', 'linewidth', 2);

filt_x2 = filter(b, a, x);

[b, a] = butter(6, 0.4);
[h, w] = freqz(b,a, 201);
plot(w/pi * 10/2,abs(h), 'color', [0.7, 0.5, 0], 'linewidth', 3);

filt_x3 = filter(b, a, x);

ylim([0 1.2]);
xlabel('Frequency Hz');

hold off;

subplot(2, 1, 2);
plot(t, filt_x, 'linewidth', 2,...
    t, filt_x2, 'linewidth', 3, 'color', 'red', ...
    t, filt_x3, 'linewidth', 4, 'color', [0.7, 0.5, 0],...
    t, x);
xlim([4 8]);
xlabel('Time');

 


GNU Octave 3.8.2フィルタの周波数特性の確認

2015-02-05 06:30:39 | アプリケーション

フィルタの周波数特性がよく分からない.

とりあえず次のようなコードを用意してみた.

% Sampling freq. 10 Hz

figure(1);
[b, a] = butter(2, 0.4);
[h, w] = freqz(b,a, 201);
plot(w/pi * 10/2,20*log10(abs(h)));

hold on;
[b, a] = butter(4, 0.4);
[h, w] = freqz(b,a, 201);
plot(w/pi * 10/2,20*log10(abs(h)), 'color', 'red', 'linewidth', 2);

[b, a] = butter(6, 0.4);
[h, w] = freqz(b,a, 201);
plot(w/pi * 10/2,20*log10(abs(h)), 'color', [0.7, 0.5, 0], 'linewidth', 3);
ylim([-60 10]);
xlabel('Frequency Hz');
ylabel('dB');

hold off;


GNU Octave filter関数の中身

2015-02-04 03:24:17 | アプリケーション

Octaveのsignal パッケージでフィルタを設計してもfilter関数の中身を知らないとOctave以外の他のところに使用できない。

中身の見当はついた。

t = 0:0.1:10; % Sampling freq. 10 Hz
x = sin(t * 2.0 * 2 * pi); %仮の信号
[b, a] = butter(2, 0.3); % 10 Hz / 2 * 0.3 (Hz)をカットオフ周波数とする2次のバターワース特性
filt_x = filter(b, a, x);

% a(1)は1と仮定しています.違う場合にはb, aをa(1)で割る必要あり.
x_buf = zeros(1, 3);
y_buf = zeros(1, 2);

y = zeros(size(x));
for k = 1:length(x);
    x_buf(2:3) = x_buf(1:2);
    x_buf(1) = x(k);
    y(k) = b(1)*x_buf(1) + b(2)*x_buf(2) + b(3)*x_buf(3) - a(2) * y_buf(1) - a(3) * y_buf(2);
    y_buf(2) = y_buf(1);
    y_buf(1) = y(k);
end
%ここまでfilter 2次の場合

figure(1);
plot(t, x, t, filt_x);

figure(2);
plot(t, x, t, filt_x, t, y);


GNU Octave 3.8.2 日本語入力は全く使えない

2015-02-01 08:27:14 | アプリケーション

GUIは非常に便利にできているが、日本語入力は使用できない。

使用できるのはコピーアンドペースト

特殊文字入りのfigureはprintで特殊文字が反映されない。
「text(1,30,'q','fontname','symbol','fontsize',24);」のようにすると入力できるが、混ぜて使用できない。

x = 0:0.1:10;
y = x.^2;

 
plot(x, y);

text(1,10,'\alpha \beta \gamma','fontname','Lucida','fontsize',20);
text(1,20,'\delta','fontname','Tahma','fontsize',24);

text(1,30,'q','fontname','symbol','fontsize',24);


set(gca, 'fontsize', 22, 'linewidth', 2);
xlabel('\mum');
ylabel('\degC');

pint('file.pdf', '-dpdf');


GNU Octave 3.8.2 プロットへの特殊文字の挿入

2015-01-28 06:27:28 | アプリケーション

例えば

x = 0:0.1:10;
y = x.^2;

plot(x, y);

text(1,10,'\alpha \beta \gamma','fontname','Lucida','fontsize',20);
text(1,20,'\delta','fontname','Tahma','fontsize',24);

set(gca, 'fontsize', 22, 'linewidth', 2);
xlabel('\mum');
ylabel('\degC');

スクリーンには特殊文字が表示されるが、「File -> Save As」でほぞんすると元に戻ってしまう。