3-Punkt-Kreisbogen in Standardform bringen

Hallo da draußen,

heute mal ein längerer Beitrag, der mal wieder aus einem Projekt entstanden ist. Aktuell entwickle ich einige Dinge am Laserteile-Shop – es wird demnächst einen sehr speziellen Ableger geben, viel mehr darf ich aber nicht verraten. Im Zuge der Umsetzung sind einige geometrische Dinge umzuwandeln. Bisher bin ich da relativ verschont geblieben, wie bereits im letzten Mathematik-Blog-Eintrag hat das in der Vergangenheit mein Kollege und guter Freund übernommen, der inzwischen aber in Rente ist und wohl nicht mehr allzu viel machen wird.
Daher musste ich mich mal etwas mit der Materie befassen.

Schluss der langen Worte, worum geht es genau?

Ich habe einen Kreisbogen gegeben, der durch 3 Punkte definiert ist:

Startpunkt: P_{0}= (x_{0}, y_{0})
beliebiger Punkt auf dem Kreisbogen: P_{1}= (x_{1}, y_{1})
Endpunkt: P_{2}= (x_{2}, y_{2})

Das entsprechende Format, was ich lesen muss, stellt einen Kreisbogen so dar. Damit ich das aber in ein Vektor-Format (speziell DXF) umgewandelt werden kann, wird der Kreisbogen mit folgenden Parametern erwartet:

  • Mittelpunkt
  • Radius
  • Anfangswinkel
  • Endwinkel

So aus der Kalten hatte ich erst mal gar keine Idee, wie ich das umsetzen kann. Glücklicherweise gibt es ja Google, dabei bin ich auf folgende nette Seite gestoßen:

https://arndt-bruenner.de/mathe/scripts/kreis3p.htm

Um aus den drei Punkten den eigentlichen Kreis zu ermitteln, wird ein Gleichungssystem aufgestellt:

A - Bx_{0} - Cy_{0} = -(x_{0}^{2} + y_{0}^{2})
A - Bx_{1} - Cy_{1} = -(x_{1}^{2} + y_{1}^{2})
A - Bx_{2} - Cy_{2} = -(x_{2}^{2} + y_{2}^{2})

Genaueres könnt ihr dem entsprechenden Link entnehmen.
Mit Hilfe des Gleichungssystems können dann Mittelpunkt und Radius ermittelt werden.
Soweit so gut. Nun muss das Gleichungssystem aber gelöst werden. Hierfür nehmen wir das Gauß-Verfahren. Ich weiß, dass ich das mal zu Studienzeiten in Turbo Pascal implementiert hatte, da das nun aber schon über 20 Jahre her ist, musste ich mir das mal selbst wieder herleiten.
Wieder einmal hat mir da Google geholfen:

https://www.matopt.de/grundlagen/gauss-algorithmus.html

Blöderweise ist dem Autor im Beispiel ein ziemlich schlimmer Fehler unterlaufen, das Tableau 1 ist definitiv falsch, weil er wohl nicht das 3fache der ersten Zeile genutzt hat, sondern das -3 fache der ersten Zeile, was definitiv falsch ist. Aber zumindest die Erklärung ist so nicht schlecht.

Hier wird genauer das Gauß-Verfahren erklärt, allerdings ist da noch keine Computer-Implementierung hinterlegt. Daher werde ich das mal hier in diesem Blog in java tun. Die entsprechenden Dinge erkläre ich dann als Kommentar. Ich habe mir die ganzen Sachen mal hergeleitet, man kann die Klasse mal als Ansatz nutzen, ggf. müsste man noch etwas mehr Fehlerbehandlung durchführen, wenn das System eben nicht lösbar ist.

package com.cutworks.btleasy.helper;
 
import java.util.Vector;
/**
 * Diese Klasse löst ein Gleichungssystem über das Gauss-Verfahren
 * Bei einem Gleichungssystem mit n Variablen werden dabei n Vektoren (Zeilen) der Größe n+1 erwartet,
 * wobei der n+1-te Wert des Vektors das Resultat (ohne unbekannte entspricht)
 * 
 * Beispiel: Für die Gleichung:
 *  2*x + 3*y + 1*z = 5 wird der Vektor: {2, 3, 1, 5} übergeben, für die nächste Zeile dann der nächste Vector, usw.
 *  
 *  
 * @author carsten
 *
 */
