sal - recorrer una imagen opencv python

python: filtro de paso de banda de una imagen

Question

Tengo una imagen de datos con un artefacto de imagen que aparece como un fondo sinusoidal, que quiero eliminar. Como es una onda sinusoidal de una sola frecuencia, parece natural que la transformada de Fourier y el filtro de paso de banda o "filtro de muesca" (donde creo que usaría un filtro gaussiano en + -omega)

Al tratar de hacer esto, noto dos cosas:

1) simplemente al realizar el pie y atrás, he reducido el componente de onda sinusoidal, que se muestra a continuación. Parece que hay un filtrado de paso alto de los datos solo por ir y volver.

import numpy as np

f = np.fft.fft2(img)                  #do the fourier transform
fshift1 = np.fft.fftshift(f)          #shift the zero to the center
f_ishift = np.fft.ifftshift(fshift1)  #inverse shift
img_back = np.fft.ifft2(f_ishift)     #inverse fourier transform
img_back = np.abs(img_back)

Esta es una imagen de img_back:

Tal vez el filtro aquí sea lo suficientemente bueno para mí, pero no tengo tanta confianza en ello ya que no entiendo bien la supresión de fondo.

2) Para estar más seguro de la supresión en las frecuencias no deseadas, hice una máscara booleana ''bandpass'' y la apliqué a los datos, pero la transformada de Fourier ignora la máscara.

a = shape(fshift1)[0]
b = shape(fshift1)[1]

ro = 8
ri = 5
y,x = np.ogrid[-a/2:a/2, -b/2:b/2] 
m1 = x*x + y*y >= ro*ro 
m2 = x*x + y*y <= ri*ri
m3=np.dstack((m1,m2))       
maskcomb =[]
for r in m3:
    maskcomb.append([any(c) for c in r])  #probably not pythonic, sorry
newma = np.invert(maskcomb)
filtdat = ma.array(fshift1,mask=newma) 
imshow(abs(filtdat))
f_ishift = np.fft.ifftshift(filtdat) 
img_back2 = np.fft.ifft2(f_ishift) 
img_back2 = np.abs(img_back2)

Aquí el resultado es el mismo que antes, porque np.fft ignora las máscaras. La solución a eso era simple:

filtdat2 = filtdat.filled(filtdat.mean())

Desafortunadamente, (pero al reflexionar también como es de esperar) el resultado se muestra aquí:

El gráfico de la izquierda es de la amplitud de la FFT, con el filtro de paso de banda aplicado. Es el anillo oscuro alrededor del componente central (DC). La fase no se muestra.

Claramente, el filtro ''brickwall'' no es la solución correcta. El fenómeno de hacer anillos a partir de este filtro está bien explicado aquí: ¿Qué sucede cuando aplica un filtro de pared de ladrillo a un conjunto de datos 1D?

Así que ahora estoy atascado. Tal vez sería mejor utilizar uno de los métodos scipy integrados, pero parecen ser para datos 1d, como en esta implementación de un filtro Butterworth . Posiblemente lo correcto es usar fftconvolve () como se hace aquí para desenfocar una imagen. Mi pregunta sobre fftconvolve es la siguiente: ¿Requiere que ambas ''imágenes'' (la imagen y el filtro) estén en el espacio real? Creo que sí, pero en el ejemplo usan un gaussiano, por lo que es ambiguo (fft (gaussian) = gaussian). Si es así, entonces parece incorrecto tratar de hacer un filtro de paso de banda de espacio real. Tal vez la estrategia correcta use convolve2d () con la imagen del espacio de Fourier y un filtro casero. Si es así, ¿sabes cómo hacer un buen filtro 2d?

Answer 1

Por lo tanto, un problema aquí es que su sinusoide de fondo tiene un período no muy diferente de los componentes de la señal que está tratando de preservar. es decir, el espaciado de los picos de la señal es aproximadamente el mismo que el período del fondo. Esto va a dificultar el filtrado.

Mi primera pregunta es si este fondo es realmente constante de un experimento a otro, o ¿depende de la muestra y la configuración experimental? Si es constante, entonces la resta del cuadro de fondo funcionaría mejor que el filtrado.

La mayoría de las funciones de filtro de señal de scipy.signal (bessel, chebychev, etc.) están, como usted dice, diseñadas para datos 1-D. Pero puede extenderlos fácilmente al filtrado isotrópico en 2-D. Cada filtro en el espacio de frecuencia es una función racional de f. Las dos representaciones son [a, b] que son los coeficientes del numerador y del denominador polinomial, o [z, p, k] que es la representación factorizada del polinomio, es decir: H(f) = k(f-z0)*(f-z1)/(f-p0)*(f-p1) Puedes tomar el polinomio de uno de los algoritmos de diseño del filtro, evaluarlo como una función de sqrt (x ^ 2 + y ^ 2) y aplicar a sus datos de dominio de frecuencia.

¿Puedes publicar un enlace a los datos de la imagen original?