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 = "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 4072 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 4081 of file processor.h.

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

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

References NAME.

04077                 {
04078                         return NAME;
04079                 }

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

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

04092                 {
04093                         TypeDict d;
04094                         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.");
04095                         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.");
04096                         d.put("niter", EMObject::INT, "how many times to apply this processing on your data.");
04097                         d.put("half_width", EMObject::INT, "processing window size = (2 * half_widthh + 1) ^ 3.");
04098                         return d;
04099                 }

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

Definition at line 4086 of file processor.h.

04087                 {
04088                         return new BilateralProcessor();
04089                 }

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

03983 {
03984         if (!image) {
03985                 LOGWARN("NULL Image");
03986                 return;
03987         }
03988 
03989         float distance_sigma = params["distance_sigma"];
03990         float value_sigma = params["value_sigma"];
03991         int max_iter = params["niter"];
03992         int half_width = params["half_width"];
03993 
03994         if (half_width < distance_sigma) {
03995                 LOGWARN("localwidth(=%d) should be larger than distance_sigma=(%f)\n",
03996                                                         half_width, distance_sigma);
03997         }
03998 
03999         distance_sigma *= distance_sigma;
04000 
04001         float image_sigma = image->get_attr("sigma");
04002         if (image_sigma > value_sigma) {
04003                 LOGWARN("image sigma(=%f) should be smaller than value_sigma=(%f)\n",
04004                                                         image_sigma, value_sigma);
04005         }
04006         value_sigma *= value_sigma;
04007 
04008         int nx = image->get_xsize();
04009         int ny = image->get_ysize();
04010         int nz = image->get_zsize();
04011 
04012         if(nz==1) {     //for 2D image
04013                 int width=nx, height=ny;
04014 
04015                 int i,j,m,n;
04016 
04017                 float tempfloat1,tempfloat2,tempfloat3;
04018                 int   index1,index2,index;
04019                 int   Iter;
04020                 int   tempint1,tempint3;
04021 
04022                 tempint1=width;
04023                 tempint3=width+2*half_width;
04024 
04025                 float* mask=(float*)calloc((2*half_width+1)*(2*half_width+1),sizeof(float));
04026                 float* OrgImg=(float*)calloc((2*half_width+width)*(2*half_width+height),sizeof(float));
04027                 float* NewImg=image->get_data();
04028 
04029                 for(m=-(half_width);m<=half_width;m++)
04030                         for(n=-(half_width);n<=half_width;n++) {
04031                    index=(m+half_width)*(2*half_width+1)+(n+half_width);
04032                    mask[index]=exp((float)(-(m*m+n*n)/distance_sigma/2.0));
04033                 }
04034 
04035                 //printf("entering bilateral filtering process \n");
04036 
04037                 Iter=0;
04038                 while(Iter<max_iter) {
04039                         for(i=0;i<height;i++)
04040                         for(j=0;j<width;j++) {
04041                                 index1=(i+half_width)*tempint3+(j+half_width);
04042                                         index2=i*tempint1+j;
04043                                 OrgImg[index1]=NewImg[index2];
04044                         }
04045 
04046                         // Mirror Padding
04047                         for(i=0;i<height;i++){
04048                                 for(j=0;j<half_width;j++) OrgImg[(i+half_width)*tempint3+(j)]=OrgImg[(i+half_width)*tempint3+(2*half_width-j)];
04049                                 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)];
04050                         }
04051                         for(i=0;i<half_width;i++){
04052                                 for(j=0;j<(width+2*half_width);j++) OrgImg[i*tempint3+j]=OrgImg[(2*half_width-i)*tempint3+j];
04053                                 for(j=0;j<(width+2*half_width);j++) OrgImg[(i+height+half_width)*tempint3+j]=OrgImg[(height+half_width-2-i)*tempint3+j];
04054                         }
04055 
04056                         //printf("finish mirror padding process \n");
04057                         //now mirror padding have been done
04058 
04059                         for(i=0;i<height;i++){
04060                                 //printf("now processing the %d th row \n",i);
04061                                 for(j=0;j<width;j++){
04062                                         tempfloat1=0.0; tempfloat2=0.0;
04063                                         for(m=-(half_width);m<=half_width;m++)
04064                                                 for(n=-(half_width);n<=half_width;n++){
04065                                                         index =(m+half_width)*(2*half_width+1)+(n+half_width);
04066                                                         index1=(i+half_width)*tempint3+(j+half_width);
04067                                                         index2=(i+half_width+m)*tempint3+(j+half_width+n);
04068                                                         tempfloat3=(OrgImg[index1]-OrgImg[index2])*(OrgImg[index1]-OrgImg[index2]);
04069 
04070                                                         tempfloat3=mask[index]*(1.0f/(1+tempfloat3/value_sigma));       // Lorentz kernel
04071                                                         //tempfloat3=mask[index]*exp(tempfloat3/Sigma2/(-2.0)); // Guassian kernel
04072                                                         tempfloat1+=tempfloat3;
04073 
04074                                                         tempfloat2+=tempfloat3*OrgImg[(i+half_width+m)*tempint3+(j+half_width+n)];
04075                                         }
04076                                         NewImg[i*width+j]=tempfloat2/tempfloat1;
04077                                 }
04078                         }
04079                         Iter++;
04080             }
04081 
04082             //printf("have finished %d  th iteration\n ",Iter);
04083 //              doneData();
04084                 free(mask);
04085                 free(OrgImg);
04086                 // end of BilaFilter routine
04087 
04088         }
04089         else {  //3D case
04090                 int width = nx;
04091                 int height = ny;
04092                 int slicenum = nz;
04093 
04094                 int slice_size = width * height;
04095                 int new_width = width + 2 * half_width;
04096                 int new_slice_size = (width + 2 * half_width) * (height + 2 * half_width);
04097 
04098                 int width1 = 2 * half_width + 1;
04099                 int mask_size = width1 * width1;
04100                 int old_img_size = (2 * half_width + width) * (2 * half_width + height);
04101 
04102                 int zstart = -half_width;
04103                 int zend = -half_width;
04104                 int is_3d = 0;
04105                 if (nz > 1) {
04106                         mask_size *= width1;
04107                         old_img_size *= (2 * half_width + slicenum);
04108                         zend = half_width;
04109                         is_3d = 1;
04110                 }
04111 
04112                 float *mask = (float *) calloc(mask_size, sizeof(float));
04113                 float *old_img = (float *) calloc(old_img_size, sizeof(float));
04114 
04115                 float *new_img = image->get_data();
04116 
04117                 for (int p = zstart; p <= zend; p++) {
04118                         int cur_p = (p + half_width) * (2 * half_width + 1) * (2 * half_width + 1);
04119 
04120                         for (int m = -half_width; m <= half_width; m++) {
04121                                 int cur_m = (m + half_width) * (2 * half_width + 1) + half_width;
04122 
04123                                 for (int n = -half_width; n <= half_width; n++) {
04124                                         int l = cur_p + cur_m + n;
04125                                         mask[l] = exp((float) (-(m * m + n * n + p * p * is_3d) / distance_sigma / 2.0f));
04126                                 }
04127                         }
04128                 }
04129 
04130                 int iter = 0;
04131                 while (iter < max_iter) {
04132                         for (int k = 0; k < slicenum; k++) {
04133                                 size_t cur_k1 = (size_t)(k + half_width) * new_slice_size * is_3d;
04134                                 int cur_k2 = k * slice_size;
04135 
04136                                 for (int i = 0; i < height; i++) {
04137                                         int cur_i1 = (i + half_width) * new_width;
04138                                         int cur_i2 = i * width;
04139 
04140                                         for (int j = 0; j < width; j++) {
04141                                                 size_t k1 = cur_k1 + cur_i1 + (j + half_width);
04142                                                 int k2 = cur_k2 + cur_i2 + j;
04143                                                 old_img[k1] = new_img[k2];
04144                                         }
04145                                 }
04146                         }
04147 
04148                         for (int k = 0; k < slicenum; k++) {
04149                                 size_t cur_k = (k + half_width) * new_slice_size * is_3d;
04150 
04151                                 for (int i = 0; i < height; i++) {
04152                                         int cur_i = (i + half_width) * new_width;
04153 
04154                                         for (int j = 0; j < half_width; j++) {
04155                                                 size_t k1 = cur_k + cur_i + j;
04156                                                 size_t k2 = cur_k + cur_i + (2 * half_width - j);
04157                                                 old_img[k1] = old_img[k2];
04158                                         }
04159 
04160                                         for (int j = 0; j < half_width; j++) {
04161                                                 size_t k1 = cur_k + cur_i + (width + half_width + j);
04162                                                 size_t k2 = cur_k + cur_i + (width + half_width - j - 2);
04163                                                 old_img[k1] = old_img[k2];
04164                                         }
04165                                 }
04166 
04167 
04168                                 for (int i = 0; i < half_width; i++) {
04169                                         int i2 = i * new_width;
04170                                         int i3 = (2 * half_width - i) * new_width;
04171                                         for (int j = 0; j < (width + 2 * half_width); j++) {
04172                                                 size_t k1 = cur_k + i2 + j;
04173                                                 size_t k2 = cur_k + i3 + j;
04174                                                 old_img[k1] = old_img[k2];
04175                                         }
04176 
04177                                         i2 = (height + half_width + i) * new_width;
04178                                         i3 = (height + half_width - 2 - i) * new_width;
04179                                         for (int j = 0; j < (width + 2 * half_width); j++) {
04180                                                 size_t k1 = cur_k + i2 + j;
04181                                                 size_t k2 = cur_k + i3 + j;
04182                                                 old_img[k1] = old_img[k2];
04183                                         }
04184                                 }
04185                         }
04186 
04187                         size_t idx;
04188                         for (int k = 0; k < slicenum; k++) {
04189                                 size_t cur_k = (k + half_width) * new_slice_size;
04190 
04191                                 for (int i = 0; i < height; i++) {
04192                                         int cur_i = (i + half_width) * new_width;
04193 
04194                                         for (int j = 0; j < width; j++) {
04195                                                 float f1 = 0;
04196                                                 float f2 = 0;
04197                                                 size_t k1 = cur_k + cur_i + (j + half_width);
04198 
04199                                                 for (int p = zstart; p <= zend; p++) {
04200                                                         size_t cur_p1 = (p + half_width) * (2 * half_width + 1) * (2 * half_width + 1);
04201                                                         size_t cur_p2 = (k + half_width + p) * new_slice_size;
04202 
04203                                                         for (int m = -half_width; m <= half_width; m++) {
04204                                                                 size_t cur_m1 = (m + half_width) * (2 * half_width + 1);
04205                                                                 size_t cur_m2 = cur_p2 + cur_i + m * new_width + j + half_width;
04206 
04207                                                                 for (int n = -half_width; n <= half_width; n++) {
04208                                                                         size_t k = cur_p1 + cur_m1 + (n + half_width);
04209                                                                         size_t k2 = cur_m2 + n;
04210                                                                         float f3 = Util::square(old_img[k1] - old_img[k2]);
04211 
04212                                                                         f3 = mask[k] * (1.0f / (1 + f3 / value_sigma));
04213                                                                         f1 += f3;
04214                                                                         size_t l1 = cur_m2 + n;
04215                                                                         f2 += f3 * old_img[l1];
04216                                                                 }
04217 
04218                                                                 idx = (size_t)k * height * width + i * width + j;
04219                                                                 new_img[idx] = f2 / f1;
04220                                                         }
04221                                                 }
04222                                         }
04223                                 }
04224                         }
04225                         iter++;
04226                 }
04227                 if( mask ) {
04228                         free(mask);
04229                         mask = 0;
04230                 }
04231 
04232                 if( old_img ) {
04233                         free(old_img);
04234                         old_img = 0;
04235                 }
04236         }
04237 
04238         image->update();
04239 }


Member Data Documentation

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

Definition at line 4101 of file processor.h.

Referenced by get_name().


The documentation for this class was generated from the following files:
Generated on Thu May 3 10:10:21 2012 for EMAN2 by  doxygen 1.4.7