4
$\begingroup$

Problem Description

  1. I have two devices which are measuring the same inputs and are using them to calculate the same (non-periodic) output signal. The process generating the inputs is deterministic.
  2. Both devices are sampling with different rates (each > 5 kHz).
  3. The capture duration of one device may be several seconds longer than the other one. Total capture duration is > 60s.
  4. Since the calculation of the output signals is done with different algorithms they are looking quite different at first sight but they do look very similar after applying some sort of linear low-pass filter.

My goal is to find out the time delay between those two output signals x and y.

Proposal for Solution

  1. Crop the longer signal such that x and y have the same capture duration.
  2. Resample x and y such that they have the same number of samples.
  3. Filter x and y with the same linear low-pass filter.
  4. Remove linear trends by applying MATLAB's detrend function to x and y.
  5. Apply MATLAB's finddelay on x and y.

Questions

  • May cropping the longer signal corrupt the result of the cross-correlation?
  • What is the best approach for resampling x and y?
  • As I mentioned, x and y are only looking very similar after filtering. Is filtering signals prior cross-correlation advisable? If yes, are there recommended filter types for that purpose?

  • Is it a must to detrend data prior cross-correlation or might it be enough to remove the mean? In my case there may be multiple linear trends and it could be hard to find all breakpoints as depicted here.


Desired Time Delay between two (filtered) Output Signals

enter image description here

$\endgroup$
  • $\begingroup$ Could you share some signals to play with? $\endgroup$ – Royi Sep 24 '17 at 19:11
  • $\begingroup$ We need both signals in the same MAT file. I can do nothing with only one of the signals. $\endgroup$ – Royi Sep 25 '17 at 6:32
  • $\begingroup$ I will generate simpler simulation for the case you presented. $\endgroup$ – Royi Sep 25 '17 at 12:08
  • 1
    $\begingroup$ @PeterK. I thought it was some sort of PWM in the beginning and I was using Octave. The signals do look like that. lR8n6i Thank you. I added one more edit along the idea of aligning those specific "landmarks", not the signal as a whole. I still cannot see how these two signals align by the way. Can you please post a subplot type of graph marking on the graph of one waveform how it corresponds with the other? ">this< should match with >this<" kind of thing, not parallel lines. $\endgroup$ – A_A Sep 25 '17 at 22:07
  • 1
    $\begingroup$ Related question Time delay estimation in low SNR $\endgroup$ – Olli Niemitalo Sep 27 '17 at 6:06
3
+50
$\begingroup$

What you are describing is called (/ has been solved via) "Dynamic Time Warping".

In DTW, you calculate a series of simple operations between two sequences that if applied on one sequence they would "bring it" directly "over" the other. Of course, there might be more than one ways to do this in which case you are selecting the one with the least cost. This "cost of distortion" is now your similarity.

For more information on DTW, please see this, this and this. More relevant to your application might be this publication. For a practical example, please see this link.

With regards to some of the rest of your questions:

  • May cropping the longer signal corrupt the result of the cross-correlation?

Yes, it would...

  • What is the best approach for resampling x and y?

...because both cropping and re-sampling assume a linear speed up. This might not be the case across the entire length of your signal.

  • As I mentioned, x and y are only looking very similar after filtering. Is filtering signals prior cross-correlation advisable?

It depends a lot on the application. However, filtering two signals $x,y$ with the same impulse response $h$ is expected to increase their cross correlation because after convolution, both $x,y$ time instances would now have the same $h$ imprint on them.

