+ static const TriCorner PT_SEG_MAP[] = {
+ P, Q,
+ Q, R,
+ R, P
+ };
+
+ const double eps = THRESHOLD_F*TransformedTriangle::MULT_PREC_F;
+ for (TriSegment seg = PQ; seg <= RP; seg = TriSegment(seg+1))
+ {
+ // Is h coordinate for both end of segment small enough?
+ int pt1 = PT_SEG_MAP[2*seg], pt2 = PT_SEG_MAP[2*seg+1];
+ if (fabs(_coords[5*pt1+3]) < eps && fabs(_coords[5*pt2+3]) < eps)
+ {
+ // If so, perturb x,y and z to reset h to exactly zero.
+ for (auto pt: {pt1, pt2}) // thx C++17
+ {
+ const double correc = _coords[pt*5+3]/3.; // this should be really small!
+ _coords[pt*5+0] += correc;
+ _coords[pt*5+1] += correc;
+ _coords[pt*5+2] += correc;
+ // And then, if x,y or z very close to 0 or 1, snap exactly to tetra corner:
+ for(int d=0; d < 3; d++)
+ {
+ if (fabs(_coords[5*pt+d]) < eps) _coords[5*pt+d] = 0.0;
+ if (fabs(_coords[5*pt+d]-1) < eps) _coords[5*pt+d] = 1.0;
+ }
+ _coords[pt*5+3] = 0.0;
+ }
+ }
+ }