implementation of HoG

hog.hpp

/*

    // prepare
    HoG hog(PATCH_SIZE, cv::Size(8, 8), cv::Size(2, 2), cv::Size(8, 8), 9);
    // ...or...
    MultiScaleHoG hog;
    hog.add(PATCH_SIZE, cv::Size(4, 4), cv::Size(2, 2), cv::Size(4, 4), 9);
    hog.add(PATCH_SIZE, cv::Size(16, 16), cv::Size(2, 2), cv::Size(8, 8), 9);
    hog.add(PATCH_SIZE, cv::Size(32, 32), cv::Size(2, 2), cv::Size(16, 16), 9);

    std::cout << cv::format("HoG feature dimensions = %d", hog.feature_dim()) << std::endl;

    // container of feature vectors & integral histograms
    std::vector<float> ihist;
    cv::Mat_<float> X(N, hog.feature_dim());

    // extract hog
    for(int i=0; i < N; i++){
        const cv::Mat& src = images[i];
        hog.prepare(src, ihist);
        hog.extract(src, ihist, cv::Point(0, 0), X.ptr<float>(i));
    }

*/

#pragma once

#include <useopencv.h>
#include <blas.hpp>
#include <cmath>
#include <vector>
#include <boost/foreach.hpp>

class HoG {
public:
    cv::Size window_size, cell_size/*pixels*/, block_size/*cells*/, stride;
    int bins;
    cv::Size block_size_pixel() const {
        return cv::Size(cell_size.width * block_size.width, cell_size.height * block_size.height);
    };
    cv::Size block_num() const {
        return
            cv::Size(
            (window_size.width  - block_size_pixel().width)  / stride.width + 1,
            (window_size.height - block_size_pixel().height) / stride.height + 1
            );
    };
    int block_dim() const {
        return block_size.width * block_size.height * bins;
    };
    int feature_dim() const {
        assert(((window_size.width  - block_size_pixel().width)  % stride.width)  == 0);
        assert(((window_size.height - block_size_pixel().height) % stride.height) == 0);
        return block_dim() * block_num().width * block_num().height;
    };
    int hist_size() const {
        return bins * window_size.width * window_size.height;
    };
public:
    HoG()
        : window_size(64, 64), cell_size(8, 8), block_size(2, 2), stride(4, 4), bins(9)
    {};
    ~HoG(){};

    HoG(const cv::Size window_size, const cv::Size cell_size, const cv::Size block_size, const cv::Size stride, const int bins)
        : window_size(window_size), cell_size(cell_size), block_size(block_size), stride(stride), bins(bins)
    {
        feature_dim(); // for assertion
    };

    static cv::Mat mo_filter(const cv::Mat& src)
    {
        assert(src.channels() == 1);
        cv::Mat mo = cv::Mat::zeros(src.size(), CV_8UC2);

        int dx, dy;
        uchar m;
        for( int y = 1; y < src.rows-1; y++ )
            for( int x = 1; x < src.cols-1; x++ ){
                int maxch = 0;
                uchar maxm = 0;
                dx = (int)( src.at<uchar>(y, x-1) - src.at<uchar>(y, x+1) );
                dy = (int)( src.at<uchar>(y-1, x) - src.at<uchar>(y+1, x) );
                m  = (uchar)std::sqrtf((float)(dx * dx + dy * dy));
                cv::Vec2b& mo_pix = mo.at<cv::Vec2b>(y, x);
                mo_pix[0] = m;
                mo_pix[1] = (uchar)(dy == 0 ? -90 : std::atanf(((float)dx/(float)dy)) * 180 / 3.14159) + 90;
            }
        return mo;
    };

