{

TriEdge();

TriEdge(int tri_, int edge_);

bool operator<(const TriEdge& other) const;

bool operator==(const TriEdge& other) const;

bool operator!=(const TriEdge& other) const;

friend std::ostream& operator<<(std::ostream& os, const TriEdge& tri_edge);

int tri, edge;

{

XY();

XY(const double& x_, const double& y_);

double angle() const; // Angle in radians with respect to x-axis.

double cross_z(const XY& other) const; // z-component of cross product.

bool is_right_of(const XY& other) const; // Compares x then y.

bool operator==(const XY& other) const;

bool operator!=(const XY& other) const;

XY operator*(const double& multiplier) const;

const XY& operator+=(const XY& other);

const XY& operator-=(const XY& other);

XY operator+(const XY& other) const;

XY operator-(const XY& other) const;

friend std::ostream& operator<<(std::ostream& os, const XY& xy);

double x, y;

{

XYZ(const double& x_, const double& y_, const double& z_);

XYZ cross(const XYZ& other) const;

double dot(const XYZ& other) const;

XYZ operator-(const XYZ& other) const;

friend std::ostream& operator<<(std::ostream& os, const XYZ& xyz);

double x, y, z;

{

BoundingBox();

void add(const XY& point);

void expand(const XY& delta);

// Consider these member variables read-only.

bool empty;

XY lower, upper;

{

ContourLine();

void push_back(const XY& point);

void write() const;

{

typedef py::array_t<double, py::array::c_style | py::array::forcecast> CoordinateArray;

typedef py::array_t<double, py::array::c_style | py::array::forcecast> TwoCoordinateArray;

typedef py::array_t<int, py::array::c_style | py::array::forcecast> TriangleArray;

typedef py::array_t<bool, py::array::c_style | py::array::forcecast> MaskArray;

typedef py::array_t<int, py::array::c_style | py::array::forcecast> EdgeArray;

typedef py::array_t<int, py::array::c_style | py::array::forcecast> NeighborArray;

/* A single boundary is a vector of the TriEdges that make up that boundary

typedef std::vector<TriEdge> Boundary;

typedef std::vector<Boundary> Boundaries;

/* Constructor with optional mask, edges and neighbors. The latter two

Triangulation(const CoordinateArray& x,

const CoordinateArray& y,

const TriangleArray& triangles,

const MaskArray& mask,

const EdgeArray& edges,

const NeighborArray& neighbors,

bool correct_triangle_orientations);

/* Calculate plane equation coefficients for all unmasked triangles from

TwoCoordinateArray calculate_plane_coefficients(const CoordinateArray& z);

// Return the boundaries collection, creating it if necessary.

const Boundaries& get_boundaries() const;

// Return which boundary and boundary edge the specified TriEdge is.

void get_boundary_edge(const TriEdge& triEdge,

int& boundary,

int& edge) const;

/* Return the edges array, creating it if necessary. */

EdgeArray& get_edges();

/* Return the triangle index of the neighbor of the specified triangle

int get_neighbor(int tri, int edge) const;

/* Return the TriEdge that is the neighbor of the specified triangle edge,

TriEdge get_neighbor_edge(int tri, int edge) const;

/* Return the neighbors array, creating it if necessary. */

NeighborArray& get_neighbors();

// Return the number of points in this triangulation.

int get_npoints() const;

// Return the number of triangles in this triangulation.

int get_ntri() const;

/* Return the index of the point that is at the start of the specified

int get_triangle_point(int tri, int edge) const;

int get_triangle_point(const TriEdge& tri_edge) const;

// Return the coordinates of the specified point index.

XY get_point_coords(int point) const;

// Indicates if the specified triangle is masked or not.

bool is_masked(int tri) const;

/* Set or clear the mask array. Clears various derived fields so they are

void set_mask(const MaskArray& mask);

// Debug function to write boundaries.

void write_boundaries() const;

// An edge of a triangulation, composed of start and end point indices.

struct Edge final

{

Edge() : start(-1), end(-1) {}

Edge(int start_, int end_) : start(start_), end(end_) {}

bool operator<(const Edge& other) const {

return start != other.start ? start < other.start : end < other.end;

}

int start, end;

};

/* An edge of a boundary of a triangulation, composed of a boundary index

struct BoundaryEdge final

{

BoundaryEdge() : boundary(-1), edge(-1) {}

BoundaryEdge(int boundary_, int edge_)

: boundary(boundary_), edge(edge_) {}

int boundary, edge;

};

/* Calculate the boundaries collection. Should normally be accessed via

void calculate_boundaries();

/* Calculate the edges array. Should normally be accessed via

void calculate_edges();

/* Calculate the neighbors array. Should normally be accessed via

void calculate_neighbors();

/* Correct each triangle so that the vertices are ordered in an

void correct_triangles();

/* Determine which edge index (0,1 or 2) the specified point index is in

int get_edge_in_triangle(int tri, int point) const;

bool has_edges() const;

bool has_mask() const;

bool has_neighbors() const;

// Variables shared with python, always set.

CoordinateArray _x, _y; // double array (npoints).

TriangleArray _triangles; // int array (ntri,3) of triangle point indices,

// ordered anticlockwise.

// Variables shared with python, may be unset (size == 0).

MaskArray _mask; // bool array (ntri).

// Derived variables shared with python, may be unset (size == 0).

// If unset, are recalculated when needed.

EdgeArray _edges; // int array (?,2) of start & end point indices.

NeighborArray _neighbors; // int array (ntri,3), neighbor triangle indices

// or -1 if no neighbor.

// Variables internal to C++ only.

Boundaries _boundaries;

// Map used to look up BoundaryEdges from TriEdges. Normally accessed via

// get_boundary_edge().

typedef std::map<TriEdge, BoundaryEdge> TriEdgeToBoundaryMap;

TriEdgeToBoundaryMap _tri_edge_to_boundary_map;

{

typedef Triangulation::CoordinateArray CoordinateArray;

typedef Triangulation::TwoCoordinateArray TwoCoordinateArray;

typedef py::array_t<unsigned char> CodeArray;

/* Constructor.

TriContourGenerator(Triangulation& triangulation,

const CoordinateArray& z);

/* Create and return a non-filled contour.

py::tuple create_contour(const double& level);

/* Create and return a filled contour.

py::tuple create_filled_contour(const double& lower_level,

const double& upper_level);

typedef Triangulation::Boundary Boundary;

typedef Triangulation::Boundaries Boundaries;

/* Clear visited flags.

void clear_visited_flags(bool include_boundaries);

/* Convert a non-filled Contour from C++ to Python.

py::tuple contour_line_to_segs_and_kinds(const Contour& contour);

/* Convert a filled Contour from C++ to Python.

py::tuple contour_to_segs_and_kinds(const Contour& contour);

/* Return the point on the specified TriEdge that intersects the specified

XY edge_interp(int tri, int edge, const double& level);

/* Find and follow non-filled contour lines that start and end on a

void find_boundary_lines(Contour& contour,

const double& level);

/* Find and follow filled contour lines at either of the specified contour

void find_boundary_lines_filled(Contour& contour,

const double& lower_level,

const double& upper_level);

/* Find and follow lines at the specified contour level that are

void find_interior_lines(Contour& contour,

const double& level,

bool on_upper);

/* Follow contour line around boundary of the Triangulation from the

bool follow_boundary(ContourLine& contour_line,

TriEdge& tri_edge,

const double& lower_level,

const double& upper_level,

bool on_upper);

/* Follow contour line across interior of Triangulation.

void follow_interior(ContourLine& contour_line,

TriEdge& tri_edge,

bool end_on_boundary,

const double& level,

bool on_upper);

// Return the Triangulation boundaries.

const Boundaries& get_boundaries() const;

/* Return the edge by which the a level leaves a particular triangle,

int get_exit_edge(int tri, const double& level, bool on_upper) const;

// Return the z-value at the specified point index.

const double& get_z(int point) const;

/* Return the point at which the a level intersects the line connecting the

XY interp(int point1, int point2, const double& level) const;

// Variables shared with python, always set.

Triangulation _triangulation;

CoordinateArray _z; // double array (npoints).

// Variables internal to C++ only.

typedef std::vector<bool> InteriorVisited; // Size 2*ntri

typedef std::vector<bool> BoundaryVisited;

typedef std::vector<BoundaryVisited> BoundariesVisited;

typedef std::vector<bool> BoundariesUsed;

InteriorVisited _interior_visited;

BoundariesVisited _boundaries_visited; // Only used for filled contours.

BoundariesUsed _boundaries_used; // Only used for filled contours.

{

typedef Triangulation::CoordinateArray CoordinateArray;

typedef py::array_t<int, py::array::c_style | py::array::forcecast> TriIndexArray;

/* Constructor. A separate call to initialize() is required to initialize

TrapezoidMapTriFinder(Triangulation& triangulation);

~TrapezoidMapTriFinder();

/* Return an array of triangle indices. Takes 1D arrays x and y of

TriIndexArray find_many(const CoordinateArray& x, const CoordinateArray& y);

/* Return a reference to a new python list containing the following

py::list get_tree_stats();

/* Initialize this object before use. May be called multiple times, if,

void initialize();

// Print the search tree as text to stdout; useful for debug purposes.

void print_tree();

/* A Point consists of x,y coordinates as well as the index of a triangle

struct Point final : XY

{

Point() : XY(), tri(-1) {}

Point(const double& x, const double& y) : XY(x,y), tri(-1) {}

explicit Point(const XY& xy) : XY(xy), tri(-1) {}

int tri;

};

/* An Edge connects two Points, left and right. It is always true that

struct Edge final

{

Edge(const Point* left_,

const Point* right_,

int triangle_below_,

int triangle_above_,

const Point* point_below_,

const Point* point_above_);

// Return -1 if point to left of edge, 0 if on edge, +1 if to right.

int get_point_orientation(const XY& xy) const;

// Return slope of edge, even if vertical (divide by zero is OK here).

double get_slope() const;

/* Return y-coordinate of point on edge with specified x-coordinate.

double get_y_at_x(const double& x) const;

// Return true if the specified point is either of the edge end points.

bool has_point(const Point* point) const;

bool operator==(const Edge& other) const;

friend std::ostream& operator<<(std::ostream& os, const Edge& edge)

{

return os << *edge.left << "->" << *edge.right;

}

void print_debug() const;

const Point* left; // Not owned.

const Point* right; // Not owned.

int triangle_below; // Index of triangle below (to right of) Edge.

int triangle_above; // Index of triangle above (to left of) Edge.

const Point* point_below; // Used only for resolving ambiguous cases;

const Point* point_above; // is 0 if corresponding triangle is -1

};

class Node; // Forward declaration.

// Helper structure used by TrapezoidMapTriFinder::get_tree_stats.

struct NodeStats final

{

NodeStats()

: node_count(0), trapezoid_count(0), max_parent_count(0),

max_depth(0), sum_trapezoid_depth(0.0)

{}

long node_count, trapezoid_count, max_parent_count, max_depth;

double sum_trapezoid_depth;

std::set<const Node*> unique_nodes, unique_trapezoid_nodes;

};

struct Trapezoid; // Forward declaration.

/* Node of the trapezoid map search tree. There are 3 possible types:

class Node final

{

public:

Node(const Point* point, Node* left, Node* right);// Type_XNode.

Node(const Edge* edge, Node* below, Node* above); // Type_YNode.

Node(Trapezoid* trapezoid); // Type_TrapezoidNode.

~Node();

void add_parent(Node* parent);

/* Recurse through the search tree and assert that everything is valid.

void assert_valid(bool tree_complete) const;

// Recurse through the tree to return statistics about it.

void get_stats(int depth, NodeStats& stats) const;

// Return the index of the triangle corresponding to this node.

int get_tri() const;

bool has_child(const Node* child) const;

bool has_no_parents() const;

bool has_parent(const Node* parent) const;

/* Recurse through the tree and print a textual representation to

void print(int depth = 0) const;

/* Remove a parent from this Node. Return true if no parents remain

bool remove_parent(Node* parent);

void replace_child(Node* old_child, Node* new_child);

// Replace this node with the specified new_node in all parents.

void replace_with(Node* new_node);

/* Recursive search through the tree to find the Node containing the

const Node* search(const XY& xy);

/* Recursive search through the tree to find the Trapezoid containing

Trapezoid* search(const Edge& edge);

/* Copy constructor and assignment operator defined but not implemented

Node(const Node& other);

Node& operator=(const Node& other);

private:

typedef enum {

Type_XNode,

Type_YNode,

Type_TrapezoidNode

} Type;

Type _type;

union {

struct {

const Point* point; // Not owned.

Node* left; // Owned.

Node* right; // Owned.

} xnode;

struct {

const Edge* edge; // Not owned.

Node* below; // Owned.

Node* above; // Owned.

} ynode;

Trapezoid* trapezoid; // Owned.

} _union;

typedef std::list<Node*> Parents;

Parents _parents; // Not owned.

};

/* A Trapezoid is bounded by Points to left and right, and Edges below and

struct Trapezoid final

{

Trapezoid(const Point* left_,

const Point* right_,

const Edge& below_,

const Edge& above_);

/* Assert that this Trapezoid is valid. Reduces to a no-op if NDEBUG

void assert_valid(bool tree_complete) const;

/* Return one of the 4 corner points of this Trapezoid. Only used for

XY get_lower_left_point() const;

XY get_lower_right_point() const;

XY get_upper_left_point() const;

XY get_upper_right_point() const;

void print_debug() const;

/* Set one of the 4 neighbor trapezoids and the corresponding reverse

void set_lower_left(Trapezoid* lower_left_);

void set_lower_right(Trapezoid* lower_right_);

void set_upper_left(Trapezoid* upper_left_);

void set_upper_right(Trapezoid* upper_right_);

/* Copy constructor and assignment operator defined but not implemented

Trapezoid(const Trapezoid& other);

Trapezoid& operator=(const Trapezoid& other);

const Point* left; // Not owned.

const Point* right; // Not owned.

const Edge& below;

const Edge& above;

// 4 neighboring trapezoids, can be 0, not owned.

Trapezoid* lower_left; // Trapezoid to left that shares below

Trapezoid* lower_right; // Trapezoid to right that shares below

Trapezoid* upper_left; // Trapezoid to left that shares above

Trapezoid* upper_right; // Trapezoid to right that shares above

Node* trapezoid_node; // Node that owns this Trapezoid.

};

// Add the specified Edge to the search tree, returning true if successful.

bool add_edge_to_tree(const Edge& edge);

// Clear all memory allocated by this object.

void clear();

// Return the triangle index at the specified point, or -1 if no triangle.

int find_one(const XY& xy);

/* Determine the trapezoids that the specified Edge intersects, returning

bool find_trapezoids_intersecting_edge(const Edge& edge,

std::vector<Trapezoid*>& trapezoids);

// Variables shared with python, always set.

Triangulation& _triangulation;

// Variables internal to C++ only.

Point* _points; // Array of all points in triangulation plus corners of

// enclosing rectangle. Owned.

typedef std::vector<Edge> Edges;

Edges _edges; // All Edges in triangulation plus bottom and top Edges of

// enclosing rectangle.

Node* _tree; // Root node of the trapezoid map search tree. Owned.