EP2732430B1

EP2732430B1 - Bestimmung der übertragungsfunktion eines signalverarbeitenden systems ohne bekanntes eingangssignal

Info

Publication number: EP2732430B1
Application number: EP12724869.8A
Authority: EP
Inventors: Andreas Thust
Original assignee: Forschungszentrum Juelich GmbH
Current assignee: Forschungszentrum Juelich GmbH
Priority date: 2011-07-14
Filing date: 2012-04-21
Publication date: 2018-09-19
Anticipated expiration: 2032-04-21
Also published as: US20140108476A1; DE102011107371A1; JP2014529386A; WO2013007226A1; JP6154810B2; DE102011107371A9; US9767072B2; CN103930922B; CN103930922A; EP2732430A1

Description

Die Erfindung betrifft Verfahren zur Bestimmung der Übertragungsfunktion eines signalverarbeitenden Systems, die kein bekanntes Eingangssignal voraussetzen.

Stand der Technik

Signalverarbeitungssysteme, insbesondere Verstärker und optische Abbildungssysteme, arbeiten im Allgemeinen nicht linear und fehlerfrei. Das in einem Raum mit beliebigen Koordinaten x dargestellte Ausgangssignal I(x) ist in der Regel die Faltung des Eingangssignals O(x) mit der Übertragungsfunktion T(x) des Systems: $I (x) = O (x) \otimes T (x)$
wobei der Operator ⊗ für die Faltung steht. Die Übertragungsfunktion T(x) ist von größter technischer Bedeutung für die Nutzung des signalverarbeitenden Systems. Ist diese Funktion bekannt, so kann aus I(x) durch Entfaltung auf das Eingangssignal O(x) und somit auf diejenige physikalische Messgröße zurückgeschlossen werden, die in I(x) steckt. Bei Bildgebungssystemen wird T(x) auch Punktantwortfunktion, englisch "point spread function", genannt. Kenntnis der Übertragungsfunktion T(x) verbessert damit das Funktionieren des ganzen Systems. Ein prominentes Beispiel hierfür ist das Hubble-Teleskop, dessen Übertragungsfunktion neben den unvermeidlichen optischen Aberrationen zusätzlich noch die Auswirkungen eines systematischen Linsenfehlers enthielt. Indem diese Übertragungsfunktion erschlossen wurde, konnten nicht nur die zunächst völlig unbrauchbaren mit dem Teleskop aufgenommenen Bilder korrigiert werden, sondern auch Rückschlüsse auf Art und Stärke des Linsenfehlers gezogen werden. Dies ermöglichte eine nachträgliche optische Korrektur am Teleskop selbst. Um an die Übertragungsfunktion T(x) zu gelangen, wird das System in der Regel mit einem bekannten Eingangssignal O(x) beaufschlagt und dieses mit dem Ausgangssignal I(x) verglichen. Bei Teleskopen wird häufig ein entfernter Stern als Objekt abgebildet, der als ideale Punktquelle mathematisch durch eine Delta-Funktion beschreibbar ist. Bei Photoapparaten werden als bekannte Objekte häufig beleuchtete Lochblenden, Schlitzblenden, oder Streifenmuster mit variablem Streifenabstand als Testobjekte herangezogen. Rauschen konnte in dem Maße als Eingangssignal O(x) verwendet werden, wie ein mittlerer Funktionsverlauf bekannt war.
Für die Kameras von Elektronenmikroskopen behilft man sich mit scharfen Kanten als Teststrukturen (R. R. Meyer, A. Kirkland, "The effects of electron and photon scattering on signal and noise transfer properties of scintillators in CCD cameras used for electron detection", Ultramicroscopy 75, 23-33 (1998); R. R. Meyer, A. I. Kirkland, "Characterisation of the Signal and Noise Transfer of CCD Cameras for Electron Detection", Microscopy Research and Technique 49, 269-280 (2000); R. R. Meyer, A. I. Kirkland, R. E. Dunin-Borkowski, J. L. Hutchison, "Experimental characterisation of CCD cameras for HREM at 300 kV", Ultramicroscopy 85, 9-13 (2000)) sowie mit durch Elektronen erzeugtem Rauschen, dessen Ensemblemittelwert bekannt ist (J. M. Zuo, "Electron Detection Characteristics of a Slow-Scan CCD Camera, Imaging Plates and Film, and Electron Image Restoration", Microscopy Research and Technique 49, 245-268 (2000); K. Du, K. von Hochmeister, F. Philipp, "Quantitative comparison of image contrast and pattern between experimental and simulated highresolution transmission electron micrographs", Ultramicroscopy 107, 281-292 (2007)).
Nachteilig versagt diese Methode, wenn kein bekanntes Testobjekt und damit kein bekanntes Eingangssignal O(x) zur Verfügung steht.

Aufgabe und Lösung

Es ist daher die Aufgabe der Erfindung, ein Verfahren zur Bestimmung der Übertragungsfunktion eines signalverarbeitenden Systems zur Verfügung zu stellen, das ohne Kenntnis des Testobjekts, und damit ohne Kenntnis des Eingangssignals O(x), auskommt.
Diese Aufgabe wird erfindungsgemäß gelöst durch Verfahren gemäß Haupt- und Nebenansprüchen. Weitere Ausgestaltungen ergeben sich jeweils aus den darauf rückbezogenen Unteransprüchen.

Gegenstand der Erfindung

