[ VIGRA Homepage | Function Index | Class Index | Namespaces | File List | Main Page ]
total_variation.cxx | ![]() |
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.
/************************************************************************/ /* */ /* Copyright 2012 by Frank Lenzen & */ /* Ullrich Koethe */ /* */ /* This file is part of the VIGRA computer vision library. */ /* The VIGRA Website is */ /* http://hci.iwr.uni-heidelberg.de/vigra/ */ /* Please direct questions, bug reports, and contributions to */ /* ullrich.koethe@iwr.uni-heidelberg.de or */ /* vigra@informatik.uni-hamburg.de */ /* */ /* Permission is hereby granted, free of charge, to any person */ /* obtaining a copy of this software and associated documentation */ /* files (the "Software"), to deal in the Software without */ /* restriction, including without limitation the rights to use, */ /* copy, modify, merge, publish, distribute, sublicense, and/or */ /* sell copies of the Software, and to permit persons to whom the */ /* Software is furnished to do so, subject to the following */ /* conditions: */ /* */ /* The above copyright notice and this permission notice shall be */ /* included in all copies or substantial portions of the */ /* Software. */ /* */ /* THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND */ /* EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES */ /* OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND */ /* NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT */ /* HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, */ /* WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING */ /* FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR */ /* OTHER DEALINGS IN THE SOFTWARE. */ /* */ /************************************************************************/ #include <iostream> #include <stdio.h> #include <string.h> #include <vigra/stdimage.hxx> #include <vigra/stdimagefunctions.hxx> #include <vigra/impex.hxx> #include <vigra/multi_array.hxx> #include <vigra/tv_filter.hxx> #ifdef _MSC_VER #define strcasecmp _stricmp #endif using namespace vigra; // Read a string from file <file>. 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; } // read a flota-value from file <file> 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: " << vigra::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"); // open parameter file if (!file){ std::cout<<"Cannot open file "<<argv[1]<<std::endl; return 1; } read_string(file,"input_file",infile); // read name of input file (an image) read_string(file,"output_file",outfile); // read name of output file (an image) mode=(int)read_value(file,"mode"); // read mode: 1- standard Total Variation // 2- Anisotropic Total Variation // 3- Anisotropic Total Variation with Higher Order term // further parameters depend on mode switch(mode){ case 1: case 2: inner_steps=1000; // read number of iteration steps alpha0=read_value(file,"alpha"); // read alpha eps=read_value(file,"epsilon"); // read eps break; case 3: outer_steps=(int)read_value(file,"outer_steps"); // read number of outer iteration steps inner_steps=(int)read_value(file,"inner_steps"); // read number of inner iteration steps alpha0=read_value(file,"alpha"); // read alpha beta0=read_value(file,"beta"); // read beta sigma=read_value(file,"sigma"); // read sigma rho=read_value(file,"rho"); // read rho K=read_value(file,"edge_factor"); // read parameter K (edge sensitivity) write_steps=(int)read_value(file,"write_outer_steps"); //read flag for writing intermediate result after each outer iteration break; case 4: outer_steps=(int)read_value(file,"outer_steps"); // read number of outer iteration steps inner_steps=(int)read_value(file,"inner_steps"); // read number of inner iteration steps alpha0=read_value(file,"alpha"); // read alpha beta0=read_value(file,"beta"); // read beta gamma0=read_value(file,"gamma"); // read gamma sigma=read_value(file,"sigma"); // read sigma rho=read_value(file,"rho"); // read rho K=read_value(file,"edge_factor");// read parameter K (edge sensitivity) write_steps=(int)read_value(file,"write_outer_steps");//read flag for writing intermediate result after each outer iteration break; default: std::cout<<"Unknown mode "<<mode<<std::endl; } int stretch=(int)read_value(file,"histogram_stretch"); // flag for histogram stretching: 0 -no, 1 -yes vigra::ImageImportInfo info(infile); // get information about input image if(info.isGrayscale()) { MultiArray<2,double> data(Shape2(info.width(),info.height())); MultiArray<2,double> out(Shape2(info.width(),info.height())); MultiArray<2,double> weight(Shape2(info.width(),info.height())); for (int y=0;y<data.shape(1);y++){ for (int x=0;x<data.shape(0);x++){ weight(x,y)=1; // set weight to 1 (needed for anisotropicTotalVariationFilter and // higherOrderTotalVariationFilter } } importImage(info, destImage(data)); // read image data data=data*(1/255.); // scale to [0,1] (smoothing parameters are tuned to this scaling) //------------------------------------ // Choose specific TV Filter //------------------------------------ switch(mode){ case 1: std::cout<<"Standard TV filter"<<std::endl; totalVariationFilter(data,out,alpha0,inner_steps,eps); break; case 2: std::cout<<"Weighted TV filter"<<std::endl; totalVariationFilter(data,out,weight,alpha0,inner_steps,eps); break; case 3:{ std::cout<<"Anisotropic TV filter"<<std::endl; MultiArray<2,double> phi(Shape2(info.width(),info.height())); MultiArray<2,double> alpha(Shape2(info.width(),info.height())); MultiArray<2,double> beta(Shape2(info.width(),info.height())); out=data; //'data' serves as initial value for (int i=0;i<outer_steps;i++){ // Outer loop to update anisotropic data std::cout<<"outer step "<<i<<"\n"; getAnisotropy(out,phi,alpha,beta,alpha0,beta0,sigma,rho,K); // get anisotropic data anisotropicTotalVariationFilter(data,weight,phi,alpha,beta,out,inner_steps); //perform smoothing if(write_steps){ char dummy[80]; sprintf(dummy,"output_step_%03d.png",i); std::cout<<"Writing temp file\n"; if (stretch) exportImage(srcImageRange(out), vigra::ImageExportInfo(dummy)); //exportImage performes histogramm stretching else{ MultiArray<2,unsigned char> buffer(Shape2(info.width(),info.height())); buffer=max(min(out*255+0.5,0.),255.); //scaling back to [0,255], clipping, cast to uint8 exportImage(srcImageRange(buffer), vigra::ImageExportInfo(dummy)); } } } break; } case 4:{ std::cout<<"Second order TV filter"<<std::endl; MultiArray<2,double> phi(Shape2(info.width(),info.height())); MultiArray<2,double> alpha(Shape2(info.width(),info.height())); MultiArray<2,double> beta(Shape2(info.width(),info.height())); MultiArray<2,double> gamma(Shape2(info.width(),info.height())); MultiArray<2,double> xedges(Shape2(info.width(),info.height())); MultiArray<2,double> yedges(Shape2(info.width(),info.height())); for (int y=0;y<data.shape(1);y++){ // set data needed for second order term for (int x=0;x<data.shape(0);x++){ gamma(x,y)=gamma0; xedges(x,y)=1; yedges(x,y)=1; } } out=data; //'data' serves as initial value for (int i=0;i<outer_steps;i++){// Outer loop to update anisotropic data std::cout<<"outer step "<<i<<"\n"; getAnisotropy(out,phi,alpha,beta,alpha0,beta0,sigma,rho,K); // get anisotropic data secondOrderTotalVariationFilter(data,weight,phi,alpha,beta,gamma,xedges,yedges,out,inner_steps);//perform smoothing if(write_steps){ char dummy[80]; sprintf(dummy,"output_step_%03d.png",i); std::cout<<"Writing temp file\n"; if (stretch) exportImage(srcImageRange(out), vigra::ImageExportInfo(dummy)); else{ MultiArray<2,unsigned char> buffer(Shape2(info.width(),info.height())); buffer=max(min(out*255+0.5,0.),255.); //scaling back to [0,255], clipping, cast to uint8 exportImage(srcImageRange(buffer), vigra::ImageExportInfo(dummy)); } } } } } //------------------------------------- if (stretch) exportImage(srcImageRange(out), vigra::ImageExportInfo(outfile)); //write result to file else{ MultiArray<2,unsigned char> buffer(Shape2(info.width(),info.height())); buffer=min(max(out*255.+0.5,0.),255.); //scaling back to [0,255], clipping, cast to uint8 exportImage(srcImageRange(buffer), vigra::ImageExportInfo(outfile)); } } else { std::cout<<"Color images are currently not supported !\n"; } } catch (vigra::StdException & e) { std::cout << e.what() << std::endl; return 1; } return 0; }
© Ullrich Köthe (ullrich.koethe@iwr.uni-heidelberg.de) |
html generated using doxygen and Python
|