求助：snake算法详解

SHASHA518 2008-03-24 04:47:44

bool termination_point
bool neighbor3by3; /* true means 3x3 neighborhood, false means 5x5 neighborhood */
int m_Cols; /* Number of Columns */
int m_Rows; /* Number of Rows */
int no_of_snake_points;
CPoint Snake_points[200]; /* Snake Points */
double *alpha,*beta, *gamma; /* Weights for energies */
double threshold_curvature;
int threshold_movement;

void Snake_algorithm()
{
bool flag_change;
int movement,i,j,iteration=0;
double avg_distance=0.0, max_curvature;
CPoint temp;
double *curvature;
alpha=new double[no_of_snake_points];
beta=new double[no_of_snake_points];
gamma=new double[no_of_snake_points];
curvature=new double[no_of_snake_points];
termination_point=false;

for(i=0;i<no_of_snake_points;i++)
{
*(alpha+i)=1.0;
*(beta+i)=1.0;
*(gamma+i)=1.2;
avg_distance+=find_distance(i, Snake_points);
}
j=no_of_snake_points;
pdc=GetDC();
while(!termination_point)
{
movement=0;
avg_distance=avg_distance/(double)no_of_snake_points;
max_curvature=-1000000000.0;
for(i=0;i<no_of_snake_points;i++)
{
temp=find_min_energy(i, Snake_points,avg_distance);
flag_change=false;
if(temp!=Snake_points&&temp!=Snake_points[(i-1+j)%j]&&temp!=Snake_points[(i+1)%j])
{
Snake_points=temp;
movement++;
}
curvature=find_curvature(i, Snake_points);
if(max_curvature<curvature)
max_curvature=curvature;
}
avg_distance=0.0;
for(i=0;i<no_of_snake_points;i++)
curvature=curvature/max_curvature;
for(i=0;i<no_of_snake_points;i++)
{
avg_distance+=find_distance(i, Snake_points);
if(curvature>threshold_curvature&&curvature>curvature[(i+1)%no_of_snake_points]&&curvature>curvature[(i-1+no_of_snake_points)%no_of_snake_points])
*(beta+i)=0;
}
if(movement<threshold_movement)
termination_point=true;
iteration++;
if(iteration>Max_Iterations)
termination_point=true;
}
delete alpha;
delete beta;
delete gamma;
delete curvature;
}

double find_distance(int no, CPoint point)
{
int x=no_of_snake_points;
point-=Snake_points[(no-1+x)%x];
return(sqrt(point.x*point.x+point.y*point.y));
}

double find_curvature(int no, CPoint point)
{
int x=no_of_snake_points;
point=Snake_points[(no-1+x)%x]-point-point+Snake_points[(no+1)%x];
return(point.x*point.x+point.y*point.y);
}

double find_continuity(int no, CPoint point, double avg_distance)
{
return(pow(avg_distance-find_distance(no,point),2));
}

CPoint find_min_energy(int no, CPoint point, double avg_distance)
{
CPoint p, min_point;
double max_curvature, max_continuity, max_internal, min_internal, min_energy, energy;
double curvatures[5][5];
double continuities[5][5];
double internal_energies[5][5];
int i,j, limit=1;
max_curvature=max_continuity=max_internal=-1000000000000.0;
min_internal=1000000000000.0;
if(!neighbor3by3)
limit++;
for(i=-limit;i<=limit;i++)
{
p.y=point.y+i;
if(p.y<0)
p.y=0;
if(p.y>=m_Rows)
p.y=m_Rows-1;
for(j=-limit;j<=limit;j++)
{
p.x=point.x+j;
if(p.x<0)
p.x=0;
if(p.x>=m_Cols)
p.x=m_Cols-1;
curvatures[i+limit][j+limit]=find_curvature(no, p); //This code can cause problem near
continuities[i+limit][j+limit]=find_continuity(no,p,avg_distance); //border of image
if(curvatures[i+limit][j+limit]>max_curvature)
max_curvature=curvatures[i+limit][j+limit];
if(continuities[i+limit][j+limit]>max_continuity)
max_continuity=continuities[i+limit][j+limit];
if(internal_energies[i+limit][j+limit]>max_internal)
max_internal=internal_energies[i+limit][j+limit];
if(internal_energies[i+limit][j+limit]<min_internal)
min_internal=internal_energies[i+limit][j+limit];
}
}
for(i=0;i<=2*limit;i++)
{
for(j=0;j<=2*limit;j++)
{
curvatures[j]=curvatures[j]/max_curvature;
continuities[j]=continuities[j]/max_continuity;
internal_energies[j]=(internal_energies[j]-min_internal)/(max_internal-min_internal);
}
}
min_point.x=-limit;
min_point.y=-limit;
min_energy=1000000000000.0;
for(i=-limit;i<=limit;i++)
{
for(j=-limit;j<=limit;j++)
{
energy=*(alpha+no)*continuities[i+limit][j+limit]+*(beta+no)*curvatures[i+limit][j+limit]-*(gamma+no)*internal_energies[i+limit][j+limit];
if(energy<min_energy||(energy==min_energy&&i==0&&j==0))
{
min_energy=energy;
min_point.x=j;
min_point.y=i;
}
}
}
min_point.x=min_point.x+point.x;
min_point.y=min_point.y+point.y;
return(min_point);
...全文
532 9 打赏 收藏 转发到动态 举报

9 条回复

TandyT 2011-03-30
• 打赏
• 举报

• 打赏
• 举报

# SHASHA518
#
#
# 等　级：
# 结帖率：0.00%
maths95 2011-03-29
• 打赏
• 举报

zljcj 2008-05-14
• 打赏
• 举报

knowledge_Is_Life 2008-04-30
• 打赏
• 举报

• 打赏
• 举报

kaizitop 2008-03-25
• 打赏
• 举报

http://blog.csdn.net/byxdaz/archive/2006/04/18/667432.aspx

SHASHA518 2008-03-25
• 打赏
• 举报

kaizitop 2008-03-24
• 打赏
• 举报

19,444

• 近7日
• 近30日
• 至今