timefrequency.jl
This unit implements time-frequency representations based on analytic signal estimations. It uses the filters.jl unit for filtering the signals in a filter bank and the hilbert.jl unit for computing the Hilbert transform.
The main object in the time-frequency domain is the TFAnalyticSignal, from which the TFAmplitude and TFPhase objects are derived. Taken together these three objects are referred to as TFobjects and vectors of them are a TFobjectsVector type.
Object | Vector type | Constructors from data |
---|---|---|
TFAnalyticSignal | TFAnalyticSignalVector | TFanalyticsignal |
TFAmplitude | TFAmplitudeVector | TFamplitude |
TFPhase | TFPhaseVector | TFphase |
TFAnalyticSignal
Categories: data objects, TFobjects.
An analytic signal object has the following structure:
struct TFAnalyticSignal
y :: Matrix{T} where T<:Complex
bandwidth :: IntOrReal
flabels :: Vector{S} where S<:Real
nonlinear :: Bool
fsmoothing :: Smoother
tsmoothing :: Smoother
end
Fields:
y
: a complex matrix holding the analytic signal in the time-frequency domain, with frequency band-pass regions (in Hz) in rows and time (samples) in columns.
bandwidth
: the bandwidth (in Hz) of the filter bank band-pass regions. See constructor TFanalyticsignal
for details. It can be an integer or a real number.
flabels
: a vector holding all center frequencies (in Hz) of the filter bank band-pass regions. Those are the frequency labels for the rows of y
.
non-linear
: a flag indicating whether the analytic signal has been normalized so as to have amplitude$1.0$ at all points. Such normalization allows non-linear univariate and bivariate estimations (see timefrequencyuni.jl and timefrequencybi.jl).
fsmoothing
: a flag of the Smoother type indicating whether the analytic signal has been smoothed across adjacent frequency band-pass regions. If no frequency smoothing has been applied, it is equal to noSmoother
, which is the default in all constructors.
tsmoothing
: a flag of the Smoother type indicating whether the analytic signal has been smoothed across adjacent samples (time). If no time smoothing has been applied, it is equal to noSmoother
, which is the default in all constructors.
Note: In Julia the fields are accessed by the usual dot notation, e.g., you may verify that for TFAnalyticSignal object Y
,
length(Y.flabels) == dim(Y.z, 1)
A vector of TFAnalyticSignal objects is of type TFAnalyticSignalVector.
Methods for TFAnalyticSignal and TFAnalyticSignalVector objects
method | TFAnalyticSignal | TFAnalyticSignalVector |
---|---|---|
amplitude | ✔ | ✔ |
phase | ✔ | ✔ |
polar | ✔ | |
extract | ✔ | ✔ |
mean | ✔ | ✔ |
smooth | ✔ | ✔ |
sameParams | ✔ | |
isLinear | ✔ | |
isNonLinear | ✔ | |
meanAmplitude | ✔ | |
concentration | ✔ | |
comodulation | ✔ | |
meanDirection | ✔ | |
coherence | ✔ |
Generic Constructors:
In order to construct a TFAnalyticSignal object from univariate data, FourierAnalysis provides two TFanalyticsignal
constuctors, which is what you will use in practice most of the time.
Manual constructors are also possible, for which you have to provide appropriate arguments. The default manual constructor of TFAnalyticSignal objects is
TFAnalyticSignal(y, bandwidth, flabels, nonlinear, fsmoothing, tsmoothing).
No other generic constructor is provided for this object.
TFAmplitude
Categories: data objects, TFobjects.
An amplitude object has the following structure:
struct TFAmplitude
y :: Matrix{T} where T<:Real
bandwidth :: IntOrReal
flabels :: Vector{S} where S<:Real
fsmoothing :: Smoother
tsmoothing :: Smoother
func :: Function
end
Fields:
y
: a real matrix holding the amplitude (modulus) of an analytic signal in the time-frequency domain, often referred to as the envelope, with frequency band-pass regions (in Hz) in rows and time (samples) in columns.
bandwidth
: the bandwidth (in Hz) of the filter bank band-pass regions. See constructor TFanalyticsignal
for details. It can be an integer or a real number.
flabels
: a vector holding all center frequencies (in Hz) of the filter bank band-pass regions. Those are the frequency labels for the rows of y
.
fsmoothing
: a flag of the Smoother type indicating whether the amplitude has been smoothed across adjacent frequency band-pass regions. If no frequency smoothing has been applied, it is equal to noSmoother
, which is the default in all constructors.
tsmoothing
: a flag of the Smoother type indicating whether the amplitude has been smoothed across adjacent samples (time). If no time smoothing has been applied, it is equal to noSmoother
, which is the default in all constructors.
Note: Smoothing flags in this object indicate that the amplitude has been smoothed, whereas in the TFAnalyticSignal object indicate that the analytic signal has been smoothed. Note that it is not equivalent to obtain amplitude from smoothed analytic signal (e.g., using the amplitude
function) or to smooth the amplitude of analytic signal, e.g., using the TFamplitude
constructor.
func
: a name of a function that has been applied element-wise to the matrix .y
holding the amplitude. All constructors from data by default use the identity
(do nothing) function.
A vector of TFAmplitude objects is of type TFAmplitudeVector.
Methods for TFAmplitude and TFAmplitudeVector objects
method | TFAmplitude | TFAmplitudeVector |
---|---|---|
extract | ✔ | ✔ |
mean | ✔ | ✔ |
smooth | ✔ | ✔ |
sameParams | ✔ | |
isLinear | ✔ | |
isNonLinear | ✔ | |
meanAmplitude | ✔ | |
comodulation | ✔ |
Generic Constructors:
In order to construct a TFAmplitude
object from univariate data, FourierAnalysis provides four TFamplitude
constuctors, which is what you will use in practice most of the time.
Manual constructors are also possible, for which you have to provide appropriate arguments. The default manual constructor of TFAmplitude objects is
TFAmplitude(y, bandwidth, flabels, fsmoothing, tsmoothing, func).
Other generic constructors are also provided:
TFAmplitude(y, bandwidth, flabels, fsmoothing, tsmoothing)
enables construction giving only y
, bandwidth
, flabels
, fsmoothing
and tsmoothing
. func
is set automatically to identity
;
TFAmplitude(y, bandwidth, flabels)
acts like the constructor above, but sets by default also both fsmoothing
and tsmoothing
to noSmoother
.
TFPhase
Categories: data objects, TFobjects.
A phase object has the following structure:
struct TFPhase
y :: Matrix{T} where T<:Real
bandwidth :: IntOrReal
flabels :: Vector{S} where S<:Real
nonlinear :: Bool
fsmoothing :: Smoother
tsmoothing :: Smoother
unwrapped :: Bool
func :: Function
end
Fields:
y
: a real matrix holding the phase (argument) of an analytic signal in the time-frequency domain, with frequency band-pass regions (in Hz) in rows and time (samples) in columns. By default the phase is represented in $[−π, π]$.
bandwidth
: the bandwidth (in Hz) of the filter bank band-pass regions. See constructor TFanalyticsignal
for details. It can be an integer or a real number.
flabels
: a vector holding all center frequencies (in Hz) of the filter bank band-pass regions. Those are the frequency labels for the rows of y
.
non-linear
: a flag indicating whether the phase has been estimated from analytic signal normalized so as to have amplitude$1.0$ at all points. Such normalization allows non-linear univariate and bivariate estimations (see timefrequencyuni.jl and timefrequencybi.jl).
fsmoothing
: a flag of Smoother indicating whether the phase has been smoothed across adjacent frequency band-pass regions. If no frequency smoothing has been applied, it is equal to noSmoother
, which is the default in all constructors.
tsmoothing
: a flag of the Smoother type indicating whether the phase has been smoothed across adjacent samples (time). If no time smoothing has been applied, it is equal to noSmoother
, which is the default in all constructors.
Note: Smoothing raw phase estimates is unappropriate since the phase is a discontinous function, however it makes sense to smooth phase if the phase is unwrapped (see below).
Note: Smoothing flags in this object indicate that the (unwrapped) phase has been smoothed, whereas in the TFAnalyticSignal object indicate that the analytic signal has been smoothed. Note that it is not equivalent to obtain unwrapped phase from smoothed analytic signal (e.g., using the phase
function) or to smooth the unwrapped phase of analytic signal, e.g., using the TFphase
constructor.
unwrapped
: a flag indicating if the phase has been unwrapped. The unwrapped phase is defined as the cumulative sum of the phase along the time dimension once the phase is represented in $[0, 2π]$.
func
: a name of a function that has been applied element-wise to the matrix .y
holding the phase. All constructors from data by default use the identity
(do nothing) function. Examples of possible functions:
func=x->x+π
return the phase in [0, 2π],func=x->x/π
return the phase in [-1, 1],func=sin
return the sine of the phase.
A vector of TFPhase objects is of type TFPhaseVector.
Methods for TFPhase and TFPhaseVector objects
method | TFAmplitude | TFAmplitudeVector |
---|---|---|
unwrapPhase | ✔ | ✔ |
isUnwrapped | ✔ | ✔ |
extract | ✔ | ✔ |
mean | ✔ | ✔ |
smooth | ✔ | ✔ |
sameParams | ✔ | |
isLinear | ✔ | |
isNonLinear | ✔ |
Generic Constructors:
In order to construct a TFPhase object from univariate data, FourierAnalysis provides four TFphase
constuctors, which is what you will use in practice most of the time.
Manual constructors are also possible, for which you have to provide appropriate arguments. The default manual constructor of TFPhase objects is
TFPhase(y, bandwidth, flabels, nonlinear,
fsmoothing, tsmoothing, unwrapped, func).
Other generic constructors are also provided:
TFPhase(y, bandwidth, flabels, nonlinear, fsmoothing, tsmoothing)
enables construction giving only y
, bandwidth
, flabels
, fsmoothing
and tsmoothing
. unwrapped
is set to false
and func
is set to identity
;
TFPhase(y, bandwidth, flabels)
acts like the constructor above, but sets by default also nonlinear
to true
and both fsmoothing
and tsmoothing
to noSmoother
.
Constructors from data:
FourierAnalysis.TFanalyticsignal
— Function(1)
function TFanalyticsignal(x :: Vector{T},
sr :: Int,
wl :: Int = 0,
bandwidth :: IntOrReal = 2;
fmin :: IntOrReal = bandwidth,
fmax :: IntOrReal = sr÷2,
filtkind :: FilterDesign = Butterworth(2),
nonlinear :: Bool = false,
fsmoothing :: Smoother = noSmoother,
tsmoothing :: Smoother = noSmoother,
planner :: Planner = getplanner,
⏩ :: Bool = true) where T<:Real
(2)
function TFanalyticsignal(𝐱 :: Vector{Vector{T}},
<same arguments as method (1)>
(1)
Given sampling rate sr
, construct a TFAnalyticSignal object from univariate data x
, that is, a time-frequency representation of real vector x
.
Call three functions in series performing the three following operations:
i. pass the data throught a bank of pass-band filters calling function filterbank
,
ii. compute the analytic signal calling function analyticsignal
(method (1) therein),
iii. If requested, smooth the analytic signal along the time and/or frequency dimension calling function smooth
.
The arguments passed to each of these functions are:
filterbank | analyticsignal | smooth |
---|---|---|
x | wl | fsmoothing |
sr | nonlinear | tsmoothing |
bandwidth | planner | |
filtkind | ⏩ | |
fmin | ||
fmax | ||
⏩ |
Refer to the documentation of each function to learn the meaning of each argument.
By default the wl
argument is set to length(x)
.
If ⏩
is true (default), filtering and computation fo the analytic signal are multi-threaded across band-pass regions as long as the number of the regions is at least twice the number of threads Julia is instructed to use. See Threads.
(2)
Given sampling rate sr
, construct a TFAnalyticSignalVector object from a vector of univariate data 𝐱
, that is, from the time-frequency representations of the vectors in 𝐱
.
This method operates as method (1) with the following exceptions:
By default the wl
argument is set to length(x)
for each vectors in 𝐱
. If another values is given, it will be used for all of them.
If ⏩
is true (default), the method is run in multi-threaded mode across the vectors in 𝐱
as long as their number is at least twice the number of threads Julia is instructed to use, otherwise this method attempts to run each analytic signal estimation in multi-threaded mode like in method (1). See Threads.
If a Planner
is not explicitly passed as an argument, the FFTW plan is computed once and applied for all analytic signal estimations.
See: filterbank
, analyticsignal
, smooth
.
Examples:
using Plots, FourierAnalysis
# generate some data
sr, t, bandwidth=128, 512, 2
h=taper(harris4, t)
x1=sinusoidal(10, 8, sr, t, 0)
x2=sinusoidal(10, 19, sr, t, 0)
x=Vector((x1+x2).*h.y+randn(t))
y1=sinusoidal(10, 6, sr, t, 0)
y2=sinusoidal(10, 16, sr, t, 0)
y=Vector((y1+y2).*h.y+randn(t))
plot([x, y])
# TFAnalyticSignal object for x (method (1))
Y = TFanalyticsignal(x, sr, sr*4; fmax=32)
# vector of TFAnalyticSignal objects for x and y (method (2))
𝒀 = TFanalyticsignal([x, y], sr, sr*4; fmax=32)
# (for shortening comments: TF=time-frequency, AS=analytic signal)
# mean AS in a TF region (8:12Hz and samples 1:128) for the first object
m=mean(𝒀[1], (8, 12), (1, 128)) # Output a complex number
# extract the AS in a TF region (8:12Hz and samples 1:128) for the first object
E=extract(𝒀[1], (8, 12), (1, 128)) # Output a complex matrix
# mean AS in a TF region (8:12Hz and samples 1:128) for the two objects
𝐦=mean(𝒀, (8, 12), (1, 128)) # Output a vector of complex numbers
# same computation without checking homogeneity of the two objects in 𝒀 (faster)
𝐦=mean(𝒀, (8, 12), (1, 128); check=false)
# extract the AS in a TF region (8:12Hz and samples 1:128) for the two objects
𝐄=extract(𝒀, (8, 12), (8, 12)) # Output a vector of complex matrices
# plot the real part of the AS of x (see unit recipes.jl)
# gather first useful attributes for the plot
using Plots.Measures
tfArgs=(right_margin = 2mm,
top_margin = 2mm,
xtickfont = font(10, "Times"),
ytickfont = font(10, "Times"))
heatmap(tfAxes(Y)..., real(Y.y);
c=:pu_or, tfArgs...)
# ...the imaginary part
heatmap(tfAxes(Y)..., imag(Y.y);
c=:bluesreds, tfArgs...)
# ...the amplitude
heatmap(tfAxes(Y)..., amplitude(Y.y);
c=:amp, tfArgs...)
# ...the amplitude of the AS smoothed in frequency and time
heatmap(tfAxes(Y)..., amplitude(smooth(hannSmoother, hannSmoother, Y).y);
c=:fire, tfArgs...)
# ...the phase
heatmap(tfAxes(Y)..., phase(Y.y);
c=:bluesreds, tfArgs...)
# or generate a TFAmplitude object from the AS
A=TFamplitude(Y)
# and plot it (with different colors)
heatmap(tfAxes(A)..., A.y;
c=:fire, tfArgs...)
# generate a TFPhase object
ϴ=TFphase(Y)
# and plot it (with custom colors)
heatmap(tfAxes(ϴ)..., ϴ.y;
c=:pu_or, tfArgs...)
# compute and plot phase in [0, 2π]
heatmap(tfAxes(Y)..., TFphase(Y; func=x->x+π).y;
c=:amp, tfArgs...)
# compute and plot unwrapped phase
heatmap(tfAxes(Y)..., TFphase(Y; unwrapped=true).y;
c=:amp, tfArgs...)
# smooth time-frequency AS: smooth frequency
Z=smooth(blackmanSmoother, noSmoother, Y)
# plot amplitude of smoothed analytic signal
heatmap(tfAxes(Z)..., amplitude(Z.y);
c=:amp, tfArgs...)
# not equivalently (!), create an amplitude object and smooth it:
# in this case the amplitude is smoothed, not the AS
A=smooth(blackmanSmoother, noSmoother, TFamplitude(Y))
heatmap(tfAxes(A)..., A.y;
c=:fire, tfArgs...)
# Smoothing raw phase estimates is unappropriate
# since the phase is a discontinous function, however it makes sense
# to smooth phase if the phase is unwrapped.
heatmap(tfAxes(Y)...,
smooth(blackmanSmoother, noSmoother, TFphase(Y; unwrapped=true)).y;
c=:amp, tfArgs...)
# smooth AS: smooth both frequency and time
E=smooth(blackmanSmoother, blackmanSmoother, Y)
# plot amplitude of smoothed analytic signal
heatmap(tfAxes(E)..., amplitude(E.y);
c=:fire, tfArgs...)
# plot phase of smoothed analytic signal
heatmap(tfAxes(E)..., phase(E.y);
c=:bluesreds, tfArgs...)
# not equivalently (!), create amplitude and phase objects and smooth them
A=smooth(blackmanSmoother, blackmanSmoother, TFamplitude(Y))
heatmap(tfAxes(A)..., A.y;
c=:fire, tfArgs...)
ϴ=smooth(blackmanSmoother, blackmanSmoother, TFphase(Y, unwrapped=true))
heatmap(tfAxes(ϴ)..., ϴ.y;
c=:pu_or, tfArgs...)
# smooth again
ϴ=smooth(blackmanSmoother, blackmanSmoother, ϴ)
heatmap(tfAxes(ϴ)..., ϴ.y;
c=:pu_or, tfArgs...)
# and again ...
# create directly smoothed AS
Y=TFanalyticsignal(x, sr, t, bandwidth;
fmax=32,
fsmoothing=hannSmoother,
tsmoothing=hannSmoother)
# plot amplitude of smoothed analytic signal
heatmap(tfAxes(Y)..., amplitude(Y.y);
c=:amp, tfArgs...)
# create directly smoothed Amplitude
A=TFamplitude(x, sr, t, bandwidth;
fmax=32,
fsmoothing=hannSmoother,
tsmoothing=hannSmoother)
# plot smoothed amplitude
heatmap(tfAxes(A)..., A.y;
c=:amp, tfArgs...)
# compute a TFAnalyticSignal object with non-linear AS
Y=TFanalyticsignal(x, sr, t, bandwidth; fmax=32, nonlinear=true)
# check that it is non-linear
Y.nonlinear
# check that the amplitude is now 1.0 everywhere
norm(amplitude(Y.y)-ones(eltype(Y.y), size(Y.y))) # must be zero
# plot non-linear phase
heatmap(tfAxes(Y)..., phase(Y.y);
c=:bkr, tfArgs...)
# get the center frequencies of TFAmplitude object A
A.flabels
# extract the amplitude in a time-frequency region
extract(A, (2, 10), (1, 256)) # output a real matrix
# extract the amplitude in a time-frequency region at only one frequency
extract(A, 10, (1, 256)) # output a row vector
# extract the amplitude at one temporal sample at one frequency
extract(A, 10, 12) # or extract(A, 10.0, 12)
# extract amplitude at one temporal sample in a frequency region
extract(A, (10, 12), 12) # or extract(A, (10.0, 12.0), 12)
# extract amplitude at one temporal sample and all frequencies
extract(A, :, 12) # output a (column) vector
# compute the mean in a time-frequency region:
mean(A, (2, 10), (1, 256)) # output a real number
# is equivalent to (but may be less efficient than)
mean(extract(A, (2, 10), (1, 256)))
# using column sign for extracting all time samples
extract(A, (2, 10), :)
# This:
extract(A, :, :)
# is equivalent to this:
A.y
# but if you don't need to extract all frequencies,
# use the extract function to control what frequencies will be extracted:
# This
extract(A, (4, 5), 10)
# is not equivalent to this
A.y[4:5, 10]
# since the `extract` function finds the correct rows corresponding
# to the sought frequencies (in Hz), while A.y[4:5, 10]
# just returns the elements [4:5, 10] in the TF amplitude object
# Although the first center frequency in A is 2Hz, its
# band-pass region is 1-3Hz, therefore the frequency range 1:10 is accepted
mean(A, (1, 10), (1, 256))
# but this result in an error (see the REPL) since 0.5 is out of range:
mean(A, (0.5, 10), (1, 256))
# using a colon sign for time range
a=mean(A, (1, 10), :)
# using an integer for time range (indicates one specific sample)
a=mean(A, (1, 10), 16)
# using a colon sign for frequency range
a=mean(A, :, (1, 16))
# using a real number for frequency range
a=mean(A, 8.5, (1, 16))
# NB: the `extract` and `mean` functions work with the same syntax
# for objects TFAnayticSignal, TFAmplitude and TFPhase.
FourierAnalysis.TFamplitude
— Function(1)
function TFamplitude(Z::TFAnalyticSignal;
func::Function=identity)
(2)
function TFamplitude(𝐙::TFAnalyticSignalVector;
func::Function=identity)
(3)
function TFamplitude(x :: Vector{T},
sr :: Int,
wl :: Int,
bandwidth :: IntOrReal = 2;
func :: Function = identity,
filtkind :: FilterDesign = Butterworth(2),
fmin :: IntOrReal = bandwidth,
fmax :: IntOrReal = sr÷2,
fsmoothing :: Smoother = noSmoother,
tsmoothing :: Smoother = noSmoother,
planner :: Planner = getplanner,
⏩ :: Bool = true) where T<:Real =
(4)
function TFamplitude(𝐱 :: Vector{Vector{T}},
< same arguments as method (3) >
(1)
Construct a TFAmplitude object computing the amplitude of TFAnalyticSignal object Z
. Optional keyword argument func
is a function to be applied element-wise to the data matrix of the output. By default, the identity
(do nothing) function is applied.
(2)
Construct a TFAmplitudeVector object from a TFAnalyticSignalVector object executing method (1) for all TFAnalyticSignal objects in 𝐙
(3)
Call TFanalyticsignal
to obtain the time-frequency analytic signal of real signal vector x
and construct a TFAmplitude object holding the time-frequency amplitude (the mudulus, often referred to as the envelope) of x
.
All arguments are used for regulating the estimation of the analytic signal, with the exception of func
, fsmoothing
and fsmoothing
:
func
is an optional function to be applied to the amplitude data matrix output.
In order to estimate the analytic signal in the time-frequency domain this function calls the TFanalyticsignal
constructor (method (1) therein), with both fsmoothing
and tsmoothing
arguments set to noSmoother
. Arguments fsmoothing
and fsmoothing
are then used to smooth the amplitude.
In order to obtain amplitude estimations on smoothed analytic signal instead, create a TFAnalyticSignal object passing a Smoother to the TFanalyticsignal
constructor and then use method (1) to obtain the amplitude. Such amplitude estimation can be further smoothed using the smooth
function, as shown in the examples.
For the meaning of all other arguments, which are passed to function TFanalyticsignal
, see the documentation therein.
(4)
Construct a TFAmplitudeVector object from a vector of real signal vectors 𝐱
, executing method (3) for all of them. In order to estimate the time-frequency analytic signal for a vector of signals, method (2) of TFanalyticsignal
is called.
See: TFanalyticsignal
, TFAmplitude.
See also: amplitude
.
Examples: see the examples of TFanalyticsignal
.
FourierAnalysis.TFphase
— Function(1)
function TFphase(Z :: TFAnalyticSignal;
func :: Function = identity,
unwrapped :: Bool = false)
(2)
function TFphase(𝐙 :: TFAnalyticSignalVector;
func ::Function = identity,
unwrapped ::Bool = false)
(3)
function TFphase(x :: Vector{T},
sr :: Int,
wl :: Int,
bandwidth :: IntOrReal = 2;
unwrapped :: Bool = false,
func :: Function = identity,
filtkind :: FilterDesign = Butterworth(2),
fmin :: IntOrReal = bandwidth,
fmax :: IntOrReal = sr÷2,
nonlinear :: Bool = false,
fsmoothing :: Smoother = noSmoother,
tsmoothing :: Smoother = noSmoother,
planner :: Planner = getplanner,
⏩ :: Bool = true) where T<:Real
(4)
function TFphase(𝐱 :: Vector{Vector{T}},
< same arguments as method (3) >
(1)
Construct a TFPhase object computing the phase of TFAnalyticSignal object Z
. By default the phase is represented in $[−π, π]$.
If optional keyword argument unwrapped
is true (false by defalut), the phase is unwrapped, that is, it holds the cumulative sum of the phase along the time dimension once this is represented in $[0, 2π]$.
Optional keyword argument func
is a function to be applied element-wise to the data matrix of the output. By default, the identity
(do nothing) function is applied. If unwrapped
is true, the function is applied on the unwrapped phase.
(2)
Construct a TFPhaseVector object from a TFAnalyticSignalVector object executing method (1) for all TFAnalyticSignal objects in 𝐙
(3)
Call TFanalyticsignal
to obtain the time-frequency analytic signal of real signal vector x
and construct a TFPhase object holding the time-frequency phase (argument) of x
.
All arguments are used for regulating the estimation of the analytic signal, with the exception of unwrapped
, func
, fsmoothing
and fsmoothing
.
unwrapped
has the same meaning as in method (1) and (2).
func
is an optional function to be applied to the phase data matrix output. If unwrapped
is true, it is applied to the unwrapped phase.
In order to estimate the analytic signal in the time-frequency domain this function calls the TFanalyticsignal
constructor (method (1) therein), with both fsmoothing
and tsmoothing
arguments set to noSmoother
. fsmoothing
and fsmoothing
are then used to smooth the phase if unwrapped
is true.
In order to obtain phase estimations on smoothed analytic signal instead, create a TFAnalyticSignal object passing a Smoother to the TFanalyticsignal
constructor and then use method (1) to obtain the phase.
For the meaning of all other arguments, which are passed to function TFanalyticsignal
, see the documentation therein.
(4)
Construct a TFPhaseVector object from a vector of real signal vectors 𝐱
, executing method (3) for all of them. In order to estimate the time-frequency analytic signal for a vector of signals, method (2) of TFanalyticsignal
is called.
See: TFanalyticsignal
, TFPhase.
See also: phase
, unwrapPhase
, polar
.
Examples: see the examples of TFanalyticsignal
.