Analyse von 3D-CT-Aufnahmen von Spänen zur Extrahierung der Segmentspanbildungsfrequenz

2.1 Vorverarbeitung

Die Daten einer CT-Aufnahme liegen als dreidimensionaler Datenwürfel vor und lassen sich als Matrix G→mit der Dimension I × J × K und den quasi-kontinuierlichen Intensitätswerten g_i,j,k interpretieren. Dabei sind i, j, und k die Indizes der Matrix. Die Umrechung in die räumlichen Koordinaten erfolgt über

Wobei Δx, Δy und Δz jeweils der Pixelabstand entlang der entsprechenden Koordinatenachse ist.

Um eine Punktwolke aus dem Datenwürfel zu erhalten, muss dieser zunächst einer Binarisierung unterzogen werden. Jeder Punkt muss einer der beiden Klassen „Hintergrund“ (0) oder „Span“ (1) zugewiesen werden. Hierzu wird die K-Means-Methode [7] verwendet. Die Grauwerte der Pixel werden als Punkte im eindimensionalen Raum aufgefasst. Zum Initiieren des K-Means-Algorithmus werden zwei zufällige Klassenschwerpunkte vorgegeben. Die wesentlichen Schritte sind nun das Zuordnen der Punkte zu den beiden Klassen. Ein Punkt wird jeweils der Klasse zugeordnet, deren Klassenschwerpunkt am nächsten liegt. Anschließend werden die Mittelwerte der Punkte einer Klasse berechnet und als neue Klassenschwerpunkte zugewiesen. Diese beiden Schritte werden so lange wiederholt, bis die Klassenschwerpunkte konvergieren. Am Ende kann der Schwellwert mittels dem Durchschnitt der beiden Klassenschwerpunkte berechnet werden.

In dieser Punktwolke kommen teilweise zufällige Punktfehler vor. Um diese Punktfehler herauszufiltern wird das größte zusammenhängende Volumen gesucht. Dazu wird ein „Connected-Component Labeling-Algorithmus“ wie er in [3] beschrieben wird, verwendet. Ausgehend von einem initialen Punkt, welcher in der Klasse „Span“ enthalten ist, werden alle angrenzenden Punkte dem ersten Volumen zugeordnet. Werden keine angrenzenden Punkte gefunden, wird ein neuer noch nicht zugeordneter initialer Punkt gewählt und alle angrenzenden Punkte dem nächsten Volumen zugeordnet. Dieser Prozess wird so lange wiederholt, bis alle Punkte einem Volumen zugeordnet sind. Es wird im folgenden das Volumen als Span definiert, welches die meisten zugeordneten Punkte hat. Die Punkte aller anderen Volumina werden der Klasse „Hintergrund“ zugeordnet.

2.2 Räumliche Transformation

Für die räumliche Transformation wird eine Trajektorie benötigt, welche den Verlauf des Spanes im dreidimensionalen Raum beschreibt. Anschließend sollen einzelne Schnitte des gekrümmten Spanes räumlich zu einem geradlinigen Span transformiert werden. In Abbildung 1 sind die Daten eines typischen Spanes dargestellt. Helixförmige Spanverläufe kommen besonders häufig vor, weswegen die Methodik auf diese Form angepasst wird.

Abb. 1

Beispiel einer binarisierten Punktwolke eines helixförmigen Spanes nach der Entfernung der Punktfehler.

Die Punkte der dreidimensionalen Matrix sind in einem äquidistanten Raster angeordnet und lassen sich in einzelne Tomogramme und damit zweidimensionale Matrizen T→kentlang der z-Achse unterteilen. Die Koordinatenachsen der Matrix werden dabei so definiert, dass die z-Achse in Richtung der Rotationsachse der Helixform zeigt.

Für jedes Tomogramm an den Stellen z_k lässt sich ein geometrischer Schwerpunkt

des Spanes bestimmen. Die Menge T_k enthält dabei alle Punkte der Matrix T→kmit dem Wert 1. Der Operator |·| bildet die Mächtigkeit der Menge. Der Vektor p→=[x y z]enthält die dreidimensionalen Koordinaten eines Punktes. Dieser Schwerpunkt ist ein einfaches Merkmal um den Spanverlauf zu bestimmen.

Da der Verlauf der Schwerpunkte über z im Allgemeinen an einzelnen Stellen Sprünge aufweist, wird der Verlauf der Schwerpunkte mit einer nichtlinearen Filterung geglättet. Hierzu wird sowohl die x-, als auch die y-Komponente über z aus mehreren benachbarten Schwerpunkten durch ein Polynom f(z) = az² + bz + c zweiter Ordnung angenähert. Die Parameter des Polynoms lassen sich durch eine Least-Square-Schätzung [6] berechnen. Als Trajektorienpunkt

