|| || || || || ||

---

LbpEval.h

Go to the documentation of this file.

#ifndef Impala_Core_VideoSet_LbpEval_h
#define Impala_Core_VideoSet_LbpEval_h

#include "Core/Array/RGB2Gray.h"
#include "Core/Feature/Computor.h"
#include "Core/VideoSet/Reporter.h"
#include "Core/Geometry/PointR.h"
#include "Core/Geometry/PointZ.h"
#include "Core/Histogram/MakeCooccurenceMatrix.h"
#include "Core/Histogram/Histogram2dTem.h"
#include "Core/Feature/MarkovStationaryFeature.h"
#include "Persistency/FeatureTableRepository.h"

namespace Impala
{
namespace Core
{
namespace VideoSet
{

namespace BitOperation
{
    int onecount(unsigned int c, int bits)
    {
        int count = 0, base = 1;
        for (int i=0;i<bits;i++)
        {
            if (c & base)
                count++;
            base <<= 1;
        }
        return count;
    }

    int transitions(unsigned int c, int bits)
    {
        int base = 1;
        int current = c & base;
        int current2;
        int changes = 0;
        for (int i=1;i<bits;i++)
        {
            base <<= 1;
            current2 = (c & base) >> i;
            if (current ^ current2)
                changes++;
            current = current2;
        }
        return changes; //(changes <= 2)? 1 : 0;
    }
}

void MakeMap(int samples, int*& map)
{
    ILOG_VAR(Impala.Core.VidSet.LbpEval.MakeMap);
    ILOG_DEBUG("creating map");
    int size = 1<<samples;
    map = new int[size];
    for (int i=0;i<(1<<samples);i++)
    {
        if (BitOperation::transitions(i,samples) <= 2)
            map[i] = BitOperation::onecount(i,samples);
        else
            map[i] = samples+1;
    }
}

void MakePoints(int samples, int radius, Geometry::PointZ*& points, double*& multipliers)
{
    points = new Geometry::PointZ[samples];
    multipliers = new double[samples*4];
    Geometry::PointR* offsets = new Geometry::PointR[samples];

    double step = 2 * M_PI / samples;
    for (unsigned int i=0;i<samples;i++)
    {
        double tmpX = radius * cos(i * step);
        double tmpY = radius * sin(i * step);
        points[i].mX = (int)tmpX;
        points[i].mY = (int)tmpY;
        offsets[i].mX = tmpX - points[i].mX;
        offsets[i].mY = tmpY - points[i].mY;
        if (offsets[i].mX < 1.0e-10 && offsets[i].mX > -1.0e-10) //rounding error
            offsets[i].mX = 0;
        if (offsets[i].mY < 1.0e-10 && offsets[i].mY > -1.0e-10) //rounding error
            offsets[i].mY = 0;
        
        if (tmpX < 0 && offsets[i].mX != 0)
        {
            points[i].mX -= 1;
            offsets[i].mX += 1;
        }
        if (tmpY < 0 && offsets[i].mY != 0)
        {
            points[i].mY -= 1;
            offsets[i].mY += 1;
        }

        double dx = 1-offsets[i].mX;
        double dy = 1-offsets[i].mY;
        multipliers[i*4+0] = dx*dy;
        multipliers[i*4+1] = offsets[i].mX*dy;
        multipliers[i*4+2] = dx*offsets[i].mY;
        multipliers[i*4+3] = offsets[i].mX*offsets[i].mY;
    }
    delete offsets;
}

template <class T>
inline double
InterpolateAtPtr(T* ptr, int i, int width, double* multipliers)
{
    if(multipliers[(i*4)] != 1.0)
    {
        return ptr[0] * multipliers[(i*4)] + 
               ptr[1] * multipliers[(i*4)+1] + 
               ptr[width] * multipliers[(i*4)+2] + 
               ptr[width+1] * multipliers[(i*4)+3] + 1e-10;
    }
    return ptr[0];
}


void
Lbp(Array::Array2dScalarUInt8*& dst,
    Array::Array2dScalarReal64* src,
    int* patternMap,
    int samples,
    int radius)
{
    ILOG_VAR(Impala.Core.VidSet.LbpEval);
    int border2 = radius * 2;
    int height = src->CH();
    int width = src->CW();
    if(dst == 0)
        dst = new Array::Array2dScalarUInt8(width - border2, height - border2, 0, 0);
    if(patternMap == 0)
    {
        ILOG_ERROR("LBP features without map not supported");
        return;
    }
    
    // by multiplying with width+1 we go down and sideways at the same time
    double* center = src->CPB() + radius * (1 + width);
    double** ptr = new double*[samples];
    Geometry::PointZ* points;
    double* multipliers;
    MakePoints(samples, radius, points, multipliers);

    for(int i=0 ; i<samples ; i++)
        ptr[i] = center + points[i].mX + points[i].mY * width;

    for(int y=0 ; y<height-border2 ; y++)
    {
        unsigned char* dstPtr = dst->CPB(0,y);
        for(int x=0 ; x<width-border2 ; x++)
        {
            unsigned int value = 0;
            unsigned int base = 1;
            for (int i=0;i<samples;i++)
            {
                if (InterpolateAtPtr(ptr[i],i,width,multipliers) >= *center)
                    value |= base;
                base <<= 1;
                ptr[i]++;
            }
            center++;
            *dstPtr = patternMap[value];
            ++dstPtr;
        }
        for(int i=0 ; i<samples ; i++)
            ptr[i] += border2;
        center += border2;
    }

    delete points;
    delete multipliers;
    delete ptr;
}

class LbpEval : public Listener
{
public:
    typedef Persistency::FeatureLocator FeatureLocator;

