DE4416363A1 - Rekonstruktion von Bildern aus MR-Signalen in inhomogenen Magnetfeldern - Google Patents

Description

Inhomogenitäten im Grundmagnetfeld führen bei den ge­ bräuchlichen MR-Bildgebungssequenzen bekanntlich zu Bildver­ zerrungen. Heute eingesetzte Pulssequenzen beruhen im allge­ meinen auf dem sogenannten "SPIN-WARP"-Verfahren, wie es beispielsweise in der US-PS 4 706 025 beschrieben ist. Dabei wird ein Kernresonanzsignal vor dem Auslesen in mindestens einer Richtung phasencodiert und während des Auslesens in einer weiteren Richtung durch einen Auslesegradienten frequenzcodiert. Inhomogenitäten des Grundmagnetfeldes in Phasencodierrichtung sind relativ unkritisch, da es nur auf Signalunterschiede zwischen den einzelnen Phasencodier­ schritten ankommt. In Richtung des Auslesegradienten führt die Überlagerung des Auslesegradienten mit Grundfeldinhomoge­ nitäten aber zu Verzerrungen.

Aus dem Artikel J.Weis, L.Budinsky "Simulation of the Influence of Magnetic Field Inhomogeneity and Distortion Correction in MR Imaging" in Magnetic Resonance Imaging, Vol. 8, pp. 483-489, 1990, ist es bekannt, Bildverzerrungen durch eine Nachverarbeitung eines auf herkömmliche Weise (d. h. mit mindestens zweidimensionaler Fourier-Transformation) ge­ wonnenen Bildes durchzuführen. Die dabei erforderliche In­ formation über die Magnetfeldinhomogenitäten, d. h. den Ver­ lauf des Grundmagnetfeldes, wird aus der Phase von separat aufgenommenen Spin-Echo-Bildern gewonnen. Bei diesem Verfah­ ren muß jedoch ein erhöhtes Untergrundrauschen in Kauf ge­ nommen werden.

Aus dem Artikel von C.-M. Lai "Reconstructing NMR Images under Magnetic Fields with Large Inhomogeneities" in J.Phys.E: Sci. Instrum., Vol. 15, 1982, pp 1093-1100 ist es bekannt, Inhomogenitäten bereits in die Bildrekonstruktion mit einzubeziehen. Diese Arbeit geht von der heute nicht mehr gebräuchlichen Projektionsrekonstruktion aus, wobei der bekannte Projektionsrekonstruktions-Algorithmus durch einen Algorithmus ersetzt wird, bei dem die vorher ermittelte Magnetfeldinhomogenität schon bei der Bildrekonstruktion berücksichtigt wird. Hierbei ist allerdings keine Intensi­ tätskorrektur vorgesehen.

Aufgabe der Erfindung ist es, ein Bildrekonstruktionsver­ fahren anzugeben, das bei Grundmagnetfeldern bekannter Inho­ mogenität weitgehend verzeichnungsfreie Bilder, also Bilder ohne geometrische Verzeichnungen und Intensitätsfehler lie­ fert.

Diese Aufgabe wird erfindungsgemäß durch die Merkmale des An­ spruchs 1 gelöst. Vorteilhafte Ausgestaltungen der Erfindung sind in den Unteransprüchen gegeben. Ein Ausführungsbeispiel der Erfindung wird nachfolgend anhand der Fig. 1 bis 11 beschrieben. Dabei zeigen:

Fig. 1 bis 5 zur Erläuterung der Problemstellung eine her­ kömmliche SPIN-WARP-Pulssequenz;

Fig. 6 die dazugehörige Rohdatenmatrix;

Fig. 7 den Zusammenhang zwischen Ort und Resonanz­ frequenz für einen linearen Magnetfeldgra­ dienten;

Fig. 8 den Zusammenhang zwischen Ort und Resonanz­ frequenz unter Berücksichtigung von Magnet­ feldinhomogenitäten;

