diff options
| -rw-r--r-- | geometry/formulars.cpp | 75 | ||||
| -rw-r--r-- | tcr.pdf | bin | 256865 -> 256525 bytes |
2 files changed, 41 insertions, 34 deletions
diff --git a/geometry/formulars.cpp b/geometry/formulars.cpp index 2e7b023..e77a0b9 100644 --- a/geometry/formulars.cpp +++ b/geometry/formulars.cpp @@ -1,8 +1,8 @@ -// Komplexe Zahlen als Darstellung für Punkte. -// Wenn immer möglich complex<int> verwenden. Achtung: Funktionen wie abs() geben dann int zurück. +// Komplexe Zahlen als Darstellung für Punkte. Wenn immer möglich +// complex<int> verwenden. Funktionen wie abs() geben dann int zurück. typedef pt complex<double>; -// Winkel zwischen Punkt und x-Achse in [0, 2 * PI), Winkel zwischen a und b. +// Winkel zwischen Punkt und x-Achse in [0, 2 * PI), bzw. zwischen a, b. double angle = arg (a), angle_a_b = arg (a - b); // Punkt rotiert um Winkel theta. @@ -12,14 +12,10 @@ pt a_rotated = a * exp (pt (0, theta)); pt centroid = (a + b + c) / 3.0; // Skalarprodukt. -double dot(pt a, pt b) { - return real(conj(a) * b); -} +double dot(pt a, pt b) { return real(conj(a) * b); } // Kreuzprodukt, 0, falls kollinear. -double cross(pt a, pt b) { - return imag(conj(a) * b); -} +double cross(pt a, pt b) { return imag(conj(a) * b); } // Flächeninhalt eines Dreicks bei bekannten Eckpunkten. double areaOfTriangle(pt a, pt b, pt c) { @@ -37,8 +33,8 @@ double areaOfTriangle(double a, double b, double c) { // zweite Zeile testet Ähnlichkeit mit unterschiedlicher Orientierung bool similar (pt a1, pt b1, pt c1, pt a2, pt b2, pt c2) { return ( - (b2 - a2) * (c1 - a1) == (b1 - a1) * (c2 - a2) || - (b2 - a2) * (conj (c1) - conj (a1)) == (conj (b1) - conj (a1)) * (c2 - a2) + (b2-a2) * (c1-a1) == (b1-a1) * (c2-a2) || + (b2-a2) * (conj(c1)-conj(a1)) == (conj(b1)-conj(a1)) * (c2-a2) ); } @@ -46,27 +42,30 @@ bool similar (pt a1, pt b1, pt c1, pt a2, pt b2, pt c2) { // Einschränken der Rückgabe auf [-1,1] ist sicherer gegen Overflows. double orientation(pt a, pt b, pt c) { double orien = cross(b - a, c - a); - if (abs(orien) < EPSILON) return 0; // Might need large EPSILON: ~1e-6 + if (abs(orien) < EPSILON) return 0; // Braucht großes EPSILON: ~1e-6 return orien < 0 ? -1 : 1; } // Test auf Streckenschnitt zwischen a-b und c-d. bool lineSegmentIntersection(pt a, pt b, pt c, pt d) { - if (orientation(a, b, c) == 0 && orientation(a, b, d) == 0) { // Falls kollinear. + if (orientation(a, b, c) == 0 && orientation(a, b, d) == 0) { double dist = abs(a - b); - return (abs(a - c) <= dist && abs(b - c) <= dist) || (abs(a - d) <= dist && abs(b - d) <= dist); + return (abs(a - c) <= dist && abs(b - c) <= dist) || + (abs(a - d) <= dist && abs(b - d) <= dist); } - return orientation(a, b, c) * orientation(a, b, d) <= 0 && orientation(c, d, a) * orientation(c, d, b) <= 0; + return orientation(a, b, c) * orientation(a, b, d) <= 0 && + orientation(c, d, a) * orientation(c, d, b) <= 0; } -// Berechnet die Schnittpunkte der Strecken a-b und c-d. -// Enthält entweder keinen Punkt, den einzigen Schnittpunkt oder die Endpunkte der Schnittstrecke. -// Achtung: operator<, min, max müssen selbst geschrieben werden! +// Berechnet die Schnittpunkte der Strecken a-b und c-d. Enthält entweder +// keinen Punkt, den einzigen Schnittpunkt oder die Endpunkte der +// Schnittstrecke. operator<, min, max müssen noch geschrieben werden! vector<pt> lineSegmentIntersection(pt a, pt b, pt c, pt d) { vector<pt> result; if (orientation(a, b, c) == 0 && orientation(a, b, d) == 0 && orientation(c, d, a) == 0 && orientation(c, d, b) == 0) { - pt minAB = min(a, b), maxAB = max(a, b), minCD = min(c, d), maxCD = max(c, d); + pt minAB = min(a, b), maxAB = max(a, b); + pt minCD = min(c, d), maxCD = max(c, d); if (minAB < minCD && maxAB < minCD) return result; if (minCD < minAB && maxCD < minAB) return result; pt start = max(minAB, minCD), end = min(maxAB, maxCD); @@ -76,8 +75,10 @@ vector<pt> lineSegmentIntersection(pt a, pt b, pt c, pt d) { } double x1 = real(b) - real(a), y1 = imag(b) - imag(a); double x2 = real(d) - real(c), y2 = imag(d) - imag(c); - double u1 = (-y1 * (real(a) - real(c)) + x1 * (imag(a) - imag(c))) / (-x2 * y1 + x1 * y2); - double u2 = (x2 * (imag(a) - imag(c)) - y2 * (real(a) - real(c))) / (-x2 * y1 + x1 * y2); + double u1 = (-y1 * (real(a) - real(c)) + x1 * (imag(a) - imag(c))) / + (-x2 * y1 + x1 * y2); + double u2 = (x2 * (imag(a) - imag(c)) - y2 * (real(a) - real(c))) / + (-x2 * y1 + x1 * y2); if (u1 >= 0 && u1 <= 1 && u2 >= 0 && u2 <= 1) { double x = real(a) + u2 * x1, y = imag(a) + u2 * y1; result.push_back(pt(x, y)); @@ -98,8 +99,10 @@ bool pointOnLine(pt a, pt b, pt p) { // Liegt p auf der Strecke a-b? bool pointOnLineSegment(pt a, pt b, pt p) { if (orientation(a, b, p) != 0) return false; - return real(p) >= min(real(a), real(b)) && real(p) <= max(real(a), real(b)) && - imag(p) >= min(imag(a), imag(b)) && imag(p) <= max(imag(a), imag(b)); + return real(p) >= min(real(a), real(b)) && + real(p) <= max(real(a), real(b)) && + imag(p) >= min(imag(a), imag(b)) && + imag(p) <= max(imag(a), imag(b)); } // Entfernung von Punkt p zur Strecke a-b. @@ -109,7 +112,8 @@ double distToSegment(pt a, pt b, pt p) { double u = ((real(p) - real(a)) * (real(b) - real(a)) + (imag(p) - imag(a)) * (imag(b) - imag(a))) / (segLength * segLength); - pt projection(real(a) + u * (real(b) - real(a)), imag(a) + u * (imag(b) - imag(a))); + pt projection(real(a) + u * (real(b) - real(a)), + imag(a) + u * (imag(b) - imag(a))); double projectionDist = abs(p - projection); if (!pointOnLineSegment(a, b, projection)) projectionDist = 1e30; return min(projectionDist, min(abs(p - a), abs(p - b))); @@ -130,24 +134,27 @@ bool isCoplanar(pt a, pt b, pt c, pt d) { } // Berechnet den Flächeninhalt eines Polygons (nicht selbstschneidend). -double areaOfPolygon(vector<pt> &polygon) { // Jeder Eckpunkt nur einmal im Vektor. +// Punkte gegen den Uhrzeigersinn: positiv, sonst negativ. +double areaOfPolygon(vector<pt> &polygon) { // Jeder Eckpunkt nur einmal. double res = 0; int n = polygon.size(); for (int i = 0; i < n; i++) - res += real(polygon[i]) * imag(polygon[(i + 1) % n]) - real(polygon[(i + 1) % n]) * imag(polygon[i]); - return 0.5 * res; // Positiv, wenn Punkte gegen den Uhrzeigersinn gegeben sind. Sonst negativ. + res += real(polygon[i]) * imag(polygon[(i + 1) % n]) - + real(polygon[(i + 1) % n]) * imag(polygon[i]); + return 0.5 * res; } -// Testet, ob sich zwei Rechtecke (p1, p2) und (p3, p4) schneiden (jeweils gegenüberliegende Ecken). +// Schneiden sich (p1, p2) und (p3, p4) (gegenüberliegende Ecken). bool rectIntersection(pt p1, pt p2, pt p3, pt p4) { - double minx12 = min(real(p1), real(p2)), maxx12 = max(real(p1), real(p2)); - double minx34 = min(real(p3), real(p4)), maxx34 = max(real(p3), real(p4)); - double miny12 = min(imag(p1), imag(p2)), maxy12 = max(imag(p1), imag(p2)); - double miny34 = min(imag(p3), imag(p4)), maxy34 = max(imag(p3), imag(p4)); - return (maxx12 >= minx34) && (maxx34 >= minx12) && (maxy12 >= miny34) && (maxy34 >= miny12); + double minx12=min(real(p1), real(p2)), maxx12=max(real(p1), real(p2)); + double minx34=min(real(p3), real(p4)), maxx34=max(real(p3), real(p4)); + double miny12=min(imag(p1), imag(p2)), maxy12=max(imag(p1), imag(p2)); + double miny34=min(imag(p3), imag(p4)), maxy34=max(imag(p3), imag(p4)); + return (maxx12 >= minx34) && (maxx34 >= minx12) && + (maxy12 >= miny34) && (maxy34 >= miny12); } // Testet, ob ein Punkt im Polygon liegt (beliebige Polygone). -bool pointInPolygon(pt p, vector<pt> &polygon) { // Jeder Eckpunkt nur einmal im Vektor. +bool pointInPolygon(pt p, vector<pt> &polygon) { // Punkte nur einmal. pt rayEnd = p + pt(1, 1000000); int counter = 0, n = polygon.size(); for (int i = 0; i < n; i++) { Binary files differ |