EMAN2
Public Member Functions | Static Public Member Functions | Static Public Attributes
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_sigmameans how large the voxel has impact on its neighbors in spatial domain. The larger it is, the more blurry the resulting image.
value_sigmaeans how large the voxel has impact on its in range domain. The larger it is, the more blurry the resulting image.
niterhow many times to apply this processing on your data.
half_widthprocessing window size = (2 * half_widthh + 1) ^ 3.

Definition at line 4268 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 4277 of file processor.h.

                {
                        return "Bilateral processing on 2D or 3D volume data. Bilateral processing does non-linear weighted averaging processing within a certain window. ";
                }
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 4272 of file processor.h.

References NAME.

                {
                        return NAME;
                }
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 4287 of file processor.h.

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

                {
                        TypeDict d;
                        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.");
                        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.");
                        d.put("niter", EMObject::INT, "how many times to apply this processing on your data.");
                        d.put("half_width", EMObject::INT, "processing window size = (2 * half_widthh + 1) ^ 3.");
                        return d;
                }
static Processor* EMAN::BilateralProcessor::NEW ( ) [inline, static]

Definition at line 4282 of file processor.h.

                {
                        return new BilateralProcessor();
                }
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:
imageThe image to be processed.

Implements EMAN::Processor.

Definition at line 4202 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, EMAN::Processor::params, square, and EMAN::EMData::update().

{
        if (!image) {
                LOGWARN("NULL Image");
                return;
        }

        float distance_sigma = params["distance_sigma"];
        float value_sigma = params["value_sigma"];
        int max_iter = params["niter"];
        int half_width = params["half_width"];

        if (half_width < distance_sigma) {
                LOGWARN("localwidth(=%d) should be larger than distance_sigma=(%f)\n",
                                                        half_width, distance_sigma);
        }

        distance_sigma *= distance_sigma;

        float image_sigma = image->get_attr("sigma");
        if (image_sigma > value_sigma) {
                LOGWARN("image sigma(=%f) should be smaller than value_sigma=(%f)\n",
                                                        image_sigma, value_sigma);
        }
        value_sigma *= value_sigma;

        int nx = image->get_xsize();
        int ny = image->get_ysize();
        int nz = image->get_zsize();

        if(nz==1) {     //for 2D image
                int width=nx, height=ny;

                int i,j,m,n;

                float tempfloat1,tempfloat2,tempfloat3;
                int   index1,index2,index;
                int   Iter;
                int   tempint1,tempint3;

                tempint1=width;
                tempint3=width+2*half_width;

                float* mask=(float*)calloc((2*half_width+1)*(2*half_width+1),sizeof(float));
                float* OrgImg=(float*)calloc((2*half_width+width)*(2*half_width+height),sizeof(float));
                float* NewImg=image->get_data();

                for(m=-(half_width);m<=half_width;m++)
                        for(n=-(half_width);n<=half_width;n++) {
                   index=(m+half_width)*(2*half_width+1)+(n+half_width);
                   mask[index]=exp((float)(-(m*m+n*n)/distance_sigma/2.0));
                }

                //printf("entering bilateral filtering process \n");

                Iter=0;
                while(Iter<max_iter) {
                        for(i=0;i<height;i++)
                        for(j=0;j<width;j++) {
                                index1=(i+half_width)*tempint3+(j+half_width);
                                        index2=i*tempint1+j;
                                OrgImg[index1]=NewImg[index2];
                        }

                        // Mirror Padding
                        for(i=0;i<height;i++){
                                for(j=0;j<half_width;j++) OrgImg[(i+half_width)*tempint3+(j)]=OrgImg[(i+half_width)*tempint3+(2*half_width-j)];
                                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)];
                        }
                        for(i=0;i<half_width;i++){
                                for(j=0;j<(width+2*half_width);j++) OrgImg[i*tempint3+j]=OrgImg[(2*half_width-i)*tempint3+j];
                                for(j=0;j<(width+2*half_width);j++) OrgImg[(i+height+half_width)*tempint3+j]=OrgImg[(height+half_width-2-i)*tempint3+j];
                        }

                        //printf("finish mirror padding process \n");
                        //now mirror padding have been done

                        for(i=0;i<height;i++){
                                //printf("now processing the %d th row \n",i);
                                for(j=0;j<width;j++){
                                        tempfloat1=0.0; tempfloat2=0.0;
                                        for(m=-(half_width);m<=half_width;m++)
                                                for(n=-(half_width);n<=half_width;n++){
                                                        index =(m+half_width)*(2*half_width+1)+(n+half_width);
                                                        index1=(i+half_width)*tempint3+(j+half_width);
                                                        index2=(i+half_width+m)*tempint3+(j+half_width+n);
                                                        tempfloat3=(OrgImg[index1]-OrgImg[index2])*(OrgImg[index1]-OrgImg[index2]);

                                                        tempfloat3=mask[index]*(1.0f/(1+tempfloat3/value_sigma));       // Lorentz kernel
                                                        //tempfloat3=mask[index]*exp(tempfloat3/Sigma2/(-2.0)); // Guassian kernel
                                                        tempfloat1+=tempfloat3;

                                                        tempfloat2+=tempfloat3*OrgImg[(i+half_width+m)*tempint3+(j+half_width+n)];
                                        }
                                        NewImg[i*width+j]=tempfloat2/tempfloat1;
                                }
                        }
                        Iter++;
            }

            //printf("have finished %d  th iteration\n ",Iter);