Im Rahmen der Erfindung wurde ein Verfahren zur Bestimmung der Übertragungsfunktion eines signalverarbeitenden Systems aus mindestens zwei Repräsentationen I₁(x) und I₂(x) eines Objekts, die das System aus unterschiedlich skalierten vom Objekt herrührenden Eingangssignalen erzeugt hat, oder aus einer Repräsentation I₁(x) eines ersten Objekts und aus einer Repräsentation I₂(x) eines hierzu bis auf einen Skalenfaktor geometrisch identischen Objekts entwickelt.
Die Koordinaten x, in denen die Repräsentationen dargestellt sind, können beispielsweise Orts- und/oder Zeitkoordinaten sein. Es kann aber auch beispielsweise die Energie eine der Koordinaten x sein, etwa bei der Energieverlustspektroskopie (electron energy loss spectroscopy, EELS), bei der ein Energieverlustspektrum auf eine CCD-Kamera geworfen wird. In den meisten Fällen können die Energie oder andere Koordinaten bei der Erzeugung der Repräsentationen I₁(x) und I₂(x) auf Orts- und Zeitkoordinaten zurückgeführt werden. So wird etwa bei EELS die Energieabhängigkeit in der Apparatur in Ortskoordinaten auf dem CCD-Chip der Kamera umgewandelt. Sofern im Folgenden von orts- und/oder zeitabhängigen Signalen die Rede ist, stehen diese Koordinaten somit ohne Einschränkung der Allgemeinheit auch stellvertretend für andere Koordinaten, die sich auf Ort und/oder Zeit zurückführen lassen.
Ein Objekt im Sinne dieser Erfindung vermag durch seine Anwesenheit eine physikalische Messgröße orts- und/oder zeitabhängig zu verändern und hierdurch ein Eingangssignal für das signalverarbeitende System zu erzeugen. Beispielsweise verändert eine Schallquelle die Schallintensität im Raum orts- und zeitabhängig. Ein photographierbares Objekt verändert die Lichtintensität im Raum. Die Schall- bzw. Lichtintensität ist jeweils das Eingangssignal O(x), das für signalverarbeitende Systeme bereitgestellt wird. Das Objekt selbst muss nicht bekannt sein.
Trifft das Eingangssignal auf ein signal verarbeitendes System, das auf diese Art von Signal sensitiv ist, so liefert dieses an seinem Ausgang I(x) gemäß Gl. (1). I(x) ist eine Repräsentation des Objekts, das das Eingangssignal O(x) ausgesendet hat, erhalten unter den Aufnahme- bzw. Abbildungsbedingungen, die in T(x) verkörpert sind.
Aus ein und demselben Objekt können unterschiedlich skalierte Repräsentationen erzeugt werden. Die Skalierung entspricht der räumlichen oder zeitlichen Abtastrate, mit der das Signal erfasst wird. Bei einer Tonaufnahme auf einem analogen Band hängt beispielsweise die zeitliche Abtastrate davon ab, wie schnell das Band läuft. Bei einer digitalen Tonaufnahme wird eine Abtastrate vorgegeben und entscheidet maßgeblich darüber, um wie viel sich die Datei pro Minute Aufnahme vergrößert. An einer Kamera entscheidet die Einstellung des Objektivs auf Tele oder Weitwinkel darüber, welcher räumliche Bildausschnitt auf die feste Anzahl Pixel des Bildsensors projiziert wird, und damit über die räumliche Abtastrate. Gleiches gilt beim Mikroskop, dessen räumliche Abtastrate von der eingestellten Vergrößerungsstufe abhängt.
Zwei unterschiedlich skalierte Repräsentationen mit gleicher Wirkung liegen auch dann vor, wenn mit dem System eine Repräsentation I₁(x) eines ersten Objekts und eine Repräsentation I₂(x) eines hierzu bis auf einen Skalenfaktor geometrisch identischen Objekts angefertigt worden sind. Ein Beispiel hierfür ist das Abphotographieren zweier unterschiedlich vergrößerter und hinreichend scharfer Ausdrucke desselben Motivs mittels eines Photoapparates unter identischen Abbildungsbedingungen.
Erfindungsgemäß werden die beiden Repräsentationen I₁(x) und I₂(x) zunächst in einen Arbeitsraum (mit Koordinaten g) transformiert, in dem sie jeweils als Produkt $I (g) = O (g) T (g)$
der Übertragungsfunktion mit dem unbekannten Eingangssignal darstellbar sind. Dieser Arbeitsraum kann beispielsweise der Frequenzraum sein, und um eine Repräsentation in ihn zu überführen, kann sie insbesondere fouriertransformiert werden. Die Transformation hat zum Ziel, durch die Umwandlung der Faltung in ein Produkt gemäß Gl. (2) das weitere Rechnen mit den Repräsentationen zu vereinfachen.
Die Darstellung der Repräsentation im Frequenzraum bezeichnet man auch als Spektrum. In diesem Raum ist die Übertragungsfunktion T(g) die Entsprechung der Punktantwortfunktion und wird Modulationstransferfunktion genannt.
Es gibt Sonderfälle, in denen die Objekte nicht bis auf einen Skalenfaktor geometrisch identisch sein müssen, sondern es ausreicht, wenn sie geometrisch ähnlich sind. Das bedeutet, dass dann Drehungen und/oder Spiegelungen des einen Objekts gegenüber dem anderen keinen zusätzlichen Unterschied zwischen den Repräsentationen I₁(x) und I₂(x) bewirken.
Beispielsweise kann die Übertragungsfunktion T(x) oder ihr Pendant T(g) im Arbeitsraum eine Symmetrie, insbesondere eine Rotationssymmetrie, aufweisen. Dann können auch Repräsentationen zweier gegeneinander gedrehter oder gespiegelter Eingangssignale O₁ und O₂ im Realraum (Koordinaten x) oder im Arbeitsraum (Koordinaten g) zur Deckung gebracht werden.
In diesem Zusammenhang ist noch ein weiteres vereinfachendes Szenario zu erörtern: Wenn eine ausreichend hohe Anzahl gedrehter und gespiegelter Varianten im Bild vorhanden sind (z.B. Hunderte oder Tausende Bakterien oder Hunderte oder Tausende Nanocluster von Atomen im Elektronenmikroskop), kann jegliche azimuthale Präferenz des Objekts von Natur aus verschwinden. In diesem Falle sind trotz gedrehter und gespiegelter Einzelobjekte (Einzelbakterium, Einzelcluster) innerhalb der Ensembles die jeweiligen Ensemblemittelwerte bis auf einen Vergrößerungsfaktor wieder identisch. In diesem Falle muss dann noch nicht einmal dasselbe Gebiet in den beiden Bildern enthalten sein (identisches Bakterium, bzw. genau dieser selbe Atomcluster), sondern es reicht sogar, wenn ein repräsentativer Ensemblemittelwert aus zwei verschiedenen Aufnahmeorten ermittelt wird.
Letzteres Prinzip macht man sich auch bei der Nutzung von unbekanntem Rauschen als Eingangssignal zunutze, da niemals identisches Rauschen bei zwei hintereinanderfolgenden Aufnahmen erfasst werden kann. Jedoch kann man davon ausgehen, dass die entsprechenden Ensemblemittelwerte der bei beiden Aufnahmen erhaltenen Eingangssignale O(x) identisch sind.
Es werden Ausschnitte der beiden Repräsentationen I₁(x) und I₂(x), die sich auf auf den gleichen Bereich bzw. gleichwertige Bereiche des Objekts bzw. der Objekte beziehen und somit auf das gleiche unbekannte Eingangssignal zurückgehen, ausgewählt und im Arbeitsraum als Funktionen I₁(g) und I₂(g) ausgedrückt. Insbesondere kann eine der Repräsentationen zur Gänze als Ausschnitt verwendet werden und aus der zweiten Repräsentation derjenige Ausschnitt, der den gleichen Objektinhalt wiedergibt wie die erste Repräsentation. Die Ausschnitte können bereits vor der Transformation in den Arbeitsraum ausgewählt werden. Eine Auswahl nach der Transformation in den Arbeitsraum ist in Sonderfällen möglich, aber deutlich schwieriger.
Um die Ausschnitte zu identifizieren, wird ein Ähnlichkeitsmaß, wie etwa die Kreuzkorrelation, optimiert, wobei Position und Ausdehnung des Ausschnitts aus der zweiten Repräsentation die freien Parameter sind. Sind die Koordinaten x im Realraum diskret, können beispielsweise alle möglichen Werte durchprobiert und das Maximum des Ähnlichkeitsmaßes bestimmt werden.
Das Verfahren versagt nicht schlagartig, wenn der Ausschnitt aus der zweiten Repräsentation nicht exakt den gleichen Objektinhalt wiedergibt wie die erste Repräsentation. Vielmehr wird die Bestimmung der Übertragungsfunktion dann progressiv ungenauer. Der Auswahlfehler bewirkt, dass sich die durch beide Repräsentationen wiedergegebenen Objektinhalte in einem Teilgebiet nicht überlappen. Der entsprechende relative Fehler ergibt sich aus dem Verhältnis der nicht überlappenden Fläche zur Gesamtfläche der ersten Repräsentation. In erster Näherung kann dieser relative Fehler auch im Arbeitsraum, beispielsweise im Frequenzraum, angenommen werden.
Werden beispielsweise in diskrete Pixel unterteilte Bilder als Repräsentationen gewählt, entsteht beim Übereinanderlegen mittels Durchprobieren für jede der beiden Dimensionen ein Fehler von maximal einem halben Pixel. Wird der gleiche Fehler für die Ausdehnung (Vergrößerung) angenommen, ergibt sich ein maximaler Gesamtfehler von einem Pixel je Dimension. Sind die Bilder quadratisch mit einer Kantenlänge von N Pixeln, enstpricht der relative Fehler dem Verhältnis der Pixelanzahl 2*N in den überstehenden, nicht überlappenden Streifen zur Gesamtzahl N*N aller Pixel, also 2*N/(N*N)=2/N. Bei typischen Bildkantenlängen der Größenordnung N=1024 Pixeln liegt der relative Fehler gemäß dieser Abschätzung somit im unteren Promillebereich.
Der Begriff des "gleichen Bereichs", auf den sich die beiden Repräsentationen beziehen, ist nicht rein geometrisch zu sehen. Vielmehr umfasst er auch gleichwertige Bereiche in dem Sinne, dass die Unterschiede zwischen diesen Bereichen keinen signifikanten Unterschied zwischen den Repräsentationen I₁(x) und I₂(x) bewirken. Dies ist insbesondere dann der Fall, wenn das Eingangssignal O(x) ein Ensemblemittelwert über sehr viele statistisch verteilte Einzelbeiträge ist, beispielsweise die Lichtintensität, die von Hunderten einzelner Bakterien emittiert wird. Analog ist der Begriff der "gleichwertigen Bereiche" auszulegen, wenn I₁(x) und I₂(x) nicht Repräsentationen desselben Objekts sind, sondern I₁(x) die Repräsentation eines ersten Objekts und I₂(x) die Repräsentation eines hierzu bis auf einen Skalenfaktor geometrisch identischen Objekts ist.
Es wird, beispielsweise punktweise im Arbeitsraum, der Quotient Q(g)=I₂(g)/I₁(g) aus den beiden Funktionen I₁(g) und I₂(g) gebildet, so dass das unbekannte Eingangssignal unterdrückt wird sowie Zähler und Nenner jeweils die gesuchte Übertragungsfunktion mit unterschiedlich skalierten Argumenten enthalten. Im Idealfall wird das unbekannte Eingangssignal O(g) vollständig eliminiert. Der Quotient Q(g) kann im Arbeitsraum nicht nur punktweise gebildet werden, sondern auch beispielsweise über Gebiete oder entlang von Wegen.
Seien die beiden Repräsentationen beispielsweise Bilder I₁ und I₂ desselben Objekts, die im Frequenzraum als Arbeitsraum die Darstellungen I₁(g) und I₂(g) haben. Ist die Vergrößerung des Bildes I₁ um den Faktor γ größer als die des Bildes I₂, zeigen aber beide Bilder den gleichen Ausschnitt O(g) aus dem Objekt, dann macht sich der Unterschied in der Vergrößerung lediglich in unterschiedlichen Abtastraten der Übertragungsfunktion T(g) bemerkbar. Für den Quotienten Q(g)=I₂(g)/I₁(g) gilt dann: $Q (g) = \frac{I_{2} (g)}{I_{1} (g)} = \frac{O (g) T (γ g)}{O (g) T (g)} = \frac{T (γ g)}{T (g)}$
Schließlich wird erfindungsgemäß aus dem Verlauf von Q(g) der gesuchte Verlauf der Übertragungsfunktion T(g) ausgewertet.
Es wurde erkannt, dass auf diese Weise die Übertragungsfunktion T(g) bestimmt werden kann, ohne dass das nominelle Eingangssignal O(g) bekannt ist. Selbst völlig zusammenhangloses Rauschen, von dem kein mittlerer Funktionsverlauf bekannt ist, kann als Eingangssignal verwendet werden. Nach dem Stand der Technik war genau diese Kenntnis erforderlich, um durch Umstellen von Gl. (2) nach T(g)=I(g)/O(g) an T(g) zu kommen; kennt man nur I(g), ist ohne jegliches zusätzliche Wissen weder eine Bestimmung von O(g) noch eine Bestimmung von T(g) möglich ("blind deconvolution"-Problem). Bei Teleskopen beispielsweise verschaffte man sich diese Kenntnis, indem ein weit entfernter Stern, der näherungsweise als Punktquelle angesehen werden kann, ins Visier genommen wurde. Bei Mikroskopen wurden Testproben mit bekannten Strukturen abgebildet. Indem nun gegenüber der "blind deconvolution"-Situation als einziges zusätzliches Wissen die unterschiedliche Skalierung von I₁ und I₂ eingesetzt wird und dafür eine Kenntnis von O(g) nicht mehr benötigt wird, wird die Bestimmung von T(g) erstmals für signalverarbeitende Systeme möglich, für die es keine geeigneten Testobjekte gibt.
Beispielsweise ist es bei einem Elektronenmikroskop, das bis in atomare Dimensionen vergrößert, nicht möglich, wohldefinierte Objekte abzubilden, da Testobjekte in atomarer Dimension nicht auf definierte Art und Weise zu präparieren sind. Man behilft sich nach dem Stand der Technik, indem man unter Umgehung der Elektronenoptik eine scharfe Kante unter Umgehung der Elektronenoptik direkt auf den Detektor legt, oder alternativ einen scharfen Schatten auf den Detektor projiziert. Da das Objektspektrum O(g) einer ideal scharfen Intensitätsstufe bekannt ist, kann man aus der gemessenen Intensitätsverteilung I(g) die Übertragungsfunktion T(g) der CCD Kamera bestimmen. Zum anderen wird die Rauschmethode (engl. noise method) angewandt, bei der man die CCD Kamera gleichmäßig mit Elektronen ausleuchtet. Aufgrund des statistischen Charakters der Anzahl der pro Pixel einfallenden Elektronen stellt das Eingangssignal O(g) weißes Rauschen dar. Durch Vergleich dieses bekannten Eingangssignals mit dem von der Kamera aufgenommenen Ausgangssignal I(g) lässt sich wiederum gemäß T(g)=I(g)/O(g) die Übertragungsfunktion T(g) bestimmen. Beide Methoden führen in der Elektronenmikroskopie zu gewissen Schwierigkeiten. Bei der Kantenmethode muss meist die evakuierte Mikroskopsäule geöffnet werden, um einen Gegenstand mit scharfer Kante auf die Kamera zu legen. Zudem muss die Ausrichtung der scharfen Kante entlang der Pixelreihen wohldefiniert sein, was bei typischerweise 2048 Pixeln der Größe 15 Mikrometer entlang einer Dimension kaum machbar ist. Zusätzlich kann es zu Problemen kommen, da die Kante eine gewisse physikalische Dicke hat und somit auch Elektronen von der Kante auf den Detektor gestreut werden können, oder es können auch Röntgenquanten an der Kante entstehen, die vom Detektor unerwünschterweise registriert werden. Die Benutzung der Rauschmethode ist ebenso problematisch, da zwar der relative funktionale Verlauf der Übertragungsfunktion T(g) bestimmt werden kann, aber eine Absolutnormierung von T(g) nicht möglich ist.
Infolgedessen können in der Elektronenmikroskopie die mit den beiden Methoden bestimmten Übertragungsfunktionen ein und derselben CCD Kamera erheblich voneinander abweichen. Aufgrund dieser Unsicherheiten bei der Bestimmung der Übertragungsfunktion wird eine exakt quantitative Bildauswertung, die häufig einen Vergleich mit numerisch berechneten Bildsimulationen verwendet, erschwert oder sogar unmöglich gemacht. Die nun erfindungsgemäß geschaffene Möglichkeit, eine dritte alternative Messmethode zur Bestimmung der Übertragungsfunktion T(g) zur Verfügung zu haben, erlaubt es, die Diskrepanzen bei der Bestimmung der Übertragungsfunktion in der Elektronenmikroskopie zu klären und die Genauigkeit bei der quantitativen Bildauswertung zu steigern.
Das erfindungsgemäße Verfahren kann nicht nur während des Betriebes des signalverarbeitenden Systems angewendet werden. Vielmehr kann T(g) auch dann noch aus alten mit dem System angefertigten Repräsentationen (etwa Bildern) bestimmt werden, wenn es das System selbst längst nicht mehr gibt, sondern nur noch die Bilder ohne weitere Aufzeichnungen hierzu. Gibt es zwei Bilder, die das gleiche Objekt mit unterschiedlichen Vergrößerungen zeigen, reicht dies zur Bestimmung von T(g) aus, und mit dieser Kenntnis können alle übrigen vorhandenen Bilder nachträglich verbessert werden. Beispielsweise können Aufnahmen von längst beendeten Raumsondenmissionen mit Hilfe des erfindungsgemäßen Verfahrens neu ausgewertet werden, so dass das Jahrzehnte alte Material auch heute noch neue Erkenntnisse liefert und im Extremfall eine neue Mission überflüssig wird.
Nullstellen von I₁(g) im Nenner von Gl. (3), die auf Nullstellen der Übertragungsfunktion T(g) zurückgehen, sind unproblematisch. Sie sind dadurch gekennzeichnet, dass I₁(g) asymptotisch auf Null abfällt. Auf Grund der Retardierung von I₂(g) gegenüber I₁(g) hat I₂(g) die Nullstellen von I₁(g) bei kleineren Werten von g. Gehen also ausgehend von g=0 die Werte von I₁(g) asymptotisch gegen Null, ist dies zuvor bereits im Zähler I₂(g) von Gl. (3) geschehen, so dass der Quotient Q(g) wohldefiniert ist.
Andere Nullstellen von I₁(g), die auf Nullstellen des Eingangssignals O(g) oder auf oszillatorisches Verhalten der Übertragungsfunktion T(g) zurückzuführen sind, lassen sich beherrschen, indem diese Singularitäten aus der Berechnung von T(g) ausgenommen werden. Insbesondere bei einem parametrisierten Ansatz für T(g) wird bereits eine gewisse Glattheit des Verlaufs von T(g) vorausgesetzt, beispielsweise in Form der Stetigkeit oder der maximalen Krümmung. Somit kann mit gutem Grund angenommen werden, dass die für den überwiegenden Teil (≥ 95 %) des Arbeitsraums zutreffende Lösung auch für die Singularitäten zutrifft. Insbesondere können Gebiete, auf denen sehr große Werte des Quotienten Q(g) bzw. abgeleiteter Größen wie z.B. ln[Q(g)] oder D(g) auftreten, aus der Berechnung von T(g) ausgenommen werden. Man bedenke dabei, dass die Singularität nur künstlich durch die Division gemäß Gl. (3) hergestellt wurde, während die Kurven T(γg) und T(g) selbst gar keine Singularität enthalten und durchweg "harmlos" sind und stetig verlaufen.
Im Folgenden wird gezeigt, dass sich aus dem Ausdruck für Q(g) in eindeutiger Weise die gesuchte Übertragungsfunktion T(g) bestimmen lässt. Ohne Beschränkung der Allgemeinheit wird nur eine Raumdimension g berücksichtigt, die sich auf das Koordinatensystem der ersten Repräsentation I₁(g) bezieht.
In einem ersten Schritt wird die Quotientenkurve Q(g) logarithmiert. Dadurch wird erreicht, dass nun aus dem Quotienten von T(γg) und T(g) eine Differenz der beteiligten Logarithmen entsteht, also $\ln [Q (g)] = \ln [T (γ g)] - \ln T [(g)] = : Δln [T (G)] .$
Dabei wurde ohne Einschränkung der Allgemeinheit davon ausgegangen, dass Q(g)>0 ist. Die Fallunterscheidung für Q(g)≤0, was bei einer oszillierenden Übertragungsfunktion T(g) auftreten kann, wurde der Übersichtlichkeit halber weggelassen.
Die Differenz der beiden Logarithmen ln[T(γg)] und ln[T(g)] wird hier mit dem Symbol Δln[T(G)] gekennzeichnet und dem zwischen den Frequenzen yg und g gelegenen arithmetischen Frequenzmittelpunkt G zugewiesen, der durch $G = g (1 + γ) / 2$
definiert ist. Die zur Ordinatendifferenz Δln[T(G)] gehörige Abszissendifferenz ΔG ergibt sich als Differenz der beteiligten Frequenzen γg und g, sodass man für die zugehörige Abszissendifferenz ΔG den Ausdruck $Δ G = g (γ - 1)$
erhält. Die Bildung der beiden Differenzen Δln[T(G)] und ΔG ist im Speziellen Beschreibungsteil in Figur 4 graphisch verdeutlicht. Mit der Ordinatendifferenz Δln[T(G)] aus Gl. (4) und der Abszissendifferenz ΔG aus Gl. (6) lässt sich schließlich ein Differenzenquotient D(G) bilden, der folgendermaßen definiert werden kann: $D (G) = \frac{Δln [T (G)]}{Δ G}$
Der Differenzenquotient D(G), der beispielhaft in Abbildung 4 dargestellt ist, kann unter Bezug auf das Referenzsystem der ersten Repräsentation bis zur Raumfrequenz g = g_N(1+1/γ)/2 gebildet werden, wobei g_N die Nyquistfrequenz der ersten Repräsentation darstellt. Entscheidend ist nun, dass der so gebildete Differenzenquotient D(G) eine finite Näherung der Ableitung des Logarithmus der Übertragungsfunktion T(G) nach der Raumfrequenz G darstellt, d.h. $D (G) \approx \frac{dln [T (G)]}{d G} = \frac{1}{T (G)} \frac{d T (G)}{d G}$
Diese Näherung lässt sich nutzen, um aus den Differenzenquotienten D(G) die Übertragungsfunktion T(G) eindeutig zu bestimmen.
Für den Lösungsweg sind zwei Fälle zu unterscheiden: Zum einen kann die Frequenz G als kontinuerliche Variable betrachtet werden, was eine direkte Integration der Gl. (8) ermöglicht und so einen einfachen Einblick in das allgemeine Lösungsprinzip gewährt. Zum anderen liegen in der Praxis der numerischen Berechnung die Werte von D(G) jedoch nur an diskreten Frequenzen G vor, was eine diskrete Summation anstatt einer kontinuierlichen Integration erforderlich macht, wobei besonderes Augenmerk auf die Metrik der erhaltenen Frequenzen G gelegt werden muss.
Betrachtet man die Frequenz G als kontinuierliche Variable, so erhält man den Logarithmus der Übertragungsfunktion ln[T(G)] näherungsweise durch Integration von D(G), d.h. $S (G) = \int_{0}^{G} D (G') d G' \approx \ln [T (G)],$
und die gesuchte Übertragungsfunktion T(G) ergibt sich schließlich durch Exponentiation, so dass $T (G) \approx Exp [S (G)] .$
Die üblicherweise anfallende Integrationskonstante in Gl. (9) wurde weggelassen, da sie exakt 0 ist. Geht man nämlich von der physikalisch vernünftigen Annahme aus, dass T(0) gleich 1 sein muss, so muss S(0) gleich 0 sein, wodurch die Integrationskonstante in Gl. (9) ebenfalls 0 sein muss.
Nach der verallgemeinerten Darstellung des Lösungsweges für kontinuierliche Frequenzen soll nun auf die Praxis der numerischen Berechnung eingegangen werden, bei denen die Eingangskurve Q(g) gemäß Gl. (3) nur an bestimmten diskreten Werten von g vorliegt. Letzterer Fall ist typisch für Frequenzspektren, die mittels FFT (Fast Fourier Transform) erhalten wurden. Seien die Repräsentationen I₁ und I₂ beispielsweise Bilder, die mit einem in diskrete Pixel unterteilten CCD-Sensor aufgenommen wurden. Dann lassen sich für ein Realraumfeld aus N Pixeln die zugehörigen diskreten Frequenzen des Fourierraums gemäß Abb. (3) mittels ganzzahliger Werte n=0, ±1, ±3 ... ± n_N indizieren, wobei der Index n_N die Nyquistfrequenz bezeichnet. Der Index n einer bestimmten Frequenz entspricht der Anzahl der Perioden der zugehörigen ebenen Welle Exp[2 i π n/N], die in dem Realraumfeld der Pixelanzahl N zu liegen kommen. Bei Wahl des Bildes 1 als Referenzbild kann man nun von Gl. (3) ausgehend g = n = 0, ±1, ±3 ... ± n_N setzen. Es ist jedoch hernach nicht mehr möglich, eine ganzzahlige Arithmetik aufrechtzuerhalten, da durch die Einführung des im Allgemeinen nicht ganzzahligen Vergrößerungsverhältnisses γ gemäß Gl. (5) Zwischenfrequenzen G auftreten, die weder der ganzzahligen Metrik des Bildes 1 noch der des Bildes 2 folgen, sondern der eines Bildes mit einer zwischen Bild 1 und 2 gelegenen Vergrößerung. Das Auftreten von Zwischenfrequenzen G ist also eine Folge der wechselweise symmetrischen Behandlung der Bilder 1 und 2. Nach dem Einsetzen dimensionsloser Indizes g = n = 0, ±1, ±3 ... ± n_N in Gl. (3) lässt sich die Rechnung jedoch bis Gl. (7) exakt wie zuvor beschrieben durchführen, wobei jetzt für G und ΔG im Allgemeinen keine ganzzahligen dimensionslosen Werte im Sinne eines Index oder einer Indexdifferenz auftreten, sondern beliebige rationale Werte. Eine von mehreren Möglichkeiten, die für die numerische Berechnung benötigte ganzzahlige Metrik des als Referenz gewählten Bildes 1 wiederherzustellen, sei nun im Folgenden beschrieben.
Da die Differenzenquotienten D(G) nach Einsetzen von g = n = 0, ±1, ±3 ... ± n_N in Gl. (3) nun an nicht- ganzzahligen dimensionslosen Werten von G zu liegen kommen, liegt es nahe, diese Werte von G als Stützstellen für ein Interpolationsverfahren zu benutzen, mittels dessen sich die Differenzenquotienten D(G) an den gewünschten Frequenzen G = n = 0, ±1, ±2 ... berechnen lassen. Diese Interpolation von D(G) lässt sich formal folgendermaßen darstellen als $D (G \in ℜ) \to D (n \in N) .$
Bezeichnet man den so an ganzzahligen Stützstellen n erhaltenen Differenzenquotienten D(n) als D_n, und analog T(k) als T_k, sowie S(k) als S_k, so wird aus der infintesimalen Integration von Gl. (9) nun eine Summation, wobei gilt: $S_{k} = \sum_{n = 0}^{k} D_{n} \approx \ln [T_{k}] .$
Schließlich erhält man analog zu Gl. (10) an der mit der ganzen Zahl k indizierten Frequenz des Referenzbildes 1 die gesuchte Lösung für die Übertragungsfunktion T_k mit $T_{k} \approx Exp [S_{k}] .$
Somit wird in einer besonders vorteilhaften Ausgestaltung der Erfindung aus dem kontinuierlich bzw. an diskreten Stützstellen erklärten Logarithmus von Q(g) ein kontinuierlich bzw. an diskreten Stützstellen erklärter Differenzenquotient für den Logarithmus von T(g) gebildet. Vorteilhaft wird dieser integriert bzw. summiert. Damit kann aus den unterschiedlich skalierten Repräsentationen direkt für jede Stützstelle im Arbeitsraum ein unabhängiger Wert für die Übertragungsfunktion T(g) als Lösung erhalten werden, ohne dass T(g) von vornherein auf eine bestimmte Klasse von Funktionen eingeschränkt wäre.
Bei hohen Raumfrequenzen g kann ΔG im Nenner von Gl. (7), das gemäß Gl. (6) mit der Raumfrequenz g ansteigt, zu groß werden. Dann wird die Näherung, dass D(G) als Ableitung von ln[T(G)] an der Stelle G angesehen wird, ungenau. Dies wiederum verschlechtert die Genauigkeit, mit der T(G) durch Summation bzw. Integration aus dem Differenzenquotienten D(G) bestimmt werden kann. Obwohl dies bei typischen Übertragungsfunktionen, die in der Regel nur langsam variieren, kaum der Fall ist, könnten bei untypischen Übertragungsfunktionen, die bei hohen Raumfrequenzen stark gekrümmt sind oder gar oszillieren, augrund der starr vorgegebenen Frequenzdifferenzen ΔG derartige Ungenauigkeiten auftreten. Je näher die beiden Funktionen I₁(g) und I₂(g) in ihrer Skalierung beieinander liegen, desto geringer sind diese Ungenauigkeiten.
Bei tiefen Ortsfrequenzen treten bei geringen Ortsfrequenzen sowohl im Zähler auch im Nenner des Differenzenquotienten nach Gl. (10) sehr kleine Werte auf, da typische Übertragungsfunktionen in der Regel nur langsam variieren. Starkes Rauschen in den Repräsentationen (Bildern) kann dann den Differenzenquotienten dominieren. Je weiter die beiden Funktionen I₁(g) und I₂(g) in ihrer Skalierung auseinanderliegen, desto größer wird ΔG im Nenner von Gl. (7) und desto mehr der Einfluss des Rauschens zurückgedrängt.
Somit ist mit dem Wert von γ ein unvermeidlicher Kompromiss zwischen der Genauigkeit bei hohen Ortsfrequenzen und der Genauigkeit bei tiefen Ortsfrequenzen einzustellen.
Daher werden in einer weiteren vorteilhaften Ausgestaltung der Erfindung für die Bildung des Differenzenquotienten für eine Position g im Arbeitsraum aus einer Mehrzahl von mindestens drei Repräsentationen zwei Ausschnitte aus verschiedenen Repräsentationen und damit zwei Funktionen I₁(g) und I₂(g) ausgewählt, die sich um einen vorgegebenen Skalierungsfaktor γ voneinander unterscheiden.
Sollte mit einem einzigen Vergrößerungsfaktor γ kein hinreichend guter Kompromiss gefunden werden können, der die Bildung des Differenzenquotienten sowohl bei niederen als auch bei hohen Frequenzen mit ausreichender Genauigkeit ermöglicht, kann auch mit mehreren Bildern, die in jeweils unterschiedlichen Vergrößerungsverhältnissen γ zueinander stehen, gearbeitet werden. Durch adaptive Wahl von γ für verschiedene Positionen g im Arbeitsraum kann der Differenzenquotient aus Gl. (7) an jeder gewünschten Stelle der Frequenzachse so gebildet werden, dass er gemäß Gl. (8) als zuverlässige und ausreichend genaue Näherung für den Differentialquotienten genutzt werden kann. Insbesondere bei oszillierenden Übertragungsfunktionen müssen die Stützpunkte für die Bildung des Differenzenquotienten so eng liegen (kleines ΔG), dass nicht eine komplette Oszillation innerhalb einer Spannbreite ΔG stattfindet. Ansonsten hat man bei höheren Frequenzen den klaren Fall einer Unterabtastung vorliegen. Deshalb ist hier in jedem Falle ein Mehrbildverfahren mit mehreren unterschiedlichen Vergrößerungen angeraten.
Dabei hängt es von γ ab, über welchen Bruchteil des Arbeitsraums die Übertragungsfunktion T(g) bestimmt werden kann (vergleiche Figur 4d). Im zweidimensionalen Fall für γ=2 kann die Übertragungsfunktion noch auf 56,25 % des Arbeitsraums bestimmt werden, im eindimensionalen Fall auf 75 % des Arbeitsraums. Insofern hat sich der Bereich 1 ≤ γ ≤ 2 als vorteilhaft erwiesen.
Die Ausnutzung des Arbeitsraums kann im zweidimensionalen Fall verbessert werden, indem die Auswertung nicht nur bis zur eindimensional definierten Nyquistfrequenz, also bis zum Rand des Spektrums, sondern bis zum √2-Fachen hiervon, also bis in die Ecken des Spektrums, fortgesetzt wird. Soll die Übertragungsfunktion T(g) ohne jeglichen Frequenzverlust bis zur eindimensionalen Nyquistfrequenz ermittelt werden, sind Werte von γ bis zu 2,41 zulässig.
Sofern I₁(g) und I₂(g) zwecks Rauschunterdrückung vor der Berechnung von Q(g) azimuthal gemittelt werden, ist hier als Nebeneffekt zu berücksichtigen, dass diese Mittelung in Richtung der Ecke des Frequenzraums in einem immer weiter eingeschränkten Winkelbereich möglich ist, der für das Eckpixel selbst dann zu Null zusammenschmilzt.
Als hohe Ortsfrequenzen sind typischerweise Frequenzen oberhalb der halben Nyquistfrequenz anzusehen. Bei diskreter Abtastung haben diese einen Index n>n_N/2, wobei n_N der Index der Nyquistfrequenz ist. Niedrige Ortsfrequenzen lassen sich über das Kriterium abgrenzen, dass der Abstand dieser Frequenzen bei diskreter Abtastung mindestens 1 Pixel betragen soll. In anderen Worten: Der Differenzenquotient soll wenigstens über den Abstand 1 Pixel gebildet werden. Dann hat man: γ n - n = n(γ-1) > 1, d.h. n > 1/(γ-1), wobei n die diskrete Indizierung der Frequenzen bedeutet. Für γ = 1.1 erhält man dann beispielsweise n > 10, für γ = 1.5 erhält man n > 2, für γ = 2 erhält man n > 1. Frequenzen die kleiner als das angegebene n sind, sind sicher als niedere Frequenzen anzusehen, also Frequenzen für die gilt: n < 1/(γ-1). Dies ist die unterste Grenze im Falle rauschfreier Repräsentationen (Bilder). Tritt zusätzlich Rauschen auf, sollte die entsprechende Unsicherheit mit berücksichtigt werden. Typischerweise sind Frequenzen mit n < 10 als niedrige Frequenzen anzusehen.
Ist der Arbeitsraum der Fourierraum, sind auf Grund der begrenzten Möglichkeiten einer Rauschmittelung die Werte von Q(g) für niedrige Frequenzen g mit besonders großen statistischen Fehlerbalken behaftet. Da die Repräsentationen I₁(g) und I₂(g) im Fourierraum für g in der Nähe von 0 typischerweise einen starken Gradienten aufweisen, können dort auch systematische Fehler verstärkt werden. Dem kann vorteilhaft entgegengewirkt werden, indem Q(g) für niedrige Frequenzen g aus Werten für benachbarte höhere Frequenzen g extrapoliert wird.
In einer weiteren vorteilhaften Ausgestaltung der Erfindung wird für die Übertragungsfunktion T(g) ein parametrisierter Ansatz gemacht und mit einem Optimierungsverfahren daraufhin optimiert, dass aus T(g) der Quotient Q(g) oder ein daraus abgeleiteter Funktionsverlauf bestmöglich reproduziert wird. Da zuvor bewiesen wurde, dass eine eindeutige Lösung für T(g) existiert, muss es sich bei einer durch Optimierung gefundenen Lösung für T(g), die Q(g) oder einen daraus abgeleiteten Funktionsverlauf exakt reproduziert, um diese eindeutige Lösung handeln. Als abgeleiteter Funktionsverlauf eignen sich z.B. der Logarithmus des Quotienten ln[Q(g)] aus Gl. (4) oder der Differenzenquotient D(G) aus Gl. (7).
In einer Ausgestaltung kann T(g) auf bestmögliche Übereinstimmung des Quotienten T(pg)/T(g) mit dem aus I₁(g) und I₂(g) erhaltenen Quotienten gemäß Gl. (3) optimiert werden, worin p ein Skalierungsfaktor ist. Als Wert für diesen Skalierungsfaktor kann vorteilhaft der vorbekannte oder bei der Anfertigung der Repräsentationen I₁(x) und I₂(x) eingestellte Skalierungsfaktor γ, um den sich diese Repräsentationen unterscheiden, eingestellt werden. Der Wert von p kann im Rahmen eines Optimierungsverfahrens bestimmt werden. Auch eine Kombination ist möglich. Voraussetzung hierfür ist, dass bei einer Änderung von p die Ausschnitte I₁(g) und I₂(g) durch passende Gebietsauswahl im Realraum entsprechend angepasst werden und eine Erkennung für eine falsche Gebietsauswahl vorgesehen ist.
Ein parametrisierter Ansatz für die Übertragungsfunktion lässt sich beispielsweise durch eine lineare Superposition aus Gaußkurven, fallenden Exponentialfunktionen, Lorentzfunktionen oder ähnlichen Funktionen herstellen, wobei die benutzte Anzahl von Parametern, mittels derer die Gewichtungsfaktoren der beteiligten Funktionen sowie deren Breite beschrieben wird, typischerweise kleiner als 10 ist. Durch den parametrisierten Ansatz wird die Klasse der möglichen Funktionen T(g), die als Lösung erhalten werden können, eingeschränkt. Dies aber bewirkt im Gegenzug, dass der Ansatz robuster gegen Rauschen in den Repräsentationen I₁(g) und I₂(g) ist. Damit kann zumindest teilweise der Tendenz eines jeden Entfaltungsverfahrens entgegengewirkt werden, dieses Rauschen zu verstärken. Meist ist nach Erhalt einer stark verrauschten Lösung die Anwendung eines Glättungsverfahrens nötig, welches benachbarte diskrete Frequenzen der Übertragungsfunktion im Sinne einer lokalen Mittelwertbildung miteinander verknüpft.
In einer weiteren vorteilhaften Ausgestaltung der Erfindung wird vor der Transformation in den Arbeitsraum die Repräsentation mit der geringeren Ausdehnung im Realraum auf die Ausdehnung der anderen Repräsentation hochinterpoliert. Dann sind beide Repräsentationen anschließend mit einer Transformation der gleichen Dimension in den Arbeitsraum überführbar. Dies wird damit erkauft, dass die Interpolation selbst wieder eine Übertragungsfunktion ist, mit der nur eine der beiden Repräsentationen behaftet ist. Ihr Einfluss auf das Endergebnis muss charakterisiert werden und hängt von der Ordnung der Interpolation ab.
Analog kann vor der Transformation in den Arbeitsraum eine der Repräsentationen im Realraum gedreht und/oder gespiegelt werden. Damit kann eine vorhandene Drehung und/oder Spiegelung zwischen zwei ansonsten identischen Repräsentationen ausgeglichen werden. Dann können in dem Fall, dass die Übertragungsfunktion T(g) rotationssymmetrisch ist, auch Repräsentationen, die sich neben einer unterschiedlichen Skalierung auch noch durch Drehung und/oder Spiegelung unterscheiden, zur Deckung gebracht werden, so dass die Quotientenbildung gemäß Gl. (3) möglich ist. Dabei ist zu berücksichtigen, dass die Drehung und die Spiegelung weitere zusätzliche Übertragungsfunktionen einbringen, mit denen nur eine der beiden Repräsentationen behaftet ist. Ist der Realraum in diskrete Pixel unterteilt, kommt das Ergebnis der Drehung und/oder Spiegelung zudem im Allgemeinen nicht auf ganzzahligen Pixelkoordinaten zu liegen; dann wird eine Interpolation, die eine weitere Übertragungsfunktion einbringt, notwendig.
Die Drehung und/oder Spiegelung kann beispielsweise ausgeführt werden, indem die Quotientenbildung gemäß Gl. (3) dahingehend abgeändert wird, dass der Quotient Q(g)=I₂(g)/I₁(Dg) gebildet wird, worin D die Transformationsmatrix der Drehung und/oder Spiegelung ist. Sofern die Funktionen I₂(g) und I₁(g) im Arbeitsraum hinreichend glatt sind, kann zumindest die zusätzliche Interpolation auf ganzzahlige Pixelkoordinaten entfallen.
In einer weiteren vorteilhaften Ausgestaltung der Erfindung werden die Repräsentationen im Arbeitsraum in Polarkoordinaten transformiert und in dieser Darstellung azimuthal gemittelt. Häufig ist die Annahme gerechtfertigt, dass es sich bei der gesuchten Übertragungsfunktion um eine rotationssymmetrische Funktion handelt. Dies ist beispielsweise bei elektronenmikroskopischen Bildern der Fall. Es wurde erkannt, dass die Division gemäß Gl. (3) instabil ist, insbesondere wenn die als Divisor dienende Repräsentation I₁(g) stark verrauscht ist. Das Rauschen kann durch die Division verstärkt werden, was wiederum die Genauigkeit beeinflusst, mit der von Q(g)=T(γg)/T(g) auf T(g) zurückgerechnet wird. Es kann daher in speziellen Fällen von Vorteil sein, die Repräsentationen durch die Mittelung vorab zu glätten. Man erhält damit eine eindimensionale Funktionen I'(g), wobei die azimuthale Mittelung der komplexwertigen Fourierkoeffizienten I(g,ϕ) beispielsweise folgendermaßen beschrieben werden kann: $I' (g) = {(\frac{1}{2 π} \int_{0}^{2 π} {| I (g, ϕ) |}^{2} d ϕ)}^{1 / 2}$
Diese azimuthale Mittelung ist ausdrücklich keine zwingende Voraussetzung dafür, dass die Berechnung der Übertragungsfunktion T(g) gemäß den Gln. (4) bis (13) in einer Raumdimension durchgeführt werden kann. Das allgemeine Problem, diese Funktion in zwei oder mehr Dimensionen Punkt für Punkt zu bestimmen, lässt sich in Polarkoordinaten durch Zerlegung in unabhängige eindimensionale radiale Sektionen durch den Koordinatenursprung immer auf eine Dimension zurückführen, auch wenn zuvor keine azimuthale Mittelung durchgeführt wird. Die azimuthale Mittelung ist lediglich eine Option im Sonderfall, dass die Übertragungsfunktion T(g) rotationssymmetrisch ist.
Die durch die azimuthale Mittelung erreichte Glättung kann jedoch zu einer systematischen Anhebung der so erhaltenen Kurve I'(g) führen, da die Betragsquadratbildung im Integranden des Integrals (14) stets einen positiven quadratischen Rauschterm enthält. Dieser mögliche systematische Rauschbeitrag kann gesondert entfernt werden, wenn man das tatsächlich vorhandene Rauschspektrum N(g) kennt. Durch quadratische Subtraktion des Rauschspektrums lässt sich die Bildintensität I'(g) korrigieren und man erhält für die rauschkorrigierte Bildintensität I(g) schließlich den Ausdruck: $I (g) = {(I' {(g)}^{2} - N {(g)}^{2})}^{1 / 2}$
Häufig ist der Rauschbeitrag N(g)² aus Gleichung (15) über beinahe den gesamten Raumfrequenzbereich sehr viel kleiner als der vom Signal dominierte Beitrag I'(g)² und die Korrektur nach Gleichung (15) macht sich nur in den äußersten hochfrequenten Bereichen eines Spektrums, wo das Objektspektrum und die Übertragungsfunktion sehr kleine Werte annehmen, überhaupt erst bemerkbar. Da die Korrektur nur über einen relativ schmalen Raumfrequenzbereich am Rande der Fouriertransformierten überhaupt wirksam ist, kann man auch annehmen, dass sich das Rauschspektrum N(g) über diesen schmalen Bereich nicht signifikant ändert und sich näherungsweise durch eine Konstante c ersetzen lässt. Man kann in einem solchen Falle näherungsweise anstelle der Gl. (15) auch den folgenden vereinfachten Ausdruck zur Rauschkorrektur benutzen: $I (g) \approx {(I' {(g)}^{2} - c^{2})}^{1 / 2} .$
Vorteilhaft wird also aus den azimuthal gemittelten Repräsentationen jeweils ein Rauschspektrum oder ersatzweise ein konstanter Rauschuntergrund herauskorrigiert.
Alternativ oder in Kombination zu der azimuthalen Mittelung können mehrere identisch skalierte Repräsentationen entweder im Realraum oder im Arbeitsraum vor der Quotientenbildung nach Gl. (3) gemittelt werden. Dieses Mittel wird dann als im weiteren Verfahren als eine Repräsentation I₁(x) bzw. I₁(g) oder I₂(x) bzw. I₂(g) verwendet. Diese alternative Möglichkeit ist besonders wichtig in dem Fall, in dem die Übertragungsfunktion T(g) nicht rotationssymmetrisch ist, weil dann eine azimuthale Mittelung nicht erlaubt ist.
Zur weiteren Rauschunterdrückung können die Repräsentationen I₁(g) und I₂(g) vor der Quotientenbildung nach Gl. (3) noch mit lokaler Faltung geglättet werden
Im Rahmen der Erfindung wurde ein weiteres Verfahren zur Bestimmung der Übertragungsfunktion eines signalverarbeitenden Systems entwickelt. Dieses Verfahren unterscheidet sich von dem Verfahren gemäß Hauptanspruch dadurch, dass noch keine Eingabedaten I₁(x) und I₂(x) vorliegen, sondern nur das signal verarbeitende System als solches. Daher werden bei diesem Verfahren mit dem System zunächst

a) mindestens zwei unterschiedlich skalierte Repräsentationen I₁(x) und I₂(x) eines Objekts oder
b) mindestens eine Repräsentation I₁(x) eines ersten Objekts und eine Repräsentation I₂(x) eines hierzu bis auf einen Skalenfaktor geometrisch identischen Objekts erzeugt.