public class GaussSolver {
 
 
	@SuppressWarnings("unchecked")
	public static Vector solve(Vector... input) throws Exception {
		if (input.length == 0 || input.length < input[0].size()-1 || input.length >= input[0].size()) {
			throw new IllegalArgumentException("Zur Lösung eines Gleichungssystems mit n unbekannten müssen auch n Vektoren übergeben werden!");
		}
		int size = input.length;
		//der erste Vektor bleibt erst einmal stehen
		//start-Wert der Normalisierung
		//wir fangen mit 1 an und zählen dann später weiter hoch!!
		for (int i= 0; i< size; i++) {
			Vector previousLine = input[i]; //Vektor des entsprechenden Schritts holen = 0. Schritt =&gt; 0. Vektor
			Double prevLineEl = previousLine.get(i); //entsprechenden Spaltenwert holen = 0. Spalte bei 0. Schritt
 
			//nun werden alle restlichen Vektoren durchgegangen
			for (int j = i+1; j<size; j++) {
				Vector line = input[j]; //entsprechenden Vektor holen
				Double elLine = line.elementAt(i); //Wert der Spalte im Vektor des zugehörigen Schrittes holen
				//Ermittlung des Faktors, mit dem die erste Zeile multipliziert wird 
				//um dann die entsprechende Zeile in der i. Spalte durch Addition nullen zu können
				double faktor = (-elLine.doubleValue()/prevLineEl.doubleValue()); 
 
 
				for (int k = i; k < line.size(); k++) {
					double element = line.get(k); //nun alle Elemente des entsprechenden Vektors durchgehen
					double result = faktor*previousLine.get(k)+element; //Faktor*Wert der i.Zeile + Wert der aktuellen Zeile
					line.set(k, result); //und im Vektor aktualisieren!!
				}
			}
		}
		//nun sollte die "Matrix" normalisiert sein 
		//jetzt gilt es noch die entsprechende Werte zu berechnen
		//dieser Vector wird dann später zurückgegeben
		Vector resultVector = new Vector<>();
		//erst mal initialisieren wir ihn mit 0-Werten in der Größe der vorhandenen Variablen
		for (int i=0; i< size; i++) { 
                        resultVector.add(0, 0.0); 
                } 
                //nun gehen wir die Matrix von hinten durch 
                for (int i=size-1; i>=0; i--) {
			Vector vect = input[i];//ermitteln des zugehörigen Vektors
			double summe = 0.0; //Zwischenelement summe initialisieren
			//nun iterieren wir über bereits alle schon ermittelten Resultate und multiplizieren diese 
			//mit den zugehörigen Werten in der Matrix
			//beim ersten Durchlauf wird hier nix iteriert, da wir ja noch kein richtiges Ergebnis haben
			for (int j=i+1; j<size;j++) {
				summe += vect.get(j)*resultVector.get(j);
			}
			//Zur Ermittlung des Wertes wird die letzte Spalte (absolutwert ohne Variable) um die Summe der
			//bereits ermittelten reduziert und durch den entsprechenden Faktor in der Matrix dividiert
			if (vect.get(i) == 0) {
				throw new Exception("Das Gleichungssystem ist nicht lösbar!!");
			}
 
			resultVector.set(i, (vect.get(size)-summe)/vect.get(i));
		}
		return resultVector;
	}
}

Sicher könnte man die Eingabe-Parameter auch als double-Array / zweidimensionales double-Array machen, aber mir erscheint die Lösung mit den Vektoren einfacher, speziell die Befüllung eines / mehrerer Vektoren ist bei weitem einfacher, als programmiertechnisch ein zweidimensionales Array zu erzeugen.

Mit Hilfe des Gauß-Algorithmus haben wir schon die halbe Miete. Wenn man folgende Vektoren übergibt:

\{1, -x_{0}, -y_{0}, -(x_{0}^{2}+y_{0}^{2})\} \{1, -x_{1}, -y_{1}, -(x_{1}^{2}+y_{1}^{2})\} \{1, -x_{2}, -y_{2}, -(x_{2}^{2}+y_{2}^{2})\}

erhält man dann die Werte A, B und C.
Folgt man dann dem Schema auf der Seite von Arndt Brünner kann man dann x_{m} und y_{m} mit folgenden Formeln ermitteln:

 x_{m} = \frac{B}{2}
y_{m} = \frac{C}{2}

und last but not least auch den Radius:
 r = \sqrt{x_{m}^{2} + y_{m}^{2} - A}

Nun haben wir ja fast alles. Bis auf die Winkel….

