Smooth an image with Total Variation (TV) regularization (standard ROF model, anisotropic, second order)
Usage: example_total_variation
reads a given parameter file (see tv.par, aniso_tv.par, secondo_oder_tv.par in the examples subfolder) and performes a smoothing with TV regularization with the parameters given in this file.
#include <iostream>
#include <stdio.h>
#include <string.h>
#include <vigra/multi_array.hxx>
#include <vigra/stdimagefunctions.hxx>
#include <vigra/tv_filter.hxx>
#ifdef _MSC_VER
#define strcasecmp _stricmp
#endif
void read_string(FILE *file,const char *name, char* out){
char dummy[80];
int s=fscanf(file,"%s %s ",dummy,out);
if (s==EOF || s<2){
std::cout<<"Could not read from file.\n";
exit(0);
}
else if (strcasecmp(name,dummy)!=0){
std::cout<<"Found parameter "<<dummy<<" instead of "<<name<<"\n";
exit(0);
}
else
std::cout<<"Parameter "<<name<<" = "<<out<<std::endl;
}
float read_value(FILE *file,const char *name){
char dummy[80];
char dummy2[80];
float f=0;
int s=fscanf(file,"%s %s ",dummy,dummy2);
if (s==EOF || s<2){
std::cout<<"Could not read from file.\n";
exit(0);
}
else if (strcasecmp(name,dummy)!=0){
std::cout<<"Found parameter "<<dummy<<" instead of "<<name<<"\n";
exit(0);
}
else{
f=atof(dummy2);
std::cout<<"Parameter "<<name<<" = "<<f<<std::endl;
}
return f;
}
int main(int argc, char ** argv)
{
using namespace multi_math;
if(argc <2)
{
std::cout << "Usage: " << argv[0] << " parameterfile" << std::endl;
std::cout <<
"(supported formats for images: " <<
impexListFormats() <<
")" << std::endl;
return 1;
}
try
{
char infile[80],outfile[80];
double alpha0=0.01*
sqrt(255.0),beta0=0.1*
sqrt(255.0),sigma=0.1,rho=.5,K=10,eps=0,gamma0=0;
int mode,inner_steps=200,outer_steps=5,write_steps=0;
FILE *file=fopen(argv[1],"r");
if (!file){
std::cout<<"Cannot open file "<<argv[1]<<std::endl;
return 1;
}
read_string(file,"input_file",infile);
read_string(file,"output_file",outfile);
mode=(int)read_value(file,"mode");
switch(mode){
case 1:
case 2:
inner_steps=1000;
alpha0=read_value(file,"alpha");
eps=read_value(file,"epsilon");
break;
case 3:
outer_steps=(int)read_value(file,"outer_steps");
inner_steps=(int)read_value(file,"inner_steps");
alpha0=read_value(file,"alpha");
beta0=read_value(file,"beta");
sigma=read_value(file,"sigma");
rho=read_value(file,"rho");
K=read_value(file,"edge_factor");
write_steps=(int)read_value(file,"write_outer_steps");
break;
case 4:
outer_steps=(int)read_value(file,"outer_steps");
inner_steps=(int)read_value(file,"inner_steps");
alpha0=read_value(file,"alpha");
beta0=read_value(file,"beta");
gamma0=read_value(file,"gamma");
sigma=read_value(file,"sigma");
rho=read_value(file,"rho");
K=read_value(file,"edge_factor");
write_steps=(int)read_value(file,"write_outer_steps");
break;
default:
std::cout<<"Unknown mode "<<mode<<std::endl;
}
int stretch=(int)read_value(file,"histogram_stretch");
ImageImportInfo info(infile);
if(info.isGrayscale())
{
MultiArray<2,double> data(info.shape());
MultiArray<2,double> out(info.shape());
MultiArray<2,double> weight(info.shape());
for (int y=0;y<data.shape(1);y++){
for (int x=0;x<data.shape(0);x++){
weight(x,y)=1;
}
}
data=data*(1/255.);
switch(mode){
case 1:
std::cout<<"Standard TV filter"<<std::endl;
break;
case 2:
std::cout<<"Weighted TV filter"<<std::endl;
break;
case 3:{
std::cout<<"Anisotropic TV filter"<<std::endl;
MultiArray<2,double> phi(info.shape());
MultiArray<2,double> alpha(info.shape());
MultiArray<2,double> beta(info.shape());
out=data;
for (int i=0;i<outer_steps;i++){
std::cout<<"outer step "<<i<<"\n";
if(write_steps){
char dummy[80];
sprintf(dummy,"output_step_%03d.png",i);
std::cout<<"Writing temp file\n";
if (stretch)
else{
MultiArray<2,unsigned char> buffer(info.shape());
buffer=max(min(out*255+0.5,0.),255.);
}
}
}
break;
}
case 4:{
std::cout<<"Second order TV filter"<<std::endl;
MultiArray<2,double> phi(info.shape());
MultiArray<2,double> alpha(info.shape());
MultiArray<2,double> beta(info.shape());
MultiArray<2,double>
gamma(info.shape());
MultiArray<2,double> xedges(info.shape());
MultiArray<2,double> yedges(info.shape());
for (int y=0;y<data.shape(1);y++){
for (int x=0;x<data.shape(0);x++){
xedges(x,y)=1;
yedges(x,y)=1;
}
}
out=data;
for (int i=0;i<outer_steps;i++){
std::cout<<"outer step "<<i<<"\n";
if(write_steps){
char dummy[80];
sprintf(dummy,"output_step_%03d.png",i);
std::cout<<"Writing temp file\n";
if (stretch)
else{
MultiArray<2,unsigned char> buffer(info.shape());
buffer=max(min(out*255+0.5,0.),255.);
}
}
}
}
}
if (stretch)
else{
MultiArray<2,unsigned char> buffer(info.shape());
buffer=min(max(out*255.+0.5,0.),255.);
}
}
else
{
std::cout<<"Color images are currently not supported !\n";
}
}
catch (std::exception & e)
{
std::cout << e.what() << std::endl;
return 1;
}
return 0;
}
void importImage(...)
Read an image from a file.
void exportImage(...)
Write an image to a file.
image import and export functions
std::string impexListFormats()
List the image formats VIGRA can read and write.
SquareRootTraits< FixedPoint< IntBits, FracBits > >::SquareRootResult sqrt(FixedPoint< IntBits, FracBits > v)
square root.
Definition: fixedpoint.hxx:616
double gamma(double x)
The gamma function.
Definition: mathutil.hxx:1587
void totalVariationFilter(...)
Performs standard Total Variation Regularization.
void secondOrderTotalVariationFilter(...)
Performs Anisotropic Total Variation Regularization.
void getAnisotropy(...)
Sets up directional data for anisotropic regularization.
void anisotropicTotalVariationFilter(...)
Performs Anisotropic Total Variation Regularization.