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:

[legend]
Collaboration diagram for EMAN::BilateralProcessor:
[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 = "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.

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::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                 }

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

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


Member Data Documentation

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

Definition at line 141 of file processor.cpp.


The documentation for this class was generated from the following files:
Generated on Fri Apr 30 15:39:24 2010 for EMAN2 by  doxygen 1.3.9.1