Fig. 9 und 10 den Vergleich einer idealisierten SPIN-Ver­ teilung und einer tatsächlich im inhomogenen Feld gemessenen SPIN-Verteilung für zwei unterschiedliche Gradientensteigungen;

Fig. 11 eine Fensterfunktion.

Die herkömmliche SPIN-WARP-Sequenz nach den Fig. 1 bis 5 soll lediglich zur Erläuterung der Problemstellung dienen. Es ist zu betonen, daß das erfindungsgemäße Verfahren für alle Pulssequenzen angewandt werden kann, bei denen Kernresonanz­ signale vor dem Auslesen phasencodiert und während des Ausle­ sens frequenzcodiert werden. Bei dem in den Fig. 1 bis 5 dargestellten Beispiel einer Pulssequenz wird zunächst unter Einwirkung eines Schichtselektionsgradienten G_S ein frequenz­ selektiver Hochfrequenzpuls RF eingestrahlt. Unter der Wir­ kung des Schichtselektionsgradienten G_S werden Kernspins nur in einer Schicht des Untersuchungsobjekts angeregt. Anschlie­ ßend wird die durch den positiven Teil-Puls des Schicht­ selektionsgradienten G_S verursachten Dephasierung durch einen negativen Teilpuls G_S⁻ wieder rückgängig gemacht. Ferner wird ein Phasencodiergradient G_P eingestrahlt. Schließlich wird im Ausführungsbeispiel noch ein negativer Auslesegradient G_R⁻ eingeschaltet.

Während der nachfolgenden Auslesephase wird lediglich ein positiver Auslesegradient G_R⁺ eingeschaltet. Das entstehende Echosignal S wird M-mal abgetastet und die so gewonnen M Meß­ werte in eine Zeile einer Rohdatenmatrix nach Fig. 6 einge­ tragen.

Die dargestellte Pulssequenz wird N-mal mit unterschiedlichen Werten des Phasencodiergradienten G_P wiederholt, so daß man insgesamt eine Meßmatrix mit N Zeilen erhält. Üblicherweise wird dabei der Phasencodiergradient von Pulssequenz zu Puls­ sequenz in gleichen Schritten vom höchsten positiven zum höchsten negativen Wert bzw. umgekehrt fortgeschaltet. Die Rohdatenmatrix RD kann man als Meßraum, bei dem im Ausfüh­ rungsbeispiel vorliegenden zweidimensionalen Fall als Meßdatenebene betrachten. Dieser Meßdatenraum wird in der Kernspintomographie als "K-Raum" bezeichnet.

Die für die Bilderzeugung notwendigen Informationen über die räumliche Herkunft der Signalbeiträge ist in den Phasenfak­ toren codiert, wobei zwischen dem Ortsraum und dem K-Raum mathematisch der Zusammenhang über eine Fourier-Transforma­ tion besteht. Es gilt:

Dabei gelten folgende Definitionen

γ = gyromagnetisches Verhältnis
G_R (t′) = Momentanwert des Auslesegradienten
G_P (t′) = Momentanwert des Phasencodiergradienten
ρ (xy) = Kernspindichte.

Bei der in Fig. 6 dargestellten Rohdatenmatrix RD entspricht jede Zeile einem einzelnen Kernresonanzsignal S. Bei schritt­ weiser Fortschaltung des Phasencodiergradienten G_PC erfolgt die Abtastung im K-Raum in aufeinanderfolgenden Zeilen. Zu Beginn jeder Teilsequenz wird jeweils vor dem ersten Kern­ resonanzsignal S₁ ein Phasencodiergradient G_p eingeschaltet, dessen Gradienten-Amplitude von Teilsequenz zu Teilsequenz kontinuierlich schrittweise ansteigt. Wenn man beispielsweise jedes Kernresonanzsignal mit 256 Meßpunkten abtastet und 256 Phasencodierschritt durchführt, so erhält man eine Rohdaten­ matrix mit 256 Zeilen und 256 Spalten, also 256×256 Meßwer­ ten im K-Raum. Die bei der Pulssequenz nach den Fig. 1 bis 5 erhaltenen analogen Meßsignale werden also auf ein Raster im K-Raum digitalisiert.

