From ef3bca7518f7946bdf904a5baa4f7d682dbf1603 Mon Sep 17 00:00:00 2001 From: Joshua Moerman Date: Sun, 10 Jul 2011 21:20:10 +0200 Subject: [PATCH] yet another version... --- kernels/UnravelHeart3D.hpp | 25 ++++++++++++++++--------- 1 file changed, 16 insertions(+), 9 deletions(-) diff --git a/kernels/UnravelHeart3D.hpp b/kernels/UnravelHeart3D.hpp index 8c1eb90..5ad0a77 100644 --- a/kernels/UnravelHeart3D.hpp +++ b/kernels/UnravelHeart3D.hpp @@ -20,20 +20,27 @@ public: vectorNew[1] = parameters[2]*(vectorOld[0] + parameters[3]); vectorNew[2] = parameters[4]*(vectorOld[1] + parameters[5]); - const double x = vectorNew[0] / parameters[6]; - const double y = vectorNew[1] / parameters[6]; - const double z = vectorNew[2] / parameters[6]; - const double dist = x*x + y*y + 2.0*z*z; - - if(parameters[7] * dist - parameters[7] > x*x*y*y*y) { - const double sqrtDist = std::sqrt(dist); - const double p = 1.0 - parameters[6] * (static_cast(sqrtDist / parameters[6]) + 1.0) / sqrtDist; + const double x = vectorNew[0]; + const double y = vectorNew[1]; + const double z = vectorNew[2]; + const double dist = x*x + y*y + 2.0*z*z - 1.0;//parameters[7]; + const double lh = dist*dist*dist; + const double rh = parameters[6]*x*x*y*y*y; + + if(lh > rh) { + const double sqrtDist = std::pow(std::abs(lh), 1.0/6.0); + const double sqrtDist2 = std::pow(std::abs(rh), 1.0/6.0) * sign(parameters[6]); + const double p = 1.0 - sqrtDist2 * (static_cast(sqrtDist / sqrtDist2) + 1.0) / sqrtDist; vectorNew[0] *= p; - vectorNew[1] *= p; + vectorNew[1] *= -p; vectorNew[2] *= p; } } + double sign(double x){ + return (x > 0) - (x < 0); + } + }; #endif // UNRAVEL_HPP