Angular Radial Transform

將 MPEG7Expr 文件內的 code 接上 OpenCV 就成了..

struct ART
{
    typedef double       BaseType;
    typedef BaseType**   F2D;
    typedef BaseType**** F4D;    

    int    nFeature;                                //特徵向量個數
    uchar* feature;                                 //特徵向量
    F2D    artR, artI;                              //ART係數
    F4D    basisR, basisI;                          //底基空間
    int    baseW,  baseH;                           //底基圖寬高, 
    int    nAngle, nRadius;                         //角度切割數, 半徑切割數
    IplImage* img;                                  //跟底基圖等大的影像
    
    ART (int w=101, int h=101, int nAngle=12, int nRadius=6)
        : baseW(w), baseH(h), nAngle(nAngle), nRadius(nRadius)
    {
        int      x, y, m, n;                                 
        int      cx = w/2, cy = h/2;
        int      sz = sizeof (BaseType);
        BaseType r, a, tmp1, tmp2;
        
        nFeature = nAngle * nRadius - 1;            //配置底基空間
        basisR   = (F4D) new_arr0 (sz, 4, nAngle, nRadius, w,h);    
        basisI   = (F4D) new_arr0 (sz, 4, nAngle, nRadius, w,h);
        artR     = (F2D) new_arr  (sz, 2, nAngle, nRadius);
        artI     = (F2D) new_arr  (sz, 2, nAngle, nRadius);
        feature  = (uchar*) new_arr  (1, 1, nFeature);
        img      = cvCreateImage (cvSize (baseW,baseH), 8, 1);

        for (y=0; y<baseH; ++y)                     //生成各種底基圖像    
            for (x=0; x<baseW; ++x) {
                r = sqrt (double((x-cx)*(x-cx) + (y-cy)*(y-cy)));
                a = atan2 (double(y-cy), double(x-cx));
                for (m=0; m<nAngle; ++m)
                    for (n=0; n<nRadius; ++n) {
                        tmp1 = cos (r * PI * n/cx);
                        tmp2 = cos (r * PI * n/cy);
                        basisR [m][n][x][y] = tmp1 * cos (a*m);
                        basisI [m][n][x][y] = tmp2 * sin (a*m);
                    }
            }
    }
   ~ART()                      
    {
        del_arr (basisR);                           //釋放底基空間
        del_arr (basisI);
        del_arr (artR);
        del_arr (artI);
        del_arr (feature);
        cvReleaseImage (&img);                      //釋放暫用影像
    }

    uchar quantize (BaseType v)
    {
        if (v < 0.003585473) return 0;
        if (v < 0.007418411) return 1;
        if (v < 0.011535520) return 2;
        if (v < 0.015982337) return 3;
        if (v < 0.020816302) return 4;
        if (v < 0.026111312) return 5;
        if (v < 0.031964674) return 6;
        if (v < 0.038508176) return 7;
        if (v < 0.045926586) return 8;
        if (v < 0.054490513) return 9;
        if (v < 0.064619488) return 10;
        if (v < 0.077016351) return 11;
        if (v < 0.092998687) return 12;
        if (v < 0.115524524) return 13;
        if (v < 0.154032694) return 14;
        return 15;
    }

    //對圖像 path 進行 ART 變換, 產生的特徵向量存於 feature[] 中

    void transform (char* path)             
    {
        IplImage* src = cvLoadImage (path, CV_LOAD_IMAGE_GRAYSCALE);
        if (!src) {    
            printf ("Can not load %s", path);
            _getch(); exit(-1);
        }        
        cvResize (src, img, CV_INTER_AREA);         //縮放源圖像
        cvReleaseImage (&src);                      //釋放原圖
        
        BaseType r, i, mean, color;
        int x, y, m, n, k;                                 
        int W = img->widthStep;                     //真實寬度 
        uchar *data = (uchar*) img->imageData;
        
        for (m=0; m<nAngle; ++m)
            for (n=0; n<nRadius; ++n) 
                artR[m][n] = artI[m][n] = 0;        //歸零

        for (y=0; y<baseH; ++y, data+=W)
            for (x=0; x<baseW; ++x) {
                color = BaseType (data[x])/ 255;
                for (m=0; m<nAngle; ++m)
                    for (n=0; n<nRadius; ++n) {
                        artR[m][n] += color * basisR[m][n][x][y];
                        artI[m][n] -= color * basisI[m][n][x][y];
                    }
            }
        mean = sqrt (artR[0][0]*artR[0][0] + artI[0][0]*artI[0][0]);

        for (k=m=0; m<nAngle; ++m)                  //將係數正規化與量化
            for (n=0; n<nRadius; ++n) 
                if (m || n) {
                    r = artR[m][n];
                    i = artI[m][n];
                    feature[k++] = quantize (sqrt (r*r + i*i)/ mean);
                }        
    }

    //取得底基圖, id= 0,1= 實數,虛數, 圖像會暫存於 img 中
    
    void get_base (int id, int m, int n, bool bShow=false)   
    {
        F4D base = id? basisI: basisR;
        int step = img->nChannels;
        int W    = img->widthStep;                  //真實寬度 
        int h    = img->height; 
        int w    = img->width;
        int w3   = w * step;
        int i, x, y;
        uchar *data = (uchar*) img->imageData;

        for (y=0; y<h; y++, data+=W)
            for (x=0; x<w3; x+=step) {
                i = int ((base[m][n][x][y]+1)*128);
                data[x] = i>255?255:i<0?0:i;
            }          
        if (bShow) {
            cvNamedWindow ("Base");
            cvShowImage ("Base", img);
            cvWaitKey (0);
            cvDestroyWindow ("Base");
        }
    }
};