Üblicherweise wird nun aus der Rohdatenmatrix durch zwei­ dimensionale Fourier-Transformation eine Bildmatrix gewonnen. Bei einem inhomogenen Grundmagnetfeld treten dabei jedoch Verzeichnungen in Richtung des Auslesegradienten G_R auf. Inhomogenitäten in Phasencodierrichtung, die über die gesamte Meßzeit konstant sind, führen jedoch nicht zu Verzeichnungen. In Phasencodierrichtung kommt es nämlich letztlich mehr auf die Signalunterschiede zwischen aufeinanderfolgenden Phasencodierschritten an. Zeitlich konstante Inhomogenitäten in Phasencodierrichtung wirken sich daher nicht aus.

In Richtung des Auslesegradienten G_R führt eine Magnetfeldin­ homogenität jedoch zu Verzeichnungen. Durch den Auslesegra­ dienten G_R sollte ein dem Idealfall linearer Zusammenhang zwischen dem Ort x in Frequenzcodierrichtung und der zugeord­ neten Resonanzfrequenz f der Kernspins erreicht werden, wie dies in Fig. 7 dargestellt ist. Magnetfeldinhomogenitäten führen jedoch zu Nicht-Linearitäten in diesem Zusammenhang, wie dies in Fig. 8 dargestellt ist. Dies rührt daher, daß dem linearen Auslesegradienten G_R Gradienten überlagert sind, die die Inhomogenität des Grundmagnetfeldes repräsentieren.

Ziel der Erfindung ist es nun, ein Verfahren zu finden, bei dem die dargestellten Verzeichnungen schon bei der Bildre­ konstruktion vermieden werden.

Dazu muß man zunächst den Feldverlauf des Grundmagnetfeldes, d. h. dessen Inhomogenität hinreichend genau kennen. Es muß also zunächst das Grundmagnetfeld vermessen werden. Hierfür gibt es eine Reihe von verschiedenen Methoden, wie sie z. B. in der DE-Al 42 18 902 oder von A.A. Maudsley et al., Siemens Forschungs- und Entwicklungsberichte, Band 8 (1979), Nr. 6, Springer Verlag 1979 "Rapid Measurement of Magnetic Field Distributions Using Nuclear Magnetic Resonance" beschrieben sind. Da für den Erfindungsgegenstand auf herkömmliche Verfahren zur Magnetfeldvermessung zurückgegriffen werden kann, werden diese hier nicht näher erläutert.

Bezüglich der Grundfeldinhomogenität müssen einige Voraus­ setzungen gegeben sein, damit das hier dargestellte Verfahren zum Erfolg führt:

• - Die Summe des linearen Auslesegradienten G_R und des inhomo­ genen Anteils des Grundmagnetfeldes, die im folgenden als effektiver Gradient bezeichnet wird, muß eine eindeutige Orts- und Frequenzzuordnung ermöglichen.
• - Die Inhomogenität muß mit hinreichender Genauigkeit bekannt sein, d. h. in der Praxis sollte es möglich sein, die Feldwerte an jedem Ort, der mit dem zu rekonstruierenden Pixel übereinstimmt, genauso gut zu kennen wie bei der konventionellen Bildgebung mit einem linearen Gradienten.
• - Wenn die Ausmessung des Grundmagnetfeldes nur einmal vor der Messung erfolgt, muß die Inhomogenität weitgehend zeitlich unabhängig sein. Diese Einschränkung kann jedoch entfallen, wenn die Ausmessung des Grundmagnetfeldes in regelmäßigen Abständen auch während der Messung erfolgt.
• - Die minimale Auflösung im Bild ist durch die geringste Steigung des effektiven Gradienten gegeben, d. h. hier liegen Optimierungsmöglichkeiten zwischen Inhomogenität und dem linearen Auslesegradienten vor. Zum Beispiel kann man Uneindeutigkeiten durch Erhöhung des linearen Auslesegra­ dienten vermeiden.

