const zero = 0.T

mixin `_`

gemm(one, kernel_col, input_shifts, zero, result)

# Now remove the unnecessary dimension of the result

result = result.squeeze()

# And remove the extra samples that sometimes are added because `array2shiftColMatrix`

# works with symmetric paddings at the start and end of input_shifts

if target_result_len < result_len:

result = result[_ ..< target_result_len]

proc convolve*[T](

t1, t2: Tensor[T],

mode = ConvolveMode.full,

down = 1): Tensor[T] {.noinit.} =

## Returns the discrete, linear convolution of two one-dimensional tensors.

##

## The convolution operator is often seen in signal processing, where it models

## the effect of a linear time-invariant system on a signal (Wikipedia,

## “Convolution”, https://en.wikipedia.org/wiki/Convolution).

##

## The convolution is defined as the integral of the product of the two tensors

## after one is reflected about the y-axis and shifted n positions, for all values

## of n in which the tensors overlap (since the integral will be zero outside of

## that window).

##

## Inputs:

## - t1, t2: Rank-1 input tensors of size N and M respectively.

## - mode: Convolution mode (`full`, `same` or `valid`):

## - `full`: This is the default mode. It returns the convolution at

## each point of overlap, with an output length of `N + M - 1`.

## At the end-points of the convolution, the signals do not

## overlap completely, and boundary effects may be seen.

## - `same`: Returns an output of length `max(M, N)`.

## Boundary effects are still visible.

## - `valid`: Returns output of length `max(M, N) - min(M, N) + 1`.

## The convolution is only given for points where the signals

## overlap completely. Values outside the signal boundary

## have no effect.

## - down: Downsample ratio applied to the result. Defaults to 1 (i.e.

## no downsampling).

##

## Result:

## - Convolution tensor of the same type as the inputs and size according

## to the mode and the selected `down` downsample ratio.

##

## Notes:

## - The API of this function is based on `numpy.convolve`, with the

## addtion of the `down` argument.

## - The `down` argument is useful for certain signal processing tasks,

## and is more efficient than applying the downsampling after the

## convolution step (which requires `down` times more operations).

# Ensure that both arrays are 1-dimensional

let f = if t1.rank > 1: t1.squeeze else: t1

let g = if t2.rank > 1: t2.squeeze else: t2

if f.rank > 1:

raise newException(ValueError,

"convolve input tensors must be 1D, but first input tensor is multi-dimensional (shape=" & $t1.shape & ")")

if g.rank > 1:

raise newException(ValueError,

"convolve input tensors must be 1D, but second input tensor is multi-dimensional (shape=" & $t2.shape & ")")

mixin `|-`

mixin `_`

correlateImpl(f, g[_|-1], mode = mode, stride = down)

type CorrelateMode* = ConvolveMode

proc correlate*[T: SomeNumber | Complex32 | Complex64](

t1, t2: Tensor[T],

mode = CorrelateMode.valid,

down = 1): Tensor[T] {.noinit.} =

## Returns the cross-correlation of two one-dimensional real tensors.

##

## The correlation is defined as the integral of the product of the two

## tensors after the second one is shifted n positions, for all values of n

## in which the tensors overlap (since the integral will be zero outside of

## that window).

##

## Inputs:

## - t1, t2: Rank-1 input tensors of size N and M respectively.

## - mode: Correlation mode (`full`, `same` or `valid`):

## - `full`: It returns the correlation at each point

## of overlap, with an output length of `N + M - 1`.

## At the end-points of the correlation, the signals do not

## overlap completely, and boundary effects may be seen.

## - `same`: Returns an output of length `max(M, N)`.

## Boundary effects are still visible.

## - `valid`: This is the default mode. Returns output of length

## `max(M, N) - min(M, N) + 1`.

## The correlation is only given for points where the signals

## overlap completely. Values outside the signal boundary

## have no effect.

## - down: Downsample ratio applied to the result. Defaults to 1 (i.e.

## no downsampling).

##

## Result:

## - Correlation tensor of the same type as the inputs and size according

## to the mode and the selected `down` downsample ratio.

##

## Notes:

## - Note that (as with np.correlate) the default correlation mode is

## `valid`, which is different than the default convolution mode (`full`).

## - The API of this function is based on `numpy.convolve`, with the

## addtion of the `down` argument.

## - The `down` argument is useful for certain signal processing tasks,

## and is more efficient than applying the downsampling after the

## correlation step (which requires `down` times more operations).

# Ensure that both arrays are 1-dimensional

let f = if t1.rank > 1: t1.squeeze else: t1

let g = if t2.rank > 1: t2.squeeze else: t2

if f.rank > 1:

raise newException(ValueError,

"correlate input tensors must be 1D, but first input tensor is multi-dimensional (shape=" & $t1.shape & ")")

if g.rank > 1:

raise newException(ValueError,

"correlate input tensors must be 1D, but second input tensor is multi-dimensional (shape=" & $t2.shape & ")")

when T is Complex:

correlateImpl(f, g.conjugate, mode=mode, stride = down)

else:

correlateImpl(f, g, mode=mode, stride = down)