How to resample frequency domain signal?
16 visualizaciones (últimos 30 días)
Mostrar comentarios más antiguos
Csanad Levente Balogh
el 22 de Mzo. de 2021
Respondida: Mathieu NOE
el 29 de Mzo. de 2021
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 comentarios
Mathieu NOE
el 23 de Mzo. de 2021
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 ?
Respuesta aceptada
Mathieu NOE
el 29 de Mzo. de 2021
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
0 comentarios
Más respuestas (0)
Ver también
Categorías
Más información sobre Spectral Measurements en Help Center y File Exchange.
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!