三次贝塞尔曲线绘制算法

Posted 2020-08-20

原文：http://www.cnblogs.com/flash3d/archive/2012/01/30/2332176.html

 

源码：http://files.cnblogs.com/flash3d/bezier.rar

 

====================================================

 

这学期学图形学，就把自己的一些粗浅的理解发上去让大家拍砖。前些天做三次贝塞尔曲线绘制的上机练习，正好将从直线扫描算法中启发得来的n次多项式批量计算用上了，自认为优化得还可以。原来写的版本是C++，为了便于网上观看特改写成AS3，对这方面内行的朋友，希望还多多指点！

在讲三次贝塞尔曲线之前，先给出n次贝塞尔曲线的一般式：

[R]( t ) = ( +..i=0,n ) ( [R]i * Bi,n( t ) ) , t属于[0,1]

Bi,n( t ) = Cni * pow( 1 - t , n - i ) * pow( t , i )

这里我把符号说明下，中括号[]括起来的表示矩阵。比如[R]就是个矩阵，如果我们讨论的贝塞尔曲线的点在平面内，那么矩阵里面有两个元素，比如我们熟知的x,y，如果是空间的曲线，就是三个元素，元素再多。。那就不知道是几次元了，反正是个数学概念。在计算的时候，将矩阵的元素单独拎出来计算也是可以的，因为公式中并没有出现矩阵和矩阵相乘，也就是说，各元素计算时其实各不相干。

( +..i=0,n ) 这个运算符，这个符号是我自己发明的，其实是累加符号，表示公式后面括号内将i分别用0,1,2,3..n代替，每个代替运算的结果全部加起来。

Cni是二项式系数，其值为 n! / ( i! * ( n - i )! )

pow是求指数函数。

其中 Bi,n( t ) = Cni * pow( 1 - t , n - i ) * pow( t , i ) 为n次二项式展开后的一项。其意义在后面的三次贝塞尔曲线里面讨论更容易理解

那么，现在先讨论3次贝塞尔曲线，我们把公式中的n全部替换成3，得到三次贝塞尔曲线的一般式如下：

[R]( t ) = ( +..i=0,3 ) ( [R]i * Bi,3( t ) ) , t属于[0,1]

Bi,3( t ) = C3i * pow( 1 - t , 3 - i ) * pow( t , i )

直接把Bi,3( t )代入吧，得到式子：

[R]( t ) = ( +..i=0,3 ) ( [R]i * C3i * pow( 1 - t , 3 - i ) * pow( t , i ) ) , t属于[0,1]

对于不熟悉线性代数或者不习惯用矩阵表示数据的朋友，我们就直接把矩阵拆开好了，[R]矩阵中总共就两个元素，分别是x和y，这里就把x和y分别计算，得到关于x和y的参数曲线：

X( t ) = ( +..i=0,3 ) ( Xi * C3i * pow( 1 - t , 3 - i ) * pow( t , i ) ) , t属于[0,1]

Y( t ) = ( +..i=0,3 ) ( Yi * C3i * pow( 1 - t , 3 - i ) * pow( t , i ) ) , t属于[0,1]

如果对参数曲线不熟悉，可以参考直线的参数方程

比如一条直线经过(x1,y1)，(x2,y2)，那么其参数方程可表示为：

X(t)=x1+t*(x2-x1)

Y(t)=y1+t*(y2-y1)

如果是线段那么t在0和1之间，如果是直线，t是实数。

言归正传，下面我们开始处理参数曲线，由于x的参数方程和y的形式差不多，我们就以x的参数方程做例子来分析

X( t ) = ( +..i=0,3 ) ( Xi * C3i * pow( 1 - t , 3 - i ) * pow( t , i ) ) , t属于[0,1]

首先我们把连加拆开，那么i分别等于0,1,2,3，得到一个式子：

X( t ) =

X0 * C30 * pow( 1 - t , 3 - 0 ) * pow( t , 0 )

+

X1 * C31 * pow( 1 - t , 3 - 1 ) * pow( t , 1 )

