计算机图形学学习系列的第二篇，我打算把Bresenham画圆算法从头到尾掰开揉碎讲一遍，并给出可以直接运行的matlab实现。这个算法在像素级绘图、嵌入式GUI、字符点阵屏这类场景里非常能打，核心思路就一句话：只用整数加法和比较，把一个圆“一笔一画”地光栅化出来。如果你上一篇跟着讲过Bresenham直线算法，那这篇对你来说基本是顺水推舟；就算没看过也问题不大，我会从判别式怎么来的开始推，确保每一步为什么这样写都清清楚楚。

1. 画圆为什么不用现成方程：Bresenham的思路起点

1.1 用参数方程画圆的三个痛处

很多人刚接触图形学时，第一个直觉是用三角函数。给定圆心(xc, yc)和半径r，让角度从0跑到2π，然后对每一个角度算一个点：

t = linspace(0, 2*pi, 1000); x = xc + r * cos(t); y = yc + r * sin(t);

如果只是为了“在电脑上画一个意思一下的圆”，这当然没问题。但如果你把它放到嵌入式的LCD屏幕上、放到一个只有几百字节RAM的MCU里、或者放到需要每帧画几百个圆的像素工具里，这个方案的三宗罪就暴露出来了。

第一，浮点运算和三角函数太重。当年的图形硬件几乎不支持浮点，就算支持，一次cos和sin的周期也远大于做几十次整数加减法。第二，角度均匀采样并不等于像素均匀。同样1度的弧长，在不同半径下对应的像素数量完全不一样，容易出现大圆点距过密、小圆点距过疏的情况。第三，用参数方程算出来的坐标是连续值，落到像素网格后需要做round，而round之后相邻点之间可能隔了两个像素，造成断断续续的虚线感。

一个自然的替代方案是把圆从方程变成增量过程：在某个八分圆弧上，每次让x固定加1，然后判断y应该保持原样还是减1。这就是增量光栅化的核心思想，也是Bresenham画圆算法能只用整数加减法工作的秘密。

1.2 光栅网格上的“选择题”

屏幕上的像素是离散的，每个像素只能用整数坐标表示。所谓光栅化，就是在一堆离散格子里挑出最接近真实线条的那一串格子。Bresenham算法的本质，其实是把“哪个像素更接近真实曲线”这个几何问题，转化成一个连除法都不用做的符号判断问题。

画直线时，我们知道斜率是固定的，x每走一步，y要么不动、要么加1，判断标准是“误差项是否超过半个像素”。画圆时情况类似，只不过斜率一直在变。在第一象限的45°弧上，x从0一直增大，y不断减小。当我们在某个像素(x, y)时，下一个候选像素只有两个：E点(x+1, y)保持当前y，或者SE点(x+1, y-1)向下走一格。关键就是要确定：到底选E还是SE，才能让这个像素离理想的圆弧最近。

这就是Bresenham画圆算法最核心的决策点。一旦你把这个问题想通了，后面的数学推导只是把这层直觉翻译成递推公式而已。

2. 核心判别式推导：一个递推式走天下

2.1 先用八分对称性把问题砍掉八分之七

这是所有画圆算法的常用起手式。圆在数学上和像素网格上都有很强的对称性：关于x轴对称、关于y轴对称、关于y=x直线对称。三个对称叠加在一起，就是八分对称。只要算出第一卦限里那一段，也就是x从0到x=y位置的45°圆弧，其余七个卦限的点全部可以通过坐标变换翻出来。

因为主循环只处理四分之一的一半，所以循环终止条件天然就是x <= y。一旦x和y相等，说明已经到了那个45°边界点，翻完它可以收工了。这个条件在实现里非常重要，写错成x < y会让圆缺最后那一个边界点，写错成x >= y则会让循环停不下来。

2.2 中点判别法：用圆的隐式方程做裁判

圆的隐式方程是：

F(x, y) = x² + y² - r²

这里有个方便的性质：如果点在圆上，F = 0；点在圆内，F < 0；点在圆外，F > 0。所以我们只需要用F的符号就能判断一个点相对于圆的位置。

假设当前画到了点(x, y)，那么下一列x+1上有两个候选点：E(x+1, y)和SE(x+1, y-1)。要判断谁更接近真实圆，可以取这两个候选点之间的中点M，也就是(x+1, y-0.5)，把M代入F：