Sodann wird mit diesen Repräsentationen das zuvor beschriebene Verfahren durchgeführt.
Die beiden Repräsentationen I₁(x) und I₂(x) können beispielsweise folgendermaßen erhalten werden:

Alternative a):

Es wird ein und dasselbe Objekt benutzt, um zwei unterschiedlich skalierte Repräsentationen I₁(x) und I₂(x) hiervon zu erzeugen. Hierfür gibt es zwei Möglichkeiten:

1. In einer weiteren vorteilhaften Ausgestaltung der Erfindung wird die unterschiedliche Skalierung der Repräsentationen I₁(x) und I₂(x) durch eine Änderung der Skalierung des vom Objekt am Eingang des Systems erzeugten Eingangssignals O(x) eingestellt. Das bedeutet, dass die unterschiedlich skalierten Repräsentationen auf dem Detektor des Systems mit Hilfe des signalverarbeitenden Systems erzeugt werden. Bei einem (Elektronen-)Mikroskop kann beispielsweise die Vergrößerung geändert werden, so dass man verschieden vergrößerte Abbildungen desselben Objekts vorliegen hat. Die Benutzung der optischen Zoom-Funktion bei einem Photoapparat hat eine analoge Wirkung, wobei sichergestellt werden muss, dass sich die Abbildungseigenschaften des Objektivs durch den Zoomvorgang nicht signifikant verändern. Bei einer EELS-Messung wird die Vergrößerung durch die am Spektrometer eingestellte Dispersion bestimmt. Je höher die Dispersion ist, desto weiter wird das Spektrum räumlich auf dem CCD-Chip der Kamera gestreckt, wenn die Energie verändert wird. Dadurch wird die Abtastrate der Kamera in Bezug auf die Energie größer, und das Spektrum kann mit höherer Energieauflösung bestimmt werden. Dies entspricht einer höheren Vergrößerungsstufe bei einer photographischen oder mikroskopischen Bildaufnahme.
2. In einer weiteren vorteilhaften Ausgestaltung der Erfindung wird die unterschiedliche Skalierung der Repräsentationen I₁(x) und I₂(x) durch eine Änderung des räumlichen Abstands zwischen dem Objekt und dem signalverarbeitenden System eingestellt. Das signalverarbeitende System selbst bleibt dabei im Gegensatz zur vorherigen Ausgestaltung unverändert und erhält bereits zwei unterschiedlich skalierte Eingangssignale O(x). Eine typische Anwendung ist die Aufnahme ein und desselben Objekts mittels eines Photoapparates aus verschiedenen Entfernungen, ohne dass die Einstellungen des Photoapparats geändert werden.

