levelsetmethods.h

图像分割算法

//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.


#pragma once
#ifndef LEVELSETMETHODS_H
#define LEVELSETMETHODS_H

#include <cmath>

#include "cimpl.h"
using namespace CIMPL;

#include "Core.h"
#include "Functions.h"
using namespace MathCore;


namespace LevelSetMethods
{

	Matrix<int> ExtractCurve(Matrix<float> phi);
	Matrix<int> ExtractCurve(Matrix<double> phi);

	Matrix<float>& ExtendConst2D(Matrix<float>& phi, int size);
	Matrix<float>& Extend2D(Matrix<float>& phi, int size);
	Matrix<double>& Extend2D(Matrix<double>& phi, int size);

	void Evolve2DNormalENO1(Matrix<float>& phi, float dx, float dy, Matrix<float>& Vn, Matrix<float>& delta, Matrix<float>& H1Abs, Matrix<float>& H2Abs);
	void Evolve2DNormalENO2(Matrix<float>& phi, float dx, float dy, Matrix<float>& Vn, Matrix<float>& delta, Matrix<float>& H1Abs, Matrix<float>& H2Abs);
	void Evolve2DNormalENO3(Matrix<float>& phi, float dx, float dy, Matrix<float>& Vn, Matrix<float>& delta, Matrix<float>& H1Abs, Matrix<float>& H2Abs);
	void Evolve2DNormalWENO(Matrix<float>& phi, float dx, float dy, Matrix<float>& Vn, Matrix<float>& delta, Matrix<float>& H1Abs, Matrix<float>& H2Abs);

	void Evolve2DNormalENO1(Matrix<double>& phi, double dx, double dy, Matrix<double>& Vn, Matrix<double>& delta, Matrix<double>& H1Abs, Matrix<double>& H2Abs);
	void Evolve2DNormalENO2(Matrix<double>& phi, double dx, double dy, Matrix<double>& Vn, Matrix<double>& delta, Matrix<double>& H1Abs, Matrix<double>& H2Abs);
	void Evolve2DNormalENO3(Matrix<double>& phi, double dx, double dy, Matrix<double>& Vn, Matrix<double>& delta, Matrix<double>& H1Abs, Matrix<double>& H2Abs);
	void Evolve2DNormalWENO(Matrix<double>& phi, double dx, double dy, Matrix<double>& Vn, Matrix<double>& delta, Matrix<double>& H1Abs, Matrix<double>& H2Abs);

	
	void Evolve2DVectorENO1(Matrix<float>& phi, float dx, float dy, Matrix<float>& u, Matrix<float>& v, Matrix<float>& delta);
	void Evolve2DVectorENO2(Matrix<float>& phi, float dx, float dy, Matrix<float>& u, Matrix<float>& v, Matrix<float>& delta);
	void Evolve2DVectorENO3(Matrix<float>& phi, float dx, float dy, Matrix<float>& u, Matrix<float>& v, Matrix<float>& delta);
	void Evolve2DVectorWENO(Matrix<float>& phi, float dx, float dy, Matrix<float>& u, Matrix<float>& v, Matrix<float>& delta);

	void Evolve2DVectorENO1(Matrix<double>& phi, double dx, double dy, Matrix<double>& u, Matrix<double>& v, Matrix<double>& delta);
	void Evolve2DVectorENO2(Matrix<double>& phi, double dx, double dy, Matrix<double>& u, Matrix<double>& v, Matrix<double>& delta);
	void Evolve2DVectorENO3(Matrix<double>& phi, double dx, double dy, Matrix<double>& u, Matrix<double>& v, Matrix<double>& delta);
	void Evolve2DVectorWENO(Matrix<double>& phi, double dx, double dy, Matrix<double>& u, Matrix<double>& v, Matrix<double>& delta);


