From 2b56e58f6d1616a8c7a3fc9e61a9d3c593534cda Mon Sep 17 00:00:00 2001
From: Honi Sanders <honi@brandeis.edu>
Date: Thu, 11 Jun 2015 16:32:23 -0400
Subject: [PATCH 01/12] ENH: add maxlag to numpy.correlate (see gh-5954)

---
 numpy/core/numeric.py                        | 172 +++++++++++++++----
 numpy/core/src/multiarray/multiarraymodule.c | 160 ++++++++++-------
 numpy/core/tests/test_numeric.py             | 156 +++++++++++++++--
 3 files changed, 384 insertions(+), 104 deletions(-)

diff --git a/numpy/core/numeric.py b/numpy/core/numeric.py
index f29b750f68b0..7bdcd1e1bb38 100644
--- a/numpy/core/numeric.py
+++ b/numpy/core/numeric.py
@@ -8,7 +8,7 @@
 from . import umath
 from .umath import (invert, sin, UFUNC_BUFSIZE_DEFAULT, ERR_IGNORE,
                     ERR_WARN, ERR_RAISE, ERR_CALL, ERR_PRINT, ERR_LOG,
-                    ERR_DEFAULT, PINF, NAN)
+                    ERR_DEFAULT, PINF, NAN, ceil)
 from . import numerictypes
 from .numerictypes import longlong, intc, int_, float_, complex_, bool_
 
@@ -398,7 +398,8 @@ def extend_all(module):
 
 
 def asarray(a, dtype=None, order=None):
-    """Convert the input to an array.
+    """
+    Convert the input to an array.
 
     Parameters
     ----------
@@ -409,9 +410,8 @@ def asarray(a, dtype=None, order=None):
     dtype : data-type, optional
         By default, the data-type is inferred from the input data.
     order : {'C', 'F'}, optional
-        Whether to use row-major (C-style) or
-        column-major (Fortran-style) memory representation.
-        Defaults to 'C'.
+        Whether to use row-major ('C') or column-major ('F' for FORTRAN)
+        memory representation.  Defaults to 'C'.
 
     Returns
     -------
@@ -468,7 +468,8 @@ class ndarray is returned.
     return array(a, dtype, copy=False, order=order)
 
 def asanyarray(a, dtype=None, order=None):
-    """Convert the input to an ndarray, but pass ndarray subclasses through.
+    """
+    Convert the input to an ndarray, but pass ndarray subclasses through.
 
     Parameters
     ----------
@@ -479,8 +480,8 @@ def asanyarray(a, dtype=None, order=None):
     dtype : data-type, optional
         By default, the data-type is inferred from the input data.
     order : {'C', 'F'}, optional
-        Whether to use row-major (C-style) or column-major
-        (Fortran-style) memory representation.  Defaults to 'C'.
+        Whether to use row-major ('C') or column-major ('F') memory
+        representation.  Defaults to 'C'.
 
     Returns
     -------
@@ -830,7 +831,34 @@ def _mode_from_name(mode):
         return _mode_from_name_dict[mode.lower()[0]]
     return mode
 
-def correlate(a, v, mode='valid'):
+def _lags_from_mode(alen, vlen, mode):
+    if type(mode) is int:       # maxlag
+        lags = (-mode+1, mode, 1)
+        mode = 3
+    elif type(mode) is tuple:     # minlag and maxlag
+        if len(mode) > 2:
+            lags = (int(mode[0]), int(mode[1]), int(mode[2]))
+        else:
+            lags = (int(mode[0]), int(mode[1]), 1)
+        mode = 3
+    elif isinstance(mode, basestring):
+        mode = _mode_from_name_dict[mode.lower()[0]]
+        if alen < vlen:
+            alen, vlen = vlen, alen
+            inverted = 1
+        else:
+            inverted = 0
+        if mode is 0:
+            lags = (0, alen-vlen+1, 1)
+        elif mode is 1:
+            lags = (-int(vlen/2), alen - int(vlen/2), 1)
+        elif mode is 2:
+            lags = (-vlen+1, alen, 1)
+        if inverted:
+            lags = (-int(ceil((lags[1]-lags[0])/float(lags[2])))*lags[2]-lags[0]+lags[2], -lags[0]+lags[2], lags[2])
+    return mode, lags
+
+def correlate(a, v, mode='valid', old_behavior=False, lagvector=False):
     """
     Cross-correlation of two 1-dimensional sequences.
 
@@ -846,22 +874,32 @@ def correlate(a, v, mode='valid'):
     ----------
     a, v : array_like
         Input sequences.
-    mode : {'valid', 'same', 'full'}, optional
+    mode : int, int tuple, or {'valid', 'same', 'full'}, optional
         Refer to the `convolve` docstring.  Note that the default
         is `valid`, unlike `convolve`, which uses `full`.
-    old_behavior : bool
-        `old_behavior` was removed in NumPy 1.10. If you need the old
-        behavior, use `multiarray.correlate`.
+    old_behavior : bool, optional
+        If True, uses the old behavior from Numeric,
+        (correlate(a,v) == correlate(v,a), and the conjugate is not taken
+        for complex arrays). If False, uses the conventional signal
+        processing definition.  False is default.
+    lagvector : bool, optional
+        If True, the function returns a lagvector array in addition to the
+        cross-correlation result.  The lagvector contains the indices of
+        the lags for which the cross-correlation was calculated.  It is
+        the same length as the return array, and corresponds one-to-one.
+        False is default.
 
     Returns
     -------
     out : ndarray
         Discrete cross-correlation of `a` and `v`.
+    lagvector : ndarray, optional
+        The indices of the lags for which the cross-correlation was calculated.
+        It is the same length as out, and corresponds one-to-one.
 
     See Also
     --------
     convolve : Discrete, linear convolution of two one-dimensional sequences.
-    multiarray.correlate : Old, no conjugate, version of correlate.
 
     Notes
     -----
@@ -876,10 +914,14 @@ def correlate(a, v, mode='valid'):
     --------
     >>> np.correlate([1, 2, 3], [0, 1, 0.5])
     array([ 3.5])
-    >>> np.correlate([1, 2, 3], [0, 1, 0.5], "same")
+    >>> np.correlate([1, 2, 3], [0, 1, 0.5], mode="same")
     array([ 2. ,  3.5,  3. ])
-    >>> np.correlate([1, 2, 3], [0, 1, 0.5], "full")
-    array([ 0.5,  2. ,  3.5,  3. ,  0. ])
+    >>> np.correlate([1, 2, 3], [0, 1, 0.5], mode="full", lagvector=True)
+    (array([ 0.5,  2. ,  3.5,  3. ,  0. ]), array([-2, -1,  0,  1,  2]))
+    >>> np.correlate([1, 2, 3], [0, 1, 0.5], mode=2)
+    array([ 2. ,  3.5,  3. ])
+    >>> np.correlate([1, 2, 3], [0, 1, 0.5], mode=(-1,2,2), lagvector=True)
+    (array([ 2.,  3.]), array([-1,  1]))
 
     Using complex sequences:
 
@@ -894,10 +936,29 @@ def correlate(a, v, mode='valid'):
     array([ 0.0+0.j ,  3.0+1.j ,  1.5+1.5j,  1.0+0.j ,  0.5+0.5j])
 
     """
-    mode = _mode_from_name(mode)
-    return multiarray.correlate2(a, v, mode)
+    mode, lags = _lags_from_mode(len(a), len(v), mode)
+
+# the old behavior should be made available under a different name, see thread
+# http://thread.gmane.org/gmane.comp.python.numeric.general/12609/focus=12630
+    if old_behavior:
+        warnings.warn("""
+The old behavior of correlate was deprecated for 1.4.0, and will be completely removed
+for NumPy 2.0.
+
+The new behavior fits the conventional definition of correlation: inputs are
+never swapped, and the second argument is conjugated for complex arrays.""",
+            DeprecationWarning)
+        if lagvector:
+            return multiarray.correlate(a, v, mode, lags[0], lags[1], lags[2]), arange(lags[0], lags[1], lags[2])
+        else:
+            return multiarray.correlate(a, v, mode, lags[0], lags[1], lags[2])
+    else:
+        if lagvector:
+            return multiarray.correlate2(a, v, mode, lags[0], lags[1], lags[2]), arange(lags[0], lags[1], lags[2])
+        else:
+            return multiarray.correlate2(a, v, mode, lags[0], lags[1], lags[2])
 
-def convolve(a,v,mode='full'):
+def convolve(a, v, mode='full', lagvector=False):
     """
     Returns the discrete, linear convolution of two one-dimensional sequences.
 
@@ -915,27 +976,53 @@ def convolve(a,v,mode='full'):
         First one-dimensional input array.
     v : (M,) array_like
         Second one-dimensional input array.
-    mode : {'full', 'valid', 'same'}, optional
+    mode : int, int tuple, or {'full', 'valid', 'same'}, optional
+        int (maxlag):
+          This calculates the convolution for all lags starting at
+          (-maxlag + 1) and ending at (maxlag - 1), with steps of size 1.
+          See the optional lagvec argument to get an array containing
+          lags corresponding to the convolution values in the return array.
+
+        tuple (minlag, maxlag) or (minlag, maxlag, lagstep):
+          This calculates the convolution for all lags starting at
+          minlag and ending at (maxlag - 1), with steps of size lagstep.
+          The lags for which the convolution will be calculated correspond
+          with the values in the vector formed by numpy.arange() with the
+          same tuple argument.
+
         'full':
           By default, mode is 'full'.  This returns the convolution
           at each point of overlap, with an output shape of (N+M-1,). At
           the end-points of the convolution, the signals do not overlap
-          completely, and boundary effects may be seen.
+          completely, and boundary effects may be seen.  This corresponds
+          with a lag tuple of (-M+1, N-1, 1) for N>M or (-N+1, M+1, 1)
+          for M>N.
 
         'same':
           Mode `same` returns output of length ``max(M, N)``.  Boundary
-          effects are still visible.
+          effects are still visible. This corresponds with a lag tuple of
+          (-M/2, N-M/2, 1) for N>M or (-M+N/2+1, N/2+1, 1) for M>N.
 
         'valid':
           Mode `valid` returns output of length
           ``max(M, N) - min(M, N) + 1``.  The convolution product is only given
           for points where the signals overlap completely.  Values outside
-          the signal boundary have no effect.
+          the signal boundary have no effect. This corresponds with a lag tuple
+          of (0, N-M+1, 1) for N>M or (-M+N, 1, 1) for M>N.
+    lagvector : bool, optional
+        If True, the function returns a lagvector array in addition to the
+        convolution result.  The lagvector contains the indices of
+        the lags for which the convolution was calculated.  It is
+        the same length as the return array, and corresponds one-to-one.
+        False is default.
 
     Returns
     -------
     out : ndarray
         Discrete, linear convolution of `a` and `v`.
+    lagvector : ndarray, optional
+        The indices of the lags for which the convolution was calculated.
+        It is the same length as out, and corresponds one-to-one.
 
     See Also
     --------
@@ -974,14 +1061,34 @@ def convolve(a,v,mode='full'):
     Contains boundary effects, where zeros are taken
     into account:
 
-    >>> np.convolve([1,2,3],[0,1,0.5], 'same')
+    >>> np.convolve([1,2,3],[0,1,0.5], mode='same')
     array([ 1. ,  2.5,  4. ])
 
     The two arrays are of the same length, so there
-    is only one position where they completely overlap:
+    is only one position where they completely overlap,
+    corresponding to a lag of 0.  lagvector=True causes
+    the function to return the lagvector corresponding
+    to the convolution in addition to the convolution
+    itself:
+
+    >>> np.convolve([1,2,3],[0,1,0.5], mode='valid', lagvector=True)
+    (array([ 2.5]), array([0]))
 
-    >>> np.convolve([1,2,3],[0,1,0.5], 'valid')
-    array([ 2.5])
+    Find the convolution for lags ranging from -1 to 1
+    (0 is the lag for which the left sides of the arrays
+    are aligned, -1 has the second vector to the left of
+    the first, and +1 has the second vector to the right
+    of the first):
+
+    >>> np.convolve([1,2,3],[0,1,0.5], mode=2, lagvector=True)
+    (array([ 1. ,  2.5,  4. ]), array([-1,  0,  1]))
+
+    Find the convolution for lags ranging from -2 to 4
+    with steps of length 2 (the maxlag member of the
+    lag range tuple is non-inclusive, similar to np.arange()):
+
+    >>> np.convolve([1,2,3,4,5],[0,1,0.5], mode=(-2,6,2), lagvector=True)
+    (array([ 0. ,  2.5,  5.5,  2.5]), array([-2,  0,  2,  4]))
 
     """
     a, v = array(a, copy=False, ndmin=1), array(v, copy=False, ndmin=1)
@@ -991,8 +1098,11 @@ def convolve(a,v,mode='full'):
         raise ValueError('a cannot be empty')
     if len(v) == 0 :
         raise ValueError('v cannot be empty')
-    mode = _mode_from_name(mode)
-    return multiarray.correlate(a, v[::-1], mode)
+    mode, lags = _lags_from_mode(len(a), len(v), mode)
+    if lagvector:
+        return multiarray.correlate2(a, v[::-1], mode, lags[0], lags[1], lags[2]), arange(lags[0], lags[1], lags[2])
+    else:
+        return multiarray.correlate2(a, v[::-1], mode, lags[0], lags[1], lags[2])
 
 def outer(a, b, out=None):
     """
@@ -1100,7 +1210,6 @@ def alterdot():
     restoredot : `restoredot` undoes the effects of `alterdot`.
 
     """
-    # 2014-08-13, 1.10
     warnings.warn("alterdot no longer does anything.", DeprecationWarning)
 
 
