TFtools: times series analysis

Created on Mon May 03 14:53:54 2010

@author: a1185872

MTpy.core.TFtools.dctrend(f)

dctrend(f) will remove a dc trend from the function f.

Inputs:

f = array to dctrend

Outputs:

fdc = array f with dc component removed

MTpy.core.TFtools.decimatef(f, m)

Will decimate a function by the factor m. First an 8th order Cheybechev type I filter with a cuttoff frequency of .8/m is applied in both directions to minimize any phase distortion and remove any aliasing. Note decimation values above 10 will typically result in bad coefficients, therefore if you decimation is more than 10 just repeat the decimation until the desired decimation is reached.

Inputs:

f = array to be decimated m = decimation factor

Outputs:

fdec = array f decimated by factor m

MTpy.core.TFtools.dwindow(window)

Calculates the derivative of the given window

Input:

window = some sort of window function

Output:

dwin = derivative of window

MTpy.core.TFtools.gausswin(winlen, alpha=2.5)

gausswin will compute a gaussian window of length winlen with a variance of alpha

Inputs:

winlen = length of desired window alpha = 1/standard deviation of window, ie full width half max of window

Outputs:

gwin = gaussian window

MTpy.core.TFtools.modifiedb(fx, tstep=2**5, nfbins=2**10, df=1.0, nh=2**8-1, beta=.2)

will calculate the modified b distribution as defined by cosh(n)^-2 beta for a function fx.

Inputs:

fx = array from which modifiedb will be calculated if computing cross

spectra input as [fx1,fx2]

tstep = time step between windows nfbins = number of frequencies df = sampling frequency (Hz) nh = length of time-domain smoothing window (needs to be odd) beta = smoothing coefficient

Outputs:

tfarray = modifiedB estimation of array fx tlst = time instances of each calculation flst = array of positive frequencies

MTpy.core.TFtools.normalizeL2(f)

normalizeL2(f) will return the function f normalized by the L2 norm -> f/(sqrt(sum(abs(x_i)^2))).

Inputs:

f = array to be normalized

Outputs:

fnorm = array f normalized in L2 sense

MTpy.core.TFtools.padzeros(f, npad=None, padpattern=None)

padzeros(f) will return a function that is padded with zeros to the next power of 2 for faster processing for fft or to length npad if given.

Inputs:

f = array to pad npad = length to pad to defaults to next power of two padpattern = pattern to pad with default is zero

Outputs:

fpad = array f padded to length npad with padpattern

MTpy.core.TFtools.plotAll(fx, tfarray, tlst, flst, fignum=1, starttime=0, timeinc='hrs', dt=1.0, title=None, vmm=None, cmap=None, aspect=None, interpolation=None, cbori=None, cbshrink=None, cbaspect=None, cbpad=None, normalize='n', scale='log')

plottf(tfarray,tlst,flst,fignum=1) will plot a calculated tfarray with limits corresponding to tlst and flst. Can have:

Inputs:

starttime = starttime measured in timeincrement timeincrement = ‘hrs’,’min’ or ‘sec’ vmm = [vmin,vmax] a list for min and max title = title string cmap = colormap scheme default is jet, type help on matplotlib.cm aspect = aspect of plot, default is auto, can be ‘equal’ or a scalar interpolation = type of color interpolation, type help on

matplotlib.pyplot.imshow

cbori = colorbar orientation ‘horizontal’ or ‘vertical’ cbshrink = percentage of 1 for shrinking colorbar cbaspect = aspect ratio of long to short dimensions cbpad = pad between colorbar and axis normalization = y or n, y for normalization n for none

Outputs:

plot

MTpy.core.TFtools.plottf(tfarray, tlst, flst, fignum=1, starttime=0, timeinc='hrs', dt=1.0, title=None, vmm=None, cmap=None, aspect=None, interpolation=None, cbori=None, cbshrink=None, cbaspect=None, cbpad=None, powscale='log', normalize='n', yscale='log', period='n')

