关于统计变换（CT/MCT/RMCT）算法的学习和实现

原文地址http://blog.sina.com.cn/s/blog_684c8d630100turx.html

刚开会每周的例会，最讨厌开会了，不过为了能顺利毕业，只能忍了。闲话不多说了，下面把上周学习的一个简单的算法总结一下，以备后面写毕业论文的时候可以参考一下。
一、Census Transform(CT)算法的学习
    Census Transform 算法是Ramin Zabih 和 John Woodfill 于1994年在他们的论文《Non-parametric LocalTransforms for Computing VisualCorrespondence》中提出的，正如他们在论文中所说，这是一种非参数变换，主要用来表征图像的局部结构特征，能够比较好的检测到图像中的边缘特征和角点特征，从这篇论文的470次的引用次数来看，CT算法用处还是挺广泛的。下面简要介绍一下CT算法的基本思想：用一个3*3或者5*5的滑动窗口遍历整幅图像，对于每次遍历的位置，以3*3为例，假设某个位置如下图所示
123 127 129
126 128 129
127 131 130
然后比较此窗口中每个像素值（中心除外）与中心像素值的大小，如果比中心像素值小，则比较结果为1，否则为0。由此，得到如下结果：

1 1 0
1   0
1 0 0

然后，把此窗口结果组成一个序列：11010100，以此二进制序列表示的值来代替原图像窗口中心点的像素。如此下去，等到窗口滑动完整幅图像，我们就得到原图像做统计变换（CT）之后的图像。

注意：只能作用于灰度图像，对于彩色图像，则需要转换为灰度图之后再操作；为了简单起见，我没有考虑图像的边缘像素值。
   下面给出一种用 matlab 实现的版本：
[plain] view plain copy
 
 在CODE上查看代码片派生到我的代码片
  1. % *************************************************************************  
  2. % Title: Function-Census Transform of a given Image  
  3. % Author: Siddhant Ahuja  
  4. % Created: May 2008  
  5. % Copyright Siddhant Ahuja, 2008  
  6. % Inputs: Image (var: inputImage), Window size assuming square window (var:  
  7. % windowSize) of 3x3 or 5x5 only.  
  8. % Outputs: Census Tranformed Image (var: censusTransformedImage),  
  9. % Time taken (var: timeTaken)  
  10. % Example Usage of Function: [a,b]=funcCensusOneImage('Img.png', 3)  
  11. % *************************************************************************  
  12. function [censusTransformedImage, timeTaken] = funcCensusOneImage(inputImage, windowSize)  
  13. % Grab the image information (metadata) using the function imfinfo  
  14. try  
  15. imageInfo = imfinfo(inputImage);  
  16. % Since Census Transform is applied on a grayscale image, determine if the  
  17. % input image is already in grayscale or color  
  18. if(getfield(imageInfo,'ColorType')=='truecolor')  
  19. % Read an image using imread function, convert from RGB color space to  
  20. % grayscale using rgb2gray function and assign it to variable inputImage  
  21. inputImage=rgb2gray(imread(inputImage));  
  22. else if(getfield(imageInfo,'ColorType')=='grayscale')  
  23. % If the image is already in grayscale, then just read it.         
  24. inputImage=imread(inputImage);  
  25. else  
  26. error('The Color Type of Input Image is not acceptable. Acceptable color types are truecolor or grayscale.');  
  27. end  
  28. end  
  29. catch  
  30. inputImage = inputImage;  
  31. end  
  32. % Find the size (columns and rows) of the image and assign the rows to  
  33. % variable nr, and columns to variable nc  
  34. [nr,nc] = size(inputImage);  
  35. % Check the size of window to see if it is an odd number.  
  36. if (mod(windowSize,2)==0)  
  37. error('The window size must be an odd number.');  
  38. end  
  39. if (windowSize==3)  
  40. bits=uint8(0);  
  41. % Create an image of size nr and nc, fill it with zeros and assign  
  42. % it to variable censusTransformedImage of type uint8  
  43. censusTransformedImage=uint8(zeros(nr,nc));  
  44. else if (windowSize==5)  
  45. bits=uint32(0);  
  46. % Create an image of size nr and nc, fill it with zeros and assign  
  47. % it to variable censusTransformedImage of type uint32         
  48. censusTransformedImage=uint32(zeros(nr,nc));  
  49. else  
  50. error('The size of the window is not acceptable. Just 3x3 and 5x5 windows are acceptable.');  
  51. end  
  52. end  
  53. % Initialize the timer to calculate the time consumed.  
  54. tic;  
  55. % Find out how many rows and columns are to the left/right/up/down of the  
  56. % central pixel  
  57. C= (windowSize-1)/2;  
  58. for j=C+1:1:nc-C % Go through all the columns in an image (minus C at the borders)  
  59. for i=C+1:1:nr-C % Go through all the rows in an image (minus C at the borders)   
  60. census = 0; % Initialize default census to 0  
  61. for a=-C:1:C % Within the square window, go through all the rows  
  62. for b=-C:1:C % Within the square window, go through all the columns  
  63. if (~(a==0 && b==0)) % Exclude the centre pixel from the calculation,原来是(C+1)，现改为0  
  64. census=bitshift(census,1); %Shift the bits to the left  by 1  
  65. % If the intensity of the neighboring pixel is less than  
  66. % that of the central pixel, then add one to the bit  
  67. % string  
  68. if (inputImage(i+a,j+b) < inputImage(i,j))  
  69. census=census+1;  
  70. end  
  71. end  
  72. end  
  73. end  
  74. % Assign the census bit string value to the pixel in imgTemp  
  75. censusTransformedImage(i,j) = census;  
  76. end  
  77. end  
  78. % Stop the timer to calculate the time consumed.  
  79. timeTaken=toc;  
   这是我在网上找到的一个实现的版本，注释比较多，除去注释的话，真正代码没有50行，比较简单，相信大家都可以看的懂。
  之前讲了CT算法的实现，下面说一下 MCT 以及 RMCT的实现。

