diff --git a/smc/filtering/resampling/ParticleFilterResamplingKDEPercent.h b/smc/filtering/resampling/ParticleFilterResamplingKDEPercent.h new file mode 100644 index 0000000..5a6945f --- /dev/null +++ b/smc/filtering/resampling/ParticleFilterResamplingKDEPercent.h @@ -0,0 +1,136 @@ +#ifndef PARTICLEFILTERRESAMPLINGKDEPERCENT_H +#define PARTICLEFILTERRESAMPLINGKDEPERCENT_H + + +#include +#include + +#include "ParticleFilterResampling.h" +#include "../../ParticleAssertions.h" + +#include "../../../math/boxkde/benchmark.h" +#include "../../../math/boxkde/DataStructures.h" +#include "../../../math/boxkde/Image3D.h" +#include "../../../math/boxkde/BoxGaus3D.h" +#include "../../../math/boxkde/Grid3D.h" +#include "../../../math/distribution/Normal.h" + +#include "../../../navMesh/NavMesh.h" +#include "../../../floorplan/v2/FloorplanHelper.h" + +namespace SMC { + + /** + * Remove the worst of particles, then calculate the KDE and reweight the remaining particles accordingly. + */ + template + class ParticleFilterResamplingKDEPercent : public ParticleFilterResampling { + + private: + + /** this is a copy of the particle-set to draw from it */ + std::vector> particlesCopy; + + /** random number generator */ + std::minstd_rand gen; + + /** boundingBox for the boxKDE */ + _BBox3 bb; + + /** histogram/grid holding the particles*/ + std::unique_ptr> grid; + + /** bandwith for KDE */ + Point3 bandwith; + + /** the current mesh */ + const NM::NavMesh* mesh; + + /** percent to remove **/ + float percent; + + public: + + /** ctor */ + ParticleFilterResamplingKDEPercent(const NM::NavMesh* mesh, const Point3 gridsize_m, const Point3 bandwith, const float percent) { + + this->mesh = mesh; + this->bandwith = bandwith; + this->bb = mesh->getBBox(); + this->bb.grow(10); + + // Create histogram + size_t nBinsX = (size_t)((this->bb.getMax().x - this->bb.getMin().x) / gridsize_m.x); + size_t nBinsY = (size_t)((this->bb.getMax().y - this->bb.getMin().y) / gridsize_m.y); + size_t nBinsZ = (size_t)((this->bb.getMax().z - this->bb.getMin().z) / gridsize_m.z); + + this->grid = std::make_unique>(bb, nBinsX, nBinsY, nBinsZ); + + this->percent = percent; + + gen.seed(1234); + } + + void resample(std::vector>& particles) override { + + // compile-time sanity checks + static_assert( HasOperatorPlusEq::value, "your state needs a += operator!" ); + static_assert( HasOperatorDivEq::value, "your state needs a /= operator!" ); + static_assert( HasOperatorMul::value, "your state needs a * operator!" ); + //static_assert( std::is_constructible::value, "your state needs a constructor with Point3!"); + //todo: static assert for getx, gety, getz, setposition + + + // comparator (highest first) + static auto comp = [] (const Particle& p1, const Particle& p2) { + return p1.weight > p2.weight; + }; + + const uint32_t cnt = (uint32_t) particles.size(); + + // sort particles by weight (highest first) + std::sort(particles.begin(), particles.end(), comp); + + // to-be-removed region + const int start = particles.size() * (1-percent); + const int end = particles.size(); + std::uniform_int_distribution dist(0, start-1); + + // remove by re-drawing + for (uint32_t i = start; i < end; ++i) { + const int rnd = dist(gen); + particles[i] = particles[rnd]; + } + + grid->clear(); + for (Particle p : particles){ + //grid.add receives position in meter! + + //if weight is to low, remove. + if((float) p.weight > 0 ){ + grid->add(p.state.getX(), p.state.getY(), p.state.getZ(), p.weight); + } + } + + // init KDE + int nFilt = 3; + float sigmaX = bandwith.x / grid->binSizeX; + float sigmaY = bandwith.y / grid->binSizeY; + float sigmaZ = bandwith.z / grid->binSizeZ; + + // process KDE + BoxGaus3D boxGaus; + boxGaus.approxGaus(grid->image(), Point3(sigmaX, sigmaY, sigmaZ), nFilt); + + + // reweight the particle using the kde. in theory this should be the same weight. + // however, as the kde "smoothes" the pdf, this values are lower. + for (Particle p : particles){ + p.weight = grid->fetch(p.state.getX(), p.state.getY(), p.state.getZ()); + } + } + }; +} + + +#endif // PARTICLEFILTERRESAMPLINGKDEPERCENT_H