@@ -1124,7 +1233,6 @@ def restoredot():
     alterdot : `restoredot` undoes the effects of `alterdot`.
 
     """
-    # 2014-08-13, 1.10
     warnings.warn("restoredot no longer does anything.", DeprecationWarning)
 
 
diff --git a/numpy/core/src/multiarray/multiarraymodule.c b/numpy/core/src/multiarray/multiarraymodule.c
index 73265c3b607f..3d335a1fe4b3 100644
--- a/numpy/core/src/multiarray/multiarraymodule.c
+++ b/numpy/core/src/multiarray/multiarraymodule.c
@@ -1142,19 +1142,24 @@ PyArray_CopyAndTranspose(PyObject *op)
  */
 static PyArrayObject*
 _pyarray_correlate(PyArrayObject *ap1, PyArrayObject *ap2, int typenum,
-                   int mode, int *inverted)
+                   int mode, int *inverted, int minlag, int maxlag, int lagstep)
 {
     PyArrayObject *ret;
     npy_intp length;
-    npy_intp i, n1, n2, n, n_left, n_right;
+    npy_intp i, n1, n2, n, n11;
+    npy_intp lag, tmplag, maxleft, maxright;
     npy_intp is1, is2, os;
     char *ip1, *ip2, *op;
     PyArray_DotFunc *dot;
 
     NPY_BEGIN_THREADS_DEF;
 
-    n1 = PyArray_DIMS(ap1)[0];
-    n2 = PyArray_DIMS(ap2)[0];
+    if (lagstep == NPY_MAX_INTP || lagstep == 0) {
+        lagstep = 1;
+    }
+
+    n1 = PyArray_DIMS(ap1)[0];      /* size of x */
+    n2 = PyArray_DIMS(ap2)[0];      /* size of y */
     if (n1 < n2) {
         ret = ap1;
         ap1 = ap2;
@@ -1163,40 +1168,67 @@ _pyarray_correlate(PyArrayObject *ap1, PyArrayObject *ap2, int typenum,
         i = n1;
         n1 = n2;
         n2 = i;
-        *inverted = 1;
-    } else {
-        *inverted = 0;
+        minlag = -minlag;
+        maxlag = -maxlag;
+        lagstep = -lagstep;
     }
 
-    length = n1;
-    n = n2;
+    if (maxlag == NPY_MAX_INTP || minlag == NPY_MAX_INTP) {
     switch(mode) {
-    case 0:
-        length = length - n + 1;
-        n_left = n_right = 0;
+    lagstep = 1;
+    case 0:         /* mode = 'valid' */
+        minlag = 0;
+        maxlag = n1 - n2 + 1;
         break;
-    case 1:
-        n_left = (npy_intp)(n/2);
-        n_right = n - n_left - 1;
+    case 1:         /* mode = 'same' */
+        minlag = -(npy_intp)(n2/2);
+        maxlag = n1 + minlag;
         break;
-    case 2:
-        n_right = n - 1;
-        n_left = n - 1;
-        length = length + n - 1;
+    case 2:         /* mode = 'full' */
+        minlag = -n2 + 1;
+        maxlag = n1;
+        break;
+    case 3:         /* mode = 'maxlag' */
+        if (maxlag == NPY_MAX_INTP) {
+            minlag = -n2 + 1;   /* nothing good to do with situation */
+            maxlag = n1;
+        }
+        else if (minlag == NPY_MAX_INTP) {
+            minlag = (maxlag > n2 ? -n2+1 : -maxlag+1);
+        }
         break;
     default:
-        PyErr_SetString(PyExc_ValueError, "mode must be 0, 1, or 2");
+        PyErr_SetString(PyExc_ValueError, "mode must be 0, 1, 2, or 3");
         return NULL;
     }
+    }
+
+    if (lagstep < 0) {
+        *inverted = 1;
+        i = minlag;
+        minlag = (int)(npy_ceil((maxlag - minlag)/(float)lagstep))*lagstep + minlag - lagstep;
+        maxlag = i - lagstep;
+        lagstep = -lagstep;
+    }
+    else {
+        *inverted = 0;
+    }
+    if (maxlag <= minlag) {
+        length = 0;
+    }
+    else {
+        length = (maxlag - minlag + lagstep - 1)/lagstep;
+    }
 
     /*
      * Need to choose an output array that can hold a sum
      * -- use priority to determine which subtype.
      */
-    ret = new_array_for_sum(ap1, ap2, NULL, 1, &length, typenum);
+    ret = new_array_for_sum(ap1, ap2, NULL, 1, &length, typenum);   /* answer */
     if (ret == NULL) {
         return NULL;
     }
+    PyArray_FILLWBYTE(ret, 0);
     dot = PyArray_DESCR(ret)->f->dotfunc;
     if (dot == NULL) {
         PyErr_SetString(PyExc_ValueError,
@@ -1205,36 +1237,44 @@ _pyarray_correlate(PyArrayObject *ap1, PyArrayObject *ap2, int typenum,
     }
 
     NPY_BEGIN_THREADS_DESCR(PyArray_DESCR(ret));
-    is1 = PyArray_STRIDES(ap1)[0];
-    is2 = PyArray_STRIDES(ap2)[0];
-    op = PyArray_DATA(ret);
-    os = PyArray_DESCR(ret)->elsize;
-    ip1 = PyArray_DATA(ap1);
-    ip2 = PyArray_BYTES(ap2) + n_left*is2;
-    n = n - n_left;
-    for (i = 0; i < n_left; i++) {
-        dot(ip1, is1, ip2, is2, op, n, ret);
-        n++;
-        ip2 -= is2;
-        op += os;
+    ip1 = PyArray_DATA(ap1);            /* x[0]             */
+    is1 = PyArray_STRIDES(ap1)[0];      /* x strides        */
+    ip2 = PyArray_DATA(ap2);            /* y[0]             */
+    is2 = PyArray_STRIDES(ap2)[0];      /* y strides        */
+    op = PyArray_DATA(ret);             /* answer data      */
+    os = PyArray_DESCR(ret)->elsize;    /* answer strides   */
+
+    lag = minlag;
+    maxleft = (0 < maxlag ? 0 : maxlag);
+    for (lag = minlag; lag < maxleft; lag+=lagstep) {      /* for lags where y is left of x */
+        n = n2 + lag;   /* overlap is length of y - lag */
+        dot(ip1, is1, ip2 - lag*is2, is2, op, n, ret);
+        op += os;       /* iterate over answer vector   */
+    }
+    if (maxlag < n1 - n2) { /* if maxlag doesn't take y all the way to the end of x */
+        n11 = maxlag + n2;      /* relevant length of x is smaller                  */
     }
-    if (small_correlate(ip1, is1, n1 - n2 + 1, PyArray_TYPE(ap1),
-                        ip2, is2, n, PyArray_TYPE(ap2),
-                        op, os)) {
-        ip1 += is1 * (n1 - n2 + 1);
-        op += os * (n1 - n2 + 1);
+    else {
+        n11 = n1;
+    }       /* starts at lag=minlag if minlag>0. Does lags where y entirely overlaps with x. */
+    if (lagstep == 1 && small_correlate(ip1 + lag*is1, is1, n11 - n2 + 1 - lag, PyArray_TYPE(ap1),
+                                        ip2, is2, n2, PyArray_TYPE(ap2),
+                                        op, os)) {
+        lag = n11 - n2 + 1;
+        op += os * (n11 - n2 + 1);
     }
     else {
-        for (i = 0; i < (n1 - n2 + 1); i++) {
-            dot(ip1, is1, ip2, is2, op, n, ret);
-            ip1 += is1;
+        tmplag = lag;
+        for (lag = tmplag; lag < (n11 - n2 + 1); lag+=lagstep) {
+            dot(ip1 + lag*is1, is1, ip2, is2, op, n2, ret);
             op += os;
         }
     }
-    for (i = 0; i < n_right; i++) {
-        n--;
-        dot(ip1, is1, ip2, is2, op, n, ret);
-        ip1 += is1;
+    maxright = (maxlag < n1  ? maxlag : n1 );
+    tmplag = lag;
+    for (lag = tmplag; lag < maxright; lag+=lagstep) {     /* for lags where y is right of x */
+        n = n1 - lag;
+        dot(ip1 + lag*is1, is1, ip2, is2, op, n, ret);
         op += os;
     }
 
@@ -1309,7 +1349,7 @@ _pyarray_revert(PyArrayObject *ret)
  * correlate(a2, a1), and conjugate the second argument for complex inputs
  */
 NPY_NO_EXPORT PyObject *
-PyArray_Correlate2(PyObject *op1, PyObject *op2, int mode)
+PyArray_Correlate2(PyObject *op1, PyObject *op2, int mode, int minlag, int maxlag, int lagstep)
 {
     PyArrayObject *ap1, *ap2, *ret = NULL;
     int typenum;
@@ -1344,7 +1384,7 @@ PyArray_Correlate2(PyObject *op1, PyObject *op2, int mode)
         ap2 = cap2;
     }
 
-    ret = _pyarray_correlate(ap1, ap2, typenum, mode, &inverted);
+    ret = _pyarray_correlate(ap1, ap2, typenum, mode, &inverted, minlag, maxlag, lagstep);
     if (ret == NULL) {
         goto clean_ap2;
     }
@@ -1377,7 +1417,7 @@ PyArray_Correlate2(PyObject *op1, PyObject *op2, int mode)
  * Numeric.correlate(a1,a2,mode)
  */
 NPY_NO_EXPORT PyObject *
-PyArray_Correlate(PyObject *op1, PyObject *op2, int mode)
+PyArray_Correlate(PyObject *op1, PyObject *op2, int mode, int minlag, int maxlag, int lagstep)
 {
     PyArrayObject *ap1, *ap2, *ret = NULL;
     int typenum;
@@ -1401,7 +1441,7 @@ PyArray_Correlate(PyObject *op1, PyObject *op2, int mode)
         goto fail;
     }
 
-    ret = _pyarray_correlate(ap1, ap2, typenum, mode, &unused);
+    ret = _pyarray_correlate(ap1, ap2, typenum, mode, &unused, minlag, maxlag, lagstep);
     if(ret == NULL) {
         goto fail;
     }
@@ -2896,28 +2936,28 @@ static PyObject *
 array_correlate(PyObject *NPY_UNUSED(dummy), PyObject *args, PyObject *kwds)
 {
     PyObject *shape, *a0;
-    int mode = 0;
-    static char *kwlist[] = {"a", "v", "mode", NULL};
+    int mode = 0, maxlag = NPY_MAX_INTP, minlag = NPY_MAX_INTP, lagstep = NPY_MAX_INTP;
+    static char *kwlist[] = {"a", "v", "mode", "minlag", "maxlag", "lagstep", NULL};
 
-    if (!PyArg_ParseTupleAndKeywords(args, kwds, "OO|i", kwlist,
-                &a0, &shape, &mode)) {
+    if (!PyArg_ParseTupleAndKeywords(args, kwds, "OO|iiii", kwlist,
+                &a0, &shape, &mode, &minlag, &maxlag, &lagstep)) {
         return NULL;
     }
-    return PyArray_Correlate(a0, shape, mode);
+    return PyArray_Correlate(a0, shape, mode, minlag, maxlag, lagstep);
 }
 
-static PyObject*
+static PyObject *
 array_correlate2(PyObject *NPY_UNUSED(dummy), PyObject *args, PyObject *kwds)
 {
     PyObject *shape, *a0;
-    int mode = 0;
-    static char *kwlist[] = {"a", "v", "mode", NULL};
+    int mode = 0, maxlag = NPY_MAX_INTP, minlag = NPY_MAX_INTP, lagstep = NPY_MAX_INTP;
+    static char *kwlist[] = {"a", "v", "mode", "minlag", "maxlag", "lagstep", NULL};
 
-    if (!PyArg_ParseTupleAndKeywords(args, kwds, "OO|i", kwlist,
-                &a0, &shape, &mode)) {
+    if (!PyArg_ParseTupleAndKeywords(args, kwds, "OO|iiii", kwlist,
+                &a0, &shape, &mode, &minlag, &maxlag, &lagstep)) {
         return NULL;
     }
-    return PyArray_Correlate2(a0, shape, mode);
+    return PyArray_Correlate2(a0, shape, mode, minlag, maxlag, lagstep);
 }
 
 static PyObject *
diff --git a/numpy/core/tests/test_numeric.py b/numpy/core/tests/test_numeric.py
index c39900ffa7bd..4adf12296d75 100644
--- a/numpy/core/tests/test_numeric.py
+++ b/numpy/core/tests/test_numeric.py
@@ -980,7 +980,6 @@ def test_negative(self):
         assert_equal(binary_repr(-1), '-1')
         assert_equal(binary_repr(-1, width=8), '11111111')
 
-
 class TestBaseRepr(TestCase):
     def test_base3(self):
         assert_equal(base_repr(3**5, 3), '100000')
@@ -1947,39 +1946,140 @@ def test_filled_like(self):
         self.check_like_function(np.full_like, 123.456, True)
         self.check_like_function(np.full_like, np.inf, True)
 
-class TestCorrelate(TestCase):
+class _TestCorrelate(TestCase):
     def _setup(self, dt):
         self.x = np.array([1, 2, 3, 4, 5], dtype=dt)
         self.xs = np.arange(1, 20)[::3]
         self.y = np.array([-1, -2, -3], dtype=dt)
         self.z1 = np.array([ -3.,  -8., -14., -20., -26., -14.,  -5.], dtype=dt)
+        self.z1s = np.array([-8., -14., -20., -26., -14.], dtype=dt)
+        self.z1v = np.array([-14., -20., -26.], dtype=dt)
         self.z1_4 = np.array([-2., -5., -8., -11., -14., -5.], dtype=dt)
         self.z1r = np.array([-15., -22., -22., -16., -10.,  -4.,  -1.], dtype=dt)
         self.z2 = np.array([-5., -14., -26., -20., -14., -8.,  -3.], dtype=dt)
+        self.z2s = np.array([-14, -26., -20., -14., -8.], dtype=dt)
+        self.z2v = np.array([-26., -20., -14.], dtype=dt)
         self.z2r = np.array([-1., -4., -10., -16., -22., -22., -15.], dtype=dt)
         self.zs = np.array([-3., -14., -30., -48., -66., -84.,
                            -102., -54., -19.], dtype=dt)
 
     def test_float(self):
         self._setup(np.float)
-        z = np.correlate(self.x, self.y, 'full')
+        lagvec = np.arange(1)
+        z, lagvec = np.correlate(self.x, self.y, 'full', old_behavior=self.old_behavior,lagvector=True)
         assert_array_almost_equal(z, self.z1)
-        z = np.correlate(self.x, self.y[:-1], 'full')
+        assert_equal(lagvec.size, z.size)
+        assert_array_equal(lagvec, np.arange(-self.y.size + 1, self.x.size))
+        z, lagvec = np.correlate(self.x, self.y, 'same', old_behavior=self.old_behavior, lagvector=True)
+        assert_array_almost_equal(z, self.z1s)
+        assert_array_almost_equal(z, self.z1[lagvec+self.y.size-1])
+        assert_equal(lagvec.size, z.size)
+        if self.x.size > self.y.size:
+            assert_array_equal(lagvec, np.arange(-int(self.y.size/2), self.x.size - int(self.y.size/2)))
+        else:
+            assert_array_equal(lagvec, np.arange(-self.y.size + int(self.x.size/2) + 1, int(self.x.size/2) + 1))
+        z, lagvec = np.correlate(self.x, self.y, 'valid', old_behavior=self.old_behavior, lagvector=True)
+        assert_array_almost_equal(z, self.z1v)
+        assert_array_almost_equal(z, self.z1[lagvec+self.y.size-1])
+        assert_equal(lagvec.size, z.size)
+        if self.x.size > self.y.size:
+            assert_array_equal(lagvec, np.arange(0, self.x.size - self.y.size + 1))
+        else:
+            assert_array_equal(lagvec, np.arange(self.x.size - self.y.size, 1))
+
+        z = np.correlate(self.x, self.y[:-1], 'full', old_behavior=self.old_behavior)
         assert_array_almost_equal(z, self.z1_4)
-        z = np.correlate(self.y, self.x, 'full')
+        z, lagvec = np.correlate(self.y, self.x, 'full', old_behavior=self.old_behavior, lagvector=True)
         assert_array_almost_equal(z, self.z2)
-        z = np.correlate(self.x[::-1], self.y, 'full')
+        assert_equal(lagvec.size, z.size)
+        assert_array_equal(lagvec, np.arange(-self.x.size + 1, self.y.size))
+        z, lagvec = np.correlate(self.y, self.x, 'same', old_behavior=self.old_behavior, lagvector=True)
+        assert_array_almost_equal(z, self.z2s)
+        assert_array_almost_equal(z, self.z2[lagvec+self.x.size-1])
+        assert_equal(lagvec.size, z.size)
+        if self.y.size > self.x.size:
+            assert_array_equal(lagvec, np.arange(-int(self.x.size/2), self.y.size - int(self.x.size/2)))
+        else:
+            assert_array_equal(lagvec, np.arange(-self.x.size + int(self.y.size/2) + 1, int(self.y.size/2) + 1))
+        z, lagvec = np.correlate(self.y, self.x, 'valid', old_behavior=self.old_behavior, lagvector=True)
+        assert_array_almost_equal(z, self.z2v)
+        assert_array_almost_equal(z, self.z2[lagvec+self.x.size-1])
+        assert_equal(lagvec.size, z.size)
+        if self.y.size > self.x.size:
+            assert_array_equal(lagvec, np.arange(0, self.y.size - self.x.size + 1))
+        else:
+            assert_array_equal(lagvec, np.arange(self.y.size - self.x.size, 1))
+        z = np.correlate(self.x[::-1], self.y, 'full', old_behavior=self.old_behavior)
         assert_array_almost_equal(z, self.z1r)
-        z = np.correlate(self.y, self.x[::-1], 'full')
+        z = np.correlate(self.y, self.x[::-1], 'full', old_behavior=self.old_behavior)
         assert_array_almost_equal(z, self.z2r)
-        z = np.correlate(self.xs, self.y, 'full')
+        z = np.correlate(self.xs, self.y, 'full', old_behavior=self.old_behavior)
         assert_array_almost_equal(z, self.zs)
 
+    def test_lags(self):
+        self._setup(np.float)
+        longlen = max(self.x.size, self.y.size)
+        lagvec = np.arange(1)
+        maxlag = 3 + longlen
+        minlag = -2
+        lagstep = 2
+        tmp = np.correlate(self.x, self.y, 'full')
+        z_full = np.zeros(tmp.size + 100)
+        z_full[:tmp.size] = tmp
+        z, lagvec = np.correlate(self.x, self.y, maxlag, lagvector=True)
+        assert_array_equal(z[0:3], np.zeros(3))
+        assert_array_equal(z[-3:], np.zeros(3))
+        assert_equal(lagvec.size, z.size)
+        assert_array_equal(lagvec, np.arange(-maxlag+1, maxlag))
+        z, lagvec = np.correlate(self.x, self.y, (minlag, maxlag), lagvector=True)
+        assert_array_almost_equal(z, z_full[lagvec+self.y.size-1])
+        assert_equal(lagvec.size, z.size)
+        assert_array_equal(lagvec, np.arange(minlag, maxlag))
+        z, lagvec = np.correlate(self.x, self.y, (minlag, maxlag, lagstep), lagvector=True)
+        assert_array_almost_equal(z, z_full[lagvec+self.y.size-1])
+        assert_equal(lagvec.size, z.size)
+        assert_array_equal(lagvec, np.arange(minlag, maxlag, lagstep))
+        lagstep = 3
+        z, lagvec = np.correlate(self.x, self.y, (minlag, maxlag, lagstep), lagvector=True)
+        assert_array_almost_equal(z, z_full[lagvec+self.y.size-1])
+        assert_equal(lagvec.size, z.size)
+        minlag = 10
+        maxlag = -2
+        lagstep = -2
+        z, lagvec = np.correlate(self.x, self.y, (minlag, maxlag, lagstep), lagvector=True)
+        assert_array_almost_equal(z, z_full[lagvec+self.y.size-1])
+        assert_equal(lagvec.size, z.size)
+        assert_array_equal(lagvec, np.arange(int(minlag), int(maxlag), lagstep))
+        maxlag = 10.2
+        minlag = -2.2
+        z, lagvec = np.correlate(self.x, self.y, (minlag, maxlag), lagvector=True)
+        assert_array_almost_equal(z, z_full[lagvec+self.y.size-1])
+        assert_equal(lagvec.size, z.size)
+        assert_array_equal(lagvec, np.arange(int(minlag), int(maxlag)))
+        z, lagvec = np.correlate(self.x, self.y, (maxlag, minlag), lagvector=True)
+        assert_array_almost_equal(z, z_full[lagvec+self.y.size-1])
+        assert_equal(lagvec.size, z.size)
+        assert_array_equal(lagvec, np.arange(int(maxlag), int(minlag)))
+        maxlag = -1
+        minlag = -2
+        z, lagvec = np.correlate(self.x, self.y, (minlag, maxlag), lagvector=True)
+        assert_array_almost_equal(z, z_full[lagvec+self.y.size-1])
+        assert_equal(lagvec.size, z.size)
+        assert_array_equal(lagvec, np.arange(int(minlag), int(maxlag)))
+        maxlag = 2
+        minlag = 4
+        z, lagvec = np.correlate(self.x, self.y, (minlag, maxlag), lagvector=True)
+        assert_array_almost_equal(z, z_full[lagvec+self.y.size-1])
+        assert_equal(lagvec.size, z.size)
+        assert_array_equal(lagvec, np.arange(int(minlag), int(maxlag)))
+        # check non-integer lagsteps?
+        # i want to be able to switch x and y and test
+
     def test_object(self):
         self._setup(Decimal)
-        z = np.correlate(self.x, self.y, 'full')
+        z = np.correlate(self.x, self.y, 'full', old_behavior=self.old_behavior)
         assert_array_almost_equal(z, self.z1)
-        z = np.correlate(self.y, self.x, 'full')
+        z = np.correlate(self.y, self.x, 'full', old_behavior=self.old_behavior)
         assert_array_almost_equal(z, self.z2)
 
     def test_no_overwrite(self):
@@ -1989,15 +2089,47 @@ def test_no_overwrite(self):
         assert_array_equal(d, np.ones(100))
         assert_array_equal(k, np.ones(3))
 
+class TestCorrelate(_TestCorrelate):
+    old_behavior = True
+    def _setup(self, dt):
+        # correlate uses an unconventional definition so that correlate(a, b)
+        # == correlate(b, a), so force the corresponding outputs to be the same
+        # as well
+        _TestCorrelate._setup(self, dt)
+        self.z2 = self.z1
+        self.z2r = self.z1r
+        self.z2s = self.z1s
+        self.z2v = self.z1v
+
+    @dec.deprecated()
+    def test_complex(self):
+        x = np.array([1, 2, 3, 4+1j], dtype=np.complex)
+        y = np.array([-1, -2j, 3+1j], dtype=np.complex)
+        r_z = np.array([3+1j, 6, 8-1j, 9+1j, -1-8j, -4-1j], dtype=np.complex)
+        z = np.correlate(x, y, 'full', old_behavior=self.old_behavior)
+        assert_array_almost_equal(z, r_z)
+
+    @dec.deprecated()
+    def test_float(self):
+        _TestCorrelate.test_float(self)
+
+    @dec.deprecated()
+    def test_object(self):
+        _TestCorrelate.test_object(self)
+
+class TestCorrelateNew(_TestCorrelate):
+    old_behavior = False
     def test_complex(self):
         x = np.array([1, 2, 3, 4+1j], dtype=np.complex)
         y = np.array([-1, -2j, 3+1j], dtype=np.complex)
         r_z = np.array([3-1j, 6, 8+1j, 11+5j, -5+8j, -4-1j], dtype=np.complex)
+        #z = np.acorrelate(x, y, 'full')
+        #assert_array_almost_equal(z, r_z)
+
         r_z = r_z[::-1].conjugate()
-        z = correlate(y, x, mode='full')
+        z = np.correlate(y, x, 'full', old_behavior=self.old_behavior)
         assert_array_almost_equal(z, r_z)
 
-
 class TestConvolve(TestCase):
     def test_object(self):
         d = [1.] * 100

From db1f0772acfcc5551ba6bc74c59ddb3a8382cce5 Mon Sep 17 00:00:00 2001
From: Honi Sanders <honi@brandeis.edu>
Date: Fri, 19 Jun 2015 17:44:07 -0400
Subject: [PATCH 02/12] MAINT: Improvements of numpy.correlate() (See #5954)

replaced int with npy_intp for lag parameters in
_pyarray_correlate and related functions
renamed lagvector input to correlate() and convolve()
to be named return_lags
---
 numpy/core/numeric.py                        | 24 ++++++++--------
 numpy/core/src/multiarray/multiarraymodule.c | 19 +++++++------
 numpy/core/tests/test_numeric.py             | 30 ++++++++++----------
 3 files changed, 38 insertions(+), 35 deletions(-)

diff --git a/numpy/core/numeric.py b/numpy/core/numeric.py
index 7bdcd1e1bb38..84ed98844e14 100644
--- a/numpy/core/numeric.py
+++ b/numpy/core/numeric.py
@@ -858,7 +858,7 @@ def _lags_from_mode(alen, vlen, mode):
             lags = (-int(ceil((lags[1]-lags[0])/float(lags[2])))*lags[2]-lags[0]+lags[2], -lags[0]+lags[2], lags[2])
     return mode, lags
 
-def correlate(a, v, mode='valid', old_behavior=False, lagvector=False):
+def correlate(a, v, mode='valid', old_behavior=False, returns_lags=False):
     """
     Cross-correlation of two 1-dimensional sequences.
 
