Calculation geometric basis of GIS Algorithms

Vector concept:

If the endpoint of a line segment is ordered, we convert this line segment into a directed line segment ).
If the start point P1 of the directed line P1P2 is at the coordinate origin, we can call it vector P2.

Vector addition and subtraction:

Set two-dimensional vector p = (x1, Y1), q = (X2, Y2 ),
Vector Addition is defined as: p + q = (X1 + X2, Y1 + y2 ),
Vector Subtraction is defined as p-q = (x1-X2, Y1-Y2 ).
Obviously, p + q = q + p, p-q =-(Q-P ).

Vector Cross Product:

Calculation of vector cross product is the core part of algorithms related to straight lines and line segments.
Set vector p = (x1, Y1), q = (X2, Y2 ),
The vector cross product is defined as the signed area of the parallelogram consisting of (0, 0), P1, P2, and P1 + p2,
That is, p × q = x1 * Y2-X2 * Y1, and the result is a scalar.
Obviously, the properties include P x q =-(q x P) and P x (-q) =-(P x q ).
The addition and subtraction of two points is the addition and subtraction of vectors, while the multiplication of vertices is regarded as the cross product of vectors.
A very important property of the cross product is that it can be used to determine the clockwise relationship between two vectors:

If p × q> 0, P is clockwise in Q.

If p × q <0, P is in the counterclockwise direction of Q.

If p × q = 0, p and q are both in the same direction, but they may be reversed.

Identify the direction of a line segment:

