1# File src/library/stats/R/fft.R 2# Part of the R package, https://www.R-project.org 3# 4# Copyright (C) 1995-2012 The R Core Team 5# 6# This program is free software; you can redistribute it and/or modify 7# it under the terms of the GNU General Public License as published by 8# the Free Software Foundation; either version 2 of the License, or 9# (at your option) any later version. 10# 11# This program is distributed in the hope that it will be useful, 12# but WITHOUT ANY WARRANTY; without even the implied warranty of 13# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 14# GNU General Public License for more details. 15# 16# A copy of the GNU General Public License is available at 17# https://www.R-project.org/Licenses/ 18 19fft <- function(z, inverse=FALSE) .Call(C_fft, z, inverse) 20 21mvfft <- function(z, inverse=FALSE) .Call(C_mvfft, z, inverse) 22 23nextn <- function(n, factors=c(2,3,5)) .Call(C_nextn, n, factors) 24 25convolve <- function(x, y, conj=TRUE, type=c("circular","open","filter")) 26{ 27 type <- match.arg(type) 28 n <- length(x) 29 ny <- length(y) 30 Real <- is.numeric(x) && is.numeric(y) 31 ## switch(type, circular = ..., ) 32 if(type == "circular") { 33 if(ny != n) 34 stop("length mismatch in convolution") 35 } 36 else { ## "open" or "filter": Pad with zeros 37 n1 <- ny - 1 38 x <- c(rep.int(0, n1), x) 39 n <- length(y <- c(y, rep.int(0, n - 1)))# n = nx+ny-1 40 } 41 x <- fft(fft(x)* (if(conj)Conj(fft(y)) else fft(y)), inverse=TRUE) 42 if(type == "filter") 43 (if(Real) Re(x) else x)[-c(1L:n1, (n-n1+1L):n)]/n 44 else 45 (if(Real) Re(x) else x)/n 46} 47 48