Narrow band to 1/3 Octave, Sum or Average? 1

[yihshawn]

Lately I have came across a vibration design criteria with 1/3 Octave Band. Been doing a lots of dynamic analysis on narrowband but this 1/3 octave band is rather new to me. Some suggest we assign narrow band amplitude (fixed interval) into octave bin and sum them, some suggest to average them.

Gut feelings tells me it has something to do with Y-axis terms. So my questions are:

When we are using velocity rms, do we sum it or average it? and why?

 

[GregLocock]

You sum it because you are looking for the total power or energy in each 1/3 octave band.

In practice if you recalibrate using a narrowband calibrator then averaging will work, but it is still wrong, if you type in a calibration factor (V/EU) it'll come out all wrong.

Bear in mind third_octave_power (or energy) is proportional to sum V^2 for every narrowband bin in the third octave. I suggest you carefully work at least one example out by hand, the relationship between calibrations for V/EU, narrowband and broadband signals is actually straightforward, but can be a little confusing to start with. So long as you think about energy all the time you should be ok.





Cheers

Greg Locock


New here? Try reading these, they might help FAQ731-376

http://eng-tips.com/market.cfm?

 

[amanuensis]

The question is more complex than one might think.

First, you have to know if you work with a random or a deterministic signal.

For random signal, things are complicated...

Listen to me very well, it's very important : We can say nothing about ONE random signal!!
This is perfectly understandable since, as the signal is random, you can't predict what will be next values of the signal.

The only thing that you can do is : To compare your signal with another signal (This another signal can be the signal itself)