The turning Judgment Method of the Line Segment can be directly introduced by the nature of the vector cross product.
For line segments p0p1 and P1P2 with common endpoints, you can determine the direction of the line segment by calculating the symbol (P2-P0) × (P1-P0:

If (P2-P0) × (P1-P0)> 0, p0p1 turns to the right side of P1 and returns P1P2.

If (P2-P0) × (P1-P0) <0, p0p1 returns P1P2 after turning to the left of P1.

If (P2-P0) × (P1-P0) = 0, P0, P1, and P2 are all in line.

For details, refer:

　　Determine whether the vertex is in the online segment:

Set the vertex to Q and the line segment to P1P2,
The basis for determining the vertex Q on this line segment is: (q-P1) × (P2-P1) = 0 and Q is in the rectangle where P1 and P2 are diagonal vertices.
The former ensures that the qpoint is on the P1P2 line,
The latter ensures that the qpoint is not in the extended or reverse extended line of the P1P2 line segment,
The following process can be used to determine this step:

If min (xp1, xp2) <= XQ <= max (xp1, xp2) and min (yp1, yp2) <= yq <= max (yp1, yp2)
{Return true };

Else
{Return false ;}

Note that horizontal and vertical line segments must be considered,
Min (XI, XJ) <= XK <= max (XI, XJ) and min (Yi, YJ) <= yk <= max (Yi, YJ) only when both conditions are met can the true value be returned.

　　Determine whether two line segments are Intersection:

We can determine whether two line segments overlap in two steps:

(1) Rapid Rejection Test

Set the rectangle with the linestring P1P2 as the diagonal line to R, and the rectangle with the linestring q1q2 as the diagonal line to T. If the R and t do not intersect, the two linestring will not.

(2) Cross-Site Test

If the two line segments overlap, the two lines must cross each other.
If P1P2 cross-sets q1q2, the vectors (P1-Q1) and (P2-Q1) are on both sides of the vector (Q2-Q1,
That is, (P1-Q1) × (Q2-Q1) * (P2-Q1) × (Q2-Q1) <0.
The preceding formula can be changed to (P1-Q1) × (Q2-Q1) * (Q2-Q1) × (P2-Q1)> 0.
When (P1-Q1) × (Q2-Q1) = 0, it indicates that (P1-Q1) and (Q2-Q1) are both linear,
However, the P1 must be on the q1q2 online segment because it has passed the fast Rejection Test;
Similarly, (Q2-Q1) × (P2-Q1) = 0 indicates that P2.
Therefore, the basis for determining the cross-q1q2 of P1P2 is: (P1-Q1) × (Q2-Q1) * (Q2-Q1) × (P2-Q1)> = 0.
Similarly, the basis for determining cross-site P1P2 in q1q2 is: (Q1-P1) × (P2-P1) * (P2-P1) × (Q2-P1)> = 0.

　　Determine whether a line segment and a straight line are Intersection:

If the P1P2 line segments and the q1q2 Line Intersect, P1P2 cross-site q1q2,
That is, (P1-Q1) × (Q2-Q1) * (Q2-Q1) × (P2-Q1)> = 0.

　　Determines whether a rectangle contains vertices.:

You only need to determine whether the horizontal and vertical coordinates of the point are between the left side of the rectangle and the upper and lower sides.

 

　　Determines whether a line segment, line segment, or polygon is in a rectangle.:

Because the rectangle is a convex set, you only need to determine whether all endpoints are in the rectangle.

　　Determines whether a rectangle is in a rectangle.:

You only need to compare the left and right boundary with the upper and lower boundary.

　　Judge whether the circle is in the rectangle:

It is easy to prove that the necessary and sufficient condition of a circle in a rectangle is that the center of the circle is in a rectangle and the radius of the circle is smaller than or equal to the minimum distance from the center to the four sides of the rectangle.

　　Determines whether a vertex is in a polygon.:

Determining whether point P is in a polygon is a very basic but important algorithm for geometric calculation. Take point P as the endpoint and take the ray L to the left. Because the polygon is bounded, the left end of the Ray l must be outside the polygon, consider moving from left to right starting from infinity along l. When we encounter the first intersection with the polygon, we enter the inside of the polygon. When we encounter the second intersection, we leave the polygon, ...... So it is easy to see that when the intersection of L and polygon number C is an odd number, P is inside the polygon, and P is outside the polygon if it is an even number.

However, some special cases should be taken into consideration. (A) (B) (c) (d. In figure (A), the vertices of L and polygon intersect. At this time, only one intersection can be calculated. In Figure (B), the intersection of L and polygon vertices should not be calculated; in the image (c) and (d), an edge of the L and the polygon overlaps, and this edge should be ignored. If l and a polygon side overlap, this side should be ignored.

For the sake of unification, when we calculate the intersection of Ray l and polygon, 1. Do not consider horizontal edges of polygon; 2. When the vertices of a polygon intersect with L, if the vertex is a vertex with a large vertical coordinate on the edge of the polygon, It is counted; otherwise, it is ignored; 3. In the case that P is located on the side of a polygon, we can directly determine that p is a multilateral row. The pseudo code of the algorithm is as follows:

Count distinct 0;

Taking P as the endpoint, as the right-to-left Ray l;

For each side of a polygon s

Do if p on edge s

Then return true;

If S is not horizontal

One endpoint of then if S is on L.

If the endpoint is the endpoint with a large vertical coordinate in the points at both ends of S

Then count distinct count + 1

Else if S and l Intersection

Then count distinct count + 1;

If Count mod 2 = 1

Then return true;

Else return false;

The method for Ray l is as follows: if the Y coordinate of P is the same as that of P, and the Y coordinate is positive infinity (a large positive number), P and P determine the ray L.

The time complexity of this algorithm for determining whether a vertex is in a polygon is O (n ).

Another algorithm is to compare the area of a signed triangle with the area of a polygon. This algorithm may cause some errors due to the use of floating point numbers. We do not recommend this algorithm.

　　Determines whether a line segment is within a polygon.:

A necessary condition for a line segment within a polygon is that the two ends of the line segment are in the polygon. However, because the polygon may be concave, this cannot be a sufficient condition for judgment. If the intersection of a line segment and a polygon side (the intersection of a two-line segment refers to the intersection of two lines and the intersection is not at the endpoint of the two lines ), because the left and right sides of a polygon edge belong to different parts inside and outside the polygon, the line segment must have a part outside the polygon (see figure ). So we get the second necessary condition of the line segment in the polygon: The line segment and all sides of the polygon are not inner.

The intersection of a line segment and a polygon at both ends of a line segment does not affect whether the line segment is in a polygon. However, if a vertex of a polygon is at the intersection of a line segment, you must also determine whether the line segments between two adjacent intersections are included in the polygon (see Figure B in the inverse example ).

Therefore, we can first find all the vertices of the polygon that intersect with the line segment, and then sort by the X-Y coordinate (the small x coordinate is in front, for the point with the same X coordinate, Y coordinates are arranged at the front. This sorting criterion is used to ensure that the horizontal and vertical conditions are correct. In this way, the adjacent two points are the adjacent two intersections on the online segment, if the midpoint of any adjacent two points is within the polygon, the line segment must be within the polygon.

The proof is as follows:

　

Proposition 1:

If the intersection of a line segment and a polygon is P1 and the point P of P2 is also within the polygon, all vertices between P1 and P2 are within the polygon.

 

Proof:

Assume that points P1 and P2 do not exist in the polygon. Set the point to Q, which is between P1 and P'. Because the polygon is a closed curve, the interior and exterior are bounded, while P1 is inside the multilateral row, Q is outside the multilateral, p' is inside the multilateral, P1-Q-P is completely continuous, so P1q and qP 'must span the boundary of polygon, so in P1, there are at least two intersection points between the line segment and the polygon. This is in conflict with P1P2, so the proposition is true. Pass.

Inferences can be drawn directly from Proposition 1:

Inference 2:

Set the intersection of polygon and PQ to P1, P2 ,...... Pn, where Pi and PI + 1 are adjacent to the two intersections. The condition that the PQ of a line segment is in the polygon is p, q is in the polygon, and for I = 1, 2, ......, The point of N-1, PI, PI + 1 is also within the polygon.

In actual programming, there is no need to calculate all the intersections. First, we should judge whether the edges of a line segment and a polygon are inner. If the inner intersection of a line segment and a certain side of a polygon is outside the polygon; if the line segment and each side of the polygon do not interwork, the intersection of the Line Segment and the polygon must be the end point of the line segment or the vertex of the polygon. You only need to judge whether the point is online.

The algorithm is as follows:

If the ends of the PQ on the wire end are not all inside the polygon.

Then return false;

Point Set pointset Initialization is empty;

For each side of a polygon s

An endpoint of the do if line is on S.

Then adds this endpoint to pointset;

An online segment of an else if s endpoint on PQ

Then adds this endpoint to pointset;

The else if s and the PQ of the line segment have already been handed in.

Then return false;

Sorts points in pointset by X-Y coordinates;

For pointset, pointset [I], pointset [I + 1]

Do if pointset [I], the midpoint of pointset [I + 1] is not in the Polygon

Then return false;

Return true;

In this process, the number of vertices must be much smaller than the number of vertices in the polygon. Therefore, it is a constant-level complexity and can be ignored. Therefore, the time complexity of the algorithm is O (n ).

　　Determines whether the line is within the polygon.:

You only need to determine whether each line segment of the broken line is within the polygon. If there are m line segments and N vertices in a polygon, the time complexity of this algorithm is O (M * n ).

　　Determines whether a polygon is inside a polygon.:

You only need to determine whether each side of a polygon is within the polygon. Determine whether a polygon with M vertices has an O (M * n) Complexity in a polygon with n vertices ).

　　Determines whether the rectangle is inside the polygon.:

Convert the rectangle into a polygon and then determine whether it is within the polygon.

　　Determines whether the circle is inside the polygon.:

You only need to calculate the shortest distance between the circle center and each side of the polygon. If the distance is greater than or equal to the circle radius, the circle is within the polygon. The algorithm for calculating the shortest distance from the center of a polygon to each side of a polygon is described later.

　　Determines whether the vertex is in the circle.:

Calculates the distance from the center of the circle to the point. If the distance is smaller than or equal to the radius, the point is within the circle.

　　Determines whether a line segment, line segment, rectangle, or polygon is in the circle.:

Because the circle is a convex set, you only need to determine whether each vertex is in the circle.

　　Judge whether the circle is inside the circle:

Set the circle to O1, O2, and radius to R1 and R2 respectively. Determine whether O2 is within O1. First, compare the sizes of R1 and R2. If R1 <R2, O2 cannot be within O1; otherwise, if the distance between the two centers is greater than R1-R2, O2 is not within O1; otherwise, O2 is in O1.

　　Calculate the closest point to a line segment:

If the line segment is parallel to the X axis (Y axis), the point is used as the vertical line of the line where the line segment is located. It is easy to obtain the vertical foot, and then calculate the vertical foot, if the vertical foot is on the online segment, the vertical foot is returned; otherwise, the endpoint close to the vertical foot is returned. If the line segment is not on the X axis and is not on the Y axis, the slope exists and is not 0. Set the two ends of a line segment to pt1 and pt2, and the slope is k = (pt2.y-pt1. y)/(pt2.x-pt1.x). The linear equation is: y = K * (X-pt1.x) + pt1.y. The vertical line slope is-1/K, and the vertical line equation is: Y = (-1/K) * (X-point. X) + point. Y.

Simultaneous two-line equations: x = (k ^ 2 * pt1.x + K * (point. y-pt1.y) + point. x)/(K ^ 2 + 1), Y = K * (X-pt1.x) + pt1.y; then, judge whether the DRDs row is in the online segment. If the DRDs row is in the online segment, the DRDs row is returned; if not, calculate the distance between the two ends and the vertical foot, and select the endpoint closer to the vertical foot to return.

　　Calculate the closest point of a point to a line, rectangle, or polygon.:

You only need to calculate the closest point to each line segment, record the closest distance, and obtain the least recent point.

　　Calculate the closest distance from a point to a circle and the coordinates of the intersection.:

If the point is in the center of the center, undefined is returned because the distance from the center to any point in the circle is equal.

The connection point P and the center O. If the Po is parallel to the X axis, the abscissa of the nearest point is calculated based on P on the left or right of O as centerpoint. X-radius or centerpoint. x + radius. If the Po is parallel to the Y axis, the ordinate of the nearest point is calculated based on the P above or below the O, which is centerpoint. Y-+ radius or centerpoint. Y-radius. If the Po is not on the X and Y axes, the PO slope exists and is not 0, and the Po slope of the straight line is k = (P. y-o. y)/(P. x-o. x ). The linear Po equation is Y = K * (X-p. x) + P. Y. Let's set the square process to: (x-o. x) ^ 2 + (y-o. y) ^ 2 = R ^ 2. The intersection of the Po and the circle in a straight line can be solved by the two equations, and the intersection point closest to the point P can be obtained.

　　Calculates the intersection of two collocated line segments.:

There are several situations in which the positional relationship between two collocated line segments is shown. There is no intersection between the two line segments in figure (A), the two line segments in figure (B) and (d) have infinite focal points, and the two line segments in figure (c) have one intersection point. Set line1 to a longer line segment and line2 to a shorter one. If line1 contains two endpoints of line2, it is the case of figure (d, the two line segments have infinite intersections. If line1 only contains one endpoint of line2, then if an endpoint of line1 is equal to the endpoint of line2 contained by line1, it is the situation in figure (c, at this time, the two line segments have only one intersection point. Otherwise, it is in the case of figure (B). The two lines also have infinite intersections. If line1 does not contain any endpoints of line2, it is in figure () in this case, there is no intersection between the two line segments.

 

　　Calculates the intersection of a line or line and a line.:

Set a line segment to L0 = P1P2 and another line segment or line to L1 = q1q2. The intersection of l0 and L1.

1. first, determine whether l0 and L1 are intersecting (the method has been discussed previously). If they do not, there is no intersection. Otherwise, it indicates that l0 and L1 must have an intersection, here we will consider both l0 and L1 as a straight line.

2. If the horizontal coordinates of P1 and P2 are the same, that is, l0 is parallel to Y axis.

A) if L1 is also parallel to the Y axis,

I. if the ordinate of P1 is the same as that of Q1, it means that l0 and L1 are collocated. If L1 is a straight line, they have infinite intersections, if L1 is a line segment, you can use the algorithm "Calculate the intersection of two collinearity lines" to calculate their intersection points (this method has been discussed previously );

Ii. Otherwise, it indicates that l0 and L1 are parallel and they have no intersection points;

B) If L1 runs unevenly on the Y axis, the abscissa of the intersection point is P1. The ordinate coordinates of the intersection point can be calculated from the linear equation of L1;

3. if the abscissa of P1 and P2 is different, but Q1 and Q2 are the same, that is, L1 is parallel to the Y axis, the abscissa of the intersection is Q1, the vertical coordinates of the intersection can be calculated from the linear equation of l0;

4. If the ordinate values of P1 and P2 are the same, that is, l0 is parallel to the X axis.

A) if L1 is parallel to the X axis,

I. if the abscissa of P1 is the same as that of Q1, it indicates that l0 and L1 are both linear. If L1 is a straight line, they have infinite intersections, if L1 is a line segment, you can use the algorithm "Calculate the intersection of two collinearity lines" to calculate their intersection points (this method has been discussed previously );

Ii. Otherwise, it indicates that l0 and L1 are parallel and they have no intersection points;

B) If L1 is not on the X axis, the vertical coordinate of the intersection point is P1. the abscissa of the intersection point can be calculated from the linear equation of L1;

5. if the vertical coordinates of P1 and P2 are different, but Q1 and Q2 are the same, that is, L1 is parallel to the X axis, the vertical coordinates of the intersection are Q1, the abscissa of the intersection point can be calculated from the linear equation of l0;

6. The rest is that the slope of L1 and l0 both exists and is not 0.

A) Calculate the slope K0 of l0 and the slope k1 of L1;