@@ -882,7 +882,7 @@ def correlate(a, v, mode='valid', old_behavior=False, lagvector=False):
         (correlate(a,v) == correlate(v,a), and the conjugate is not taken
         for complex arrays). If False, uses the conventional signal
         processing definition.  False is default.
-    lagvector : bool, optional
+    returns_lags : bool, optional
         If True, the function returns a lagvector array in addition to the
         cross-correlation result.  The lagvector contains the indices of
         the lags for which the cross-correlation was calculated.  It is
@@ -916,11 +916,11 @@ def correlate(a, v, mode='valid', old_behavior=False, lagvector=False):
     array([ 3.5])
     >>> np.correlate([1, 2, 3], [0, 1, 0.5], mode="same")
     array([ 2. ,  3.5,  3. ])
-    >>> np.correlate([1, 2, 3], [0, 1, 0.5], mode="full", lagvector=True)
+    >>> np.correlate([1, 2, 3], [0, 1, 0.5], mode="full", returns_lags=True)
     (array([ 0.5,  2. ,  3.5,  3. ,  0. ]), array([-2, -1,  0,  1,  2]))
     >>> np.correlate([1, 2, 3], [0, 1, 0.5], mode=2)
     array([ 2. ,  3.5,  3. ])
-    >>> np.correlate([1, 2, 3], [0, 1, 0.5], mode=(-1,2,2), lagvector=True)
+    >>> np.correlate([1, 2, 3], [0, 1, 0.5], mode=(-1,2,2), returns_lags=True)
     (array([ 2.,  3.]), array([-1,  1]))
 
     Using complex sequences:
@@ -948,17 +948,17 @@ def correlate(a, v, mode='valid', old_behavior=False, lagvector=False):
 The new behavior fits the conventional definition of correlation: inputs are
 never swapped, and the second argument is conjugated for complex arrays.""",
             DeprecationWarning)
-        if lagvector:
+        if returns_lags:
             return multiarray.correlate(a, v, mode, lags[0], lags[1], lags[2]), arange(lags[0], lags[1], lags[2])
         else:
             return multiarray.correlate(a, v, mode, lags[0], lags[1], lags[2])
     else:
-        if lagvector:
+        if returns_lags:
             return multiarray.correlate2(a, v, mode, lags[0], lags[1], lags[2]), arange(lags[0], lags[1], lags[2])
         else:
             return multiarray.correlate2(a, v, mode, lags[0], lags[1], lags[2])
 
-def convolve(a, v, mode='full', lagvector=False):
+def convolve(a, v, mode='full', returns_lags=False):
     """
     Returns the discrete, linear convolution of two one-dimensional sequences.
 
@@ -1009,7 +1009,7 @@ def convolve(a, v, mode='full', lagvector=False):
           for points where the signals overlap completely.  Values outside
           the signal boundary have no effect. This corresponds with a lag tuple
           of (0, N-M+1, 1) for N>M or (-M+N, 1, 1) for M>N.
-    lagvector : bool, optional
+    returns_lags : bool, optional
         If True, the function returns a lagvector array in addition to the
         convolution result.  The lagvector contains the indices of
         the lags for which the convolution was calculated.  It is
@@ -1071,7 +1071,7 @@ def convolve(a, v, mode='full', lagvector=False):
     to the convolution in addition to the convolution
     itself:
 
-    >>> np.convolve([1,2,3],[0,1,0.5], mode='valid', lagvector=True)
+    >>> np.convolve([1,2,3],[0,1,0.5], mode='valid', returns_lags=True)
     (array([ 2.5]), array([0]))
 
     Find the convolution for lags ranging from -1 to 1
@@ -1080,14 +1080,14 @@ def convolve(a, v, mode='full', lagvector=False):
     the first, and +1 has the second vector to the right
     of the first):
 
-    >>> np.convolve([1,2,3],[0,1,0.5], mode=2, lagvector=True)
+    >>> np.convolve([1,2,3],[0,1,0.5], mode=2, returns_lags=True)
     (array([ 1. ,  2.5,  4. ]), array([-1,  0,  1]))
 
     Find the convolution for lags ranging from -2 to 4
     with steps of length 2 (the maxlag member of the
     lag range tuple is non-inclusive, similar to np.arange()):
 
-    >>> np.convolve([1,2,3,4,5],[0,1,0.5], mode=(-2,6,2), lagvector=True)
+    >>> np.convolve([1,2,3,4,5],[0,1,0.5], mode=(-2,6,2), returns_lags=True)
     (array([ 0. ,  2.5,  5.5,  2.5]), array([-2,  0,  2,  4]))
 
     """
@@ -1099,7 +1099,7 @@ def convolve(a, v, mode='full', lagvector=False):
     if len(v) == 0 :
         raise ValueError('v cannot be empty')
     mode, lags = _lags_from_mode(len(a), len(v), mode)
-    if lagvector:
+    if returns_lags:
         return multiarray.correlate2(a, v[::-1], mode, lags[0], lags[1], lags[2]), arange(lags[0], lags[1], lags[2])
     else:
         return multiarray.correlate2(a, v[::-1], mode, lags[0], lags[1], lags[2])
diff --git a/numpy/core/src/multiarray/multiarraymodule.c b/numpy/core/src/multiarray/multiarraymodule.c
index 3d335a1fe4b3..9b08c46f7ed3 100644
--- a/numpy/core/src/multiarray/multiarraymodule.c
+++ b/numpy/core/src/multiarray/multiarraymodule.c
@@ -1142,7 +1142,7 @@ PyArray_CopyAndTranspose(PyObject *op)
  */
 static PyArrayObject*
 _pyarray_correlate(PyArrayObject *ap1, PyArrayObject *ap2, int typenum,
-                   int mode, int *inverted, int minlag, int maxlag, int lagstep)
+                   int mode, int *inverted, npy_intp minlag, npy_intp maxlag, npy_intp lagstep)
 {
     PyArrayObject *ret;
     npy_intp length;
@@ -1206,7 +1206,7 @@ _pyarray_correlate(PyArrayObject *ap1, PyArrayObject *ap2, int typenum,
     if (lagstep < 0) {
         *inverted = 1;
         i = minlag;
-        minlag = (int)(npy_ceil((maxlag - minlag)/(float)lagstep))*lagstep + minlag - lagstep;
+        minlag = (npy_intp)(npy_ceil((maxlag - minlag)/(float)lagstep))*lagstep + minlag - lagstep;
         maxlag = i - lagstep;
         lagstep = -lagstep;
     }
@@ -1349,7 +1349,7 @@ _pyarray_revert(PyArrayObject *ret)
  * correlate(a2, a1), and conjugate the second argument for complex inputs
  */
 NPY_NO_EXPORT PyObject *
-PyArray_Correlate2(PyObject *op1, PyObject *op2, int mode, int minlag, int maxlag, int lagstep)
+PyArray_Correlate2(PyObject *op1, PyObject *op2, int mode, npy_intp minlag, npy_intp maxlag, npy_intp lagstep)
 {
     PyArrayObject *ap1, *ap2, *ret = NULL;
     int typenum;
@@ -1417,7 +1417,7 @@ PyArray_Correlate2(PyObject *op1, PyObject *op2, int mode, int minlag, int maxla
  * Numeric.correlate(a1,a2,mode)
  */
 NPY_NO_EXPORT PyObject *
-PyArray_Correlate(PyObject *op1, PyObject *op2, int mode, int minlag, int maxlag, int lagstep)
+PyArray_Correlate(PyObject *op1, PyObject *op2, int mode, npy_intp minlag, npy_intp maxlag, npy_intp lagstep)
 {
     PyArrayObject *ap1, *ap2, *ret = NULL;
     int typenum;
@@ -2936,10 +2936,11 @@ static PyObject *
 array_correlate(PyObject *NPY_UNUSED(dummy), PyObject *args, PyObject *kwds)
 {
     PyObject *shape, *a0;
-    int mode = 0, maxlag = NPY_MAX_INTP, minlag = NPY_MAX_INTP, lagstep = NPY_MAX_INTP;
+    int mode = 0;
+    npy_intp maxlag = NPY_MAX_INTP, minlag = NPY_MAX_INTP, lagstep = NPY_MAX_INTP;
     static char *kwlist[] = {"a", "v", "mode", "minlag", "maxlag", "lagstep", NULL};
 
-    if (!PyArg_ParseTupleAndKeywords(args, kwds, "OO|iiii", kwlist,
+    if (!PyArg_ParseTupleAndKeywords(args, kwds, "OO|innn", kwlist,
                 &a0, &shape, &mode, &minlag, &maxlag, &lagstep)) {
         return NULL;
     }
@@ -2950,10 +2951,12 @@ static PyObject *
 array_correlate2(PyObject *NPY_UNUSED(dummy), PyObject *args, PyObject *kwds)
 {
     PyObject *shape, *a0;
-    int mode = 0, maxlag = NPY_MAX_INTP, minlag = NPY_MAX_INTP, lagstep = NPY_MAX_INTP;
+    int mode = 0;
+    npy_intp maxlag = NPY_MAX_INTP, minlag = NPY_MAX_INTP, lagstep = NPY_MAX_INTP;
+    PyObject *tmp;
     static char *kwlist[] = {"a", "v", "mode", "minlag", "maxlag", "lagstep", NULL};
 
-    if (!PyArg_ParseTupleAndKeywords(args, kwds, "OO|iiii", kwlist,
+    if (!PyArg_ParseTupleAndKeywords(args, kwds, "OO|innn", kwlist,
                 &a0, &shape, &mode, &minlag, &maxlag, &lagstep)) {
         return NULL;
     }
diff --git a/numpy/core/tests/test_numeric.py b/numpy/core/tests/test_numeric.py
index 4adf12296d75..8df19265ee55 100644
--- a/numpy/core/tests/test_numeric.py
+++ b/numpy/core/tests/test_numeric.py
@@ -1966,11 +1966,11 @@ def _setup(self, dt):
     def test_float(self):
         self._setup(np.float)
         lagvec = np.arange(1)
-        z, lagvec = np.correlate(self.x, self.y, 'full', old_behavior=self.old_behavior,lagvector=True)
+        z, lagvec = np.correlate(self.x, self.y, 'full', old_behavior=self.old_behavior,returns_lags=True)
         assert_array_almost_equal(z, self.z1)
         assert_equal(lagvec.size, z.size)
         assert_array_equal(lagvec, np.arange(-self.y.size + 1, self.x.size))
-        z, lagvec = np.correlate(self.x, self.y, 'same', old_behavior=self.old_behavior, lagvector=True)
+        z, lagvec = np.correlate(self.x, self.y, 'same', old_behavior=self.old_behavior, returns_lags=True)
         assert_array_almost_equal(z, self.z1s)
         assert_array_almost_equal(z, self.z1[lagvec+self.y.size-1])
         assert_equal(lagvec.size, z.size)
@@ -1978,7 +1978,7 @@ def test_float(self):
             assert_array_equal(lagvec, np.arange(-int(self.y.size/2), self.x.size - int(self.y.size/2)))
         else:
             assert_array_equal(lagvec, np.arange(-self.y.size + int(self.x.size/2) + 1, int(self.x.size/2) + 1))
-        z, lagvec = np.correlate(self.x, self.y, 'valid', old_behavior=self.old_behavior, lagvector=True)
+        z, lagvec = np.correlate(self.x, self.y, 'valid', old_behavior=self.old_behavior, returns_lags=True)
         assert_array_almost_equal(z, self.z1v)
         assert_array_almost_equal(z, self.z1[lagvec+self.y.size-1])
         assert_equal(lagvec.size, z.size)
@@ -1989,11 +1989,11 @@ def test_float(self):
 
         z = np.correlate(self.x, self.y[:-1], 'full', old_behavior=self.old_behavior)
         assert_array_almost_equal(z, self.z1_4)
-        z, lagvec = np.correlate(self.y, self.x, 'full', old_behavior=self.old_behavior, lagvector=True)
+        z, lagvec = np.correlate(self.y, self.x, 'full', old_behavior=self.old_behavior, returns_lags=True)
         assert_array_almost_equal(z, self.z2)
         assert_equal(lagvec.size, z.size)
         assert_array_equal(lagvec, np.arange(-self.x.size + 1, self.y.size))
-        z, lagvec = np.correlate(self.y, self.x, 'same', old_behavior=self.old_behavior, lagvector=True)
+        z, lagvec = np.correlate(self.y, self.x, 'same', old_behavior=self.old_behavior, returns_lags=True)
         assert_array_almost_equal(z, self.z2s)
         assert_array_almost_equal(z, self.z2[lagvec+self.x.size-1])
         assert_equal(lagvec.size, z.size)
@@ -2001,7 +2001,7 @@ def test_float(self):
             assert_array_equal(lagvec, np.arange(-int(self.x.size/2), self.y.size - int(self.x.size/2)))
         else:
             assert_array_equal(lagvec, np.arange(-self.x.size + int(self.y.size/2) + 1, int(self.y.size/2) + 1))
-        z, lagvec = np.correlate(self.y, self.x, 'valid', old_behavior=self.old_behavior, lagvector=True)
+        z, lagvec = np.correlate(self.y, self.x, 'valid', old_behavior=self.old_behavior, returns_lags=True)
         assert_array_almost_equal(z, self.z2v)
         assert_array_almost_equal(z, self.z2[lagvec+self.x.size-1])
         assert_equal(lagvec.size, z.size)
@@ -2026,49 +2026,49 @@ def test_lags(self):
         tmp = np.correlate(self.x, self.y, 'full')
         z_full = np.zeros(tmp.size + 100)
         z_full[:tmp.size] = tmp
-        z, lagvec = np.correlate(self.x, self.y, maxlag, lagvector=True)
+        z, lagvec = np.correlate(self.x, self.y, maxlag, returns_lags=True)
         assert_array_equal(z[0:3], np.zeros(3))
         assert_array_equal(z[-3:], np.zeros(3))
         assert_equal(lagvec.size, z.size)
         assert_array_equal(lagvec, np.arange(-maxlag+1, maxlag))
-        z, lagvec = np.correlate(self.x, self.y, (minlag, maxlag), lagvector=True)
+        z, lagvec = np.correlate(self.x, self.y, (minlag, maxlag), returns_lags=True)
         assert_array_almost_equal(z, z_full[lagvec+self.y.size-1])
         assert_equal(lagvec.size, z.size)
         assert_array_equal(lagvec, np.arange(minlag, maxlag))
-        z, lagvec = np.correlate(self.x, self.y, (minlag, maxlag, lagstep), lagvector=True)
+        z, lagvec = np.correlate(self.x, self.y, (minlag, maxlag, lagstep), returns_lags=True)
         assert_array_almost_equal(z, z_full[lagvec+self.y.size-1])
         assert_equal(lagvec.size, z.size)
         assert_array_equal(lagvec, np.arange(minlag, maxlag, lagstep))
         lagstep = 3
-        z, lagvec = np.correlate(self.x, self.y, (minlag, maxlag, lagstep), lagvector=True)
+        z, lagvec = np.correlate(self.x, self.y, (minlag, maxlag, lagstep), returns_lags=True)
         assert_array_almost_equal(z, z_full[lagvec+self.y.size-1])
         assert_equal(lagvec.size, z.size)
         minlag = 10
         maxlag = -2
         lagstep = -2
-        z, lagvec = np.correlate(self.x, self.y, (minlag, maxlag, lagstep), lagvector=True)
+        z, lagvec = np.correlate(self.x, self.y, (minlag, maxlag, lagstep), returns_lags=True)
         assert_array_almost_equal(z, z_full[lagvec+self.y.size-1])
         assert_equal(lagvec.size, z.size)
         assert_array_equal(lagvec, np.arange(int(minlag), int(maxlag), lagstep))
         maxlag = 10.2
         minlag = -2.2
-        z, lagvec = np.correlate(self.x, self.y, (minlag, maxlag), lagvector=True)
+        z, lagvec = np.correlate(self.x, self.y, (minlag, maxlag), returns_lags=True)
         assert_array_almost_equal(z, z_full[lagvec+self.y.size-1])
         assert_equal(lagvec.size, z.size)
         assert_array_equal(lagvec, np.arange(int(minlag), int(maxlag)))
-        z, lagvec = np.correlate(self.x, self.y, (maxlag, minlag), lagvector=True)
+        z, lagvec = np.correlate(self.x, self.y, (maxlag, minlag), returns_lags=True)
         assert_array_almost_equal(z, z_full[lagvec+self.y.size-1])
         assert_equal(lagvec.size, z.size)
         assert_array_equal(lagvec, np.arange(int(maxlag), int(minlag)))
         maxlag = -1
         minlag = -2
-        z, lagvec = np.correlate(self.x, self.y, (minlag, maxlag), lagvector=True)
+        z, lagvec = np.correlate(self.x, self.y, (minlag, maxlag), returns_lags=True)
         assert_array_almost_equal(z, z_full[lagvec+self.y.size-1])
         assert_equal(lagvec.size, z.size)
         assert_array_equal(lagvec, np.arange(int(minlag), int(maxlag)))
         maxlag = 2
         minlag = 4
-        z, lagvec = np.correlate(self.x, self.y, (minlag, maxlag), lagvector=True)
+        z, lagvec = np.correlate(self.x, self.y, (minlag, maxlag), returns_lags=True)
         assert_array_almost_equal(z, z_full[lagvec+self.y.size-1])
         assert_equal(lagvec.size, z.size)
         assert_array_equal(lagvec, np.arange(int(minlag), int(maxlag)))

From 8cd6e3bd782ca86f24b470f38a9d74b8b23c1ece Mon Sep 17 00:00:00 2001
From: Honi Sanders <honi@brandeis.edu>
Date: Wed, 24 Jun 2015 15:11:29 -0400
Subject: [PATCH 03/12] STY: shortened lines in multiarraymodule.c

---
 numpy/core/numeric.py                        |  2 +-
 numpy/core/src/multiarray/multiarraymodule.c | 45 +++++++++++++-------
 2 files changed, 31 insertions(+), 16 deletions(-)

diff --git a/numpy/core/numeric.py b/numpy/core/numeric.py
index 84ed98844e14..5a1a5e48e015 100644
--- a/numpy/core/numeric.py
+++ b/numpy/core/numeric.py
@@ -995,7 +995,7 @@ def convolve(a, v, mode='full', returns_lags=False):
           at each point of overlap, with an output shape of (N+M-1,). At
           the end-points of the convolution, the signals do not overlap
           completely, and boundary effects may be seen.  This corresponds
-          with a lag tuple of (-M+1, N-1, 1) for N>M or (-N+1, M+1, 1)
+          with a lag tuple of (-M+1, N, 1) for N>M or (-N+1, M, 1)
           for M>N.
 
         'same':
diff --git a/numpy/core/src/multiarray/multiarraymodule.c b/numpy/core/src/multiarray/multiarraymodule.c
index 9b08c46f7ed3..c57e0ee9e003 100644
--- a/numpy/core/src/multiarray/multiarraymodule.c
+++ b/numpy/core/src/multiarray/multiarraymodule.c
@@ -1142,11 +1142,12 @@ PyArray_CopyAndTranspose(PyObject *op)
  */
 static PyArrayObject*
 _pyarray_correlate(PyArrayObject *ap1, PyArrayObject *ap2, int typenum,
-                   int mode, int *inverted, npy_intp minlag, npy_intp maxlag, npy_intp lagstep)
+                   int mode, int *inverted,
+                   npy_intp minlag, npy_intp maxlag, npy_intp lagstep)
 {
     PyArrayObject *ret;
     npy_intp length;
-    npy_intp i, n1, n2, n, n11;
+    npy_intp i, i1, n1, n2, n, n11;
     npy_intp lag, tmplag, maxleft, maxright;
     npy_intp is1, is2, os;
     char *ip1, *ip2, *op;
@@ -1176,21 +1177,26 @@ _pyarray_correlate(PyArrayObject *ap1, PyArrayObject *ap2, int typenum,
     if (maxlag == NPY_MAX_INTP || minlag == NPY_MAX_INTP) {
     switch(mode) {
     lagstep = 1;
-    case 0:         /* mode = 'valid' */
+    case 0:
+        /* mode = 'valid' */
         minlag = 0;
         maxlag = n1 - n2 + 1;
         break;
-    case 1:         /* mode = 'same' */
+    case 1:
+        /* mode = 'same' */
         minlag = -(npy_intp)(n2/2);
         maxlag = n1 + minlag;
         break;
-    case 2:         /* mode = 'full' */
+    case 2:
+        /* mode = 'full' */
         minlag = -n2 + 1;
         maxlag = n1;
         break;
-    case 3:         /* mode = 'maxlag' */
+    case 3:
+        /* mode = 'maxlag' */
         if (maxlag == NPY_MAX_INTP) {
-            minlag = -n2 + 1;   /* nothing good to do with situation */
+            /* nothing good to do with situation */
+            minlag = -n2 + 1;
             maxlag = n1;
         }
         else if (minlag == NPY_MAX_INTP) {
@@ -1206,7 +1212,8 @@ _pyarray_correlate(PyArrayObject *ap1, PyArrayObject *ap2, int typenum,
     if (lagstep < 0) {
         *inverted = 1;
         i = minlag;
-        minlag = (npy_intp)(npy_ceil((maxlag - minlag)/(float)lagstep))*lagstep + minlag - lagstep;
+        i1 = (npy_intp)(npy_ceil((maxlag - minlag)/(float)lagstep))*lagstep;
+        minlag =  i1 + minlag - lagstep;
         maxlag = i - lagstep;
         lagstep = -lagstep;
     }
@@ -1223,8 +1230,9 @@ _pyarray_correlate(PyArrayObject *ap1, PyArrayObject *ap2, int typenum,
     /*
      * Need to choose an output array that can hold a sum
      * -- use priority to determine which subtype.
+     * ret is the array that will be returned as the answer
      */
-    ret = new_array_for_sum(ap1, ap2, NULL, 1, &length, typenum);   /* answer */
+    ret = new_array_for_sum(ap1, ap2, NULL, 1, &length, typenum);
     if (ret == NULL) {
         return NULL;
     }
@@ -1246,18 +1254,24 @@ _pyarray_correlate(PyArrayObject *ap1, PyArrayObject *ap2, int typenum,
 
     lag = minlag;
     maxleft = (0 < maxlag ? 0 : maxlag);
-    for (lag = minlag; lag < maxleft; lag+=lagstep) {      /* for lags where y is left of x */
+    for (lag = minlag; lag < maxleft; lag+=lagstep) {
+        /* for lags where y is left of x */
         n = n2 + lag;   /* overlap is length of y - lag */
         dot(ip1, is1, ip2 - lag*is2, is2, op, n, ret);
         op += os;       /* iterate over answer vector   */
     }
-    if (maxlag < n1 - n2) { /* if maxlag doesn't take y all the way to the end of x */
-        n11 = maxlag + n2;      /* relevant length of x is smaller                  */
+    if (maxlag < n1 - n2) {
+        /* if maxlag doesn't take y all the way to the end of x */
+        n11 = maxlag + n2;      /* relevant length of x is smaller */
     }
     else {
         n11 = n1;
-    }       /* starts at lag=minlag if minlag>0. Does lags where y entirely overlaps with x. */
-    if (lagstep == 1 && small_correlate(ip1 + lag*is1, is1, n11 - n2 + 1 - lag, PyArray_TYPE(ap1),
+    }
+    /* starts at lag=minlag if minlag>0.
+     * Does lags where y entirely overlaps with x.
+     */
+    if (lagstep == 1 && small_correlate(ip1 + lag*is1, is1,
+                                        n11 - n2 + 1 - lag, PyArray_TYPE(ap1),
                                         ip2, is2, n2, PyArray_TYPE(ap2),
                                         op, os)) {
         lag = n11 - n2 + 1;
@@ -1272,7 +1286,8 @@ _pyarray_correlate(PyArrayObject *ap1, PyArrayObject *ap2, int typenum,
     }
     maxright = (maxlag < n1  ? maxlag : n1 );
     tmplag = lag;
-    for (lag = tmplag; lag < maxright; lag+=lagstep) {     /* for lags where y is right of x */
+    for (lag = tmplag; lag < maxright; lag+=lagstep) {
+        /* for lags where y is right of x */
         n = n1 - lag;
         dot(ip1 + lag*is1, is1, ip2, is2, op, n, ret);
         op += os;

From 85b4e48968dca463d81368510190fcbbb8eae661 Mon Sep 17 00:00:00 2001
From: Honi Sanders <honi@brandeis.edu>
Date: Wed, 24 Jun 2015 18:43:59 -0400
Subject: [PATCH 04/12] DEP,MAINT: Remove old_behavior keyword from
 numeric.correlate.

this is to keep in line with commit bac2fdf
---
 numpy/core/numeric.py            | 33 ++++-------------
 numpy/core/tests/test_numeric.py | 63 +++++++-------------------------
 2 files changed, 22 insertions(+), 74 deletions(-)

diff --git a/numpy/core/numeric.py b/numpy/core/numeric.py
index 5a1a5e48e015..e43a334ee87e 100644
--- a/numpy/core/numeric.py
+++ b/numpy/core/numeric.py
@@ -858,7 +858,7 @@ def _lags_from_mode(alen, vlen, mode):
             lags = (-int(ceil((lags[1]-lags[0])/float(lags[2])))*lags[2]-lags[0]+lags[2], -lags[0]+lags[2], lags[2])
     return mode, lags
 
-def correlate(a, v, mode='valid', old_behavior=False, returns_lags=False):
+def correlate(a, v, mode='valid', returns_lags=False):
     """
     Cross-correlation of two 1-dimensional sequences.
 