wird jeweils der Punkt z_k des approximierenden Polynoms verwendet. Damit ist eine Trajektorie gegeben, die einen glatten Verlauf entlang des Spanes beschreibt.

Für die Berechnung der Segmentlänge ist lediglich die Oberfläche des Spans von Interesse. Damit die Datenmenge reduziert werden kann, wird eine Kantendetektion durchgeführt. Dazu wird jedes Tomogramm T→kmit einem Sobel-Filter gefiltert und über einen Schwellwert erneut binarisiert.

Der Span kann nun transformiert werden, indem er in mehrere Schnitte, orthogonal zur Trajektorie, unterteilt wird, welche anschließend aneinandergereiht einen geradlinigen Span ergeben. Ein Schnitt wird realisiert, indem an jedem Trajektorienpunkt t→⌈k⌉eine Ebene orthogonal zur Trajektorientangente aufgespannt wird. Die Trajektorie liegt nicht in analytischer Form vor, weswegen die Tangente und damit die Normale der Ebene durch

an jedem Trajektorienpunkt z_k geschätzt wird. Für einen Punkt, welcher auf der Seite der Ebene liegt, in die der Normalenvektor zeigt, gilt

Damit können alle Punkte, die zwischen zwei Ebenen liegen, der Menge eines Schnittes zugeordnet werden. Bei helixförmigen Spänen kann es vorkommen, dass in einem Schnitt der beiden Ebenen Punkte aus anderen Spanwindungen mit eingeschlossen werden. Dies kann verhindert werden, indem ein maximaler Radius r_max,z um t→[zk]festgelegt wird, welcher die Menge der Punkte eines Schnittes begrenzt. Dieser Radius wird berechnet, indem er iterativ immer größer gewählt wird, bis die Menge der Punkte des Schnittes größer als Null ist und nicht mehr weiter anwächst. Dadurch wird die kleinstmögliche zusammenhängende Oberfläche gefunden. Die Punkte eines Schnittes an der Stelle z_k erfüllen die Bedingungen

Für jeden Schnitt soll eine räumliche Transformation durchgeführt werden. Der x-y-z-Raum wird in den x'-y'- z'-Raum überführt. Dazu wird jeder Punkt eines Schnittes an der Stelle z_k durch

mithilfe der Rotationsmatrizen R→kund K→kund der Translationen t→[zk]und l→ktransformiert. Die erste Translation um den Trajektorienpunkt t→[zk]bildet die Trajektorie auf die neue z'-Achse ab. Der zweite Translationsvektor

verschiebt jeden Schnitt in z'-Richtung an das Ende des Vorherigen, da sonst alle Schnitte im Ursprung bei z' = 0 liegen würden.

Die Rotationsmatrix

eliminiert die Verdrehung der einzelnen Schnitte, die durch die Helixform des Spanes hervorgerufen wird, und besteht aus dem Normalenvektor und den beiden Basisvektoren e→1[zk]und e→2[zk]der Ebene durch den Trajektorienpunkt t→zk.Für die zwei Basisvektoren verbleibt durch die notwendige Orthogonalitätsbedingung ein Freiheitsgrad. Es wird für k = 0 ein zufälliges Basisvektorenpaar gewählt und für k > 1 sukzessiv eine möglichst geringe Aktualisierung vorgenommen. Das bedeutet, dass die beiden anderen Basisvektoren durch

und

berechnet werden können. In der Rotationsmatrix werden die normierten Basisvektoren verwendet. Das Nachführen der beiden Basisvektoren hat den Vorteil, dass in den transformierten Schnitten keine Sprünge in der Rotation vorkommen. Es kann lediglich vorkommen, dass sie eine Rotation entlang des Spanverlaufes aufweisen.

Um diese Rotation zu eliminieren kann die typische Form eines Schnittes ausgenutzt werden. Die Oberflächenpunkte eines Schnittes beschreiben meist einen elliptischen Verlauf. Von dieser sollen die Eckpunkte eines jeden Schnittes bestimmt werden. Diese lassen sich ermitteln, indem in jedem Schnitt der euklidische Abstand zwischen allen Extrempunkten xk, min',  xk,max',  yk,min'und yk,max',den äußersten Punkten eines Schnittes in x- und y-Richtung, berechnet werden. Es werden die beiden Extrempunkte als Eckpunkte κ→0,kund κ→1,kdefiniert, die den größten Abstand zueinander aufweisen. Um auch hier mögliche Sprünge zwischen zwei Schnitten zu vermeiden wird dem Eckpunkt der Index 0 zugewiesen, welcher den geringeren euklidischen Abstand zu κ→0,k−1aufweist.

