# anyone have an answer to this?

Started by April 7, 2008
```On Apr 8, 2:02&#2013266080;pm, maury <maury...@core.com> wrote:
> On Apr 7, 4:31&#2013266080;pm, "bronx" <branko_blagoje...@hotmail.com> wrote:
>
>
>
>
>
> > >You should post the equations for the frequencies, length of signals,
> > >size of FFT, ...
>
> > >I suspect this method is going to run into problems where there are a
> > >lot of FFT bins where there is virtually no signal, such that the
> > >phase is really based more on numerical noise than any signal
> > >present. &#2013266080;In that case the phase difference between these bins would
> > >not be dependent on the delay and would not contribute to your
> > >result. &#2013266080;I don't know if you have this problem or not, but it would
> > >somethng to check for.
>
> > >Dirk
>
> > here is my code in detail
>
> > freq1 = 0.5e3; &#2013266080; &#2013266080; %Frequency of incoming source in Hz
> > time = 20e-3; &#2013266080; %Total time of simulation in s
> > tau1 = 0.375e-3; &#2013266080; &#2013266080; %time difference in arrival of 2 signals
>
> > N = pow2(11);
> > %t = (0:N-1)*time;
> > t = linspace(0,time,N);
>
> > % Signals arriving at input of microphones...with different phase
>
> > % MIC A
> > x1 = sin(2*pi*freq1*t);
>
> > % MIC B
> > x2 = sin(2*pi*freq1*(t + tau1));
>
> > figure(1);
> > plot(t, x1, "r");
> > hold on;
> > plot(t, x2, "g");
>
> > x1pad = [(x1) zeros(1,N-1)];
> > x2pad = [(x2) zeros(1,N-1)];
>
> > XR_CC_num = X1.*conj(X2);
> > XR_CC_den = abs(XR_CC_num);
> > XR = (XR_CC_num)./(XR_CC_den);
> > %Phase Correlation
> > PhaseCorr = fftshift((ifft(XR)));
> > %Cross Correlation
> > CrossCorr = fftshift(ifft(XR_CC_num));
> > CorrTime = time*linspace(-1,1, NPC);
>
> > figure(2);
> > %Plot Phase correlation in ms
> > plot(1e3*CorrTime,real(PhaseCorr),'r');
>
> > figure(3);
> > %Plot Cross Correlation in ms
> > plot(1e3*CorrTime,real(CrossCorr),'g');
>
> > %the location of the peak of each result is the time delay between signal
> > %%x1 and signal x2
>
> > let me know if you know why the peaks don't show the same time delay
> > index.
>
> > Thanks
>
> x1 is a delayed version of x2 by tau1. &#2013266080;However, the question I think
> you should be asking yourself is why you think the FFT of x1 would be
> different than the FFT of x2.
>
> Maurice Givens- Hide quoted text -
>
> - Show quoted text -

Do you think they should be the same?

You will notice that the OP never takes the FFT of x1 and x2.

Dirk
```
```bronx wrote:
>> You should post the equations for the frequencies, length of signals,
>> size of FFT, ...
>>
>> I suspect this method is going to run into problems where there are a
>> lot of FFT bins where there is virtually no signal, such that the
>> phase is really based more on numerical noise than any signal
>> present.  In that case the phase difference between these bins would
>> not be dependent on the delay and would not contribute to your
>> result.  I don't know if you have this problem or not, but it would
>> somethng to check for.
>>
>> Dirk
>>
>
> here is my code in detail
>
> freq1 = 0.5e3;     %Frequency of incoming source in Hz

.
.
.

>
> let me know if you know why the peaks don't show the same time delay
> index.
>
> Thanks

bronx,

I had seen this in a few years ago in a toy radar problem. So I thought
I'd stab at this one. Two main points: noise and choice of signal.

=Noise=
At first I thought it might be due to some numerical issue as Dirk
suggested so I added some noise, but to no avail. I left the noise code
in but notice the magnitude is set to a very small level compare to the
signal. Some noise is necessary to avoid any dividing by too small of
numbers. Try it with zero and see the difference. Also experiment with
larger noise to see how the accuracy decreases.

=Signal selection=
I believe the main issue with your example is the signal.  From my
understanding you need a somewhat richer signal than a perfect cosine.

What's a simple signal with rich frequency content?
A pulse, which is what some simple systems use (I've seen it in cheap
ultrasound). Better than that would be a chirp signal. If you are
looking at audio source location, most audio should be rich enough.

=Misc=
On a side note, for efficiency it would be ideal if you made the FFT
lengths a power of two.  Currently they are an odd length, although it
looks you were aiming at that with the definition of N.

Cheers,
J Elms

PS: Maybe I shouldn't just give you the answer, but I'm feeling generous.

Modified code:

freq1 = 0.5e3;     %Frequency of incoming source in Hz
time = 20e-3;   %Total time of simulation in s
tau1 = 0.375e-3;     %time difference in arrival of 2 signals

N = pow2(11);
%t = (0:N-1)*time;
t = linspace(0,time,N);

% Signals arriving at input of microphones...with different phase
NoiseMag  = 1e-5;
tPulseLen = 1e-3; % 1ms length pulse
nPulseLen = round(tPulseLen/t(2)); % convert to samples

BaseSig1 = ( (t>0) & (t<0+tPulseLen) );
BaseSig2 = ( (t>tau1) & (t<tau1+tPulseLen) );

% MIC A
x1 = BaseSig1 + NoiseMag*rand(1,N);

% MIC B
x2 = BaseSig2 + NoiseMag*rand(1,N);  %delay

% figure(1);
% plot(t, x1, "r");
% hold on;
% plot(t, x2, "g");

XR_CC_num = X1.*conj(X2);
XR_CC_den = abs(XR_CC_num);
XR = (XR_CC_num)./(XR_CC_den);
%Phase Correlation
PhaseCorr = real(fftshift((ifft(XR))));
%Cross Correlation
CrossCorr = real(fftshift(ifft(XR_CC_num)));
CorrTime = time*linspace(-1,1, NPC);

% figure(2);
% %Plot Phase correlation in ms
% plot(1e3*CorrTime,real(PhaseCorr),'r');

% figure(3);
% %Plot Cross Correlation in ms
% plot(1e3*CorrTime,real(CrossCorr),'g');

%the location of the peak of each result is the time delay between signal
%%x1 and signal x2

[tmp,indPhase] = max(PhaseCorr);
[tmp,indCross] = max(CrossCorr);
tPhase = CorrTime(indPhase)
tCross = CorrTime(indCross)

figure;

subplot(3,1,1)
hold on

subplot(3,1,2)
plot(abs(X1),':')
hold on
plot(abs(X2),'--')

subplot(3,1,3)
plot(CrossCorr/max(CrossCorr) )
hold on
plot(PhaseCorr/max(PhaseCorr) )
legend('Normalized CrossCorr','Normalized PhaseCorr')
```
```On Apr 9, 4:23&#2013266080;pm, J Elms <jeffrey.e...@gatech.edu> wrote:
> bronx wrote:
> >> You should post the equations for the frequencies, length of signals,
> >> size of FFT, ...
>
> >> I suspect this method is going to run into problems where there are a
> >> lot of FFT bins where there is virtually no signal, such that the
> >> phase is really based more on numerical noise than any signal
> >> present. &#2013266080;In that case the phase difference between these bins would
> >> not be dependent on the delay and would not contribute to your
> >> result. &#2013266080;I don't know if you have this problem or not, but it would
> >> somethng to check for.
>
> >> Dirk
>
> > here is my code in detail
>
> > freq1 = 0.5e3; &#2013266080; &#2013266080; %Frequency of incoming source in Hz
>
> .
> .
> .
>
>
>
> > let me know if you know why the peaks don't show the same time delay
> > index.
>
> > Thanks
>
> bronx,
>
> I had seen this in a few years ago in a toy radar problem. So I thought
> I'd stab at this one. Two main points: noise and choice of signal.
>
> =Noise=
> At first I thought it might be due to some numerical issue as Dirk
> suggested so I added some noise, but to no avail. I left the noise code
> in but notice the magnitude is set to a very small level compare to the
> signal. Some noise is necessary to avoid any dividing by too small of
> numbers. Try it with zero and see the difference. Also experiment with
> larger noise to see how the accuracy decreases.
>
> =Signal selection=
> I believe the main issue with your example is the signal. &#2013266080;From my
> understanding you need a somewhat richer signal than a perfect cosine.
>
> What's a simple signal with rich frequency content?
> A pulse, which is what some simple systems use (I've seen it in cheap
> ultrasound). Better than that would be a chirp signal. If you are
> looking at audio source location, most audio should be rich enough.
>
> =Misc=
> On a side note, for efficiency it would be ideal if you made the FFT
> lengths a power of two. &#2013266080;Currently they are an odd length, although it
> looks you were aiming at that with the definition of N.
>
> Cheers,
> J Elms
>
> PS: Maybe I shouldn't just give you the answer, but I'm feeling generous.
>
> Modified code:
>
> freq1 = 0.5e3; &#2013266080; &#2013266080; %Frequency of incoming source in Hz
> time = 20e-3; &#2013266080; %Total time of simulation in s
> tau1 = 0.375e-3; &#2013266080; &#2013266080; %time difference in arrival of 2 signals
>
> N = pow2(11);
> %t = (0:N-1)*time;
> t = linspace(0,time,N);
>
> % Signals arriving at input of microphones...with different phase
> NoiseMag &#2013266080;= 1e-5;
> tPulseLen = 1e-3; % 1ms length pulse
> nPulseLen = round(tPulseLen/t(2)); % convert to samples
>
> BaseSig1 = ( (t>0) & (t<0+tPulseLen) );
> BaseSig2 = ( (t>tau1) & (t<tau1+tPulseLen) );
>
> % MIC A
> x1 = BaseSig1 + NoiseMag*rand(1,N);
>
> % MIC B
> x2 = BaseSig2 + NoiseMag*rand(1,N); &#2013266080;%delay
>
> % figure(1);
> % plot(t, x1, "r");
> % hold on;
> % plot(t, x2, "g");
>
>
> XR_CC_num = X1.*conj(X2);
> XR_CC_den = abs(XR_CC_num);
> XR = (XR_CC_num)./(XR_CC_den);
> %Phase Correlation
> PhaseCorr = real(fftshift((ifft(XR))));
> %Cross Correlation
> CrossCorr = real(fftshift(ifft(XR_CC_num)));
> CorrTime = time*linspace(-1,1, NPC);
>
> % figure(2);
> % %Plot Phase correlation in ms
> % plot(1e3*CorrTime,real(PhaseCorr),'r');
>
> % figure(3);
> % %Plot Cross Correlation in ms
> % plot(1e3*CorrTime,real(CrossCorr),'g');
>
> %the location of the peak of each result is the time delay between signal
> %%x1 and signal x2
>
> [tmp,indPhase] = max(PhaseCorr);
> [tmp,indCross] = max(CrossCorr);
> tPhase = CorrTime(indPhase)
> tCross = CorrTime(indCross)
>
> figure;
>
> subplot(3,1,1)
> hold on
>
> subplot(3,1,2)
> plot(abs(X1),':')
> hold on
> plot(abs(X2),'--')
>
> subplot(3,1,3)
> plot(CrossCorr/max(CrossCorr) )
> hold on
> plot(PhaseCorr/max(PhaseCorr) )
> legend('Normalized CrossCorr','Normalized PhaseCorr')- Hide quoted text -
>
> - Show quoted text -

Not only is the signal richer in frequency content, but over the
interval they are defined, your delayed signal is really the original
signal delayed, which was not true for the OP's signal. That makes the
FFTs the same except for an approximately linear phase term which is
important for the method to work. It is not quite linear phase due to
the essentially random phase near zeros in the spectrum, but most of
the phase values should be correct.

Dirk
```
```dbell wrote:
> On Apr 9, 4:23 pm, J Elms <jeffrey.e...@gatech.edu> wrote:
>> bronx wrote:
>>>> You should post the equations for the frequencies, length of signals,
>>>> size of FFT, ...
>>>> I suspect this method is going to run into problems where there are a
>>>> lot of FFT bins where there is virtually no signal, such that the
>>>> phase is really based more on numerical noise than any signal
>>>> present.  In that case the phase difference between these bins would
>>>> not be dependent on the delay and would not contribute to your
>>>> result.  I don't know if you have this problem or not, but it would
>>>> somethng to check for.
>>>> Dirk
>>> here is my code in detail
>>> freq1 = 0.5e3;     %Frequency of incoming source in Hz
>> .
>> .
>> .
>>
>>
>>
>>> let me know if you know why the peaks don't show the same time delay
>>> index.
>>> Thanks
>> bronx,
>>
>> I had seen this in a few years ago in a toy radar problem. So I thought
>> I'd stab at this one. Two main points: noise and choice of signal.
>>
>> =Noise=
>> At first I thought it might be due to some numerical issue as Dirk
>> suggested so I added some noise, but to no avail. I left the noise code
>> in but notice the magnitude is set to a very small level compare to the
>> signal. Some noise is necessary to avoid any dividing by too small of
>> numbers. Try it with zero and see the difference. Also experiment with
>> larger noise to see how the accuracy decreases.
>>
>> =Signal selection=
>> I believe the main issue with your example is the signal.  From my
>> understanding you need a somewhat richer signal than a perfect cosine.
>>
>> What's a simple signal with rich frequency content?
>> A pulse, which is what some simple systems use (I've seen it in cheap
>> ultrasound). Better than that would be a chirp signal. If you are
>> looking at audio source location, most audio should be rich enough.
>>
>> =Misc=
>> On a side note, for efficiency it would be ideal if you made the FFT
>> lengths a power of two.  Currently they are an odd length, although it
>> looks you were aiming at that with the definition of N.
>>
>> Cheers,
>> J Elms
>>
>> PS: Maybe I shouldn't just give you the answer, but I'm feeling generous.
>>
>> Modified code:
>>
>> freq1 = 0.5e3;     %Frequency of incoming source in Hz
>> time = 20e-3;   %Total time of simulation in s
>> tau1 = 0.375e-3;     %time difference in arrival of 2 signals
>>
>> N = pow2(11);
>> %t = (0:N-1)*time;
>> t = linspace(0,time,N);
>>
>> % Signals arriving at input of microphones...with different phase
>> NoiseMag  = 1e-5;
>> tPulseLen = 1e-3; % 1ms length pulse
>> nPulseLen = round(tPulseLen/t(2)); % convert to samples
>>
>> BaseSig1 = ( (t>0) & (t<0+tPulseLen) );
>> BaseSig2 = ( (t>tau1) & (t<tau1+tPulseLen) );
>>
>> % MIC A
>> x1 = BaseSig1 + NoiseMag*rand(1,N);
>>
>> % MIC B
>> x2 = BaseSig2 + NoiseMag*rand(1,N);  %delay
>>
>> % figure(1);
>> % plot(t, x1, "r");
>> % hold on;
>> % plot(t, x2, "g");
>>
>>
>> XR_CC_num = X1.*conj(X2);
>> XR_CC_den = abs(XR_CC_num);
>> XR = (XR_CC_num)./(XR_CC_den);
>> %Phase Correlation
>> PhaseCorr = real(fftshift((ifft(XR))));
>> %Cross Correlation
>> CrossCorr = real(fftshift(ifft(XR_CC_num)));
>> CorrTime = time*linspace(-1,1, NPC);
>>
>> % figure(2);
>> % %Plot Phase correlation in ms
>> % plot(1e3*CorrTime,real(PhaseCorr),'r');
>>
>> % figure(3);
>> % %Plot Cross Correlation in ms
>> % plot(1e3*CorrTime,real(CrossCorr),'g');
>>
>> %the location of the peak of each result is the time delay between signal
>> %%x1 and signal x2
>>
>> [tmp,indPhase] = max(PhaseCorr);
>> [tmp,indCross] = max(CrossCorr);
>> tPhase = CorrTime(indPhase)
>> tCross = CorrTime(indCross)
>>
>> figure;
>>
>> subplot(3,1,1)
>> hold on
>>
>> subplot(3,1,2)
>> plot(abs(X1),':')
>> hold on
>> plot(abs(X2),'--')
>>
>> subplot(3,1,3)
>> plot(CrossCorr/max(CrossCorr) )
>> hold on
>> plot(PhaseCorr/max(PhaseCorr) )
>> legend('Normalized CrossCorr','Normalized PhaseCorr')- Hide quoted text -
>>
>> - Show quoted text -
>
> Not only is the signal richer in frequency content, but over the
> interval they are defined, your delayed signal is really the original
> signal delayed, which was not true for the OP's signal. That makes the
> FFTs the same except for an approximately linear phase term which is
> important for the method to work. It is not quite linear phase due to
> the essentially random phase near zeros in the spectrum, but most of
> the phase values should be correct.
>
> Dirk

That's a good point, so I tried it with a real audio clip. For
simulation purposes, I picked samples in an interlaced fashion so no
samples were exactly overlapped between the signals (highly correlated,
but not identical). I got accuracy within .25 ms and within 1 sample
agreement between the two methods.

I did encounter some anomalies with it detecting a zero delay when
longer delays were used. But this begins to break the original
assumption of this algorithm, which is that for the windows selected the
data is highly correlated.

However, even in these long delay cases there was also a significant
peak at the correct location.  An algorithm for peak detection could be
used to improve results if these cases are encountered.

Notes: With a real signal I did not have to add any noise. I also tried
without padding the signals and with windowing the signal.  These only
had minor effects on the results.

J Elms

PS: Sorry if I ramble.
```
```But why does a cross correlation plot contain information of the proper
delay between the signals, while the phase correlation does not?? If both
didn't give me anything, I would agree with you guys, but the fact that
one is giving something reasonable and the other is not makes me think
that there is something i am overlooking in the matlab coding, not in the
signals used.
```
```>But why does a cross correlation plot contain information of the proper
>delay between the signals, while the phase correlation does not?? If
both
>didn't give me anything, I would agree with you guys, but the fact that
>one is giving something reasonable and the other is not makes me think
>that there is something i am overlooking in the matlab coding, not in
the
>signals used.
>

whoops, sorry guys, i didnt realize these responses came to 2 pages...this
is actually a response to the last response on the first page....disregard
it, and i will respond with more questions after i try your mod...

thanks very much,

bronx
```