    LbpEval(Reporter* reporter, CmdOptions& options)
    {
        ILOG_DEBUG("LbpEval c'tor");
        int featureLengthFactor = 1;
        mDoMSF = options.GetBool("msf-feature");
        if(mDoMSF)
        {
            featureLengthFactor = 2;
            ILOG_INFO("computing msf feature");
        }
        else
        {
            ILOG_INFO("making histogram feature");
        }
        mTableSet = new Feature::FeatureTableSet();
        int tableSize = 1000;
        int regions = 6;
        mMaps = new int*[3];
        mTableSet->Add(new Feature::FeatureTable(GetDefinition(8, 1.0), tableSize,
                                        (8+2)*featureLengthFactor*regions));
        MakeMap(8, mMaps[0]);
        mTableSet->Add(new Feature::FeatureTable(GetDefinition(16, 2.0), tableSize,
                                        (16+2)*featureLengthFactor*regions));
        MakeMap(16, mMaps[1]);
        mTableSet->Add(new Feature::FeatureTable(GetDefinition(24, 3.0), tableSize,
                                        (24+2)*featureLengthFactor*regions));
        MakeMap(24, mMaps[2]);
    }

    virtual ~LbpEval()
    {
        delete mMaps[0];
        delete mMaps[1];
        delete mMaps[2];
        delete mTableSet;
    }

    virtual void
    HandleNewWalkPartial(VideoSet* vs, int startFrame, int stepSize,
                         int numberFrames)
    {
        mLoc.SetNumberFrames(numberFrames);
        mLoc.SetStartFrame(startFrame);
    }

    virtual void
    HandleNewWalk(VideoSet* vs, String walkType)
    {
        mLoc = FeatureLocator(vs->GetLocator(), false, false, walkType,
                              "unknown", "unknown");
    }

