0
私は受信したFT信号を評価し、信号のためのモートルウェーブレットを計算するためにmatlab関数（バージョン7.10.0.499（R2010a））を書いています。私は似たようなプログラムを持っていますが、私はそれをより読みやすくし、数学的な用語に近づける必要がありました。出力プロットは、周波数の強度を示す色の2Dプロットであると想定されます。私のプロットは時間ごとにすべての周波数を同じように持っているようです。プログラムは各周波数の時間ごとにfftを作成するので、それを見るもう1つの方法は、同じ行がforループのステップごとに繰り返されるということです。問題は正しいプロットを返す元のプログラムでチェックしたことです。値の名前とコードの構成方法を超えて差異を見つけることはできません。モーレットウェーブレット変換関数は無意味なプロットを返します
function[msg] = mile01_wlt(FT_y, f_mn, f_mx, K, N, F_s)
%{
Fucntion to perform a full wlt of a morlet wavelett.
optimization of the number of frequencies to be included.
FT_y satisfies the FT(x) of 1 envelope and is our ft signal.
f min and max enter into the analysis and are decided from
the f-image for optimal values.
While performing the transformation there are different scalings
on the resulting "intensity".
Plot is made with a 2D array and a colour code for intensity.
version 05.05.2016
%}
%--------------------------------------------------------------%
%{
tableofcontents:
1: determining nr. of analysis f, prints and readies f's to be used.
2: ensuring correct orientation of FT_y
3:defining arrays
4: declaring waveletdiagram and storage of frequencies
5: for-loop over all frequencies:
6: reducing file to manageable size by truncating time.
7: marking plot to highlight ("randproblemer")
8: plotting waveletdiagram
%}
%--------------------------------------------------------------%
%1: determining nr. of analysis f, prints and readies f's to be used.
DF = floor(log(f_mx/f_mn)/log(1+(1/(8*K)))) + 1;% f-spectre analysed
nr_f_analysed = DF %output to commandline
f_step = (f_mx/f_mn)^(1/(DF-1)); % multiplicative step for new f_a
f_a = f_mn; %[Hz] frequency of analysis
T = N/F_s; %[s] total time sampled
C = 2.0; % factor to scale Psi
%--------------------------------------------------------------%
%2: ensuring correct orientation of FT_y
siz = size(FT_y);
if (siz(2)>siz(1))
FT_y = transpose(FT_y);
end;
%--------------------------------------------------------------%
%3:defining arrays
t = linspace(0, T*(N-1)/N, N); %[s] timespan
f = linspace(0, F_s*(N-1)/N, N); %[Hz] f-specter
%--------------------------------------------------------------%
%4: declaring waveletdiagram and storage of frequencies
WLd = zeros(DF,N); % matrix of DF rows and N collumns for storing our wlt
f_store = zeros(1,DF); % horizontal array for storing DF frequencies
%--------------------------------------------------------------%
%5: for-loop over all frequencies:
for jj = 1:DF
o = (K/f_a)*(K/f_a); %factor sigma
Psi = exp(- 0*(f-f_a).*(f-f_a)); % FT(\psi) for 1 envelope
Psi = Psi - exp(-K*K)*exp(- o*(f.*f)); % correctional element
Psi = C*Psi; %factor. not set in stone
%next step fits 1 row in the WLd (3 alternatives)
%WLd(jj,:) = abs(ifft(Psi.*transpose(FT_y)));
WLd(jj,:) = sqrt(abs(ifft(Psi.*transpose(FT_y))));
%WLd(jj,:) = sqrt(abs(ifft(Psi.*FT_y))); %for different array sizes
%and emphasizes weaker parts.
%prep for next round
f_store (jj) = f_a; % storing used frequencies
f_a = f_a*f_step; % determines the next step
end;
%--------------------------------------------------------------%
%6: reducing file to manageable size by truncating time.
P = floor((K*F_s)/(24*f_mx));%24 not set in stone
using_every_P_point = P %printout to cmdline for monitoring
N_P = floor(N/P);
points_in_time = N_P %printout to cmdline for monitoring
% truncating WLd and time
WLd2 = zeros(DF,N_P);
for jj = 1:DF
for ii = 1:N_P
WLd2(jj,ii) = WLd(jj,ii*P);
end
end
t_P = zeros(1,N_P);
for ii = 1:N_P % set outside the initial loop to reduce redundancy
t_P(ii) = t(ii*P);
end
%--------------------------------------------------------------%
%7: marking plot to highlight boundary value problems
maxval = max(WLd2);%setting an intensity
mxv = max(maxval);
% marks in wl matrix
for jj= 1:DF
m = floor(K*F_s/(P*pi*f_store(jj))); %finding edges of envelope
WLd2(jj,m) = mxv/2; % lower limit
WLd2(jj,N_P-m) = mxv/2;% upper limit
end
%--------------------------------------------------------------%
%8: plotting waveletdiagram
figure;
imagesc(t_P, log10(f_store), WLd2, 'Ydata', [1 size(WLd2,1)]);
set(gca, 'Ydir', 'normal');
xlabel('Time [s]');
ylabel('log10(frequency [Hz])');
%title('wavelet power spectrum'); % for non-sqrt inensities
title('sqrt(wavelet power spectrum)'); %when calculating using sqrt
colorbar('location', 'southoutside');
msg = 'done.';
エラーメッセージはありません。正確に何が間違っているのかは不明です。 私はすべてのガイドラインに従っていきたいと考えています。そうでなければ、私は謝罪します。
編集： 私の呼び出しプログラム： ％設定パラメータ N = 2 ^（16）; ％|サンプリングする点の数 F_s = 3.2e6; ％Hz |サンプリング周波数： T_t = N/F_s; ％s |サンプル時間の秒単位の長さ f_c = 2.0e5; ％Hz |搬送波周波数 f_m = 8./T_t; ％Hz |変調波周波数 w_c = 2 * pi * f_c; ％Hz |搬送波の角周波数（「ω」） w_m = 2 * pi * f_m; ％Hz |波
% establishing parameter arrays
t = linspace(0, T_t, N);
% function variables
T_h = 2*f_m.*t; % dimless | 1/2 of the period for square signal
% combined carry and modulated wave
% y(t) eq. 1):
y_t = 0.5.*cos(w_c.*t).*(1+cos(w_m.*t));
% y(t) eq. 2):
% y_t = 0.5.*cos(w_c.*t)+0.25*cos((w_c+w_m).*t)+0.25*cos((w_c-w_m).*t);
%square wave
sq_t = cos(w_c.*t).*(1 - mod(floor(t./T_h), 2)); % sq(t)
% the following can be exchanged between sq(t) and y(t)
plot(t, y_t)
% plot(t, sq_t)
xlabel('time [s]');
ylabel('signal amplitude');
title('plot of harmonically modulated signal with carrying wave');
% title('plot of square modulated signal with carrying wave');
figure()
hold on
% Fourier transform and plot of freq-image
FT_y = mile01_fftplot(y_t, N, F_s);
% FT_sq = mile01_fftplot(sq_t, N, F_s);
% Morlet wavelet transform and plot of WLdiagram
%determining K, check t-image
K_h = 57*4; % approximation based on 1/4 of an envelope, harmonious
%determining f min and max, from f-image
f_m = 1.995e5; % minimum frequency. chosen to showcase all relevant f
f_M = 2.005e5; % maximum frequency. chosen to showcase all relevant f
%calling wlt function.
name = 'mile'
msg = mile01_wlt(FT_y, f_m, f_M, K_h, N, F_s)
siz = size(FT_y);
if (siz(2)>siz(1))
FT_y = transpose(FT_y);
end;
name = 'arnt'
msg = arnt_wltransf(FT_y, f_m, f_M, K_h, N, F_s)
を変調角周波数（「オメガ」）は、時間画像は、一定の周波数を有しているが、振幅はガウス曲線をresempling発振します。私のコードは時間の経過と共に鮮明に分割された画像を返します。各画像には1つの周波数しか保持されません。これは、時間の経過とともにスペクトル全体にわたる強度の変化を反映するはずである。 ご協力いただきありがとうございます！
プロットに入力したいデータを見ましたか？もしそうなら、あなたの理解に正しいデータがありますか？もしそうでなければ、プロットは問題ではありませんが、もしそうなら、私たちにデータを与えて、どのようにプロットを以下のようにしたいかを教えてください。 – tim
なぜFFT出力でウェーブレット変換をしていますか？ –
morletウェーブレットは、単一の "エンベロープ"に分解され、これらは本質的に何かに分割されます。 ifft（fft（x）* fft（psi）） psiは私が実行する分析ソリューションです各ステップとfft（x） は私が送る信号です – Anders