二、Modified Census Transform (MCT)算法
    MCT 算法是 CT 算法的一个修改版本，它是由Bernhard Froba 在做人脸检测的时候提出来的，在他2004年发表的论文《Face Detection with theModified Census Transform》中，Bernhard Froba将 CT算法中“滑动窗口中每个像素值与中心位置像素做比较”改为“滑动窗口中每个像素值与整个窗口中像素的均值做比较”，这样，原有的每个3*3的窗口可能产生256种序列（因为没有算中心像素），现在变为可能产生512种序列（其实全0和全1的序列表示的是同样的信息，可以排除一个），也就是做完MCT之后，图像的每个像素值的范围为0——511，这样就能够比较充分的利用3*3的核（至于为什么这么说，可以看看前面提到的那篇论文）。如此来计算的话，则前面例子产生的结果窗口应该为：
1 1 0
1 0 0
1 0 0
对应的二进制序列为：110100100。然后以此作为中心像素点的像素值，循环完毕之后便得到MCT之后的图像。需要注意的一点是，如果要使变换之后的图像得到显示，应该对像素值做一下归一化，使其在0——255之间。
 
三、Revised Modified Census Transform (RMCT)算法
    RMCT 算法其实又是对 MCT的又一次修改，它与 MCT的不同之处仅仅在于一个微小的△m，即：在滑动窗口像素均值上加上一个微小的变量△m=1或者2。其他都是完全一样的。
 
    下面附上这两种修改版统计变换的 C++代码，代码是我自己编的，是基于VS2008和OpenCV2.0的，仅供参考：
