The TriangularInterpolation class. More...

Collaboration diagram for TriangularInterpolation:

[legend]

Public Member Functions
const vector< xy_t > &	getPoints () const
	returns list of verticies

const vector< triangle_t > &	getTriangles () const
	returns list of triangles

	TriangularInterpolation ()=default
	Default constructor.

	TriangularInterpolation (vector< xy_t > &pc, vector< triangle_t > &ts, double d)
	More complex constructor.

	~TriangularInterpolation ()=default
	Destructor.

void	init (vector< xy_t > &points, vector< triangle_t > &triangles, double d)
	Calculate extents of a triangular mesh and build spatial index.

xy_t	triangleCenter (const triangle_t &p) const
	Calculate triangle center.

void	makeIndex (int nx, double xmin, double xmax, int ny, double ymin, double ymax)
	Make spatial index.

short int	sideCross (short int prev, short int curr, const xy_t &r, const xy_t &v0) const
	Determine which triangle side is crossed by a line segment defined by r and v0 points.

void	weights (short int i, const xy_t &r, double &w0, double &w1, double &w2) const
	Calculate barycentric coordinates of a point inside triangle.

short int	findTriangle (const xy_t &r0) const
	Find the triangle which contain the point.

Protected Attributes
vector< triangle_t >	m_triangles
	Triangle list.

vector< xy_t >	m_points
	Vertex list.

vector< xy_t >	m_triangleCenters
	Triangle centers.

vector< double >	m_triangleAreas
	Triangle areas.

vector< short int >	m_spatialIndex
	Spatial index.

double	m_xmin
	Border of the region where the spatial index is constructed.

double	m_xmax
	Border of the region where the spatial index is constructed.

double	m_ymin
	Border of the region where the spatial index is constructed.

double	m_ymax
	Border of the region where the spatial index is constructed.

unsigned int	m_nx
	Spatial index grid size.

unsigned int	m_ny
	Spatial index grid size.

double	m_ixnorm {1}
	Reciprocals to speedup the index calculation.

double	m_iynorm {1}
	Reciprocals to speedup the index calculation.

Detailed Description

The TriangularInterpolation class.

This class travers triangular meshes which satisfies the Delaunay condition. In other meshes it may give a wrong result. The mesh is represented by a list of triangles coupled to a list of xy-points. To speed up the traverse a sort of spatial index is constructed where on a regular cartesian grid the closeses triangle center is memorized.

Definition at line 61 of file BFieldComponentBeamline.cc.

Constructor & Destructor Documentation

◆ TriangularInterpolation()

TriangularInterpolation	(	vector< xy_t > &	pc,
		vector< triangle_t > &	ts,
		double	d
	)

inline

More complex constructor.

Definition at line 72 of file BFieldComponentBeamline.cc.

    {
      init(pc, ts, d);
    }

Member Function Documentation

◆ findTriangle()

short int findTriangle ( const xy_t & r0 ) const

inline

Find the triangle which contain the point.

If not returns the closest one. First using the spatial index locate triangle close to the point and then traverse the mesh using Delaunay triangulation properties.

Parameters

r0	2d cartesian point

Returns: Triangle index in the list

Definition at line 232 of file BFieldComponentBeamline.cc.

    {
      unsigned int ix = lrint((r0.x - m_xmin) * m_ixnorm), iy = lrint((r0.y - m_ymin) * m_iynorm);
      short int curr = (ix < m_nx && iy < m_ny) ? m_spatialIndex[iy + m_ny * ix] : 0;
      xy_t r = m_triangleCenters[curr];
      const short int end = m_triangles.size();
      short int prev = end;
      do {
        short int next = sideCross(prev, curr, r, r0);
        if (next == end) break;
        prev = curr;
        curr = next;
      } while (true);
      return curr;
    }

◆ getPoints()

const vector< xy_t > & getPoints ( ) const

inline

returns list of verticies

Definition at line 64 of file BFieldComponentBeamline.cc.

64{ return m_points;}

Belle2::TriangularInterpolation::m_points

vector< xy_t > m_points

Vertex list.

Definition: BFieldComponentBeamline.cc:251

◆ getTriangles()

const vector< triangle_t > & getTriangles ( ) const

inline

returns list of triangles

Definition at line 66 of file BFieldComponentBeamline.cc.

66{ return m_triangles;}

◆ init()

void init	(	vector< xy_t > &	points,
		vector< triangle_t > &	triangles,
		double	d
	)

inline

Calculate extents of a triangular mesh and build spatial index.

Moves vector contents inside the class.

Parameters

points	List of verticies
triangles	List of triangles
d	Hint how close spatial index should be built

