How to resample frequency domain signal?

Question

Csanad Levente Balogh 2021-3-22

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/780242-how-to-resample-frequency-domain-signal

回答： Mathieu NOE 2021-3-29

Hi!

I have a spectrum with freqency resolution 0.01875 Hz, and I want to have 1Hz. Is there a way to convert the signal to have 1Hz frequency resolution without deforming and ruining either the phase or the amplitude of the signal? (I can not change the length or the sampling frequency of the time domain signal).

4 个评论
显示 2更早的评论隐藏 2更早的评论

Mathieu NOE 2021-3-23

ok - so you would a kind of energy (rms) averaging over the original samples within each 1 Hz band ?

are you willing to share your data ?

Csanad Levente Balogh 2021-3-29

在 MATLAB Online 中打开

data.mat

In data you can find a vector containing the frequencies and another vector containing the spectrum values. Single sided spectrum is calculated like:

x time domain signal
L length of time domain signal
fs sampling frequency
Y = fft(x);
Y = Y/L;
Y = Y(1:L/2+1);
Y(2:end-1) = 2*Y(2:end-1);

everywhere in the code, and the frequency is calculated like:

f = (fs/L)*(0:(L/2));

My problem is, I want to reduce the number of frequency components in the spectrum, without ruinning the amplitude of phase, and without lowering the energy stored in the signal. (altough I know that I necessarily loose information by reducing the number of frequeny components).

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

Mathieu NOE 2021-3-29

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/780242-how-to-resample-frequency-domain-signal#answer_661689

在 MATLAB Online 中打开

hello

first code below : this is a script showing true averaging (even on complex Y data) over N samples

your actual frequency resolution is 0.0625 Hz, so to get results with 1 Hz spacing you need to do N = 16 samples averaging

as Y (the fft data ) is complex , you can see that summing complex data can lead to "apparent" loss of magnitude, this is because the phase changes very abruptly from one frequency to the next

so summation of 2 complex with almost opposite phase can lead to much lower magnitude in the summation;

but this preserves the best the phase information

the second code below acts like a "autospectrum" averaging , we average the magnitude and the phase separately, so it better preserve the magnitude , but the resulting angle is not really meaningfull to me ; by the way "autospectrum" averaging is never displayed for phase - only magnitude

Code 1

load('data.mat')
N = 16;             % Number of elements to create the mean over
%% slow for loop method 
% zero overlap mean averaging
[m,n] = size(f);
for ci=1:fix(m/ N)
    start_index = 1+(ci-1)*N;
    stop_index = min(start_index+ N -1,m);
    avg_f(ci) =mean(f(start_index:stop_index));  % 
    avg_Y(ci) =mean(Y(start_index:stop_index));  % 
end
%% faster method 
f_new = my_mean_N_samples(f,N);
Y_new = my_mean_N_samples(Y,N);
%% plot
subplot(2,1,1),semilogy(f, abs(Y),'b',f_new, abs(Y_new),'dr',avg_f, abs(avg_Y),'g');
legend('original data','averaged with slow for loop','averaged with faster method');
xlim([3950 4050]);
xlabel('Hz');
ylabel('Magnitude');
subplot(2,1,2),plot(f, angle(Y),'b',f_new, angle(Y_new),'dr',avg_f, angle(avg_Y),'g');
legend('original data','averaged with slow for loop','averaged with faster method');
xlim([3950 4050]);
xlabel('Hz');
ylabel('Phase');
function out = my_mean_N_samples(x,n)
    x = x(:);             % make sure the data are column oriented
    s1 = size(x, 1);      % Find the next smaller multiple of n
    m  = s1 - mod(s1, n);
    y  = reshape(x(1:m), n, []);     % Reshape x to a [n, m/n] matrix
    out = transpose(sum(y, 1) / n);  % Calculate the mean over the 1st dimension
end

code 2 :

load('data.mat')
N = 16;             % Number of elements to create the mean over
Ymod = abs(Y);
Yangl = angle(Y);
%% slow for loop method 
% zero overlap mean averaging
[m,n] = size(f);
for ci=1:fix(m/ N)
    start_index = 1+(ci-1)*N;
    stop_index = min(start_index+ N -1,m);
    avg_f(ci) =mean(f(start_index:stop_index));  % 
    avg_Ymod(ci) =mean(Ymod(start_index:stop_index));  % 
    avg_Yangl(ci) =mean(Yangl(start_index:stop_index));  % 
end
%% faster method 
f_new = my_mean_N_samples(f,N);
Ymod_new = my_mean_N_samples(Ymod,N);
Yangl_new = my_mean_N_samples(Yangl,N);
%% plot
subplot(2,1,1),semilogy(f, Ymod,'b',f_new, Ymod_new,'dr',avg_f, avg_Ymod,'g');
legend('original data','averaged with slow for loop','averaged with faster method');
xlim([3950 4050]);
xlabel('Hz');
ylabel('Magnitude');
subplot(2,1,2),plot(f, Yangl,'b',f_new, Yangl_new,'dr',avg_f, avg_Yangl,'g');
legend('original data','averaged with slow for loop','averaged with faster method');
xlim([3950 4050]);
xlabel('Hz');
ylabel('Phase');
function out = my_mean_N_samples(x,n)
    x = x(:);             % make sure the data are column oriented
    s1 = size(x, 1);      % Find the next smaller multiple of n
    m  = s1 - mod(s1, n);
    y  = reshape(x(1:m), n, []);     % Reshape x to a [n, m/n] matrix
    out = transpose(sum(y, 1) / n);  % Calculate the mean over the 1st dimension
end