Bei selektiver Schichtanregung führen Grundfeldinhomogeni­ täten auch dazu, daß das angeregte Schichtprofil verzerrt ist. In vielen Fällen ist dies weniger störend als die Ver­ zeichnung in Ausleserichtung. Das Problem kann aber auch da­ durch umgangen werden, daß man auf eine Schichtselektion verzichtet und ein 3-D-Verfahren verwendet, d. h., Signale in zwei Richtungen phasencodiert.

Wie bereits erläutert, verursachen Grundfeldinhomogenitäten in Phasencodierrichtung keine Verzeichnungen, so daß in dieser Richtung in herkömmlicher Weise eine Fourier-Trans­ formation durchgeführt werden kann, wobei sich hier der übliche FFT-Algorithmus (Fast Fourier Transformation) an­ bietet. Dies ist in Fig. 6 durch den Hinweis "FFT" für eine Spalte der Rohdatenmatrix schematisch dargestellt. Bei 3- dimensionaler Bildgebung kann entsprechend in zwei Phasen­ codierrichtungen eine Fourier-Transformation durchgeführt werden.

Eine spezifische Methode zur Vermeidung von Verzeichnungen muß daher nur in Ausleserichtung durchgeführt werden. Es ge­ nügt daher, diese Methode eindimensional zu beschreiben. Unter dem Begriff Verzeichnungen werden im folgenden nicht nur geometrische Verschiebungen, sondern auch Intensitäts­ fehler verstanden.

Vorausgesetzt wird, daß die Abhängigkeit des inhomogenen Magnetfelds B (x) von der Ortskoordinate x in Ausleserichtung bekannt ist. Das Grundmagnetfeld B_gesamt (x) setzt sich zu­ sammen aus dem homogenen Grundmagnetfeld B₀ und einem Anteil der Inhomogenität B (x). Die Inhomogenität kann man auch als Gradient betrachten, der zusammen mit den linearen Auslese­ gradienten g einen effektiven Gradienten g_eff ergibt:

Diese Ortsabhängigkeit des effektiven Gradienten bedeutet insbesondere, daß Auflösung und Bandbreite nicht notwendiger­ weise für das gesamte Bild konstant sind. Die Wirkung von Inhomogenitäten ist in den Fig. 9 und 10 nochmals dargestellt. Dabei ist mit B_S jeweils die Intensitätsverteilung im Bild bei einer eindimensionalen, stufenförmigen Spinverteilung ohne Einwirkung einer Inhomogenität dargestellt. Wenn auf­ grund von Inhomogenitäten der effektive Gradient nicht mehr linear ist, sondern eine quadratische Ortsabhängigkeit auf­ weist, wird das tatsächlich erhaltene Bild B_I je nach Art der Abhängigkeit entweder wie in Fig. 9 nach links (|g_eff| < |g|) oder wie in Fig. 10 nach rechts (|g_eff| < |g|) verschoben. In beiden Fällen treten auch Intensitätsfehler auf.

Entsprechend der Erfindung wird zur Vermeidung von Orts- und Intensitätsfehlern in der Richtung x des Auslesegradienten keine Fourier-Transformation, sondern eine Transformation durchgeführt, die man als generalisierte Fresnel-Transfor­ mation (GFT) bezeichnen kann. Das Signal (x_m) für jedes Bildpixel x_m erhält man nach folgender Formel aus den vorher in Phasencodierrichtung fouriertransformierten Meßdaten für eine Zeile der Rohdatenmatrix:

Dabei ist n die Spaltennummer innerhalb der insgesamt N Spal­ ten der Rohdatenmatrix und k_n der entsprechende k-Faktor, g der Gradient in Ausleserichtung und B (x_m) die am Ort x tat­ sächlich vorhandene Feldstärke des Grundmagnetfeldes. I(x) ist ein Intensitätskorrekturfaktor, mit dem die obengenannten Intensitätsfehler korrigiert werden sollen. Dieser Inten­ sitätskorrekturfaktor I(x) wird wie folgt bestimmt:

Dabei ist es nicht notwendig, das Grundmagnetfeld schon bei seiner Berechnung in kartesischen Koordinaten zu differen­ zieren, sondern es genügt, bei der Rekonstruktion den Differentialquotienten aus den Grundmagnetfeldwerten zweier Pixel zu verwenden.

Im folgenden wird gezeigt, daß man bei Verwendung der oben angegebenen Rekonstruktion ein verzeichnungsfreies Bild er­ hält:
Bekanntlich erhält man bei einer eindimensionalen Spin-Ver­ teilung s(x) im inhomogenen Magnetfeld bei einer Spin-Echo- Sequenz folgendes Meßsignal s(t). Dabei ist mit t_n das Zeitraster der Abtastung und mit g der Gradienten in x-Rich­ tung bezeichnet.

bzw. nach der Digitalisierung auf ein Raster im k-Raum:

Wenn man das Signal aus Gleichung 4 in Gleichung 1 einsetzt, so erhält man:

Um diese Gleichung weiter auszuwerten, wird die Summe durch ein Integral ersetzt und es soll N gegen unendlich gehen. Dem entspricht physikalisch, daß das Signal unendlich lang und beliebig fein aufgenommen wird.

Unter Verwendung der Beziehungen

und unter der Annahme, daß der effektive Gradient eindeutig ist, kann der Ausdruck nach Gleichung 6 vereinfacht werden:

Damit ist nachgewiesen, daß das nach der oben angegebenen Re­ konstruktion gewonnene Signal verzeichnungsfrei sowohl hin­ sichtlich des Orts als auch hinsichtlich der Intensität ist, d. h., daß im Idealfall unendlich langer Meßzeit und beliebig hoher Bandbreite das Bild s(x) mit der originalen SPIN- Dichte s(x) übereinstimmt. Im realen Experiment hat sich gezeigt, daß sich durch das hier vorgestellte Verfahren keine wesentlich anderen Diskretisierungseffekte ergeben, als man sie von der Fourier-Rekonstruktion her gewöhnt ist.

Bei der Bildrekonstruktion im inhomogenen Feld liegt für je­ des Pixel einer Zeile ein anderer effektiver Gradient zu­ grunde.

Zwischen der räumlichen Auflösung Δx und dem effektiven Gradienten g_eff (x) besteht bei gegebenem Ausleseintervall T folgender Zusammenhang:

Nun ist es aber eigentlich nicht notwendig, die hohe Auf­ lösung, die der stellenweise sehr starke effektive Gradient möglicherweise liefert, tatsächlich abzubilden, da durch einen hohen Gradienten auch eine hohe Bandbreite und damit ein Verlust von Signal-Rausch-Verhältnis in Kauf genommen werden muß. Ein Objekt sieht unter einem größeren Gradienten ein breiteres Frequenzband und hat folglich im Zeitbereich ein schmaleres Echo. Daher kann das Zeitfenster bei der Re­ konstruktion verkürzt werden bzw. in der Rekonstruktion die Datenpunkte großer k-Werte nicht mitgenommen werden. Dies kann durch Gewichtung der Meßdaten mit einer Fensterfunktion w(k) erfolgen, die adaptiv für jedes Pixel die Auflösung genau auf das Maß reduziert, welches unter gleichen Bedin­ gungen in einem homogenen Magnetfeld erhalten worden wäre.