plottf(tfarray,tlst,flst,fignum=1) will plot a calculated tfarray with limits corresponding to tlst and flst.

Inputs:

starttime = starttime measured in timeincrement tinc = ‘hrs’,’min’ or ‘sec’ vmm = [vmin,vmax] a list for min and max title = title string cmap = colormap scheme default is jet, type help on matplotlib.cm aspect = aspect of plot, default is auto, can be ‘equal’ or a scalar interpolation = type of color interpolation, type help on

matplotlib.pyplot.imshow

cbori = colorbar orientation ‘horizontal’ or ‘vertical’ cbshrink = percentage of 1 for shrinking colorbar cbaspect = aspect ratio of long to short dimensions cbpad = pad between colorbar and axis powscale = linear or log for power normalize = y or n, yes for normalization, n for no yscale = linear or log plot yscale period = ‘y’ or ‘n’ to plot in period instead of frequency

Outputs:

plot

MTpy.core.TFtools.reassignedSmethod(fx, nh=127, tstep=16, nfbins=512, df=1.0, alpha=4, thresh=0.01, L=5)

reassignedSmethod(fx,nh=2**7-2,tstep=2**4,nfbins=2**9,df=1.0,alpha=4,

thresh=.05,L=5)

will calulate the reassigned S-method as described by Djurovic[1999] by using the spectrogram to estimate the reassignment

Inputs:

fx = 1-d array to be processed nh = window length for each time instance tstep = step between time instances nfbins = number of frequency bins, note output will be nfbins/2 due to

symmetry of the FFT

df = sampling rate (Hz) alpha = inverse of full-width half max of gaussian window, smaller

numbers mean broader windows

thresh = threshold for reassignment, lower numbers more points

reassigned, higer numbers less points reassigned

L = length of window for S-method calculation, higher numbers tend

tend toward WVD

Outputs:

rtfarray = reassigned S-method shape of (nfbins/2,len(fx)/tstep) tlst = list of time instances where rtfarray was calculated flst = positive frequencies sm = S-method array

MTpy.core.TFtools.reassignedstft(fx, nh=63, tstep=32, nfbins=1024, df=1.0, alpha=4, threshold=None)

reassignedstft(fx,nh=2**5-1,tstep=2**8,nfbins=2**10,df=1.0,alpha=20) will compute the reassigned spectrogram by estimating the center of gravity of the signal and condensing dispersed energy back to that location.

Inputs:

fx = time series to be analyzed nh = length of gaussian window, should be odd tstep = time step for each window calculation nfbins = number of frequency bins to calculate, note result will be

length nfbins/2

df = sampling frequency (Hz) alpha = reciprocal of full width half max of gaussian window threshold = threshold value for reassignment

Outputs:

rtfarray = reassigned spectrogram in units of amplitude tlst = array of time instances where windows were calculated for ploting flst = array of frequencies for plotting stft = standard spectrogram in units of amplitude

MTpy.core.TFtools.robustSmethod(fx, L=5, nh=128, tstep=32, nfbins=1024, df=1.0, robusttype='median', sigmal=None)

robustSmethod(fx,L=15,nh=2**7,tstep=2**5,nfbins=2**10,df=1.0) computes the robust Smethod via the robust spectrogram.

Inputs:

fx = array of data, if computing cross-spectra input as [fa,fb] L = frequency smoothing window if robusttype=’median’ nh = window length for STFT tstep = time step for each STFT to be computed nfbins = number of frequency bins to be calculate df = sampling frequency robusttype = type of robust STFT calculation can be ‘median’ or ‘L’ simgal = full-width half max of gaussian window applied in frequency

Outputs:

tfarray = robust S-method estimation of array fx tlst = time instances of each calculation flst = array of positive frequencies

MTpy.core.TFtools.robuststftL(fx, alpha=0.325, nh=256, tstep=32, df=1.0, nfbins=1024)

robuststftL(fx,nh=2**8,tstep=2**5,ng=1,df=1.0) will output an array of the time-frequency robust spectrogram by estimating the vector median and summing terms estimated by alpha coefficients.