	void Evolve2DNormalVectorENO1SignedDistance(Matrix<float>& phi, float dx, float dy, Matrix<float>& u, Matrix<float>& v, Matrix<float>& Vn, Matrix<float>& delta, Matrix<float>& H1Abs, Matrix<float>& H2Abs);
	void Evolve2DNormalVectorENO2SignedDistance(Matrix<float>& phi, float dx, float dy, Matrix<float>& u, Matrix<float>& v, Matrix<float>& Vn, Matrix<float>& delta, Matrix<float>& H1Abs, Matrix<float>& H2Abs);
	void Evolve2DNormalVectorENO3SignedDistance(Matrix<float>& phi, float dx, float dy, Matrix<float>& u, Matrix<float>& v, Matrix<float>& Vn, Matrix<float>& delta, Matrix<float>& H1Abs, Matrix<float>& H2Abs);
	void Evolve2DNormalVectorWENOSignedDistance(Matrix<float>& phi, float dx, float dy, Matrix<float>& u, Matrix<float>& v, Matrix<float>& Vn, Matrix<float>& delta, Matrix<float>& H1Abs, Matrix<float>& H2Abs);

	void Evolve2DNormalVectorENO1SignedDistance(Matrix<double>& phi, double dx, double dy, Matrix<double>& u, Matrix<double>& v, Matrix<double>& Vn, Matrix<double>& delta, Matrix<double>& H1Abs, Matrix<double>& H2Abs);
	void Evolve2DNormalVectorENO2SignedDistance(Matrix<double>& phi, double dx, double dy, Matrix<double>& u, Matrix<double>& v, Matrix<double>& Vn, Matrix<double>& delta, Matrix<double>& H1Abs, Matrix<double>& H2Abs);
	void Evolve2DNormalVectorENO3SignedDistance(Matrix<double>& phi, double dx, double dy, Matrix<double>& u, Matrix<double>& v, Matrix<double>& Vn, Matrix<double>& delta, Matrix<double>& H1Abs, Matrix<double>& H2Abs);
	void Evolve2DNormalVectorWENOSignedDistance(Matrix<double>& phi, double dx, double dy, Matrix<double>& u, Matrix<double>& v, Matrix<double>& Vn, Matrix<double>& delta, Matrix<double>& H1Abs, Matrix<double>& H2Abs);


	void Evolve2DKappa(Matrix<float>& phi, float dx, float dy, float dx2, float dy2, Matrix<float>& b, Matrix<float>& delta);
	void Evolve2DKappa(Matrix<double>& phi, double dx, double dy, double dx2, double dy2, Matrix<double>& b, Matrix<double>& delta);


	float Getdt2DNormal(float alpha, float dx, float dy, Matrix<float>& H1Abs, Matrix<float>& H2Abs);
	float Getdt2DNormalVector(float alpha, float dx, float dy, Matrix<float>& H1Abs, Matrix<float>& H2Abs);
	float Getdt2DVector(float alpha, float dx, float dy, Matrix<float>& u, Matrix<float>& v);
	float Getdt2DKappa(float alpha, float dx2, float dy2, Matrix<float>& b);
	float Getdt2DVectorKappa(float alpha, float dx, float dy, float dx2, float dy2, Matrix<float>& u, Matrix<float>& v, Matrix<float>& b);
	float Getdt2DNormalKappa(float alpha, float dx, float dy, float dx2, float dy2, Matrix<float>& H1Abs, Matrix<float>& H2Abs, Matrix<float>& b);
	float Getdt2DNormalVectorKappa(float alpha, float dx, float dy, float dx2, float dy2, Matrix<float>& H1Abs, Matrix<float>& H2Abs, Matrix<float>& b);

	double Getdt2DNormal(double alpha, double dx, double dy, Matrix<double>& H1Abs, Matrix<double>& H2Abs);
	double Getdt2DNormalVector(double alpha, double dx, double dy, Matrix<double>& H1Abs, Matrix<double>& H2Abs);
	double Getdt2DVector(double alpha, double dx, double dy, Matrix<double>& u, Matrix<double>& v);
	double Getdt2DKappa(double alpha, double dx2, double dy2, Matrix<double>& b);
	double Getdt2DVectorKappa(double alpha, double dx, double dy, double dx2, double dy2, Matrix<double>& u, Matrix<double>& v, Matrix<double>& b);
	double Getdt2DNormalKappa(double alpha, double dx, double dy, double dx2, double dy2, Matrix<double>& H1Abs, Matrix<double>& H2Abs, Matrix<double>& b);
	double Getdt2DNormalVectorKappa(double alpha, double dx, double dy, double dx2, double dy2, Matrix<double>& H1Abs, Matrix<double>& H2Abs, Matrix<double>& b);



