1#!/usr/bin/python
2
3# Copyright (C) 2012 The Android Open Source Project
4#
5# Licensed under the Apache License, Version 2.0 (the "License");
6# you may not use this file except in compliance with the License.
7# You may obtain a copy of the License at
8#
9#       http://www.apache.org/licenses/LICENSE-2.0
10#
11# Unless required by applicable law or agreed to in writing, software
12# distributed under the License is distributed on an "AS IS" BASIS,
13# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14# See the License for the specific language governing permissions and
15# limitations under the License.
16
17import numpy as np
18import numpy.linalg
19import scipy as sp
20import scipy.fftpack
21import scipy.signal
22import math
23import sys
24from multiprocessing import Pool
25
26def convolution(data0, data1reversed, n):
27    """calculate convolution part of data0 with data1 from pos n"""
28    N = len(data1reversed)
29    return np.dot(data0[n:N+n], data1reversed)
30
31
32def convolutionstar(args):
33    return convolution(*args)
34
35def calc_delay(data0, data1):
36    """Calcuate delay between two data. data0 is assumed to be recorded first,
37       and will have longer length than data1
38       returns delay between data0 and data1 in number of samples in data0's point of view"""
39
40    len0 = len(data0)
41    len1 = len(data1)
42    if len1 > len0:
43        print "data1 longer than data0"
44        return -1
45    searchLen = len0 - len1
46    data1reverse = data1[::-1]
47
48    # This is faster than signal.correlate as there is no need to process
49    # full data, but still it is slow. about 18 secs for data0 of 4 secs with data1 of 1 secs
50    print "***Caluclating delay, may take some time***"
51    gData0 = data0
52    gData1 = data1reverse
53    pool = Pool(processes = 4)
54    TASK = [(data0, data1reverse, i) for i in range(searchLen)]
55    result = pool.map(convolutionstar, TASK)
56
57    return np.argmax(result)
58
59
60# test code
61if __name__=="__main__":
62    samplingRate = 44100
63    durationInSec = 0.001
64    if len(sys.argv) > 1:
65        durationInSec = float(sys.argv[1])
66    signalFrequency = 1000
67    samples = float(samplingRate) * float(durationInSec)
68    index = np.linspace(0.0, samples, num=samples, endpoint=False)
69    time = index / samplingRate
70    multiplier = 2.0 * np.pi * signalFrequency / float(samplingRate)
71    data0 = np.sin(index * multiplier)
72    DELAY = durationInSec / 2.0 * samplingRate
73    data1 = data0[DELAY:]
74    delay = calc_delay(data0, data1)
75    print "calc_delay returned", delay, " while expecting ", DELAY
76
77
78