function [t, y, w, amp, phase] = sqrfunc(freq_on, a, b, t0, c, w0);

t = -4.00001:0.01:3.9999;  % prevent division by 0
ImpResp = 0.01*sin(w0*t)./(pi*t);

if (freq_on==0),
   T = 1; 
   t = -4:0.01:4;
   y = a/c*[ustep(b/c*t+T-t0) - ustep(b/c*t-T-t0)];
 
   w = -6*pi:0.05:6*pi;
   Y = 2*a/b*sin(c*w/b*T)./(c*w/b).*exp(-(j*w*(t0+0.001)));

   phase = angle(Y);
   amp = abs(Y);

else
   T = 1;
   t = -4:0.01:4;
   y = a/c*(ustep(b/c*t+T-t0) - ustep(b/c*t-T-t0));
   yc = conv(y, ImpResp);
   y = yc(400:400+length(y)-1);

   w = -6*pi:0.05:6*pi;
   Y = 2*a/b*sin(c*w/b*T)./(c*w/b).*exp(-(j*w*t0)).*(ustep(w+w0) - ustep(w-w0));

   phase = angle(Y);
   amp = abs(Y);
end