c++ - Where is my kd tree traversal code wrong? -

i optimizing c++ raytracer. i'm tracing single rays through kdtrees. far using havran's recursive algorithm 'b', seems antique , overblown oop. new code short possible (and more optimized compiler , more maintained):

    struct stackelement{         kdtreenode<pt>* node;         float tmax;         array<float, 3> origin;     };      //initializing explicit stack     stack<stackelement> mystack;      //initialize local variables     kdtreenode<pt>* node = tree.root;     array<float, 3> origin {ray.origin[0], ray.origin[1], ray.origin[2]};     const array<float, 3> direction {ray.direction[0], ray.direction[1], ray.direction[2]};     const array<float, 3> invdirection {1.0f / ray.direction[0], 1.0f / ray.direction[1], 1.0f / ray.direction[2]};     float tmax = numeric_limits<float>::max();     float tclosestintersection = numeric_limits<float>::max();     bool notfullytraversed = true;      while(notfullytraversed) {         if (node->isleaf()) {              //test primitives inside leaf             (auto p : node->primitives()) {                 p->intersect(ray, tclosestintersection, intersection, tmax);             }              //test if leaf + empty stack => return             if (nodestack.empty()) {                 notfullytraversed = false;             } else {                 //pop stack                 origin = mystack.top().origin;                 tmax = mystack.top().tmax;                 node = mystack.top().node;                 mystack.pop();             }         } else {             //get axis of node , split plane             const int axis = node->axis();             const float plane = node->splitposition();              //test if ray not parallel plane             if ((fabs(direction[axis]) > epsilon)) {                 const float t = (plane - origin[axis]) * invdirection[axis];                  //case of ray intersecting plane, test both childs                 if (0.0f < t && t < tmax) {                     //traverse near first, far. set tmax = t near                     tmax = t;                     //push far child onto stack                     mystack.push({                                  (origin[axis] > plane )  ? node->leftchild() : node->rightchild() ,                                  tmax - t,                                  {origin[0] + direction[0] * t, origin[1] + direction[1] * t, origin[2] + direction[2] * t}                                  });                 }             }             //in every case: traverse near child first             node = (origin[axis] > plane) ? node->rightchild() : node->leftchild();          }     }      return intersection.found; 

it's not traversing far child enough. miss related case?

one problem small (original, wrong code):

      //traverse near first, far. set tmax = t near       tmax = t;       //push far child onto stack       mystack.push({ ... , tmax - t,  ...    }); 

it pushes 0.0f onto stack exit distance far node, meaning no positive t accepted intersections.

swapping both lines fixes problem.

my resurcive stack trace / decisions still different (havran takes ~25% more iterations) output picture having 99,5% of same pixels. that's inside floating point rounding issues, still not answer question: case not recognized simplified implementation? or operation not (numerically) stable enough in version?