求散点包围的闭合平面的面积

一、问题描述

求从点A1 到 点Ai 所包围的平面面积 <从A1到Ai 为闭合平面的一笔画顺序所得点>

二、原理

假设有任意的从 点A 到 点J 的闭合平面，如图：

图1

 

对该闭合图形以步长step切割成 面积为 step*step 的正方形。数该平面的格子即可得到近似的面积，精确度随step而改变，如图：

图2

三、程序实现原理

在图2中，我们可以看到闭合图形最低点为：点H，最高点为：点C， 最左点为：点B，最右点是：点F，由此可知：该闭合图形长为9*step，宽为 8*step。那么可以构建一个一个7*8的 全0二维数组。

闭合平面映射到二维数组方法：从点A 到 点B，可以从图2可得 在纵轴上升高3*step，在横轴偏移0*step。即在二维数组上可以映射成：area[0][2] = 1,area[0][3] = 1,area[0][4] = 1。如图3：

图3

继续重复上诉步骤，即可将散点包围的平面的边映射到二维数组，如图：

图4

从平面映射到二维数组后可轻易得到闭合曲面包围的格子数。

四、缺点

上诉的程序实现方法有一个致命缺点：当step减小 或者 平面图形的面积很大，有可能导致数组长度过长，导致堆栈溢出。

五、改进

使用二维数组完全记录闭合图形的空间位置是不必要的。

可知：闭合图形交y = ki 点A、B，蚂蚁从边的任意一点不回头移动总会经过A、B。因此：仅需得到交点A、B就能得到 y=ki 一行的格子数。在程序上仅需要一个长度为闭合图形的 宽/step 长度的一维数组就能满足需求 <数组元素保存与数组下标相对应的闭合图形的行的格子数>。在图2中的宽为8*step。如下图：

图5

 

六、程序实现

#include <iostream>


double find_maxPoint(float *X, int len)
{
    double max;
    double now = X[0];
    if(len==0){
        return 0;
    }
    for(int i=1; i<len; i++){
        if(now < X[i]){
            now = X[i];
        }
    }
    max = now;
    return max;
}

double find_minPoint(float *X, int len)
{
    double min;
    double now = X[0];
    if(len==0){
        return 0;
    }
    for(int i=1; i<len; i++){
        if(now > X[i]){
            now = X[i];
        }
    }
    min = now;
    return min;
}

/*
int abs(int num)
{
    return num>0 ? num : -num;
}
*/

//限制：只能计算任意方向的线过闭合曲面的交点 为偶数，奇数会导致该行总计算出一段格子 或0
double Area(float *X, float *Y, int len, double step = 0.01)
{
    //double step = 0.01; //步长，也指示小格子的面积
    //double littleLattice = step*step;
    double area = 0;
    double now_pointX, now_pointY;
    int adjust_nowPointX, adjust_nowPointY;
    double next_pointX, next_pointY;
    int adjust_nextPointX, adjust_nextPointY;
    double offset;
    int count = len;
    double highPoint, lowPoint, rightPoint, leftPoint;
    int matrixWidth, matrixLong;
    //double areaMatrix = 0;
    short int *link;
    short int jump; //跳跃
    short int i, j;


    if(count < 3){ //包围图形的最少点数是3点
        return 0;
    }

    //建立面积矩阵
    highPoint = find_maxPoint(Y, len);
    lowPoint = find_minPoint(Y, len);
    rightPoint = find_maxPoint(X, len);
    leftPoint = find_minPoint(X, len);
    area = ( (highPoint-lowPoint)/step + 1.0), matrixLong = (int)area;
    area = ( (rightPoint-leftPoint)/step  + 1.0), matrixWidth = (int)area;
    link = new short int [ matrixLong ];
    for (i=0; i<matrixLong; i++)
    {
        link[i] = 0x7fff;
    }area = 0;

    while(count--){
        now_pointX = X[count]/step;
        now_pointY = Y[count]/step;
        adjust_nowPointX = (int)(now_pointX);
        adjust_nowPointY = (int)(now_pointY);
        adjust_nowPointY -= (int)(lowPoint/step);
        adjust_nowPointY += 1;

        next_pointX = X[count-1]/step;
        next_pointY = Y[count-1]/step;
        adjust_nextPointX = (int)(next_pointX);
        adjust_nextPointY = (int)(next_pointY);
        adjust_nextPointY -= (int)(lowPoint/step);
        adjust_nextPointY += 1;

        if(count == 0){ // 为了在循环的最后一布回归到最后一点与最初一点闭合
            now_pointX = X[count]/step;
            now_pointY = Y[count]/step;
            adjust_nowPointX = (int)(now_pointX);
            adjust_nowPointY = (int)(now_pointY);
            adjust_nowPointY -= (int)(lowPoint/step);
            adjust_nowPointY += 1;

            next_pointX = X[len-1]/step;
            next_pointY = Y[len-1]/step;
            adjust_nextPointX = (int)(next_pointX);
            adjust_nextPointY = (int)(next_pointY);
            adjust_nextPointY -= (int)(lowPoint/step);
            adjust_nextPointY += 1;
        }

        jump = adjust_nextPointY-adjust_nowPointY;
        jump = abs(jump)>=1 ? jump>0 ? jump-1 : jump+1 : jump;

        if(adjust_nowPointX > 0x7ff0){
            continue;
        }

        for(i=0, j=0; i<=abs(jump); i++){
            if(jump > 0){
                j = i;
                offset = 1.0*(adjust_nextPointX - adjust_nowPointX)/jump;
            }else if(jump < 0){
                j = -i;
                offset = 1.0*(adjust_nextPointX - adjust_nowPointX)/jump;
            }else if(jump == 0){
                offset = 1.0*(adjust_nextPointX - adjust_nowPointX);
            }

            if(adjust_nowPointY+j >= matrixLong || adjust_nowPointY+j < 0){
                continue;
            }

            if(link[ adjust_nowPointY+j ] == 0x7fff){
                 link[ adjust_nowPointY+j ] = adjust_nowPointX + (int)( offset*i );
            }else if( link[ adjust_nowPointY+j ] < 0x7fff){
                area += abs( link[ adjust_nowPointY+j ] - (adjust_nowPointX + (int)( offset*i )) );
                link[ adjust_nowPointY+j ] = 0x7fff;
            }
        }
    }
    delete [] link;

    area = area*step*step;

    return area;
}


int main()
{
//    float X[] ={0,1,2,3,4,5,6};
//    float Y[] ={0,1,0,1,0,1,-1};
    float X[] ={1,1,0,0};
    float Y[] ={0,1,1,0};
    double c = Area(X,Y,4);
    printf("num is:%f",c);
    return 0;
}

七、二次改进

使用该方法进行计算散点包围得图形面积，要求散点是“顺序”的。那么能不能无序？可以的。只要先对散点进行空间排序就可以了。如下图：

对于上诉要求的空间排序方法，换种说法可转变为：规划一条路线从散点A1到Ai，其路线最短。

对于该算法网上可以找到。有一种最简单粗暴的方式是：穷举所有可能找最短。