+

X2 * C32 * pow( 1 - t , 3 - 2 ) * pow( t , 2 )

+

X3 * C33 * pow( 1 - t , 3 - 3 ) * pow( t , 3 )

t属于[0,1]

看这个式子，我们可以发现，C30 * pow( 1 - t , 3 - 0 ) * pow( t , 0 )+C31 * pow( 1 - t , 3 - 1 ) * pow( t , 1 )+C32 * pow( 1 - t , 3 - 2 ) * pow( t , 2 )+C33 * pow( 1 - t , 3 - 3 ) * pow( t , 3 )=pow( 1 - t + t , 3)=1

那就是说，x0，x1，x2，x3相加前的系数和为1，我们称这些系数为各点的权重。当然，点的权重是随着t的值连续变化的，这也是为什么三次贝塞尔曲线能光滑地通过四个点控制的原因。

下面我们进一步对式子进行化简。

其中C30等于1，C31等于3，C32等于3，C33等于1，把能合并的都合并，能计算值的都计算出来，将式子化简：

X( t ) =

X0 * pow( 1 - t , 3 )

+

X1 * 3 * pow( 1 - t , 2 ) * t

+

X2 * 3 * ( 1 - t ) * pow( t , 2 )

+

X3 * pow( t , 3 )

t属于[0,1]

这样的式子应该很好求了，他是一个最高次为3次的多项式，在《Alchemy的使用和多项式批量计算的优化》中提到，a*pow(x,n)+b*pow(x,n-1)+..这样形式的式子，都能通过多项式批量计算方法来快速求得结果。

对于最高三次的式子来说，我们需要前三个原始数据，以及式子求导三次的常数，算法细节就不讲了，细节参考上面那篇博文。

 

