-
Notifications
You must be signed in to change notification settings - Fork 1
/
_dtw.c
159 lines (132 loc) · 4.46 KB
/
_dtw.c
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
#include <Python.h>
#include <numpy/arrayobject.h>
#include "dtw.h"
/* Docstrings */
static char module_docstring[] =
"Dynamic time warping (DTW) distance C extension for Python";
static char dtw_docstring[] =
"Compute the dynamic time warping (DTW) distance between two sequences.\n"
"\n"
"Parameters\n"
"----------\n"
"x : numpy array of floats\n"
" First sequence\n"
"y : numpy array of floats\n"
" Second sequence\n"
"window_frac: float\n"
" Locality constraint, given as a fraction from 0 to 1 of the size of\n"
" the larger sequence.\n"
"\n"
"Returns\n"
"-------\n"
"float\n"
" The DTW distance between x and y";
static char dtw_path_docstring[] =
"Determine the optimal warping between two sequences.\n"
"\n"
"Parameters\n"
"----------\n"
"x : numpy array of floats\n"
" First sequence\n"
"y : numpy array of floats\n"
" Second sequence\n"
"window_frac: float\n"
" Locality constraint, given as a fraction from 0 to 1 of the size of\n"
" the larger sequence.\n"
"\n"
"Returns\n"
"-------\n"
"numpy.ndarray\n"
" The pairings between the two sequences that provide the optimal\n"
" warping path.";
/* Available functions */
static PyObject *dtw_dtw(PyObject *self, PyObject *args);
static PyObject *dtw_dtw_path(PyObject *self, PyObject *args);
/* Module specification */
static PyMethodDef module_methods[] = {
{"dtw", dtw_dtw, METH_VARARGS, dtw_docstring},
{"dtw_path", dtw_dtw_path, METH_VARARGS, dtw_path_docstring},
{NULL, NULL, 0, NULL}
};
/* Initialize the module */
PyMODINIT_FUNC init_dtw(void)
{
PyObject *m = Py_InitModule3("_dtw", module_methods, module_docstring);
if (m == NULL)
return;
import_array();
}
static PyObject *dtw_dtw(PyObject *self, PyObject *args)
{
double window_frac;
PyObject *x_obj, *y_obj;
/* Parse the input tuple */
if(!PyArg_ParseTuple(args, "OOd", &x_obj, &y_obj, &window_frac))
return NULL;
/* Interpret the input objects as numpy arrays. */
PyObject *x_array = PyArray_FROM_OTF(x_obj, NPY_DOUBLE, NPY_IN_ARRAY);
PyObject *y_array = PyArray_FROM_OTF(y_obj, NPY_DOUBLE, NPY_IN_ARRAY);
/* If that didn't work, throw an exception. */
if(x_array == NULL || y_array == NULL)
{
Py_XDECREF(x_array);
Py_XDECREF(y_array);
return NULL;
}
/* Get the number of data points in each array */
int xsize = (int)PyArray_DIM(x_array, 0);
int ysize = (int)PyArray_DIM(y_array, 0);
/* Get pointers to the data as C-types. */
double *x = (double*)PyArray_DATA(x_array);
double *y = (double*)PyArray_DATA(y_array);
/* Call the external C function to compute the DTW distance. */
double value = dtw(x, y, xsize, ysize, window_frac);
/* Clean up. */
Py_DECREF(x_array);
Py_DECREF(y_array);
if(value < 0.0)
{
PyErr_SetString(PyExc_RuntimeError,
"dtw returned an impossible value.");
return NULL;
}
/* Build the output tuple */
PyObject *ret = Py_BuildValue("d", value);
return ret;
}
static PyObject *dtw_dtw_path(PyObject *self, PyObject *args)
{
double window_frac;
PyObject *x_obj, *y_obj;
/* Parse the input tuple */
if(!PyArg_ParseTuple(args, "OOd", &x_obj, &y_obj, &window_frac))
return NULL;
/* Interpret the input objects as numpy arrays. */
PyObject *x_array = PyArray_FROM_OTF(x_obj, NPY_DOUBLE, NPY_IN_ARRAY);
PyObject *y_array = PyArray_FROM_OTF(y_obj, NPY_DOUBLE, NPY_IN_ARRAY);
/* If that didn't work, throw an exception. */
if(x_array == NULL || y_array == NULL)
{
Py_XDECREF(x_array);
Py_XDECREF(y_array);
return NULL;
}
/* Get the number of data points in each array */
int xsize = (int)PyArray_DIM(x_array, 0);
int ysize = (int)PyArray_DIM(y_array, 0);
/* Get pointers to the data as C-types. */
double *x = (double*)PyArray_DATA(x_array);
double *y = (double*)PyArray_DATA(y_array);
/* Call the external C function to compute the DTW distance. */
int *path = dtw_path(x, y, xsize, ysize, window_frac);
/* Clean up. */
Py_DECREF(x_array);
Py_DECREF(y_array);
/* Copy to numpy array. */
int path_size = path[0];
npy_intp dims[2] = {path_size, 2};
PyObject *ret = PyArray_SimpleNew(2, dims, NPY_INT);
memcpy(PyArray_DATA(ret), &path[1], (path_size)*sizeof(int)*2);
free(path);
return ret;
}