In die Rekonstruktionsformel nach Gleichung 1 kann man dafür eine Gewichtungsfunktion w (k_n) einführen.

In Fig. 11 ist die zeitliche Abhängigkeit von Kernresonanz­ signalen s(t) schematisch dargestellt, wobei die Kurve a ein Signal unter einem kleineren Gradienten und Kurve b ein Signal unter einem größeren Gradienten repräsentiert. Die Gewichtungsfunktion w(k_n) könnte beispielsweise als Fenster­ funktion in Form eines Rechtecks F wie in Fig. 11 ausgeführt sein, wobei die Breite des Fensters an die Breite des Echosignals s(t) angepaßt wird. Beispielsweise könnte die Gewichtungsfunktion w(t) folgende Form haben:

wobei x_m = Pixel Nummer m und Δx = x_m-x_m-1.

Die Gewichtungsfunktion muß aber nicht notwendigerweise einen rechteckförmigen Verlauf haben. Durch die Fensterfunktion erzielt man eine Verbesserung des Signal-Rausch-Verhält­ nisses. Ferner werden Artefakte vermieden, die sonst durch die aufgrund des örtlich unterschiedlichen effektiven Gra­ dienten variierende Auflösung erzielt werden.

Mit dem dargestellten Rekonstruktionsverfahren gelingt es, auch in inhomogenen Feldern Bilder zu rekonstruieren, die sowohl örtlich als auch bezüglich der Intensität weitgehend verzeichnungsfrei sind. Aus physikalischen Gründen wird das erzielbare Signal-Rausch-Verhältnis bei inhomogenen Feldern stets etwas schlechter sein als bei homogenen Feldern. Die Verschlechterung des Signal-Rausch-Verhältnisses, bezogen auf Signalwert zu Untergrund, fällt jedoch geringer aus als bei Methoden, die auf der Bild-Nachverarbeitung beruhen.

Das beschriebene Verfahren benötigt zwar eine deutlich höhere Rechenzeit als herkömmliche FFT-Algorithmen. Da die zur Ver­ fügung stehenden Rechenleistungen jedoch ständig steigen, fällt dieser Nachteil zunehmend weniger ins Gewicht.

Claims (3)

1. Verfahren zur Rekonstruktion von Bildern aus MR-Signalen in inhomogenen Magnetfeldern mit folgenden Schritten:

• a) Die Magnetfeldverteilung des Grundmagnetfeldes (B(x)) wird bestimmt.
• b) Es wird eine Anzahl von in Phasencodierrichtung unter­ schiedlich phasencodierten Kernresonanzsignalen (S) gewon­ nen, die durch Einwirkung eines Gradienten (g_R) in Aus­ leserichtung frequenzcodiert werden.
• c) Die Kernresonanzsignale (S) werden abgetastet, auf ein Raster im K-Raum digitalisiert und im K-Raum in eine Roh­ datenmatrix (RD) eingetragen.
• d) In der Rohdatenmatrix wird in Phasencodierrichtung eine Fourier-Transformation durchgeführt.
• e) Für jedes Bildsignal wir ein Signal S(x_m) aus den Signalen S(k_n) im K-Raum wie folgt gewonnen: wobei I ein Intensitätskorrekturfaktor ist, der aufgrund folgender Gleichung bestimmt wird:

2. Verfahren nach Anspruch 1, dadurch ge­ kennzeichnet, daß in der Gleichung nach Anspruch l noch eine Gewichtungsfunktion W(k_n) eingeführt wird, die entsprechend der variablen Steilheit eines effektiven Gradienten bei der Rekonstruktion ein mit zunehmender Steilheit engeres Zeitfenster bei der Datenerfassung bewirkt.

3. Verfahren nach Anspruch 2, dadurch ge­ kennzeichnet, daß die Gewichtungsfunktion (W(k_n) eine Rechteckfunktion ist.