ADD the following function call to your script:
circle((S1+S2)/2, 0, Tm) % function call
Then create a new function and save to the same (current) folder as your script:
    function h = circle(x,y,r)
    hold on
    th = 0:pi/50:2*pi;
    xunit = r * cos(th) + x;
    yunit = r * sin(th) + y;
    h = plot(xunit, yunit);
    grid on
    hold off
    end