B) If k1 = k2

I. if Q1 is on l0, it indicates that l0 and L1 are both linear. If L1 is a straight line, there is an infinite point of intersection, if L1 is a line segment, you can use the algorithm "Calculate the intersection of two collinearity lines" to calculate their intersection points (this method has been discussed previously );

II. If Q1 is not on l0, it means that l0 is parallel to L1, and there is no intersection between them.

C) A two-line equations can solve the intersection.

This algorithm is not complex, but it should be discussed clearly in different situations. Especially when the two line segments are collocated, we need to consider them separately. Therefore, we will write the algorithm for finding two collocated lines separately in the previous article. In addition, the vector cross-multiplication is used to determine whether a line segment and a line segment (or line) are intersecting at the beginning. If the result is an intersection, all the lines can be considered as a straight line at the end. It should be noted that we can rewrite the linear or line segment equation to the form of AX + by + c = 0, so that part of the steps in the above process can be merged, shortening the code length, however, because a parameter is required first, this algorithm takes more time.

　　Calculates the intersection of a line or line with a line, rectangle, or polygon.:

Calculate the intersection point with each edge separately.

　　Intersection of a line or line and a circle:

Set the center to O, the circle radius to R, and the two points on the straight line (OR Line Segment) L are P1 and P2.

1. If l is a line segment and P1 and P2 are included in the circle o, there is no intersection. Otherwise, proceed to the next step.