EDIT: Following a brief review of the data provided with the question, it seems that these are some sort of counts (?) or other digital / impulsive data. Therefore, the question "...is filtering signals prior to cross-correlation advisable?" is ill posed. The fact remains that if you filter two signals through the same filter, their cross-correlation will increase. But here, it seems that the signals do not exist in their final form yet. To clarify: We have to agree on which one is the signal. If you were dealing (for example) with a pulse width modulation (PWM) derived sound signal, we would not try to correlate the two PWM pulse trains that produce it but their integrated versions because THAT is what we consider as our signal. Consequently, in this particular case we are dealing with here, we have to understand what is considered as "signal". The only thing we know about it is that it comes out of some process. That's fine but, at the same time, if you have to integrate the data to PRODUCE the signal of interest then the process itself or its context would tell you what sort of integration you have to apply. It might not even be called "low pass filtering". In conclusion: If your signal has to be integrated first, then don't try to compare the "pulse sequences". To get back to the PWM-Audio case, it is like asking "Should I integrate the PWM pulses before cross correlation or would this spoil my signal?". Of course you have to apply the integration because otherwise, you don't have the signal to work with. Decide on the parameters of integration, apply them to your original signal and then, if you still have drifts and compressions and dilations of features in the signal, see if you can use DWT to establish the time delay.

  • Is it a must to detrend data prior cross-correlation or might it be enough to remove the mean?

Again, depends on the application and the amount of drift expected in the data but in general, it is a good idea to remove any long term drifts to avoid the cross-correlation being thrown off by the fact that the trends seem to match despite those little "details" on top of them. Where "details", of course, might be your desired signal. This is illustrative but a hypothetical scenario, your data might not even include such a case, it simply shows how detrending might help.

Hope this helps.

EDIT (2):

In the light of more information, I have to say that this does not look like a task for DTW neither cross-correlation the way it is described.

What is attempted here is an estimation of the delay between two waveforms by aligning "landmark" points. NOT the waveforms themselves. Personally, I am struggling to see why the marked points in the signal should align at all. Therefore, what I would suggest is this:

  1. Pick the higher frequency waveform
  2. Mark the peak you are after
  3. Extract $\pm m$ points (or ms) around it. For example, if the peak is at sample 4500, extract x[4500-200:4500+200]. Call this your template $q$.
  4. Subtract the mean of $q$ from $q$
  5. Downsample $q$ so that it matches the sampling rate of the lower frequency signal. Call this $q_d$.

This is now your template. It is shorter than the signal you are looking for and you can "slide" it over the lower rate signal and estimate the cross correlation. At the point that the cross correlation peaks is where this second waveform is, in the lower rate signal. Prior to applying cross correlation, you might want to be subtracting the mean of the signal portion you are about to compare with your template or instead of removing the mean, you might want to be applying detrending. Because the duration of the signal is small, simple linear detrending might work good enough.

This is an attempt to align specific points of interest in the signal. Not the waveform as a whole

Hope this helps.

EDIT (3): Ideally, a large amount of the information provided in the comments here, should appear in the original question and then perhaps we can delete all of these comments.

If you can control the load of the inverter, maybe you can insert a "sync event". Apply a known load on and off after the inverter has reached its nominal operation to create a "click" with a predictable response. You can use that "click" later to sync the two signals. If you have this capability, then you can average the responses from a number of "clicks" so that you have a good idea of what you are looking for, via cross correlation, as the "sync pulse".

Of course, this covers what happens when the two signals are supposed to have a constant time invariant delay.

