11.1 Spectrum analysis (SPA) - how window width impacts frequency resolution
The file ex_current.m illustrates spectral analysis of a current signal, which was measured when a transformer was started. The signal has a sampling interval of 1ms. One would suspect that the signal should contain several harmonics of the 50Hz base frequency. This is however not seen using the default parameters of the spa command. Find better parameters M and w for the command Phi = spa(z,M,w) where
- Phi = estimated spectral density
- z = iddata(y,) contains output y and input u from an experiment,
- M = number of lags in impulse response
- w = frequency vector for which Phi should be estimated
Use spectrumplot command to plot resulting spectrum
11.2 Spectrum Analysis using the System Identification GUI (SPA and ETFE)
Download the file vibrations.mat containing data from two experiments of a mechanical system having several vibrational modes. Import the data zh and zs into the GUI, and the model Gd (the true transfer function) and estimate spectral models using the menus "Estimate -> SPA" and "Estimate-> ETFE". Experiment with parameters "Frequencies" (the number of frequencies estimated) and "Frequency Resolution", a larger number gives better resolution, until you think you have found all resonances (located somewhere in the range [2 60] rad/s).
- How many resonant modes are there, and where are they located ?
- Which data set seems to give the best results, zh or zs?
- You should also learn how to save resulting figures from the GUI. For future use.
11.3 Empirical Transfer Function Estimation by FFT (EFTE)
The file ex_etfe.m estimates the 2nd order transfer function G(z) from noisy data . The sample rate is fs = 100Hz.
a) Use a chirp signal as input u(t) with frequency increasing logarithmically from f0=0.01 Hz to f1 = 50Hz. Study how the resulting ETFE depends on noise level . Make sure you understand how the ETFE is calculated using FFT of y and u.
b) Experiment with other input signals u such as
- u(t) = a1sin(w1t) + a2cos(w2 t)
- u(t) = step function
- u(t) = white noise
- u(t) = random binary sequence, using idinput
11.4 System Identification by correlation analysis (CRA)
The command g = cra(z,M,NA,plot) identifies a system using correlation analysis, where
- g = impulse response
- z = iddata(y,u) contains output y and input u from an experiment,
- M = number of lags in impulse response
- NA = order of an ARX whitening filter
- plot = 0,1 or 2 controls what should be plotted
The analysis is based on the fact that where * denotes convolution, and that when e is white noise.
a) Show that when the input is white noise then . (there was a typo here earlier, that has been corrected)
b) In matlab type "doc cra" and use the documentation to explain how first finding an ARX model of the input on the form makes it possible to use the result from a) to find the impulse response g (Hint: Use the whitening filter on the equation y(t) = G(q) u(t)
c) Study the file ex_cra.m that uses CRA. The code contains data from two experiments on the system . One with =white noise, the other with a signal which is not white.
- Is the true impulse response identified ?
- Vary the number of data points N. Does it improve accuracy ?
- Change NA = 0 to NA=1,2 or 5. This introduces the whitening filter A(q). Does the experiment with work better now ?
11.5 Acoustics models of two rooms (CRA)
Hint: Do NOT use the GUI for this problem.
The file ex_sound.m illustrates how to find a dynamical model of the acoustics of two rooms of different sizes. The audio source u(t) is sampled at 16kHz and the sound is recorded in two different rooms, giving the data y1 and y2. The sound data is available in file akustikdata.mat.
a) Plot y1 an y2. You might also be able to play the sound signal y running the command soundsc(y,16e3). Which data is recorded in the larger room (giving longer reverbarations) ?
b) Guess suitable parameters M (lags) and NA (order of whitening filter) for use in the CRA command and insert it in the file ex_sound.m. (You will need quite large values.) Plot the estimated impulse response g1hat and g2hat and compare with the true g1 and g2
c) Use your estimated impulse responses to simulate the two rooms' acoustics. This can be done with y1hat = filter(g1hat,1,u); How large must M and NA be for the simulated sound to be realistic ?
Note: Some problems and data on this exercise were taken from similar problems at LiU. Thanks.