I am getting wrong results for the Continuous wavelet transform (CWT) on converting the frequency axis (Y Axis) from log scale to linear scale for the same frequency limits.
11 次查看(过去 30 天)
显示 更早的评论
Hello All,
I have been trying to apply CWT on a signal. The function cwt(signal, Fs) gives the correct frequency when the Frequency scale is in the log scale but the incorrect frequency on converting the axis to linear scale.
Please suggest how to correct this problem. Will upgrading to newer version of MATLAB help?
Code -
fs = 20000; % Sampling frequency in Hz
t1 = 0:1/fs:10-1/fs; % Time vector for the first 10 seconds
t2 = 10:1/fs:20-1/fs; % Time vector for 10 to 20 seconds
t3 = 20:1/fs:30-1/fs; % Time vector for 20 to 30 seconds
% Create Signal
% 1. Sine wave of amplitude 10 and frequency 1 kHz for the first 10 seconds
signal1 = 10 * sin(2 * pi * 1000 * t1) + 3 * sin(2 * pi * 3000 * t1);
% 2. Combination of 2 sine waves:
signal2 = 4 * sin(2 * pi * 2000 * t2) + 6 * sin(2 * pi * 500 * t2);
% 3. Sine wave of amplitude 20 and frequency 300 Hz for time 20-30 seconds
signal3 = 8 * sin(2 * pi * 300 * t3) + 5 * sin(2 * pi * 1500 * t3);
% Concatenate all parts of the signal to form the complete signal
completeSignal = [signal1, signal2, signal3];
% Time vector for the complete signal
t = [t1, t2, t3];
% Wavelet Analysis
frequencyLimits = [0 10000]; % Frequency range from 0 Hz to 10,000 Hz
voicesPerOctave = 48; % Choosing 48 voices per octave for finer frequency resolution
% Perform Continuous Wavelet Transform
% [wt, f] = cwt(completeSignal, 'amor', fs, 'FrequencyLimits', frequencyLimits, 'VoicesPerOctave', voicesPerOctave);
cwt(completeSignal, 'amor', fs);
% Plot the results
% figure (2);
% imagesc(t, f, abs(wt));
% axis xy;
% xlabel('Time (s)');
% ylabel('Frequency (Hz)');
% title(['Continuous Wavelet Transform (CWT) of the Signal (0 to 10 kHz) with ', num2str(voicesPerOctave), ' Voices/Octave']);
% colorbar;
Plot on using log scale for the Frequency Axis (Y Axis) -
Plot on using linear scale for the Frequency Axis (Y Axis) -
2 个评论
Voss
2024-9-18
@Sumanta Kumar Dutta: Looks like your images accidentally got deleted when I edited the question to apply code formatting. Sorry about that; add them again if necessary.
采纳的回答
Voss
2024-9-18
Use a surface rather than an image. A surface allows you to use arbitrary x and y values, whereas an image requires linearly spaced x and y, which is not the case with your f vector since it is logarithmically spaced.
See below for how to (more-or-less) reproduce the plot created by cwt() on both a linear and log y-scale. (The code wouldn't run here in the forum environment when specifying 48 voicesPerOctave or when plotting all columns of t and wt, so I made modifications to get it to run.)
fs = 20000; % Sampling frequency in Hz
t1 = 0:1/fs:10-1/fs; % Time vector for the first 10 seconds
t2 = 10:1/fs:20-1/fs; % Time vector for 10 to 20 seconds
t3 = 20:1/fs:30-1/fs; % Time vector for 20 to 30 seconds
% Create Signal
% 1. Sine wave of amplitude 10 and frequency 1 kHz for the first 10 seconds
signal1 = 10 * sin(2 * pi * 1000 * t1) + 3 * sin(2 * pi * 3000 * t1);
% 2. Combination of 2 sine waves:
signal2 = 4 * sin(2 * pi * 2000 * t2) + 6 * sin(2 * pi * 500 * t2);
% 3. Sine wave of amplitude 20 and frequency 300 Hz for time 20-30 seconds
signal3 = 8 * sin(2 * pi * 300 * t3) + 5 * sin(2 * pi * 1500 * t3);
% Concatenate all parts of the signal to form the complete signal
completeSignal = [signal1, signal2, signal3];
% Time vector for the complete signal
t = [t1, t2, t3];
% Wavelet Analysis
frequencyLimits = [0 10000]; % Frequency range from 0 Hz to 10,000 Hz
voicesPerOctave = 48; % Choosing 48 voices per octave for finer frequency resolution
% Perform Continuous Wavelet Transform
[wt, f] = cwt(completeSignal, 'amor', fs, 'FrequencyLimits', frequencyLimits);%, 'VoicesPerOctave', voicesPerOctave);
cwt(completeSignal, 'amor', fs);
yl = ylim();
% Plot the results
figure();
ax = gca();
% surface(ax,t, f/1000, abs(wt), 'EdgeColor', 'none');
surface(ax,t(1:10:end), f/1000, abs(wt(:,1:10:end)), 'EdgeColor', 'none');
ax.YLim = yl;
ax.YScale = 'log';
ax.YTick = 10.^(-4:0);
xlabel(ax,'Time (s)');
ylabel(ax,'Frequency (kHz)');
title(ax,'reproduction of cwt() plot with log y-scale');
colorbar(ax)
% Plot the results
figure()
ax = gca();
% surface(ax,t, f/1000, abs(wt), 'EdgeColor', 'none');
surface(ax,t(1:10:end), f/1000, abs(wt(:,1:10:end)), 'EdgeColor', 'none');
ax.YLim = yl;
xlabel(ax,'Time (s)');
ylabel(ax,'Frequency (kHz)');
title(ax,'reproduction of cwt() plot with linear y-scale');
colorbar(ax)
2 个评论
Voss
2024-9-19
You're welcome!
The cone of influence is the third output from cwt, as described here:
See also:
更多回答(1 个)
Jonas
2024-9-18
编辑:Jonas
2024-9-19
you could interpolate the image pixels, but this makes the upper frequency values barely visible
fs = 10000; % Sampling frequency in Hz
t1 = 0:1/fs:3-1/fs; % Time vector for the first 10 seconds
t2 = 3:1/fs:6-1/fs; % Time vector for 10 to 20 seconds
t3 = 6:1/fs:9-1/fs; % Time vector for 20 to 30 seconds
% Create Signal
% 1. Sine wave of amplitude 10 and frequency 1 kHz for the first 10 seconds
signal1 = 10 * sin(2 * pi * 1000 * t1) + 3 * sin(2 * pi * 3000 * t1);
% 2. Combination of 2 sine waves:
signal2 = 4 * sin(2 * pi * 2000 * t2) + 6 * sin(2 * pi * 500 * t2);
% 3. Sine wave of amplitude 20 and frequency 300 Hz for time 20-30 seconds
signal3 = 8 * sin(2 * pi * 300 * t3) + 5 * sin(2 * pi * 1500 * t3);
% Concatenate all parts of the signal to form the complete signal
completeSignal = [signal1, signal2, signal3];
% Time vector for the complete signal
t = [t1, t2, t3];
% Wavelet Analysis
frequencyLimits = [0 fs/2];
voicesPerOctave = 48; % Choosing 48 voices per octave for finer frequency resolution
% Perform Continuous Wavelet Transform
[wt, f] = cwt(completeSignal, 'amor', fs, 'FrequencyLimits', frequencyLimits, 'VoicesPerOctave', voicesPerOctave);
wt=abs(wt);
% Plot the results
figure;
cwt(completeSignal, 'amor', fs, 'FrequencyLimits', frequencyLimits, 'VoicesPerOctave', voicesPerOctave);
f_linear = min(f):10:max(f); % to reduce memory a bit we use step of 10 here
% Interpolate the wavelet transform to the linear frequency scale
wt_linear = interp1(f, abs(wt), f_linear, 'makima');
t = (0:length(completeSignal)-1)/fs;
% Plot the scalogram with linear frequency scale
imagesc(t, f_linear, wt_linear);
set(gca,'YDir','normal')
2 个评论
Jonas
2024-9-19
you can get the cone of influence y coordinates from the third output of the cwt function
另请参阅
类别
在 Help Center 和 File Exchange 中查找有关 Continuous Wavelet Transforms 的更多信息
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!