2012-04-04 8 views
5

He creado un código python para suavizar una señal dada usando la transformada Weierstrass, que es básicamente la convolución de un gaussiano normalizado con una señal.¿Cómo eliminar los efectos límite que surgen debido al relleno cero en fft scipy/numpy?

El código es el siguiente:


#Importing relevant libraries 
from __future__ import division 
from scipy.signal import fftconvolve 
import numpy as np 

def smooth_func(sig, x, t= 0.002): 
    N = len(x) 
    x1 = x[-1] 
    x0 = x[0]  


# defining a new array y which is symmetric around zero, to make the gaussian symmetric. 
    y = np.linspace(-(x1-x0)/2, (x1-x0)/2, N) 
    #gaussian centered around zero. 
    gaus = np.exp(-y**(2)/t)  

#using fftconvolve to speed up the convolution; gaus.sum() is the normalization constant. 
    return fftconvolve(sig, gaus/gaus.sum(), mode='same') 

Si ejecuta este código para decir una función escalonada, se alisa la esquina, pero en el límite interpreta otra esquina y suaviza que también, como resultado, dando un comportamiento innecesario en el límite. Explico esto con una figura que se muestra en el siguiente enlace.
Boundary effects

Este problema no aparece si lo integramos directamente para encontrar la convolución. Por lo tanto, el problema no está en la transformación de Weierstrass, y de ahí que el problema esté en la función fftconvolucionar de scipy.

Para entender por qué surge este problema, primero debemos entender el funcionamiento de fftconvolve in scipy.
La función fftconvolve básicamente utiliza el teorema de convolución para acelerar el cálculo.
En resumen dice:
convolución (INT1, INT2) = IFFT (FFT (INT1) * FFT (INT2))
Si aplicamos directamente este teorema que no te dan el resultado deseado. Para obtener el resultado deseado, necesitamos tomar fft en una matriz del doble del tamaño de max (int1, int2). Pero esto conduce a los efectos de frontera no deseados. Esto se debe a que en el código fft, si el tamaño (int) es mayor que el tamaño (sobre el cual tomar fft) cero rellena la entrada y luego toma el fft. Este relleno cero es exactamente lo que es responsable de los efectos de frontera no deseados.

¿Puede sugerir una forma de eliminar estos efectos límite?

He intentado eliminarlo por un simple truco. Después de suavizar la función, estoy comparando el valor de la señal suavizada con la señal original cerca de los límites y, si no coinciden, reemplace el valor del func suavizado con la señal de entrada en ese punto.
Es como sigue:


i = 0 
eps=1e-3 
while abs(smooth[i]-sig[i])> eps: #compairing the signals on the left boundary 
    smooth[i] = sig[i] 
    i = i + 1 
j = -1 

while abs(smooth[j]-sig[j])> eps: # compairing on the right boundary. 
    smooth[j] = sig[j] 
    j = j - 1 

Hay un problema con este método, debido al uso de un epsilon hay pequeños saltos en la función smoothened, como se muestra a continuación:
jumps in the smooth func

¿Se pueden realizar cambios en el método anterior para resolver este problema de límite?

+0

Dupe of http://math.stackexchange.com/q/127875/2206 – endolith

Respuesta

3

Lo que produce un núcleo de filtro simétrico en los extremos depende de lo que suponga que los datos están más allá de los extremos.

Si no le gusta el aspecto del resultado actual, que supone ceros más allá de ambos extremos, intente extender los datos con otra suposición, por ejemplo, un reflejo de los datos, o la continuación de regresión polinomial. Extienda los datos en ambos extremos al menos a la mitad de la longitud del kernel de filtro (excepto si su extensión es cero, que viene de forma gratuita con el relleno de cero existente requerido para la convolución no circular). A continuación, elimine las extensiones de los extremos añadidos después de filtrar, y vea si le gusta el aspecto de su suposición. Si no, intente con otra suposición. O mejor aún, use datos reales más allá de los extremos si los tiene.

+0

Gracias, su respuesta abre una amplia gama de posibilidades. – Omkar

6

mejor enfoque es probablemente usar mode = 'valid':

The output consists only of those elements that do not rely on the zero-padding.

A menos que usted puede envolver su señal, o la señal que se está procesando es un extracto de una señal más grande (en cuyo caso: Proceso de señal completa luego recorta la región de interés) siempre vas a tener efectos de borde al hacer la convolución.Tienes que elegir cómo quieres lidiar con ellos. Usar mode = valid simplemente los corta, lo cual es una solución bastante buena. Si sabe que la señal es siempre "similar a un paso", puede extender el frente y el final de la señal procesada según corresponda.

+0

Gracias su solución funcionó. – Omkar

Cuestiones relacionadas