stft
Short-time Fourier transform
Syntax
Description
returns the
Short-Time Fourier Transform (STFT) of
s
= stft(x
)x
.
specifies additional options using name-value pair arguments. Options include the FFT
window and length. These arguments can be added to any of the previous input
syntaxes.s
= stft(___,Name,Value
)
stft(___)
with no output arguments plots the
magnitude of the STFT in the current figure window.
Examples
3D STFT Visualization
Generate two seconds of a voltage controlled oscillator output, controlled by a sinusoid sampled at 10 kHz.
fs = 10e3; t = 0:1/fs:2; x = vco(sin(2*pi*t),[0.1 0.4]*fs,fs);
Compute and plot the STFT of the signal. Use a Kaiser window of length 256 and shape parameter . Specify the length of overlap as 220 samples and DFT length as 512 points. Plot the STFT with default colormap and view.
stft(x,fs,'Window',kaiser(256,5),'OverlapLength',220,'FFTLength',512)
Change the view to display the STFT as a waterfall plot. Set the colormap to jet
.
view(-45,65)
colormap jet
STFT of Quadratic Chirp
Generate a quadratic chirp sampled at 1 kHz for 2 seconds. The instantaneous frequency is 100 Hz at and crosses 200 Hz at second.
ts = 0:1/1e3:2; f0 = 100; f1 = 200; x = chirp(ts,f0,1,f1,'quadratic',[],'concave');
Compute and display the STFT of the quadratic chirp with a duration of 1 ms.
d = seconds(1e-3); win = hamming(100,'periodic'); stft(x,d,'Window',win,'OverlapLength',98,'FFTLength',128);
STFT Frequency Ranges
Generate a signal sampled at 5 kHz for 4 seconds. The signal consists of a set of pulses of decreasing duration separated by regions of oscillating amplitude and fluctuating frequency with an increasing trend. Plot the signal.
fs = 5000; t = 0:1/fs:4-1/fs; x = besselj(0,600*(sin(2*pi*(t+1).^3/30).^5)); plot(t,x)
Compute the one-sided, two-sided, and centered short-time Fourier transforms of the signal. In all cases, use a 202-sample Kaiser window with shape factor to window the signal segments. Display the frequency range used to compute each transform.
rngs = ["onesided" "twosided" "centered"]; for kj = 1:length(rngs) opts = {'Window',kaiser(202,10),'FrequencyRange',rngs(kj)}; [~,f] = stft(x,fs,opts{:}); subplot(length(rngs),1,kj) stft(x,fs,opts{:}) title(sprintf('''%s'': [%5.3f, %5.3f] kHz',rngs(kj),[f(1) f(end)]/1000)) end
Repeat the computation, but now change the length of the Kaiser window to 203, an odd number. The 'twosided'
frequency interval does not change. The other two frequency intervals become open at the higher end.
for kj = 1:length(rngs) opts = {'Window',kaiser(203,10),'FrequencyRange',rngs(kj)}; [~,f] = stft(x,fs,opts{:}); subplot(length(rngs),1,kj) stft(x,fs,opts{:}) title(sprintf('''%s'': [%5.3f, %5.3f] kHz',rngs(kj),[f(1) f(end)]/1000)) end
STFT of Multichannel Signals
Generate a three-channel signal consisting of three different chirps sampled at 1 kHz for one second.
The first channel consists of a concave quadratic chirp with instantaneous frequency 100 Hz at t = 0 and crosses 300 Hz at t = 1 second. It has an initial phase equal to 45 degrees.
The second channel consists of a convex quadratic chirp with instantaneous frequency 100 Hz at t = 0 and crosses 500 Hz at t = 1 second.
The third channel consists of a logarithmic chirp with instantaneous frequency 300 Hz at t = 0 and crosses 500 Hz at t = 1 second.
Compute the STFT of the multichannel signal using a periodic Hamming window of length 128 and an overlap length of 50 samples.
fs = 1e3; t = 0:1/fs:1-1/fs; x = [chirp(t,100,1,300,'quadratic',45,'concave'); chirp(t,100,1,500,'quadratic',[],'convex'); chirp(t,300,1,500,'logarithmic')]'; [S,F,T] = stft(x,fs,'Window',hamming(128,'periodic'),'OverlapLength',50);
Visualize the STFT of each channel as a waterfall plot. Control the behavior of the axes using the helper function helperGraphicsOpt
.
waterfall(F,T,abs(S(:,:,1))') helperGraphicsOpt(1)
waterfall(F,T,abs(S(:,:,2))') helperGraphicsOpt(2)
waterfall(F,T,abs(S(:,:,3))') helperGraphicsOpt(3)
This helper function sets the appearance and behavior of the current axes.
function helperGraphicsOpt(ChannelId) ax = gca; ax.XDir = 'reverse'; ax.ZLim = [0 30]; ax.Title.String = ['Input Channel: ' num2str(ChannelId)]; ax.XLabel.String = 'Frequency (Hz)'; ax.YLabel.String = 'Time (seconds)'; ax.View = [30 45]; end
Input Arguments
x
— Input signal
vector | matrix | timetable
Input signal, specified as a vector, a matrix, or a MATLAB®
timetable
.
Note
If you invert s
using istft
and want the result to be the same length as x
, the value of
(length(x)-noverlap)/(length(window)-noverlap)
must be an
integer.
If the input has multiple channels, specify
x
as a matrix where each column corresponds to a channel.For timetable input,
x
must contain uniformly increasing finite row times. If a timetable has missing or duplicate time points, you can fix it using the tips in Clean Timetable with Missing, Duplicate, or Nonuniform Times.For multichannel timetable input, specify
x
as a timetable with a single variable containing a matrix or a timetable with multiple variables each containing a column vector. All variables must have the same precision.
Each channel of x
must have a length greater than or
equal to the window length.
Example: chirp(0:1/4e3:2,250,1,500,'quadratic')
specifies a
single-channel chirp.
Example: timetable(rand(5,2),'SampleRate',1)
specifies a
two-channel random variable sampled at 1 Hz for 4 seconds.
Example: timetable(rand(5,1),rand(5,1),'SampleRate',1)
specifies a
two-channel random variable sampled at 1 Hz for 4 seconds.
Data Types: double
| single
Complex Number Support: Yes
fs
— Sample rate
2π
(default) | positive scalar
Sample rate, specified as a positive scalar. This argument applies only when
x
is a vector or a matrix.
Data Types: double
| single
Name-Value Arguments
Specify optional pairs of arguments as
Name1=Value1,...,NameN=ValueN
, where Name
is
the argument name and Value
is the corresponding value.
Name-value arguments must appear after other arguments, but the order of the
pairs does not matter.
Before R2021a, use commas to separate each name and value, and enclose
Name
in quotes.
Example: 'Window',hamming(100),'OverlapLength',50,'FFTLength',128
windows the data using a 100-sample Hamming window, with 50 samples of overlap between
adjoining segments and a 128-point FFT.
Window
— Spectral window
hann(128,'periodic')
(default) | vector
Spectral window, specified as a vector. If you do not specify the window or
specify it as empty, the function uses a Hann window of length 128. The length of
'Window'
must be greater than or equal to 2.
For a list of available windows, see Windows.
Example: hann(N+1)
and
(1-cos(2*pi*(0:N)'/N))/2
both specify a Hann window of length
N
+ 1.
Data Types: double
| single
OverlapLength
— Number of overlapped samples
75%
of window length (default) | nonnegative integer
Number of overlapped samples, specified as a positive integer smaller than the
length of 'Window'
. If you omit
'OverlapLength'
or specify it as empty, it is set to the
largest integer less than 75% of the window length, which is 96 samples for the
default Hann window.
Data Types: double
| single
FFTLength
— Number of DFT points
128
(default) | positive integer
Number of DFT points, specified as a positive integer. The value must be greater than or equal to the window length. If the length of the input signal is less than the DFT length, the data is padded with zeros.
Data Types: double
| single
FrequencyRange
— STFT frequency range
'centered'
(default) | 'twosided'
| 'onesided'
STFT frequency range, specified as 'centered'
,
'twosided'
, or 'onesided'
.
'centered'
— Compute a two-sided, centered STFT. If'FFTLength'
is even, thens
is computed over the interval (–π, π] rad/sample. If'FFTLength'
is odd, thens
is computed over the interval (–π, π) rad/sample. If you specify time information, then the intervals are (–fs, fs/2] cycles/unit time and (–fs, fs/2) cycles/unit time, respectively, where fs is the effective sample rate.'twosided'
— Compute a two-sided STFT over the interval [0, 2π) rad/sample. If you specify time information, then the interval is [0, fs) cycles/unit time.'onesided'
— Compute a one-sided STFT. If'FFTLength'
is even, thens
is computed over the interval [0, π] rad/sample. If'FFTLength'
is odd, thens
is computed over the interval [0, π) rad/sample. If you specify time information, then the intervals are [0, fs/2] cycles/unit time and [0, fs/2) cycles/unit time, respectively, where fs is the effective sample rate. This option is valid only for real signals.Note
When this argument is set to
'onesided'
,stft
outputs the values in the positive Nyquist range and does not conserve the total power.
For an example, see STFT Frequency Ranges.
Data Types: char
| string
OutputTimeDimension
— Output time dimension
'acrosscolumns'
(default) | 'downrows'
Output time dimension, specified as 'acrosscolumns'
or
'downrows'
. Set this value to 'downrows'
if
you want the time dimension of s
down the rows and the frequency
dimension across the columns. Set this value to 'acrosscolumns'
if
you want the time dimension of s
across the columns and the
frequency dimension down the rows. This input is ignored if the function is called
without output arguments.
Output Arguments
s
— Short-time Fourier transform
matrix | 3-D array
Short-time Fourier transform, returned as a matrix or a 3-D array. Time increases
across the columns of s
and frequency increases down the rows. The
third dimension, if present, corresponds to the input channels.
If the signal
x
has Nx time samples, thens
has k columns, where k = ⌊(Nx–L)/(M–L)⌋, M is the length of'Window'
, L is the'OverlapLength'
, and the ⌊ ⌋ symbols denote the floor function.The number of rows in
s
is equal to the value specified in'FFTLength'
.
Data Types: double
| single
f
— Frequencies
vector
Frequencies at which the STFT is evaluated, returned as a vector.
Data Types: double
| single
t
— Time instants
vector
Time instants, returned as a vector. t
contains the time values
corresponding to the centers of the data segments used to compute short-time power
spectrum estimates.
Data Types: double
| single
More About
Short-Time Fourier Transform
The short-time Fourier transform (STFT) is used to analyze how the frequency content of a nonstationary signal changes over time.
The STFT of a signal is calculated by sliding an analysis window of length over the signal and calculating the discrete Fourier transform of the windowed data. The window hops over the original signal at intervals of samples. Most window functions taper off at the edges to avoid spectral ringing. If a nonzero overlap length is specified, overlap-adding the windowed segments compensates for the signal attenuation at the window edges. The DFT of each windowed segment is added to a matrix that contains the magnitude and phase for each point in time and frequency. The number of columns in the STFT matrix is given by
where is the length of the original signal and the ⌊⌋ symbols denote the floor function. The number of rows in the matrix equals NDFT, the number of DFT points, for centered and two-sided transforms and ⌊NDFT/2⌋ + 1 for one-sided transforms.
The STFT matrix is given by such that the th element of this matrix is
where
— Window function of length .
— DFT of windowed data centered about time .
— Hop size between successive DFTs. The hop size is the difference between the window length and the overlap length .
The magnitude squared of the STFT yields the spectrogram
representation of the power spectral density of the function.
Perfect Reconstruction
In general, computing the STFT of an input signal and inverting it does not result in perfect reconstruction. If you want the output of ISTFT to match the original input signal as closely as possible, the signal and the window must satisfy the following conditions:
Input size — If you invert the output of
stft
usingistft
and want the result to be the same length as the input signalx
, the value of must be an integer.COLA compliance — Use COLA-compliant windows, assuming that you have not modified the short-time Fourier transform of the signal.
Padding — If the length of the input signal is such that the value of k is not an integer, zero-pad the signal before computing the short-time Fourier transform. Remove the extra zeros after inverting the signal.
References
[1] Mitra, Sanjit K. Digital Signal Processing: A Computer-Based Approach. 2nd Ed. New York: McGraw-Hill, 2001.
[2] Sharpe, Bruce. Invertibility of Overlap-Add Processing. https://gauss256.github.io/blog/cola.html, accessed July 2019.
[3] Smith, Julius Orion. Spectral Audio Signal Processing. https://ccrma.stanford.edu/~jos/sasp/, online book, 2011 edition, accessed Nov 2018.
Extended Capabilities
Tall Arrays
Calculate with arrays that have more rows than fit in memory.
Usage notes and limitations:
'OutputTimeDimension'
must always be specified and set to
'downrows'
.
For more information, see Tall Arrays.
C/C++ Code Generation
Generate C and C++ code using MATLAB® Coder™.
Usage notes and limitations:
Timetables are not supported for code generation.
GPU Code Generation
Generate CUDA® code for NVIDIA® GPUs using GPU Coder™.
Usage notes and limitations:
Timetables are not supported for code generation.
Thread-Based Environment
Run code in the background using MATLAB® backgroundPool
or accelerate code with Parallel Computing Toolbox™ ThreadPool
.
Usage notes and limitations:
The syntax with no output arguments is not supported.
For more information, see Run MATLAB Functions in Thread-Based Environment.
GPU Arrays
Accelerate code by running on a graphics processing unit (GPU) using Parallel Computing Toolbox™.
This function fully supports GPU arrays. For more information, see Run MATLAB Functions on a GPU (Parallel Computing Toolbox).
Version History
Beispiel öffnen
Sie haben eine geänderte Version dieses Beispiels. Möchten Sie dieses Beispiel mit Ihren Änderungen öffnen?
MATLAB-Befehl
Sie haben auf einen Link geklickt, der diesem MATLAB-Befehl entspricht:
Führen Sie den Befehl durch Eingabe in das MATLAB-Befehlsfenster aus. Webbrowser unterstützen keine MATLAB-Befehle.
Select a Web Site
Choose a web site to get translated content where available and see local events and offers. Based on your location, we recommend that you select: .
You can also select a web site from the following list:
How to Get Best Site Performance
Select the China site (in Chinese or English) for best site performance. Other MathWorks country sites are not optimized for visits from your location.
Americas
- América Latina (Español)
- Canada (English)
- United States (English)
Europe
- Belgium (English)
- Denmark (English)
- Deutschland (Deutsch)
- España (Español)
- Finland (English)
- France (Français)
- Ireland (English)
- Italia (Italiano)
- Luxembourg (English)
- Netherlands (English)
- Norway (English)
- Österreich (Deutsch)
- Portugal (English)
- Sweden (English)
- Switzerland
- United Kingdom (English)