Definition at line 88 of file BFieldComponentBeamline.cc.

    {
      std::swap(points, m_points);
      std::swap(triangles, m_triangles);
      const double inf = numeric_limits<double>::infinity();
      double xmin = inf, xmax = -inf;
      double ymin = inf, ymax = -inf;
      auto limit = [&xmin, &xmax, &ymin, &ymax](const xy_t & p) -> void {
        xmin = std::min(xmin, p.x); xmax = std::max(xmax, p.x);
        ymin = std::min(ymin, p.y); ymax = std::max(ymax, p.y);
      };
      for (const auto& t : m_triangles) {
        limit(m_points[t.j0]);
        limit(m_points[t.j1]);
        limit(m_points[t.j2]);
      }
      double xw = xmax - xmin;
      double yw = ymax - ymin;
 
      xmin -= xw * (1. / 256), xmax += xw * (1. / 256);
      ymin -= yw * (1. / 256), ymax += yw * (1. / 256);
      int nx = lrint(xw * (1 + 1. / 128) / d), ny = lrint(yw * (1 + 1. / 128) / d);
      nx = std::max(nx, 2);
      ny = std::max(ny, 2);
      makeIndex(nx, xmin, xmax, ny, ymin, ymax);
    }

◆ makeIndex()

void makeIndex	(	int	nx,
		double	xmin,
		double	xmax,
		int	ny,
		double	ymin,
		double	ymax
	)

inline

Make spatial index.

Parameters

nx	Grid size in X axis
xmin	Low border on in X
xmax	Upper border on in X
ny	Grid size in Y axis
ymin	Low border on in Y
ymax	Upper border on in Y

Definition at line 136 of file BFieldComponentBeamline.cc.

    {
      m_nx = nx, m_ny = ny, m_xmin = xmin, m_xmax = xmax, m_ymin = ymin, m_ymax = ymax;
      m_ixnorm = (nx - 1) / (xmax - xmin), m_iynorm = (ny - 1) / (ymax - ymin);
      double dx = (xmax - xmin) / (nx - 1), dy = (ymax - ymin) / (ny - 1);
      //    cout<<nx<<" "<<xmin<<" "<<xmax<<" "<<ny<<" "<<ymin<<" "<<ymax<<endl;
      m_spatialIndex.resize(nx * ny);
      m_triangleCenters.resize(m_triangles.size());
      for (unsigned int i = 0; i < m_triangles.size(); i++) m_triangleCenters[i] = triangleCenter(m_triangles[i]);
 
      m_triangleAreas.resize(m_triangles.size());
      auto getTriangleArea = [this](const triangle_t& p) ->double{
        const xy_t& p0 = m_points[p.j0], &p1 = m_points[p.j1], &p2 = m_points[p.j2];
        double d21x = p2.x - p1.x, d21y = p2.y - p1.y;
        double d01x = p0.x - p1.x, d01y = p0.y - p1.y;
        return 1 / (d21x * d01y - d21y * d01x);
      };
      for (unsigned int i = 0; i < m_triangles.size(); i++) m_triangleAreas[i] = getTriangleArea(m_triangles[i]);
 
 
      for (int ix = 0; ix < nx; ix++) {
        double x = xmin + ix * dx;
        for (int iy = 0; iy < ny; iy++) {
          double y = ymin + iy * dy;
          int imin = -1; double dmin = 1e100;
          for (unsigned int i = 0; i < m_triangleCenters.size(); i++) {
            const xy_t& p = m_triangleCenters[i];
            double d = pow(p.x - x, 2) + pow(p.y - y, 2);
            if (d < dmin) {imin = i; dmin = d;}
          }
          int k = iy + ny * ix;
          m_spatialIndex[k] = imin;
        }
      }
    }

◆ sideCross()

short int sideCross	(	short int	prev,
		short int	curr,
		const xy_t &	r,
		const xy_t &	v0
	)		const

inline

Determine which triangle side is crossed by a line segment defined by r and v0 points.

Parameters

prev	Triangle number which has been already checked
curr	Triangle number which is being checked
r	Starting point of the line segment
v0	Ending point of the line segment

Returns: Next triangle index in the list if nothing found returns the total number of triangles in the list

Definition at line 181 of file BFieldComponentBeamline.cc.

    {
      const double vx = r.x - v0.x, vy = r.y - v0.y;
      auto isCrossed = [&vx, &vy, &v0](const xy_t & p1, const xy_t & p0) -> bool {
        double u0x = p0.x, u0y = p0.y;
        double ux = p1.x - u0x, uy = p1.y - u0y;
        double dx = u0x - v0.x, dy = u0y - v0.y;
        double D = uy * vx - ux * vy;
        double t = dx * vy - dy * vx;
        double s = dx * uy - dy * ux;
        return ((t < D) != (t < 0)) && ((s < D) != (s < 0));
      };
 
      const triangle_t& p = m_triangles[curr];
      const xy_t& p0 = m_points[p.j0], &p1 = m_points[p.j1], &p2 = m_points[p.j2];
      if (p.n0 != prev && isCrossed(p1, p2)) return p.n0;
      if (p.n1 != prev && isCrossed(p2, p0)) return p.n1;
      if (p.n2 != prev && isCrossed(p0, p1)) return p.n2;
      return m_triangles.size();
    }