p = F(x+1, y-0.5) = (x+1)² + (y-0.5)² - r²

如果p < 0，说明中点落在圆内部。这时候真实圆弧在这个中点的上方，所以E点离真实弧线更近，选择E，也就是y保持不变。如果p ≥ 0，说明中点落在圆外或圆上，选择SE，也就是y减1。

这个p就是我们的决策参数。光有初始式子还不够，因为每走一步都要重新算一遍平方的话，等于又回到了浮点运算。Bresenham的精华是递推：用上一步的p值，加上一个增量，就能得到下一步的p值。

推导过程不复杂，直接展开看。选E时，下一步的判别式是：

p_next = F(x+2, y-0.5) = (x+2)² + (y-0.5)² - r²

把p的表达式减去它，正好得到增量：

p_next = p + 2x + 3

选SE时，y变成了y-1，所以新的判别式是：

p_next = F(x+2, y-1.5) = (x+2)² + (y-1.5)² - r²

同样减一下：

p_next = p + 2x + 3 - 2y + 2 = p + 2(x - y) + 5

这两个增量式子，每一步都只涉及加法和当前坐标，完全没有乘方、浮点、三角函数。

2.3 初始值整数化：一个小细节

初始位置是(0, r)，所以第一个判别式是：

p0 = F(1, r-0.5) = 1 + (r-0.5)² - r² = 5/4 - r

这里有个0.5，看起来不太“整数”。但注意一点：我们只关心p的符号，不关心它的绝对值。既然4和1都是正数，把p乘以4或者整体平移，并不会改变p的正负判断。于是工程实践里常用的初始值就出来了：

p = 1 - r

它等价于原始判别式的符号判断，只是去掉了分数。很多教材里看到p初始值是3-2r的版本，其实也是把5/4-r乘以4后取整化简得到的，原理完全一样，别再被不同版本的符号吓到。

至此，算法框架非常简洁：

x = 0 y = r p = 1 - r while x <= y: 生成八个对称点(x, y)并平移圆心 if p < 0: p = p + 2*x + 3 # 选择E else: p = p + 2*(x - y) + 5 # 选择SE y = y - 1 x = x + 1

3. 用Matlab写一个干净的实现

3.1 主算法函数逐行解析

下面这个函数接受圆心坐标和半径，返回圆周上所有点的x坐标和y坐标。我在注释里把每一block在做什么都标清楚了。

function [px, py] = bresenhamCircle(xc, yc, r) % BRESENHAMCIRCLE 用Bresenham画圆算法生成圆周像素点 % 输入: % xc, yc - 圆心坐标 % r - 半径 % 输出: % px, py - 圆上所有像素点的坐标向量 x = 0; y = r; p = 1 - r; px = []; py = []; while x <= y % 把当前点的八个对称点全部加到结果里 [px, py] = addSymmetryPoints(px, py, xc, yc, x, y); % 根据判别式决定下一个像素选E还是SE if p < 0 p = p + 2*x + 3; % 选E，y不变 else p = p + 2*(x - y) + 5; % 选SE，y减1 y = y - 1; end x = x + 1; end end

在写这个函数时，我特意把“更新判别式”放在“x=x+1”之前，因为递推式里的x和y就应该是当前列的坐标。这个顺序一旦搞反，画出来的圆虽然大体是圆的，但某些半径下会在45°边界附近多出一个或少掉一个点。

3.2 八个对称点一次性生成

对称点函数是整段代码里比较容易写乱的地方。八个对称点分别是(x,y)的四种正负组合，以及(y,x)的四种正负组合。

function [px, py] = addSymmetryPoints(px, py, xc, yc, x, y) pts = [ x y; -x y; x -y; -x -y; y x; -y x; y -x; -y -x]; px = [px; xc + pts(:, 1)]; py = [py; yc + pts(:, 2)]; end

这里要注意，pt矩阵里每一行是一个对称点，列方向分别是x坐标和y坐标。如果只写了前四个，你只会得到一个x轴对称的四分圆；如果只写后四个，得到的则是y=x直线对称后的另外四条弧。八分对称要求两者都有，否则圆形会在45°方向“缺肉”。