@@ -877,17 +877,15 @@ def correlate(a, v, mode='valid', old_behavior=False, returns_lags=False):
     mode : int, int tuple, or {'valid', 'same', 'full'}, optional
         Refer to the `convolve` docstring.  Note that the default
         is `valid`, unlike `convolve`, which uses `full`.
-    old_behavior : bool, optional
-        If True, uses the old behavior from Numeric,
-        (correlate(a,v) == correlate(v,a), and the conjugate is not taken
-        for complex arrays). If False, uses the conventional signal
-        processing definition.  False is default.
     returns_lags : bool, optional
         If True, the function returns a lagvector array in addition to the
         cross-correlation result.  The lagvector contains the indices of
         the lags for which the cross-correlation was calculated.  It is
         the same length as the return array, and corresponds one-to-one.
         False is default.
+    old_behavior : bool
+        `old_behavior` was removed in NumPy 1.10. If you need the old
+        behavior, use `multiarray.correlate`.
 
     Returns
     -------
@@ -900,6 +898,7 @@ def correlate(a, v, mode='valid', old_behavior=False, returns_lags=False):
     See Also
     --------
     convolve : Discrete, linear convolution of two one-dimensional sequences.
+    multiarray.correlate : Old, no conjugate, version of correlate.
 
     Notes
     -----
@@ -937,26 +936,10 @@ def correlate(a, v, mode='valid', old_behavior=False, returns_lags=False):
 
     """
     mode, lags = _lags_from_mode(len(a), len(v), mode)
-
-# the old behavior should be made available under a different name, see thread
-# http://thread.gmane.org/gmane.comp.python.numeric.general/12609/focus=12630
-    if old_behavior:
-        warnings.warn("""
-The old behavior of correlate was deprecated for 1.4.0, and will be completely removed
-for NumPy 2.0.
-
-The new behavior fits the conventional definition of correlation: inputs are
-never swapped, and the second argument is conjugated for complex arrays.""",
-            DeprecationWarning)
-        if returns_lags:
-            return multiarray.correlate(a, v, mode, lags[0], lags[1], lags[2]), arange(lags[0], lags[1], lags[2])
-        else:
-            return multiarray.correlate(a, v, mode, lags[0], lags[1], lags[2])
+    if returns_lags:
+        return multiarray.correlate2(a, v, mode, lags[0], lags[1], lags[2]), arange(lags[0], lags[1], lags[2])
     else:
-        if returns_lags:
-            return multiarray.correlate2(a, v, mode, lags[0], lags[1], lags[2]), arange(lags[0], lags[1], lags[2])
-        else:
-            return multiarray.correlate2(a, v, mode, lags[0], lags[1], lags[2])
+        return multiarray.correlate2(a, v, mode, lags[0], lags[1], lags[2])
 
 def convolve(a, v, mode='full', returns_lags=False):
     """
diff --git a/numpy/core/tests/test_numeric.py b/numpy/core/tests/test_numeric.py
index 8df19265ee55..8df9fa42c0c6 100644
--- a/numpy/core/tests/test_numeric.py
+++ b/numpy/core/tests/test_numeric.py
@@ -1946,7 +1946,7 @@ def test_filled_like(self):
         self.check_like_function(np.full_like, 123.456, True)
         self.check_like_function(np.full_like, np.inf, True)
 
-class _TestCorrelate(TestCase):
+class TestCorrelate(TestCase):
     def _setup(self, dt):
         self.x = np.array([1, 2, 3, 4, 5], dtype=dt)
         self.xs = np.arange(1, 20)[::3]
@@ -1965,12 +1965,11 @@ def _setup(self, dt):
 
     def test_float(self):
         self._setup(np.float)
-        lagvec = np.arange(1)
-        z, lagvec = np.correlate(self.x, self.y, 'full', old_behavior=self.old_behavior,returns_lags=True)
+        z, lagvec = np.correlate(self.x, self.y, 'full', returns_lags=True)
         assert_array_almost_equal(z, self.z1)
         assert_equal(lagvec.size, z.size)
         assert_array_equal(lagvec, np.arange(-self.y.size + 1, self.x.size))
-        z, lagvec = np.correlate(self.x, self.y, 'same', old_behavior=self.old_behavior, returns_lags=True)
+        z, lagvec = np.correlate(self.x, self.y, 'same', returns_lags=True)
         assert_array_almost_equal(z, self.z1s)
         assert_array_almost_equal(z, self.z1[lagvec+self.y.size-1])
         assert_equal(lagvec.size, z.size)
@@ -1978,7 +1977,7 @@ def test_float(self):
             assert_array_equal(lagvec, np.arange(-int(self.y.size/2), self.x.size - int(self.y.size/2)))
         else:
             assert_array_equal(lagvec, np.arange(-self.y.size + int(self.x.size/2) + 1, int(self.x.size/2) + 1))
-        z, lagvec = np.correlate(self.x, self.y, 'valid', old_behavior=self.old_behavior, returns_lags=True)
+        z, lagvec = np.correlate(self.x, self.y, 'valid', returns_lags=True)
         assert_array_almost_equal(z, self.z1v)
         assert_array_almost_equal(z, self.z1[lagvec+self.y.size-1])
         assert_equal(lagvec.size, z.size)
@@ -1987,13 +1986,13 @@ def test_float(self):
         else:
             assert_array_equal(lagvec, np.arange(self.x.size - self.y.size, 1))
 
-        z = np.correlate(self.x, self.y[:-1], 'full', old_behavior=self.old_behavior)
+        z = np.correlate(self.x, self.y[:-1], 'full')
         assert_array_almost_equal(z, self.z1_4)
-        z, lagvec = np.correlate(self.y, self.x, 'full', old_behavior=self.old_behavior, returns_lags=True)
+        z, lagvec = np.correlate(self.y, self.x, 'full', returns_lags=True)
         assert_array_almost_equal(z, self.z2)
         assert_equal(lagvec.size, z.size)
         assert_array_equal(lagvec, np.arange(-self.x.size + 1, self.y.size))
-        z, lagvec = np.correlate(self.y, self.x, 'same', old_behavior=self.old_behavior, returns_lags=True)
+        z, lagvec = np.correlate(self.y, self.x, 'same', returns_lags=True)
         assert_array_almost_equal(z, self.z2s)
         assert_array_almost_equal(z, self.z2[lagvec+self.x.size-1])
         assert_equal(lagvec.size, z.size)
@@ -2001,7 +2000,7 @@ def test_float(self):
             assert_array_equal(lagvec, np.arange(-int(self.x.size/2), self.y.size - int(self.x.size/2)))
         else:
             assert_array_equal(lagvec, np.arange(-self.x.size + int(self.y.size/2) + 1, int(self.y.size/2) + 1))
-        z, lagvec = np.correlate(self.y, self.x, 'valid', old_behavior=self.old_behavior, returns_lags=True)
+        z, lagvec = np.correlate(self.y, self.x, 'valid', returns_lags=True)
         assert_array_almost_equal(z, self.z2v)
         assert_array_almost_equal(z, self.z2[lagvec+self.x.size-1])
         assert_equal(lagvec.size, z.size)
@@ -2009,17 +2008,16 @@ def test_float(self):
             assert_array_equal(lagvec, np.arange(0, self.y.size - self.x.size + 1))
         else:
             assert_array_equal(lagvec, np.arange(self.y.size - self.x.size, 1))
-        z = np.correlate(self.x[::-1], self.y, 'full', old_behavior=self.old_behavior)
+        z = np.correlate(self.x[::-1], self.y, 'full')
         assert_array_almost_equal(z, self.z1r)
-        z = np.correlate(self.y, self.x[::-1], 'full', old_behavior=self.old_behavior)
+        z = np.correlate(self.y, self.x[::-1], 'full')
         assert_array_almost_equal(z, self.z2r)
-        z = np.correlate(self.xs, self.y, 'full', old_behavior=self.old_behavior)
+        z = np.correlate(self.xs, self.y, 'full')
         assert_array_almost_equal(z, self.zs)
 
     def test_lags(self):
         self._setup(np.float)
         longlen = max(self.x.size, self.y.size)
-        lagvec = np.arange(1)
         maxlag = 3 + longlen
         minlag = -2
         lagstep = 2
@@ -2077,9 +2075,9 @@ def test_lags(self):
 
     def test_object(self):
         self._setup(Decimal)
-        z = np.correlate(self.x, self.y, 'full', old_behavior=self.old_behavior)
+        z = np.correlate(self.x, self.y, 'full')
         assert_array_almost_equal(z, self.z1)
-        z = np.correlate(self.y, self.x, 'full', old_behavior=self.old_behavior)
+        z = np.correlate(self.y, self.x, 'full')
         assert_array_almost_equal(z, self.z2)
 
     def test_no_overwrite(self):
@@ -2089,45 +2087,12 @@ def test_no_overwrite(self):
         assert_array_equal(d, np.ones(100))
         assert_array_equal(k, np.ones(3))
 