$\endgroup$
  • $\begingroup$ Thank you for your answer. However, at least MATLAB's dtw function seems not to be able to handle my signals due to their lengths. Do you think that cross-correlation isn't applicable at all in my scenario? $\endgroup$ – lR8n6i Sep 24 '17 at 22:05
  • $\begingroup$ At the moment, I cannot access the provided signals. Perhaps down-sampling until satisfying MATLAB's restrictions will provide some indicative results in the short term (?). Difficult to answer, about the suitability of cross-correlation, without further knowledge of signals / process that generates them. Might amend response once I have the opportunity to review the provided signals. $\endgroup$ – A_A Sep 24 '17 at 22:13
  • $\begingroup$ Downsampling the signal with a factor of 1000 lets me execute dtw. However it seems to return no meaningful results if I add 100 leading zeros and compare the original signal with the shifted one. $\endgroup$ – lR8n6i Sep 24 '17 at 22:43
  • 1
    $\begingroup$ @lR8n6i I had a look at your signals. Are these the exact same ones that you want to filter prior to cross-correlation? You might need to pre-process them a little bit more before you do anything else with them. If these are meant to be pulses, you might have to apply integration. If not, there might be some overflow happening there. Difficult to say for sure without knowing more about the signal. If you can, could you please update the question with a couple of plots (but not full extent, e.g. x(1:1500)) and a little bit more about the signal itself?Agree with PeterK,1000 might be too much. $\endgroup$ – A_A Sep 25 '17 at 10:48
  • 1
    $\begingroup$ @lR8n6i I have to say that putting the signals up on a plot, I would not call them "very similar". So, it is not that there is a shortage of solutions but I feel that we are reaching the limits of "you could try this". Can you share a bit of background for the signals? No matter how "special interest" you might think it is. Try us :D There is a very wide diversity of talent here, it's wonderful. But maybe the interest is not immediate yet. Except of course if there are any IP issues at risk that prohibit you from revealing more (?) $\endgroup$ – A_A Sep 26 '17 at 14:51
2
$\begingroup$

A_A's answer is probably closest to what you want, but I thought I'd try to simulate estimating the delay using cross-correlation.

The picture below shows five plots

  • The top plot is an undelayed, corrupted, smoothed, and noisy pulse.

  • The second plot is a delayed version with a different noise realization.

  • The third plot is the cross-correlation of the undelayed and delayed noiseless signals. You can see that the peak is at the right place (100) but it's VERY wide, which means it will be hard to localize.

  • The fourth plot is if I apply a DC blocker before taking the cross-correlation. This improves localization.

  • The final plot is the cross-correlation of the DC blocked noisy undelayed and the DC blocked noisy delayed signals. It's a mess.

The bottom line is that I don't think cross-correlation, as-is, will help. I'll think about it a bit more and see if something else comes to mind.

Various attempts at using cross correlation


R Code Only Below

delayed_smoothed_step <- function(delay)
{
  N <-  10000
  C1 <- 1050
  C2 <- 850
  step_function <- c(rep(C1,N+delay), rep(C2,N-delay))

  alpha <- 0.999

  smoothed_step_function <- step_function

  for (t in seq(2,length(smoothed_step_function)))
  {
    smoothed_step_function[t] <- alpha*smoothed_step_function[t-1] + (1-alpha)*step_function[t-1]
  }

  return(smoothed_step_function)
}


dc_blocker <- function(x)
{
  alpha <- 0.5
  y <- rep(0,length(x))
  for (t in seq(2,length(x)))
  {
    y[t] <- alpha*y[t-1] + x[t] - x[t-1]
  }

  return(y)
}

s1 <- delayed_smoothed_step(0)
s2 <- delayed_smoothed_step(100)

n1 <- s1 + 100*rnorm(2*N)
n2 <- s2 + 100*rnorm(2*N)

ds1 <- dc_blocker(s1)
ds2 <- dc_blocker(s2)
dn1 <- dc_blocker(n1)
dn2 <- dc_blocker(n2)

par(mfrow=c(6,1))
plot(n1, col='red', type='l')
lines(s1)

plot(n2, col='red', type='l')
lines(s2)

ccf(s1[2:(2*N)]-s1[1:(2*N-1)],s2[2:(2*N)]-s2[1:(2*N-1)], lag.max=200)
ccf(ds1[2:(2*N)]-ds1[1:(2*N-1)],ds2[2:(2*N)]-ds2[1:(2*N-1)], lag.max=200)
ccf(ds1,ds2, lag.max=200)
ccf(n1[2:(2*N)]-n1[1:(2*N-1)],n2[2:(2*N)]-n2[1:(2*N-1)], lag.max=200)
#ccf(dn1[2:(2*N)]-dn1[1:(2*N-1)],dn2[2:(2*N)]-dn2[1:(2*N-1)], lag.max=200)
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.