How to calculate the phase spectrum of the square wave with fft?

there is an article said that the phase spectrum is π or 0. I don't know it's right or not. so i used fft function to check it. the result is not same as the article.the matlab code is below, is there something wrong?thank you very much.

if true
  Fs=1e2;
  N=100;
  t=0:1/Fs:(N-1)/Fs;             %计算时间范围 
   S=[zeros(1,25) ones(1,50) zeros(1,25)];% %     
  subplot(3,1,1);plot(t,S);           
   N=length(S);
  plot(S);
  title('original signal');
  figure;
  Y = fft(S,N); %做FFT变换
   Ayy = (abs(Y)); %取模
   figure;
  Ayy=Ayy/(N/2); %换算成实际的幅度
   Ayy(1)=Ayy(1)/2;
  F=Fs/N*([1:N]-1); %换算成实际的频率值
   plot(F(1:N/2),Ayy(1:N/2)); %显示换算后的FFT模值结果
   title('magnitude spectra ');
  figure;
  Pyy=angle(Y);
  plot(F(1:N/2),Pyy(1:N/2)); %显示相位图
   title(' phase spectrum');
end