Since you need 2 signals, this is the reason why the unit of random signal is EU²/Hz (I don't explain here why).

So if the unit is EU²/Hz, you have to sum up over the whole 1/3 Octave band.

If this is a deterministic signal, you have to average, else unit wil not be consistent.

 

[yihshawn]

Greg,
By energy point of view, total energy in octave equals to total energy in narrowband.
Both my fft amplitude term and octave band amplitude need to be in rms velocity term (mms-1).

So:
1. square the fft lines (vi->vi^2)
2. assign them to bin and sum them, Sum(v1^2+ v2^2+...)
3. square root the bin value, sqrt(Sum (v1^2+ v2^2+...))
Is this the right track to go?


Ama,

I thought since we are using FFT, we always assume that the time domain signal repeats hence it is deterministic?
I am working on an environment assessment, there are deterministic harmonic sources and rather random impact sources (determined amplitude, random occurrence).
To work out the worst case, it would be easier to treat those random impact sources as deterministic.
However, I still don't understand why do we have to do average in converting narrowband deterministic signal into octave band.

 

[GregLocock]

1. square the fft lines (vi->vi^2) (now energy, remember to ignore phase)
2. assign them to bin and sum them, Sum(v1^2+ v2^2+...) (yes, still energy)
3. square root the bin value, sqrt(Sum (v1^2+ v2^2+...)) (back to rms, if that is what you want)

Either 2 or 3 could be your desired endpoints depending on whether you want to stay in the energy domain or back in linear units.








Cheers

Greg Locock


New here? Try reading these, they might help FAQ731-376

http://eng-tips.com/market.cfm?

 

[SomptingGuy]

Can't agree more with doing at least one calc by hand if that's possible.

If you can start with a (time domain) signal containing a few sine waves of known amplitudes at different frequencies (within the same 1/3 octave band), you can satisfy yourself that the RMS values you get from the different methods (time domain summation, frequency domain summation) are the same.

Also very important is getting your head around peak vs RMS, which you can also do if you are messing with a time series in a tool such as Matlab or similar. No real need for maths, you can just do actual means of squares of time values.

Steve

 

[yihshawn]

Did tests on a signal composed of sinusoidal waves of different phase in the same octave band (center 25Hz)

x(t)=2*sin(2pi*23*t)
y(t)=5*cos(2pi*24*t)
z(t)=15*cos(2pi*25*t+pi/5)

The sum of fft line values are 260, very close compared to hand calc, 2^2+5^2+15^2=254, a discrepancy ~ 2.4% in energy.

x(t)=2*sin(2pi*23*t)
y(t)=5*cos(2pi*24*t)
z(t)=8*cos(2pi*25*t+pi/5)

Gives a sum of 95, compared to hand calc, 93, ~2.5% error.
Gain some confidence here, thanks for the help given guys.


The error seems to be coming from FFT. Is this magnitude of error acceptable?
This is going far from the topic, but it is very interesting to know what is spurring the errors. I have tested sample size with fixed sampling internal (varying fft resolution), it does not affect the sum of square. Could it be inherent in FFT algorithm? Does it stay constant?

 

[SomptingGuy]

Here's my Matlab working. And it's the kind of thing I still do in Matlab, when I need to get back to basics. Note that my time vector is just the right length so that all signals are exactly periodic within the window. The n-1 is important here, don't want the start point duplicated.

n=1024; % Power of 2

t=(0:n-1)'/n; % Sample length is 1s, FFT will have 1Hz resolution

x=2*sin(2*pi*23*t);
y=5*cos(2*pi*24*t);
z=8*sin(2*pi*25*t+pi/5);

s=x+y+z;

fs=abs(fft(s)); % Lose phase info
famp=fs(1:n/2)*2/n; % Scale for amplitude
frms=fs(1:n/2)*2/n/sqrt(2); % Scale for RMS

famp(24:26) % amplitudes at 23,24,25 Hz
frms(24:26) % rms at 23,24,25 Hz

sqrt(mean(s.^2)) % RMS of whole signal in time domain
sqrt(sum(frms.^2)) % RMS of whole signal in frequency domain


= = = = =

The console output:

>> tmp

ans =

2.0000
5.0000
8.0000


ans =

1.4142
3.5355
5.6569


ans =

6.8191


ans =

6.8191

 

[GregLocock]

yihshawn - that's a bigger error than I'd have expected, what happens when you IFFT the signal?

Incidentally, for these tests, don't use windowing.

Cheers

Greg Locock


New here? Try reading these, they might help FAQ731-376

http://eng-tips.com/market.cfm?

 

[yihshawn]

No window is used.

I have been using this code modified from Embedded Lab. I found out If we choose the sample size in the power of 2 it gives very small error (~0.1%). If we not, eg. 1000, then there is still 1024 bins and ifft go flat after the 1000th point, the error grows.

The other weird stuff is that the frequency vector f1 (or amplitude vector Xsig) is (N/2)+1 x 1 vector unlike Somptinguys's, which is a N/2 x 1 vector. Should I remove the +1 in the colored code?



%% Author :- Embedded Laboratory

fSampling = 1000; %Sampling Frequency
tSampling = 1/fSampling; %Sampling Time
L = 32768; %Sample size N
t = (0:L-1)*tSampling; %Time Vector

xsig = xlsread('Octavetest.xls','Sheet1','E1:E32768');

%% Original signal
plot(t,xsig);
grid on;
axis([0 5 -250 250])
pause(2);

%%Frequency Transform of above Signal
subplot(2,1,2)
NFFT = 2^nextpow2(L);
Xsig = fft(xsig,NFFT)/L;
AXsig = 2*abs(Xsig(1:NFFT/2+1));
f1 = fSampling/2*(linspace(0,1,NFFT/2+1));
bar(f1,2*abs(Xsig(1:NFFT/2+1)),'r');
grid on;
axis([15 40 0 40])
pause(2);

figure(2);
inv=ifft(Xsig);
plot(inv)

 

[GregLocock]

FFT relies on 2^n datapoints. If you want to use a different number then use some other fourier transform. Matlab is being a bit too clever letting you use 1000, the results will vary depending on the assumptions they made(zero padding I suspect).

I don't know why your spectrum would have a different number of bins, Matlab's documentation is usually thorough if not always clear.

I would /always/ generate a sine signal and run it through an unknown fourier transform, to make sure I understood what amplitude scaling they had used, and what they have used in terms of bins and frequency.





Cheers

Greg Locock


New here? Try reading these, they might help FAQ731-376

http://eng-tips.com/market.cfm?

 

[SomptingGuy]

... there is an implied +1 in my code snippet for the frequency axis:

"famp(24:26) % amplitudes at 23,24,25 Hz"

First bin is 0Hz (but need to scale it only by 1/N, not 2/N). I ignored it, since it was zero in your test signal. I always forget any magic formula for the frequency axis and just remember that the first (zeroth) bin is DC, the next bin corresponds to the sample length and so on (easier than juggling with sampling rate and number of points). How many frequency values you compute really depends on what your FFT function does around and after the Nyquist frequency.


Steve

 

[yihshawn]

Thanks for your brainpower Steve and Greg!

Looking back at your comments after a few months of learning and practice I found them extremely helpful.

The error mentioned previously could be caused by spectral leakage where the signal block is not an integer period of its frequency component contained. When a perfect sample length is designed as what Steve did there is no discrepancy in total energy.

However in a case of monitoring environment vibration, which the frequency component are more random, it is almost impossible to look for a perfect sample length. As I have read in many textbook or lecture materials, windowing is often used to reduce spectral leakage in harmonic signal case. I am wondering if it can be used before octave band analysis.

*Updated: No window has better accuracy then Hanning window.

 

[SomptingGuy]

So in the Real World, we don't have perfectly periodic test signals and we need to use windowing to minimise the effects.

Realise that a window will remove some of your signal, so the spectrum of the windowed signal will need to be multiplied by some "correction factor" to restore what's been lost. This is where you need to make a choice based on what you are going to do with the spectrum. If you want to present it with the peaks close to the amplitudes in your input signal, use "Amplitude" correction. If you want the energy to be preserved (e.g. if you are summing across frequency bands), use "Energy" correction.

It's quite easy to derive the correction factors of a window: Amplitude correction is 1/mean(w), Energy correction is 1/rms(w), where w is the window. For a Hanning window, this gives 2.0 for amplitude and sqrt(8/3) for energy.

Also be sure to choose the best window for your application. I've always used flattop with amplitude correction for grabbing the amplitude of a calibration signal, hanning with energy correction for 1/3 octave summation.

HTH

Steve

 

[GregLocock]

As a guru once observed, if you know so much about your signal that you can select a better window than Hanning, then you probably don't need a window at all.

One of his more inscrutable comments at first sight!


Cheers

Greg Locock


New here? Try reading these, they might help FAQ731-376

http://eng-tips.com/market.cfm?

 

[yihshawn]

Can't agree more! As I learn more on this topic I found I know less than I imagined.

Upon testing on a signal composed of 1, 10 , 50 Hz, using hanning window with the correct energy correction give an excellent match on high freq band (10, 50), where leakage on low frequency band (1) persists.

If only comparing very low frequencies band, it appears that the side lobes can still slip to the neighbor band. No window at these low frequencies band is still doing better than hanning.

fsampling = 1000
no. of sample = 4096

Perhaps the right way to do it is to increase the sampling rate in the next measurement.

 

[SomptingGuy]

... you'll probably find that increasing the length of the measurement block will have an effect on the lower frequencies, whereas changing the sampling rate will not.

Steve

 

[yihshawn]

You are right Steve, in this problem the time block is constraint between 2s-8s. Have a limited parameter to play with...