Faster Algorithm to convolve/correlate two sparse 1-D signals in python (or any language)

user_6910902 0
Viewed 264

I have two signals which I need to correlate or convolve. Each signal is sampled non-uniformly and the values of the signal I have with me are the timestamp and the magnitude of the signal at that timestamp. The values of the signal at all the other times can be assumed to be zero. The timestamps of the signal have a resolution running in microseconds.

An example of how the signal looks is shown below:

enter image description here

As can be seen, the resolution of the signal is in microseconds and signal is mostly sparse.

If I were to the convolve two signals of this type, I would first have to pad the signals with zeros (since I would have to discretise the signal). While the padding can be done with resolution of microseconds, The number of values to be multiplied becomes too big and the operation becomes increasingly slow. Most of the multiplications in this convolution would be multiplication of zeros(which are pretty much useless). I have therefore chosen a round off value of 2 places (0.xxxxxx becomes 0.xx),since I have to perform 40,000 similar convolutions. I have written my resampling function as shown below.

import numpy as np
import math

def resampled_signal_lists_with_zeros(signal_dict, xlimits):
  '''
      resamples the given signal with precision determined by the round function. 
      signal_dict is a dictionary with timestamp as key and signal magnitude is the value of the key. 
      xlimits is an array containing the start and stop time of the signal.
  '''
  t_stamps_list = list(signal_dict.keys())
  t_list = list(np.arange(int(math.floor(xlimits[0])), int(math.ceil(xlimits[1])), 0.005))
  t_list = [round(t, 2) for t in t_list]
  s_list = list()
  time_keys = [round(t, 2) for t in t_stamps_list]

  i = 0
  for t in t_list:
    if i < len(t_stamps_list):
      if t==time_keys[i]:
        s_list.append(signal_dict[t_stamps_list[i]])
        i+=1

      else:
        s_list.append(0)
    else:
      s_list.append(0)
  return t_list, s_list

The correlation of two signals padded in the above manner is done using scipy as follows:

from scipy.signal import correlate
output = correlate(s_1, s_2, mode='same')

The output calculated in the above manner is pretty slow .Since the signal is pretty sparse and most of the multiplications in the signal are multiplications of zero, I think there should be a better way to do the same operations. Is there a way to get the result of the convolutions of the two sparse signals faster?

convolution signal-processing python
1 Answers

Let f(x) and g(x) be "sparse" pulse functions defined with support in [0,inf). Discretize them over some linear mesh, so that f = [f0,f1,...], and g = [g0,g1,...]. Let l_f be the length of f, and l_g be the length of g.

If the signals are very sparse, like point measures, then this could be done combinatorically. Let indices f_sup = u1,u2,... be the sparse support of f, of length k_f, and g_sup = v1,v2,... be the sparse support of g, of length k_g. Then for convolution, the support points simply add, requiring k_f*k_g operations:

out_function_by_index = []
# (used for loops, but this could be done via comprehension for speed)
for ii in f_sup:
    for jj in g_sup:
        out_function_by_index.append((ii+jj,f[ii]*g[jj]))

To reconstitute the output to a function over the same linear mesh:

out_function = [0.0 for x in range(l_f+l_g-1)]
for tup in out_function:
    out_function[tup[0]] = tup[1]

Do a sign flip on f or g for the respective (non-commutative) cross-correlation.

If on the other hand the sparse signals are "piecewise sparse," you could pre-process the signals to identify sub-intervals that constitute sufficiently dense support, and do the combinatoric approach with convolution per sub-interval--ie just ordinary convolution between these dense-support sub-intervals, and add the result of all these inter-function sub-interval convolutions.

user_141794920
Related