当x等于y时，比如循环到了(5,5)，前后两组点里会重复出现同一个点，因为(x,y)和(y,x)变成了一模一样的坐标。这会在结果里产生一个重复点，但它不影响绘制，也不影响后续填充和连通性分析。如果你强迫症犯了，可以在函数末尾用unique去重。

3.3 绘制验证：看看像素圆和理论圆差多少

算法写完后，第一件事就是画出来跟理想圆对比。在matlab命令窗口里可以这样跑：

[x, y] = bresenhamCircle(0, 0, 30); plot(x, y, 'r.', 'MarkerSize', 10); hold on; t = linspace(0, 2*pi, 500); plot(30*cos(t), 30*sin(t), 'b-'); axis equal; legend('Bresenham像素', '理论圆');

画完你会看到红色点阵基本都落在蓝色理论圆的边界附近，最大误差不超过一个像素，这正是光栅化算法该有的表现。axis equal这行一定不能少。如果省略了，matlab默认的横纵坐标比例并不相等，一个正圆会被显示成椭圆，第一眼很容易误以为是算法错了。

如果你想更像“真实屏幕”那样一格一格地显示像素，可以把结果写进一个零矩阵，再用imshow看。这里有一个坑：矩阵索引是(行, 列)，也就是(y, x)，跟plot用的(x, y)相反，必须互换。

N = 101; img = zeros(N, N); [x, y] = bresenhamCircle(50, 50, 40); % 注意sub2ind里传入的是[y, x]，不是[x, y] idx = sub2ind([N, N], y + 1, x + 1); img(idx) = 1; imshow(img, 'InitialMagnification', 400);

圆心坐标50对应矩阵索引51，所以全部坐标都加1。半径40的圆在101×101的图像里刚好不会越界。用imshow看的时候，圆的台阶感会非常明显，尤其是半径比较小的情况下，那种“像素风”反而让算法特性一目了然。

4. 实操中会踩的坑和排查清单

4.1 圆形缺角或出现断点

最常见的问题是45°边界处少点。不少人会把循环条件写成while x < y，想着到x=y就多算了一个点。这种写法确实避免了重复点，但也丢掉了一个关键像素：那个刚好在45°对角线上的点。由于后续所有点都是从这一点对称展开的，缺了它，整个圆看起来就像在两条45°线上被各切了一刀。

第二个常见问题是判别式更新顺序错乱。比如写成先x=x+1再更新p，那么p公式里的x就不是当前列，误差项会整体偏移一位，导致某些半径下一段弧偏到圆外。实测下来，半径越大这种偏移越明显。

4.2 圆看起来是扁的

这个坑很多人会纠结半天。代码没改，圆就是扁的，原因通常不是算法，而是显示环境。在matlab里plot，如果没有axis equal，横纵坐标的缩放比例不一样，正圆当然变成椭圆。在LCD屏上，像素本身是长方形的，那么再正确的算法画出来的也是“物理上的椭圆”。这时候不是去改算法，而是要根据屏幕的物理宽高比做一次坐标缩放，或者把半径按比例换算。

另外，图像矩阵用imshow显示时，如果矩阵本身是正方形并且坐标系单位一致，椭圆问题一般不会出现。如果用了imagesc并且没有指定坐标轴范围，也可能被拉伸。

4.3 小半径的圆“不成圆”

半径小于等于3的时候，画出来的东西与其说是圆，不如说是一个对称的方块。这不是算法的bug，而是光栅网格的固有分辨率限制。你想在5×5网格里表达一个半径2的圆，总共就那么多像素，无论怎么选，看起来都不够圆。

遇到这种场景，建议考虑分辨率和半径的比例。如果半径在10以下还要追求视觉圆润，就需要上抗锯齿或者灰度加权。我在某个像素画工具里试过，半径小于5的形状直接用预置点位表，效果更好，速度也更快。

4.4 常见问题速查表

4.5 越界与裁剪

当圆心靠近屏幕边缘，或者半径大到超出显示范围，直接往图像矩阵里写就会越界报错。工程上一般先做一次简单的边界判断：如果xc+r、yc+r超过画布大小，就需要做裁剪。比如在图像矩阵里写入像素前，对坐标做一次范围筛选：

valid = (px >= 0) & (py >= 0) & (px < N) & (py < N); px = px(valid); py = py(valid);