Alternative b):

Man benutzt zwei oder mehrere tatsächlich physikalisch vorliegende Objekte bzw. Signalformen als Eingangssignal, die wechselseitig identische Form (Funktionsverlauf) besitzen, jedoch von unterschiedlicher Ausdehnung (Skalierung) sind. Ein Beispiel hierfür ist das Abphotographieren zweier unterschiedlich vergrößerter und hinreichend scharfer Ausdrucke desselben Motivs mittels eines Photoapparates unter identischen Abbildungsbedingungen. Der Abstand zum signalverarbeitenden System sowie die Einstellungen des Systems selbst können unverändert bleiben.
Zusätzlich kann man insbesondere im Zusammenhang mit Ausgestaltung 1 von Alternative a) anstatt eines herkömmlichen Objektes auch Rauschen als Eingangssignal benutzen, wobei der funktionale Verlauf des Rauschspektrums aber nun nicht mehr bekannt sein muss. Es kann also im Gegensatz zum Stand der Technik ein beliebiges Rauschsignal als Eingangssignal benutzt werden, indem dieses in mindestens zwei verschieden skalierten Varianten registriert wird.
In einer weiteren vorteilhaften Ausgestaltung der Erfindung wird mindestens eine der Repräsentationen I₁(x) und I₂(x) als Zusammenfassung von mindestens zwei Einzelrepräsentationen erzeugt. Beispielsweise können mit einem Photoapparat oder einem (Elektronen-)Mikroskop bei identischen Einstellungen mehrere Bilder nacheinander aufgenommen werden. Die Einzelrepräsentationen können direkt im Realraum addiert werden. Alternativ können Ausschnitte aus den Einzelrepräsentationen ausgewählt und deren Entsprechungen im Arbeitsraum zu den Funktionen I₁(g) bzw. I₂(g) zusammengefasst werden.
Hinsichtlich eines Elektronenmikroskops als signalverarbeitendes System ist anzumerken, dass die Vergrößerungsänderung durch die Hintereinanderschaltung der sogenannten Zwischenlinsen und der sogenannten Projektorlinse bewirkt wird, welche als fehlerfrei angesehen werden. Diese Linsen deshalb optisch "neutral" und wirken im Sinne einer mathematisch idealen Vergrößerung; es wird lediglich die Objektfunktion bezüglich des Detektors umskaliert. Aufgrund der optischen Neutralität dieser Linsen entspricht eine mit dem hier dargestellten Verfahren gemessene Übertragungsfunktion nur der Übertragungsfunktion des Detektors und ist nicht mehr vom vorgeschalteten Linsensystem beeinflusst (siehe Figur 1 im Speziellen Beschreibungsteil).
Hinsichtlich eines Photoapparats als signalverarbeitendes System ist anzumerken, dass sich die Übertragungsfunktion aus einem Beitrag des Kameraobjektivs und einem Beitrag des Films bzw. des CCD-Sensors zusammensetzt (siehe Figur 1 im Speziellen Beschreibungsteil). Hier ist es wichtig dass sich der Beitrag des Kameraobjektivs zwischen den einzelnen Aufnahmen nicht ändert. Bei der zweiten Möglichkeit von Alternative a) kann dies dadurch erreicht werden, dass sich das Objekt bei allen Aufnahmen jenseits der sogenannten hyperfokalen Entfernung befindet, womit ein individuelles "Scharfstellen" der einzelnen Aufnahmen nicht mehr nötig ist. Auf diese Weise kann man mit der Entfernungseinstellung "Unendlich" durch reine Abstandsänderung unterschiedlich vergrößerte und ausreichend scharfe Aufnahmen ein und desselben weit entfernten Objekts erhalten. Bei der Alternative b) kann der Beitrag des Kameraobjektivs konstant gehalten werden, indem bei festem Objektabstand zwei verschieden vergrößerte und hinreichend scharfe Ausdrucke eines beliebigen Motivs als Objekt benutzt werden.
Im Rahmen der Erfindung wurde ein weiteres Verfahren zur Bestimmung der Übertragungsfunktion T(x) eines signalverarbeitenden Systems und des Eingangssignals O(x), das ein unbekanntes Objekt am Eingang dieses Systems erzeugt, aus mindestens zwei Repräsentationen I₁(x) und I₂(x) des Objekts, die das System aus unterschiedlichen Skalierungen dieses Eingangssignals O(x) erzeugt hat, entwickelt.
Dabei können die unterschiedlichen Skalierungen des Eingangssignals O(x) auf jede zuvor beschriebene Weise erzeugt werden. Insbesondere kann das System selbst dieses Eingangssignal skalieren, oder es können als Repräsentationen I₁(x) und I₂(x) eine Repräsentation I₁(x) eines ersten Objekts und eine Repräsentation I₂(x) eines hierzu bis auf einen Skalenfaktor geometrisch identischen Objekts verwendet werden. Dabei können die beiden Repräsentationen analog dem zuvor Gesagten wie bei dem Verfahren gemäß Hauptanspruch als Eingabedaten vorgegeben sein oder aber zu Beginn des Verfahrens wie bei dem Verfahren gemäß Nebenanspruch erzeugt werden.
Erfindungsgemäß werden parametrisierte Ansätze für die Objektfunktion O(x) und für die Übertragungsfunktion T(x) selbstkonsistent dahingehend optimiert, dass T(x) angewendet auf O(x) die beiden Repräsentationen I₁(x) und I₂(x) reproduziert.
Da es nach der Ermittlung der Übertragungsfunktion T(x) möglich ist, durch Entfaltung beispielsweise des Bildes 1 auch die Objektfunktion O(x) zu erhalten, steckt also die komplette Information sowohl über die Objektfunktion O(x) als auch über die Übertragungsfunktion T(x) in dem aus den Repräsentationen I₁(x) und I₂(x) bestehenden Eingabedaten. Mit dem erfindungsgemäß aufgestellten Gleichungssystem $I_{1} (x) = O (x) \otimes T (x)$
$I_{2} (x) = O (x) \otimes T (γ x)$
können beide gesuchten Funktionen O(x) und T(x) synchron bestimmt werden. Bei diesem Verfahren wird gar kein expliziter Gebrauch von der Quotientenkurve nach Gl. (7) und den daraus abgeleiteten funktionalen Zusammenhängen der Gln. (8) bis (14) gemacht. Ist man durch Probieren oder durch Anwendung eines beliebigen numerischen Optimierungsverfahrens in der Lage, eine Objektfunktion und eine Übertragungsfunktion zu finden, welche die aus Bild 1 und 2 bestehenden Eingabedaten reproduzieren, dann ist die gefundene Lösung auch eindeutig.
Dieses Verfahren ist eine logische Verallgemeinerung der vorherigen Verfahren. Alle für die vorherigen Verfahren offenbarten Maßnahmen sind daher auch auf dieses Verfahren anwendbar. Insbesondere kann auch der Vergrößerungsfaktor γ in die Optimierung einbezogen werden. Die Optimierung kann im Realraum oder in einem Arbeitsraum stattfinden, wobei der Arbeitsraum insbesondere der Frequenzraum sein kann.
Durch Einsatz von Kompressionsmethoden ist es möglich, die Parameteranzahl für das Objekt in die Größenordnung von 100 zu bringen. In diesem Falle sind Optimierungsverfahren wie z.B. Gradientenverfahren, Simulated Annealing oder genetische Algorithmen durchaus in der Lage, zuverlässig eine Lösung zu finden.

