DE102008050049A1

DE102008050049A1 - Verfahren zur Bestimmung einer die Objektgrenze eines Objekts in einem digitalen 3D-Datensatz annähernden 3D-Oberfläche

Info

Publication number: DE102008050049A1
Application number: DE102008050049A
Authority: DE
Inventors: Karl Dr. Barth
Original assignee: Siemens AG
Current assignee: Siemens Healthcare GmbH
Priority date: 2008-10-01
Filing date: 2008-10-01
Publication date: 2010-04-08
Anticipated expiration: 2028-10-02
Also published as: US20100079456A1; DE102008050049B4

Abstract

Ein Verfahren zur Bestimmung einer die Objektgrenze (6) eines Objekts (4) in einem digitalen 3-D-Datensatzes (2) annähernden 3-D-Oberfläche (14) weist folgende Schritte auf: a) Der 3-D-Datensatz (2) wird in zusammenhängende 3-D-Teilbereiche (10) von Voxeln (8) zerlegt, b) in jedem 3-D-Teilbereich (10) wird ein Wahrscheinlichkeitsmaß (21) für das Vorhandensein der Objektgrenze (6) in diesem 3-D-Teilbereich (10) sowie eine potentielle Lage (26) und Orientierung (24) der Objektgrenze (6) ermittelt, c) in jedem 3-D-Teilbereich (10) wird das Wahrscheinlichkeitsmaß (21) anhand eines Kantenkriteriums (22) dahingehend untersucht, ob die Objektgrenze (6) im 3-D-Teilbereich (10) verläuft, d) aus einem Satz von für den 3-D-Bereich (10) gegebenen Flächenstücken (F1-9) wird dasjenige ausgewählt, das der Lage (26) und Orientierung (24) der Objektgrenze (6) am nächsten kommt, e) die 3-D-Oberfläche (14) wird als Vereinigungsmenge der ausgewählten Flächenstücke (F1-9) gebildet.

Description

