曲率计算

仁者 · 发表于 2010-5-15 16:57:53

下面是第三届东北大学三点法算曲率的代码

研究了一下，基本看不懂，只知道他是求的三角形外接圆半径作为曲率半径，大部分代码都是做处理的，确保取得点是正确的，里面有一点我的注释，不知正确与否（原文无一点注释），那位强人能帮忙讲下不???大家一起研究下

//东北大学猎豹队
/*
void CountCurvature(uchar i1,uchar i2) // i1起始行，i2结束行
{
uchar mid; //临时值，辅助作用
int dx1,dx2,dx3,dy1,dy2,dy3; //三点坐标，用起始行和结束行表示的？？？
//应该是确保三点是稳定的点的
if(i1>i2)
{ //i1>i2，则交换
mid=i1;
i1=i2;
i2=mid;
}
//如果取得行数超过7行，则取中间行作为另一个点的位置
if(i1+6<=i2)
{
mid=(uchar)((i1+i2)/2); //取始末点行数的一半为之间点的行数,为uchar型的
dx1=(mid-i1)*3; //
dx2=(i2-i1-1)*3;
dx3=(i2-mid-1)*3;
dy1=(BlackDrift[mid]+BlackDrift[mid+1]-BlackDrift[i1]-BlackDrift[i1+1])/2; //BlackDrift 黑线趋势
dy2=(BlackDrift[i2]+BlackDrift[i2-1]-BlackDrift[i1]-BlackDrift[i1+1])/2;
dy3=(BlackDrift[i2]+BlackDrift[i2-1]-BlackDrift[mid]-BlackDrift[mid+1])/2;
PreCurvatureNum=dx1*dy2-dx2*dy1;
dy1/=2;
dy2/=2;
dy3/=2;
//Curvature 曲率
PreCurvatureDen=(long)m_sqrt(dx1*dx1+dy1*dy1)*m_sqrt(dx2*dx2+dy2*dy2)*m_sqrt(dx3*dx3+dy3*dy3);
if(PreCurvatureNum==0)
Curvature=0;
else
{
if(PreCurvatureDen>10000)
{
PreCurvatureNum=PreCurvatureNum*90/(PreCurvatureDen/100);
Curvature=(int)PreCurvatureNum;
}
else
{
PreCurvatureNum=PreCurvatureNum*9000/PreCurvatureDen;
Curvature=(int)PreCurvatureNum;
}
}
if(Curvature<-99)
Curvature=-99;
else if(Curvature>99)
Curvature=99;
CurvatureCache[3]=CurvatureCache[2];
CurvatureCache[2]=CurvatureCache[1];
CurvatureCache[1]=CurvatureCache[0];
if(i2-i1<=6)
CurvatureCache[0]=(Curvature*5+CurvatureCache[1]*5)/10;
else if(i2-i1<=10)
CurvatureCache[0]=(Curvature*6+CurvatureCache[1]*4)/10;
else if(i2-i1<=14)
CurvatureCache[0]=(Curvature*7+CurvatureCache[1]*3)/10;
else if(i2-i1<=18)
CurvatureCache[0]=(Curvature*8+CurvatureCache[1]*2)/10;
else
CurvatureCache[0]=Curvature;
if(CurvatureCache[0]>99)
CurvatureCache[0]=99;
else if(CurvatureCache[0]<-99)
CurvatureCache[0]=-99;
}
DiffCurvature=(2*CurvatureCache[0]-CurvatureCache[1]-CurvatureCache[2])/3;
}