Phase spectrum of the DTFT of a rect pulse using FFT

Question

I'm trying to use the FFT to plot the phase spectrum of the DTFT of a sampled rect impulse function. Here's the Matlab code I've written.

fs=10; % sampling frequency Hz
Ts=1/fs;
T=4; % length of time vector in seconds
t=-T/2:Ts:T/2-Ts; %vector of 40 time samples -2:0.1:1.9
y=rectpuls(t); %unit length rect impulse
z=fft(y,10000).*Ts; 
%fft multiplied by Ts to get rid of the 1/Ts factor in the spectral replicas
%of the DTFT of the signal
f=linspace(0, (length(z)-1)/length(z),length(z));
figure
plot(f,abs(z));
grid on

%plot of the phase of z as it is, so with the time sequence
%considered to be shifted by length(y)/2 samples by fft.
figure
plot(f,unwrap(angle(z)),'r');
grid on

%I try to get rid of the linear phase delay component to get the phase
%spectrum of the DTFT of the signal in its frequency period.
%I divide each fft sample by its own exponential delay factor.
delay=length(y)/(2*length(z));
expd=zeros(1,length(z));
for index=1:1:length(z)
expd(index)=exp(-1i*2*pi*(index-1)*delay);
end
z2=z./expd;
figure
plot(f,unwrap(angle(z2))./pi,'r');
grid on

Now the problem is that this phase plot is wrong, the DTFT is the sum of real sinc functions so its phase spectrum con only assume k*pi values. There's still a linear component in the phase plot and the slope is positive so the exponentials that I'm using are over compensanting the delay. I have found out that the problem is that the fft function introduces a linear phase delay as if the sequence was shifted by 19.5 samples instead of 20 samples as it should be in my example. Infact if I set in my code delay=19.5/10000 I get the correct phase plot (or at least one that makes sense :D):

So what's going on here, is it a precision problem? Are there mistakes in my Matlab code and/or my dsp theory formulas? Thank you!!

Unfortunately I can't link more than two images in the post :/

Answer 1

In order for the result of an even number N length FFT to be strictly real, the input vector data vector has to be exactly symmetric around sample N/2, which is the sample point before t=0, given your sample window, which is offset slightly to the left.

Answer 2

I've found out the reason of the problem. The problem is that to get a real FFT (as it should be in my case) the sequence must be symmetric relatively to x[N/2] which in my case would be x[20]. This doesn't happen (notwithstading the fact that the rect signal is even so its countinous Fourier Transform is real) because of how the rectpuls function I used is defined. Rectpuls(-0.5)=1 and rectpuls(0.5)=0 If instead I create a symmetric sequence using the value 1/2 at times t=-0.5s and t=0.5s (also according to the Dirichlet conditions), I get a real FFT and the results I wanted to obtain.

1)Phase spectrum of real DTFT with linear delay factor introduced by fft shifting

2)Phase spectrum of real DTFT I wanted to plot