Main Page | Modules | Namespace List | Class Hierarchy | Alphabetical List | Class List | Directories | File List | Namespace Members | Class Members | File Members

EMAN::BilateralProcessor Class Reference

Bilateral processing on 2D or 3D volume data. More...

#include <processor.h>

Inheritance diagram for EMAN::BilateralProcessor:

Inheritance graph
[legend]
Collaboration diagram for EMAN::BilateralProcessor:

Collaboration graph
[legend]
List of all members.

Public Member Functions

void process_inplace (EMData *image)
 To process an image in-place.
string get_name () const
 Get the processor's name.
string get_desc () const
 Get the descrition of this specific processor.
TypeDict get_param_types () const
 Get processor parameter information in a dictionary.

Static Public Member Functions

static ProcessorNEW ()

Static Public Attributes

static const string NAME = "bilateral"

Detailed Description

Bilateral processing on 2D or 3D volume data.

Bilateral processing does non-linear weighted averaging processing within a certain window.

Parameters:
distance_sigma means how large the voxel has impact on its neighbors in spatial domain. The larger it is, the more blurry the resulting image.
value_sigma eans how large the voxel has impact on its in range domain. The larger it is, the more blurry the resulting image.
niter how many times to apply this processing on your data.
half_width processing window size = (2 * half_widthh + 1) ^ 3.

Definition at line 3993 of file processor.h.


Member Function Documentation

string EMAN::BilateralProcessor::get_desc  )  const [inline, virtual]
 

Get the descrition of this specific processor.

This function must be overwritten by a subclass.

Returns:
The description of this processor.

Implements EMAN::Processor.

Definition at line 4002 of file processor.h.

04003                 {
04004                         return "Bilateral processing on 2D or 3D volume data. Bilateral processing does non-linear weighted averaging processing within a certain window. ";
04005                 }

string EMAN::BilateralProcessor::get_name  )  const [inline, virtual]
 

Get the processor's name.

Each processor is identified by a unique name.

Returns:
The processor's name.

Implements EMAN::Processor.

Definition at line 3997 of file processor.h.

References NAME.

03998                 {
03999                         return NAME;
04000                 }

TypeDict EMAN::BilateralProcessor::get_param_types  )  const [inline, virtual]
 

Get processor parameter information in a dictionary.

Each parameter has one record in the dictionary. Each record contains its name, data-type, and description.

Returns:
A dictionary containing the parameter info.

Reimplemented from EMAN::Processor.

Definition at line 4012 of file processor.h.

References EMAN::EMObject::FLOAT, EMAN::EMObject::INT, and EMAN::TypeDict::put().

04013                 {
04014                         TypeDict d;
04015                         d.put("distance_sigma", EMObject::FLOAT, "means how large the voxel has impact on its neighbors in spatial domain. The larger it is, the more blurry the resulting image.");
04016                         d.put("value_sigma", EMObject::FLOAT, "means how large the voxel has impact on its in  range domain. The larger it is, the more blurry the resulting image.");
04017                         d.put("niter", EMObject::INT, "how many times to apply this processing on your data.");
04018                         d.put("half_width", EMObject::INT, "processing window size = (2 * half_widthh + 1) ^ 3.");
04019                         return d;
04020                 }

static Processor* EMAN::BilateralProcessor::NEW  )  [inline, static]
 

Definition at line 4007 of file processor.h.

04008                 {
04009                         return new BilateralProcessor();
04010                 }

void BilateralProcessor::process_inplace EMData image  )  [virtual]
 

To process an image in-place.

For those processors which can only be processed out-of-place, override this function to just print out some error message to remind user call the out-of-place version.

Parameters:
image The image to be processed.

Implements EMAN::Processor.

Definition at line 3841 of file processor.cpp.

References EMAN::EMData::get_attr(), EMAN::EMData::get_data(), EMAN::EMData::get_xsize(), EMAN::EMData::get_ysize(), EMAN::EMData::get_zsize(), LOGWARN, EMAN::Processor::params, square, and EMAN::EMData::update().