裁剪虽然让圆变成残圆，但总比程序崩溃好得多。嵌入式中由于没有异常机制，这种边界检查几乎必须做。

5. 扩展玩法与我的使用心得

5.1 中点法、Bresenham法、整数判别式：同一个东西

我在这篇文章里用的是整数判别式版本，不少教材把它叫中点圆算法，又有一些书直接叫Bresenham画圆算法。这两者的关系经常把人搞晕。严格追根溯源，Bresenham在1962年提出的是基于误差项的直线算法；画圆算法的完整推导是后来基于隐式方程的中点法体系。但工程社区里早就习惯把这一整套“只比较符号、只用整数加法”的画圆实现统称为Bresenham画圆。理解了这一点，你看到不同教材里初始p值不一样就不会慌，它们只是对同一个符号判断做了不同化简。

如果你喜欢用“误差项”而不是“中点”去理解，可以换个角度：假设当前像素到理想圆弧的平方距离误差是ε，每次x前进一格后，比较“y不动”和“y减1”两种误差哪个小。这样推出来的递推式跟中点法完全等价，但会更贴近“误差累积”这个直觉。

5.2 从画圆到填充圆

只要圆周能画出来，填充圆就顺理成章。最简单的扫描线填充法是：对每一行，找到这一行上圆周左、右两端的像素，然后把两点之间的所有像素全部点亮。用上面这个算法生成圆周点后，可以按y坐标分组，对每个y取最小的x和最大的x作为左右边界，再画一条水平线段。

matlab示例：

function filledCircle = fillFromOutline(px, py, imgSize) filledCircle = zeros(imgSize, imgSize); for i = 1:length(px) xp = px(i) + 1; yp = py(i) + 1; if xp > 0 && yp > 0 && yp <= imgSize filledCircle(yp, xp) = 1; end end for y = 1:imgSize cols = find(filledCircle(y, :) > 0); if ~isempty(cols) filledCircle(y, min(cols):max(cols)) = 1; end end end

这个方法对凸图形非常稳，圆正好是凸的，所以不用处理复杂的边界交错问题。如果你要画很多同心圆，更聪明的做法是只算一次45°弧线点表，之后不同的半径直接查表平移缩放，速度能快出来一个量级。

5.3 抗锯齿：追求视觉圆润的下一步

纯Bresenham产生的像素只有0和1，所以边缘必然有锯齿。想更平滑，可以把每个像素看成一个小方块，计算这个方块被圆覆盖的面积比例，按比例给灰度。这就是“像素覆盖率抗锯齿”的基础思路。matlab里画图本身有抗锯齿渲染，但在图像矩阵上模拟逐像素光栅时并不生效，需要自己实现超采样或者面积加权。

我在实际使用中其实很少为静态图像做抗锯齿，因为成本高、效果也不如在更高分辨率下直接画完再缩小。但在字体渲染、地图应用那种需要长期凝视的界面里，边缘平滑度对观感影响非常大，值得多花一倍计算量。

5.4 我在实际项目里的一点点体会

这个算法我前前后后写了不知道多少遍，最早是在某块开发板上驱动一块点阵屏，那时候还没有浮点单元，所有绘图函数都得老老实实走整数路线。Bresenham画圆算法最大的好处不是快，而是“可预测”——你永远知道它下一列点亮的是哪个像素，这对做逐帧动画和碰撞检测特别友好。

后面做图像处理demo时，我用它往矩阵里画圆形掩膜，就踩过sub2ind行列颠倒的坑，也踩过imshow显示时坐标轴没有等比例导致圆形变椭圆的坑。这些坑都不是算法本身的，而是环境对坐标系的定义差异。我的经验是：凡是光栅化代码，第一步先把坐标系约定写在注释第一行，比如“这里的x是列方向，矩阵索引用(y+1, x+1)”，能帮你省掉无数排查时间。

最后给你一个小建议：把半径依次设成1、2、5、10、30、100，在matlab里多跑几遍，看看不同半径下的像素分布差异。这比我在这里说一百句“光栅化是离散化过程”都更直观。当你亲眼看到半径从1到100的变化轨迹，你就真正理解为什么Bresenham算法能在计算机图形学里存活几十年还不退休了。