-class TestCorrelate(_TestCorrelate):
-    old_behavior = True
-    def _setup(self, dt):
-        # correlate uses an unconventional definition so that correlate(a, b)
-        # == correlate(b, a), so force the corresponding outputs to be the same
-        # as well
-        _TestCorrelate._setup(self, dt)
-        self.z2 = self.z1
-        self.z2r = self.z1r
-        self.z2s = self.z1s
-        self.z2v = self.z1v
-
-    @dec.deprecated()
-    def test_complex(self):
-        x = np.array([1, 2, 3, 4+1j], dtype=np.complex)
-        y = np.array([-1, -2j, 3+1j], dtype=np.complex)
-        r_z = np.array([3+1j, 6, 8-1j, 9+1j, -1-8j, -4-1j], dtype=np.complex)
-        z = np.correlate(x, y, 'full', old_behavior=self.old_behavior)
-        assert_array_almost_equal(z, r_z)
-
-    @dec.deprecated()
-    def test_float(self):
-        _TestCorrelate.test_float(self)
-
-    @dec.deprecated()
-    def test_object(self):
-        _TestCorrelate.test_object(self)
-
-class TestCorrelateNew(_TestCorrelate):
-    old_behavior = False
     def test_complex(self):
         x = np.array([1, 2, 3, 4+1j], dtype=np.complex)
         y = np.array([-1, -2j, 3+1j], dtype=np.complex)
         r_z = np.array([3-1j, 6, 8+1j, 11+5j, -5+8j, -4-1j], dtype=np.complex)
-        #z = np.acorrelate(x, y, 'full')
-        #assert_array_almost_equal(z, r_z)
-
         r_z = r_z[::-1].conjugate()
-        z = np.correlate(y, x, 'full', old_behavior=self.old_behavior)
+        z = np.correlate(y, x, 'full')
         assert_array_almost_equal(z, r_z)
 
 class TestConvolve(TestCase):

