ndft {cooltools} | R Documentation |
Non-uniform Discrete Fourier Transform
Description
Compute the one-dimensional Non-uniform Discrete Fourier Transform (NDFT). This is needed if the data are sampled at irregularly spaced times.
Usage
ndft(
f,
x = seq(0, length(f) - 1)/length(f),
nu = seq(0, length(f) - 1),
inverse = FALSE,
weighing = TRUE,
simplify = TRUE
)
Arguments
f |
vector of real or complex function values |
x |
vector of points in direct space, typically time or position coordinates. If |
nu |
vector of frequencies in units of [1/units of x]. If |
inverse |
logical flag; if TRUE, the inverse Fourier transform is performed. |
weighing |
logical flag; if TRUE, irregularly spaced evaluations of |
simplify |
logical flag; if TRUE, the complex output array will be simplified to a real array, if it is real within the floating point accuracy. |
Details
The one-dimensional NDFT of a vector f=(f_1,...,f_N)
is defined as
F_j=\sum_i w_i f_i exp(-2\pi i*x_i*\nu_j)
where w_i
are optional weights, proportional to the interval around x_i
, only used if weighing=TRUE
. Likewise, the inverse NDFT is defined as
f_i=\sum_j w_j F_j exp(+2\pi i*x_i*nu_j)
where w_j
are optional weights, proportional to the interval around \nu_j
. In this implementation NDFTs are computed using a brute force algorithm, scaling as O(N*N)
, which is considerably worse than the O(N)*log(N)
scaling of FFT algorithms. It is therefore important to pick the required frequencies wisely to minimise computing times.
Value
Returns a vector of the same length as x
(if inverse=FALSE
) or nu
(if inverse=TRUE
).
Author(s)
Danail Obreschkow
See Also
Examples
# Define an example signal
nu1 = 1 # [Hz] first frequency
nu2 = 8 # [Hz] second frequency in the signal
s = function(t) sin(2*pi*nu1*t)+0.7*cos(2*pi*nu2*t+5)
# Discretize signal
N = 50 # number of samples
t.uniform = seq(0,N-1)/N
t.nonuniform = t.uniform^1.3
s.uniform = s(t.uniform)
s.nonuniform = s(t.nonuniform)
# Plot signal
oldpar = par(mfrow = c(1, 2))
curve(s,0,1,500,xaxs='i',main='Time signal',xlab='Time t',ylab='s(t)',col='grey')
points(t.uniform,s.uniform,pch=16,cex=0.8)
points(t.nonuniform,s.nonuniform,pch=4,col='blue')
legend('topright',c('Continuous signal','Uniform sample','Non-uniform sample'),
lwd=c(1,NA,NA),pch=c(NA,16,4),col=c('grey','black','blue'),pt.cex=c(1,0.8,1))
# Uniform and non-uniform DFT
nu = seq(0,N-1) # discrete frequencies
spectrum.uniform = stats::fft(s.uniform)
spectrum.nonuniform = ndft(s.nonuniform,t.nonuniform,nu)
spectrum.wrong = stats::fft(s.nonuniform)
# Evaluate power
power.uniform = Mod(spectrum.uniform)^2
power.nonuniform = Mod(spectrum.nonuniform)^2
power.wrong = Mod(spectrum.wrong)^2
# Plot DFT and NDFT up to Nyquist frequency
plot(nu,power.uniform,pch=16,cex=0.8,xlim=c(0,N/2),xaxs='i',
main='Power spectrum',xlab=expression('Frequency'~nu~'[Hz]'),ylab='Power')
points(nu,power.nonuniform,pch=4,col='blue')
points(nu,power.wrong,pch=1,col='red')
abline(v=c(nu1,nu2),col='grey',lty=2)
legend('topright',c('DFT of uniform sample','NDFT of non-uniform sample',
'DFT of non-uniform sample (wrong)','Input frequencies'),
lwd=c(NA,NA,NA,1),lty=c(NA,NA,NA,2),pch=c(16,4,1,NA),
col=c('black','blue','red','grey'),pt.cex=c(0.8,1,1,NA))
par(oldpar)