Spezieller Beschreibungsteil

Nachfolgend wird der Gegenstand der Erfindung an Hand von Figuren erläutert, ohne dass der Gegenstand der Erfindung hierdurch beschränkt wird. Es ist gezeigt:

Figur 1: Optisches Abbildungssystem als Beispiel für ein signalverarbeitendes System.
Figur 2: Skizze des erfindungsgemäßen Verfahrens.
Figur 3: Auswahl des Ausschnitts aus der zweiten Repräsentation, so dass beide Repräsentationen sich auf den gleichen Bereich des Objekts beziehen.
Figur 4: Graphische Erläuterung zur Berechnung des Differenzenquotienten D(G).
Figur 5: Ausfiihrungsbeispiel des Verfahrens an einem Elektronenmikroskop.

Figur 1 verdeutlicht, inwiefern ein optisches System OS als Signalverarbeitungssystem anzusehen ist. Das System umfasst ein Transfersystem T und einen Detektor D. Das Transfersystem T projiziert das vom Objekt O ausgehende Licht auf den Detektor D, so dass dort ein scharfes Bild entsteht. Die Übertragungsfunktion des Systems ist aus einem Beitrag des Transfersystems T und einem Beitrag des Detektors D zusammengesetzt.
Figur 2 illustriert den Gang des erfindungsgemäßen Verfahrens. Ohne Beschränkung der Allgemeinheit wird hier von einem quadratisch bemessenen Detektor ausgegangen, wobei anders bemessene Detektoren vom Prinzip her nicht ausgeschlossen sind. Der Detektor bestehe aus MxM Pixeln und habe entlang einer Richtung die physikalische Abmessung d_D=b Pixel (Detektorreferenz DR). Gegeben seien zwei mittels dieses Detektors gemachte Aufnahmen I1 und I2' als Repräsentationen desselben Objekts, die bei jeweils verschiedenen Vergrößerungen gemacht wurden. Die Aufnahme I1 sei bei einer beliebigen Vergrößerung gemacht und enthalte ein Objektgebiet der Ausdehnung d_O=a Nanometer (Objektreferenz OR). Die Aufnahme I2' sei bei einer Vergrößerung gemacht, die sich im Bezug auf I1 um den Faktor 1/γ unterscheidet. Ohne Beschränkung der Allgemeinheit soll hier von γ>1 ausgegangen werden, was im vorliegenden Fall bedeutet, dass die Vergrößerung von I2' geringer ist als die von I1. Aufgrund der geringeren Vergrößerung von I2' enthält I2' dann ein Objektgebiet, das im Vergleich zu I1 entlang einer Dimension um den Faktor γ größer ist, d.h. I2' enthält ein Objektgebiet der Ausdehnung d_O=γa.
Es ist nun möglich, ein Gebiet auf I2' auszuwählen, das bezüglich der Abmessung des dargestellten Objektgebiets genau den Abmessungen des auf I1 erfassten Objektgebiets entspricht.
Die Positionierung des in I2' ausgewählten Gebietes kann dann so gewählt werden, dass das gesamte I1 und der ausgewählte Bereich I2 von I2' genau dasselbe Objektgebiet wiedergeben. Die so vorgenommene Auswahl eines Bildausschnittes bei I2' entspricht einer fiktiven physikalischen Verkleinerung des Detektors von seiner tatsächlichen physikalischen Abmessung d_D=b auf die Abmessung d_D=b/γ. Es ist wichtig zu erwähnen, dass diese Verkleinerung nicht durch Skalierung der Pixel, sondern durch das Weglassen von Randpixeln bei gleicher Pixelgröße erfolgt. Da ein Detektor meist aus diskreten Pixeln aufgebaut ist, ergibt sich für die so entstandene kleinere Pixelzahl N des fiktiven Detekors N = NINT(M/γ), wobei die Funktion NINT die nächstliegende ganze Zahl bedeutet. Bei einem hinreichend großen Betrag von N, der bei gängigen Bilddetektoren typischerweise größer als 1000 ist, kann der bei beliebigem γ durch die ganzzahlige NINT-Rundung entstandene Fehler in Bezug auf eine exakt deckungsgleiche Ausschnittswahl in der Regel vernachlässigt werden.
Im nächsten Schritt werden die Bilder vom Realraum R in den Fourierraum F transformiert, wobei das aus MxM Pixeln bestehende I1 einer diskreten MxM Fouriertransformation unterzogen wird. Das auf die Größe NxN reduzierte deckungsgleiche I2 wird einer NxN Fouriertransformation unterzogen. Die mittlere Zeile von Figur 2 enthält von links nach rechts die Spektren der Objektfunktion O(g) und der Übertragungsfunktion T(g) für die Bilder I1, I2' und I2.
Die Fouriertransformation kann mit dem FFT (Fast-Fourier-Transform) Algorithmus berechnet werden. Da insbesondere N nicht unbedingt eine Potenz von 2 ist, eignen sich jedoch viele gängige FFT Programme nicht, da diese häufig auf einer Relation M,N=2ⁿ basieren (Radix 2 Algorithmus). Es kann jedoch stets ein sogenannter mixed-radix FFT Algorithmus herangezogen werden, der eine Zerlegung von N in allgemeine Primzahlen, die nicht notwendigerweise 2 sind, ermöglicht. Bei günstiger Zerlegbarkeit von N in viele kleine Primzahlen kann der mixed-radix Algorithmus eine Recheneffizienz erreichen, die sehr nahe an der des Radix-2 FFT Algorithmus liegt. Ist N im ungünstigsten Fall selbst eine Primzahl, so erniedrigt sich die Effizienz des mixed-radix Algorithmus auf die Effizienz einer direkten Fouriertransformation. Letztlich kann jedoch immer wenigstens eine direkte Fouriertransformation für beliebige Zahlen M,N verwendet werden. Bei einer typischen Größenordnung von N ≈ 10³ ist selbst der ungünstigste Fall einer direkten Fouriertransformation auf zeitgemäßen Computern kein Problem mehr. Bei Größenordnungen von N ≈ 10⁴ könnte man durch eine zusätzliche künstliche Erhöhung oder Erniedrigung von N um typischerweise ± 1 die Primzahlzerlegbarkeit verbessern und die Rechengeschwindigkeit dadurch erhöhen, wobei der durch die künstliche Rundung eintretende Skalierungsfehler von etwa 10^-4 bei den meisten praktischen Anwendungen ebenfalls vernachlässigbar ist.
Zum besseren Verständnis des Verfahrens werden im Folgenden anhand von Figur 2 zwei Szenarien verglichen, wobei das erste Szenario keine Gebietsauswahl bei I2' vorsieht, während beim zweiten Szenario die eben erläuterte Gebietsauswahl vorgenommen wird.
Hätte man das gesamte I2' mit seiner Objektausdehnung d_O=γa und seiner ursprünglichen Detektorausdehnung d_D=b ohne weiteres Zutun ebenso wie I1 einer diskreten MxM Fouriertransformation unterworfen, so würden die zu den I1 und I2' gehörigen Objektspektren voneinander abweichen, zum einen weil sie verschieden skaliert sind, zum anderen weil die Objektgebiete, die zur jeweiligen Transformierten beitragen, nicht identisch sind. Im Gegensatz zum Objektspektrum wäre jedoch die allein durch den Detektor gegebene Übertragungsfunktion für beide Bilder identisch, da diese nicht vom Bildinhalt, sondern nur von dem in beiden Fällen auf identische Weise genutzten Detektor abhängt (siehe mittlere Spalte von Figur 2).
Umgekehrt verhält es sich mit dem Vergleich der Fouriertransformierten von I1 und I2 nach der Reduktion von I2' auf die Objektausdehnung d_O=a mit der zugehörigen fiktiven Detektorausdehnung d_D=b/γ und der zugehörigen Pixelanzahl N. Nun kommen trotz der ursprünglich verschieden gewählten Vergrößerung wieder identische Objektfrequenzen auf identischen Frequenzen der Fouriertransformierten zu liegen, während die zugehörige Koeffizienten der Übertragungsfunktion jetzt nicht mehr auf identischen Frequenzen liegen. Das Objektspektrum und das Spektrum der Übertragungsfunktion haben also in den beiden besprochenen Szenarien bezüglich ihrer gegenseitigen Passung jeweils die Rolle vertauscht (siehe rechte Spalte von Figur 2).
Figur 3 verdeutlich das Prinzip der Passung von Objektspektren aus verschieden vergrößerten Aufnahmen des gleichen Objekts mittels einer diskreten Fouriertransformation ist für ein diskretes Objektspektrum, das aus nur zwei Kosinuswellen verschiedener Frequenz besteht, im Detail. Hierbei wurde der Übersichtlichkeit halber auf die Darstellung einer Übertragungsfunktion verzichtet. Die linke Spalte von Figur 3 enthält Darstellungen im Realraum R, die rechte Spalte von Figur 3 enthält Darstellungen im Fourierraum F. Die in Figur 3 auf der rechten Seite zur Indizierung der Frequenzachse benutzten Zahlen bedeuten die Anzahl der Perioden pro Bildausschnitt.
Eine Kosinuswelle, die im Realraum beispielsweise n Perioden innerhalb eines Bildausschnitts von d_D=b aufweist, lässt sich in zwei ebene Wellen der Form Exp[2 π i n/b x] und Exp [-2 π i n/b x] zerlegen. Teilbild A von Figur 3 zeigt eine solche Überlagerung zweier Kosinuswellen im Realraum, die auf dem Detektor M Pixel belegt. Nach der Fouriertransformation (Teilbild a von Figur 3) entstehen für jede Kosinuswelle zu den diskreten Frequenzen n/b und -n/b Fourierkoeffizienten, die jeweils über der Anzahl n der im Bild enthaltenen Perioden dargestellt sind. Das Fourierspektrum weist für zwei Kosinuswellen also vier von Null verschiedene Fourierkoeffizienten auf.
Teilbild B von Figur 3 zeigt die gleiche Überlagerung der beiden Kosinuswellen, die ebenfalls einen Ausschnitt von d_D=b enthält, der auf dem Detektor M Pixel breit ist. Die Überlagerung wurde hier jedoch bei einer anderen Vergrößerung mit γ=3/2 aufgenommen. In der Fouriertransformation (Teilbild b von Figur 3) verschieben sich dementsprechend die Fourierkoeffizienten gegenüber Teilbild a zu anderen Frequenzen.
Wird aus Teilbild B von Figur 3 nun der in Teilbild C gezeigte Ausschnitt von d_D=b/γ ausgewählt, der M/γ Pixeln auf dem Detektor entspricht, dann enthält dieser Ausschnitt für beide Kosinuswellen die gleiche Anzahl Perioden wie der in Teilbild A gezeigte Ausschnitt. Durch diese identische Gebietswahl wird also erreicht, dass abgesehen von der gewählten Abtastrate exakt identische Szenarien vorliegen. Dementsprechend enthält das in Teilbild c gezeigte Fourierspektrum für die gleichen Frequenzen Beiträge wie das in Teilbild a gezeigte, aus Teilbild A erzeugte Fourierspektrum. Die Möglichkeit, eine wechselseitige Passung der Objektfrequenzen bei gleichzeitiger wechselseitiger Identität der Fourierkoeffizienten herstellen zu können, beruht auf der Tatsache dass das Ergebnis einer Fouriertransformation bei hinreichend großer Anzahl der Stützpunkte M, N unabhängig von der Anzahl der Stützpunkte und somit unabhängig von der gewählten Abtastrate ist. Dadurch können sowohl die Objektfrequenzen der beiden Bilder als auch die zu den jeweiligen Objektfrequenzen gehörigen Fourierkoeffizienten des Objektspektrums zur Passung gebracht werden.
Die hier beschriebene Passung der Objektspektren zweier unterschiedlich vergrößerter Aufnahmen über die Fouriertransformation in unterschiedlichen Dimensionen M und N könnte alternativ auch auf andere Art erreicht werden: Wenn man in Figur 2 das aus I2' ausgewählte Untergebiet I2 der Ausdehnung d_O=a und der Pixelzahl NxN mittels eines Interpolationsverfahrens im Realraum auf die Pixelzahl MxM von I1 bringt, kann man danach bei I2 eine diskrete Fouriertransformation der Dimension MxM vornehmen, um die Passung der Objektfrequenzen zwischen I1 und I2 herzustellen. Allerdings bringt diese Methode mit sich, dass jede Interpolation wiederum eine Übertragungsfunktion besitzt, die zusätzlich charakterisiert werden müsste, da I2 dann mit dieser Übertragungsfunktion behaftet ist, während I1 frei von dieser zusätzlichen Übertragungsfunktion ist. Eine Steigerung der polynomialen Ordnung der Interpolation verringert den Fehler der Interpolation, welcher dann bei der höchstmöglichen Ordnung N dann minimal würde. Jedoch entspricht der vorhin aufgezeigte Weg mittels der direkten Fouriertransformation der Dimension N×N genau dieser alternativen Realrauminterpolation der Ordnung N. Aufgrund der deutlich höheren Recheneffizienz und der absolut symmetrischen Behandlung der beiden Bilder, die es nicht notwendig macht, ein Bild zu interpolieren, während man das andere Bild unbehandelt lässt, ist der vorhin beschriebene Weg der direkten Fouriertransformation mit verschiedenen Dimensionen M und N jedoch eindeutig vorzuziehen.
Unabhängig von der Art der Bildprozessierung ist die Tatsache, dass nach dem Erreichen einer Passung der Objektfrequenzen die jeweils höchste darstellbare Objektfrequenz der Bilder I1 und I2 unterschiedlich ist. Bei der hier gewählten Darstellung, bei der ohne Beschränkung der Allgemeinheit I2 eine geringere Vergrößerung aufweist, ist auch die höchste nutzbare Objektfrequenz (Nyquistfrequenz) bei I2 niedriger als bei I1. Dies liegt daran, dass man bei der niedrigeren Abtastrate des Objektes, wie sie beim geringer vergrößerten I2 vorliegt, entsprechend weniger feine Objektdetails auflösen kann. Wie man aus Figur 2 leicht erkennen kann, besteht zwischen der Nyquistfrequenz g^N ₂ des frequenzgepassten zweiten Bildes und der Nyquistfrequenz g^N ₁ des ersten Bildes die Beziehung: g^N ₂ = 1/γ g^N ₁.
Stellt man ein zweidimensionales Bild in Polarkoordinaten dar, so betrifft die vorhin beschriebene Passung der Objektfrequenzen nur die radiale Koordinate g, während die azimuthale Koordinate ϕ davon unberührt bleibt. Weiterhin setzt der in der Beschreibung zum Hauptanspruch aufgeführte Beweis, dass es eine eindeutige Lösung für die Übertragungsfunktion T(g) gibt, keine Abhängigkeiten zwischen verschiedenen Azimuthwinkeln voraus und erzeugt diese auch nicht. Daher kann jede azimuthale Richtung ϕ eines Spektrums für sich gesondert behandelt werden. Obwohl hier zweidimensionale Bilder behandelt werden, genügt es, eine einzige exemplarische Raumrichtung ϕ zu behandeln. Der vorgenannte Beweis gilt also bei voller Aufrechterhaltung der Zweidimensionalität für jede beliebige Raumrichtung ϕ einzeln.
In Figur 4a ist beispielhaft eine Quotientenkurve Q(g) dargestellt, die sich bis zur Grenzfrequenz 1/γ g_N bilden lässt, wobei g_N die Nyquistfrequenz des ersten Bildes ist, und γ > 1 vorausgesetzt wurde. Es sind beispielhaft drei Stützstellen g_A, g_B und g_C eingezeichnet. Figur 4b zeigt beispielhaft die noch zu erzielende Dekomposition der dargestellten Quotientenkurve Q(g) in Zähler T(γg) und Nenner T(g). In Figur 4b ist ersichtlich, wie durch die Quotientenbildung stets jeweils zwei Punkte g und yg auf ein und derselben Kurve T(g) in eine paarweise Relation gesetzt werden: Zwischen den Punkten g_A und γg_A fällt die Kurve von T(g) gerade bis auf den Wert ab, den die Kurve T(γg) an der Stelle g_A annimmt.
Diese Relation bleibt erhalten, wenn beide Kurven logarithmiert werden (Figur 4c). Die Differenz zwischen den beiden Kurven ln[T(g)] und ln[T(γg)] ist nun genau der Logarithmus der gegebenen Quotientenkurve Q(g). Das bedeutet, dass der Ordinatenabschnitt des Differenzenquotienten für ln[T(g)] bekannt ist. Somit steht für alle Frequenzen bis g_N(1+1/γ)/2 der Differenzenquotient Δln[T(G)]/ΔG zur Verfügung (Figur 4d). Damit aber kann T(G) eindeutig durch Integration bestimmt werden.
Figur 5 zeigt ein Ausführungsbeispiel des erfindungsgemäßen Verfahrens zur Bestimmung der Übertragungsfunktion der CCD Kamera eines Transmissionselektronenmikroskops. In Figur 5a ist die elektronenmikroskopische Aufnahme einer handelsüblichen Testprobe zu sehen, die aus einem dünnen, auf einem Kupfernetz liegenden Kohlenstofffilm besteht. Die nominale Vergrößerung dieser als Bild 1 ("Image 1") bezeichneten Aufnahme beträgt 13000. Die zur Aufnahme benutzte CCD Kamera besteht aus 2048 x 2048 Pixeln der Größe 15 Mikrometer. In Figur 5b ist ein Ausschnitt einer zweiten Aufnahme ("Image 2") zu sehen, die vom selben Objektgebiet mit 10000-facher Nominalvergrößerung gemacht wurde. Dieser zu Bild 1 deckungsgleiche Ausschnitt besteht aus 1556 x 1556 Pixeln und wird als Bild 2 bezeichnet. Das genaue Vergrößerungsverhältnis zwischen den Bildern 1 und 2 beträgt γ=1.316. Das zweidimensionale Frequenzspektrum von Bild 1 wurde mittels des Algorithmus der schnellen Fouriertransformation (Fast Fourier Transform, FFT) auf 2048 x 2048 Pixeln ermittelt, das entsprechende Frequenzspektrum des Bildes 2 mittels einer mixed-radix FFT auf 1556 x 1556 Pixeln. Anschließend wurden die beiden zweidimensionalen Spektren zur Rauschminderung nach Gl. (4) azimuthal gemittelt. Zur weiteren Rauschminderung wurden für Bild 1 und für Bild 2 jeweils vier solcher unabhängig gewonnener Spektren aus jeweils vier Aufnahmen gemittelt. Die so erhaltenen eindimensionalen Spektren wurden anschließend gemäß Gl. (6) bezüglich einer Rauschkonstante korrigiert. Dieser Schritt ist notwendig, da ein geringer additiver Sockelbetrag der Intensitätspektren, der auf Rauschen zurückzuführen ist, bei hohen Raumfrequenzen dominant werden kann. Da dieser Sockelbetrag weder vom Objekt noch von der Übertragungsfunktion abhängt, lässt er sich auch nicht nach Gl. (2) behandeln und würde das Ergebnis bei hohen Raumfrequenzen systematisch verfälschen. Die um den Sockelbetrag reduzierten Intensitätsspektren I₁(g) und I₂(g) sind in Figur 5c dargestellt. Die in Figur 5 gewählte diskrete Indizierung der Frequenzachse entspricht der in Figur 3 gewählten Darstellung, wobei die Frequenz n=1024 die Nyquistfrequenz des aus 2048 x 2048 Pixeln bestehenden Bildes 1 bezeichnet.
In Figur 5d ist die nach Gl. (3) aus den beiden Intensitätskurven I₁(g) und I₂(g) gebildete Quotientenkurve Q(g) dargestellt. Aus der Quotientenkurve Q(g) wird danach unter Zuhilfenahme der Gln. (4)-(7) der Differenzenquotient D(G) gebildet. Die mittels der Interpolation nach Gl. (11) gewonnene Darstellung des Differenzenquotienten D_n = D(n), welche die ursprünglichen Metrik des Bildes 1 wieder herstellt, ist in Figur 5e dargestellt. Bei dieser Darstellung gehört ein mit dem Index n versehener Differenzenquotient D_n = D(n) exakt zu einer mit dem Index n versehenen Frequenz des Bildes 1. Die gesuchte Übertragungsfunktion T(k) der CCD Kamera, die sich durch die Summation nach Gl. (12) und die anschließende Exponentiation nach Gl. (13) ergibt, ist in Figur 5f dargestellt. Die offensichtlich glatte Kurvenform der erhaltenen Übertragungsfunktion T(k) ist nicht - wie man vielleicht meinen könnte - durch einen zusätzlichen Glättungsschritt an dieser Stelle entstanden, sondern ist ein intrinsischer glättender Effekt der Summation nach Gl. (12), bzw. der dazu analogen allgemeineren Integration nach Gl. (9).

