Convolution y(t) = a(t − s) · x(s) ds is one of the main techniques in digital signal processing. A straightforward computation of the convolution y(t) requires O(n2) steps, where n is the number of observations x(t0), . . . , x(tn−1). It is well known that by using the Fast Fourier Transform (FFT) algorithm, we can compute convolution much faster, with computation time O(n · log(n)). In practice, we only know the signal x(t) and the function a(t) with uncertainty. Sometimes, we know them with interval uncertainty, i.e., we know intervals [x(t), x(t)] and [a(t), a(t)] that contain the actual (unknown) functions x(t) and a(t). In such situations, it is desirable, for every t, to compute the range [y(t), y(t)] of possible values of y(t). Of course, it is possible to use straightforward interval computations to compute this range, i.e., replace every computational step in FFT by the corresponding operations of interval arithmetic. However, the resulting enPreprint submitted to Else...