2. If l is parallel to the Y axis,

A) calculate the distance from the center of the circle to the distance dis of L;

B) If dis> r, there is no intersection between l and circle;

C) the coordinate of the two intersections can be obtained by using the alignment theorem, but attention should be paid to the tangent of the L and the circle.

3. If l is parallel to the X axis, it is similar to L's parallel to the Y axis;

4. If l is neither parallel to the X axis nor parallel to the Y axis, the slope k of l can be obtained, and then the point-oblique equation of l can be listed. The intersection of L and the circle can be solved by combining with the circle equation;

5. If l is a line segment, for the intersection points obtained in 2, 3, and 4, you must determine whether it belongs to the range of the line segment.

　　Concept of convex hull:

Convex Hull refers to a smallest Convex Polygon that satisfies the vertices in Q or on or within the polygon. The polygon represented by a red line segment is the convex hull of the Point Set q = {P0, P1,... p12.

 

　

　　Optimization of convex hull:

It has been proved that the lower bound of the time complexity of the convex packet algorithm is O (n * logn). However, when the number of vertices H in the convex packet is also taken into account, the pruning search algorithms of krikpatrick and Seidel can reach O (N * logh), and achieve the optimum in a progressive sense. The most common convex hull algorithm is Graham scan and Jarvis step. This article only briefly introduces the Graham scanning method. For more information about its correctness and the process of Jarvis step-by-step method, see Introduction to algorithms.

For a set of Q with three or more vertices, the process of Graham scanning is as follows:

Make P0 the smallest point in the order of Y-X coordinates in Q

Set <P1, P2 ,... PM> This is the point set that is obtained by sorting the rest of the points alphabetically by the Polar Angle centered on P0 (if multiple points have the same polar angle, all except the Farthest Points from P0 are removed

Press P0 into Stack s

Press P1 into Stack s

Press P2 into Stack s

For I need 3 to m

Do While consists of the next element of S's stack top element, S's stack top element, and PI line segment without turning to the left

Stack s

Press PI into Stack s

Return S;

　

After this process is executed, the elements in stack s from the bottom to the top are the point sequence of the convex packet vertices of Q arranged counterclockwise. Note that we do not need to find the polar angle in the order of any two points when we sort the points alphabetically by the polar angle. This step can be achieved by the aforementioned vector cross product.