4 Correction to the some dodgy english.
source | link

The intention is to show how the cross-correlation can be performed for non-integer $\tau$, which is defined by the closure correlate_point. This uses the omega array, which describes the rotation of the complex phasor at each discrete frequency corresponding to a time-shift $\tau=1$. $\tau$ then scales this for each shift. It should be apparent that to maintain a real time signal, the rotations of the negative frequencies are just the negation of$-1$ times the rotations of the positive frequencies (for corresponding frequency pairs).

The intention is to show how the cross-correlation can be performed for non-integer $\tau$, which is defined by the closure correlate_point. This uses the omega array, which describes the rotation of the complex phasor at each discrete frequency corresponding to a time-shift $\tau=1$. $\tau$ then scales this for each shift. It should be apparent that to maintain a real time signal, the rotations of the negative frequencies are just the negation of the rotations of the positive frequencies.

The intention is to show how the cross-correlation can be performed for non-integer $\tau$, which is defined by the closure correlate_point. This uses the omega array, which describes the rotation of the complex phasor at each discrete frequency corresponding to a time-shift $\tau=1$. $\tau$ then scales this for each shift. It should be apparent that to maintain a real time signal, the rotations of the negative frequencies are just $-1$ times the rotations of the positive frequencies (for corresponding frequency pairs).

3 small correction to the period stuff.
source | link

The input variables a and b describe $f\left(n\right)$ and $g\left(n\right)$ for $n = \{0, 1, ... , N-1\}$ and are both assumed to be band limited and periodic with period $N-1$$N$ (the shift is implemented in the discrete Fourier domain). $\tau$ is in range $[-N+1, N-1]$.

The input variables a and b describe $f\left(n\right)$ and $g\left(n\right)$ for $n = \{0, 1, ... , N-1\}$ and are both assumed to be band limited and periodic with period $N-1$ (the shift is implemented in the discrete Fourier domain). $\tau$ is in range $[-N+1, N-1]$.

The input variables a and b describe $f\left(n\right)$ and $g\left(n\right)$ for $n = \{0, 1, ... , N-1\}$ and are both assumed to be band limited and periodic with period $N$ (the shift is implemented in the discrete Fourier domain). $\tau$ is in range $[-N+1, N-1]$.

2 Changed a variable name in the code to make it more intuitive.
source | link
import numpy
from numpy import fft
from scipy import optimize

def arg_max_corr(a, b):

    if len(a.shape) > 1:
        raise ValueError('Needs a 1-dimensional array.')

    length = len(a)
    if not length % 2 == 0:
        raise ValueError('Needs an even length array.')

    if not a.shape == b.shape:
        raise ValueError('The 2 arrays need to be the same shape')

    # Start by finding the coarse discretised arg_max
    coarse_max = numpy.argmax(numpy.correlate(a, b, mode='full')) - length+1

    omega = numpy.zeros(length)
    omega[0:length/2] = (2*numpy.pi*numpy.arange(length/2))/length
    omega[length/2+1:] = (2*numpy.pi*
            (numpy.arange(length/2+1, length)-length))/length

    fft_a = fft.fft(a)

    def correlate_point(tau):
        rotate_vec = numpy.exp(1j*tau*omega)
        rotate_vec[length/2] = numpy.cos(numpy.pi*tau)

        return numpy.sum((fft.ifft(fft_a*rotate_vec)).real*b)

    start_arg, end_arg = (float(coarse_max)-1, float(coarse_max)+1)

    min_argmax_arg = optimize.fminbound(lambda tau: -correlate_point(tau), 
            start_arg, end_arg)

    return min_argmax_arg
import numpy
from numpy import fft
from scipy import optimize

def arg_max_corr(a, b):

    if len(a.shape) > 1:
        raise ValueError('Needs a 1-dimensional array.')

    length = len(a)
    if not length % 2 == 0:
        raise ValueError('Needs an even length array.')

    if not a.shape == b.shape:
        raise ValueError('The 2 arrays need to be the same shape')

    # Start by finding the coarse discretised arg_max
    coarse_max = numpy.argmax(numpy.correlate(a, b, mode='full')) - length+1

    omega = numpy.zeros(length)
    omega[0:length/2] = (2*numpy.pi*numpy.arange(length/2))/length
    omega[length/2+1:] = (2*numpy.pi*
            (numpy.arange(length/2+1, length)-length))/length

    fft_a = fft.fft(a)

    def correlate_point(tau):
        rotate_vec = numpy.exp(1j*tau*omega)
        rotate_vec[length/2] = numpy.cos(numpy.pi*tau)

        return numpy.sum((fft.ifft(fft_a*rotate_vec)).real*b)

    start_arg, end_arg = (float(coarse_max)-1, float(coarse_max)+1)

    min_arg = optimize.fminbound(lambda tau: -correlate_point(tau), 
            start_arg, end_arg)

    return min_arg
import numpy
from numpy import fft
from scipy import optimize

def arg_max_corr(a, b):

    if len(a.shape) > 1:
        raise ValueError('Needs a 1-dimensional array.')

    length = len(a)
    if not length % 2 == 0:
        raise ValueError('Needs an even length array.')

    if not a.shape == b.shape:
        raise ValueError('The 2 arrays need to be the same shape')

    # Start by finding the coarse discretised arg_max
    coarse_max = numpy.argmax(numpy.correlate(a, b, mode='full')) - length+1

    omega = numpy.zeros(length)
    omega[0:length/2] = (2*numpy.pi*numpy.arange(length/2))/length
    omega[length/2+1:] = (2*numpy.pi*
            (numpy.arange(length/2+1, length)-length))/length

    fft_a = fft.fft(a)

    def correlate_point(tau):
        rotate_vec = numpy.exp(1j*tau*omega)
        rotate_vec[length/2] = numpy.cos(numpy.pi*tau)

        return numpy.sum((fft.ifft(fft_a*rotate_vec)).real*b)

    start_arg, end_arg = (float(coarse_max)-1, float(coarse_max)+1)

    max_arg = optimize.fminbound(lambda tau: -correlate_point(tau), 
            start_arg, end_arg)

    return max_arg
1
source | link