Claims

Verfahren zur Bestimmung der Übertragungsfunktion eines signalverarbeitenden Systems aus mindestens zwei Repräsentationen I₁(x) und I₂(x) eines Objekts, die das System aus unterschiedlich skalierten vom Objekt herrührenden Eingangssignalen erzeugt hat, oder aus einer Repräsentation I₁(x) eines ersten Objekts und aus einer Repräsentation I₂(x) eines hierzu bis auf einen Skalenfaktor geometrisch identischen Objekts, wobei das Objekt durch seine Anwesenheit eine physikalische Messgröße orts- und/oder zeitabhängig verändert und hierdurch das Eingangssignal erzeugt, wobei das Verfahren die Übertragungsfunktion ohne Kenntnis eines nominellen Eingangssignals bestimmt indem
• die beiden Repräsentationen in einen Arbeitsraum transformiert werden, in dem sie jeweils als Produkt der Übertragungsfunktion mit dem unbekannten Eingangssignal darstellbar sind;
• Ausschnitte der beiden Repräsentationen I₁(x) und I₂(x), die sich auf den gleichen Bereich bzw. gleichwertige Bereiche des Objekts bzw. der Objekte beziehen und somit auf das gleiche unbekannte Eingangssignal zurückgehen, ausgewählt werden;