	void SubtractPhi2D(Matrix<float>& phi, Matrix<float>& delta, float dt);
	void AddPhi2D(Matrix<float>& phi, Matrix<float>& delta, float dt);

	void SubtractPhi2D(Matrix<double>& phi, Matrix<double>& delta, double dt);
	void AddPhi2D(Matrix<double>& phi, Matrix<double>& delta, double dt);



	float KappaAbsPhi2D(float phi_i_j, float phi_ip1_j, float phi_ip1_jp1, 
		float phi_i_jp1, float phi_im1_jp1, float phi_im1_j, float phi_im1_jm1, 
		float phi_i_jm1, float phi_ip1_jm1, float dx, float dy, float dx2, float dy2);

	double KappaAbsPhi2D(double phi_i_j, double phi_ip1_j, double phi_ip1_jp1, 
		double phi_i_jp1, double phi_im1_jp1, double phi_im1_j, double phi_im1_jm1, 
		double phi_i_jm1, double phi_ip1_jm1, double dx, double dy, double dx2, double dy2);

	
	void ReinitSignedDistanceENO1(Matrix<float>& phi, float dx, float dy, float alpha, int iterations);
	void ReinitSignedDistanceENO2(Matrix<float>& phi, float dx, float dy, float alpha, int iterations);
	void ReinitSignedDistanceENO3(Matrix<float>& phi, float dx, float dy, float alpha, int iterations);
	void ReinitSignedDistanceWENO(Matrix<float>& phi, float dx, float dy, float alpha, int iterations);

	void ReinitSignedDistanceENO1(Matrix<double>& phi, double dx, double dy, double alpha, int iterations);
	void ReinitSignedDistanceENO2(Matrix<double>& phi, double dx, double dy, double alpha, int iterations);
	void ReinitSignedDistanceENO3(Matrix<double>& phi, double dx, double dy, double alpha, int iterations);
	void ReinitSignedDistanceWENO(Matrix<double>& phi, double dx, double dy, double alpha, int iterations);



	float SelectDerNormal(float phi_x_m, float phi_x_p, float Vn);
	float SelectDerNormalVectorSD(float phi_x_m, float phi_x_p, float u, float Vn);

	double SelectDerNormal(double phi_x_m, double phi_x_p, double Vn);
	double SelectDerNormalVectorSD(double phi_x_m, double phi_x_p, double u, double Vn);


	float DerENO1Minus(float v1, float v2, float dx);
	float DerENO1Plus(float v1, float v2, float dx);
	float DerENO2Minus(float v1, float v2, float v3, float v4, float dx);
	float DerENO2Plus(float v1, float v2, float v3, float v4, float dx);
	float DerENO3Minus(float v1, float v2, float v3, float v4, float v5, float v6, float dx);
	float DerENO3Plus(float v1, float v2, float v3, float v4, float v5, float v6, float dx);
	float DerWENOMinus(float v1, float v2, float v3, float v4, float v5, float v6, float dx);
	float DerWENOPlus(float v1, float v2, float v3, float v4, float v5, float v6, float dx);

	double DerENO1Minus(double v1, double v2, double dx);
	double DerENO1Plus(double v1, double v2, double dx);
	double DerENO2Minus(double v1, double v2, double v3, double v4, double dx);
	double DerENO2Plus(double v1, double v2, double v3, double v4, double dx);
	double DerENO3Minus(double v1, double v2, double v3, double v4, double v5, double v6, double dx);
	double DerENO3Plus(double v1, double v2, double v3, double v4, double v5, double v6, double dx);
	double DerWENOMinus(double v1, double v2, double v3, double v4, double v5, double v6, double dx);
	double DerWENOPlus(double v1, double v2, double v3, double v4, double v5, double v6, double dx);


	//inline int Sign(float f)
	//{
	//	if(f>0)
	//	{
	//		return 1;
	//	}
	//	else if(f < 0)
	//	{
	//		return -1;
	//	}
	//	else
	//	{
	//		return 0;
	//	}
	//}

	











#include "./LevelSetMethods.inl"

}; //namespace



#endif