Inputs:

fx = the function to have a spectrogram computed for alpha = robust parameter [0,.5] -> 0 gives spectrogram, .5 gives median stft nh = window length for each time step tstep = time step between short windows df = sampling frequency nfbins = number of frequency bins

Outputs:

tfarray = robust L-estimation of array fx tlst = time instances of each calculation flst = array of positive frequencies

MTpy.core.TFtools.robuststftMedian(fx, nh=256, tstep=32, df=1.0, nfbins=1024)

robuststftMedian(fx,nh=2**8,tstep=2**5,ng=1,df=1.0) will output an array of the time-frequency robust spectrogram calculated using the vector median simplification.

Inputs:

fx = the function to have a spectrogram computed for nh = window length for each time step tstep = time step between short windows df = sampling frequency nfbins = number of frequency bins

Outputs:

tfarray = WVD estimation of array fx tlst = time instances of each calculation flst = array of positive frequencies

MTpy.core.TFtools.robustwvd(fx, nh=127, ng=15, tstep=16, nfbins=256, df=1.0, sigmanh=None, sigmang=None)

robustwvd(fx,tstep=2**5,nfbins=2**10,df=1.0,nh=2**8-1,ng=2**5-1,

sigmanh=None,sigmang=None)

will calculate the smoothed pseudo Wigner-Ville distribution for a function fx. smoothed with Gaussians windows to get best localization.

Inputs:

fx = array to estimate spwvd, input as [fx1,fx2] if computing cross

spectra

tstep = time step between windows nfbins = number of frequencies df = sampling frequency (Hz) ng = length of time-domain smoothing window (needs to be odd) nh = length of frequency-domain smoothing window (needs to be odd) sigmanh = std of window h, ie full width half max of gaussian sigmang = std of window g, ie full width half max of gaussian

Outputs:

tfarray = SPWVD estimation of array fx tlst = time instances of each calculation flst = array of positive frequencies

MTpy.core.TFtools.sfilter(f, fcutoff=10.0, w=10.0, dt=0.001)

Will apply a sinc filter of width w to the function f by multipling in the frequency domain. Returns filtered function

Inputs:

f = array to filter fcuttoff = cutoff frequency w = length of filter dt = sampling time (s)

Outputs:

filtfunc = filtered function

MTpy.core.TFtools.smethod(fx, L=11, nh=256, tstep=128, ng=1, df=1.0, nfbins=1024, sigmaL=None)

smethod(fx,L=11,nh=2**8,tstep=2**7,ng=1,df=1.0,nfbins=2**10) will calculate the smethod by estimating the STFT first and computing the WV of window length L in the frequency domain. For larger L more of WV estimation, if L=0 get back STFT

Inputs:

fx = the function to have a S-methoc computed for, if computing cross

spectra input as [fx1,fx2]

L = window length in frequency domain nh = window length for each time step tstep = time step between short windows ng = smoothing window along frequency plane should be odd df = sampling frequency nfbins = number of frequency bins

Outputs:

tfarray = S-method estimation of array fx tlst = time instances of each calculation flst = array of positive frequencies

MTpy.core.TFtools.specwv(fx, tstep=32, nfbins=1024, nhs=256, nhwv=511, ngwv=7, df=1.0)

specwv(f,tstep=2**5,nfbins=2**10,nh=2**8-1,ng=1,df=1.0) will calculate the Wigner-Ville distribution mulitplied by the STFT windowed by the common gaussian window h for a function f.

Inputs:

fx = array to compute the specwv tstep = time step between windows nfbins = number of frequencies nhs = length of time-domain smoothing window for STFT should be even nhwv = length of time-domain smoothing window for WV (needs to be odd) ngwv = lenght of frequency-domain smoothing window (needs to be odd) df = sampling frequency (Hz)

Outputs:

tfarray = SPECWV estimation of array fx tlst = time instances of each calculation flst = array of positive frequencies