//              doneData();
                free(mask);
                free(OrgImg);
                // end of BilaFilter routine

        }
        else {  //3D case
                int width = nx;
                int height = ny;
                int slicenum = nz;

                int slice_size = width * height;
                int new_width = width + 2 * half_width;
                int new_slice_size = (width + 2 * half_width) * (height + 2 * half_width);

                int width1 = 2 * half_width + 1;
                int mask_size = width1 * width1;
                int old_img_size = (2 * half_width + width) * (2 * half_width + height);

                int zstart = -half_width;
                int zend = -half_width;
                int is_3d = 0;
                if (nz > 1) {
                        mask_size *= width1;
                        old_img_size *= (2 * half_width + slicenum);
                        zend = half_width;
                        is_3d = 1;
                }

                float *mask = (float *) calloc(mask_size, sizeof(float));
                float *old_img = (float *) calloc(old_img_size, sizeof(float));

                float *new_img = image->get_data();

                for (int p = zstart; p <= zend; p++) {
                        int cur_p = (p + half_width) * (2 * half_width + 1) * (2 * half_width + 1);

                        for (int m = -half_width; m <= half_width; m++) {
                                int cur_m = (m + half_width) * (2 * half_width + 1) + half_width;

                                for (int n = -half_width; n <= half_width; n++) {
                                        int l = cur_p + cur_m + n;
                                        mask[l] = exp((float) (-(m * m + n * n + p * p * is_3d) / distance_sigma / 2.0f));
                                }
                        }
                }

                int iter = 0;
                while (iter < max_iter) {
                        for (int k = 0; k < slicenum; k++) {
                                size_t cur_k1 = (size_t)(k + half_width) * new_slice_size * is_3d;
                                int cur_k2 = k * slice_size;

                                for (int i = 0; i < height; i++) {
                                        int cur_i1 = (i + half_width) * new_width;
                                        int cur_i2 = i * width;

                                        for (int j = 0; j < width; j++) {
                                                size_t k1 = cur_k1 + cur_i1 + (j + half_width);
                                                int k2 = cur_k2 + cur_i2 + j;
                                                old_img[k1] = new_img[k2];
                                        }
                                }
                        }

                        for (int k = 0; k < slicenum; k++) {
                                size_t cur_k = (k + half_width) * new_slice_size * is_3d;

                                for (int i = 0; i < height; i++) {
                                        int cur_i = (i + half_width) * new_width;

                                        for (int j = 0; j < half_width; j++) {
                                                size_t k1 = cur_k + cur_i + j;
                                                size_t k2 = cur_k + cur_i + (2 * half_width - j);
                                                old_img[k1] = old_img[k2];
                                        }

                                        for (int j = 0; j < half_width; j++) {
                                                size_t k1 = cur_k + cur_i + (width + half_width + j);
                                                size_t k2 = cur_k + cur_i + (width + half_width - j - 2);
                                                old_img[k1] = old_img[k2];
                                        }
                                }


                                for (int i = 0; i < half_width; i++) {
                                        int i2 = i * new_width;
                                        int i3 = (2 * half_width - i) * new_width;
                                        for (int j = 0; j < (width + 2 * half_width); j++) {
                                                size_t k1 = cur_k + i2 + j;
                                                size_t k2 = cur_k + i3 + j;
                                                old_img[k1] = old_img[k2];
                                        }

                                        i2 = (height + half_width + i) * new_width;
                                        i3 = (height + half_width - 2 - i) * new_width;
                                        for (int j = 0; j < (width + 2 * half_width); j++) {
                                                size_t k1 = cur_k + i2 + j;
                                                size_t k2 = cur_k + i3 + j;
                                                old_img[k1] = old_img[k2];
                                        }
                                }
                        }

                        size_t idx;
                        for (int k = 0; k < slicenum; k++) {
                                size_t cur_k = (k + half_width) * new_slice_size;

                                for (int i = 0; i < height; i++) {
                                        int cur_i = (i + half_width) * new_width;

                                        for (int j = 0; j < width; j++) {
                                                float f1 = 0;
                                                float f2 = 0;
                                                size_t k1 = cur_k + cur_i + (j + half_width);

                                                for (int p = zstart; p <= zend; p++) {
                                                        size_t cur_p1 = (p + half_width) * (2 * half_width + 1) * (2 * half_width + 1);
                                                        size_t cur_p2 = (k + half_width + p) * new_slice_size;

                                                        for (int m = -half_width; m <= half_width; m++) {
                                                                size_t cur_m1 = (m + half_width) * (2 * half_width + 1);
                                                                size_t cur_m2 = cur_p2 + cur_i + m * new_width + j + half_width;

                                                                for (int n = -half_width; n <= half_width; n++) {
                                                                        size_t k = cur_p1 + cur_m1 + (n + half_width);
                                                                        size_t k2 = cur_m2 + n;
                                                                        float f3 = Util::square(old_img[k1] - old_img[k2]);

                                                                        f3 = mask[k] * (1.0f / (1 + f3 / value_sigma));
                                                                        f1 += f3;
                                                                        size_t l1 = cur_m2 + n;
                                                                        f2 += f3 * old_img[l1];
                                                                }

                                                                idx = (size_t)k * height * width + i * width + j;
                                                                new_img[idx] = f2 / f1;
                                                        }
                                                }
                                        }
                                }
                        }
                        iter++;
                }
                if( mask ) {
                        free(mask);
                        mask = 0;
                }

                if( old_img ) {
                        free(old_img);
                        old_img = 0;
                }
        }

        image->update();
}

Member Data Documentation

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

Definition at line 4297 of file processor.h.

Referenced by get_name().


The documentation for this class was generated from the following files: