I am trying to use matched filters for signal detection and advanced processing - to get additional information from the signal after the convolution with various templates.
The task I have is pretty simple - to find patterns with known forms in the signal I have as an input - but I need to be able to distinguish the results - e.g. detect signals with various width so the signals of similar forms but with another width will not be detected or result of filtering will have much lower amplitude comparing it to the result where we have the exact match of signal & pattern we are trying to detect.
The problem is that after the convolution of the signal & template, I have similar amplitudes of the resulting convolution for different shapes of input signals - so it's pretty hard to interpret - which signal is more similar to the template using convolution results.
The code I use in Python to do the matched filtering:
def matched_filter(template, signal):
# removing the DC from signal and template
# template = template - template.mean()
# signal = signal - signal.mean()
# Compute the conjugate of the template
conjugate_template = np.conjugate(template)
# Flip the template for zero-phase delay
flip_template = np.flip(conjugate_template)
# Compute the convolution of the signal and flipped template
convolution = np.convolve(signal, flip_template, mode='same')
return convolution
The code for toy-demonstration of a problem.
# Define the clean Gaussian template
# template 1 - Gaussian width (1,1)
template_size=10
template_mean1 = 5
template_std1 = 1
template_A1 = 1
template = template_A1 * np.exp(-(np.arange(template_size) - template_mean1)**2 / (2 * template_std1**2))
# template 2 - Gaussian width (2, 1) - high amplitude the same width
template_mean2 = 5
template_std2 = 1
template_A2 = 2
template2 = template_A2 * np.exp(-(np.arange(template_size) - template_mean2)**2 / (2 * template_std2**2))
# template 3 - Gaussian width (1, 2) - low amplitude, wider
template_mean3 = 5
template_std3 = 2
template_A3 = 3
template3 = template_A3 * np.exp(-(np.arange(template_size) - template_mean3)**2 / (2 * template_std3**2))
# Generate a sum of Gaussians with various amplitudes as the input signal
# first
input_mean1 = 10
input_std1 = 1
A1 = 1
# second
input_mean2 = 20
input_std2 = 1
A2 = 2
# third
input_mean3 = 30
input_std3 = 2
A3 = 3
# Generate a not-noised sum of Gaussians signal
not_noised_signal = A1 * np.exp(-(np.arange(40) - input_mean1)**2 / (2 * input_std1**2)) + \
A2 * np.exp(-(np.arange(40) - input_mean2)**2 / (2 * input_std2**2)) + \
A3 * np.exp(-(np.arange(40) - input_mean3)**2 / (2 * input_std3**2))
# add some noise
#input_signal = not_noised_signal + np.random.normal(0, 0.5, 40)
input_signal = not_noised_signal
# Perform matched filtering with zero-phase delay
filtered_signal = matched_filter(template, input_signal)
filtered_signal2 = matched_filter(template2, input_signal)
filtered_signal3 = matched_filter(template3, input_signal)
# Plot the results
#signal 1 & template 1
plt.cla()
plt.plot(input_signal, label='Input Signal')
plt.plot(template, label='template A='+str(template_A1)+', s='+str(template_std1), linestyle='--')
plt.plot(filtered_signal, label='Filter result')
#plt.plot(not_noised_signal, label='Not Noised Signal')
plt.legend()
plt.show()
Convolving the modelled gaussians with the gaussians of various templates (various width and amplitudes) - I have very similar results - the resulting amplitude after the convolution is larger for parts of a signal with the largest amplitude, but not for the parts of a signal where you have exact match of the template form with the signal form. See the figures for explanation:
Why is that?
I suppose I am wrong in understanding the results of convolution. Does the convolution result (matched filtering) show how similar is current part of the signal under test with the template used? Or it's not like that?
In my case - I need to find the best matches - to show the positions where the signal has the most similarity with the template used. So the amplitude and form of a Gaussian in this case must be taken into account - the max. of filtering must be on the positions where the pattern really is. I need to get the positions of possible peaks, along with the information about some properties of a signal - like the parameters of the detected signal using the convolution/filtering result.
Gaussian is taken only as an example. A real application is very noised and the normal distribution of a noise is assumed,
Signals, that can be detected have the form of Lorentzians with various parameters - so in real application various amplitudes/widths are present in the signal - that is why I am trying to build dome logic to distinguish one peak type from another.
P.S. Another thing I need to understand - is how to deal with the nonuniform grid in this case? I have non-equidistant samples, dx <> const
x(i+1) - x(i) <> x(i) - x(i-1)
While for the template I have only a model - how can I apply matched filtering in this case?
UPDATE 1. Using proposed approach,
def matched_filter(template, signal):
# removing the DC from signal and template
template = template - template.mean()
signal = signal - signal.mean()
temp_std = np.std(template)
signal_std = np.std(signal)
# Compute the conjugate of the template
conjugate_template = np.conjugate(template)
# Flip the template for zero-phase delay
flip_template = np.flip(conjugate_template)
# Compute the convolution of the signal and flipped template
convolution = np.convolve(signal, flip_template, mode='same')
convolution = convolution / (temp_std * signal_std)
return convolution
But in this case I get "rescaled" version of what I had in the input signal, not the magnitude where we have the most similarity.