Die Erfindung betrifft ein Verfahren zur Bestimmung einer die Objektgrenze eines Objekts in einem digitalen 3D-Datensatz annähernden 3D-Oberfläche.
Auf vielen Gebieten der zerstörungsfreien nichtinvasiven Mess- und Prüftechnik, z. B. auch der medizinischen Bildgebung im Rahmen einer Diagnostik, ist es heute üblich, von zu untersuchenden Objekten, z. B. einem Patienten, einen digitalen dreidimensionalen Bilddatensatz zu erzeugen. Ein solcher Bilddatensatz, auch Volumenscan genannt, ist ein geschlossener, also vollständig mit Voxeln ausgefüllter dreidimensionaler Voxelblock. Bei der bildlichen Darstellung eines derartigen Voxelblockes soll jedoch ein Einblick auf Punkte innerhalb des Volumens möglich sein.
Ein bekanntes Verfahren hierzu ist das Volumenrendern (VR, volume rendering). Das Volumenrendern ermöglicht die direkte Visualisierung, auch direct volume rendering genannt, und liefert in der Regel eine transparente bzw. halbtransparente Darstellung des gesamten Voxelblocks mit gleitenden, also in der Regel unscharfen Struktur- bzw. Objektübergängen. Mit anderen Worten wird so eine „Milchglas”-ähnliche Einsicht in die Tiefe des Volumenblocks möglich. Objektgrenzen eines im Volumenblock abgebildeten Objekts werden mit der Methode des Volumenrenderns jedoch nicht klar dargestellt. Diese Aufgabe ist speziell der Zweck des so genannten „Shaded Surface Display (SSD)”-Verfahrens. SSD setzt voraus, dass innerhalb des 3D-Volumens Objektgrenzen bzw. Punkte (Voxels) an der 3D-Oberfläche von Teilobjekten ermittelt und markiert sind und sodann nur diese Voxels dargestellt werden. Die Ermittlung solcher Oberflächenpunkte (Oberflächenvoxels) innerhalb eines Volumenblocks wird gemeinhin als Segmentierung bezeichnet.
In einer ersten bekannten Variante werden hierzu in einem ersten Schritt jeweils zusammengehörige Volumenelemente, also Voxel, welche gemäß ihrem Messwert ein und derselben Objektart (Knochen) zugeordnet werden sollen, identifiziert und in jeweils gleicher Weise markiert. So entstehen markierte Teilvolumina im 3D-Datensatz, welche im Wesentlichen dem Volumen der Objekte entsprechen. In medizinischen tomographischen Bildern kann der Messwert z. B. der knochenspezifische Dichtewert sein, wenn knöcherne Objekte segmentiert werden sollen. In einem zweiten Schritt werden dann sämtliche Randpunkte jedes Teilvolumens bzw. Objekts bestimmt. Solche Randpunkte sind jene markierte Punkte, welche mindestens einen Nachbarpunkt besitzen, der einem anderen Objekt oder der Objektumgebung (sog. Hindergrundvoxel) zugeordnet ist.
In einer alternativen Variante werden alle Randpunkte direkt, also ohne Bestimmung der zum Volumen des Objekts gehörenden Voxel, bestimmt. Das geschieht z. B. durch Anwendung eines Gradientenoperators in Verbindung mit einem Schwellwert, da ein geeigneter Gradientenoperator an den Stellen hohe Signalwerte liefert, an denen ein Übergang vom Objekt zur Umgebung/Hintegrund vorliegt. Das segmentierte Objektvolumen kann dann umgekehrt zur vorigen Variante als Gesamtheit aller von den Randpunkten eingeschlossenen Voxel ermittelt werden.
Mit der reinen Segmentierung ist aber noch kein Zusammenhang zwischen den Voxeln hergestellt, d. h. es ist nicht klar, welche markierten Voxels jeweils zueinander benachbart sind, und damit ist z. B. auch nicht bekannt, wieviele jeweils in sich geschlossene Objekte sich im Voxelblock befinden. Auch ist bei der SSD-Darstellung immer eine gewisse Stufigkeit der Oberfläche gegeben, besonders wenn die Darstellung in reduzierter Auflösung erfolgt, z. B. wenn jeweils kleine Voxelblöcke von z. B. 8 (2×2×2) oder 27 (3×3×3) Voxeln zur Darstellung zusammengefasst werden.
Für viele Anforderungen und Anwendungen ist daher z. B. eine Oberflächendarstellung als Maschen eines Netzes (Meshes), welche sich an die segmentierten Oberflächenvoxel anschmiegen, besser geeignet. Vorzugsweise sind es die Knoten dieser Maschen, welche die lokalen Oberflächenpositionen definieren und die sich aus den Positionen der in der Segmentierung identifizierten Oberflächenvoxel ergeben. Die Netz-Maschen können beliebige Vielecke sein (Dreiecke, Vierecke, ..., Polygone).
Im Falle von Dreiecken als Maschenelementen spricht man bei der Berechnung der Netzdarstellung von „Triangulierung”. In einem letzten Schritt erfolgt die Visualisierung des Oberflächennetzes, z. B. unter Berücksichtigung geeigneter Farbgebung, Beleuchtung und Orientierung der jeweiligen Oberflächenelemente bezüglich Betrachter und Lichtquellen.
Die Oberflächendarstellung bzw. -beschreibung in Netzform (Maschen) hat generell gegenüber VR und SSD den Vorteil, dass die Nachbarschaft ihrer Elemente klar ist und die so dargestellte Oberfläche einen kontinuierlichen Verlauf hat. Damit können die dargestellten Objekte in Verbindung mit geeigneter Kolorierung und Schattierung klar konturiert und mit optimaler Oberflächenqualität abgebildet werden. Dies ist möglich, da für die Flächenbearbeitung bzw. Darstellung im Wesentlichen nur Aufwand anfällt, welcher mit dem Flächenquadrat steigt. Bei einer gegebenen Rechenkapazität kann so relativ viel Aufwand für eine feine Modellierung jedes Flächenelementes betrieben werden. Ist die Oberfläche des Objektes klar definiert bzw. ermittelt, können die Berechnungen, z. B. um den Kontakt zwischen verschiedenen Objekten herzustellen, um in anderen Anwendungen Kollisionen zu verhindern, Objekte zu verändern oder Objekte zusammenzufügen wesentlich einfacher, zuverlässiger und schneller durchgeführt werden. Auch ist der Speicherbedarf innerhalb eines Rechners, z. B. einer medizinischen Workstation für die zweidimensionale Oberfläche bzw. deren Bearbeitung oder Behandlung in 3D-Bildern in der Netzdarstellung deutlich geringer als bei der oben genannten direkten Volumendarstellung wie z. B. beim VR-Volumenrendern oder beim Shaded Surface Display.
Ein Teil der bekannten Segmentierungsverfahren orientiert sich zunächst an einzelnen 2D-Schichten, welche aus dem 3D-Volumendatensatz nacheinander herausgegriffen werden, und bestimmt zunächst in den 2D-Schichten 2D-Umrisse der gesuchten Objekte. Die gefundenen 2D-Konturen werden anschließend über alle betroffenen Schichten des Voxelblocks hinweg zu einer 3D-Kontur zusammengesetzt.
Direkte 3D-Segmentierungsverfahren verhalten sich zwar besser isotrop als die o. g. Schichtverfahren. Bisher bekannte, direkte 3D-Segmentierungsverfahren, welche direkt am Voxelblock in Originalauflösung arbeiten, also der gesamten Datenmenge, weisen hierbei aber das Problem auf, dass der Rechen- bzw. Suchaufwand, welcher während der Segmentierung in allen 3 Raumrichtungen durchgeführt werden muss, erheblich ist.
In vielen Fällen wird daher zunächst eine besser strukturierte Beschreibung des ursprünglichen 3D-Bilddatensatzes generiert, um den Berechnungsaufwand zu verringern. Oft wird hierzu eine hierarchische Strukturierung des Volumens benutzt, z. B. werden – evtl. in mehreren Stufen – sogenannte Octrees oder Octtrees gebildet. Hierbei werden jeweils 8 Voxel, die eine gemeinsame Ecke besitzen, zusammengefasst, so dass sich nach einem Verfahrensdurchlauf das gesamte Volumen der 3D-Bilddaten mit halber Auflösung pro Achse darstellt. Dies führt zu einem Achtel an Speicherbedarf und entsprechend verringertem Rechenaufwand, wenn zunächst ein Segmentierungsverfahren an einem 3D-Bilddatensatz mit gröberer Auflösung vorgenommen wird. Die volle Auflösung ist jedoch wieder herstellbar, wenn die in der groben Auflösung vorgenommene „Grob-„ Segmentierung schlussendlich zur benötigten höheren Auflösung hin verfeinert wird.
Bei vielen der bekannten Methoden werden die Objekte, Teilobjekte oder Oberflächen bzw. die hierzu gehörenden Voxel anhand von Schwellwerten ermittelt. Hierzu wird zunächst ein Schwellwert definiert: Voxel mit Werten, welche genau diesem Schwellwert entsprechen, sollen den Übergang des Objekts zur Umgebung darstellen, sind also Oberflächenpunkte. Voxel mit höheren Werten liegen z. B. innerhalb, solche mit niedrigeren Werten außerhalb des Objektes. Diese ausschließliche Schwellwertabhängigkeit ist ein großer Nachteil dieses Verfahrens speziell bei medizinischen Bilddaten. In vielen Fällen kann nicht von einem konstanten Signalpegel innerhalb des 3D-Volumens und auch nicht innerhalb des Objektes ausgegangen werden. Z. B. ändern sich Voxelwerte bei der Magnetresonanztomographie durch zunehmendem Abstand von der Empfangsspule. Voxel, die real zu ein und demselben Objekt gehören, können also unterschiedliche Werte aufweisen. Des Weiteren sind sämtliche Schwellwertverfahren anfällig für Störungen bzw. Rauschen in den Volumendaten. Störungen und Rauschen können insbesondere bei medizinischen Datensätzen z. B. aufgrund des Röntgenrauschens oder allgemein des Signalsrauschens – bei Magnetresonanztomographie oder Ultraschall – häufig auftreten. Derartiges Rauschen wirkt sich bei einer schwellwertbasierten Definition unmittelbar aus. Störungen zeigen sich z. B. darin, dass die resultierenden Oberflächen lückenhaft sind, weil es zu vielen Unterbrechungen der Nachbarschaft der Oberflächenelemente, also „Löchern” in der Oberfläche kommt.
Aufgabe der Erfindung ist es, ein verbessertes Verfahren zur Bestimmung von 3D-Oberflächen anzugeben, welche die Objektgrenzen von Objekten in einem digitalen 3D-Datensatz annähern.
Die Aufgabe wird gelöst durch ein Verfahren zur Bestimmung einer die Objektgrenze eines Objekts in einem digitalen 3D-Datensatz annähernden 3D-Oberfläche, mit folgenden Schritten: In einem Schritt a) wird der 3D-Datensatz in zusammenhängende 3D-Teilbereiche von Voxeln zerlegt. Jeder 3D-Teilbereich ist also ein lückenlos zusammenhängender Teil-voxelblock. Die Größe des Teilvoxelblockes umfasst hierbei mindestens zwei Voxel in jeder der drei Raumrichtungen und ist dabei wesentlich kleiner – also z. B. ein Hundertstel bis ein Tausendstel – als die Ausdehnung des 3D-Datensatzes (3D-Voxelblocks) in der entsprechenden Raumrichtung.
In einem Schritt b) wird in jedem 3D-Teilbereich ein Wahrscheinlichkeitsmaß für das Vorhandensein der Objektgrenze in diesem 3D-Teilbereich sowie eine potentielle Lage und Orientierung der Objektgrenze ermittelt. Es werden als Wahrscheinlichkeitsmaß zusammengesetzte Werte, z. B. in Form von Matrizen oder vektoriellen Größen aus den Voxelwerten ermittelt, die also ein Maß für die Abschätzung des Vorhandenseins einer Objektgrenze erlauben, z. B. die Stärke der Voxel-Wertänderung. In jedem Fall wird auch vorsorglich die potentielle (falls also eine solche vorhanden ist) Lage und Orientierung der potentiellen Objektgrenze bestimmt.
In einem Schritt c) werden in jedem 3D-Teilbereich die Werte der Voxel in ihrer Gesamtheit in Form des Wahrscheinlichkeitsmaßes anhand eines Kantenkriteriums dahingehend untersucht, ob die Objektgrenze bzw. ein Teil dieser tatsächlich im vorliegenden 3D-Teilbereich verläuft. Mit anderen Worten wird also der 3D-Datensatz nicht Voxel für Voxel, sondern bereichsweise nach der Objektgrenze durchsucht. Dazu bzw. zur Ermittlung des Wahrscheinlichkeitsmaßes und von Lage und Orientierung der Objektgrenze werden numerische Berechnungen angestellt, die im Wesentlichen kontinuierliche Werte liefern für die Eigenschaften der Objektgrenze, z. B. Voxelwerte, Wertestatistik, Oberflächenorientierung, Farbübergang und/oder Stärke der Ausprägung der Objektgrenze.
Im Schritt c) wird also entschieden, ob der gerade untersuchte 3D-Teilbereich tatsächlich einen flächenhaften Abschnitt der Objektgrenze enthält, und es werden diesem 3D-Teilbereich die im Schritt b) ermittelte Lage und Orientierung zugeordnet. Andernfalls wird dieser 3D-Teilbereich nicht weiter betrachtet.
Falls die Objektgrenze tatsächlich in diesem 3D-Teilbereich verläuft, können in einer optionalen Ausgestaltung des Ver fahrens hier nochmals verfeinerte Werte zur Beschreibung und Darstellung der Objektgrenze im 3D-Teilbereich ermittelt werden, z. B. nochmals deren Lage, lokaler Verlauf und Ausrichtung präzisiert werden.
In einem Schritt d) wird in diesem Falle aus einem Satz von für den 3D-Teilbereich gegebenen Flächenstücken dasjenige ausgewählt, das der Lage und/oder Orientierung der Objektgrenze am nächsten kommt. Mit anderen Worten wird also kein beliebiges individuelles Flächenstück in den 3D-Bereich, einen Teil dessen oder seiner Umgebung eingepasst, sondern es steht nur ein begrenzter Satz von „Musterflächen” zu Verfügung, aus denen dann ein passendes Flächenstück ausgewählt wird.
In einem Schritt e) wird dann die 3D-Oberfläche als Vereinigungsmenge aller wie oben erläutert aufgefundenen bzw. ausgewählten Flächenstücke gebildet.
Durch die teilbereichsweise anstelle der voxelweisen Suche erfolgt bereits bei der Detektion der Objektgrenze gegenüber bekannten Verfahren eine Glättung, Filterung bzw. der Ausgleich von lokalen, einzelnen Voxelfehlern. Die Objektgrenze wird insgesamt störsicherer ermittelt. Im erfindungsgemäßen Verfahren wird also nicht punktweise – Voxel für Voxel – nach der Objektgrenze gesucht, sondern in den 3D-Teilbereichen, also in lokalen Umgebungen der eigentlich zu untersuchenden einzelnen Voxel. Die Umgebung hat mit anderen Worten einen Durchmesser von n Punkten mit 2 ≤ n << N, wenn N der Durchmesser des 3D-Datensatzes ist.
Durch die Benutzung eines Satzes gegebener Flächenelemente erfolgt eine weitere Glättung der Objektoberfläche. Der Satz von Flächenelementen kann z. B. so gewählt werden, dass die entstehende Fläche nur natürlich vorkommenden Krümmungen des Objekts folgen kann, Voxelfehlern, Artefakten o. ä. aber nicht folgen kann und diese somit ausglättet.
Durch das erfindungsgemäße Verfahren wird also das bekannte Schwellwertverfahren, welches zur Ermittlung der Objektgrenze dient, verbessert, was vor allem bei verrauschten Eingangsdaten und/oder nicht verlässlichen Signalpegeln einen Gewinn darstellt. Probleme in bekannten Segmentierungsverfahren, welche an dieser Stelle – durch eine falsche Ermittlung der Objektgrenze – entstehen, werden im erfindungsgemäßen Verfahren eliminiert.
Die Schritte a) bis c) des erfindungsgemäßen Verfahrens stellen damit eine rauscharme und zuverlässige Bestimmung des Verlaufs der Objektgrenze im 3D-Datensatz dar. Die Schritte d) und e) bilden ebenfalls ein neues bzw. verbessertes Verfahren zur Bildung der 3D-Oberfläche oder mit anderen Worten zur Belegung der Objektgrenze mit besonderen Flächenstücken. Die zur Verfügung gestellten Flächenstücke, welche hierzu zur Auswahl stehen, besitzen jeweils einen eigenen Normalenvektor. Entsprechend der ermittelten Orientierung des Normalenvektors der Objektgrenze im entsprechenden 3D-Teilbereich wird also aus der beschränkten Menge von Flächenstücken dasjenige ausgewählt, dessen Normalenvektor dem ermittelten Normalenvektor der tatsächlichen Objektgrenze am nächsten kommt. Der Satz von Flächenstücken ist hierbei so gewählt, dass er es erlaubt, eine zu erwartenden Oberfläche geschlossen abzudecken.
Durch das Verfahren ist es möglich, auch mehrere 3D-(Teil-)Oberflächen auf einmal zu bestimmen. So können durch das Verfahren in einem medizinischen Bilddatensatz z. B. mehrere Knochenfragmente und Organe gleichzeitig segmentiert werden.
In einer bevorzugten Ausführungsform des Verfahrens wird in Schritt b) im 3D-Teilbereich als Wahrscheinlichkeitsmaß z. B. ein 3D-Gradient auf Basis der Werte der Voxel gebildet. Mit anderen Worten wird also in der lokalen Umgebung – dem 3D-Bereich – eines eigentlich zu untersuchenden Voxels im Kern des 3D-Bereiches ein dreidimensionaler Gradient bestimmt, welcher daher auch Nachbarschaftspunkte, nämlich alle Voxel im 3D-Bereich, mit berücksichtigt und somit ein mehr tendenzielles, geglättetes bzw. gefiltertes Maß für das Vorhandensein der Objektgrenze liefert. Ein Kantenkriterium wird dann auf den 3D-Gradienten angewendet, womit entschieden wird, ob durch den 3D-Teilbereich tatsächlich eine 3D-Oberfläche (Kante) verläuft oder nicht. Wenn nicht, wird dieser 3D-Teilbereich nicht weiter berücksichtigt.
In einer alternativen Variante des Verfahrens wird in Schritt c) als Kantenkriterium der Betrag des 3D-Gradienten anhand eines Schwellwertes geprüft. Aufgrund der oben beschriebenen Filter- bzw. Glättungseffekte bzw. -mechanismen ist in diesem Fall die Benutzung eines Schwellwertes deutlich unkritischer als bei den bekannten Verfahren der unmittelbar auf die Voxelwerte angewandte Schwellwert und führt zu zuverlässigeren und damit genaueren Aussagen über den Verlauf der Objektgrenze bzw. deren Vorhandensein in einem 3D-Teilbereich.
Eine besonders zuverlässige Glättung bei der Bestimmung des Gradienten erfolgt in einer Variante des Verfahrens, bei welcher in Schritt b) zunächst im 3D-Teilbereich in jeder der drei Raumrichtungen für jede Voxelebene quer zur jeweiligen Raumrichtung ein Mittelwert der Werte der Voxel gebildet wird. Der Gradient wird dann aus den Mittelwerten gebildet. Mit anderen Worten wird also z. B. zur Bestimmung des Gradienten in x-Richtung der n Voxel messende 3D-Teilbereich in n y-z-Voxelebenen aufgeteilt, welche in x-Richtung schichtweise hintereinander liegen. In jeder y-z-Ebene wird einer von n Mittelwerten der Voxelwerte gebildet. Der x-Gradient berechnet sich dann aus den Werten des n-Tupels der sich in x-Richtung erstreckenden n Mittelwerte. Die so ermittelten Gradientenwerte in allen drei Raumrichtungen werden zum 3D-Gradienten in Form eines Vektors zusammengesetzt.
Im Falle der Ermittlung eines entsprechenden Gradientenvektors als 3D-Gradient kann in einer weiteren Variante des Ver fahrens in Schritt c) stufenlos die Orientierung der Objektgrenze als Richtung des 3D-Gradienten ermittelt werden.
In einer alternativen Ausführungsform des Verfahrens kann in Schritt b) auch für den zu untersuchenden 3D-Teilbereich ein Kantenmuster für die Objektgrenze festgelegt werden. Das Kantenmuster entspricht hierbei dem relativen Verlauf der Voxelwerte, welcher für einen Kantenübergang vom Objekt zum Hindergrund erwartet wird. Es wird dann anschließend als Wahrscheinlichkeitsmaß ein Ähnlichkeitsmaß zwischen den – eventuell gemittelten – Werten der Voxel im 3D-Teilbereich und dem Kantenmuster ermittelt. Das Kantenkriterium wird dann in Schritt c) auf das Ähnlichkeitsmaß angewandt.
Das Kantenmuster wird also vorab so definiert, wie der relative Verlauf von Voxelwerten senkrecht zur gesuchten Objektgrenze typischerweise zu erwarten ist. Mit anderen Worten erfolgt also eine geeignete Verrechnung aller Nachbarschaftspunkte eines zentralen Voxels im 3D-Teilbereich, d. h. aller Voxel des 3D-Teilbereiches, mit einem typischen gesuchten Kantenmuster. Hierbei erfolgt eine gute Störunterdrückung, also Rauschunterdrückung für alle 3D-Oberflächenpunkte bzw. 3D-Kantenpunkte, welche die Objektgrenze im 3D-Datensatz beschreiben.
Auch hier kann in einer Variante der Erfindung als Kantenkriterium das Ähnlichkeitsmaß anhand eines Schwellwertes geprüft werden.
In einer bevorzugten Variante des Verfahrens werden in Schritt a) als 3D-Teilbereiche Quader mit Kantenlängen von mindestens 2 Voxeln gewählt. Durch die Wahl von Quadern oder insbesondere Würfeln kann der gesamte 3D-Datensatz leicht in einzelne 3D-Bereiche unterteilt werden, welche den 3D-Datensatz lückenlos ausfüllen.
In einer bevorzugten Ausführungsform der Erfindung überlagern sich dann benachbarte Quader in jeder der drei Raumrichtungen um jeweils mindestens eine Voxelebene. Mit anderen Worten gehört dann die eine Seitenfläche eines Quaders bildende Voxelebene des einen 3D-Bereiches auch gleichzeitig zum benachbarten 3D-Bereich bzw. bildet auch dessen Seitenfläche. Bei einer entsprechend voxelgenauen Anpassung von Flächenelementen in die 3D-Bereiche entstehen so sich lückenlos aneinanderschließende Gesamtflächen, wenn die Schnittlinien der Flächenstücke mit den Quadergrenzen an zwei benachbarten 3D-Bereichen deckungsgleich sind.
Beim erfindungsgemäßen Verfahren kann es dennoch vorkommen, dass Lücken bzw. ein Versatz zwischen Flächenstücken benachbarter 3D-Teilbereiche in der 3D-Oberfläche bestehen. In einer bevorzugten Ausführungsform des Verfahrens werden daher nach Schritt e) nicht aneinandergrenzende Flächenstücke benachbarter 3D-Teilbereiche durch Interpolation aneinander angepasst. Die Lückenschließung erfolgt z. B. am Rand oder an Eckpunktpositionen der 3D-Teilbereiche senkrecht zur Hauptnormalenrichtung der unpassenden Flächenstücke. Die Hauptnormalenrichtung ist die Richtung der größten Komponente des lokalen Richtungs- bzw. Normalenvektors der Flächenelemente.
An Stellen, wo trotz aller Verbesserungen des Verfahrens Diskontinuitäten (Sprünge) verblieben sind, wird auf mindestens einer Seite eines solchen Sprunges ein Flächenstück in seiner Lage und/oder Orientierung modifiziert oder ein anderes Flächenstück ausgewählt. In einer Tabelle ist festgehalten, welche Kombination von zwei Flächenstücken dabei jeweils zulässig ist.
In einer Variante des Verfahrens wird auf mindestens einer Seite des Versatzes und der angrenzenden Flächenstücke ein anderes als das (die) vorhandene(n) Flächenstück(e) ausgewählt. In einer Tabelle ist festgehalten, welche Kombination von zwei Flächenstücken dabei jeweils zulässig ist.
Auftretende Lücken in der 3D-Oberfläche können in einer weiteren Variante der Erfindung außerdem nach Schritt e) durch Einfügen zusätzlicher interpolierter Flächenstücke zwischen umgrenzenden Flächenstücken geschlossen werden.
Für eine weitere Beschreibung der Erfindung wird auf die Ausführungsbeispiele der Zeichnungen verwiesen. Es zeigen jeweils in einer schematischen Prinzipskizze:
1 einen 3D-Datensatz eines Patienten mit Zerlegung in 3D-Bereiche,
2 einen 3D-Bereich aus 1 im Detail,
3 die Berechnung des Gradienten im 3D-Bereich aus 2,
4 den 3D-Bereich aus 2 mit einem Satz von Flächenstücken,
5 die Schließung eines Versatzes oder einer Lücke in der 3D-Oberfläche durch Interpolation.
1 zeigt einen 3D-Datensatz 2, in welchem als Objekt 4 die Leber eines Patienten abgebildet ist. Beim 3D-Datensatz 2 handelt es sich um einen Dicom-Volumenblock, welcher im Rahmen einer medizinischen Bildgebung aufgenommen wurde. Im 3D-Datensatz 2 soll die Oberfläche, also die Objektgrenze 6 des Objektes 4 segmentiert werden. Der 3D-Datensatz 2 besteht aus einer Vielzahl von Voxeln 8.
In einem ersten Schritt wird der 3D-Datensatz 2 in 3D-Teilbereiche 10 zerlegt. Die Gesamtheit aller 3D-Teilbereiche 10 bildet also den gesamten 3D-Datensatz 2. Jeder 3D-Teilbereich 10 enthält hierbei eine Vielzahl von Voxeln 8. Für die weiteren Erklärungen wird dem 3D-Datensatz 2 ein Koordinatensystem 12 mit Achsen in den Raumrichtungen x, y, z zugeordnet.
2 zeigt einen der 3D-Teilbereiche 10 aus 1 im Detail, welcher sowohl in x, y, und z-Richtung jeweils n = 4 Voxel 8 umfasst. Der 3D-Teilbereich 10 ist also ein Würfel aus insgesamt 4 × 4 × 4 = 64 Voxeln 8. Der 3D-Teilbereich 10 besitzt außerdem 8 Eckpunkte A bis H.
Anschließend wird nun sukzessive jeder der 3D-Teilbereiche 10 auf das Vorhandensein der Objektgrenze 6 im jeweiligen 3D-Bereich 10 hin untersucht. 3 zeigt, wie hierzu in einem ersten Teilschritt der 3D-Teilbereich 10 in eine Anzahl von n_x = 4 sich in y-z-Ebene erstreckende Voxelebenen 20a aufgeteilt wird. In jeder Voxelebene 20a von n_x = 1–4 werden nun die jeweiligen 16 Voxel der betreffenden Voxelebene 20a bzw. Schicht zu einer Summe S_X1-4 summiert und dann zu einem Mittelwert Mx_1-4 gemittelt. Aus den jeweiligen Mittelwerten M_x1-4 wird anschließend ein Gradient g_x gebildet, welcher somit eine Aussage über den Verlauf der gemittelten Voxelwerte innerhalb des 3D-Teilbereiches 10 entlang der Raumrichtung x liefert.
Entsprechendes wird gemäß 3 in weiteren Teilschritten für die Raumrichtungen y und z im 3D-Teilbereich 10 durchgeführt. Aus jeweils vier Voxelebenen 20b und 20c in den Raumrichtungen y und z werden entsprechend Summen Sy1–4 und Sz1–4, Mittelwerte My1–4 und Mz1–4 und Gradienten g_y und g_z gebildet. Schließlich wird aus den Gradienten g_x-z als Komponenten der 3D-Gradient
als Vektor gebildet. Der 3D-Gradient g ⇀ stellt ein Wahrscheinlichkeitsmaß 21 für das Vorhandensein der Objektgrenze (6) in 3D-Teilbereichen dar.
Der Betrag des 3D-Gradienten
stellt ein Maß dafür dar, ob im 3D-Bereich 10 ein Teil der Objektgrenze 6 verläuft. In Form eines Kantenkriteriums 22 wird daher der Betrag des 3D-Gradienten
mit einem Schwellwert S1 verglichen.
Für die Bildung der Gradienten g_x-z existieren manigfaltige Möglichkeiten, z. B. wird g_x zu M_x3 + M_x4 – M_x1 – M_x2 als Differenz bestimmt. Alternativ kann der Gradient gy z. B. auch in Form einer Korrelation zu COR(M_y1-4, K) bestimmt werden. K ist hierbei ein Kantenmuster, welches einen zu erwartenden Verlauf von Voxelwerten bzw. deren Änderung wiederspiegelt, wenn im 3D-Teilbereich 10 die Objektgrenze 6 verläuft. Im vorlie genden Fall ist K z. B. (100, 50, –50, –100) entsprechend den vier zu erwartenden (normierten) Mittelwerten My1–4 bei Vorhandensein einer Objektgrenze 6 im 3D-Teilbereich 10. Der Operator COR bildet den Korrelationskoeffizienten. Alternativ kann auch z. B. der Kovarianzoperator COV verwendet werden. In den letzten beiden Fällen bilden die Operatoren ein Ähnlichkeitsmaß 28 zwischen dem Kantenmuster K und den Werten der Voxel 8 im 3D-Bereich 10. Das Ähnlichkeitsmaß 28 ist ein alternatives Wahrscheinlichkeitsmaß 21.
Hier wird dann anschließend das Ähnlichkeitsmaß 28 gegen einen Schwellwert S₂ geprüft, um zu ermitteln, ob ein Teil der Objektgrenze 6 im 3D-Teilbereich 10 verläuft.
Werden die Schwellwerte S1 oder S2 überschritten, ist also im 3D-Teilbereich 10 lokal die Objektgrenze 6 detektiert, ansonsten nicht. Im Falle des 3D-Teilbereichs von 2 ergibt sich eine positive Detektion. Die Vektorrichtung des 3D-Gradienten
gibt dann an, welche Orientierung 24 die Objektgrenze im 3D-Teilbereich 10 besitzt. Die Lage 26 des 3D-Teilbereiches 10 (bzw. dessen Zentrums) im 3D-Datensatz 2 gibt den Ort der lokalen Objektgrenze 6 an.
In allen 3D-Teilbereichen 10, in denen eine Objektgrenze 6 detektiert ist, soll nun die Objektgrenze 6 durch Flächenelemente angenähert werden.
4 zeigt hierzu nochmals den 3D-Teilbereich 10 mit seinen Eckpunkten A–H. Für den 3D-Teilbereich 10 wird nun ein Satz von insgesamt 9 Flächenelementen F_1-9 generiert, von denen jedes einen Normalenvektor
besitzt. Z. B. besitzt das Flächenelement F1 zwischen den Eckpunkten ACGE den normierten Normalenvektor (1, 0, 0). Das Flächenelement F2 zwischen BCGF den normierten Normalenvektor
und das Flächenelement F₃ zwischen ADHE den normierten Normalenvektor
Die drei Flächenelemente F_1-3 werden auch x-dominierte Flächenelemente genannt, da deren Hauptorientie rung in Richtung der x-Achse verläuft. 4 zeigt außerdem entsprechend zu oben die y-dominierten Flächenelemente F_4-6 und die z-dominierten Flächenelemente F_7-9.
Für den Fall, dass nun in einem der 3D-Teilbereiche 10 der Schwellwert S1 oder S2 überschritten ist, wird aus dem Satz der neun Flächenelement F_1-9 dasjenige ausgewählt, dessen Normalenvektor
dem 3D-Gradienten
bezüglich seiner Orientierung 24 am nächsten kommt. Der Ort des Flächenelements ist durch die Lage 26 bestimmt.
Nachdem das Verfahren für sämtliche 3D-Teilbereiche 10 des 3D-Datensatzes 2 durchgeführt wurde, bilden die entsprechend ausgewählten Flächenelemente F_1-9 zusammen eine 3D-Oberfläche 14, welche die Objektgrenze 6 in etwa annähert. In 1 ist symbolisch ein Ausschnitt der 3D-Oberfläche 14 gezeigt.
5a zeigt eine Situation nach Komplettierung der 3D-Oberfläche 14, wobei diese aufgrund von Signalrauschen nicht geschlossen ist, da in den beiden dargestellten 3D-Teilbereichen 10 die gefundenen Flächenelemente F_i und F_j nicht aneinander stoßen. Zum Schließen der 3D-Oberfläche 14 wird daher zwischen den korrespondierenden Punkten A und C und B und D der 3D-Teilbereiche 10 eine Interpolation der Flächenelemente F_i und F_j und damit der 3D-Oberfläche 14 durchgeführt.
5b zeigt einen weiteren Fall, der auftritt, wenn zwischen zwei 3D-Teilbereichen 10 ein dritter (markiert durch den Pfeil 18) liegt, welchem fehlerhafterweise überhaupt kein Flächenelement zugeordnet wurde. Im betreffenden 3D-Teilbereich 10 wird daher nach Abschluss des eigentlichen Segmentierungsverfahrens ein passendes der Flächenelemente F_k ausgewählt, welches die schräg schraffierten bereites ermittelten Flächenelemente F_i und F_j durch Interpolation verbindet, mit anderen Worten die entstandene Lücke schließt.

