The IMSL_CONVOL1D function computes the discrete convolution of two one-dimensional arrays.

This routine requires an IDL Advanced Math and Stats license. For more information, contact your sales or technical support representative.

Let n_x be the length of x, and n_y denote the length of y. If the keyword PERIODIC is set, then n_z = max{n_x, n_y}, otherwise n_z is set to the smallest whole number, n_z ≥ n_x + n_y – 1, of the form:

The arrays x and y are then zero-padded to a length n_z. Then, we compute:

where the index on x is interpreted as a nonnegative number between 0 and n_z – 1.

The technique used to compute the z_i’s is based on the fact that the (complex discrete) Fourier transform maps convolution into multiplication. Thus, the Fourier transform of z is given by:

where the following equation is true:

The technique used here to compute the convolution is to take the discrete Fourier transform of x and y, multiply the results together component-wise, and then take the inverse transform of this product. It is very important to make sure that nz is the product of small primes if PERIODIC is set. If n_z is a product of small primes, then the computational effort will be proportional to n_z log (n_z). If PERIODIC is not set, then nz is chosen to be a product of small primes.

We point out that if x and y are not complex, then no complex transforms of x or y are taken, since a real transforms can simulate the complex transform above. Such a strategy is six times faster and requires less space than when using the complex transform.

Example

---

This example computes simple moving-average digital filter plots of 5-point and 25- point moving-average filters of noisy data. Results are shown in the figures that follow the example.

PRO Convol1d_ex1

  IMSL_RANDOMOPT, SET = 1234579L

  ; Set the random number seed. 

  ny = 100

  t = FINDGEN(ny)/(ny-1)

  y = SIN(2*!PI*t) + .5*IMSL_RANDOM(ny, /Uniform) -.25

  ; Define a 1-period sine wave with added noise.

  win=0

  FOR nfltr = 5, 25, 20 DO BEGIN

    nfltr_str = strcompress(nfltr,/Remove_All)

    fltr = fltarr(nfltr)

    fltr(*) = 1./nfltr

    ; Define the NFLTR-point moving average array.

    z = IMSL_CONVOL1D(fltr, y, /Periodic)

    ; Convolve filter and signal, using keyword Periodic.

    WINDOW, win++

    PLOT, y, LINESTYLE = 1, TITLE = nfltr_str + $

      '-point Moving Average'

    OPLOT, shift(z, -nfltr/2)

  ENDFOR

END

Syntax

---

Result = IMSL_CONVOL1D(X, Y [, DIRECT=value] [, PERIODIC=value])

Return Value

---

A one-dimensional array containing the discrete convolution of x and y.

Arguments

---

X

One-dimensional array.

Y

One-dimensional array.

Keywords

---

DIRECT (optional)

If present and nonzero, causes the computations to be done by the direct method instead of the FFT method regardless of the size of the vectors passed in.

PERIODIC (optional)

If present and nonzero, then a circular convolution is computed.

Version History

---

6.4	Introduced