f0 = 20000;
a0 = 1;
alfa = 0;
T1 = f0/a0;
sigma1 = T1*cos(deg2rad(alfa));
tau1 = T1*sin(deg2rad(alfa));
for teta=0:1:180
tau(teta+1) = (sigma1/2)*sin(deg2rad(2*teta));
sigma(teta+1) = (sigma1/2)*(1+cos(deg2rad(2*teta)));
end
plot(sigma,tau);