• diese Ausschnitte im Arbeitsraum als Funktionen I₁(g) und I₂(g) ausgedrückt werden;

• der Quotient Q(g)=I₂(g)/I₁(g) aus diesen beiden Funktionen gebildet wird, so dass das unbekannte Eingangssignal unterdrückt wird sowie Zähler und Nenner jeweils die gesuchte Übertragungsfunktion mit unterschiedlich skalierten Argumenten enthalten;

• aus dem Verlauf des Quotienten der gesuchte Verlauf der Übertragungsfunktion ausgewertet wird,
dadurch gekennzeichnet, dass
die Ausschnitte der beiden Repräsentationen ausgewählt werden, indem ein Ähnlichkeitsmaß optimiert wird, wobei Position und Ausdehnung des Ausschnitts aus der zweiten Repräsentation die freien Parameter der Optimierung sind.
Verfahren nach Anspruch 1, dadurch gekennzeichnet, dass die beiden Repräsentationen in den Frequenzraum als Arbeitsraum transformiert, insbesondere fouriertransformiert, werden.
Verfahren nach einem der Ansprüche 1 bis 2, dadurch gekennzeichnet, dass aus dem kontinuierlich bzw. an diskreten Stützstellen erklärten Logarithmus von Q(g) ein kontinuierlich bzw. an diskreten Stützstellen erklärter Differenzenquotient für den Logarithmus von T(g) gebildet wird, wobei dieser Differenzenquotient insbesondere integriert bzw. summiert werden kann.
Verfahren nach Anspruch 3, dadurch gekennzeichnet, dass für die Bildung des Differenzenquotienten für eine Position g im Arbeitsraum aus einer Mehrzahl von mindestens drei Repräsentationen zwei Ausschnitte aus verschiedenen Repräsentationen und damit zwei Funktionen I₁(g) und I₂(g) ausgewählt werden, die sich um einen vorgegebenen Skalierungsfaktor γ voneinander unterscheiden, wobei der Skalierungsfaktor γ insbesondere zwischen 1 und 2 gewählt werden kann.
Verfahren nach einem der Ansprüche 1 bis 4, dadurch gekennzeichnet, dass für die Übertragungsfunktion T(g) ein parametrisierter Ansatz gemacht und mit einem Optimierungsverfahren daraufhin optimiert wird, dass aus T(g) der Quotient Q(g) oder ein daraus abgeleiteter Funktionsverlauf bestmöglich reproduziert wird.
Verfahren nach Anspruch 5, dadurch gekennzeichnet, dass T(g) auf bestmögliche Übereinstimmung des Quotienten T(pg)/T(g) mit dem aus I₁(g) und I₂(g) erhaltenen Quotienten Q(g)=I₂(g)/I₁(g) optimiert wird, worin p ein Skalierungsfaktor ist.
Verfahren nach einem der Ansprüche 5 bis 6, dadurch gekennzeichnet, dass der vorbekannte oder bei der Anfertigung der Repräsentationen I₁(x) und I₂(x) eingestellte Skalierungsfaktor γ, um den sich diese Repräsentationen unterscheiden, als Skalierungsfaktor p gewählt wird.
Verfahren nach einem der Ansprüche 5 bis 7, dadurch gekennzeichnet, dass der Wert von p im Rahmen eines Optimierungsverfahrens bestimmt wird.
Verfahren nach einem der Ansprüche 1 bis 8, dadurch gekennzeichnet, dass vor der Transformation in den Arbeitsraum die Repräsentation mit der geringeren Ausdehnung im Realraum auf die Ausdehnung der anderen Repräsentation hochinterpoliert wird.
Verfahren nach einem der Ansprüche 1 bis 9, dadurch gekennzeichnet, dass die Repräsentationen im Arbeitsraum in Polarkoordinaten transformiert und in dieser Darstellung azimuthal gemittelt werden.
Verfahren nach Anspruch 10, dadurch gekennzeichnet, dass aus den azimuthal gemittelten Repräsentationen jeweils ein Rauschspektrum oder ein konstanter Rauschuntergrund herauskorrigiert wird.
Verfahren nach einem der Ansprüche 1 bis 11, dadurch gekennzeichnet, dass die beiden Repräsentationen I₁(x) und I₂(x) mit einem in diskrete Pixel unterteilten CCD-Sensor aufgenommene Bilder sind.
Verfahren zur Bestimmung der Übertragungsfunktion eines signalverarbeitenden Systems mit den Schritten:
• es werden mit dem System mindestens zwei unterschiedlich skalierte Repräsentationen I₁(x) und I₂(x) eines Objekts oder mindestens eine Repräsentation I₁(x) eines ersten Objekts und eine Repräsentation I₂(x) eines hierzu bis auf einen Skalenfaktor geometrisch identischen Objekts erzeugt;

