Use of PSD functions to obtain FRFs

Hi,

I am doing some modal analysis on a steel structure and have encountered some troubling issues while using Matlab's Power Spectral Density functions to obtain the Frequency Response Functions of the structure.

I have obtained a force vector (in Newtons) from an impact hammer, denoted 'y', and an acceleration vector (in g's) from an accelerometer, denoted 'x'. These vectors have been trimmed to a length 'L' of 153600, equating to 30 seconds of data at a sampling rate of 5120 Hz. Having already applied my own windowing function to the data, I am using the rectwin function in Matlab. The data has been padded with zeros to a length of 'NFFT' = 262144 as per instructions in the Matlab help file.

In order to obtain the FRFs, I have used the periodogram and cpsd functions as follows:

[Pxx,f]=periodogram(x,rectwin(L),NFFT,Fs); % autospectrum of the response

[Pff,f]=periodogram(y,rectwin(L),NFFT,Fs); % autospectrum of the impact

[Pfx,f]=cpsd(y,x,rectwin(L),0,NFFT,Fs); % cross spectrum

[Pxf,f]=cpsd(x,y,rectwin(L),0,NFFT,Fs); % cross spectrum

H1=Pfx./Pff; % FRF estimate 1

H2=Pxx./Pxf; % FRF estimate 2

Coherence=H1./H2; % This is the coherence function

I have noticed two problems with these FRF estimates that do not match with the theory. Firstly, the coherence function obtained by this method is exactly equal to unity (i.e., 1.00000) across the entire frequency spectrum. This shouldn't be possible! I would expect to see the coherence values drop in the vicinity of anti-resonance. The second problem I have noticed is that the FRF estimates appear to be rotated 90 degrees in phase. This is particularly noticeable when plotting the Real vs. Imaginary portions of the FRF near resonance. Theory tells me that I should see a circle centred on the imaginary axis, however, all my plots are showing a circle centred on the real axis.

This is all very annoying because I don't know what I have done wrong.

Please advise,

Craig Cowled.

Answer 1

[Wayne King]

Hi Craig, one problem is that you use one window for your coherence estimate that is the entire length of your data. This will always result in a coherence estimate of unity at all frequencies.

That is not just a MATLAB issue, that is a known result. You always need more than one bivariate sample to get a meaningful coherence estimate.

If you use mscohere (or if you want to keep it doing it your way with cpsd), then you must use overlapped segments. And, with these overlapped, windowed segments, I would recommend against using a rectangular window. You can use the default Hamming, that would be fine.

[Craig Cowled]

Wayne, thank you very much for your guidance. As you have recommended, I played around with various window lengths and number of overlaps. The periodogram function does not allow overlaps, so I have used the pwelch function for autospectra. Because the frequencies of interest are between 5 - 25Hz, I need to maximise the window length as much as possible. As I mentioned previously, my experiment involved impact testing, so I have to use the rectwin function as the hamming window eliminates most of my useful data (i.e., at the beginning of the signal). After some trial and error, I have adopted a window length that of 76800, which is half the data length. The number of overlaps does not appear to affect the results very much. Using the following syntax, I have been able to obtain meaningful coherence plots:

[Pxx,f]=pwelch(x,rectwin(76800),16,NFFT,Fs);

[Pff,f]=pwelch(y,rectwin(76800),16,NFFT,Fs);

[Pfx,f]=cpsd(y,x,rectwin(76800),16,NFFT,Fs);

[Pxf,f]=cpsd(x,y,rectwin(76800),16,NFFT,Fs);

Thank you very much for helping me to solve the problem of the coherence plots. I am, however, still getting the same problem with the Nyquist plots of the FRFs near resonance. The circles are still centred on the real axis. Any thoughts on this issue?

Regards,

Craig Cowled.

Answer 2

[Craig Cowled]

I thought I should tie up a loose end with this question. I discovered that the data acquisition from my impact hammer was faulty. The hammer had an IEPE sensor which requires a constant current of 2mA. In my initial testing, I was taking a passive reading of voltage from the sensor, which explains why my FRFs were out by 90 degrees. Once I applied the correct settings in my data acquisition, the FRFs matched with the theory.