?? imagesegmentation.cpp
字號(hào):
//Copyright (c) 2004-2005, Baris Sumengen
//All rights reserved.
//
// CIMPL Matrix Performance Library
//
//Redistribution and use in source and binary
//forms, with or without modification, are
//permitted provided that the following
//conditions are met:
//
// * No commercial use is allowed.
// This software can only be used
// for non-commercial purposes. This
// distribution is mainly intended for
// academic research and teaching.
// * Redistributions of source code must
// retain the above copyright notice, this
// list of conditions and the following
// disclaimer.
// * Redistributions of binary form must
// mention the above copyright notice, this
// list of conditions and the following
// disclaimer in a clearly visible part
// in associated product manual,
// readme, and web site of the redistributed
// software.
// * Redistributions in binary form must
// reproduce the above copyright notice,
// this list of conditions and the
// following disclaimer in the
// documentation and/or other materials
// provided with the distribution.
// * The name of Baris Sumengen may not be
// used to endorse or promote products
// derived from this software without
// specific prior written permission.
//
//THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT
//HOLDERS AND CONTRIBUTORS "AS IS" AND ANY
//EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT
//NOT LIMITED TO, THE IMPLIED WARRANTIES OF
//MERCHANTABILITY AND FITNESS FOR A PARTICULAR
//PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE
//CONTRIBUTORS BE LIABLE FOR ANY
//DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
//EXEMPLARY, OR CONSEQUENTIAL DAMAGES
//(INCLUDING, BUT NOT LIMITED TO, PROCUREMENT
//OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
//DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
//HOWEVER CAUSED AND ON ANY THEORY OF
//LIABILITY, WHETHER IN CONTRACT, STRICT
//LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR
//OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
//OF THIS SOFTWARE, EVEN IF ADVISED OF THE
//POSSIBILITY OF SUCH DAMAGE.
#include "./ImageSegmentation.h"
#include <queue>
#include <vector>
namespace ImageProcessing
{
// Gaussian Filter
Matrix<float> Gaussian2D(int side, float sigma_x, float angle, float ratio)
{
Matrix<float> temp(2*side+1, 2*side+1);
float sigma_y = sigma_x*ratio;
double sin_angle = sin(angle*PI/180);
double cos_angle = cos(angle*PI/180);
int sample;
double d;
if (sigma_x < 2.0 || sigma_y < 2.0)
{
// use denser sample for better filter quality
sample = 5;
d = 1.0/(2.0*sample+1.0);
}
else
{
sample = 0;
d = 1.0;
}
double T = (double) (2*sample+1)*(2*sample+1);
for (int i=0;i<2*side+1;i++)
{
for (int j=0;j<2*side+1;j++)
{
double sum_G = 0.0;
double y = (double) (j-side)-d*sample-d;
for (int sx=0;sx<2*sample+1;sx++)
{
y += d;
double x = (double) (i-side)-d*sample-d;
for (int sy=0;sy<2*sample+1;sy++)
{
x += d;
double x_new = x*cos_angle + y*sin_angle;
double y_new = -x*sin_angle + y*cos_angle;
double g = 1.0/(2.0*PI*sigma_x*sigma_y)*exp(-(x_new*x_new/(sigma_x*sigma_x)+y_new*y_new/(sigma_y*sigma_y))/2.0);
sum_G += g;
}
}
temp.ElemNC(j,i) = (float)(sum_G/T);
}
}
return temp;
}
Matrix<double> Gaussian2D(int side, double sigma_x, double angle, double ratio)
{
Matrix<double> temp(2*side+1, 2*side+1);
double sigma_y = sigma_x*ratio;
double sin_angle = sin(angle*PI/180);
double cos_angle = cos(angle*PI/180);
int sample;
double d;
if (sigma_x < 2.0 || sigma_y < 2.0)
{
// use denser sample for better filter quality
sample = 5;
d = 1.0/(2.0*sample+1.0);
}
else
{
sample = 0;
d = 1.0;
}
double T = (double) (2*sample+1)*(2*sample+1);
for (int i=0;i<2*side+1;i++)
{
for (int j=0;j<2*side+1;j++)
{
double sum_G = 0.0;
double y = (double) (j-side)-d*sample-d;
for (int sx=0;sx<2*sample+1;sx++)
{
y += d;
double x = (double) (i-side)-d*sample-d;
for (int sy=0;sy<2*sample+1;sy++)
{
x += d;
double x_new = x*cos_angle + y*sin_angle;
double y_new = -x*sin_angle + y*cos_angle;
double g = 1.0/(2.0*PI*sigma_x*sigma_y)*exp(-(x_new*x_new/(sigma_x*sigma_x)+y_new*y_new/(sigma_y*sigma_y))/2.0);
sum_G += g;
}
}
temp.ElemNC(j,i) = sum_G/T;
}
}
return temp;
}
// First Directional Derivative of Gaussian Filter
Matrix<float> FDGaussian2D(int side, float sigma_x, float angle, float ratio)
{
Matrix<float> temp(2*side+1, 2*side+1);
float sigma_y = sigma_x*ratio;
double sin_angle = sin(angle*PI/180);
double cos_angle = cos(angle*PI/180);
int sample;
double d;
if (sigma_x < 2.0 || sigma_y < 2.0)
{
// use denser sample for better filter quality
sample = 5;
d = 1.0/(2.0*sample+1.0);
}
else
{
sample = 0;
d = 1.0;
}
double T = (double) (2*sample+1)*(2*sample+1);
for (int i=0;i<2*side+1;i++)
{
for (int j=0;j<2*side+1;j++)
{
double sum_G = 0.0;
double y = (double) (j-side)-d*sample-d;
for (int sx=0;sx<2*sample+1;sx++)
{
y += d;
double x = (double) (i-side)-d*sample-d;
for (int sy=0;sy<2*sample+1;sy++)
{
x += d;
double x_new = x*cos_angle + y*sin_angle;
double y_new = -x*sin_angle + y*cos_angle;
double g = -x_new/(2.0*PI*sigma_x*sigma_y*sigma_x*sigma_x)*exp(-(x_new*x_new/(sigma_x*sigma_x)+y_new*y_new/(sigma_y*sigma_y))/2.0);
sum_G += g;
}
}
temp.ElemNC(j,i) = (float)(sum_G/T);
}
}
return temp;
}
Matrix<double> FDGaussian2D(int side, double sigma_x, double angle, double ratio)
{
Matrix<double> temp(2*side+1, 2*side+1);
double sigma_y = sigma_x*ratio;
double sin_angle = sin(angle*PI/180);
double cos_angle = cos(angle*PI/180);
int sample;
double d;
if (sigma_x < 2.0 || sigma_y < 2.0)
{
// use denser sample for better filter quality
sample = 5;
d = 1.0/(2.0*sample+1.0);
}
else
{
sample = 0;
d = 1.0;
}
double T = (double) (2*sample+1)*(2*sample+1);
for (int i=0;i<2*side+1;i++)
{
for (int j=0;j<2*side+1;j++)
{
double sum_G = 0.0;
double y = (double) (j-side)-d*sample-d;
for (int sx=0;sx<2*sample+1;sx++)
{
y += d;
double x = (double) (i-side)-d*sample-d;
for (int sy=0;sy<2*sample+1;sy++)
{
x += d;
double x_new = x*cos_angle + y*sin_angle;
double y_new = -x*sin_angle + y*cos_angle;
double g = -x_new/(2.0*PI*sigma_x*sigma_y*sigma_x*sigma_x)*exp(-(x_new*x_new/(sigma_x*sigma_x)+y_new*y_new/(sigma_y*sigma_y))/2.0);
sum_G += g;
}
}
temp.ElemNC(j,i) = sum_G/T;
}
}
return temp;
}
// Second Directional Derivative of Gaussian Filter
Matrix<float> SDGaussian2D(int side, float sigma_x, float angle, float ratio)
{
Matrix<float> temp(2*side+1, 2*side+1);
float sigma_y = sigma_x*ratio;
double sin_angle = sin(angle*PI/180);
double cos_angle = cos(angle*PI/180);
int sample;
double d;
if (sigma_x < 2.0 || sigma_y < 2.0)
{
// use denser sample for better filter quality
sample = 5;
d = 1.0/(2.0*sample+1.0);
}
else
{
sample = 0;
d = 1.0;
}
double T = (double) (2*sample+1)*(2*sample+1);
for (int i=0;i<2*side+1;i++)
{
for (int j=0;j<2*side+1;j++)
{
double sum_G = 0.0;
double y = (double) (j-side)-d*sample-d;
for (int sx=0;sx<2*sample+1;sx++)
{
y += d;
double x = (double) (i-side)-d*sample-d;
for (int sy=0;sy<2*sample+1;sy++)
{
x += d;
double x_new = x*cos_angle + y*sin_angle;
double y_new = -x*sin_angle + y*cos_angle;
double g = (x_new*x_new-sigma_x*sigma_x)/(2.0*PI*sigma_x*sigma_y*sigma_x*sigma_x*sigma_x*sigma_x)*exp(-(x_new*x_new/(sigma_x*sigma_x)+y_new*y_new/(sigma_y*sigma_y))/2.0);
sum_G += g;
}
}
temp.ElemNC(j,i) = (float)(sum_G/T);
}
}
return temp;
}
Matrix<double> SDGaussian2D(int side, double sigma_x, double angle, double ratio)
{
Matrix<double> temp(2*side+1, 2*side+1);
double sigma_y = sigma_x*ratio;
double sin_angle = sin(angle*PI/180);
double cos_angle = cos(angle*PI/180);
int sample;
double d;
if (sigma_x < 2.0 || sigma_y < 2.0)
{
// use denser sample for better filter quality
sample = 5;
d = 1.0/(2.0*sample+1.0);
}
else
{
sample = 0;
d = 1.0;
}
double T = (double) (2*sample+1)*(2*sample+1);
for (int i=0;i<2*side+1;i++)
{
for (int j=0;j<2*side+1;j++)
{
double sum_G = 0.0;
double y = (double) (j-side)-d*sample-d;
for (int sx=0;sx<2*sample+1;sx++)
{
y += d;
double x = (double) (i-side)-d*sample-d;
for (int sy=0;sy<2*sample+1;sy++)
{
x += d;
double x_new = x*cos_angle + y*sin_angle;
double y_new = -x*sin_angle + y*cos_angle;
double g = (x_new*x_new-sigma_x*sigma_x)/(2.0*PI*sigma_x*sigma_y*sigma_x*sigma_x*sigma_x*sigma_x)*exp(-(x_new*x_new/(sigma_x*sigma_x)+y_new*y_new/(sigma_y*sigma_y))/2.0);
sum_G += g;
}
}
temp.ElemNC(j,i) = sum_G/T;
}
}
return temp;
}
// Laplacian of Gaussian Filter
Matrix<float> LOG(int side, float sigma_x, float angle, float ratio)
{
Matrix<float> temp(2*side+1, 2*side+1);
float sigma_y = sigma_x*ratio;
double sin_angle = sin(angle*PI/180);
double cos_angle = cos(angle*PI/180);
int sample;
double d;
if (sigma_x < 2.0 || sigma_y < 2.0)
{
// use denser sample for better filter quality
sample = 5;
d = 1.0/(2.0*sample+1.0);
}
else
{
sample = 0;
d = 1.0;
}
double T = (double) (2*sample+1)*(2*sample+1);
for (int i=0;i<2*side+1;i++)
{
for (int j=0;j<2*side+1;j++)
{
double sum_G = 0.0;
double y = (double) (j-side)-d*sample-d;
for (int sx=0;sx<2*sample+1;sx++)
{
y += d;
double x = (double) (i-side)-d*sample-d;
for (int sy=0;sy<2*sample+1;sy++)
{
x += d;
double x_new = x*cos_angle + y*sin_angle;
double y_new = -x*sin_angle + y*cos_angle;
double g = (x_new*x_new+y_new*y_new-sigma_x*sigma_x-sigma_y*sigma_y)
/(2.0*PI*sigma_x*sigma_y*(sigma_x*sigma_x*sigma_x*sigma_x+sigma_y*sigma_y*sigma_y*sigma_y))
*exp(-(x_new*x_new/(sigma_x*sigma_x)+y_new*y_new/(sigma_y*sigma_y))/2.0);
sum_G += g;
}
}
temp.ElemNC(j,i) = (float)(sum_G/T);
}
}
return temp;
}
Matrix<double> LOG(int side, double sigma_x, double angle, double ratio)
{
Matrix<double> temp(2*side+1, 2*side+1);
double sigma_y = sigma_x*ratio;
double sin_angle = sin(angle*PI/180);
double cos_angle = cos(angle*PI/180);
int sample;
double d;
if (sigma_x < 2.0 || sigma_y < 2.0)
{
// use denser sample for better filter quality
sample = 5;
d = 1.0/(2.0*sample+1.0);
}
else
{
sample = 0;
d = 1.0;
}
double T = (double) (2*sample+1)*(2*sample+1);
for (int i=0;i<2*side+1;i++)
{
for (int j=0;j<2*side+1;j++)
{
double sum_G = 0.0;
double y = (double) (j-side)-d*sample-d;
for (int sx=0;sx<2*sample+1;sx++)
{
y += d;
double x = (double) (i-side)-d*sample-d;
for (int sy=0;sy<2*sample+1;sy++)
{
x += d;
double x_new = x*cos_angle + y*sin_angle;
double y_new = -x*sin_angle + y*cos_angle;
double g = (x_new*x_new+y_new*y_new-sigma_x*sigma_x-sigma_y*sigma_y)
/(2.0*PI*sigma_x*sigma_y*(sigma_x*sigma_x*sigma_x*sigma_x+sigma_y*sigma_y*sigma_y*sigma_y))
*exp(-(x_new*x_new/(sigma_x*sigma_x)+y_new*y_new/(sigma_y*sigma_y))/2.0);
sum_G += g;
}
}
temp.ElemNC(j,i) = sum_G/T;
}
}
return temp;
}
// Difference of offset Gaussians filter
Matrix<float> DOOG2D(int side, float sigma_x, float offset, float angle, float ratio)
{
Matrix<float> temp(2*side+1, 2*side+1);
float sigma_y = sigma_x*ratio;
double sin_angle = sin(angle*PI/180);
double cos_angle = cos(angle*PI/180);
int sample;
double d;
if (sigma_x < 2.0 || sigma_y < 2.0)
{
// use denser sample for better filter quality
sample = 5;
d = 1.0/(2.0*sample+1.0);
}
else
{
sample = 0;
d = 1.0;
}
?? 快捷鍵說明
復(fù)制代碼
Ctrl + C
搜索代碼
Ctrl + F
全屏模式
F11
切換主題
Ctrl + Shift + D
顯示快捷鍵
?
增大字號(hào)
Ctrl + =
減小字號(hào)
Ctrl + -