How to find out the amplitude and duration of a sine wave of different amplitudes?
Show older comments
Hi How to find out the amplitude and duration of a sine wave of different amplitudes?
Consider a voltage sag ie, a 50hz signal of amplitude less than the normal one.I have to find out the amplitude and duration at wich this disturbance occurs.
I have read some articles which gives some information that we can extract time and frequency information using wavelets db4 using multi resolution analysis .As i am a beginner i have just a rough idea about the work.Can anyone please direct me to accomplish my works?
Thanks
1 Comment
A Lumbini
on 20 Feb 2024
Hi have you got the answer
Answers (4)
I would not use wavelets. , I would use stft().
Make a signal which is 3 seconds long. The frequency is 50 Hz. The amplitude is 230 V rms (i.e. peak amplitude=325 V). The peak amplitude dips to 275 V for 1 second in the middle of the 3 second long signal.
fs=1000; % sampling rate (Hz)
t=0:1/fs:3-1/fs; % time (s)
f=50; % frequency (Hz)
a1=325; a2=275; % amplitude (V)
x=a1*sin(2*pi*f*t); % signal
x(t>=1 & t<2)=a2*sin(2*pi*f*t(t>=1 & t<2));
Plot the signal versus time.
plot(t,x,'-r');
grid on; xlabel('Time (s)'); ylabel('Voltage (V)');
The plot above shows that the amplitude is lower from t=1 to t=2.
Compute the time-varying Fourier transform, and plot it:
nfft=100; % number of points to use for FFT
[s,f1,t1]=stft(x,fs,Window=ones(1,nfft),OverlapLength=nfft/2,FFTLength=nfft,FrequencyRange="onesided");
% the STFT scales like nfft/2, so we adjust the plotted value accordingly
surf(t1,f1,abs(s)/(nfft/2));
colorbar; ylim([0,100])
xlabel('Time (s)'); ylabel('Frequency (Hz)'); zlabel('Amplitude')
Plot the ampltude of the STFT for frequency=50, as a function of time. This means plotting the amplitude along the top ridge of the 3D plot above.
i50=find(f1==50); % frequency index for 50 Hz
y=abs(s(i50,:))/(nfft/2);
figure; plot(t1,y,'-rx')
grid on; xlabel('Time (s)'); ylabel('Amplitude'); title('Amplitude at 50 Hz')
Now let us find the time period when the signal amplitude is abnormal, and let us find the mean amplitude during that time period. The normal voltage is 325 V (=230V rms x sqrt(2)). Any amplitude that deviates by more than 3 V from the normal value will be considered abnormal.
t_ab=t1(abs(y-325)>3); % time vector for abnormal amplitude
tabs=min(t_ab); tabe=max(t_ab); % start, end times for abnormal amplitude
a_ab=mean(y(abs(y-325)>3)); % mean amplitude during abnormal time period
fprintf('Ampltude is abnormal from t=%.2f to %.2f.\n',tabs,tabe)
fprintf('Mean abnormal amplitude=%.1f.\n',a_ab)
The calculated mean abnormal ampltude = 277.4 V, not 275 V, because the first and last abnormal values are halfway between 325 V (normal) and 275 V (abnormal).
2 Comments
William Rose
on 20 Feb 2024
There are simpler ways to extract the amplitude, as a function of time, for the signal above, if you do not care about frequency. For example, you could compute the maximum value of the absolute value of the signal, using a moving window that you slide across the signal.
A Lumbini
on 21 Feb 2024
@William Rose Thank you for giving me this information
But, i am working on Wavelets So, i need to do it using wavelet transform
Alexander
on 20 Feb 2024
0 votes
Very easy approach with "movmax":
clear;clf;
fs=1000; % sampling rate [Hz]
ts=1/fs; % sampling rate [Hz]
t=0:0.001:5; % time [s]
tlength = length(t);
fl=50; % frequency [Hz] to look for
tl=1/50; % period time s
nl=round(tl/ts);
a=(1:tlength)/100; % amplitude
y=a.*sin(2*pi*fl*t); % signal
subplot(211)
plot(t,y);grid minor;
yAmp = movmax(y,nl); % look for maxima
subplot(211)
plot(t,y);grid minor;
subplot(212)
plot(t,yAmp);grid minor;
I think it can be done much more elegant by performing a FFT at a single frequency, in this case "fl". But I have to think about it and if I get a solution I will come back later.
3 Comments
William Rose
on 20 Feb 2024
I agree that movmax() is good, if you do not care about a specific frequency. The original poster specified interst in a 50 Hz signal. If there is broadband noise, movmax() will not do a good job of detecting the amplitude, or the change in amplitude, of the 50 Hz component of the signal. But stft() will do a good job - see answer above.
Alexander
on 21 Feb 2024
@William Rose you are right, but @A Lumbini was investigating the jitter of the amplitude, so I think he can work with this. but I will post another possible solution.
A Lumbini
on 21 Feb 2024
Thank you, your inputs are really helpful
But,i want to do it using wavelets
Alexander
on 21 Feb 2024
0 votes
Here is another approach as promissed above:
clear;clf;commandwindow;
fs=1000; % sampling rate [Hz]
ts=1/fs; % sampling rate [Hz]
t=0:ts:5-ts; % time 0-5 s
tlength = length(t);
fl=50; % frequency [Hz] to look for
tl=1/50; % period time s
N=round(tl/ts); % number of elements in period
a=ones(1,tlength)*4.5+rand(1,tlength); % amplitude 5 and +-0.5 random noise
y=a.*sin(2*pi*fl*t); % signal
yRes = reshape(y,N,[]); % Preparing the time slots
tN = (1:N)*ts;
Real = 2*fl/fs*sum(yRes.*sin(2*pi*fl*tN)');
Simag = 2*fl/fs*sum(yRes.*cos(2*pi*fl*tN)');
yARes = sqrt(Simag.^2+Real.^2); % vector of amplitudes
tRes = reshape(t,N,[]); % Preparing the time slots
tRes = mean(tRes); % t-axis according to the amplitudes in the middle of a time slot
plot(tRes,yARes);grid minor;
ylabel('Amplitude');xlabel('Time [sec]');
The advantage is, that the frequency fl is easy to adapt if used in (e.g.) order analysis of an engine runup.
Just an remainder: the equations Real and Simag are handmade (to be honest: I didn't found something appropriate in the net). Hence, check your results carefully.
@A Lumbini, If you must use wavelets, do a continuous wavelet transform (cwt). THe Matlab documentation is good.
Make a 50 Hz signal, 3 seconds long, with an amplitude dip in the middle:
fs=1000; % sampling rate (Hz)
t=0:1/fs:3-1/fs; % time (s)
f=50; % frequency (Hz)
a1=325; a2=275; % amplitude (V)
x=a1*sin(2*pi*f*t); % signal
x(t>=1 & t<2)=a2*sin(2*pi*f*t(t>=1 & t<2));
Do the cwt:
[cfs,f] = cwt(x,fs);
Plot the cwt:
surf(t,f,abs(cfs),Edgecolor='none')
xlabel("Time (s)"); ylabel("Frequency (Hz)");
ylim([0 200]); grid on; colorbar
title('C.W.T. of x(t)')
The surface plot, above, shows the dip in 50 Hz amplitude, from t=1 to t=2.
Make a 2D plot of amplitude versus time, for the 50 Hz component of the CWT:
[~,i50]=min(abs(f-50)); % frequency index closest to 50 Hz
y=abs(cfs(i50,:));
figure; plot(t,y,'-rx')
grid on; xlabel('Time (s)'); ylabel('Amplitude'); title('Amplitude at 50 Hz')
The plot above shows the ampltude of the 50 Hz component versus time. It shows the dip in amplitude from t=1 to t=2.
6 Comments
A Lumbini
on 23 Feb 2024
@William RoseThank you so much
Can we do it using Discrete Wavelet Transform?
William Rose
on 23 Feb 2024
The discrete wavelet transform is less suitable than the cwt for time-frequency analysis. See the discussion here. It includes this: "If your primary goal is a detailed time-frequency (scale) analysis or precise localization of signal transients, use cwt. "
The original goal was to detect the time when the 50 Hz amplitude droops. That requires time-frequency analysis. Use the cwt to do it.
A Lumbini
on 23 Feb 2024
@William Rose Mam, For example If we have a signal with fundamental frequency as 50Hz and harmonics in it.
If we want to find the amplitude of 50Hz freq i.e; fundamental component.
In this scenario is it possible to do it with DWT
William Rose
on 23 Feb 2024
Try it and see.
I have not tried it, but I doubt that the dwt will be great for identifying harmonics. I assume, based on your earlier comments, what you are looking for time-resolved harmonics, i.e. harmonics that change over time. I especially think the dwt won't be good for that, for the reasons mentioned above.
In the dwt, the scale factor, which correspnds more or less to frequency resolution, changes by factors of 2 only. In other words, you stretch and shrink the mother wavelet by factors of 2, only, when you do the dwt. That limits the frequency resolution of the dwt. That is why the dwt is said to have "one voice per octave" (an octave being a factor of 2 change in period or frequency).
A Lumbini
on 23 Feb 2024
William Rose
on 23 Feb 2024
@A Lumbini, you are welcome. Good luck with your work.
Categories
Find more on Continuous Wavelet Transforms in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!