From e1acde80311d55c326af6bfb46f03a5ef4797096 Mon Sep 17 00:00:00 2001
From: Honi Sanders <honi@brandeis.edu>
Date: Mon, 26 Oct 2015 20:13:48 -0400
Subject: [PATCH 05/12] BUG: fixed default arguments in correlate pipeline (See
 #5954) _pyarray_correlate ignores minlag, maxlag, and lagstep for modes other
 than 3

---
 numpy/core/numeric.py                        |  8 +++---
 numpy/core/src/multiarray/multiarraymodule.c | 30 +++++++-------------
 2 files changed, 15 insertions(+), 23 deletions(-)

diff --git a/numpy/core/numeric.py b/numpy/core/numeric.py
index e43a334ee87e..0d8b24f6b0a1 100644
--- a/numpy/core/numeric.py
+++ b/numpy/core/numeric.py
@@ -937,9 +937,9 @@ def correlate(a, v, mode='valid', returns_lags=False):
     """
     mode, lags = _lags_from_mode(len(a), len(v), mode)
     if returns_lags:
-        return multiarray.correlate2(a, v, mode, lags[0], lags[1], lags[2]), arange(lags[0], lags[1], lags[2])
+        return multiarray.correlate2(a, v, 3, lags[0], lags[1], lags[2]), arange(lags[0], lags[1], lags[2])
     else:
-        return multiarray.correlate2(a, v, mode, lags[0], lags[1], lags[2])
+        return multiarray.correlate2(a, v, 3, lags[0], lags[1], lags[2])
 
 def convolve(a, v, mode='full', returns_lags=False):
     """
@@ -1083,9 +1083,9 @@ def convolve(a, v, mode='full', returns_lags=False):
         raise ValueError('v cannot be empty')
     mode, lags = _lags_from_mode(len(a), len(v), mode)
     if returns_lags:
-        return multiarray.correlate2(a, v[::-1], mode, lags[0], lags[1], lags[2]), arange(lags[0], lags[1], lags[2])
+        return multiarray.correlate2(a, v[::-1], 3, lags[0], lags[1], lags[2]), arange(lags[0], lags[1], lags[2])
     else:
-        return multiarray.correlate2(a, v[::-1], mode, lags[0], lags[1], lags[2])
+        return multiarray.correlate2(a, v[::-1], 3, lags[0], lags[1], lags[2])
 
 def outer(a, b, out=None):
     """
diff --git a/numpy/core/src/multiarray/multiarraymodule.c b/numpy/core/src/multiarray/multiarraymodule.c
index c57e0ee9e003..a5f4d74637c2 100644
--- a/numpy/core/src/multiarray/multiarraymodule.c
+++ b/numpy/core/src/multiarray/multiarraymodule.c
@@ -1155,7 +1155,7 @@ _pyarray_correlate(PyArrayObject *ap1, PyArrayObject *ap2, int typenum,
 
     NPY_BEGIN_THREADS_DEF;
 
-    if (lagstep == NPY_MAX_INTP || lagstep == 0) {
+    if (lagstep == 0 || mode != 3) {
         lagstep = 1;
     }
 
@@ -1174,9 +1174,7 @@ _pyarray_correlate(PyArrayObject *ap1, PyArrayObject *ap2, int typenum,
         lagstep = -lagstep;
     }
 
-    if (maxlag == NPY_MAX_INTP || minlag == NPY_MAX_INTP) {
     switch(mode) {
-    lagstep = 1;
     case 0:
         /* mode = 'valid' */
         minlag = 0;
@@ -1194,20 +1192,13 @@ _pyarray_correlate(PyArrayObject *ap1, PyArrayObject *ap2, int typenum,
         break;
     case 3:
         /* mode = 'maxlag' */
-        if (maxlag == NPY_MAX_INTP) {
-            /* nothing good to do with situation */
-            minlag = -n2 + 1;
-            maxlag = n1;
-        }
-        else if (minlag == NPY_MAX_INTP) {
-            minlag = (maxlag > n2 ? -n2+1 : -maxlag+1);
-        }
+        /* use minlag, maxlag, and lagstep that were passed in */
         break;
     default:
         PyErr_SetString(PyExc_ValueError, "mode must be 0, 1, 2, or 3");
         return NULL;
     }
-    }
+
 
     if (lagstep < 0) {
         *inverted = 1;
@@ -1270,10 +1261,11 @@ _pyarray_correlate(PyArrayObject *ap1, PyArrayObject *ap2, int typenum,
     /* starts at lag=minlag if minlag>0.
      * Does lags where y entirely overlaps with x.
      */
-    if (lagstep == 1 && small_correlate(ip1 + lag*is1, is1,
-                                        n11 - n2 + 1 - lag, PyArray_TYPE(ap1),
-                                        ip2, is2, n2, PyArray_TYPE(ap2),
-                                        op, os)) {
+    if (lagstep == 1 && lag < maxlag &&
+            small_correlate(ip1 + lag*is1, is1,
+                            n11 - n2 + 1 - lag, PyArray_TYPE(ap1),
+                            ip2, is2, n2, PyArray_TYPE(ap2),
+                            op, os)) {
         lag = n11 - n2 + 1;
         op += os * (n11 - n2 + 1);
     }
@@ -2951,8 +2943,8 @@ static PyObject *
 array_correlate(PyObject *NPY_UNUSED(dummy), PyObject *args, PyObject *kwds)
 {
     PyObject *shape, *a0;
-    int mode = 0;
-    npy_intp maxlag = NPY_MAX_INTP, minlag = NPY_MAX_INTP, lagstep = NPY_MAX_INTP;
+    int mode = 0; /* input modes other than 3 cause minlag, maxlag, and lagstep to be ignored */
+    npy_intp maxlag = 0, minlag = 0, lagstep = 0;
     static char *kwlist[] = {"a", "v", "mode", "minlag", "maxlag", "lagstep", NULL};
 
     if (!PyArg_ParseTupleAndKeywords(args, kwds, "OO|innn", kwlist,
@@ -2967,7 +2959,7 @@ array_correlate2(PyObject *NPY_UNUSED(dummy), PyObject *args, PyObject *kwds)
 {
     PyObject *shape, *a0;
     int mode = 0;
-    npy_intp maxlag = NPY_MAX_INTP, minlag = NPY_MAX_INTP, lagstep = NPY_MAX_INTP;
+    npy_intp maxlag = 0, minlag = 0, lagstep = 0;
     PyObject *tmp;
     static char *kwlist[] = {"a", "v", "mode", "minlag", "maxlag", "lagstep", NULL};
 

From eb0dda5531c2e200d1fadd6f7b5573111b304125 Mon Sep 17 00:00:00 2001
From: Honi Sanders <honi@brandeis.edu>
Date: Tue, 27 Oct 2015 18:21:26 -0400
Subject: [PATCH 06/12] API: added new C API function PyArray_CorrelateLags to
 multiarraymodule.c

Removed correlate and corresponding array_correlate from multiarraymodule.c
Reason: old_behavior has been deprecated from correlate and convolve in numeric.py

Also:
- restore PyArray_Correlate and PyArray_Correlate2 signatures
- add bounds checking to maxlag and minlag in _pyarray_correlate.
- minor cleanups in correlate/convolve code.
---
 numpy/core/code_generators/numpy_api.py      |  1 +
 numpy/core/src/multiarray/multiarraymodule.c | 70 ++++++++++++--------
 2 files changed, 42 insertions(+), 29 deletions(-)

diff --git a/numpy/core/code_generators/numpy_api.py b/numpy/core/code_generators/numpy_api.py
index 972966627b99..5ebe93ff0b52 100644
--- a/numpy/core/code_generators/numpy_api.py
+++ b/numpy/core/code_generators/numpy_api.py
@@ -344,6 +344,7 @@
     # End 1.9 API
     'PyArray_CheckAnyScalarExact':          (300, NonNull(1)),
     # End 1.10 API
+    'PyArray_CorrelateLags':                (301,),
 }
 
 ufunc_types_api = {
diff --git a/numpy/core/src/multiarray/multiarraymodule.c b/numpy/core/src/multiarray/multiarraymodule.c
index a5f4d74637c2..a101922202f8 100644
--- a/numpy/core/src/multiarray/multiarraymodule.c
+++ b/numpy/core/src/multiarray/multiarraymodule.c
@@ -1134,6 +1134,8 @@ PyArray_CopyAndTranspose(PyObject *op)
     return (PyObject *)ret;
 }
 
+
+
 /*
  * Implementation which is common between PyArray_Correlate
  * and PyArray_Correlate2.
@@ -1227,7 +1229,9 @@ _pyarray_correlate(PyArrayObject *ap1, PyArrayObject *ap2, int typenum,
     if (ret == NULL) {
         return NULL;
     }
-    PyArray_FILLWBYTE(ret, 0);
+    if (length > 0) {
+        PyArray_FILLWBYTE(ret, 0);
+    }
     dot = PyArray_DESCR(ret)->f->dotfunc;
     if (dot == NULL) {
         PyErr_SetString(PyExc_ValueError,
@@ -1244,8 +1248,16 @@ _pyarray_correlate(PyArrayObject *ap1, PyArrayObject *ap2, int typenum,
     os = PyArray_DESCR(ret)->elsize;    /* answer strides   */
 
     lag = minlag;
+    if (lag < -n2+1) {
+        /* if minlag is before any overlap between the vectors,
+         * then skip to first relevant lag
+         */
+        op += os*((-n2+1) - lag);
+        lag = -n2+1;
+    }
     maxleft = (0 < maxlag ? 0 : maxlag);
-    for (lag = minlag; lag < maxleft; lag+=lagstep) {
+    tmplag = lag;
+    for (lag = tmplag; lag < maxleft; lag+=lagstep) {
         /* for lags where y is left of x */
         n = n2 + lag;   /* overlap is length of y - lag */
         dot(ip1, is1, ip2 - lag*is2, is2, op, n, ret);
@@ -1269,7 +1281,7 @@ _pyarray_correlate(PyArrayObject *ap1, PyArrayObject *ap2, int typenum,
         lag = n11 - n2 + 1;
         op += os * (n11 - n2 + 1);
     }
-    else {
+    else if (lag < maxlag) {
         tmplag = lag;
         for (lag = tmplag; lag < (n11 - n2 + 1); lag+=lagstep) {
             dot(ip1 + lag*is1, is1, ip2, is2, op, n2, ret);
@@ -1350,13 +1362,13 @@ _pyarray_revert(PyArrayObject *ret)
 }
 
 /*NUMPY_API
- * correlate(a1,a2,mode)
+ * correlate(a1,a2,mode,minlag,maxlag,lagstep)
  *
- * This function computes the usual correlation (correlate(a1, a2) !=
- * correlate(a2, a1), and conjugate the second argument for complex inputs
+ * This function correlates a1 with a2.  mode=0,1,2 retains functionality of PyArray_Correlate2
+ * mode=3 allows for specification of minlag, maxlag, and lagstep.
  */
 NPY_NO_EXPORT PyObject *
-PyArray_Correlate2(PyObject *op1, PyObject *op2, int mode, npy_intp minlag, npy_intp maxlag, npy_intp lagstep)
+PyArray_CorrelateLags(PyObject *op1, PyObject *op2, int mode, npy_intp minlag, npy_intp maxlag, npy_intp lagstep)
 {
     PyArrayObject *ap1, *ap2, *ret = NULL;
     int typenum;
@@ -1420,17 +1432,36 @@ PyArray_Correlate2(PyObject *op1, PyObject *op2, int mode, npy_intp minlag, npy_
     return NULL;
 }
 
+/*NUMPY_API
+ * correlate(a1,a2,mode)
+ *
+ * This function computes the usual correlation (correlate(a1, a2) !=
+ * correlate(a2, a1), and conjugate the second argument for complex inputs
+ */
+NPY_NO_EXPORT PyObject *
+PyArray_Correlate2(PyObject *op1, PyObject *op2, int mode)
+{
+    npy_intp z = 0;
+    /* For modes other than 3 (which was not a mode before maxlags were added),
+    minlag, maxlag, and lagstep (the last three arguments) will be ignored */
+    return PyArray_CorrelateLags(op1, op2, mode, z, z, z);
+}
+
 /*NUMPY_API
  * Numeric.correlate(a1,a2,mode)
  */
 NPY_NO_EXPORT PyObject *
-PyArray_Correlate(PyObject *op1, PyObject *op2, int mode, npy_intp minlag, npy_intp maxlag, npy_intp lagstep)
+PyArray_Correlate(PyObject *op1, PyObject *op2, int mode)
 {
     PyArrayObject *ap1, *ap2, *ret = NULL;
     int typenum;
     int unused;
     PyArray_Descr *typec;
 
+    npy_intp z = 0;
+    /* For modes other than 3 (which was not a mode before maxlags were added),
+    minlag, maxlag, and lagstep (the last three arguments of _pyarray_correlate) will be ignored */
+
     typenum = PyArray_ObjectType(op1, 0);
     typenum = PyArray_ObjectType(op2, typenum);
 
@@ -1448,7 +1479,7 @@ PyArray_Correlate(PyObject *op1, PyObject *op2, int mode, npy_intp minlag, npy_i
         goto fail;
     }
 
-    ret = _pyarray_correlate(ap1, ap2, typenum, mode, &unused, minlag, maxlag, lagstep);
+    ret = _pyarray_correlate(ap1, ap2, typenum, mode, &unused, z, z, z);
     if(ret == NULL) {
         goto fail;
     }
@@ -2939,35 +2970,19 @@ array_fastCopyAndTranspose(PyObject *NPY_UNUSED(dummy), PyObject *args)
     return PyArray_Return((PyArrayObject *)PyArray_CopyAndTranspose(a0));
 }
 
-static PyObject *
-array_correlate(PyObject *NPY_UNUSED(dummy), PyObject *args, PyObject *kwds)
-{
-    PyObject *shape, *a0;
-    int mode = 0; /* input modes other than 3 cause minlag, maxlag, and lagstep to be ignored */
-    npy_intp maxlag = 0, minlag = 0, lagstep = 0;
-    static char *kwlist[] = {"a", "v", "mode", "minlag", "maxlag", "lagstep", NULL};
-
-    if (!PyArg_ParseTupleAndKeywords(args, kwds, "OO|innn", kwlist,
-                &a0, &shape, &mode, &minlag, &maxlag, &lagstep)) {
-        return NULL;
-    }
-    return PyArray_Correlate(a0, shape, mode, minlag, maxlag, lagstep);
-}
-
 static PyObject *
 array_correlate2(PyObject *NPY_UNUSED(dummy), PyObject *args, PyObject *kwds)
 {
     PyObject *shape, *a0;
     int mode = 0;
     npy_intp maxlag = 0, minlag = 0, lagstep = 0;
-    PyObject *tmp;
     static char *kwlist[] = {"a", "v", "mode", "minlag", "maxlag", "lagstep", NULL};
 
     if (!PyArg_ParseTupleAndKeywords(args, kwds, "OO|innn", kwlist,
                 &a0, &shape, &mode, &minlag, &maxlag, &lagstep)) {
         return NULL;
     }
-    return PyArray_Correlate2(a0, shape, mode, minlag, maxlag, lagstep);
+    return PyArray_CorrelateLags(a0, shape, mode, minlag, maxlag, lagstep);
 }
 
 static PyObject *
@@ -4154,9 +4169,6 @@ static struct PyMethodDef array_module_methods[] = {
     {"_fastCopyAndTranspose",
         (PyCFunction)array_fastCopyAndTranspose,
         METH_VARARGS, NULL},
-    {"correlate",
-        (PyCFunction)array_correlate,
-        METH_VARARGS | METH_KEYWORDS, NULL},
     {"correlate2",
         (PyCFunction)array_correlate2,
         METH_VARARGS | METH_KEYWORDS, NULL},

From 380c3fa52f3cc2ad25788688fa4ecfed387f3d93 Mon Sep 17 00:00:00 2001
From: Honi Sanders <honi@brandeis.edu>
Date: Tue, 10 Nov 2015 23:00:22 -0500
Subject: [PATCH 07/12] DOC: add documentation of new C API function
 PyArray_CorrelateLags

---
 doc/source/reference/c-api.array.rst    | 30 +++++++++++++++++++++++++
 numpy/core/code_generators/numpy_api.py |  1 +
 2 files changed, 31 insertions(+)

diff --git a/doc/source/reference/c-api.array.rst b/doc/source/reference/c-api.array.rst
index f5f753292098..d49b54eb19d0 100644
--- a/doc/source/reference/c-api.array.rst
+++ b/doc/source/reference/c-api.array.rst
@@ -2161,6 +2161,36 @@ Array Functions
 
       z[k] = sum_n op1[n] * conj(op2[n+k])
 
+.. cfunction:: PyObject* PyArray_CorrelateLags(PyObject* op1, PyObject* op2, int mode, npy_intp minlag, npy_intp maxlag, npy_intp lagstep)
+
+    Updated version of PyArray_Correlate2, which allows for custom lags in calculation.
+    The correlation is computed at each output point by multiplying *op1* by
+    a shifted version of *op2* and summing the result.
+    As a result of the shift, needed values outside of the defined range of
+    *op1* and *op2* are interpreted as zero. The mode determines how many
+    shifts to return: 0 - return only shifts that did not need to assume zero-
+    values; 1 - return an object that is the same size as *op1*, 2 - return all
+    possible shifts (any overlap at all is accepted).
+
+    This version differs from the other PyArray_Correlate* functions in the addition of a
+    third mode, mode = 3.  In this mode, PyArray_CorrelateLags calculates the correlation
+    for lags starting at minlag and ending at (maxlag - 1), with steps of size lagstep.
+    The lags for which the correlation will be calculated correspond with the values in
+    the vector formed by numpy.arange((minlag, maxlag, lagstep)).
+
+    .. rubric:: Notes
+
+    If the mode argument is 0, 1, or 2, the values of minlag, maxlag, and lagstep will be
+    ignored.  Instead, the following lag values will be used, corresponding to the
+    description of modes above (n1 and n2 are the vector lengths):
+       mode   |   minlag   |   maxlag   |   lagstep
+         0    |     0      |  n1-n2+1   |      1
+         1    |   -n2/2    |  n1-n2/2   |      1
+         2    |   -n2+1    |     n1     |      1
+
+    PyArray_CorrelateLags uses the same definition of correlation as PyArray_Correlate2
+    in that the conjugate is taken and the argument order matters.
+
 .. cfunction:: PyObject* PyArray_Where(PyObject* condition, PyObject* x, PyObject* y)
 
     If both ``x`` and ``y`` are ``NULL``, then return
diff --git a/numpy/core/code_generators/numpy_api.py b/numpy/core/code_generators/numpy_api.py
index 5ebe93ff0b52..d6a0fbdd70e6 100644
--- a/numpy/core/code_generators/numpy_api.py
+++ b/numpy/core/code_generators/numpy_api.py
@@ -345,6 +345,7 @@
     'PyArray_CheckAnyScalarExact':          (300, NonNull(1)),
     # End 1.10 API
     'PyArray_CorrelateLags':                (301,),
+    # End 1.11 API
 }
 
 ufunc_types_api = {

From 4b0a9e05221551c9191bdab5b6d28ca4e64dfb3b Mon Sep 17 00:00:00 2001
From: Honi Sanders <honi@brandeis.edu>
Date: Tue, 21 Jun 2016 12:29:30 -0400
Subject: [PATCH 08/12] STY: fixed comments and line lengths in
 multiarraymodule.c to conform to C-style guide

---
 numpy/core/src/multiarray/multiarraymodule.c | 58 +++++++++++++-------
 1 file changed, 37 insertions(+), 21 deletions(-)

diff --git a/numpy/core/src/multiarray/multiarraymodule.c b/numpy/core/src/multiarray/multiarraymodule.c
index 0a7cf97f4b24..9fef139d65bc 100644
--- a/numpy/core/src/multiarray/multiarraymodule.c
+++ b/numpy/core/src/multiarray/multiarraymodule.c
@@ -1111,8 +1111,9 @@ _pyarray_correlate(PyArrayObject *ap1, PyArrayObject *ap2, int typenum,
         lagstep = 1;
     }
 
-    n1 = PyArray_DIMS(ap1)[0];      /* size of x */
-    n2 = PyArray_DIMS(ap2)[0];      /* size of y */
+    /* size of x (n1) and y (n2) */
+    n1 = PyArray_DIMS(ap1)[0];
+    n2 = PyArray_DIMS(ap2)[0];
     if (n1 < n2) {
         ret = ap1;
         ap1 = ap2;
@@ -1190,12 +1191,14 @@ _pyarray_correlate(PyArrayObject *ap1, PyArrayObject *ap2, int typenum,
     }
 
     NPY_BEGIN_THREADS_DESCR(PyArray_DESCR(ret));
-    ip1 = PyArray_DATA(ap1);            /* x[0]             */
-    is1 = PyArray_STRIDES(ap1)[0];      /* x strides        */
-    ip2 = PyArray_DATA(ap2);            /* y[0]             */
-    is2 = PyArray_STRIDES(ap2)[0];      /* y strides        */
-    op = PyArray_DATA(ret);             /* answer data      */
-    os = PyArray_DESCR(ret)->elsize;    /* answer strides   */
+    /* p is pointer to array start and s is stride */
+    /* i1 is x, i2, is y, o is answer */
+    ip1 = PyArray_DATA(ap1);
+    is1 = PyArray_STRIDES(ap1)[0];
+    ip2 = PyArray_DATA(ap2);
+    is2 = PyArray_STRIDES(ap2)[0];
+    op = PyArray_DATA(ret);
+    os = PyArray_DESCR(ret)->elsize;
 
     lag = minlag;
     if (lag < -n2+1) {
@@ -1208,14 +1211,19 @@ _pyarray_correlate(PyArrayObject *ap1, PyArrayObject *ap2, int typenum,
     maxleft = (0 < maxlag ? 0 : maxlag);
     tmplag = lag;
     for (lag = tmplag; lag < maxleft; lag+=lagstep) {
-        /* for lags where y is left of x */
-        n = n2 + lag;   /* overlap is length of y - lag */
+        /* for lags where y is left of x,
+         * overlap is length of y - lag
+         */
+        n = n2 + lag;
         dot(ip1, is1, ip2 - lag*is2, is2, op, n, ret);
-        op += os;       /* iterate over answer vector   */
+        /* iterate over answer vector */
+        op += os;
     }
     if (maxlag < n1 - n2) {
-        /* if maxlag doesn't take y all the way to the end of x */
-        n11 = maxlag + n2;      /* relevant length of x is smaller */
+        /* if maxlag doesn't take y all the way to the end of x,
+         * relevant length of x is smaller
+         */
+        n11 = maxlag + n2;
     }
     else {
         n11 = n1;
@@ -1308,11 +1316,13 @@ _pyarray_revert(PyArrayObject *ret)
 /*NUMPY_API
  * correlate(a1,a2,mode,minlag,maxlag,lagstep)
  *
- * This function correlates a1 with a2.  mode=0,1,2 retains functionality of PyArray_Correlate2
+ * This function correlates a1 with a2.
+ * mode=0,1,2 retains functionality of PyArray_Correlate2.
  * mode=3 allows for specification of minlag, maxlag, and lagstep.
  */
 NPY_NO_EXPORT PyObject *
-PyArray_CorrelateLags(PyObject *op1, PyObject *op2, int mode, npy_intp minlag, npy_intp maxlag, npy_intp lagstep)
+PyArray_CorrelateLags(PyObject *op1, PyObject *op2, int mode,
+                      npy_intp minlag, npy_intp maxlag, npy_intp lagstep)
 {
     PyArrayObject *ap1, *ap2, *ret = NULL;
     int typenum;
@@ -1347,14 +1357,15 @@ PyArray_CorrelateLags(PyObject *op1, PyObject *op2, int mode, npy_intp minlag, n
         ap2 = cap2;
     }
 
-    ret = _pyarray_correlate(ap1, ap2, typenum, mode, &inverted, minlag, maxlag, lagstep);
+    ret = _pyarray_correlate(ap1, ap2, typenum, mode, &inverted,
+                             minlag, maxlag, lagstep);
     if (ret == NULL) {
         goto clean_ap2;
     }
 
     /*
-     * If we inverted input orders, we need to reverse the output array (i.e.
-     * ret = ret[::-1])
+     * If we inverted input orders, we need to reverse the output array
+     * (i.e. ret = ret[::-1])
      */
     if (inverted) {
         st = _pyarray_revert(ret);
@@ -1387,7 +1398,9 @@ PyArray_Correlate2(PyObject *op1, PyObject *op2, int mode)
 {
     npy_intp z = 0;
     /* For modes other than 3 (which was not a mode before maxlags were added),
-    minlag, maxlag, and lagstep (the last three arguments) will be ignored */
+     * minlag, maxlag, and lagstep (the last three args of _pyarray_correlate)
+     * will be ignored
+     */
     return PyArray_CorrelateLags(op1, op2, mode, z, z, z);
 }
 
@@ -1404,7 +1417,9 @@ PyArray_Correlate(PyObject *op1, PyObject *op2, int mode)
 
     npy_intp z = 0;
     /* For modes other than 3 (which was not a mode before maxlags were added),
-    minlag, maxlag, and lagstep (the last three arguments of _pyarray_correlate) will be ignored */
+     * minlag, maxlag, and lagstep (the last three args of _pyarray_correlate)
+     * will be ignored
+     */
 
     typenum = PyArray_ObjectType(op1, 0);
     typenum = PyArray_ObjectType(op2, typenum);
@@ -2914,7 +2929,8 @@ array_correlate2(PyObject *NPY_UNUSED(dummy), PyObject *args, PyObject *kwds)
     PyObject *shape, *a0;
     int mode = 0;
     npy_intp maxlag = 0, minlag = 0, lagstep = 0;
-    static char *kwlist[] = {"a", "v", "mode", "minlag", "maxlag", "lagstep", NULL};
+    static char *kwlist[] = {"a", "v", "mode", "minlag", "maxlag", "lagstep",
+                             NULL};
 
     if (!PyArg_ParseTupleAndKeywords(args, kwds, "OO|innn", kwlist,
                 &a0, &shape, &mode, &minlag, &maxlag, &lagstep)) {

From 03c6a341e3be1eeb5edaaa8cbeef7891690cc2c2 Mon Sep 17 00:00:00 2001
From: Honi Sanders <honi@brandeis.edu>
Date: Tue, 21 Jun 2016 17:39:32 -0400
Subject: [PATCH 09/12] now no mode argument in PyArray_CorrelateLags() or
 _pyarray_correlate() in multiarraymodule.c, created _lags_from_mode() to
 compensate

---
 doc/source/reference/c-api.array.rst         |  26 +--
 numpy/core/src/multiarray/multiarraymodule.c | 206 ++++++++++++++-----
 2 files changed, 169 insertions(+), 63 deletions(-)

diff --git a/doc/source/reference/c-api.array.rst b/doc/source/reference/c-api.array.rst
index 4aecd010db35..fc3845f87cd3 100644
--- a/doc/source/reference/c-api.array.rst
+++ b/doc/source/reference/c-api.array.rst
@@ -2161,32 +2161,28 @@ Array Functions
 
       z[k] = sum_n op1[n] * conj(op2[n+k])
 
-.. c:function:: PyObject* PyArray_CorrelateLags(PyObject* op1, PyObject* op2, int mode, npy_intp minlag, npy_intp maxlag, npy_intp lagstep)
+.. c:function:: PyObject* PyArray_CorrelateLags(PyObject* op1, PyObject* op2, npy_intp minlag, npy_intp maxlag, npy_intp lagstep)
 
     Updated version of PyArray_Correlate2, which allows for custom lags in calculation.
     The correlation is computed at each output point by multiplying *op1* by
     a shifted version of *op2* and summing the result.
     As a result of the shift, needed values outside of the defined range of
-    *op1* and *op2* are interpreted as zero. The mode determines how many
-    shifts to return: 0 - return only shifts that did not need to assume zero-
-    values; 1 - return an object that is the same size as *op1*, 2 - return all
-    possible shifts (any overlap at all is accepted).
+    *op1* and *op2* are interpreted as zero.
 
-    This version differs from the other PyArray_Correlate* functions in the addition of a
-    third mode, mode = 3.  In this mode, PyArray_CorrelateLags calculates the correlation
-    for lags starting at minlag and ending at (maxlag - 1), with steps of size lagstep.
+    This version differs from the other PyArray_Correlate* functions in that
+    PyArray_CorrelateLags calculates the correlation for lags starting at minlag
+    and ending at (maxlag - 1), with steps of size lagstep.
     The lags for which the correlation will be calculated correspond with the values in
     the vector formed by numpy.arange((minlag, maxlag, lagstep)).
 
     .. rubric:: Notes
 
-    If the mode argument is 0, 1, or 2, the values of minlag, maxlag, and lagstep will be
-    ignored.  Instead, the following lag values will be used, corresponding to the
-    description of modes above (n1 and n2 are the vector lengths):
-       mode   |   minlag   |   maxlag   |   lagstep
-         0    |     0      |  n1-n2+1   |      1
-         1    |   -n2/2    |  n1-n2/2   |      1
-         2    |   -n2+1    |     n1     |      1
+    The following lag values should be used to reproduce the mode parameter used by
+    PyArray_Correlate2 (n1 and n2 are the vector lengths):
+       mode     |   minlag   |   maxlag   |   lagstep
+    0 ("valid") |     0      |  n1-n2+1   |      1
+     1 ("same") |   -n2/2    |  n1-n2/2   |      1
+     2 ("full") |   -n2+1    |     n1     |      1
 
     PyArray_CorrelateLags uses the same definition of correlation as PyArray_Correlate2
     in that the conjugate is taken and the argument order matters.
diff --git a/numpy/core/src/multiarray/multiarraymodule.c b/numpy/core/src/multiarray/multiarraymodule.c
index 9fef139d65bc..0aebf565b187 100644
--- a/numpy/core/src/multiarray/multiarraymodule.c
+++ b/numpy/core/src/multiarray/multiarraymodule.c
@@ -1093,8 +1093,8 @@ PyArray_CopyAndTranspose(PyObject *op)
  * inverted is set to 1 if computed correlate(ap2, ap1), 0 otherwise
  */
 static PyArrayObject*
-_pyarray_correlate(PyArrayObject *ap1, PyArrayObject *ap2, int typenum,
-                   int mode, int *inverted,
+_pyarray_correlate(PyArrayObject *ap1, PyArrayObject *ap2,
+                   int typenum, int *inverted,
                    npy_intp minlag, npy_intp maxlag, npy_intp lagstep)
 {
     PyArrayObject *ret;
@@ -1107,7 +1107,7 @@ _pyarray_correlate(PyArrayObject *ap1, PyArrayObject *ap2, int typenum,
 
     NPY_BEGIN_THREADS_DEF;
 
-    if (lagstep == 0 || mode != 3) {
+    if (lagstep == 0) {
         lagstep = 1;
     }
 
@@ -1127,32 +1127,6 @@ _pyarray_correlate(PyArrayObject *ap1, PyArrayObject *ap2, int typenum,
         lagstep = -lagstep;
     }
 
-    switch(mode) {
-    case 0:
-        /* mode = 'valid' */
-        minlag = 0;
-        maxlag = n1 - n2 + 1;
-        break;
-    case 1:
-        /* mode = 'same' */
-        minlag = -(npy_intp)(n2/2);
-        maxlag = n1 + minlag;
-        break;
-    case 2:
-        /* mode = 'full' */
-        minlag = -n2 + 1;
-        maxlag = n1;
-        break;
-    case 3:
-        /* mode = 'maxlag' */
-        /* use minlag, maxlag, and lagstep that were passed in */
-        break;
-    default:
-        PyErr_SetString(PyExc_ValueError, "mode must be 0, 1, 2, or 3");
-        return NULL;
-    }
-
-
     if (lagstep < 0) {
         *inverted = 1;
         i = minlag;
@@ -1313,15 +1287,71 @@ _pyarray_revert(PyArrayObject *ret)
     return 0;
 }
 
+/*
+ * Generate the lags corresponding to each mode
+ * n1 and n2 are sizes of two vectors
+ * minlag, maxlag, and lagstep are edited in-place
+ */
+void
+_lags_from_mode(int mode, npy_intp n1, npy_intp n2,
+                npy_intp *minlag, npy_intp *maxlag, npy_intp *lagstep) {
+
+    int i, i1, inverted = 0;
+    /* place holders for minlag, maxlag, and lagstep respectively */
+    npy_intp m0, m1, s;
+
+
+
+    if (n1 < n2) {
+        i = n1;
+        n1 = n2;
+        n2 = i;
+        inverted = 1;
+    }
+
+    s = 1;
+    switch(mode) {
+    case 0:
+        /* mode = 'valid' */
+        m0 = 0;
+        m1 = n1 - n2 + 1;
+        break;
+    case 1:
+        /* mode = 'same' */
+        m0 = -(npy_intp)(n2/2);
+        m1 = n1 + m0;
+        break;
+    case 2:
+        /* mode = 'full' */
+        m0 = -n2 + 1;
+        m1 = n1;
+        break;
+    default:
+        PyErr_SetString(PyExc_ValueError, "mode must be 0, 1, or 2");
+    }
+
+    if (inverted) {
+        i = m0;
+        i1 = (npy_intp)(npy_ceil((m1 - m0)/(float)s))*s;
+        m0 =  -(i1 + m0 - s);
+        m1 = -(i - s);
+    }
+
+    /* set minlag, maxlag, and lagstep values */
+    *minlag = m0;
+    *maxlag = m1;
+    *lagstep = s;
+}
+
 /*NUMPY_API
- * correlate(a1,a2,mode,minlag,maxlag,lagstep)
+ * correlate(a1,a2,minlag,maxlag,lagstep)
  *
  * This function correlates a1 with a2.
- * mode=0,1,2 retains functionality of PyArray_Correlate2.
- * mode=3 allows for specification of minlag, maxlag, and lagstep.
+ * This function does not accept modes.
+ * See _lags_from_mode() to convert modes to lags.
  */
 NPY_NO_EXPORT PyObject *
-PyArray_CorrelateLags(PyObject *op1, PyObject *op2, int mode,
+PyArray_CorrelateLags(PyObject *op1, PyObject *op2,
                       npy_intp minlag, npy_intp maxlag, npy_intp lagstep)
 {
     PyArrayObject *ap1, *ap2, *ret = NULL;
@@ -1357,7 +1387,7 @@ PyArray_CorrelateLags(PyObject *op1, PyObject *op2, int mode,
         ap2 = cap2;
     }
 
-    ret = _pyarray_correlate(ap1, ap2, typenum, mode, &inverted,
+    ret = _pyarray_correlate(ap1, ap2, typenum, &inverted,
                              minlag, maxlag, lagstep);
     if (ret == NULL) {
         goto clean_ap2;
@@ -1396,12 +1426,39 @@ PyArray_CorrelateLags(PyObject *op1, PyObject *op2, int mode,
 NPY_NO_EXPORT PyObject *
 PyArray_Correlate2(PyObject *op1, PyObject *op2, int mode)
 {
-    npy_intp z = 0;
-    /* For modes other than 3 (which was not a mode before maxlags were added),
-     * minlag, maxlag, and lagstep (the last three args of _pyarray_correlate)
-     * will be ignored
-     */
-    return PyArray_CorrelateLags(op1, op2, mode, z, z, z);
+    PyArrayObject *ap1, *ap2;
+    int typenum;
+    PyArray_Descr *typec;
+    npy_intp minlag, maxlag, lagstep;
+    npy_intp n1, n2;
+
+    typenum = PyArray_ObjectType(op1, 0);
+    typenum = PyArray_ObjectType(op2, typenum);
+
+    typec = PyArray_DescrFromType(typenum);
+    Py_INCREF(typec);
+    ap1 = (PyArrayObject *)PyArray_FromAny(op1, typec, 1, 1,
+                                            NPY_ARRAY_DEFAULT, NULL);
+    if (ap1 == NULL) {
+        Py_DECREF(typec);
+        return NULL;
+    }
+    ap2 = (PyArrayObject *)PyArray_FromAny(op2, typec, 1, 1,
+                                           NPY_ARRAY_DEFAULT, NULL);
+    if (ap2 == NULL) {
+        goto fail;
+    }
+
+    /* size of x (n1) and y (n2) */
+    n1 = PyArray_DIMS(ap1)[0];
+    n2 = PyArray_DIMS(ap2)[0];
+    _lags_from_mode(mode, n1, n2, &minlag, &maxlag, &lagstep);
+    return PyArray_CorrelateLags(op1, op2, minlag, maxlag, lagstep);
+
+fail:
+    Py_XDECREF(ap1);
+    Py_XDECREF(ap2);
+    return NULL;
 }
 
 /*NUMPY_API
@@ -1414,12 +1471,8 @@ PyArray_Correlate(PyObject *op1, PyObject *op2, int mode)
     int typenum;
     int unused;
     PyArray_Descr *typec;
-
-    npy_intp z = 0;
-    /* For modes other than 3 (which was not a mode before maxlags were added),
-     * minlag, maxlag, and lagstep (the last three args of _pyarray_correlate)
-     * will be ignored
-     */
+    npy_intp minlag, maxlag, lagstep;
+    npy_intp n1, n2;
 
     typenum = PyArray_ObjectType(op1, 0);
     typenum = PyArray_ObjectType(op2, typenum);
@@ -1438,7 +1491,12 @@ PyArray_Correlate(PyObject *op1, PyObject *op2, int mode)
         goto fail;
     }
 
-    ret = _pyarray_correlate(ap1, ap2, typenum, mode, &unused, z, z, z);
+    /* size of x (n1) and y (n2) */
+    n1 = PyArray_DIMS(ap1)[0];
+    n2 = PyArray_DIMS(ap2)[0];
+    _lags_from_mode(mode, n1, n2, &minlag, &maxlag, &lagstep);
+    ret = _pyarray_correlate(ap1, ap2, typenum, &unused,
+                             minlag, maxlag, lagstep);
     if(ret == NULL) {
         goto fail;
     }
@@ -2936,7 +2994,59 @@ array_correlate2(PyObject *NPY_UNUSED(dummy), PyObject *args, PyObject *kwds)
                 &a0, &shape, &mode, &minlag, &maxlag, &lagstep)) {
         return NULL;
     }
-    return PyArray_CorrelateLags(a0, shape, mode, minlag, maxlag, lagstep);
+    if (mode == -1) {
+        if (minlag == 0 && maxlag == 0 && lagstep == 0) {
+            /* if no lag parameters passed, use default: mode = 'valid' */
+            mode = 0;
+        }
+        else {
+            /* if lag parameters were passed, use them */
+            mode = 3;
+        }
+    }
+    if (mode != 3) {
+        PyArrayObject *ap1, *ap2;
+        int typenum;
+        PyArray_Descr *typec;
+        npy_intp n1, n2;
+        npy_intp m0, m1, s;
+
+        typenum = PyArray_ObjectType(a0, 0);
+        typenum = PyArray_ObjectType(shape, typenum);
+
+        typec = PyArray_DescrFromType(typenum);
+        Py_INCREF(typec);
+        ap1 = (PyArrayObject *)PyArray_FromAny(a0, typec, 1, 1,
+                                                NPY_ARRAY_DEFAULT, NULL);
+        if (ap1 == NULL) {
+            Py_DECREF(typec);
+            return NULL;
+        }
+        ap2 = (PyArrayObject *)PyArray_FromAny(shape, typec, 1, 1,
+                                               NPY_ARRAY_DEFAULT, NULL);
+        if (ap2 == NULL) {
+            Py_XDECREF(ap1);
+            Py_XDECREF(ap2);
+            return NULL;
+        }
+        /* size of x (n1) and y (n2) */
+        n1 = PyArray_DIMS(ap1)[0];
+        n2 = PyArray_DIMS(ap2)[0];
+        _lags_from_mode(mode, n1, n2, &m0, &m1, &s);
+        if (minlag != 0) {
+            assert(m0 == minlag);
+        }
+        if (maxlag != 0) {
+            assert(m1 == maxlag);
+        }
+        if (lagstep != 0) {
+            assert(s == lagstep);
+        }
+        minlag = m0;
+        maxlag = m1;
+        lagstep = s;
+    }
+    return PyArray_CorrelateLags(a0, shape, minlag, maxlag, lagstep);
 }
 
 static PyObject *

From 2f96454e78ef93210df543ca3c2dd98e25d91966 Mon Sep 17 00:00:00 2001
From: Honi Sanders <honi@brandeis.edu>
Date: Tue, 21 Jun 2016 17:40:14 -0400
Subject: [PATCH 10/12] added array_correlate() and multiarray.correlate() back
 to multiarraymodule.c

---
 numpy/core/src/multiarray/multiarraymodule.c | 21 ++++++++++++++++++--
 1 file changed, 19 insertions(+), 2 deletions(-)

diff --git a/numpy/core/src/multiarray/multiarraymodule.c b/numpy/core/src/multiarray/multiarraymodule.c
index 0aebf565b187..41c3ee9e4333 100644
--- a/numpy/core/src/multiarray/multiarraymodule.c
+++ b/numpy/core/src/multiarray/multiarraymodule.c
@@ -2982,11 +2982,25 @@ array_fastCopyAndTranspose(PyObject *NPY_UNUSED(dummy), PyObject *args)
 }
 
 static PyObject *
-array_correlate2(PyObject *NPY_UNUSED(dummy), PyObject *args, PyObject *kwds)
+array_correlate(PyObject *NPY_UNUSED(dummy), PyObject *args, PyObject *kwds)
 {
     PyObject *shape, *a0;
     int mode = 0;
-    npy_intp maxlag = 0, minlag = 0, lagstep = 0;
+    static char *kwlist[] = {"a", "v", "mode", NULL};
+
+    if (!PyArg_ParseTupleAndKeywords(args, kwds, "OO|i", kwlist,
+                &a0, &shape, &mode)) {
+        return NULL;
+    }
+    return PyArray_Correlate(a0, shape, mode);
+}
+
+static PyObject *
+array_correlate2(PyObject *NPY_UNUSED(dummy), PyObject *args, PyObject *kwds)
+{
+    PyObject *shape, *a0;
+    int mode = -1;
+    npy_intp minlag = 0, maxlag = 0, lagstep = 0;
     static char *kwlist[] = {"a", "v", "mode", "minlag", "maxlag", "lagstep",
                              NULL};
 
@@ -4344,6 +4358,9 @@ static struct PyMethodDef array_module_methods[] = {
     {"_fastCopyAndTranspose",
         (PyCFunction)array_fastCopyAndTranspose,
         METH_VARARGS, NULL},
+    {"correlate",
+        (PyCFunction)array_correlate,
+        METH_VARARGS | METH_KEYWORDS, NULL},
     {"correlate2",
         (PyCFunction)array_correlate2,
         METH_VARARGS | METH_KEYWORDS, NULL},

From bb87262fa85d4285ebf6131c7d3405520ef83eb1 Mon Sep 17 00:00:00 2001
From: Honi Sanders <honi@brandeis.edu>
Date: Tue, 21 Jun 2016 19:09:13 -0400
Subject: [PATCH 11/12] now correlate() and convolve() in numeric.py have a
 'lags' mode, along with a new lags argument.

---
 numpy/core/numeric.py            | 158 +++++++++++++++++++------------
 numpy/core/tests/test_numeric.py |  57 +++++++----
 2 files changed, 136 insertions(+), 79 deletions(-)

diff --git a/numpy/core/numeric.py b/numpy/core/numeric.py
index 647b90b3ed2e..886179534534 100644
--- a/numpy/core/numeric.py
+++ b/numpy/core/numeric.py
@@ -815,43 +815,56 @@ def flatnonzero(a):
 
 _mode_from_name_dict = {'v': 0,
                         's': 1,
-                        'f': 2}
-
+                        'f': 2,
+                        'l': 3,
+                        }
 
 def _mode_from_name(mode):
-    if isinstance(mode, basestring):
+    if type(mode) is int:
+        assert(mode in _mode_from_name_dict.values())
+        return mode
+    elif isinstance(mode, basestring):
         return _mode_from_name_dict[mode.lower()[0]]
     return mode
 
-
-def _lags_from_mode(alen, vlen, mode):
-    if type(mode) is int:       # maxlag
-        lags = (-mode+1, mode, 1)
-        mode = 3
-    elif type(mode) is tuple:     # minlag and maxlag
-        if len(mode) > 2:
-            lags = (int(mode[0]), int(mode[1]), int(mode[2]))
-        else:
-            lags = (int(mode[0]), int(mode[1]), 1)
-        mode = 3
-    elif isinstance(mode, basestring):
-        mode = _mode_from_name_dict[mode.lower()[0]]
-        if alen < vlen:
-            alen, vlen = vlen, alen
-            inverted = 1
+def _lags_from_lags(l):
+    if type(l) is int:          # maxlag
+        lags = (-l+1, l, 1)
+    elif type(l) is tuple:      # minlag and maxlag
+        if len(l) > 2:
+            lags = (int(l[0]), int(l[1]), int(l[2]))
         else:
-            inverted = 0
-        if mode is 0:
-            lags = (0, alen-vlen+1, 1)
-        elif mode is 1:
-            lags = (-int(vlen/2), alen - int(vlen/2), 1)
-        elif mode is 2:
-            lags = (-vlen+1, alen, 1)
-        if inverted:
-            lags = (-int(ceil((lags[1]-lags[0])/float(lags[2])))*lags[2]-lags[0]+lags[2], -lags[0]+lags[2], lags[2])
-    return mode, lags
-
-def correlate(a, v, mode='valid', returns_lags=False):
+            lags = (int(l[0]), int(l[1]), 1)
+    return lags
+
+def _lags_from_mode(alen, vlen, mode, lags=()):
+    inverted = 0
+    if alen < vlen:
+        alen, vlen = vlen, alen
+        inverted = 1
+
+    mode = _mode_from_name(mode)
+    if mode is 0:
+        mode_lags = (0, alen-vlen+1, 1)
+    elif mode is 1:
+        mode_lags = (-int(vlen/2), alen - int(vlen/2), 1)
+    elif mode is 2:
+        mode_lags = (-vlen+1, alen, 1)
+    elif mode is 3:
+        mode_lags = _lags_from_lags(lags)
+
+    if inverted and mode != 3:
+        mode_lags = (-int(ceil((mode_lags[1]-mode_lags[0])/float(mode_lags[2])))
+                     *mode_lags[2]-mode_lags[0]+mode_lags[2],
+                     -mode_lags[0]+mode_lags[2], mode_lags[2])
+
+    if lags:
+        lags = _lags_from_lags(lags)
+        assert lags == mode_lags, \
+               'Custom lags should be passed to correlate() with mode="lags".'
+    return mode, mode_lags
+
+def correlate(a, v, mode='default', lags=(), returns_lagvector=False):
     """
     Cross-correlation of two 1-dimensional sequences.
 
@@ -867,10 +880,13 @@ def correlate(a, v, mode='valid', returns_lags=False):
     ----------
     a, v : array_like
         Input sequences.
-    mode : int, int tuple, or {'valid', 'same', 'full'}, optional
+    mode : {'valid', 'same', 'full', 'lags'}, optional
         Refer to the `convolve` docstring.  Note that the default
         is `valid`, unlike `convolve`, which uses `full`.
-    returns_lags : bool, optional
+    lags : int or int tuple, optional
+        Refer to the `convolve` docstring.
+        mode should be unset or set to 'lags' to use the lags argument
+    returns_lagvector : bool, optional
         If True, the function returns a lagvector array in addition to the
         cross-correlation result.  The lagvector contains the indices of
         the lags for which the cross-correlation was calculated.  It is
@@ -908,11 +924,11 @@ def correlate(a, v, mode='valid', returns_lags=False):
     array([ 3.5])
     >>> np.correlate([1, 2, 3], [0, 1, 0.5], mode="same")
     array([ 2. ,  3.5,  3. ])
-    >>> np.correlate([1, 2, 3], [0, 1, 0.5], mode="full", returns_lags=True)
+    >>> np.correlate([1, 2, 3], [0, 1, 0.5], mode="full", returns_lagvector=True)
     (array([ 0.5,  2. ,  3.5,  3. ,  0. ]), array([-2, -1,  0,  1,  2]))
-    >>> np.correlate([1, 2, 3], [0, 1, 0.5], mode=2)
+    >>> np.correlate([1, 2, 3], [0, 1, 0.5], mode="lags", lags=2)
     array([ 2. ,  3.5,  3. ])
-    >>> np.correlate([1, 2, 3], [0, 1, 0.5], mode=(-1,2,2), returns_lags=True)
+    >>> np.correlate([1, 2, 3], [0, 1, 0.5], mode="lags", lags=(-1,2,2), returns_lagvector=True)
     (array([ 2.,  3.]), array([-1,  1]))
 
     Using complex sequences:
@@ -928,14 +944,24 @@ def correlate(a, v, mode='valid', returns_lags=False):
     array([ 0.0+0.j ,  3.0+1.j ,  1.5+1.5j,  1.0+0.j ,  0.5+0.5j])
 
     """
-    mode, lags = _lags_from_mode(len(a), len(v), mode)
-    if returns_lags:
-        return multiarray.correlate2(a, v, 3, lags[0], lags[1], lags[2]), arange(lags[0], lags[1], lags[2])
+    if len(a) == 0:
+        raise ValueError('a cannot be empty')
+    if len(v) == 0:
+        raise ValueError('v cannot be empty')
+    if mode == 'default':
+        if lags:
+            mode = 'lags'
+        else:
+            mode = 'valid'
+    mode, lags = _lags_from_mode(len(a), len(v), mode, lags)
+    if returns_lagvector:
+        return multiarray.correlate2(a, v, 3, lags[0], lags[1], lags[2]), \
+                                     arange(lags[0], lags[1], lags[2])
     else:
         return multiarray.correlate2(a, v, 3, lags[0], lags[1], lags[2])
 
 
-def convolve(a, v, mode='full', returns_lags=False):
+def convolve(a, v, mode='full', lags=(), returns_lagvector=False):
     """
     Returns the discrete, linear convolution of two one-dimensional sequences.
 
@@ -954,19 +980,6 @@ def convolve(a, v, mode='full', returns_lags=False):
     v : (M,) array_like
         Second one-dimensional input array.
     mode : int, int tuple, or {'full', 'valid', 'same'}, optional
-        int (maxlag):
-          This calculates the convolution for all lags starting at
-          (-maxlag + 1) and ending at (maxlag - 1), with steps of size 1.
-          See the optional lagvec argument to get an array containing
-          lags corresponding to the convolution values in the return array.
-
-        tuple (minlag, maxlag) or (minlag, maxlag, lagstep):
-          This calculates the convolution for all lags starting at
-          minlag and ending at (maxlag - 1), with steps of size lagstep.
-          The lags for which the convolution will be calculated correspond
-          with the values in the vector formed by numpy.arange() with the
-          same tuple argument.
-
         'full':
           By default, mode is 'full'.  This returns the convolution
           at each point of overlap, with an output shape of (N+M-1,). At
@@ -986,7 +999,26 @@ def convolve(a, v, mode='full', returns_lags=False):
           for points where the signals overlap completely.  Values outside
           the signal boundary have no effect. This corresponds with a lag tuple
           of (0, N-M+1, 1) for N>M or (-M+N, 1, 1) for M>N.
-    returns_lags : bool, optional
+
+        'lags':
+          Mode 'lags' uses the lags argument to define the lags for which
+          to perform the convolution.
+    lags : int or int tuple, optional
+        mode should be unset or set to 'lags' to use the lags argument
+
+        int (maxlag):
+          This calculates the convolution for all lags starting at
+          (-maxlag + 1) and ending at (maxlag - 1), with steps of size 1.
+          See the optional returns_lagvec argument to get an array containing
+          lags corresponding to the convolution values in the return array.
+
+        tuple (minlag, maxlag) or (minlag, maxlag, lagstep):
+          This calculates the convolution for all lags starting at
+          minlag and ending at (maxlag - 1), with steps of size lagstep.
+          The lags for which the convolution will be calculated correspond
+          with the values in the vector formed by numpy.arange() with the
+          same tuple argument.
+    returns_lagvector : bool, optional
         If True, the function returns a lagvector array in addition to the
         convolution result.  The lagvector contains the indices of
         the lags for which the convolution was calculated.  It is
@@ -1048,7 +1080,7 @@ def convolve(a, v, mode='full', returns_lags=False):
     to the convolution in addition to the convolution
     itself:
 
-    >>> np.convolve([1,2,3],[0,1,0.5], mode='valid', returns_lags=True)
+    >>> np.convolve([1,2,3],[0,1,0.5], mode='valid', returns_lagvector=True)
     (array([ 2.5]), array([0]))
 
     Find the convolution for lags ranging from -1 to 1
@@ -1057,14 +1089,14 @@ def convolve(a, v, mode='full', returns_lags=False):
     the first, and +1 has the second vector to the right
     of the first):
 