2: 3D-Datensatz
4: Objekt
6: Objektgrenze
8: Voxel
10: 3D-Bereich
12: Koordinatensystem
14: 3D-Oberfläche
18: Pfeil
20a–c: Voxelebene
22: Kantenkriterium
24: Orientierung
26: Lage
x, y, z: Raumrichtung
n_x,y,z: Anzahl
A–H: Eckpunkt
S_x-z,1-4: Summe
M_x-z,1-4: Mittelwert
g_x-z: Gradient
g ⇀: 3D-Gradient
S1, 2: Schwellwert
F_1-9,i,j,k: Flächenelement
n ⇀_1-9: Normalenvektor
K: Kantenmuster

Claims

Verfahren zur Bestimmung einer die Objektgrenze (6) eines Objekts (4) in einem digitalen 3D-Datensatz (2) annähernden 3D-Oberfläche (14), mit folgenden Schritten: a) der 3D-Datensatz (2) wird in zusammenhängende 3D-Teilbereiche (10) von Voxeln (8) zerlegt, b) in jedem 3D-Teilbereich (10) wird ein Wahrscheinlichkeitsmaß (21) für das Vorhandensein der Objektgrenze (6) in diesem 3D-Teilbereich (10) sowie eine potentielle Lage (26) und Orientierung (24) der Objektgrenze (6) ermittelt, c) in jedem 3D-Teilbereich (10) wird das Wahrscheinlichkeitsmaß (21) anhand eines Kantenkriteriums (22) dahingehend untersucht, ob die Objektgrenze (6) im 3D-Teilbereich (10) verläuft, d) aus einem Satz von für den 3D-Bereich (10) gegebenen Flächenstücken (F_1-9) wird dasjenige ausgewählt, das der Lage (26) und Orientierung (24) der Objektgrenze (6) am nächsten kommt, e) die 3D-Oberfläche (14) wird als Vereinigungsmenge der ausgewählten Flächenstücke (F_1-9) gebildet.
Verfahren nach Anspruch 1, bei dem in Schritt b): – im 3D-Teilbereich (10) als Wahrscheinlichkeitsmaß (21) ein 3D-Gradient
der Werte der Voxel (8) gebildet wird, – das Kantenkriterium (22) auf den 3D-Gradienten
angewendet wird.
Verfahren nach Anspruch 2, bei dem in Schritt c) als Kantenkriterium (22) der Betrag des 3D-Gradienten
anhand eines Schwellwertes (S₁) geprüft wird.
Verfahren nach Anspruch 2 oder 3, bei dem in Schritt b): – im 3D-Bereich (10) in jeder der drei Raumrichtungen (x, y, z) für jede Voxelebene (20a–c) quer zur jeweiligen Raumrichtung (x, y, z) ein Mittelwert (M_x,y,z) der Werte der Voxel (8) gebildet wird, – der 3D-Gradient
aus den Mittelwerten (M_x,y,z) gebildet wird.
Verfahren nach einem der Ansprüche 2 bis 4, bei dem in Schritt b) die Orientierung (24) der Objektgrenze (6) als Richtung des 3D-Gradienten
ermittelt wird.
Verfahren nach einem der vorhergehenden Ansprüche, bei dem in Schritt b): – als Wahrscheinlichkeitsmaß (21) für den 3D-Bereich (10) ein Kantenmuster (K) für die Objektgrenze (6) festgelegt wird, – ein Ähnlichkeitsmaß (28) zwischen den Werten der Voxel (8) im 3D-Bereich (10) und dem Kantenmuster (K) ermittelt wird, – in Schritt c) das Kantenkriterium (22) auf das Ähnlichkeitsmaß (28) angewandt wird.
Verfahren nach Anspruch 6, bei dem in Schritt c) als Kantenkriterium (22) das Ähnlichkeitsmaß (28) anhand eines Schwellwertes (S₂) geprüft wird.
Verfahren nach einem der vorhergehenden Ansprüche, bei dem in Schritt a) als 3D-Bereiche (10) Quader mit Kantenlängen von mindestens zwei Voxeln (8) gewählt werden.
Verfahren nach Anspruch 8, bei dem in Schritt a) sich benachbarte 3D-Bereich (10) in jeder der drei Raumrichtungen (x, y, z) um jeweils mindestens eine Voxelebene (20a–c) überlagern.
Verfahren nach einem der vorhergehenden Ansprüche, bei dem nach Schritt e) nicht aneinander schließende Flächenstücke (F_i,j) benachbarter 3D-Bereiche (10) durch Interpolation aneinander angepasst werden.
Verfahren nach einem der vorhergehenden Ansprüche, bei dem nach Schritt e) eine Lücke (18) in der 3D-Oberfläche (14) durch ein zusätzlich eingesetztes und interpoliertes Flächenstück (F_k) geschlossen wird.
Verfahren nach einem der vorhergehenden Ansprüche, bei dem nach Schritt e) eine Lücke in der 3D-Oberfläche (14) durch Austausch von Flächenstücken (F_i,j) gemäß einer Tabelle zulässiger Nachbarschaften von Flächenstücken (F_i) mit anderen Flächenstücken (F_j) aus der Menge der zulässigen Nachbarn geschlossen wird.