Skip to content

Optimisation of Surface::intersect(...)

Andreas Salzburger requested to merge 666-speed-optimisation-intersect into master

Closes #666 (closed)

Triggered by the following remark in the code:

  // Evaluate (if necessary in terms of boundaries)
  // @todo: speed up isOnSurface - we know that it is on surface

This MR adds a benchmark for Surface::intersect(...) and introduces some optimisations.

  • PlaneSurface during the intersectEstimate(...) asking for insideBounds(...) did a full globalToLocal transformation. This is not necessary, as part of the information is already available in this regard.

This is a test of 100 million intersections with PlaneSurface in master: start

When shifting the inverse into the method it becomes this: medium

And optimising the usage of Eigen a bit, we arrive at: wow

As a cross-check, the final implementation of this MR is faster than a by hand calculation, which is certainly also uglier:

    double lx = tMatrix(0,0)*intersection.position.x()+tMatrix(0,1)*intersection.position.y()+tMatrix(0,2)*intersection.position.z();
    double ly = tMatrix(1,0)*intersection.position.x()+tMatrix(1,1)*intersection.position.y()+tMatrix(1,2)*intersection.position.z();
    if (not insideBounds(Vector2D(lx,ly), bcheck)) {
      intersection.status = Intersection::Status::missed;
    }
  • DiscSurface - following the implementation of the PlaneSurface, we start from:

DiscSurface_-_start

And then re-using the information available inside the intersectionEstimate(...) we arrive at:

DiscSurface_-_medium

The gain is not as impressive, as we spend a significant amount of time on atan2 for the cartesian to global transformation.

  • CylinderSurface:

Here, the last MR introduced some overhead via the second solution being stored in a std::vector - I reverted that to be a simple object, and allowed for maximum of two intersections with a surface for the moment.

Start time lies at >20s for 100 million tests: CylinderSurface_-_start

While after optimisation, the current time is ~4s: CylinderSurface_-_final

One trick was to introduce a bool fullAzimuth() method to the CylinderBounds, which indicates full azimuthal coverage and hence saves us from checking the radial compatibility (when knowing that we are compatible).

  • StrawSurface or LineSurface:

Same play here, start: Straw_-start

And with the same modifications: Straw_-end

Additionally, the boundaryCheck for PerigeeSurface : public LineSurface is always omitted, as a closest approach as a perigee does not have bounds.

Edited by Andreas Salzburger

Merge request reports