◆ triangleCenter()

xy_t triangleCenter ( const triangle_t & p ) const

inline

Calculate triangle center.

Parameters

p Triangle

Definition at line 120 of file BFieldComponentBeamline.cc.

    {
      const xy_t& p0 = m_points[p.j0], &p1 = m_points[p.j1], &p2 = m_points[p.j2];
      return {(p0.x + p1.x + p2.x)* (1. / 3), (p0.y + p1.y + p2.y)* (1. / 3)};
    }

◆ weights()

void weights	(	short int	i,
		const xy_t &	r,
		double &	w0,
		double &	w1,
		double &	w2
	)		const

inline

Calculate barycentric coordinates of a point inside triangle.

Parameters

i	Triangle index in the list
r	2d cartesian point
w0	Weight of 0 vertex
w1	Weight of 1 vertex
w2	Weight of 2 vertex

Definition at line 211 of file BFieldComponentBeamline.cc.

    {
      const triangle_t& p = m_triangles[i];
      const xy_t& p0 = m_points[p.j0], &p1 = m_points[p.j1], &p2 = m_points[p.j2];
      double dx2  = p2.x -  r.x, dy2  = p2.y -  r.y;
      double d21x = p2.x - p1.x, d21y = p2.y - p1.y;
      double d02x = p0.x - p2.x, d02y = p0.y - p2.y;
      w0 = (dx2 * d21y - dy2 * d21x) * m_triangleAreas[i];
      w1 = (dx2 * d02y - dy2 * d02x) * m_triangleAreas[i];
      w2 = 1 - (w0 + w1);
    }

Member Data Documentation

◆ m_ixnorm

double m_ixnorm {1}

protected

Reciprocals to speedup the index calculation.

Definition at line 271 of file BFieldComponentBeamline.cc.

◆ m_iynorm

double m_iynorm {1}

protected

Reciprocals to speedup the index calculation.

Definition at line 273 of file BFieldComponentBeamline.cc.

◆ m_nx

unsigned int m_nx

protected

Spatial index grid size.

Definition at line 267 of file BFieldComponentBeamline.cc.

◆ m_ny

unsigned int m_ny

protected

Spatial index grid size.

Definition at line 269 of file BFieldComponentBeamline.cc.

◆ m_points

vector<xy_t> m_points

protected

Vertex list.

Definition at line 251 of file BFieldComponentBeamline.cc.

◆ m_spatialIndex

vector<short int> m_spatialIndex

protected

Spatial index.

Definition at line 257 of file BFieldComponentBeamline.cc.

◆ m_triangleAreas

vector<double> m_triangleAreas

protected

Triangle areas.

Definition at line 255 of file BFieldComponentBeamline.cc.

◆ m_triangleCenters

vector<xy_t> m_triangleCenters

protected

Triangle centers.

Definition at line 253 of file BFieldComponentBeamline.cc.

◆ m_triangles

vector<triangle_t> m_triangles

protected

Triangle list.

Definition at line 249 of file BFieldComponentBeamline.cc.

◆ m_xmax

double m_xmax

protected

Border of the region where the spatial index is constructed.

Definition at line 261 of file BFieldComponentBeamline.cc.

◆ m_xmin

double m_xmin

protected

Border of the region where the spatial index is constructed.

Definition at line 259 of file BFieldComponentBeamline.cc.

◆ m_ymax

double m_ymax

protected

Border of the region where the spatial index is constructed.

Definition at line 265 of file BFieldComponentBeamline.cc.

◆ m_ymin

double m_ymin

protected

Border of the region where the spatial index is constructed.

Definition at line 263 of file BFieldComponentBeamline.cc.

The documentation for this class was generated from the following file:

geometry/bfieldmap/src/BFieldComponentBeamline.cc

Public Member Functions

Protected Attributes

Detailed Description

Constructor & Destructor Documentation

◆ TriangularInterpolation()

Member Function Documentation

◆ findTriangle()

◆ getPoints()

◆ getTriangles()

◆ init()

◆ makeIndex()

◆ sideCross()

◆ triangleCenter()

◆ weights()

Member Data Documentation

◆ m_ixnorm

◆ m_iynorm

◆ m_nx

◆ m_ny

◆ m_points

◆ m_spatialIndex

◆ m_triangleAreas

◆ m_triangleCenters

◆ m_triangles

◆ m_xmax

◆ m_xmin

◆ m_ymax

◆ m_ymin