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

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

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

04001                 {
04002                         return NAME;
04003                 }

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

References EMAN::TypeDict::put().

04016                 {
04017                         TypeDict d;
04018                         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.");
04019                         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.");
04020                         d.put("niter", EMObject::INT, "how many times to apply this processing on your data.");
04021                         d.put("half_width", EMObject::INT, "processing window size = (2 * half_widthh + 1) ^ 3.");
04022                         return d;
04023                 }

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

Definition at line 4010 of file processor.h.

04011                 {
04012                         return new BilateralProcessor();
04013                 }

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 3895 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().

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


Member Data Documentation

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

Definition at line 147 of file processor.cpp.


The documentation for this class was generated from the following files:
Generated on Tue Jul 12 13:52:15 2011 for EMAN2 by  doxygen 1.3.9.1