package
{

    import flash.display.Bitmap;

    import flash.display.BitmapData;

    import flash.display.Shape;

    import flash.display.Sprite;

    import flash.events.Event;

    import flash.events.MouseEvent;

    import flash.geom.Point;

    import flash.geom.Rectangle;

    /**

        *贝塞尔曲线绘制算法

        * @author boycy815

        */

    public class Main extends Sprite

    {

        /*C30到C33的值总共才4个直接枚举*/

        private const C30:int = 1;

        private const C31:int = 3;

        private const C32:int = 3;

        private const C33:int = 1;

        /*定义曲线上点的个数（包括两端）*/

        private const LENGTH:int = 10000;

        /*btC30到btC33分别是R0 R1 R2 R3四个点当t取前三个值时的权重值*/

        private var btC30:Vector.<Number> = new Vector.<Number>(3);

        private var btC31:Vector.<Number> = new Vector.<Number>(3);

        private var btC32:Vector.<Number> = new Vector.<Number>(3);

        private var btC33:Vector.<Number> = new Vector.<Number>(3);

        /*画点的位图*/

        private var data:BitmapData;

        private var map:Bitmap;

        /*用于连接四个点的四条直线*/

        private var lin:Shape;

        /*分别是四个点的显示对象*/

        private var ct0:Sprite;

        private var ct1:Sprite;

        private var ct2:Sprite;

        private var ct3:Sprite;

        public function Main():void

        {

            if (stage) init();

            else addEventListener(Event.ADDED_TO_STAGE, init);

        }

        

        private function init(e:Event = null):void

        {

            removeEventListener(Event.ADDED_TO_STAGE, init);

            /*初始化位图和一些显示对象*/

            data = new BitmapData(800, 600,false);

            map = new Bitmap(data);

            this.addChild(map);

            lin = new Shape();

            this.addChild(lin);

            ct0 = new Sprite();

            ct0.buttonMode = true;

            ct0.graphics.lineStyle(0, 0x00000);

            ct0.graphics.beginFill(0xffffff);

            ct0.graphics.drawCircle(0, 0, 5);

            ct0.graphics.endFill();

            ct0.addEventListener(MouseEvent.MOUSE_DOWN, mouseDown);

            ct0.addEventListener(MouseEvent.MOUSE_UP, mouseUp);

            this.addChild(ct0);

            ct1 = new Sprite();

            ct1.buttonMode = true;

            ct1.graphics.lineStyle(0, 0x00000);

            ct1.graphics.beginFill(0xffffff);

            ct1.graphics.drawCircle(0, 0, 5);

            ct1.graphics.endFill();

            ct1.addEventListener(MouseEvent.MOUSE_DOWN, mouseDown);

            ct1.addEventListener(MouseEvent.MOUSE_UP, mouseUp);

            this.addChild(ct1);

            ct2 = new Sprite();

            ct2.buttonMode = true;

            ct2.graphics.lineStyle(0, 0x00000);

            ct2.graphics.beginFill(0xffffff);

            ct2.graphics.drawCircle(0, 0, 5);

            ct2.graphics.endFill();

            ct2.addEventListener(MouseEvent.MOUSE_DOWN, mouseDown);

            ct2.addEventListener(MouseEvent.MOUSE_UP, mouseUp);

            this.addChild(ct2);

            ct3 = new Sprite();

            ct3.buttonMode = true;

            ct3.graphics.lineStyle(0, 0x00000);

            ct3.graphics.beginFill(0xffffff);

            ct3.graphics.drawCircle(0, 0, 5);

            ct3.graphics.endFill();

            ct3.addEventListener(MouseEvent.MOUSE_DOWN, mouseDown);

            ct3.addEventListener(MouseEvent.MOUSE_UP, mouseUp);

            this.addChild(ct3);
            /*初始化计算的一些准备数据*/

            btInit();

            /*默认进行第一次定位和绘制*/

            ct0.x = 100;

            ct0.y = 100;

            ct1.x = 700;

            ct1.y = 100;

            ct2.x = 700;

            ct2.y = 500;

            ct3.x = 100;

            ct3.y = 500;

            lin.graphics.lineStyle(0, 0x000000);

            lin.graphics.moveTo(ct0.x, ct0.y);

            lin.graphics.lineTo(ct1.x, ct1.y);

            lin.graphics.lineTo(ct2.x, ct2.y);

            lin.graphics.lineTo(ct3.x, ct3.y);

            drawBezier3(new Point(ct0.x, ct0.y), new Point(ct1.x, ct1.y), new Point(ct2.x, ct2.y), new Point(ct3.x, ct3.y));

        }

        /*按钮按下便开始可以拖动点*/

        private function mouseDown(e:MouseEvent):void

        {

            e.target.startDrag();

            stage.addEventListener(MouseEvent.MOUSE_MOVE, mouseMove);

        }

        /*鼠标移动便不断绘制曲线*/

        private function mouseMove(e:MouseEvent):void

        {

            /*绘制连接四个点的直线*/

            lin.graphics.clear();

            lin.graphics.lineStyle(0, 0x000000);

            lin.graphics.moveTo(ct0.x, ct0.y);

            lin.graphics.lineTo(ct1.x, ct1.y);

            lin.graphics.lineTo(ct2.x, ct2.y);

            lin.graphics.lineTo(ct3.x, ct3.y);

            /*绘制曲线*/

            data.fillRect(new Rectangle(0, 0, 800, 600), 0xffffffff);

            drawBezier3(new Point(ct0.x, ct0.y), new Point(ct1.x, ct1.y), new Point(ct2.x, ct2.y), new Point(ct3.x, ct3.y));

        }

        /*释放按钮停止拖动*/

        private function mouseUp(e:MouseEvent):void

        {

            var n:int;

            e.target.stopDrag();

            stage.removeEventListener(MouseEvent.MOUSE_MOVE, mouseMove);

        }

        /***********算法部分**********/

        

        /*初始化计算所需数据，主要是btC30到btC33*/

        private function btInit():void

        {

            var add:Number;

            /*由于t的值域为[0,1]，要绘制LENGTH个点，则必须计算出每次t的增量*/

            add = 1.0 / (LENGTH - 1.0);

            /*三次白赛尔曲线Ri的权重公式为：C3i * pow( t - 1 , 3 - i ) * pow( t , i )

                    *前三次t的值分别为0，add，add*2

                    *则代入公式就能求出R0 R1 R2 R3的前三次的权重值

                    * */

            btC30[0] = 1;

            btC30[1] = C30 * Math.pow(1 - add, 3);

            btC30[2] = C30 * Math.pow(1 - add * 2, 3);

            btC31[0] = 0;

            btC31[1] = C31 * Math.pow(1 - add, 2) * add;

            btC31[2] = C31 * Math.pow(1 - add * 2, 2) * add * 2;

            btC32[0] = 0;

            btC32[1] = C32 * (Math.pow(add, 2) - Math.pow(add, 3));

            btC32[2] = C32 * (Math.pow(add * 2, 2) - Math.pow(add * 2, 3));

            btC33[0] = 0;

            btC33[1] = C33 * Math.pow(add, 3);

            btC33[2] = C33 * Math.pow(add * 2, 3);

        }

        /*利用准备好的初始数据，通过递推方法求出四个点对应的贝塞尔曲线上每个点的坐标，并绘制

            *计算中用到多项式批量计算优化策略，上一篇博文中有提及具体的方法和思想

            * */

        private function drawBezier3(p0:Point,p1:Point,p2:Point,p3:Point):void

        {

            /*qun_x和qun_x本别为x坐标和y坐标的临时数据数组*/

            var qun_x:Vector.<Number> = new Vector.<Number>(3);

            var qun_y:Vector.<Number> = new Vector.<Number>(3);

            /*dd_x和dd_y分别为x坐标的多项式和y坐标式求导三次后得到的常数*/

            var dd_x:Number;

            var dd_y:Number;

            var n:int;

            /*计算前三次种子数据*/

            dd_x = 6 * (3 * p1.x - 3 * p2.x - p0.x + p3.x) / Math.pow(LENGTH - 1.0, 3);

            qun_x[0] = p0.x * btC30[0] + p1.x * btC31[0] + p2.x * btC32[0] + p3.x * btC33[0];

            dd_y = 6 * (3 * p1.y - 3 * p2.y - p0.y + p3.y) / Math.pow(LENGTH - 1.0, 3);

            qun_y[0] = p0.y * btC30[0] + p1.y * btC31[0] + p2.y * btC32[0] + p3.y * btC33[0];

            data.setPixel(qun_x[0], qun_y[0], 0x000000);

            qun_x[1] = p0.x * btC30[1] + p1.x * btC31[1] + p2.x * btC32[1] + p3.x * btC33[1];

            qun_x[0] = qun_x[1] - qun_x[0];

            qun_y[1] = p0.y * btC30[1] + p1.y * btC31[1] + p2.y * btC32[1] + p3.y * btC33[1];

            qun_y[0] = qun_y[1] - qun_y[0];

            data.setPixel(qun_x[1], qun_y[1], 0x000000);

            qun_x[2] = p0.x * btC30[2] + p1.x * btC31[2] + p2.x * btC32[2] + p3.x * btC33[2];

            qun_x[1] = qun_x[2] - qun_x[1];

            qun_x[0] = qun_x[1] - qun_x[0];

            qun_y[2] = p0.y * btC30[2] + p1.y * btC31[2] + p2.y * btC32[2] + p3.y * btC33[2];

            qun_y[1] = qun_y[2] - qun_y[1];

            qun_y[0] = qun_y[1] - qun_y[0];

            data.setPixel(qun_x[2], qun_y[2], 0x000000);

            /*进行递推计算的循环*/

            for (n = 3; n < LENGTH; n++)

            {

                qun_x[0] += dd_x;

                qun_x[1] += qun_x[0];

                qun_x[2] += qun_x[1];

                qun_y[0] += dd_y;

                qun_y[1] += qun_y[0];

                qun_y[2] += qun_y[1];

                data.setPixel(qun_x[2], qun_y[2], 0x000000);

            }

        }

    }
}