Свертка для цифровой обработки сигналов в R

У меня есть простая цифровая система, которая имеет вход x (n) = u (n) - u (n-4). Система 1

Я пытаюсь найти выход y (n) с помощью функции conv () из пакета «сигнал» или функции свертки () из пакета «статистика» и построить график зависимости y (n) от n для -10 ≤ n ≤ 10.

Пока у меня есть следующий код:

library(signal)

n <- c(-10:10)                           # Time index
x <- c(rep(0, 10), rep(1, 4), rep(0, 7)) # Input Signal
h1 <- c(rep(0, 11), 0.5, rep(0, 9))      # Filter 1
h2 <- 0.8^n                              # Filter 2
h2[0:11] <- 0                            #

system <- data.frame(n, x, h1, h2)


y <- conv(x + conv(x, h1), h2)           # Output Signal

system <- transform(system, y=y[1:21]) 

plot(system$n, system$y)

Я проверил этот сюжет, и он очень неправильный. Я думаю, что есть некоторая переработка векторов, когда я выполняю свертку, и выходные данные функции conv(), похоже, не совпадают с исходным временным индексом. Я просто не могу понять, как исправить мою логику здесь. Я понимаю, что функция conv(n, m) возвращает вектор длины (m+n)-1, есть ли хороший способ легко сопоставить этот вектор с вектором индекса времени?

Это потребует некоторых знаний в области цифровой обработки сигналов, а также кодирования в R, и было бы здорово, если бы кто-то имел опыт использования R для этой цели и мог бы дать несколько советов. Заранее спасибо.


r vector convolution signal-processing digital
person N8TRO    schedule 02.02.2013    source источник

Ответы (1)


arrow_upward
4
arrow_downward

Я понял это. Центр вывода функции conv() совпадает с центром вектора индекса времени. Как таковой:

library(signal)

n <- c(-10:10)                           # Time index
x <- c(rep(0, 10), rep(1, 4), rep(0, 7)) # Input Signal, square pulse
h1 <- c(rep(0, 11), 0.5, rep(0, 9))      # Filter 1
h2 <- 0.8^n                              # Filter 2
h2[1:10] <- 0                            #

system <- data.frame(n, x, h1, h2)

y <- conv(x + conv(x, h1)[11:31], h2)    # Output Signal

system <- transform(system, y=y[11:31]) 

plot(system$n, system$y)

Я буду работать над общей формой, чтобы добиться этого, так как я буду делать это регулярно и не хотел бы делать это каждый раз вручную. Если кто-то опередит меня в этом, пожалуйста, поделитесь. :)

ОБНОВИТЬ

Создал общую форму функции conv() для автоматического выравнивания индексов входных и выходных векторов. Это достигается за счет отсутствия полной свертки, поэтому вам придется сначала настроить ввод так, чтобы он отображал всю интересующую область.

library(signal) # Should this be inside the func. with attach(), detach()?

conv2 <- function(x, y){
    conv(x, y)[ceiling(length(x)/2):(length(x)+floor(length(x)/2))]
}

# so 
y <- conv2(x + conv2(x, h1), h2)

ОБНОВЛЕНИЕ 2

Мне нужна функция для сравнения с БПФ. Я не совсем доволен этой версией, я хотел использовать sapply(), но она работает. Пока сойдет.. Буду работать над улучшениями.

conv3 <- function(x, h){
m <- length(x)
n <- length(h)
X <- c(x, rep(floor(n/2), 0, floor(n/2)))   
H <- c(h, rep(floor(m/2), 0, floor(m/2)))
   Y <- vector()

for(i in 1:n+m-1){
    Y[i] <- 0 
    for(j in 1:m){
        Y[i] <- ifelse(i-j+1>0, Y[i] + X[j]*H[i-j+1], 0)
    }
}
Y[is.na(Y)] <- 0
Y[ceiling(m/2):(m+floor(m/2))]
}

Далее, я думаю, мне нужно поработать над тем, чтобы сделать его многомерным.

person Community    schedule 02.02.2013