@@ -276,82 +276,73 @@ void InsertionAlgorithm::insertionPhase()
276276 const type::Vec3 tipToLastCP = lastCP - tipProx->getPosition ();
277277
278278 // Only add a new coupling point if the needle tip has advanced far enough
279- if (tipToLastCP.norm () > tipDistThreshold)
279+ if (tipToLastCP.norm () < tipDistThreshold) return ;
280+
281+ // Prepare operations before entering loop
282+ auto createShaftProximity =
283+ Operations::CreateCenterProximity::Operation::get (l_shaftGeom->getTypeInfo ());
284+ auto projectOnShaft = Operations::Project::Operation::get (l_shaftGeom);
285+ auto findClosestProxOnVol = Operations::FindClosestProximity::Operation::get (l_volGeom);
286+ auto projectOnVol = Operations::Project::Operation::get (l_volGeom);
287+ auto containsPointInVol = Operations::ContainsPointInProximity::Operation::get (l_volGeom);
288+
289+ // Iterate over shaft segments to find which one contains the next candidate CP
290+ for (auto itShaft = l_shaftGeom->begin (); itShaft != l_shaftGeom->end (); itShaft++)
280291 {
281- // Prepare operations before entering loop
282- auto createShaftProximity =
283- Operations::CreateCenterProximity::Operation::get (l_shaftGeom->getTypeInfo ());
284- auto projectOnShaft = Operations::Project::Operation::get (l_shaftGeom);
285- auto findClosestProxOnVol = Operations::FindClosestProximity::Operation::get (l_volGeom);
286- auto projectOnVol = Operations::Project::Operation::get (l_volGeom);
287- auto containsPointInVol = Operations::ContainsPointInProximity::Operation::get (l_volGeom);
288-
289- // Iterate over shaft segments to find which one contains the next candidate CP
290- for (auto itShaft = l_shaftGeom->begin (); itShaft != l_shaftGeom->end (); itShaft++)
291- {
292- BaseProximity::SPtr shaftProx = createShaftProximity (itShaft->element ());
292+ BaseProximity::SPtr shaftProx = createShaftProximity (itShaft->element ());
293293
294- if (!shaftProx) continue ;
294+ if (!shaftProx) continue ;
295295
296- const EdgeProximity::SPtr edgeProx = dynamic_pointer_cast<EdgeProximity>(shaftProx);
296+ const EdgeProximity::SPtr edgeProx = dynamic_pointer_cast<EdgeProximity>(shaftProx);
297297
298- if (!edgeProx) continue ;
298+ if (!edgeProx) continue ;
299299
300- const type::Vec3 p0 = edgeProx->element ()->getP0 ()->getPosition ();
301- const type::Vec3 p1 = edgeProx->element ()->getP1 ()->getPosition ();
302- const type::Vec3 shaftEdgeDir = (p1 - p0).normalized ();
303- const type::Vec3 lastCPToP1 = p1 - lastCP;
300+ const type::Vec3 p0 = edgeProx->element ()->getP0 ()->getPosition ();
301+ const type::Vec3 p1 = edgeProx->element ()->getP1 ()->getPosition ();
302+ const type::Vec3 shaftEdgeDir = (p1 - p0).normalized ();
303+ const type::Vec3 lastCPToP1 = p1 - lastCP;
304304
305- // Skip if last CP lies after edge end point
306- if (dot (shaftEdgeDir, lastCPToP1) < 0_sreal) continue ;
305+ // Skip if last CP lies after edge end point
306+ if (dot (shaftEdgeDir, lastCPToP1) < 0_sreal) continue ;
307307
308- const int numCPs = floor (lastCPToP1.norm () / tipDistThreshold);
308+ const int numCPs = floor (lastCPToP1.norm () / tipDistThreshold);
309309
310- for (int idCP = 0 ; idCP < numCPs; idCP++)
311- {
312- // Candidate coupling point along shaft segment
313- const type::Vec3 candidateCP = lastCP + tipDistThreshold * shaftEdgeDir;
310+ for (int idCP = 0 ; idCP < numCPs; idCP++)
311+ {
312+ // Candidate coupling point along shaft segment
313+ const type::Vec3 candidateCP = lastCP + tipDistThreshold * shaftEdgeDir;
314314
315- // Project candidate CP onto the edge element and compute scalar coordinate
316- // along segment
317- const SReal edgeSegmentLength = (p1 - p0).norm ();
318- const type::Vec3 p0ToCandidateCP = candidateCP - p0;
319- const SReal projPtOnEdge = dot (p0ToCandidateCP, shaftEdgeDir);
315+ // Project candidate CP onto the edge element and compute scalar coordinate
316+ // along segment
317+ const SReal edgeSegmentLength = (p1 - p0).norm ();
318+ const type::Vec3 p0ToCandidateCP = candidateCP - p0;
319+ const SReal projPtOnEdge = dot (p0ToCandidateCP, shaftEdgeDir);
320320
321- // Skip if candidate CP is outside current edge segment
322- if (projPtOnEdge < 0_sreal || projPtOnEdge > edgeSegmentLength) break ;
321+ // Skip if candidate CP is outside current edge segment
322+ if (projPtOnEdge < 0_sreal || projPtOnEdge > edgeSegmentLength) break ;
323323
324- // Project candidate CP onto shaft geometry ...
325- shaftProx = projectOnShaft (candidateCP, itShaft->element ()).prox ;
324+ // Project candidate CP onto shaft geometry ...
325+ shaftProx = projectOnShaft (candidateCP, itShaft->element ()).prox ;
326326
327- if (!shaftProx) continue ;
327+ if (!shaftProx) continue ;
328328
329- // ... then find nearest volume proximity
330- const BaseProximity::SPtr volProx =
331- findClosestProxOnVol (shaftProx, l_volGeom.get (), projectOnVol, getFilterFunc ());
329+ // ... then find nearest volume proximity
330+ const BaseProximity::SPtr volProx =
331+ findClosestProxOnVol (shaftProx, l_volGeom.get (), projectOnVol, getFilterFunc ());
332332
333- if (!volProx) continue ;
333+ if (!volProx) continue ;
334334
335- // Proximity can be detected before the tip enters the tetra (e.g. near a
336- // boundary face) Only accept proximities if the tip is inside the tetra
337- // during insertion
338- if (containsPointInVol (shaftProx->getPosition (), volProx))
339- {
340- volProx->normalize ();
341- m_couplingPts.push_back (volProx);
342- lastCP = volProx->getPosition ();
343- }
335+ // Proximity can be detected before the tip enters the tetra (e.g. near a
336+ // boundary face) Only accept proximities if the tip is inside the tetra
337+ // during insertion
338+ if (containsPointInVol (shaftProx->getPosition (), volProx))
339+ {
340+ volProx->normalize ();
341+ m_couplingPts.push_back (volProx);
342+ lastCP = volProx->getPosition ();
344343 }
345344 }
346345 }
347- else // Don't bother with removing the point that was just added
348- {
349- // Remove coupling points that are ahead of the tip in the insertion direction
350- ElementIterator::SPtr itShaft = l_shaftGeom->begin (l_shaftGeom->getSize () - 2 );
351- auto prunePointsAheadOfTip =
352- Operations::Needle::PrunePointsAheadOfTip::get (itShaft->getTypeInfo ());
353- prunePointsAheadOfTip (m_couplingPts, itShaft->element ());
354- }
355346}
356347
357348InsertionAlgorithm::AlgorithmOutput InsertionAlgorithm::reprojectCouplingPoints ()
0 commit comments