File Exchange

image thumbnail

Spectral proper orthogonal decomposition (SPOD)

version (78.3 MB) by Oliver Schmidt
Implements the frequency domain form of proper orthogonal decomposition (POD)


Updated 21 May 2020

From GitHub

View Version History

View license on GitHub

SPOD() is a Matlab implementation of the frequency domain form of proper orthogonal decomposition (POD, also known as principle component analysis or Karhunen-Loève decomposition) called spectral proper orthogonal decomposition (SPOD). SPOD is derived from a space-time POD problem for stationary flows [1,2] and leads to modes that each oscillate at a single frequency. SPOD modes represent dynamic structures that optimally account for the statistical variability of stationary random processes.
The large-eddy simulation data provided along with this example is a subset of the database of a Mach 0.9 turbulent jet described in [3] and was calculated using the unstructured flow solver Charles developed at Cascade Technologies. If you are using the database in your research or teaching, please include explicit mention of Brès et al. [3]. The test database consists of 5000 snapshots of the symmetric component (m=0) of a round turbulent jet.
spod.m is a stand-alone Matlab function with no toolbox dependencies. All other Matlab files contained in this repository are related to the six examples that demonstrate the functionality of the code (see file descriptions below). A physical interpretation of the results obtained from the examples can be found in [4].

Repository zip file with examples (81.5 MB):
Matlab function only (15 KB):

File Description
spod.m - Spectral proper orthogonal decomposition in Matlab
example_1.m - Inspect data and plot SPOD spectrum
example_2.m - Plot SPOD spectrum and inspect SPOD modes
example_3.m - Specify spectral estimation parameters and use weighted inner product
example_4.m - Calculate the SPOD of large data and save results on hard drive
example_5.m - Calculate full SPOD spectrum of large data
example_6.m - Calculate and plot confidence intervals for SPOD eigenvalues
jet_data/getjet.m - Interfaces external data source with SPOD() (examples 4-5)
utils/trapzWeightsPolar.m - Integration weight matrix for cylindrical coordinates (examples 3-6)
utils/jetLES.mat - Mach 0.9 turbulent jet test database
LICENSE.txt - License

[1] Towne, A., Schmidt, O. T., Colonius, T., Spectral proper orthogonal decomposition and its relationship to dynamic mode decomposition and resolvent analysis, J. of Fluid Mech. 847, 821–867, 2018
[2] Lumley, J. L., Stochastic tools in turbulence, Academic Press, 1970
[3] G. A. Brès, P. Jordan, M. Le Rallic, V. Jaunet, A. V. G. Cavalieri, A. Towne, S. K. Lele, T. Colonius, O. T. Schmidt, Importance of the nozzle-exit boundary-layer state in subsonic turbulent jets, J. of Fluid Mech. 851, 83-124, 2018
[4] Schmidt, O. T., Towne, A., Rigas, G., Colonius, T., Bres, G. A., Spectral analysis of jet turbulence, J. of Fluid Mech. 855, 953–982, 2018

Cite As

Towne, A., Schmidt, O. T., Colonius, T., Spectral proper orthogonal decomposition and its relationship to dynamic mode decomposition and resolvent analysis, J. of Fluid Mech. 847, 821–867, 2018

Comments and Ratings (9)

Lican Wang

Chi Young Moon

Oliver Schmidt

Hi Christian,
Regarding your first question: increasing the frequency bins by reducing the number of snapshots per block is one option, as you mentioned (you're sacrificing frequency resolution doing that though). A better option is to use more data, if available/computable, to get more realizations of the Fourier transform. Regarding your second question: assuming that the sampling frequency is fixed, your only option here is to increase the number of snapshots per block. It really sounds like you'll need more data to get nicely converged modes...

Hi Oliver. Do you know any form to reduce the noise that could appear on the figures? I'm trying to calculate the first SPOD mode from some frequencies and I can obtain a cleaner figure from low frequencies instead from highest. I reduced the number of snapshots per block, reaching some clean figures but only for certain frequencies.

Also I'm trying to replicate some results that I have previously obtained using a DMD on my data, so I want to calculate the SPOD modes using an array with some frequencies that I founded relevant after calculate the DMD modes. It is possible to introduce manually some frequencies instead of those that the program calculate itself?

Thanks, Christian.

Oliver Schmidt

Have you tried out the examples? Each of the 6 examples builds on the previous one and introduces new features. The examples also show you how to visualize the results.

Hi Oliver, Do you have any tutorial to work with this program? I appreciate you if you can prepare a tutorial on youtube or anywhere.


Oliver Schmidt

Below this message is a code snippet that you can append at the end of example_2.m to animate the modes on screen.
Cheers, Oliver

%% Animate the same modes.
% Note how all wavepackets travel at approximately the same phase
% speed c_ph. The reason is that their streamwise wavenumber k_x changes
% with frequency such that c_ph = omega/k_x is approximately constant.
nt = 30;
T = 1/f(10); % period of the 10th frequency
time = linspace(0,T,nt); % animate over one period
count = 1;
for ti = 1:nt
for fi = [10 15 20]
for mi = [1 2]
pcolor(x,r,real(squeeze(P(fi,:,:,mi)*exp(2i*pi*f(fi)*time(ti))))), shading interp, axis equal tight, caxis(max(abs(caxis))*[-1 1])
xlabel('x'), ylabel('r'), title(['$f=' num2str(f(fi),'%.2f$') ', mode ' num2str(mi) ', $\lambda=' num2str(L(fi,mi),'%.2g$')])
xlim([0 10]); ylim([0 2])
count = count + 1;
hold on
hold off
count = 1;

Yonatan Cadavid

I wonder if is it possible to see the time evolution or oscilation of one coherent structure or several coherent structures with SPOD.

MATLAB Release Compatibility
Created with R2017b
Compatible with any release
Platform Compatibility
Windows macOS Linux

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!