[cpp] view plain copy
 
 在CODE上查看代码片派生到我的代码片
  1. #include "stdafx.h"   
  2. #include "MCT.h"   
  3. #include "highgui.h"   
  4.    
  5. MCT::MCT()   
  6. {   
  7.     window_size = 0;   
  8. }   
  9.    
  10. MCT::~MCT()   
  11. {   
  12. }   
  13.    
  14. void MCT::ModifiedCensusTransform(IplImage *input_image, IplImage *mct_image, const int window_size, const int delta )   
  15. {   
  16.     CvSize image_size = cvGetSize(input_image);   
  17.     int image_width = image_size.width;   
  18.     int image_height = image_size.height;   
  19.    
  20.     IplImage *gray_image = cvCreateImage(cvGetSize(input_image), input_image->depth, 1);   
  21.     cvSetZero(gray_image);   
  22.    
  23.     if(input_image->nChannels != 1)   
  24.     {   
  25.         cvCvtColor(input_image, gray_image, CV_RGB2GRAY);   
  26.     }   
  27.     else   
  28.     {   
  29.         cvCopy(input_image, gray_image);   
  30.     }   
  31.    
  32.    
  33.     IplImage *modified_image = NULL;   
  34.     switch (window_size)   
  35.     {   
  36.     case 3:   
  37.         modified_image = cvCreateImage(image_size, IPL_DEPTH_16U, 1);   
  38.         cvSetZero(modified_image);   
  39.         break;   
  40.     case 5:   
  41.         modified_image = cvCreateImage(image_size, IPL_DEPTH_32S, 1);   
  42.         cvSetZero(modified_image);   
  43.         break;   
  44.     default:   
  45.         printf("window size must be 3 or 5! ");   
  46.         exit(EXIT_FAILURE);   
  47.     }   
  48.    
  49.     CvMat window;   
  50.     for(int i = 0; i < image_height - window_size; i++)   
  51.     {   
  52.         for(int j = 0; j < image_width - window_size; j++)   
  53.         {   
  54.             unsigned long census = 0;   
  55.             CvRect roi = cvRect(j, i, window_size, window_size);   
  56.             cvGetSubRect(gray_image, &window, roi);   
  57.             CvScalar m = cvAvg(&window, NULL);   
  58.    
  59.             for(int w = 0; w < window_size; w++)   
  60.             {   
  61.                 for(int h = 0; h < window_size; h++)   
  62.                 {   
  63.                     census = census << 1; //左移1位   
  64.                     double tempvalue = cvGetReal2D(&window, w, h);   
  65.                     if(tempvalue < m.val[0] + delta)   
  66.                         census += 1;   
  67.                 }   
  68.             }   
  69.             cvSetReal2D(modified_image, i, j, census);   
  70.         }   
  71.     }   
  72.     //cvConvertScaleAbs(modified_image, mct_image, 1, 0);   
  73.     Normalize(modified_image, mct_image);   
  74.     cvReleaseImage(&gray_image);   
  75.     cvReleaseImage(&modified_image);   
  76. }   
  77.    
  78. void MCT::Normalize(IplImage *mct_image, IplImage *nor_image)   
  79. {   
  80.     double minv, maxv;   
  81.     cvMinMaxLoc(mct_image, &minv, &maxv);   
  82.     for (int i = 0; i < mct_image->height; i++)   
  83.     {   
  84.         for (int j = 0; j < mct_image->width; j++)   
  85.         {   
  86.             double tempv = cvGetReal2D(mct_image, i, j);   
  87.             tempv = (tempv - minv) / (maxv - minv) * 255;   
  88.             cvSetReal2D(nor_image, i, j, tempv);   
  89.         }   
  90.     }   
  91. }   
由于算法比较简单，所以没有写注释，应该比较容易理解。
四、CT/MCT/RMCT的应用
   目前我读到的几篇论文里面，他们主要用于人脸检测或者面部伪装检测，用于做图像的预处理，这可能是因为CT算法对光照不敏感，可以比较好的排除光照对图像的影响。这里给出用到该算法的几篇paper：
1、《Face Detection with the Modified CensusTransform》（前面提到的那篇）
2、《Adaboost Based Disguised Face Discrimination on EmbeddedDevices》
3、《Disguised-Face Discriminator for Embedded Systems》
 
OVER！
【推广】 免费学中医，健康全家人
原文地址：https://www.cnblogs.com/liuchaogege/p/5160828.html