14 #ifndef ESYS_LSMSPHEREFITTER_H
15 #define ESYS_LSMSPHEREFITTER_H
17 #include "Geometry/Sphere3d.h"
18 #include "Geometry/Sphere2d.h"
19 #include "Foundation/BoundingSphere.h"
28 template <
typename TmplFitTraits>
32 typedef TmplFitTraits FitTraits;
33 typedef typename FitTraits::Validator Validator;
34 typedef typename FitTraits::Particle Particle;
35 typedef typename Validator::ParticleVector ParticleVector;
36 typedef typename FitTraits::Plane Plane;
38 SphereFitter(
const std::string &name, Validator &validator)
39 : m_pValidator(&validator),
40 m_successfulFitCount(0),
51 virtual Particle getFitParticle(
52 const Particle &particle,
53 const ParticleVector &neighbours,
57 static Particle getInvalidParticle()
59 return Particle::INVALID;
72 void incrSuccessfulFit()
74 m_successfulFitCount++;
77 const std::string &getName()
const
82 void write(std::ostream &oStream)
const
86 <<
"fit count = " << m_getFitCount
87 <<
", success = " << m_successfulFitCount
88 <<
", fail = " << m_failedFitCount;
91 std::string toString()
const
93 std::stringstream sStream;
98 bool particleIsValid(
const Particle &particle)
const
100 return getValidator().particleIsValid(particle);
104 Validator &getValidator()
106 return *m_pValidator;
109 const Validator &getValidator()
const
111 return *m_pValidator;
114 Validator *m_pValidator;
115 int m_successfulFitCount;
117 int m_failedFitCount;
121 template <
typename TmplFitTraits>
122 inline std::ostream &operator<<(
123 std::ostream &oStream,
127 fitter.write(oStream);
131 template <
typename TmplFitTraits>
132 class MoveToSurfaceFitter :
public SphereFitter<TmplFitTraits>
135 typedef SphereFitter<TmplFitTraits> Inherited;
136 typedef typename Inherited::Validator Validator;
137 typedef typename Inherited::Particle Particle;
138 typedef typename Inherited::ParticleVector ParticleVector;
139 typedef typename Inherited::Plane
Plane;
140 MoveToSurfaceFitter(Validator &validator)
141 : Inherited(
"Mv to Surface", validator)
145 Particle moveToSurface(
146 const Particle &stationary,
147 const Particle &movable
150 Particle moved = movable;
151 const Vec3 centreDiff = movable.getPos() - stationary.getPos();
152 const double centreDist = centreDiff.norm();
153 if (centreDist > 0.0) {
154 const Vec3 newCentrePos =
157 (centreDiff/(centreDist))*(stationary.getRad() + movable.getRad());
158 moved.moveTo(newCentrePos);
163 virtual Particle getFitParticle(
164 const Particle &particle,
165 const ParticleVector &neighbours,
170 Particle newParticle = particle;
171 if (neighbours.size() > 0) {
173 (particle.getPos() - neighbours[0]->getPos()).norm()
175 (particle.getRad() + neighbours[0]->getRad())
177 newParticle = moveToSurface(*(neighbours[0]), particle);
180 if (newParticle.isValid() && !this->particleIsValid(newParticle)) {
181 newParticle = this->getInvalidParticle();
182 this->incrFailedFit();
183 }
else if (newParticle.isValid()) {
184 this->incrSuccessfulFit();
190 template <
typename TmplFitTraits>
195 typedef typename Inherited::Validator Validator;
196 typedef typename Inherited::Particle Particle;
197 typedef typename Inherited::ParticleVector ParticleVector;
198 typedef typename Inherited::Plane Plane;
207 const ParticleVector &particleVector
210 Particle fittedParticle = this->getInvalidParticle();
215 if (particleVector.size() > 3) {
216 const bool foundFit =
218 particleVector[0]->getPos(),
219 particleVector[1]->getPos(),
220 particleVector[2]->getPos(),
221 particleVector[3]->getPos(),
222 particleVector[0]->getRad(),
223 particleVector[1]->getRad(),
224 particleVector[2]->getRad(),
225 particleVector[3]->getRad(),
230 fittedParticle = Particle(M, r,
id);
233 throw std::runtime_error(
"findAFit: particleVector argument has fewer than 4 elements.");
235 return fittedParticle;
238 virtual Particle getFitParticle(
239 const Particle &particle,
240 const ParticleVector &neighbours,
245 Particle newParticle = this->getInvalidParticle();
246 if(neighbours.size() > 3){
247 const Particle &Pi = *(neighbours[0]);
248 const double ndist=(particle.getPos()-Pi.getPos()).norm();
250 newParticle = particle;
251 if (ndist < Pi.getRad()){
253 Pi.getPos()+((particle.getPos()-Pi.getPos())*(Pi.getRad()/ndist))
256 const double dist_p = plane.sep(newParticle.getPos());
257 const double dist_3 = (newParticle.getPos()-neighbours[3]->getPos()).norm()-neighbours[3]->getRad();
258 if (dist_p > dist_3) {
259 newParticle = findAFit(newParticle, neighbours);
262 newParticle = this->getInvalidParticle();
266 if (newParticle.isValid() && !this->particleIsValid(newParticle)) {
267 newParticle = this->getInvalidParticle();
268 this->incrFailedFit();
269 }
else if (newParticle.isValid()) {
270 this->incrSuccessfulFit();
276 template <
typename TmplFitTraits>
281 typedef typename Inherited::Validator Validator;
282 typedef typename Inherited::Particle Particle;
283 typedef typename Inherited::ParticleVector ParticleVector;
284 typedef typename Inherited::Plane Plane;
293 const ParticleVector &particleVector
296 Particle fittedParticle = this->getInvalidParticle();
301 if (particleVector.size() > 2) {
302 const bool foundFit =
304 particleVector[0]->getPos(),
305 particleVector[1]->getPos(),
306 particleVector[2]->getPos(),
307 particleVector[0]->getRad(),
308 particleVector[1]->getRad(),
309 particleVector[2]->getRad(),
314 fittedParticle = Particle(M, r,
id);
317 throw std::runtime_error(
"findAFit: particleVector argument has fewer than 3 elements.");
319 return fittedParticle;
322 virtual Particle getFitParticle(
323 const Particle &particle,
324 const ParticleVector &neighbours,
329 Particle newParticle = this->getInvalidParticle();
330 if(neighbours.size() > 2){
331 const Particle &Pi = *(neighbours[0]);
332 const double ndist=(particle.getPos()-Pi.getPos()).norm();
334 newParticle = particle;
335 if (ndist < Pi.getRad()){
337 Pi.getPos()+((particle.getPos()-Pi.getPos())*(Pi.getRad()/ndist))
340 const double dist_p = plane.sep(newParticle.getPos());
341 const double dist_2 = (newParticle.getPos()-neighbours[2]->getPos()).norm()-neighbours[2]->getRad();
342 if (dist_p > dist_2) {
343 newParticle = findAFit(newParticle, neighbours);
346 newParticle = this->getInvalidParticle();
350 if (newParticle.isValid() && !this->particleIsValid(newParticle)) {
351 newParticle = this->getInvalidParticle();
352 this->incrFailedFit();
353 }
else if (newParticle.isValid()) {
354 this->incrSuccessfulFit();
360 template <
typename TmplFitTraits>
365 typedef typename Inherited::Validator Validator;
366 typedef typename Inherited::Particle Particle;
367 typedef typename Inherited::ParticleVector ParticleVector;
368 typedef typename Inherited::Plane Plane;
371 Validator &validator,
374 :
Inherited(
" 2D Sph Fitter", validator),
381 const ParticleVector &particleVector
384 Particle fittedParticle = this->getInvalidParticle();
389 if (particleVector.size() >= 2) {
390 const bool foundFit =
392 m_sphere.getCentre(),
393 particleVector[0]->getPos(),
394 particleVector[1]->getPos(),
395 -m_sphere.getRadius(),
396 particleVector[0]->getRad(),
397 particleVector[1]->getRad(),
402 fittedParticle = Particle(M, r,
id);
405 throw std::runtime_error(
406 "TwoDSphereSphereFitter::findAFit: particleVector "
407 "argument has fewer than 2 elements."
410 return fittedParticle;
413 virtual Particle getFitParticle(
414 const Particle &particle,
415 const ParticleVector &neighbours,
420 Particle newParticle = this->getInvalidParticle();
421 if (neighbours.size() >= 2) {
422 const Particle &Pi = *(neighbours[0]);
423 const double ndist=(particle.getPos()-Pi.getPos()).norm();
428 (neighbours.size() == 2)
431 fabs((m_sphere.getCentre()-particle.getPos()).norm()-m_sphere.getRadius())
433 (particle.getPos()-neighbours[2]->getPos()).norm()-neighbours[2]->getRad()
437 newParticle = particle;
438 if (ndist < Pi.getRad()) {
440 Pi.getPos() + ((particle.getPos() - Pi.getPos()) * (Pi.getRad()/ndist))
443 newParticle = findAFit(newParticle, neighbours);
446 if (newParticle.isValid() && !this->particleIsValid(newParticle)) {
447 newParticle = this->getInvalidParticle();
448 this->incrFailedFit();
449 }
else if (newParticle.isValid()) {
450 this->incrSuccessfulFit();
459 template <
typename TmplFitTraits>
464 typedef typename Inherited::Validator Validator;
465 typedef typename Inherited::Particle Particle;
466 typedef typename Inherited::ParticleVector ParticleVector;
467 typedef typename Inherited::Plane Plane;
470 Validator &validator,
473 :
Inherited(
" 3D Sph Fitter", validator),
480 const ParticleVector &particleVector
483 Particle fittedParticle = this->getInvalidParticle();
488 if (particleVector.size() > 2) {
489 const bool foundFit =
491 m_sphere.getCentre(),
492 particleVector[0]->getPos(),
493 particleVector[1]->getPos(),
494 particleVector[2]->getPos(),
495 -m_sphere.getRadius(),
496 particleVector[0]->getRad(),
497 particleVector[1]->getRad(),
498 particleVector[2]->getRad(),
503 fittedParticle = Particle(M, r,
id);
508 "ThreeDSphereSphereFitter::findAFit: particleVector argument"
509 " has fewer than 3 elements."
512 return fittedParticle;
515 virtual Particle getFitParticle(
516 const Particle &particle,
517 const ParticleVector &neighbours,
522 Particle newParticle = this->getInvalidParticle();
523 if (neighbours.size() >= 3) {
524 const Particle &Pi = *(neighbours[0]);
525 const double ndist=(particle.getPos()-Pi.getPos()).norm();
530 (neighbours.size() == 3)
533 fabs((m_sphere.getCentre()-particle.getPos()).norm()-m_sphere.getRadius())
535 (particle.getPos()-neighbours[3]->getPos()).norm()-neighbours[3]->getRad()
539 newParticle = particle;
540 if (ndist < Pi.getRad()) {
542 Pi.getPos() + ((particle.getPos() - Pi.getPos()) * (Pi.getRad()/ndist))
545 newParticle = findAFit(newParticle, neighbours);
548 if (newParticle.isValid() && !this->particleIsValid(newParticle)) {
549 newParticle = this->getInvalidParticle();
550 this->incrFailedFit();
551 }
else if (newParticle.isValid()) {
552 this->incrSuccessfulFit();
561 template <
typename TmplFitTraits>
566 typedef typename Inherited::Validator Validator;
567 typedef typename Inherited::Particle Particle;
568 typedef typename Inherited::ParticleVector ParticleVector;
569 typedef typename Inherited::Plane Plane;
577 const Particle& particle,
578 const ParticleVector &particleVector,
582 Particle fittedParticle = this->getInvalidParticle();
585 const int id = particle.getID();
587 if (particleVector.size() > 1) {
588 const bool foundFit =
590 particleVector[0]->getPos(),
591 particleVector[1]->getPos(),
594 particleVector[0]->getRad(),
595 particleVector[1]->getRad(),
600 fittedParticle = Particle(M, r,
id);
603 throw std::runtime_error(
"findAFit: particleVector vector contains less than 2 particles.");
606 return fittedParticle;
609 virtual Particle getFitParticle(
610 const Particle &particle,
611 const ParticleVector &neighbours,
616 Particle newParticle = this->getInvalidParticle();
617 if (neighbours.size() >= 2) {
618 const Particle &Pi = *(neighbours[0]);
619 const double ndist=(particle.getPos()-Pi.getPos()).norm();
624 (neighbours.size() == 2)
627 plane.sep(particle.getPos())
629 (particle.getPos()-neighbours[2]->getPos()).norm()-neighbours[2]->getRad()
633 newParticle = particle;
634 if (ndist < Pi.getRad()) {
636 Pi.getPos() + ((particle.getPos() - Pi.getPos()) * (Pi.getRad()/ndist))
639 newParticle = findAFit(newParticle, neighbours, plane);
642 if (newParticle.isValid() && !this->particleIsValid(newParticle)) {
643 newParticle = this->getInvalidParticle();
644 this->incrFailedFit();
645 }
else if (newParticle.isValid()) {
646 this->incrSuccessfulFit();
652 template <
typename TmplFitTraits>
657 typedef typename Inherited::Validator Validator;
658 typedef typename Inherited::Particle Particle;
659 typedef typename Inherited::ParticleVector ParticleVector;
660 typedef typename Inherited::Plane Plane;
668 const Particle& particle,
669 const ParticleVector &particleVector,
673 Particle fittedParticle = this->getInvalidParticle();
676 const int id = particle.getID();
678 if (particleVector.size() > 2) {
679 const bool foundFit =
681 particleVector[0]->getPos(),
682 particleVector[1]->getPos(),
683 particleVector[2]->getPos(),
686 particleVector[0]->getRad(),
687 particleVector[1]->getRad(),
688 particleVector[2]->getRad(),
693 fittedParticle = Particle(M, r,
id);
696 throw std::runtime_error(
"findAFit: particleVector vector contains less than 3 particles.");
699 return fittedParticle;
702 virtual Particle getFitParticle(
703 const Particle &particle,
704 const ParticleVector &neighbours,
709 Particle newParticle = this->getInvalidParticle();
710 if (neighbours.size() >= 3) {
711 const Particle &Pi = *(neighbours[0]);
712 const double ndist=(particle.getPos()-Pi.getPos()).norm();
717 (neighbours.size() == 3)
720 plane.sep(particle.getPos())
722 (particle.getPos()-neighbours[3]->getPos()).norm()-neighbours[3]->getRad()
726 newParticle = particle;
727 if (ndist < Pi.getRad()) {
729 Pi.getPos() + ((particle.getPos() - Pi.getPos()) * (Pi.getRad()/ndist))
732 newParticle = findAFit(newParticle, neighbours, plane);
735 if (newParticle.isValid() && !this->particleIsValid(newParticle)) {
736 newParticle = this->getInvalidParticle();
737 this->incrFailedFit();
738 }
else if (newParticle.isValid()) {
739 this->incrSuccessfulFit();
747 #include "Geometry/SphereFitter.hpp"