    template <class T>
    void prepare(const cv::Mat& src, std::vector<T>& ihist) const {
        std::vector<T> col_sum(this->bins);
        if(ihist.size() < this->bins * src.rows * src.cols){
            ihist.resize(this->bins * src.rows * src.cols, 0);
        }
        cv::Mat mo = mo_filter(src);

        // image size (width * height * bin)
        for( int y = 0; y < src.rows; y++ ){
            for(int i=0; i < col_sum.size(); i++) col_sum[i] = 0;
            for( int x = 0; x < src.cols; x++ ){
                /* compute bin */
                cv::Vec2b& mo_pix = mo.at<cv::Vec2b>(y, x);
                int bin = (int)(mo_pix[1] / (180 / this->bins));
                T a = (T)( mo_pix[1] % (180 / this->bins) ) / (T)(180 / this->bins);

                // HOG
                assert(bin < bins);
                col_sum[bin] += mo_pix[0] * (1-a);
                col_sum[(bin == this->bins - 1) ? 0 : bin] += mo_pix[0] * a;
                //col_sum[bin] += PIX(moc, x, y, 0);
            
                blas::copy(this->bins, &col_sum[0], 1, &ihist[ (y * src.cols + x) * this->bins ], 1);
                if( y > 0 )
                    blas::axpy(this->bins, 1, &ihist[ ((y-1) * src.cols + x) * this->bins ], 1,
                                              &ihist[ (y     * src.cols + x) * this->bins ], 1);
            }
        }
    };

    template <class T>
    void extract(const cv::Mat& src, const std::vector<T>& ihist, const cv::Point offset, T* F) const
    {
        T* tF = F;
        auto block_num = this->block_num();
        for( int by = 0; by < block_num.height; by++ ){
            for( int bx = 0; bx < block_num.width; bx++ ){
                int x = bx * stride.width  + offset.x;
                int y = by * stride.height + offset.y;
                int xp, yp, xm, ym;

                for( int cy = 0; cy < block_size.height; cy++ ){
                    for( int cx = 0; cx < block_size.width; cx++ ){
                        xm = x + cx * cell_size.width;
                        ym = y + cy * cell_size.height;
                        xp = xm + cell_size.width  - 1;
                        yp = ym + cell_size.height - 1;
                        blas::copy(bins,     &ihist[(xp + yp * src.cols) * bins], 1, tF, 1);
                        blas::axpy(bins, -1, &ihist[(xm + yp * src.cols) * bins], 1, tF, 1);
                        blas::axpy(bins, -1, &ihist[(xp + ym * src.cols) * bins], 1, tF, 1);
                        blas::axpy(bins, +1, &ihist[(xm + ym * src.cols) * bins], 1, tF, 1);
                        tF += bins;
                    }
                }

                int block_dim = bins * block_size.width * block_size.height;
                T nrm = blas::nrm2(block_dim, F, 1);
                blas::scal(block_dim, (T)( 1 / sqrt(nrm * nrm + 1) ), F, 1);
                F += block_dim;
                //printf("%g, ", nrm);
            }
        }
    }
};

class MultiScaleHoG {
    std::vector<HoG> hogs;
public:
    MultiScaleHoG()
        : hogs()
    {};
    ~MultiScaleHoG(){};

    void add(const cv::Size window_size, const cv::Size cell_size, const cv::Size block_size, const cv::Size stride, const int bins)
    {
        hogs.push_back( HoG(window_size, cell_size, block_size, stride, bins) );
    };

    int feature_dim(){
        int d = 0;
        BOOST_FOREACH(auto& hog, hogs){
            d += hog.feature_dim();
        }
        return d;
    };

    template <class T>
    void prepare(const cv::Mat& src, std::vector<T>& ihist) const {
        assert(hogs.size() > 0);
        hogs[0].prepare(src, ihist);
    };

    template <class T>
    void extract(const cv::Mat& src, const std::vector<T>& ihist, const cv::Point offset, T* F) const
    {
        BOOST_FOREACH(const auto& hog, hogs){
            hog.extract(src, ihist, offset, F);
            F += hog.feature_dim();
        }
    };
};