22

それぞれに位相シフトがある 3 つの同一の波を生成しました。例えば：

t = 1:10800; % generate time vector
fs = 1; % sampling frequency (seconds)
A = 2; % amplitude
P = 1000; % period (seconds), the time it takes for the signal to repeat itself
f1 = 1/P; % number of cycles per second (i.e. how often the signal repeats itself every second).
y1 = A*sin(2*pi*f1*t); % signal 1
phi = 10; % phase shift
y2 = A*sin(2*pi*f1*t + phi); % signal 2
phi = 15; % phase shift
y3 = A*sin(2*pi*f1*t + phi); % signal 3

YY = [y1',y2',y3'];

plot(t,YY)

ここに画像の説明を入力

ここで、波間のこの位相シフトを検出する方法を使用したいと思います。これを行うポイントは、最終的にこの方法を実際のデータに適用し、信号間の位相シフトを特定できるようにすることです。

これまでのところ、各波と最初の波の間のクロススペクトルを計算することを考えていました (つまり、位相シフトなし):

for i = 1:3;
    [Pxy,Freq] = cpsd(YY(:,1),YY(:,i));
    coP = real(Pxy);
    quadP = imag(Pxy);
    phase(:,i) = atan2(coP,quadP);
end

しかし、これが意味をなすかどうかはわかりません。

他の誰かがこれに似たようなことをしましたか? 望ましい結果は、波 2 と 3 に対してそれぞれ 10 と 15 の位相シフトを示すはずです。

アドバイスをいただければ幸いです。

4

7 に答える 7

5

2つの正弦波信号の場合、複素相関係数の位相が必要なものを提供します。テストするためのmatlabがないため、pythonの例（scipyを使用）のみを提供できます。

x1 = sin( 0.1*arange(1024) )
x2 = sin( 0.1*arange(1024) + 0.456)
x1h = hilbert(x1)
x2h = hilbert(x2)
c = inner( x1h, conj(x2h) ) / sqrt( inner(x1h,conj(x1h)) * inner(x2h,conj(x2h)) )
phase_diff = angle(c)

matlab にも関数 corrcoeff があり、これも機能するはずです (python は虚数部を破棄します)。つまり、 c = corrcoeff(x1h,x2h) は matlab で動作するはずです。

于 2014-12-19T09:22:57.840 に答える
1

これはあなたのコードを少し変更したものです:phi = 10実際には度数であり、正弦関数では、位相情報はほとんどラジアンで表されるためdeg2rad(phi)、次のように変更する必要があります:

t = 1:10800; % generate time vector
fs = 1; % sampling frequency (seconds)
A = 2; % amplitude
P = 1000; % period (seconds), the time it takes for the signal to repeat itself
f1 = 1/P; % number of cycles per second (i.e. how often the signal repeats itself every second).
y1 = A*sin(2*pi*f1*t); % signal 1
phi = deg2rad(10); % phase shift 
y2 = A*sin(2*pi*f1*t + phi); % signal 2
phi = deg2rad(15); % phase shift
y3 = A*sin(2*pi*f1*t + phi); % signal 3

YY = [y1',y2',y3'];

plot(t,YY)

次に、前述のチポデットのように周波数領域法を使用します

fft_y1 = fft(y1);
fft_y2 = fft(y2);
phase_rad = angle(fft_y1(1:end/2)/fft_y2(1:end/2));
phase_deg = rad2deg(angle(fft_y1(1:end/2)/fft_y2(1:end/2)));

これで位相シフトの推定値が得られますerror = +-0.2145

于 2016-01-14T02:56:30.077 に答える
1

正しい信号で:

t = 1:10800; % generate time vector
fs = 1; % sampling frequency (seconds)
A = 2; % amplitude
P = 1000; % period (seconds), the time it takes for the signal to repeat itself
f1 = 1/P; % number of cycles per second (i.e. how often the signal repeats itself every second).
y1 = A*sin(2*pi*f1*t); % signal 1
phi = 10*pi/180; % phase shift in radians
y2 = A*sin(2*pi*f1*t + phi); % signal 2
phi = 15*pi/180; % phase shift in radians
y3 = A*sin(2*pi*f1*t + phi); % signal 3

以下が機能するはずです。

>> acos(dot(y1,y2)/(norm(y1)*norm(y2)))
>> ans*180/pi
ans =  9.9332
>> acos(dot(y1,y3)/(norm(y1)*norm(y3)))
ans =  0.25980
>> ans*180/pi
ans =  14.885

それがあなたの「本当の」信号にとって十分かどうかは、あなただけが知ることができます.

于 2014-12-18T11:42:06.730 に答える
0

周波数が分かっていて、完全な FFT を使用するのではなく、位相だけを見つけたい場合は、Goertzel アルゴリズムを検討することをお勧めします。これは、単一周波数の DFT を計算するより効率的な方法です (FFT は、すべての周波数）。

適切な実装については、https: //www.mathworks.com/matlabcentral/fileexchange/35103-generalized-goertzel-algorithmおよびhttps://asp-eurasipjournals.springeropen.com/track/pdf/10.1186/1687-6180を参照してください。 -2012-56.pdf

于 2020-11-20T09:15:04.980 に答える