function [y, Ph] = freqres(A1, M, shift)
if nargin < 3
    shift = 0;
end
A1 = A1(:);
Ph = ((0:M - 1) - floor(M/2))*2*pi/M;
N = length(A1);
num = (0:N - 1) - floor(N/2);

phi = exp(1j*Ph.'*num);
A1 = repmat(A1, 1, length(shift));
y = zeros(M, length(shift));
for k = 1:length(shift)
    A1(:, k) = A1(:, k).*exp(1j*num.'*shift(k));
    y(:, k) = sum(phi*diag(A1(:, k)), 2);
end
