Построение графика функции (Спираль Ферма) в полярных координатах

procedure TForm1.FormCreate(Sender: TObject);
begin
with Form1 do
 begin
  Position:=poScreenCenter;
  clientwidth:=500;
  clientheight:=550;
 end;
with Image1 do
 begin
  left:=10;
  top:=10;
  width:=Form1.ClientWidth-20;
  height:=Form1.ClientHeight-70;
  canvas.Brush.Color:=clWhite;
  canvas.Rectangle(0,0,width,height);
 end;
with Button1 do
 begin
  top:=Form1.ClientHeight-50;
  caption:='Нарисовать';
 end;
end;
 
procedure TForm1.Button1Click(Sender: TObject);
var x0,y0,i,a,q,r,x1,y1,dx,dy:integer;
    m,r1,f,mx:real;
begin
x0:=Image1.Width div 2;
y0:=Image1.Height div 2;
a:=2;
q:=2;
mx:=a*sqrt(q*2*pi);{максимальное значение радиуса}
m:=(y0-80)/mx;{масштаб координатной сетки}
//КООРДИНАТНАЯ  СЕТКА
with Image1.Canvas do
 begin
  pen.Color:=clBlue;
  brush.Style:=bsClear;
  for i:=1 to trunc(mx)+1 do
   begin
    r:=round(i*m);
    ellipse(x0-r,y0-r,x0+r,y0+r);//рисуем окружности
    textout(x0+r,y0+10,inttostr(i));//пишем шкалу
   end;
  f:=30;//шаг по кругу=30 град
  for i:=0 to 12 do//делим на 12 частей
   begin
   //рисуем радиусы пунктиром
    pen.Style:=psDot;
    x1:=x0+round((r+20)*cos(i*f*pi/180));
    y1:=y0-round((r+20)*sin(i*f*pi/180));
    moveto(x0,y0); lineto(x1,y1);
   //пишем подписи
    if x1>=x0 then dx:=5 else dx:=-10;
    if y1<=y0 then dy:=-15 else dy:=10;
    if i*30<360 then textout(x1+dx,y1+dy,inttostr(i*30));
   end;
  //рисуем осевые линии сплошной линией
  pen.Style:=psSolid;
  moveto(10,y0);lineto(Image1.Width-10,y0);
  moveto(x0,10);lineto(x0,Image1.Height-10);
  font.Color:=clBlue;
  font.Style:=[fsBold];
  textout(Image1.Width-30,y0+10,'R');
  //ГРАФИК
  f:=0;//начальный угол
  pen.Width:=2;
  pen.Color:=clRed;
  moveto(x0,y0);
  while f<=2*q*pi do//полный оборот
   begin
    r1:=a*sqrt(f);//уравнение в полярных координатах
    x1:=x0+round(m*r1*cos(f));//приведение к экранным
    y1:=y0-round(m*r1*sin(f));
    lineto(x1,y1);
    f:=f+0.001;
   end;
  moveto(x0,y0);
  f:=0;
  while f<=2*q*pi do//полный оборот
   begin
    r1:=a*sqrt(f);//уравнение в полярных координатах
    x1:=x0-round(m*r1*cos(f));//приведение к экранным
    y1:=y0+round(m*r1*sin(f));
    lineto(x1,y1);
    f:=f+0.001;
   end;
 end;
end;