@@ -262,27 +262,30 @@ std::pair<double, array<int, 3>> RectLattice::distance(
262262 double y0 {copysign (0.5 * pitch_[1 ], u.y )};
263263 double z0;
264264
265- double d = std::min (
266- u.x != 0.0 ? (x0 - x) / u.x : INFTY, u.y != 0.0 ? (y0 - y) / u.y : INFTY);
267- if (is_3d_) {
268- z0 = copysign (0.5 * pitch_[2 ], u.z );
269- d = std::min (d, u.z != 0.0 ? (z0 - z) / u.z : INFTY);
265+ // Evaluate axial distances.
266+ double dx = u.x != 0.0 ? (x0 - x) / u.x : INFTY;
267+ double dy = u.y != 0.0 ? (y0 - y) / u.y : INFTY;
268+ double dz = is_3d_ ? (u.z != 0.0 ? (z0 - z) / u.z : INFTY) : INFTY;
269+
270+ double d = dx;
271+ array<int , 3 > lattice_trans = {0 ,0 ,0 };
272+ lattice_trans[0 ] = copysign (1 , u.x );
273+
274+ if (dy < d) {
275+ d = dy;
276+ lattice_trans[0 ] = 0 ;
277+ lattice_trans[1 ] = copysign (1 , u.y );
270278 }
271-
272- // Determine which lattice boundaries are being crossed
273- array<int , 3 > lattice_trans = {0 , 0 , 0 };
274- if (u.x != 0.0 && std::abs (x + u.x * d - x0) < FP_PRECISION)
275- lattice_trans[0 ] = copysign (1 , u.x );
276- if (u.y != 0.0 && std::abs (y + u.y * d - y0) < FP_PRECISION)
277- lattice_trans[1 ] = copysign (1 , u.y );
278- if (is_3d_) {
279- if (u.z != 0.0 && std::abs (z + u.z * d - z0) < FP_PRECISION)
280- lattice_trans[2 ] = copysign (1 , u.z );
279+
280+ if (is_3d_ && dz < d) {
281+ d = dz;
282+ lattice_trans[0 ] = 0 ;
283+ lattice_trans[1 ] = 0 ;
284+ lattice_trans[2 ] = copysign (1 , u.z );
281285 }
282286
283287 return {d, lattice_trans};
284288}
285-
286289// ==============================================================================
287290
288291void RectLattice::get_indices (
0 commit comments