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

ProcessorNEW ()

Static Public Attributes

const string NAME = "filter.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 4027 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 4036 of file processor.h.

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

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 4031 of file processor.h.

04032                 {
04033                         return NAME;
04034                 }

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 4046 of file processor.h.

References EMAN::TypeDict::put().

04047                 {
04048                         TypeDict d;
04049                         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.");
04050                         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.");
04051                         d.put("niter", EMObject::INT, "how many times to apply this processing on your data.");
04052                         d.put("half_width", EMObject::INT, "processing window size = (2 * half_widthh + 1) ^ 3.");
04053                         return d;
04054                 }

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

Definition at line 4041 of file processor.h.

04042                 {
04043                         return new BilateralProcessor();
04044                 }

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 3910 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, nx, ny, square, and EMAN::EMData::update().

03911 {
03912         if (!image) {
03913                 LOGWARN("NULL Image");
03914                 return;
03915         }
03916 
03917         float distance_sigma = params["distance_sigma"];
03918         float value_sigma = params["value_sigma"];
03919         int max_iter = params["niter"];
03920         int half_width = params["half_width"];
03921 
03922         if (half_width < distance_sigma) {
03923                 LOGWARN("localwidth(=%d) should be larger than distance_sigma=(%f)\n",
03924                                                         half_width, distance_sigma);
03925         }
03926 
03927         distance_sigma *= distance_sigma;
03928 
03929         float image_sigma = image->get_attr("sigma");
03930         if (image_sigma > value_sigma) {
03931                 LOGWARN("image sigma(=%f) should be smaller than value_sigma=(%f)\n",
03932                                                         image_sigma, value_sigma);
03933         }
03934         value_sigma *= value_sigma;
03935 
03936         int nx = image->get_xsize();
03937         int ny = image->get_ysize();
03938         int nz = image->get_zsize();
03939 
03940         if(nz==1) {     //for 2D image
03941                 int width=nx, height=ny;
03942 
03943                 int i,j,m,n;
03944 
03945                 float tempfloat1,tempfloat2,tempfloat3;
03946                 int   index1,index2,index;
03947                 int   Iter;
03948                 int   tempint1,tempint3;
03949 
03950                 tempint1=width;
03951                 tempint3=width+2*half_width;
03952 
03953                 float* mask=(float*)calloc((2*half_width+1)*(2*half_width+1),sizeof(float));
03954                 float* OrgImg=(float*)calloc((2*half_width+width)*(2*half_width+height),sizeof(float));
03955                 float* NewImg=image->get_data();
03956 
03957                 for(m=-(half_width);m<=half_width;m++)
03958                         for(n=-(half_width);n<=half_width;n++) {
03959                    index=(m+half_width)*(2*half_width+1)+(n+half_width);
03960                    mask[index]=exp((float)(-(m*m+n*n)/distance_sigma/2.0));
03961                 }
03962 
03963                 //printf("entering bilateral filtering process \n");
03964 
03965                 Iter=0;
03966                 while(Iter<max_iter) {
03967                         for(i=0;i<height;i++)
03968                         for(j=0;j<width;j++) {
03969                                 index1=(i+half_width)*tempint3+(j+half_width);
03970                                         index2=i*tempint1+j;
03971                                 OrgImg[index1]=NewImg[index2];
03972                         }
03973 
03974                         // Mirror Padding
03975                         for(i=0;i<height;i++){
03976                                 for(j=0;j<half_width;j++) OrgImg[(i+half_width)*tempint3+(j)]=OrgImg[(i+half_width)*tempint3+(2*half_width-j)];
03977                                 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)];
03978                         }
03979                         for(i=0;i<half_width;i++){
03980                                 for(j=0;j<(width+2*half_width);j++) OrgImg[i*tempint3+j]=OrgImg[(2*half_width-i)*tempint3+j];
03981                                 for(j=0;j<(width+2*half_width);j++) OrgImg[(i+height+half_width)*tempint3+j]=OrgImg[(height+half_width-2-i)*tempint3+j];
03982                         }
03983 
03984                         //printf("finish mirror padding process \n");
03985                         //now mirror padding have been done
03986 
03987                         for(i=0;i<height;i++){
03988                                 //printf("now processing the %d th row \n",i);
03989                                 for(j=0;j<width;j++){
03990                                         tempfloat1=0.0; tempfloat2=0.0;
03991                                         for(m=-(half_width);m<=half_width;m++)
03992                                                 for(n=-(half_width);n<=half_width;n++){
03993                                                         index =(m+half_width)*(2*half_width+1)+(n+half_width);
03994                                                         index1=(i+half_width)*tempint3+(j+half_width);
03995                                                         index2=(i+half_width+m)*tempint3+(j+half_width+n);
03996                                                         tempfloat3=(OrgImg[index1]-OrgImg[index2])*(OrgImg[index1]-OrgImg[index2]);
03997 
03998                                                         tempfloat3=mask[index]*(1.0f/(1+tempfloat3/value_sigma));       // Lorentz kernel
03999                                                         //tempfloat3=mask[index]*exp(tempfloat3/Sigma2/(-2.0)); // Guassian kernel
04000                                                         tempfloat1+=tempfloat3;
04001 
04002                                                         tempfloat2+=tempfloat3*OrgImg[(i+half_width+m)*tempint3+(j+half_width+n)];
04003                                         }
04004                                         NewImg[i*width+j]=tempfloat2/tempfloat1;
04005                                 }
04006                         }
04007                         Iter++;
04008             }
04009 
04010             //printf("have finished %d  th iteration\n ",Iter);
04011 //              doneData();
04012                 free(mask);
04013                 free(OrgImg);
04014                 // end of BilaFilter routine
04015 
04016         }
04017         else {  //3D case
04018                 int width = nx;
04019                 int height = ny;
04020                 int slicenum = nz;
04021 
04022                 int slice_size = width * height;
04023                 int new_width = width + 2 * half_width;
04024                 int new_slice_size = (width + 2 * half_width) * (height + 2 * half_width);
04025 
04026                 int width1 = 2 * half_width + 1;
04027                 int mask_size = width1 * width1;
04028                 int old_img_size = (2 * half_width + width) * (2 * half_width + height);
04029 
04030                 int zstart = -half_width;
04031                 int zend = -half_width;
04032                 int is_3d = 0;
04033                 if (nz > 1) {
04034                         mask_size *= width1;
04035                         old_img_size *= (2 * half_width + slicenum);
04036                         zend = half_width;
04037                         is_3d = 1;
04038                 }
04039 
04040                 float *mask = (float *) calloc(mask_size, sizeof(float));
04041                 float *old_img = (float *) calloc(old_img_size, sizeof(float));
04042 
04043                 float *new_img = image->get_data();
04044 
04045                 for (int p = zstart; p <= zend; p++) {
04046                         int cur_p = (p + half_width) * (2 * half_width + 1) * (2 * half_width + 1);
04047 
04048                         for (int m = -half_width; m <= half_width; m++) {
04049                                 int cur_m = (m + half_width) * (2 * half_width + 1) + half_width;
04050 
04051                                 for (int n = -half_width; n <= half_width; n++) {
04052                                         int l = cur_p + cur_m + n;
04053                                         mask[l] = exp((float) (-(m * m + n * n + p * p * is_3d) / distance_sigma / 2.0f));
04054                                 }
04055                         }
04056                 }
04057 
04058                 int iter = 0;
04059                 while (iter < max_iter) {
04060                         for (int k = 0; k < slicenum; k++) {
04061                                 size_t cur_k1 = (size_t)(k + half_width) * new_slice_size * is_3d;
04062                                 int cur_k2 = k * slice_size;
04063 
04064                                 for (int i = 0; i < height; i++) {
04065                                         int cur_i1 = (i + half_width) * new_width;
04066                                         int cur_i2 = i * width;
04067 
04068                                         for (int j = 0; j < width; j++) {
04069                                                 size_t k1 = cur_k1 + cur_i1 + (j + half_width);
04070                                                 int k2 = cur_k2 + cur_i2 + j;
04071                                                 old_img[k1] = new_img[k2];
04072                                         }
04073                                 }
04074                         }
04075 
04076                         for (int k = 0; k < slicenum; k++) {
04077                                 size_t cur_k = (k + half_width) * new_slice_size * is_3d;
04078 
04079                                 for (int i = 0; i < height; i++) {
04080                                         int cur_i = (i + half_width) * new_width;
04081 
04082                                         for (int j = 0; j < half_width; j++) {
04083                                                 size_t k1 = cur_k + cur_i + j;
04084                                                 size_t k2 = cur_k + cur_i + (2 * half_width - j);
04085                                                 old_img[k1] = old_img[k2];
04086                                         }
04087 
04088                                         for (int j = 0; j < half_width; j++) {
04089                                                 size_t k1 = cur_k + cur_i + (width + half_width + j);
04090                                                 size_t k2 = cur_k + cur_i + (width + half_width - j - 2);
04091                                                 old_img[k1] = old_img[k2];
04092                                         }
04093                                 }
04094 
04095 
04096                                 for (int i = 0; i < half_width; i++) {
04097                                         int i2 = i * new_width;
04098                                         int i3 = (2 * half_width - i) * new_width;
04099                                         for (int j = 0; j < (width + 2 * half_width); j++) {
04100                                                 size_t k1 = cur_k + i2 + j;
04101                                                 size_t k2 = cur_k + i3 + j;
04102                                                 old_img[k1] = old_img[k2];
04103                                         }
04104 
04105                                         i2 = (height + half_width + i) * new_width;
04106                                         i3 = (height + half_width - 2 - i) * new_width;
04107                                         for (int j = 0; j < (width + 2 * half_width); j++) {
04108                                                 size_t k1 = cur_k + i2 + j;
04109                                                 size_t k2 = cur_k + i3 + j;
04110                                                 old_img[k1] = old_img[k2];
04111                                         }
04112                                 }
04113                         }
04114 
04115                         size_t idx;
04116                         for (int k = 0; k < slicenum; k++) {
04117                                 size_t cur_k = (k + half_width) * new_slice_size;
04118 
04119                                 for (int i = 0; i < height; i++) {
04120                                         int cur_i = (i + half_width) * new_width;
04121 
04122                                         for (int j = 0; j < width; j++) {
04123                                                 float f1 = 0;
04124                                                 float f2 = 0;
04125                                                 size_t k1 = cur_k + cur_i + (j + half_width);
04126 
04127                                                 for (int p = zstart; p <= zend; p++) {
04128                                                         size_t cur_p1 = (p + half_width) * (2 * half_width + 1) * (2 * half_width + 1);
04129                                                         size_t cur_p2 = (k + half_width + p) * new_slice_size;
04130 
04131                                                         for (int m = -half_width; m <= half_width; m++) {
04132                                                                 size_t cur_m1 = (m + half_width) * (2 * half_width + 1);
04133                                                                 size_t cur_m2 = cur_p2 + cur_i + m * new_width + j + half_width;
04134 
04135                                                                 for (int n = -half_width; n <= half_width; n++) {
04136                                                                         size_t k = cur_p1 + cur_m1 + (n + half_width);
04137                                                                         size_t k2 = cur_m2 + n;
04138                                                                         float f3 = Util::square(old_img[k1] - old_img[k2]);
04139 
04140                                                                         f3 = mask[k] * (1.0f / (1 + f3 / value_sigma));
04141                                                                         f1 += f3;
04142                                                                         size_t l1 = cur_m2 + n;
04143                                                                         f2 += f3 * old_img[l1];
04144                                                                 }
04145 
04146                                                                 idx = (size_t)k * height * width + i * width + j;
04147                                                                 new_img[idx] = f2 / f1;
04148                                                         }
04149                                                 }
04150                                         }
04151                                 }
04152                         }
04153                         iter++;
04154                 }
04155                 if( mask ) {
04156                         free(mask);
04157                         mask = 0;
04158                 }
04159 
04160                 if( old_img ) {
04161                         free(old_img);
04162                         old_img = 0;
04163                 }
04164         }
04165 
04166         image->update();
04167 }


Member Data Documentation

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

Definition at line 142 of file processor.cpp.


The documentation for this class was generated from the following files:
Generated on Thu Nov 17 12:46:17 2011 for EMAN2 by  doxygen 1.3.9.1