-    >>> np.convolve([1,2,3],[0,1,0.5], mode=2, returns_lags=True)
+    >>> np.convolve([1,2,3],[0,1,0.5], mode='lags', lags=2, returns_lagvector=True)
     (array([ 1. ,  2.5,  4. ]), array([-1,  0,  1]))
 
     Find the convolution for lags ranging from -2 to 4
     with steps of length 2 (the maxlag member of the
     lag range tuple is non-inclusive, similar to np.arange()):
 
-    >>> np.convolve([1,2,3,4,5],[0,1,0.5], mode=(-2,6,2), returns_lags=True)
+    >>> np.convolve([1,2,3,4,5],[0,1,0.5], mode='lags', lags=(-2,6,2), returns_lagvector=True)
     (array([ 0. ,  2.5,  5.5,  2.5]), array([-2,  0,  2,  4]))
 
     """
@@ -1075,9 +1107,15 @@ def convolve(a, v, mode='full', returns_lags=False):
         raise ValueError('a cannot be empty')
     if len(v) == 0:
         raise ValueError('v cannot be empty')
-    mode, lags = _lags_from_mode(len(a), len(v), mode)
-    if returns_lags:
-        return multiarray.correlate2(a, v[::-1], 3, lags[0], lags[1], lags[2]), arange(lags[0], lags[1], lags[2])
+    if mode == 'default':
+        if lags:
+            mode = 'lags'
+        else:
+            mode = 'full'
+    mode, lags = _lags_from_mode(len(a), len(v), mode, lags)
+    if returns_lagvector:
+        return multiarray.correlate2(a, v[::-1], 3, lags[0], lags[1], lags[2]),\
+                                     arange(lags[0], lags[1], lags[2])
     else:
         return multiarray.correlate2(a, v[::-1], 3, lags[0], lags[1], lags[2])
 
diff --git a/numpy/core/tests/test_numeric.py b/numpy/core/tests/test_numeric.py
index 1d0cc60f72d8..e7f4def24939 100644
--- a/numpy/core/tests/test_numeric.py
+++ b/numpy/core/tests/test_numeric.py
@@ -2053,19 +2053,21 @@ def _setup(self, dt):
 
     def test_float(self):
         self._setup(np.float)
-        z, lagvec = np.correlate(self.x, self.y, 'full', returns_lags=True)
+        z, lagvec = np.correlate(self.x, self.y, 'full', returns_lagvector=True)
         assert_array_almost_equal(z, self.z1)
         assert_equal(lagvec.size, z.size)
         assert_array_equal(lagvec, np.arange(-self.y.size + 1, self.x.size))
-        z, lagvec = np.correlate(self.x, self.y, 'same', returns_lags=True)
+        z, lagvec = np.correlate(self.x, self.y, 'same', returns_lagvector=True)
         assert_array_almost_equal(z, self.z1s)
         assert_array_almost_equal(z, self.z1[lagvec+self.y.size-1])
         assert_equal(lagvec.size, z.size)
         if self.x.size > self.y.size:
-            assert_array_equal(lagvec, np.arange(-int(self.y.size/2), self.x.size - int(self.y.size/2)))
+            assert_array_equal(lagvec,
+                np.arange(-int(self.y.size/2), self.x.size - int(self.y.size/2)))
         else:
-            assert_array_equal(lagvec, np.arange(-self.y.size + int(self.x.size/2) + 1, int(self.x.size/2) + 1))
-        z, lagvec = np.correlate(self.x, self.y, 'valid', returns_lags=True)
+            assert_array_equal(lagvec,
+                np.arange(-self.y.size + int(self.x.size/2) + 1, int(self.x.size/2) + 1))
+        z, lagvec = np.correlate(self.x, self.y, 'valid', returns_lagvector=True)
         assert_array_almost_equal(z, self.z1v)
         assert_array_almost_equal(z, self.z1[lagvec+self.y.size-1])
         assert_equal(lagvec.size, z.size)
@@ -2076,19 +2078,21 @@ def test_float(self):
 
         z = np.correlate(self.x, self.y[:-1], 'full')
         assert_array_almost_equal(z, self.z1_4)
-        z, lagvec = np.correlate(self.y, self.x, 'full', returns_lags=True)
+        z, lagvec = np.correlate(self.y, self.x, 'full', returns_lagvector=True)
         assert_array_almost_equal(z, self.z2)
         assert_equal(lagvec.size, z.size)
         assert_array_equal(lagvec, np.arange(-self.x.size + 1, self.y.size))
-        z, lagvec = np.correlate(self.y, self.x, 'same', returns_lags=True)
+        z, lagvec = np.correlate(self.y, self.x, 'same', returns_lagvector=True)
         assert_array_almost_equal(z, self.z2s)
         assert_array_almost_equal(z, self.z2[lagvec+self.x.size-1])
         assert_equal(lagvec.size, z.size)
         if self.y.size > self.x.size:
-            assert_array_equal(lagvec, np.arange(-int(self.x.size/2), self.y.size - int(self.x.size/2)))
+            assert_array_equal(lagvec,
+                np.arange(-int(self.x.size/2), self.y.size - int(self.x.size/2)))
         else:
-            assert_array_equal(lagvec, np.arange(-self.x.size + int(self.y.size/2) + 1, int(self.y.size/2) + 1))
-        z, lagvec = np.correlate(self.y, self.x, 'valid', returns_lags=True)
+            assert_array_equal(lagvec,
+                np.arange(-self.x.size + int(self.y.size/2) + 1, int(self.y.size/2) + 1))
+        z, lagvec = np.correlate(self.y, self.x, 'valid', returns_lagvector=True)
         assert_array_almost_equal(z, self.z2v)
         assert_array_almost_equal(z, self.z2[lagvec+self.x.size-1])
         assert_equal(lagvec.size, z.size)
@@ -2112,49 +2116,64 @@ def test_lags(self):
         tmp = np.correlate(self.x, self.y, 'full')
         z_full = np.zeros(tmp.size + 100)
         z_full[:tmp.size] = tmp
-        z, lagvec = np.correlate(self.x, self.y, maxlag, returns_lags=True)
+        z, lagvec = np.correlate(self.x, self.y, lags=maxlag,
+                                 returns_lagvector=True)
         assert_array_equal(z[0:3], np.zeros(3))
         assert_array_equal(z[-3:], np.zeros(3))
         assert_equal(lagvec.size, z.size)
         assert_array_equal(lagvec, np.arange(-maxlag+1, maxlag))
-        z, lagvec = np.correlate(self.x, self.y, (minlag, maxlag), returns_lags=True)
+        z, lagvec = np.correlate(self.x, self.y, mode='lags',
+                                 lags=(minlag, maxlag), returns_lagvector=True)
         assert_array_almost_equal(z, z_full[lagvec+self.y.size-1])
         assert_equal(lagvec.size, z.size)
         assert_array_equal(lagvec, np.arange(minlag, maxlag))
-        z, lagvec = np.correlate(self.x, self.y, (minlag, maxlag, lagstep), returns_lags=True)
+        z, lagvec = np.correlate(self.x, self.y, lags=(minlag, maxlag, lagstep),
+                                 returns_lagvector=True)
         assert_array_almost_equal(z, z_full[lagvec+self.y.size-1])
         assert_equal(lagvec.size, z.size)
         assert_array_equal(lagvec, np.arange(minlag, maxlag, lagstep))
         lagstep = 3
-        z, lagvec = np.correlate(self.x, self.y, (minlag, maxlag, lagstep), returns_lags=True)
+        z, lagvec = np.correlate(self.x, self.y, mode='lags',
+                                 lags=(minlag, maxlag, lagstep),
+                                 returns_lagvector=True)
         assert_array_almost_equal(z, z_full[lagvec+self.y.size-1])
         assert_equal(lagvec.size, z.size)
         minlag = 10
         maxlag = -2
         lagstep = -2
-        z, lagvec = np.correlate(self.x, self.y, (minlag, maxlag, lagstep), returns_lags=True)
+        z, lagvec = np.correlate(self.x, self.y, mode='lags',
+                                 lags=(minlag, maxlag, lagstep),
+                                 returns_lagvector=True)
         assert_array_almost_equal(z, z_full[lagvec+self.y.size-1])
         assert_equal(lagvec.size, z.size)
         assert_array_equal(lagvec, np.arange(int(minlag), int(maxlag), lagstep))
         maxlag = 10.2
         minlag = -2.2
-        z, lagvec = np.correlate(self.x, self.y, (minlag, maxlag), returns_lags=True)
+        z, lagvec = np.correlate(self.x, self.y, mode='lags',
+                                 lags=(minlag, maxlag),
+                                 returns_lagvector=True)
         assert_array_almost_equal(z, z_full[lagvec+self.y.size-1])
         assert_equal(lagvec.size, z.size)
         assert_array_equal(lagvec, np.arange(int(minlag), int(maxlag)))
-        z, lagvec = np.correlate(self.x, self.y, (maxlag, minlag), returns_lags=True)
+        z, lagvec = np.correlate(self.x, self.y, mode='lags',
+                                 lags=(maxlag, minlag),
+                                 returns_lagvector=True)
         assert_array_almost_equal(z, z_full[lagvec+self.y.size-1])
         assert_equal(lagvec.size, z.size)
         assert_array_equal(lagvec, np.arange(int(maxlag), int(minlag)))
         maxlag = -1
         minlag = -2
-        z, lagvec = np.correlate(self.x, self.y, (minlag, maxlag), returns_lags=True)
+        z, lagvec = np.correlate(self.x, self.y, mode='lags',
+                                 lags=(minlag, maxlag),
+                                 returns_lagvector=True)
         assert_array_almost_equal(z, z_full[lagvec+self.y.size-1])
         assert_equal(lagvec.size, z.size)
         assert_array_equal(lagvec, np.arange(int(minlag), int(maxlag)))
         maxlag = 2
         minlag = 4
-        z, lagvec = np.correlate(self.x, self.y, (minlag, maxlag), returns_lags=True)
+        z, lagvec = np.correlate(self.x, self.y, mode='lags',
+                                 lags=(minlag, maxlag),
+                                 returns_lagvector=True)
         assert_array_almost_equal(z, z_full[lagvec+self.y.size-1])
         assert_equal(lagvec.size, z.size)
         assert_array_equal(lagvec, np.arange(int(minlag), int(maxlag)))

From ecd85dd01a09514d8f39e4485c3e2cf0de4f5776 Mon Sep 17 00:00:00 2001
From: Honi Sanders <honi@brandeis.edu>
Date: Wed, 22 Jun 2016 13:26:43 -0400
Subject: [PATCH 12/12] working on benchmarks

---
 numpy/core/numeric.py                        | 23 ++---
 numpy/core/src/multiarray/multiarraymodule.c | 92 +++++++++-----------
 2 files changed, 56 insertions(+), 59 deletions(-)

diff --git a/numpy/core/numeric.py b/numpy/core/numeric.py
index 886179534534..8fed76d6d6c3 100644
--- a/numpy/core/numeric.py
+++ b/numpy/core/numeric.py
@@ -843,7 +843,6 @@ def _lags_from_mode(alen, vlen, mode, lags=()):
         alen, vlen = vlen, alen
         inverted = 1
 
-    mode = _mode_from_name(mode)
     if mode is 0:
         mode_lags = (0, alen-vlen+1, 1)
     elif mode is 1:
@@ -862,7 +861,7 @@ def _lags_from_mode(alen, vlen, mode, lags=()):
         lags = _lags_from_lags(lags)
         assert lags == mode_lags, \
                'Custom lags should be passed to correlate() with mode="lags".'
-    return mode, mode_lags
+    return mode_lags
 
 def correlate(a, v, mode='default', lags=(), returns_lagvector=False):
     """