MTpy.core.TFtools.spwvd(fx, tstep=32, nfbins=1024, df=1.0, nh=None, ng=None, sigmat=None, sigmaf=None)

spwvd(fx,tstep=2**5,nfbins=2**10,df=1.0,nh=2**8-1,ng=2**5-1,sigmat=None,

sigmaf=None)

will calculate the smoothed pseudo Wigner-Ville distribution for a function fx. smoothed with Gaussians windows to get best localization.

Inputs:

fx = array to estimate spwvd, input as [fx1,fx2] if computing cross

spectra

tstep = time step between windows nfbins = number of frequencies df = sampling frequency (Hz) ng = length of time-domain smoothing window (needs to be odd) nh = length of frequency-domain smoothing window (needs to be odd) sigmat = std of window h, ie full width half max of gaussian sigmaf = std of window g, ie full width half max of gaussian

Outputs:

tfarray = SPWVD estimation of array fx tlst = time instances of each calculation flst = array of positive frequencies

MTpy.core.TFtools.stfbss(X, nsources=5, ng=31, nh=511, tstep=63, df=1.0, nfbins=1024, tftol=1e-08, normalize=True)

btfssX,nsources=5,ng=2**5-1,nh=2**9-1,tstep=2**6-1,df=1.0,nfbins=2**10,

tftol=1.E-8,normalize=True)

estimates sources using a blind source algorithm based on spatial time-frequency distributions. At the moment this algorithm uses the SPWVD to estimate TF distributions.

Inputs:

X = m x n array of time series, where m is number of time series and n

is length of each time series

nsources = number of estimated sources ng = frequency window length nh = time window length tstep = time step increment df = sampling frequency (Hz) nfbins = number of frequencies tftol = tolerance for a time-frequency point to be estimated as a cross

term or as an auto term, the higher the number the more auto terms.

normalization = True or False, True to normalize, False if already

normalized

Outputs:

Se = estimated individual signals up to a permutation and scale Ae = estimated mixing matrix as X=A*S

MTpy.core.TFtools.stft(fx, nh=256, tstep=128, ng=1, df=1.0, nfbins=1024)

stft(fx,nh=2**8,tstep=2**7,ng=1,df=1.0) will calculate the spectrogam of the given function by calculating the fft of a window of length nh at each time instance with an interval of tstep. The frequency resolution is nfbins Can compute the cross STFT by inputting fx as [fx1,fx2]

Inputs:

fx = the function to have a spectrogram computed for can be two functions

input as [fx1,fx2]

nh = window length for each time step tstep = time step between short windows ng = smoothing window along frequency plane should be odd df = sampling frequency nfbins = number of frequency bins

Outputs:

tfarray = spectrogram in units of amplitude tlst = time instance array where each window was calculated flst = frequency array containing only positive frequencies

MTpy.core.TFtools.wvd(fx, nh=255, tstep=32, nfbins=1024, df=1.0)

wvd(f,nh=2**8-1,tstep=2**5,nfbins=2**10,df=1.0) will calculate the Wigner-Ville distribution for a function f. Can compute the cross spectra by inputting fx as [fx1,fx2]

Inputs:

fx = array for which WVD will be calculated, input as [fx1,fx2] for

cross-spectra calculation

nh = window length, needs to be odd so centered on zero tstep = time step between windows nfbins = number of frequencies df = sampling frequency (Hz)

Outputs:

tfarray = WVD estimation of array fx tlst = time instances of each calculation flst = array of positive frequencies

MTpy.core.TFtools.wvdas(fx)

wvdas(fx) will compute the analytic signal for WVVD as defined by J. M. O’ Toole, M. Mesbah, and B. Boashash, (2008), “A New Discrete Analytic Signal for Reducing Aliasing in the Discrete Wigner-Ville Distribution”, IEEE Trans. on Signal Processing,

Inputs:

fx = signal to compute anlytic signal for with length N

Outputs:

fxa = analytic signal of fx with length 2*N