synsq_cwt_fw | R Documentation |
This function calculates the synchrosqueezing transform.
This code is translated from MATLAB Synchrosqueezing Toolbox, version 1.1 developed by Eugene Brevdo (http://www.math.princeton.edu/~ebrevdo/).
synsq_cwt_fw(tt, x, nv=16, opt=NULL)
tt |
vector of times samples are taken (e.g. |
x |
vector of signal samples (e.g. |
nv |
number of voices (default: 16, recommended: 32 or 64) |
opt |
list of options. |
This function
calculates the synchrosqueezing transform of vector x
, with
samples taken at times given in vector tt
. Uses nv
voices. opt
is a struct of options for the way synchrosqueezing is done, the
type of wavelet used, and the wavelet parameters. This
implements the algorithm described in Sec. III of [1].
Note that Wx
itself is not thresholded by opt$gamma
.
The instantaneoues frequency w
is calculated using Wx
by cwt_freq_direct
and
cwt_freq
. After the calculation, instantaneoues frequency w
is treated as NA
where abs(Wx) < opt$gamma
.
Tx |
synchrosqueezed output of |
fs |
frequencies associated with rows of |
Wx |
wavelet transform of |
asc |
scales associated with rows of |
w |
instantaneous frequency estimates for each element of |
[1] Thakur, G., Brevdo, E., Fuckar, N. S. and Wu, H-T. (2013) The Synchrosqueezing algorithm for time-varying spectral analysis: Robustness properties and new paleoclimate applications. Signal Processing, 93, 1079–1094.
cwt_fw
, est_riskshrink_thresh
, wfiltfn
.
tt <- seq(0, 10, , 1024) nv <- 32 f0 <- (1+0.6*cos(2*tt))*cos(4*pi*tt+1.2*tt^2) sigma <- 0.5 f <- f0 + sigma*rnorm(length(tt)) # Continuous wavelet transform opt <- list(type = "bump") cwtfit <- cwt_fw(f, opt$type, nv, tt[2]-tt[1], opt) # Hard thresholing thresh <- est_riskshrink_thresh(cwtfit$Wx, nv) cwtfit$Wx[which(abs(cwtfit$Wx) < thresh)] <- 0.0 # Reconstruction opt$gamma <- thresh #[1] 0.0593984 #opt$gamma <- 10^-5 cwtrec <- cwt_iw(cwtfit$Wx, opt$type, opt) par(mfrow=c(1,1)) plot(tt, f, type="p", lty=2, xlab="time", ylab="f", col="red", cex=0.1) lines(tt, f0, col="blue") lines(tt, cwtrec) # Synchrosqueezed wavelet transform sstfit <- synsq_cwt_fw(tt, f, nv, opt) #par(mfrow=c(2,2)) #plot(tt, f, type="p", lty=2, xlab="time", ylab="f", col="red", cex=0.1) #lines(tt, f0, col="blue") #image.plot(list(x=tt, y=sstfit$asc, z=t(abs(sstfit$Wx))), log="y", # xlab="Time", ylab="Scale", main="Time-Scale Representation by CWT", # col=designer.colors(64, c("azure", "cyan", "blue", "darkblue")), ylim=rev(range(sstfit$asc))) #image.plot(list(x=tt, y=sstfit$fs, z=t(abs(sstfit$Tx))), log="y", # xlab="Time", ylab="Frequency", main="Time-Frequency Representation by SST", # col=designer.colors(64, c("azure", "cyan", "blue", "darkblue")), ylim=c(0.5, 25)) #image.plot(list(x=tt, y=sstfit$asc, z=t(sstfit$w)), log="y", # xlab="Time", ylab="Scale", main="Instantaneous Frequency", # col=designer.colors(64, c("azure", "cyan", "blue", "darkblue")), ylim=rev(range(sstfit$asc)))
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.