@@ -944,16 +943,16 @@ def correlate(a, v, mode='default', lags=(), returns_lagvector=False):
     array([ 0.0+0.j ,  3.0+1.j ,  1.5+1.5j,  1.0+0.j ,  0.5+0.5j])
 
     """
-    if len(a) == 0:
-        raise ValueError('a cannot be empty')
-    if len(v) == 0:
-        raise ValueError('v cannot be empty')
     if mode == 'default':
         if lags:
             mode = 'lags'
         else:
             mode = 'valid'
-    mode, lags = _lags_from_mode(len(a), len(v), mode, lags)
+    mode = _mode_from_name(mode)
+    if mode in (0, 1, 2) and not returns_lagvector:
+        return multiarray.correlate2(a, v, mode)
+
+    lags = _lags_from_mode(len(a), len(v), mode, lags)
     if returns_lagvector:
         return multiarray.correlate2(a, v, 3, lags[0], lags[1], lags[2]), \
                                      arange(lags[0], lags[1], lags[2])
@@ -1101,18 +1100,22 @@ def convolve(a, v, mode='full', lags=(), returns_lagvector=False):
 
     """
     a, v = array(a, copy=False, ndmin=1), array(v, copy=False, ndmin=1)
-    if (len(v) > len(a)):
-        a, v = v, a
     if len(a) == 0:
         raise ValueError('a cannot be empty')
     if len(v) == 0:
         raise ValueError('v cannot be empty')
+    if (len(v) > len(a)):
+        a, v = v, a
     if mode == 'default':
         if lags:
             mode = 'lags'
         else:
             mode = 'full'
-    mode, lags = _lags_from_mode(len(a), len(v), mode, lags)
+    mode = _mode_from_name(mode)
+    if mode in (0, 1, 2) and not returns_lagvector:
+        return multiarray.correlate2(a, v[::-1], mode)
+
+    lags = _lags_from_mode(len(a), len(v), mode, lags)
     if returns_lagvector:
         return multiarray.correlate2(a, v[::-1], 3, lags[0], lags[1], lags[2]),\
                                      arange(lags[0], lags[1], lags[2])
diff --git a/numpy/core/src/multiarray/multiarraymodule.c b/numpy/core/src/multiarray/multiarraymodule.c
index 41c3ee9e4333..50ed082a39fb 100644
--- a/numpy/core/src/multiarray/multiarraymodule.c
+++ b/numpy/core/src/multiarray/multiarraymodule.c
@@ -1300,8 +1300,6 @@ _lags_from_mode(int mode, npy_intp n1, npy_intp n2,
     /* place holders for minlag, maxlag, and lagstep respectively */
     npy_intp m0, m1, s;
 
-
-
     if (n1 < n2) {
         i = n1;
         n1 = n2;
@@ -1422,15 +1420,18 @@ PyArray_CorrelateLags(PyObject *op1, PyObject *op2,
  *
  * This function computes the usual correlation (correlate(a1, a2) !=
  * correlate(a2, a1), and conjugate the second argument for complex inputs
+ *
+ * This function requires a mode.
  */
 NPY_NO_EXPORT PyObject *
 PyArray_Correlate2(PyObject *op1, PyObject *op2, int mode)
 {
-    PyArrayObject *ap1, *ap2;
+    PyArrayObject *ap1, *ap2, *ret = NULL;
     int typenum;
     PyArray_Descr *typec;
     npy_intp minlag, maxlag, lagstep;
     npy_intp n1, n2;
+    int inverted, st;
 
     typenum = PyArray_ObjectType(op1, 0);
     typenum = PyArray_ObjectType(op2, typenum);
@@ -1438,7 +1439,7 @@ PyArray_Correlate2(PyObject *op1, PyObject *op2, int mode)
     typec = PyArray_DescrFromType(typenum);
     Py_INCREF(typec);
     ap1 = (PyArrayObject *)PyArray_FromAny(op1, typec, 1, 1,
-                                            NPY_ARRAY_DEFAULT, NULL);
+                                           NPY_ARRAY_DEFAULT, NULL);
     if (ap1 == NULL) {
         Py_DECREF(typec);
         return NULL;
@@ -1446,18 +1447,50 @@ PyArray_Correlate2(PyObject *op1, PyObject *op2, int mode)
     ap2 = (PyArrayObject *)PyArray_FromAny(op2, typec, 1, 1,
                                            NPY_ARRAY_DEFAULT, NULL);
     if (ap2 == NULL) {
-        goto fail;
+        goto clean_ap1;
+    }
+
+    if (PyArray_ISCOMPLEX(ap2)) {
+        PyArrayObject *cap2;
+        cap2 = (PyArrayObject *)PyArray_Conjugate(ap2, NULL);
+        if (cap2 == NULL) {
+            goto clean_ap2;
+        }
+        Py_DECREF(ap2);
+        ap2 = cap2;
     }
 
     /* size of x (n1) and y (n2) */
     n1 = PyArray_DIMS(ap1)[0];
     n2 = PyArray_DIMS(ap2)[0];
     _lags_from_mode(mode, n1, n2, &minlag, &maxlag, &lagstep);
-    return PyArray_CorrelateLags(op1, op2, minlag, maxlag, lagstep);
+    ret = _pyarray_correlate(ap1, ap2, typenum, &inverted,
+                             minlag, maxlag, lagstep);
+    if (ret == NULL) {
+        goto clean_ap2;
+    }
 
-fail:
-    Py_XDECREF(ap1);
-    Py_XDECREF(ap2);
+    /*
+     * If we inverted input orders, we need to reverse the output array
+     * (i.e. ret = ret[::-1])
+     */
+    if (inverted) {
+        st = _pyarray_revert(ret);
+        if(st) {
+            goto clean_ret;
+        }
+    }
+
+    Py_DECREF(ap1);
+    Py_DECREF(ap2);
+    return (PyObject *)ret;
+
+clean_ret:
+    Py_DECREF(ret);
+clean_ap2:
+    Py_DECREF(ap2);
+clean_ap1:
+    Py_DECREF(ap1);
     return NULL;
 }
 
@@ -3019,46 +3052,7 @@ array_correlate2(PyObject *NPY_UNUSED(dummy), PyObject *args, PyObject *kwds)
         }
     }
     if (mode != 3) {
-        PyArrayObject *ap1, *ap2;
-        int typenum;
-        PyArray_Descr *typec;
-        npy_intp n1, n2;
-        npy_intp m0, m1, s;
-
-        typenum = PyArray_ObjectType(a0, 0);
-        typenum = PyArray_ObjectType(shape, typenum);
-
-        typec = PyArray_DescrFromType(typenum);
-        Py_INCREF(typec);
-        ap1 = (PyArrayObject *)PyArray_FromAny(a0, typec, 1, 1,
-                                                NPY_ARRAY_DEFAULT, NULL);
-        if (ap1 == NULL) {
-            Py_DECREF(typec);
-            return NULL;
-        }
-        ap2 = (PyArrayObject *)PyArray_FromAny(shape, typec, 1, 1,
-                                               NPY_ARRAY_DEFAULT, NULL);
-        if (ap2 == NULL) {
-            Py_XDECREF(ap1);
-            Py_XDECREF(ap2);
-            return NULL;
-        }
-        /* size of x (n1) and y (n2) */
-        n1 = PyArray_DIMS(ap1)[0];
-        n2 = PyArray_DIMS(ap2)[0];
-        _lags_from_mode(mode, n1, n2, &m0, &m1, &s);
-        if (minlag != 0) {
-            assert(m0 == minlag);
-        }
-        if (maxlag != 0) {
-            assert(m1 == maxlag);
-        }
-        if (lagstep != 0) {
-            assert(s == lagstep);
-        }
-        minlag = m0;
-        maxlag = m1;
-        lagstep = s;
+        return PyArray_Correlate2(a0, shape, mode);
     }
     return PyArray_CorrelateLags(a0, shape, minlag, maxlag, lagstep);
 }