03842 {
03843         if (!image) {
03844                 LOGWARN("NULL Image");
03845                 return;
03846         }
03847 
03848         float distance_sigma = params["distance_sigma"];
03849         float value_sigma = params["value_sigma"];
03850         int max_iter = params["niter"];
03851         int half_width = params["half_width"];
03852 
03853         if (half_width < distance_sigma) {
03854                 LOGWARN("localwidth(=%d) should be larger than distance_sigma=(%f)\n",
03855                                                         half_width, distance_sigma);
03856         }
03857 
03858         distance_sigma *= distance_sigma;
03859 
03860         float image_sigma = image->get_attr("sigma");
03861         if (image_sigma > value_sigma) {
03862                 LOGWARN("image sigma(=%f) should be smaller than value_sigma=(%f)\n",
03863                                                         image_sigma, value_sigma);
03864         }
03865         value_sigma *= value_sigma;
03866 
03867         int nx = image->get_xsize();
03868         int ny = image->get_ysize();
03869         int nz = image->get_zsize();
03870 
03871         if(nz==1) {     //for 2D image
03872                 int width=nx, height=ny;
03873 
03874                 int i,j,m,n;
03875 
03876                 float tempfloat1,tempfloat2,tempfloat3;
03877                 int   index1,index2,index;
03878                 int   Iter;
03879                 int   tempint1,tempint3;
03880 
03881                 tempint1=width;
03882                 tempint3=width+2*half_width;
03883 
03884                 float* mask=(float*)calloc((2*half_width+1)*(2*half_width+1),sizeof(float));
03885                 float* OrgImg=(float*)calloc((2*half_width+width)*(2*half_width+height),sizeof(float));
03886                 float* NewImg=image->get_data();
03887 
03888                 for(m=-(half_width);m<=half_width;m++)
03889                         for(n=-(half_width);n<=half_width;n++) {
03890                    index=(m+half_width)*(2*half_width+1)+(n+half_width);
03891                    mask[index]=exp((float)(-(m*m+n*n)/distance_sigma/2.0));
03892                 }
03893 
03894                 //printf("entering bilateral filtering process \n");
03895 
03896                 Iter=0;
03897                 while(Iter<max_iter) {
03898                         for(i=0;i<height;i++)
03899                         for(j=0;j<width;j++) {
03900                                 index1=(i+half_width)*tempint3+(j+half_width);
03901                                         index2=i*tempint1+j;
03902                                 OrgImg[index1]=NewImg[index2];
03903                         }
03904 
03905                         // Mirror Padding
03906                         for(i=0;i<height;i++){
03907                                 for(j=0;j<half_width;j++) OrgImg[(i+half_width)*tempint3+(j)]=OrgImg[(i+half_width)*tempint3+(2*half_width-j)];
03908                                 for(j=0;j<half_width;j++) OrgImg[(i+half_width)*tempint3+(j+width+half_width)]=OrgImg[(i+half_width)*tempint3+(width+half_width-j-2)];
03909                         }
03910                         for(i=0;i<half_width;i++){
03911                                 for(j=0;j<(width+2*half_width);j++) OrgImg[i*tempint3+j]=OrgImg[(2*half_width-i)*tempint3+j];
03912                                 for(j=0;j<(width+2*half_width);j++) OrgImg[(i+height+half_width)*tempint3+j]=OrgImg[(height+half_width-2-i)*tempint3+j];
03913                         }
03914 
03915                         //printf("finish mirror padding process \n");
03916                         //now mirror padding have been done
03917 
03918                         for(i=0;i<height;i++){
03919                                 //printf("now processing the %d th row \n",i);
03920                                 for(j=0;j<width;j++){
03921                                         tempfloat1=0.0; tempfloat2=0.0;
03922                                         for(m=-(half_width);m<=half_width;m++)
03923                                                 for(n=-(half_width);n<=half_width;n++){
03924                                                         index =(m+half_width)*(2*half_width+1)+(n+half_width);
03925                                                         index1=(i+half_width)*tempint3+(j+half_width);
03926                                                         index2=(i+half_width+m)*tempint3+(j+half_width+n);
03927                                                         tempfloat3=(OrgImg[index1]-OrgImg[index2])*(OrgImg[index1]-OrgImg[index2]);
03928 
03929                                                         tempfloat3=mask[index]*(1.0f/(1+tempfloat3/value_sigma));       // Lorentz kernel
03930                                                         //tempfloat3=mask[index]*exp(tempfloat3/Sigma2/(-2.0)); // Guassian kernel
03931                                                         tempfloat1+=tempfloat3;
03932 
03933                                                         tempfloat2+=tempfloat3*OrgImg[(i+half_width+m)*tempint3+(j+half_width+n)];
03934                                         }
03935                                         NewImg[i*width+j]=tempfloat2/tempfloat1;
03936                                 }
03937                         }
03938                         Iter++;
03939             }
03940 
03941             //printf("have finished %d  th iteration\n ",Iter);
03942 //              doneData();
03943                 free(mask);
03944                 free(OrgImg);
03945                 // end of BilaFilter routine
03946 
03947         }
03948         else {  //3D case
03949                 int width = nx;
03950                 int height = ny;
03951                 int slicenum = nz;
03952 
03953                 int slice_size = width * height;
03954                 int new_width = width + 2 * half_width;
03955                 int new_slice_size = (width + 2 * half_width) * (height + 2 * half_width);
03956 
03957                 int width1 = 2 * half_width + 1;
03958                 int mask_size = width1 * width1;
03959                 int old_img_size = (2 * half_width + width) * (2 * half_width + height);
03960 
03961                 int zstart = -half_width;
03962                 int zend = -half_width;
03963                 int is_3d = 0;
03964                 if (nz > 1) {
03965                         mask_size *= width1;
03966                         old_img_size *= (2 * half_width + slicenum);
03967                         zend = half_width;
03968                         is_3d = 1;
03969                 }
03970 
03971                 float *mask = (float *) calloc(mask_size, sizeof(float));
03972                 float *old_img = (float *) calloc(old_img_size, sizeof(float));
03973 
03974                 float *new_img = image->get_data();
03975 
03976                 for (int p = zstart; p <= zend; p++) {
03977                         int cur_p = (p + half_width) * (2 * half_width + 1) * (2 * half_width + 1);
03978 
03979                         for (int m = -half_width; m <= half_width; m++) {
03980                                 int cur_m = (m + half_width) * (2 * half_width + 1) + half_width;
03981 
03982                                 for (int n = -half_width; n <= half_width; n++) {
03983                                         int l = cur_p + cur_m + n;
03984                                         mask[l] = exp((float) (-(m * m + n * n + p * p * is_3d) / distance_sigma / 2.0f));
03985                                 }
03986                         }
03987                 }
03988 
03989                 int iter = 0;
03990                 while (iter < max_iter) {
03991                         for (int k = 0; k < slicenum; k++) {
03992                                 int cur_k1 = (k + half_width) * new_slice_size * is_3d;
03993                                 int cur_k2 = k * slice_size;
03994 
03995                                 for (int i = 0; i < height; i++) {
03996                                         int cur_i1 = (i + half_width) * new_width;
03997                                         int cur_i2 = i * width;
03998 
03999                                         for (int j = 0; j < width; j++) {
04000                                                 int k1 = cur_k1 + cur_i1 + (j + half_width);
04001                                                 int k2 = cur_k2 + cur_i2 + j;
04002                                                 old_img[k1] = new_img[k2];
04003                                         }
04004                                 }
04005                         }
04006 
04007                         for (int k = 0; k < slicenum; k++) {
04008                                 int cur_k = (k + half_width) * new_slice_size * is_3d;
04009 
04010                                 for (int i = 0; i < height; i++) {
04011                                         int cur_i = (i + half_width) * new_width;
04012 
04013                                         for (int j = 0; j < half_width; j++) {
04014                                                 int k1 = cur_k + cur_i + j;
04015                                                 int k2 = cur_k + cur_i + (2 * half_width - j);
04016                                                 old_img[k1] = old_img[k2];
04017                                         }
04018 
04019                                         for (int j = 0; j < half_width; j++) {
04020                                                 int k1 = cur_k + cur_i + (width + half_width + j);
04021                                                 int k2 = cur_k + cur_i + (width + half_width - j - 2);
04022                                                 old_img[k1] = old_img[k2];
04023                                         }
04024                                 }
04025 
04026 
04027                                 for (int i = 0; i < half_width; i++) {
04028                                         int i2 = i * new_width;
04029                                         int i3 = (2 * half_width - i) * new_width;
04030                                         for (int j = 0; j < (width + 2 * half_width); j++) {
04031                                                 int k1 = cur_k + i2 + j;
04032                                                 int k2 = cur_k + i3 + j;
04033                                                 old_img[k1] = old_img[k2];
04034                                         }
04035 
04036                                         i2 = (height + half_width + i) * new_width;
04037                                         i3 = (height + half_width - 2 - i) * new_width;
04038                                         for (int j = 0; j < (width + 2 * half_width); j++) {
04039                                                 int k1 = cur_k + i2 + j;
04040                                                 int k2 = cur_k + i3 + j;
04041                                                 old_img[k1] = old_img[k2];
04042                                         }
04043                                 }
04044                         }
04045 
04046                         size_t idx;
04047                         for (int k = 0; k < slicenum; k++) {
04048                                 int cur_k = (k + half_width) * new_slice_size;
04049 
04050                                 for (int i = 0; i < height; i++) {
04051                                         int cur_i = (i + half_width) * new_width;
04052 
04053                                         for (int j = 0; j < width; j++) {
04054                                                 float f1 = 0;
04055                                                 float f2 = 0;
04056                                                 int k1 = cur_k + cur_i + (j + half_width);
04057 
04058                                                 for (int p = zstart; p <= zend; p++) {
04059                                                         int cur_p1 = (p + half_width) * (2 * half_width + 1) * (2 * half_width + 1);
04060                                                         int cur_p2 = (k + half_width + p) * new_slice_size;
04061 
04062                                                         for (int m = -half_width; m <= half_width; m++) {
04063                                                                 int cur_m1 = (m + half_width) * (2 * half_width + 1);
04064                                                                 int cur_m2 = cur_p2 + cur_i + m * new_width + j + half_width;
04065 
04066                                                                 for (int n = -half_width; n <= half_width; n++) {
04067                                                                         int k = cur_p1 + cur_m1 + (n + half_width);
04068                                                                         int k2 = cur_m2 + n;
04069                                                                         float f3 = Util::square(old_img[k1] - old_img[k2]);
04070 
04071                                                                         f3 = mask[k] * (1.0f / (1 + f3 / value_sigma));
04072                                                                         f1 += f3;
04073                                                                         int l1 = cur_m2 + n;
04074                                                                         f2 += f3 * old_img[l1];
04075                                                                 }
04076 
04077                                                                 idx = k * height * width + i * width + j;
04078                                                                 new_img[idx] = f2 / f1;
04079                                                         }
04080                                                 }
04081                                         }
04082                                 }
04083                         }
04084                         iter++;
04085                 }
04086                 if( mask ) {
04087                         free(mask);
04088                         mask = 0;
04089                 }
04090 
04091                 if( old_img ) {
04092                         free(old_img);
04093                         old_img = 0;
04094                 }
04095         }
04096 
04097         image->update();
04098 }


Member Data Documentation

const string BilateralProcessor::NAME = "bilateral" [static]
 

Definition at line 4022 of file processor.h.

Referenced by get_name().


The documentation for this class was generated from the following files:
Generated on Tue May 25 17:37:24 2010 for EMAN2 by  doxygen 1.4.4