Attenuation and Dispersion correction

9 次查看(过去 30 天)
Hello,
I am currently working on correcting the attenuation and dispersion in my measurement signal, but I am not getting the desired result.
Briefly about the measurement signal: These are pressure signals, some of which pass through a material and are measured at a piezo sensor at time t1 = 1 us and t2 = 1.45 us. The pressure signal (signal1) at time t1 can be assumed to be undamped and non-dispersive, as this is formed directly in front of the piezo sensor and is measured without passing through the material. For the pressure signal (signal2), which arrives at the sensor at time t2, an attenuation and dispersion correction must be carried out, as this signal2 has passed through the material, but originally also had the same absolute value and pulse width as signal1 (just with a different sign).
I performed the calculation using the usual formulae, but did not get the expected and desired correction.
I have the feeling that it's not such a big problem, but I can't find the solution. Perhaps one of you will spot my mistake(s) or have a new or different approach.
I am grateful for any help.
%% Correction of Dispersion and Attenuation
clear;
clc;
close all;
dt = 5e-10;
fs = 1/dt;
time = 0:dt:2.000500000000000e-06;
%Variables for gaussian-functions
GRD_Ampl = -1.02251282295196;
HS_Ampl = 0.497967501946497;
a1 = 6.28704925677956e+15;
a2 = 3.38523088091301e+14;
GRD_time = 2e-07;
HS_time = 1.453e-06;
% Variables of the specimen
d = 0.002;
v = 1600;
% Building of the gaussian-functions
Ideal_GRD = GRD_Ampl*exp(-a1 * (time - GRD_time).^2);
Ideal_HS = HS_Ampl*exp(-a2 * (time - HS_time).^2);
Signal = Ideal_HS + Ideal_GRD;
% Fourier transformation of the pressure signals
freq = fs/length(Signal)*(-length(Signal)/2+1:length(Signal)/2);
F_Ideal_GRD = fft(Ideal_GRD);
F_Ideal_HS = fft(Ideal_HS);
F_Signal = fft(Signal);
% Calculation of the attenuation and dispersion factors
alpha = -(1/d)*log(abs(F_Ideal_HS./F_Ideal_GRD)); % attenuation factor
beta = -(1/d)*(phase(F_Ideal_HS)-phase(F_Ideal_GRD)); % dispersion factor
figure;
tiledlayout(2,1);
nexttile(1)
plot(freq, fftshift(alpha), 'LineWidth', 2);
title('attenuation factor alpha');
xlabel('freqency');
ylabel('Amplitude');
grid on;
nexttile(2)
plot(freq, fftshift(beta), 'LineWidth', 2);
title('dispersion factor beta');
xlabel('freqency');
ylabel('Amplitude');
grid on;
% Calculation of the corrected signal
G = exp((-alpha*d)+(-1i*beta*d));
F_Signal_corrected = F_Signal.*G;
Signal_corrected = ifft(F_Signal_corrected);
figure;
plot(time, Signal_corrected, 'LineWidth', 2);
Warning: Imaginary parts of complex X and/or Y arguments ignored.
title('Ideal corrected pressure pulse signal');
xlabel('Time in s');
ylabel('Amplitude');
grid on

回答(1 个)

Sanju
Sanju 2024-1-4
I understand that you are trying to correct the attenuation and dispersion in your measurement signal. However, you mentioned that you are not getting the desired result. It would be helpful if you could provide more information about the expected and actual results, as well as any specific issues or error messages you encountered.
In the provided code, you are calculating the attenuation and dispersion factors using the Fourier transformation of the pressure signals. Then, you are applying these factors to correct the signal. One potential issue could be calculation of the dispersion factor, as the phase difference between the Fourier transforms might not be accurate. You can try using the “unwrap” function to improve the accuracy of the phase calculation.
Here is an example of how you can use the unwrap function to calculate the dispersion factor:
beta = -(1/d)*unwrap(phase(F_Ideal_HS) - phase(F_Ideal_GRD));
Additionally, there’s a warning about ignoring the imaginary parts of the complex arguments when plotting the corrected signal. This warning suggests that the corrected signal might have complex values. To visualize the real part of the corrected signal, you can modify the plotting code as follows:
plot(time, real(Signal_corrected), 'LineWidth', 2);
You can refer to the below documentation on “unwrap” function if required,
Hope this Helps!
Thanks.

类别

Help CenterFile Exchange 中查找有关 Filter Analysis 的更多信息

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by