00001 #ifndef Impala_Core_VideoSet_LbpEval_h
00002 #define Impala_Core_VideoSet_LbpEval_h
00003
00004 #include "Core/Array/RGB2Gray.h"
00005 #include "Core/Feature/Computor.h"
00006 #include "Core/VideoSet/Reporter.h"
00007 #include "Core/Geometry/PointR.h"
00008 #include "Core/Geometry/PointZ.h"
00009 #include "Core/Histogram/MakeCooccurenceMatrix.h"
00010 #include "Core/Histogram/Histogram2dTem.h"
00011 #include "Core/Feature/MarkovStationaryFeature.h"
00012 #include "Persistency/FeatureTableRepository.h"
00013
00014 namespace Impala
00015 {
00016 namespace Core
00017 {
00018 namespace VideoSet
00019 {
00020
00021 namespace BitOperation
00022 {
00023 int onecount(unsigned int c, int bits)
00024 {
00025 int count = 0, base = 1;
00026 for (int i=0;i<bits;i++)
00027 {
00028 if (c & base)
00029 count++;
00030 base <<= 1;
00031 }
00032 return count;
00033 }
00034
00035 int transitions(unsigned int c, int bits)
00036 {
00037 int base = 1;
00038 int current = c & base;
00039 int current2;
00040 int changes = 0;
00041 for (int i=1;i<bits;i++)
00042 {
00043 base <<= 1;
00044 current2 = (c & base) >> i;
00045 if (current ^ current2)
00046 changes++;
00047 current = current2;
00048 }
00049 return changes;
00050 }
00051 }
00052
00053 void MakeMap(int samples, int*& map)
00054 {
00055 ILOG_VAR(Impala.Core.VidSet.LbpEval.MakeMap);
00056 ILOG_DEBUG("creating map");
00057 int size = 1<<samples;
00058 map = new int[size];
00059 for (int i=0;i<(1<<samples);i++)
00060 {
00061 if (BitOperation::transitions(i,samples) <= 2)
00062 map[i] = BitOperation::onecount(i,samples);
00063 else
00064 map[i] = samples+1;
00065 }
00066 }
00067
00068 void MakePoints(int samples, int radius, Geometry::PointZ*& points, double*& multipliers)
00069 {
00070 points = new Geometry::PointZ[samples];
00071 multipliers = new double[samples*4];
00072 Geometry::PointR* offsets = new Geometry::PointR[samples];
00073
00074 double step = 2 * M_PI / samples;
00075 for (unsigned int i=0;i<samples;i++)
00076 {
00077 double tmpX = radius * cos(i * step);
00078 double tmpY = radius * sin(i * step);
00079 points[i].mX = (int)tmpX;
00080 points[i].mY = (int)tmpY;
00081 offsets[i].mX = tmpX - points[i].mX;
00082 offsets[i].mY = tmpY - points[i].mY;
00083 if (offsets[i].mX < 1.0e-10 && offsets[i].mX > -1.0e-10)
00084 offsets[i].mX = 0;
00085 if (offsets[i].mY < 1.0e-10 && offsets[i].mY > -1.0e-10)
00086 offsets[i].mY = 0;
00087
00088 if (tmpX < 0 && offsets[i].mX != 0)
00089 {
00090 points[i].mX -= 1;
00091 offsets[i].mX += 1;
00092 }
00093 if (tmpY < 0 && offsets[i].mY != 0)
00094 {
00095 points[i].mY -= 1;
00096 offsets[i].mY += 1;
00097 }
00098
00099 double dx = 1-offsets[i].mX;
00100 double dy = 1-offsets[i].mY;
00101 multipliers[i*4+0] = dx*dy;
00102 multipliers[i*4+1] = offsets[i].mX*dy;
00103 multipliers[i*4+2] = dx*offsets[i].mY;
00104 multipliers[i*4+3] = offsets[i].mX*offsets[i].mY;
00105 }
00106 delete offsets;
00107 }
00108
00109 template <class T>
00110 inline double
00111 InterpolateAtPtr(T* ptr, int i, int width, double* multipliers)
00112 {
00113 if(multipliers[(i*4)] != 1.0)
00114 {
00115 return ptr[0] * multipliers[(i*4)] +
00116 ptr[1] * multipliers[(i*4)+1] +
00117 ptr[width] * multipliers[(i*4)+2] +
00118 ptr[width+1] * multipliers[(i*4)+3] + 1e-10;
00119 }
00120 return ptr[0];
00121 }
00122
00123
00124 void
00125 Lbp(Array::Array2dScalarUInt8*& dst,
00126 Array::Array2dScalarReal64* src,
00127 int* patternMap,
00128 int samples,
00129 int radius)
00130 {
00131 ILOG_VAR(Impala.Core.VidSet.LbpEval);
00132 int border2 = radius * 2;
00133 int height = src->CH();
00134 int width = src->CW();
00135 if(dst == 0)
00136 dst = new Array::Array2dScalarUInt8(width - border2, height - border2, 0, 0);
00137 if(patternMap == 0)
00138 {
00139 ILOG_ERROR("LBP features without map not supported");
00140 return;
00141 }
00142
00143
00144 double* center = src->CPB() + radius * (1 + width);
00145 double** ptr = new double*[samples];
00146 Geometry::PointZ* points;
00147 double* multipliers;
00148 MakePoints(samples, radius, points, multipliers);
00149
00150 for(int i=0 ; i<samples ; i++)
00151 ptr[i] = center + points[i].mX + points[i].mY * width;
00152
00153 for(int y=0 ; y<height-border2 ; y++)
00154 {
00155 unsigned char* dstPtr = dst->CPB(0,y);
00156 for(int x=0 ; x<width-border2 ; x++)
00157 {
00158 unsigned int value = 0;
00159 unsigned int base = 1;
00160 for (int i=0;i<samples;i++)
00161 {
00162 if (InterpolateAtPtr(ptr[i],i,width,multipliers) >= *center)
00163 value |= base;
00164 base <<= 1;
00165 ptr[i]++;
00166 }
00167 center++;
00168 *dstPtr = patternMap[value];
00169 ++dstPtr;
00170 }
00171 for(int i=0 ; i<samples ; i++)
00172 ptr[i] += border2;
00173 center += border2;
00174 }
00175
00176 delete points;
00177 delete multipliers;
00178 delete ptr;
00179 }
00180
00186 class LbpEval : public Listener
00187 {
00188 public:
00189 typedef Persistency::FeatureLocator FeatureLocator;
00190
00191 LbpEval(Reporter* reporter, CmdOptions& options)
00192 {
00193 ILOG_DEBUG("LbpEval c'tor");
00194 int featureLengthFactor = 1;
00195 mDoMSF = options.GetBool("msf-feature");
00196 if(mDoMSF)
00197 {
00198 featureLengthFactor = 2;
00199 ILOG_INFO("computing msf feature");
00200 }
00201 else
00202 {
00203 ILOG_INFO("making histogram feature");
00204 }
00205 mTableSet = new Feature::FeatureTableSet();
00206 int tableSize = 1000;
00207 int regions = 6;
00208 mMaps = new int*[3];
00209 mTableSet->Add(new Feature::FeatureTable(GetDefinition(8, 1.0), tableSize,
00210 (8+2)*featureLengthFactor*regions));
00211 MakeMap(8, mMaps[0]);
00212 mTableSet->Add(new Feature::FeatureTable(GetDefinition(16, 2.0), tableSize,
00213 (16+2)*featureLengthFactor*regions));
00214 MakeMap(16, mMaps[1]);
00215 mTableSet->Add(new Feature::FeatureTable(GetDefinition(24, 3.0), tableSize,
00216 (24+2)*featureLengthFactor*regions));
00217 MakeMap(24, mMaps[2]);
00218 }
00219
00220 virtual ~LbpEval()
00221 {
00222 delete mMaps[0];
00223 delete mMaps[1];
00224 delete mMaps[2];
00225 delete mTableSet;
00226 }
00227
00228 virtual void
00229 HandleNewWalkPartial(VideoSet* vs, int startFrame, int stepSize,
00230 int numberFrames)
00231 {
00232 mLoc.SetNumberFrames(numberFrames);
00233 mLoc.SetStartFrame(startFrame);
00234 }
00235
00236 virtual void
00237 HandleNewWalk(VideoSet* vs, String walkType)
00238 {
00239 mLoc = FeatureLocator(vs->GetLocator(), false, false, walkType,
00240 "unknown", "unknown");
00241 }
00242
00243 virtual void
00244 HandleNewFrame(VideoSet* vs, int fileId, Stream::RgbDataSrc* src)
00245 {
00246 Array::Array2dVec3UInt8* im = Array::ArrayCreate<Array::Array2dVec3UInt8>
00247 (src->FrameWidth(), src->FrameHeight(), 0, 0, src->DataPtr(), true);
00248 Array::Array2dScalarReal64* greyim = 0;
00249 Array::RGB2Gray(greyim, im);
00250 delete im;
00251 Quid quid = vs->GetQuidFrame(fileId, src->FrameNr());
00252 for(int i=0 ; i<mTableSet->Size() ; ++i)
00253 {
00254 Feature::FeatureDefinition fd = mTableSet->GetFeatureDefinition(i);
00255
00256 int samples = atoi(fd.GetParameterValue("samples"));
00257 int radius = atoi(fd.GetParameterValue("radius"));
00258 Array::Array2dScalarUInt8* lbpIm = 0;
00259 Lbp(lbpIm, greyim, mMaps[i], samples, radius);
00260
00261
00262 Vector::VectorTem<Real64> vector;
00263 Vector::VectorTem<Real64>& accu = vector;
00264 std::vector<Geometry::Rectangle> rects;
00265 ILOG_DEBUG("greyim w:"<<lbpIm->CW()<<" h:"<<lbpIm->CH());
00266 GetRectangles(rects, lbpIm->CW(), lbpIm->CH(), radius);
00267 for(int r=0 ; r<6 ; ++r)
00268 {
00269 ILOG_DEBUG("rect "<<r<<" -> "<<rects[r]);
00270 }
00271 for(int r=0 ; r<6 ; ++r)
00272 {
00273
00274 if(mDoMSF)
00275 {
00276 Histogram::Histogram2dTem<double> hist(-1.5, samples+1.5, samples+2, -1.5, samples+1.5, samples+2);
00277 Histogram::MakeCooccurenceMatrix(&hist, lbpIm, rects[i]);
00278 const Vector::VectorTem<Real64>& v =
00279 Feature::MarkovStationaryFeature(&hist);
00280 accu = Append(accu, v);
00281 }
00282 else
00283 {
00284 double dummy=0;
00285 Histogram::Histogram1dTem<double> hist(-0.5, samples+1.5, samples+2, &dummy);
00286 MakeHistogram1d(&hist, lbpIm, rects[i]);
00287 hist.Normalize();
00288 if(dummy != 0)
00289 ILOG_ERROR("lbp with map should not create outliers in histogram");
00290 accu = Append(accu, hist);
00291 }
00292 }
00293 ILOG_DEBUG("feature length: " << accu.Size());
00294 Feature::FeatureTable* t = mTableSet->GetTable(i);
00295 t->Add(quid, accu);
00296 delete lbpIm;
00297 }
00298 delete greyim;
00299 }
00300
00301 virtual void
00302 HandleDoneFile(VideoSet* vs, int fileId, Stream::RgbDataSrc* src)
00303 {
00304 for (int i=0 ; i<mTableSet->Size() ; i++)
00305 {
00306 if (mLoc.GetIsPartial() && (mTableSet->GetTable(i)->Size() < 2))
00307 continue;
00308 Feature::FeatureDefinition def = mTableSet->GetFeatureDefinition(i);
00309 mLoc.SetFeatureDef(def);
00310 mLoc.SetContainer(vs->GetContainerFile(fileId));
00311 Persistency::FeatureTableRepository().Add(mLoc,
00312 mTableSet->GetTable(i));
00313 mTableSet->GetTable(i)->SetSize(0);
00314 }
00315 }
00316
00317 virtual void
00318 HandleDoneWalk(VideoSet* vs)
00319 {
00320
00321
00322
00323 }
00324
00325 private:
00326 void
00327 GetRectangles(std::vector<Geometry::Rectangle>& rects, int width, int height, int border)
00328 {
00329 rects.resize(6);
00330 int halfWidth = width/2;
00331 int halfHeight = height/2;
00332 rects[0] = Geometry::Rectangle(border, border, width-border-1, height-border-1);
00333 rects[1] = Geometry::Rectangle(border, border, halfWidth-1, halfHeight-1);
00334 rects[2] = Geometry::Rectangle(halfWidth, border, width-border-1, halfHeight-1);
00335 rects[3] = Geometry::Rectangle(border, halfHeight, halfWidth-1, height-border-1);
00336 rects[4] = Geometry::Rectangle(halfWidth, halfHeight, width-border-1, height-border-1);
00337 int width4 = width/4;
00338 int height4 = height/4;
00339 rects[5] = Geometry::Rectangle(width4, height4, height4+halfWidth-1, height4+halfHeight-1);
00340 }
00341
00342 Feature::FeatureDefinition
00343 GetDefinition(int samples, double radius)
00344 {
00345 Feature::FeatureDefinition def("lbp");
00346 def.AddParameter("samples", samples);
00347 def.AddParameter("radius", radius);
00348 if(mDoMSF)
00349 def.AddParameter("vector", "msf");
00350 else
00351 def.AddParameter("vector", "hist");
00352 return def;
00353 }
00354
00355 int mGridHeight;
00356 int mGridWidth;
00357 bool mDoMSF;
00358 int** mMaps;
00359 Feature::FeatureTableSet* mTableSet;
00360 FeatureLocator mLoc;
00361
00362 ILOG_VAR_DEC;
00363 };
00364
00365 ILOG_VAR_INIT(LbpEval, Impala.Core.VideoSet);
00366
00367
00368 }
00369 }
00370 }
00371
00372 #endif