Gut, dass ich mich für den letzten Mathematik-Beitrag mal wieder etwas mit den Winkelfunktionen beschäftigt hatte. Mit Hilfe dieser kann man dann die Winkel relativ einfach ausrechnen.

Für den Anfangswinkel nutzen wir einfach den Sinus-Satz. Wir erinnern uns:

\sin \alpha = \frac{Gegenkathete}{Hypothenuse}

Die Hypothenuse ist in unserem Fall dann der Radius, die Gegenkathete der Abstand zwischen x_{0} und x_{m}.
Hier ergibt sich dann also:

\sin \alpha = \frac{x_{0}-x_{m}}{r}

Um auf den Anfangs-Winkel zu kommen nutzt man die Umkehrfunktion des Sinus, also den arcussin (asin)

\alpha_{a} = \arcsin(\frac{x_{0}-x_{m}}{r})

Der Endwinkel kann man dann analog ermitteln:

\alpha_{e} = \arcsin(\frac{x_{2}-x_{m}}{r})

Hier ist darauf zu achten, wenn man die Winkel wirklich für eine DXF-Generierung einsetzen möchte, dass an dieser Stelle komischerweise Anfangs- und Endwinkel vertauscht sind.  Erzeugt man die, so wie man sich das eigentlich vorstellen würde, wird der Bogen lustigerweise umgedreht. Wenn man das weiß, ist das aber kein Problem.

Mit Hilfe dieser Sachen konnte ich dann in meinem Projekt aus einem XML-Dokument Polygonzüge auslesen, mit Hilfe der Formeln umwandeln und mit ein paar Freemarker-Templates für ein DXF-Dokument und die entsprechenden DXF-Entities ein im CAD und im Shop lesbares DXF erzeugen.

Ich hätte nie gedacht, mal irgendwann was mathematisches aus meiner Studienzeit wieder auskramen zu müssen. Aber ist schon interessant, wie  das Leben so spielt.

Update 14.05.:

So einfach ist dann das mit den Start- und Endwinkeln wohl doch nicht. Das Problem ist, dass die Vertauschung von Start- und Endwinkeln nur in bestimmten Situationen anwendbar ist, außerdem ist in bestimmten Fällen auch die Subtraktion der Winkel von 360° notwendig.

Durch mehrfaches Trial & Error konnte ich folgende Gesetzmäßigkeiten erkennen:

Vertauscht werden Anfangs- und Endwinkel bei folgenden Situationen:

x_{m} < x_{1}, y_{m} >= y_{1}, x_{2} - x_{0} >= 0, y_{2} - y_{0} >= 0
x_{m} >= x_{1}, y_{m} < y_{1}, x_{2} - x_{0} < 0, y_{2} - y_{0} < 0
x_{m} >= x_{1}, y_{m} >= y_{1}, x_{2} - x_{0} >= 0, y_{2} - y_{0} < 0
x_{m} < x_{1}, y_{m} < y_{1}, x_{2} - x_{0} < 0, y_{2} - y{0} >= 0
x_{m} >= x_{1}, y_{m} < y_{1}, x_{2} - x_{0} < 0, y_{2} - y{0} >= 0
x_{m} < x_{1}, y_{m} >= y_{1}, x_{2} - x_{0} < 0, y_{2} - y{0} >= 0

In diesem Fall werden Anfangs- und Endwinkel getauscht. Zumindest in meinen Beispielen hat das soweit funktioniert.

Abschließend gibt es noch Fälle, in denen der Winkel von 360° subtrahiert wird:

Der Anfangswinkel wird von 360° abgezogen, wenn die Tauschung von Anfangs- und Endwinkel vorgenommen wurde und :
 y_{m} > y_{0}
oder wenn keine Tauschung vorgenommen wurde:
 y_{m} > y_{2}

Der Endwinkel wird von 360° abgezogen, wenn die Tauschung vorgenommen wurde und:
 y_{m} > y_{2}
oder wenn keine Tauschung vorgenommen wurde:
 y_{m} > y_{0}

Nach Einbau dieser Sonderlocken konnte ich dann aus vorgegebenen Konturen auch eine passende DXF erzeugen.

3-Punkt-Kreisbogen in Standardform bringen

Ein Kommentar zu „3-Punkt-Kreisbogen in Standardform bringen

Schreibe einen Kommentar

Deine E-Mail-Adresse wird nicht veröffentlicht. Erforderliche Felder sind mit * markiert.