Mithilfe dieser Eckpunkte wird jeder Schnitt um den neuen Ursprung rotiert, sodass der Vektor zwischen den beiden Eckpunkten

entlang der y-Achse liegt. Dazu wird die Rotationsmatrix

verwendet. Somit ist der transformierte Span geradlinig und weist keine Drehung mehr auf. In Abbildung 2 ist ein beispielhafter transformierter Span dargestellt.

Abb. 2

Ausschnitt aus dem in Abbildung 1 gezeigten Span nach der räumlichen Transformation.

2.3 Analyse der Ortsfrequenz

Um die Segmentspanbildungsfrequenz zu berechnen, muss der transformierte Span einer Ortsfrequenzanalyse unterzogen werden. Ein Span besitzt in der Regel eine glatte Seite, die an der Werkzeugschneide entstanden ist, und eine Seite, auf der sich die Spansegmente ausbilden. Daher soll die Ortsfrequenzanalyse nur auf einer Seite durchgeführt werden. Hier kann die Beschaffenheit des Spanes ausgenutzt werden, da die Oberfläche von Eckpunkt zu Eckpunkt auf der glatten Seite in der Regel einen kürzeren Verlauf aufweist als auf der Spansegmentenseite. Die Punkte eines Schnittes haben allerdings keine Ordnung. Es ist unklar, ob ein Punkt zu der glatten Seite oder der Spansegmentseite gehört. Dieses Problem kann gelöst werden, indem der Winkel eines Punktes zu dem Durchschnitt μ→kaller Punkte eines Schnittes

berechnet wird. Durch die Kenntnis der Eckpunkte können die beiden Seiten voneinander getrennt werden, indem die Winkel der Eckpunkte als Klassengrenzen verwendet werden. Über die geordneten Punkte können die Pfadlängen beider Seite bestimmt werden. Die längere Seite ist in der Regel die Seite mit den Spansegmenten. Da diese Schätzung jedoch teilweise fehlerhaft ist, wird eine Entscheidung über den gesamten Span getroffen. Es wird die Seite ausgewählt, für welche die Mehrheit der einzelnen Schätzungen spricht.

Für die Ortsfrequenzanalyse wird nun ein Signal benötigt, welches den Oberflächenverlauf entlang einer Geraden in z'-Richtung beschreibt. Die transformierten Punkte im x'-y'-z'-Raum sind allerdings nicht in einem äquidistanten Raster. Daher wird ein schmales Intervall Δy′in y'-Richtung definiert, welches alle Punkte einschließt, die verwendet werden, um das Signal zu erzeugen. Diese werden im Anschluss entsprechend ihrer z'-Komponente sortiert. Es kommt bei wenigen Spänen durch die Form der Schnitte vor, dass einzelne Punkte der glatten Seite im Signal auftauchen. Diese lassen sich herausfiltern, indem nur Punkte dem Signal zugeordnet werden, die eine x'-Komponente größer Null aufweisen. Wenn das Intervall Δy′hinreichend klein gewählt wird, ist der Fehler, welcher durch die nicht identische y'-Komponente entsteht, sehr gering. Um ein Signal mit äquidistanter Abtastung zu erzeugen, wird eine lineare Interpolation durchgeführt.

Um die Segmentlänge zu schätzen, sollen die Segmente in einem bestimmten Bereich gezählt werden und im Anschluss der mittlere Abstand zwischen diesen Segmenten gemessen werden. Die Segmentanzahl eines Bereiches lässt sich bestimmen, indem die Anzahl der Nulldurchgänge N₀ des um den Mittelwert befreiten Teilsignales berechnet wird. Da es an einzelnen Stellen des Signals zu hochfrequenten Störungen kommen kann, wird das Signal an äquidistanten Abtastpunkten linear interpoliert und anschließend mit einer gleitenden Mittelung geglättet. Der mittlere Abstand eines Teilsignals ΔL lässt sich durch die Differenz des ersten und letzten Nulldurchganges schätzen. Die mittlere Segmentlänge kann mit

bestimmt werden. Die Segmentspanbildungsfrequenz lässt sich mit Gleichung (1) berechnen.