• mit diesen Repräsentationen wird das Verfahren nach einem der Ansprüche 1 bis 12 durchgeführt.
Verfahren nach Anspruch 13, dadurch gekennzeichnet, dass die unterschiedliche Skalierung der Repräsentationen I₁(x) und I₂(x) durch eine Änderung der Skalierung des vom Objekt am Eingang des Systems erzeugten Eingangssignals O(x) eingestellt wird.
Verfahren nach Anspruch 14, dadurch gekennzeichnet, dass ein Elektronenmikroskop, insbesondere ein Transmissionselektronenmikroskop, als signalverarbeitendes System gewählt und die Vergrößerungsänderung zwischen den Repräsentationen I₁(x) und I₂(x) durch die Hintereinanderschaltung der Zwischenlinse und der Projektorlinse eingestellt wird.
Verfahren nach einem der Ansprüche 13 bis 15, dadurch gekennzeichnet, dass die unterschiedliche Skalierung der Repräsentationen I₁(x) und I₂(x) durch eine Änderung des räumlichen Abstands zwischen dem Objekt und dem signalverarbeitenden System eingestellt wird.
Verfahren nach einem der Ansprüche 13 bis 16, dadurch gekennzeichnet, dass mindestens eine der Repräsentationen I₁(x) und I₂(x) als Zusammenfassung von mindestens zwei Einzelrepräsentationen erzeugt wird.
Verfahren nach einem der Ansprüche 1 bis 17, dadurch gekennzeichnet, dass Rauschen als Objekt gewählt wird.