Computes Short Time Fourier Transforms
Processes a dataset, creating an object contained processed time-frequency analyses. These can then be plotted.
stft(X, start=0, end=1, length=NULL, time.format = c("auto"),
type = c("mv", "svm", "sum"), mv.indices,
date.col, reassign = TRUE, plot.it = FALSE,...)X |
The dataset to be processed. |
start, end, length, time.format |
A specification for the segment to process, as in get.intervals. |
type |
The type of STFT to compute. |
mv.indices |
For type = "mv" or type = "sum", the indices to process and the order to process them in. |
date.col |
logical. Whether the first column should be ignored and treated as a timestamp. If unset, is automatically chosen. |
reassign |
logical. If TRUE, compute the time-reassigned STFT. For type c("mv", "sum"), this is done with the first coordinate in mv.indices. |
plot.it |
logical. Whether to plot the STFT immediately when processing is complete, using the default plot.stft options. |
... |
Additional optional arguments to control the STFT computation. These are:
|
This function accepts input in a variety of forms and computes short time fourier transforms to extract frequency structure from the data.X may be an array, a vector, or an AccData object. If date.col is TRUE, the first column of an array X would be used to determine timestamps. Otherwise indices would be used. If date.col is not set, the function will attempt to determine whether the first column is timestamp-like. The timestamp column is removed from X (and so not included in consideration of mv.indices, for instance). With vectors, the basic method is to compute the STFT by creating windows of size win seconds every inc seconds, and computing the fourier transform. With multi-dimensional data and AccData, processing is done on the dimensions that are in mv.indices, or the first three non-date columns if that is unavailable. Three methods are possible:
1. type = "mv": The one dimensional method is first applied to each of the chosen column indices. These are then collated by choosing, for each time-frequency combination, the maximum such value across each of the indices.
2. type = "svm": The SVM is computed first for each time step by computing the square rooted sum of squares. This is then dealt with using the one dimensional method.
3. type = "sum": As in "mv", the 1d method is applied. The square of the modulus of the result is then summed and square rooted.
If reassign is set, the time reassigned stft is also computed for the first element of mv.indices or the svm as appropriate, by using finite differencing. This gives potentially better resolution results for data with a clear signal component.
A "stft" class object - a list with the following components:
call: The function call.
type: Type of STFT computed.
values: Mod of FFT computed, with each row corresponding to a specific time increment.
increment,windowsize,centre,sampling.frequency: Various control parameters used in the computation.
null.logmean,null.logsd: Log of the square rooted mean and standard deviation of the Mod FFT squared for the randomised data, if calc.null = TRUE.
p.values: Wilcoxian pvalues, if pvalues = TRUE.
principals: Principal frequencies.
frequency: Frequencies at which FFT is computed.
time: Timestamps for FFT windows.
LGD: Local group delay matrix for reassigned STFT.
CIF: Channelized instantaneous frequency matrix for reassigned STFT.
Fulop, S.A. & Fitz, K. (2006). Algorithms for computing the time-corrected instantaneous frequency (reassigned) spectrogram, with applications J Acoustical Society of America 119(1), 360–371. Nelson. D.J. (2001). Cross-spectral methods for processing speech J Acoustical Society of America 110(1), 2575-2592.
## Not run:
#Some artificial data
time = 1:5000
#sum of two sine curves at 0.3 Hz and 0.05 Hz
f1 = 0.3; f2 = 0.05
sin1 = sin(time * f1 * 2*pi)
sin2 = sin(time * f2 * 2*pi)
#add a bit of noise
signal = sin1 + sin2 + 1*rnorm(5000)
#non-reassigned
stft(signal, plot = TRUE, reassign = FALSE, win = 100)
#reassigned
stft(signal, plot = TRUE, reassign = TRUE, win = 100)
#add a third component: varying frequency.
stft(signal + sin(cumsum(seq(f2, f1, length = 5000))*2*pi),
plot = TRUE, reassign = TRUE, win = 100)
# Real data
binfile = system.file("binfile/TESTfile.bin", package = "GENEAread")[1]
# Read in the entire file, calibrated
procfile<-read.bin(binfile)
# Default is mv
stft(procfile, plot.it = TRUE)
# Try sum?
stft(procfile, plot.it = TRUE, type = "sum", reassign = FALSE)
# Just look at the last 50% of the data
stft(procfile, start = 0.5, plot.it = TRUE)
# not reassigned, svm
stft(procfile, type = "svm", reassign = FALSE, plot.it = TRUE)
# a narrower 5 second window means better time resolution
stft(procfile, type = "svm", reassign = FALSE, plot.it = TRUE, win = 5)
# choose increments so as not to overlap
stft(procfile, type = "svm", reassign = FALSE, plot.it = TRUE, win = 5, inc = 5)
# uniform windows
stft(procfile, type = "svm", reassign = FALSE, plot.it = TRUE, wtype = "uniform.window")
# Svm, reassigned, quietly
obj = stft(procfile, type = "svm", quiet = TRUE)
plot(obj, cex = 3, showmax = FALSE, mode = "pval")
#example code
plot(stft(subs(mag, 0.94,0.96), win = 1024, plot = F, coef = 512), zlog = T, log="y")
plot(stft(subs(mag, 0.7,8), win = 1024, plot = F, coef = 512), zlog = T, log="y")
plot(stft(subs(mag, 0.0001,0.005), win = 1024, plot = F, coef = 512), zlog = T)
plot(stft(subs(mag, 0.7,0.8), win = 1024, plot = F), zlog = T, log = "y")
plot(stft(rep(1, 1000) +
c(sin(1:500/ 10 * 2*pi), rep(0, 500)) +
c(rep(0, 300),sin(1:500/ 20 * 2*pi), rep(0, 200)),
freq = 1, plot.it = F), log="x")
stft(sin(1:1000 / (1 +sqrt(1000:1)) * 2 * pi), freq = 1)
stft(rep(1, 1000) + sin(1:1000/ 10 * 2*pi), freq = 1)
## End(Not run)Please choose more modern alternatives, such as Google Chrome or Mozilla Firefox.