    virtual void
    HandleNewFrame(VideoSet* vs, int fileId, Stream::RgbDataSrc* src)
    {
        Array::Array2dVec3UInt8* im = Array::ArrayCreate<Array::Array2dVec3UInt8>
            (src->FrameWidth(), src->FrameHeight(), 0, 0, src->DataPtr(), true);
        Array::Array2dScalarReal64* greyim = 0;
        Array::RGB2Gray(greyim, im);
        delete im;
        Quid quid = vs->GetQuidFrame(fileId, src->FrameNr());
        for(int i=0 ; i<mTableSet->Size() ; ++i)
        {
            Feature::FeatureDefinition fd = mTableSet->GetFeatureDefinition(i);
            //ILOG_DEBUG("processing " << fd.AsString());
            int samples = atoi(fd.GetParameterValue("samples"));
            int radius = atoi(fd.GetParameterValue("radius"));
            Array::Array2dScalarUInt8* lbpIm = 0;
            Lbp(lbpIm, greyim, mMaps[i], samples, radius);

            //here we could loop over the image parts (patatsnijden)
            Vector::VectorTem<Real64> vector;
            Vector::VectorTem<Real64>& accu = vector;
            std::vector<Geometry::Rectangle> rects;
            ILOG_DEBUG("greyim w:"<<lbpIm->CW()<<" h:"<<lbpIm->CH());
            GetRectangles(rects, lbpIm->CW(), lbpIm->CH(), radius);
            for(int r=0 ; r<6 ; ++r)
            {
                ILOG_DEBUG("rect "<<r<<" -> "<<rects[r]);
            }
            for(int r=0 ; r<6 ; ++r)
            {
                //either simple histogram, or cooccurence matrix
                if(mDoMSF)
                {
                    Histogram::Histogram2dTem<double> hist(-1.5, samples+1.5, samples+2, -1.5, samples+1.5, samples+2);
                    Histogram::MakeCooccurenceMatrix(&hist, lbpIm, rects[i]);
                    const Vector::VectorTem<Real64>& v =
                        Feature::MarkovStationaryFeature(&hist);
                    accu = Append(accu, v);
                }
                else
                {
                    double dummy=0;
                    Histogram::Histogram1dTem<double> hist(-0.5, samples+1.5, samples+2, &dummy);
                    MakeHistogram1d(&hist, lbpIm, rects[i]);
                    hist.Normalize();
                    if(dummy != 0)
                        ILOG_ERROR("lbp with map should not create outliers in histogram");
                    accu = Append(accu, hist);
                }
            }
            ILOG_DEBUG("feature length: " << accu.Size());
            Feature::FeatureTable* t = mTableSet->GetTable(i);
            t->Add(quid, accu);
            delete lbpIm;
        }
        delete greyim;
    }

    virtual void
    HandleDoneFile(VideoSet* vs, int fileId, Stream::RgbDataSrc* src)
    {
        for (int i=0 ; i<mTableSet->Size() ; i++)
        {
            if (mLoc.GetIsPartial() && (mTableSet->GetTable(i)->Size() < 2))
                continue;
            Feature::FeatureDefinition def = mTableSet->GetFeatureDefinition(i);
            mLoc.SetFeatureDef(def);
            mLoc.SetContainer(vs->GetContainerFile(fileId));
            Persistency::FeatureTableRepository().Add(mLoc,
                                                      mTableSet->GetTable(i));
            mTableSet->GetTable(i)->SetSize(0);
        }
    }

    virtual void
    HandleDoneWalk(VideoSet* vs)
    {
//        std::cout << "Feature definitions: " << std::endl;
//        for (int i=0 ; i<mTableSet->Size() ; i++)
//            std::cout << mTableSet->GetDefinitionAsString(i) << std::endl;
    }

private:
    void
    GetRectangles(std::vector<Geometry::Rectangle>& rects, int width, int height, int border)
    {
        rects.resize(6);
        int halfWidth = width/2;
        int halfHeight = height/2;
        rects[0] = Geometry::Rectangle(border, border, width-border-1, height-border-1);
        rects[1] = Geometry::Rectangle(border, border, halfWidth-1, halfHeight-1);
        rects[2] = Geometry::Rectangle(halfWidth, border, width-border-1, halfHeight-1);
        rects[3] = Geometry::Rectangle(border, halfHeight, halfWidth-1, height-border-1);
        rects[4] = Geometry::Rectangle(halfWidth, halfHeight, width-border-1, height-border-1);
        int width4 = width/4;
        int height4 = height/4;
        rects[5] = Geometry::Rectangle(width4, height4, height4+halfWidth-1, height4+halfHeight-1);
    }

    Feature::FeatureDefinition
    GetDefinition(int samples, double radius)
    {
        Feature::FeatureDefinition def("lbp");
        def.AddParameter("samples", samples);
        def.AddParameter("radius", radius);
        if(mDoMSF)
            def.AddParameter("vector", "msf");
        else
            def.AddParameter("vector", "hist");
        return def;
    }

    int mGridHeight;
    int mGridWidth;
    bool               mDoMSF;
    int**              mMaps;
    Feature::FeatureTableSet* mTableSet; // a table for each LBP variation
    FeatureLocator     mLoc;

    ILOG_VAR_DEC;
};

ILOG_VAR_INIT(LbpEval, Impala.Core.VideoSet);


} // namespace VideoSet
} // namespace Core
} // namespace Impala

#endif

---