プレゼンハムで円を描画する

プレゼンハムのアルゴリズム

一般に乗除算は加減算に比べて数十倍の実行時間を要します。
プレゼンハムのアルゴリズムの凄い所は、乗除算を一切使わずに円周を8分割して単純な加減算で描画してしまうのです。
その方法を説明しましょう。
  1. 中心座標を x,y 半径を r とします。
  2. 次の計算をします。
    f= r*2+3
    cx= r
    dx= 0
  3. cx>=dx の間、以下の手続きを繰り返します。
これで見事に円が描画されるのです。
乗除算を一切使わないと言いながら、乗算を使っているではないか? と思うような人はいませんよね。
cx とか dx とか変な名前を使っていますが、元ネタは数十年前にアセンブラで組んだプログラムの名残です。
それでは実際に乗除算を一切使わずに円を描画して下さい。

全ソースコード

DOS Graph Mode で作成した全ソースコードを提供します。 (^_^;)
そのままでは、VC++ の現在のバージョンでは動かないかも知れません。 (;_;)

/*★ プレゼンハムで円を描く ★*/
#include <stdio.h>
#include <stdlib.h>
#include <graph.h>

/*プレゼンハム関数*/
void    en(short x, short y, short r)
{   short   f,dx,cx;

    f= r+r+3;
    cx= r;
    dx= 0;
    while(cx>=dx)
    {   _setpixel(x+dx,y+cx);       //5時の方向
        _setpixel(x+dx,y-cx);       //1時の方向
        _setpixel(x-dx,y+cx);       //7時の方向
        _setpixel(x-dx,y-cx);       //11時の方向
        _setpixel(x+cx,y+dx);       //4時の方向
        _setpixel(x+cx,y-dx);       //2時の方向
        _setpixel(x-cx,y+dx);       //8時の方向
        _setpixel(x-cx,y-dx);       //10時の方向
        if (f>=0)
        {   cx--;
            f-= cx*4;
            //f-= cx<<2;
            //f-= (cx+cx+cx+cx);
        }
        dx++;
        f+= dx*4+2;
    }
}

/*★ MAIN ★*/
void    main(void)
{   short   i,x,y,color;

    //スクリーンの表示モードの設定
    if (!_setvideomode(_MAXRESMODE))    exit(1);
    //円の描画
    for(color=0; color<16; color++)
    {   _setcolor(color);
        en(320,240,color*5);
    }
    //どれかキーを押すと終了
    _getch();
    _setvideomode(_DEFAULTMODE);
}


Assembler Source Code

16ビットアセンブラの円を描くサブルーチンのソースコードです。
そのうち忘却のかなたに消えてしまいそうなので掲載してみました。
古き良き時代に思いをはせて、何かの参考にして下さい。 (^_^;)

;★SI(行),DI(列)から、オフセットとビットパターンを求める★
;★BX←SI(行),DI(列)のオフセット  AH←ビットパターン★
setoff  PROC    NEAR
        push    cx
        mov     bx,di
        mov     cl,bl
        and     cl,7            ;CL=列の下3ビット
        mov     ah,80H          ;AH=ビットパターン
        shr     ah,cl
        shr     bx,3            ;BX=列/8
        imul    cx,si,80        ;CX=行*80
        add     bx,cx           ;BX=行*80+列/8
        pop     cx
        ret
setoff  ENDP
;★ALの指示する色領域の[BX]に、AHのパターンをorする★
;★AL(色):[BX]←AH    AL=1:青,2:赤,4:緑 ★
setdot  PROC    NEAR
        push    ax
        push    dx
        and     al,7
        mov     dx,0A800H       ;DX=青色セグメント
sl1:    test    al,1            ;該当色=0
        jz      @f
        mov     es,dx           ;色の種類毎にセット
        or      es:[bx],ah
@@:     add     dx,800H
        shr     al,1            ;次の色に設定
        jnz     sl1
        pop     dx
        pop     ax
        ret
setdot  ENDP
;★SI(行),DI(列)←AL(色)のドットを表示★
dot     PROC    NEAR
        call    setoff
        call    setdot
        ret
dot     ENDP
.CODE
;★SI(行),DI(列)を中心にして、半径BP(最大200)の円を表示★
en      PROC    NEAR
        pusha
        cmp     cx,0            ;(半径≦0)の時、中心点に表示
        jg      @f
        call    dot
        jmp     owari
@@:     mov     y0,si           ;中心点をセーブ
        mov     x0,di
        mov     f,cx
        shl     f,1
        neg     f
        add     f,3             ;F=-(R*2)+3
        mov     dx,0
pset:   cmp     cx,dx           ;(CX<DX)の時終了
        jae     @f
        jmp     owari
;5時の方向(Y0+CX,X0+DX)の座標にドットを表示する
@@:     mov     si,y0
        add     si,cx
        mov     di,x0
        add     di,dx
        call    dot
;7時の方向(Y0+CX,X0-DX)の座標にドットを表示する
        sub     di,dx
        sub     di,dx
        call    dot
;11時の方向(Y0-CX,X0-DX)の座標にドットを表示する
        sub     si,cx
        sub     si,cx
        call    dot
;1時の方向(Y0-CX,X0+DX)の座標にドットを表示する
        add     di,dx
        add     di,dx
        call    dot
;4時の方向(Y0+DX,X0+CX)の座標にドットを表示する
        mov     si,y0
        add     si,dx
        mov     di,x0
        add     di,cx
        call    dot
;8時の方向(Y0+DX,X0-CX)の座標にドットを表示する
        sub     di,cx
        sub     di,cx
        call    dot
;10時の方向(Y0-DX,X0-CX)の座標にドットを表示する
        sub     si,dx
        sub     si,dx
        call    dot
;2時の方向(Y0-DX,X0+CX)の座標にドットを表示する
        add     di,cx
        add     di,cx
        call    dot
        cmp     cx,2
        jl      owari           ;(半径<2)の時終了
        cmp     f,0
        jl      @f
        dec     cx              ;CX=CX+1
        mov     bx,cx
        shl     bx,2
        sub     f,bx            ;F=F-CX*4
;        jmp     pset
@@:     inc     dx              ;DX=DX+1
        mov     bx,dx
        shl     bx,2
        add     bx,2
        add     f,bx            ;F=F+DX*4+2
        jmp     pset
owari:  popa
        ret
.DATA
x0      DW      ?
y0      DW      ?
f       DW      ?
en      ENDP
        END     start

超初心者のプログラム入門(Win32API C++)