01 Mai

Lightmodulation: Gerchberg-Saxton-Algorithmus (GSA)


Hier will ich kurz das Vorgehen bei dem Phasensuchalgorithmus GSA erläutern. Ich habe den Algorithmus in Python übertragen und eine Diskretisierungsfunktion (Link zum Code) hinzugefügt. 

Der GSA wurde in dem wegweisenden Paper „A Practical Algorithm for the Determination of Phase from Image and Diffraction Plane Pictures“ von R. W. Gerchberg und W. O. Saxton 1972 in Optik veröffentlicht und wird wegen seiner Einfachheit heute angewandt.
Der Wikipedia-Artikel Gerchberg-Saxton algorithm braucht unbedingt eine deutsche Übersetzung.

Phasenverschub für die Herstellung eines homogenen quadratischen Amplitudenverlauf visualisiert [0,2\pi] \rightarrow[scharz,weiß]

In der Optik ist oft die Zielintensitätsverteilung entscheidend und der Phasenversatz \Delta\phi zur Herstellung der Intensitätsverteilung unbekannt. Iterative Lösungsverfahren, wie der GSA können in Verbindung mit der Eingangsintensitätsverteilung den benötigten Phasenversatz approximieren. Dafür werden die Funktionen der Intensitätsverteilungen \sqrt{I} auf diskrete Werte 2-dimensionalen Matrizen reduziert. Der GSA ist nicht auf 2-dimensionale Intensitätsverteilungen beschränkt. Der GSA berechnet aus der Intensitätsverteilung in der Beugungsebene und der Intensitätsverteilung in der Bildebene, den erforderlichen Phasenversatz.

Der GSA benötigt eine Eingangsintensitätsverteilung I_{ein}, wobei \phi die Eingangsphase der Eingangsintensität ist. Im Normalfall wird von einem zeitlich und räumlich kohärenten Strahl ausgegangen, \phi=const.

(1)   \begin{equation*} \sqrt{I_{ein}}=A_{ein} \cdot e^{i\phi}=A_{ein} \cdot e^{i \cdot const}  \end{equation*}

Die Ausgangsintensität I_{aus} beinhaltet die Zielintensitätsverteilung und den jeweiligen Phasenversatz. Bei einer Phasenmodulation soll primär eine vorgegebene Zielintensitätsverteilung erzeugt werden, ein Phasenversatz der Ausgangsintensität ist nicht gewollt. Daher wird die Phase der Ausgangsintensität \phi=const gesetzt.

(2)   \begin{equation*} \sqrt{I_{aus}}=A_{aus} \cdot e^{i\phi}=A_{aus} \cdot e^{i \cdot const} \end{equation*}

Die Eingangs- und Ausgangsintenstitätsverteilung wird für ein 2-dimensionales DOE in einer 2-dimensionalen Matrix dargestellt, wobei der Eintrag der diskreten räumlichen Anordnung entspricht.
Im ersten Schritt des Algorithmus wird die Eingangsintensität I_{ein} mit dem Phasenversatz \Delta \phi kombiniert, welcher später den Phasenversatz des DOE erzeugen soll.

Der Phasenversatz wird im ersten Zyklus aus Gleichung (1) entnommen:

(3)   \begin{equation*} A_{ein} \cdot e^{i \cdot 0} \rightarrow A_{ein} \cdot e^{i\phi} \end{equation*}

Diese modifizierte Ausgangsverteilung wird über die 2-dimensionale Fast-Fourier-Trans\-formation (FFT) auf die Ausgangsseite übersetzt.

(4)   \begin{equation*} [ A_{ein} \cdot e^{i\phi} ] \xrightarrow[]{FFT} [ \bar A_{aus} \cdot e^{i\bar\phi} ]  \end{equation*}

Die Intensität der fourier-transformierten Eingangsintensität | \bar A_{aus} \cdot e^{i\bar \phi} |^2 entspricht der Intensitätsverteilung des Strahls im Fernfeld nach dem Einbringen des Phasenversatz durch das DOE. Der Phasenversatz \bar\phi wird aus der transformierten Eingangsintensität genommen und mit der Wurzel der Zielintensitätsverteilung \sqrt{I_{aus}}=A_{aus} getauscht.

(5)   \begin{equation*} A_{aus} \cdot e^{i \cdot 0} \rightarrow A_{aus} \cdot e^{i\bar\phi} \end{equation*}

Die \bar\phi-modifizierte Zielintensitätsverteilung wird über die 2-dimensionale Inverse Fast-Fourier-Transformation zurück auf die Eingangsseite transformiert.

(6)   \begin{equation*} [ \bar A_{ein} \cdot e^{i\phi} ] \xleftarrow[]{IFFT} [ A_{aus} \cdot e^{i\bar\phi} ] \end{equation*}

Die neue Phase \phi aus der Inversen Fast-Fourier-Transformation stellt den korrigierten Phasenversatz \Delta \phi für die Frauenhoferbeugung dar und kann als Phasenversatz für das DOE verwendet werden. Die erste Iteration ist abgeschlossen.
Nach der ersten Iteration ist der Abbildungsfehler und das Signal-Rausch-Verhältnis noch groß, daher wird hier die Eingangsintensität I_{ein} wieder mit dem neuen Phasenversatz \Delta \phi kombiniert mit Gleichung (1) und erneut transformiert. Der Algorithmus wird erneut durchlaufen. Ein Abbruch des Algorithmus erfolgt, wenn die gewollte Güte des DOEs erreicht wurde oder wenn einer Anzahl N an Iterationen erreicht ist.

Amplitudenverlauf bei 20 Interationen ohne Diskretisierung


Die Ausgabe des GSA ist eine 2-dimensionale Matrix, welche die kontinuierlichen Phasenversätze \Delta \phi \in [0, 2\pi] des jeweiligen Bereiches beinhaltet. Ein kontinuierlicher Abtrag nach Maß des Phasenversatz kann technisch nur schwer umgesetzt werden. Daher wird in der Anwendung der Phasenversatz diskretisiert. Die Diskretisierung kann innerhalb des Algorithmus oder nach Abbruch durchgeführt werden und bestimmt die theoretische Effizienz maßgeblich. Der GSA diskretisiert damit 2-mal den Phasenversatz: Die Positionsmatrix der Phasenversätze hat eine begrenzte Anzahl an Einträgen, ein kontinuierlicher Verlauf ist daher hier nicht möglich. Zusätzlich werden die Werte der Einträge, also der Phasenversatz, diskretisiert, es kommt erneut zu einer Verringerung der Ausgabemöglichkeiten.

2 Gedanken zu „Lightmodulation: Gerchberg-Saxton-Algorithmus (GSA)

  1. Hey Dominik, danke für den tollen Artikel. Ich wundere mich ein wenig, dass die Wellenlänge des Eingangsstrahls keine Rolle im Algorithmus spielt. Beugt man eine Welle an einem Amplitudengitter ist das Beugungsbild von der Wellenlänge abhängig. Dies sollte doch auch bei einem Phasengitter der Fall sein?

    • Hallo Maximilian, vielen Dank für deinen Kommentar.
      Genau, die Größe des Beugungsbildes ist abhängig von der Wellenlänge. Die Wellenlänge $\lambda$ bestimmt zusammen mit der Breite der eingebrachten Phasen im DOE den Beugungswinkel (genauso wie beim Gitter). Damit meine ich die ‚Vergrößerung‘ der einzelnen Ordnungen des gebeugten Strahls.
      Das relative Intensitätsmuster (an dem wir interessiert sind) ändert sich bei unterschiedlichen Wellenlängen nicht.

      Monochromatisches kohärentes Licht (bspw. ein Laserstrahl) kann man lokal wie eine ebene Welle betrachten.
      Wenn man nun eine Störung (bspw. eine Blende oder ein DOE) in den Strahl stellt, wird das Licht -genauso wie bei der Wasserwelle- das Licht gebeugt.
      Allgemein ist das schwierig zu beschreiben, deshalb haben sich zwei Näherungen historisch durchgesetzt:
      Die Beugung, die nahe nach der Störung eintritt: Fresnel-Beugung.
      Und die Beugung, die weit entfernt von der Störung auftritt: Frauenhofer-Beugung.

      Letztere wird in dem obigen Artikel behandelt.
      Hier handelt es sich also um die Näherung, wie das Beugungsbild in weiter Ferne aussieht.
      Um das zu berechnen nähert man das fresnel-kirchhoffsche Beugungsintegral für lange Strecken und stellt fest, dass es -zufälligerweise- die Form einer Fouriertransformation hat.

      Im Beugungsintegral für die Fraunhofer-Näherung sieht man, wie die Wellenlänge das Resultat in der Fouriertransformation beeinflusst: Es wird die resultierende Funktion in der Bildebene nur gestreckt oder gestaucht. (bspw. http://www.physik.kit.edu/Studium/Lehramt/Fourieroptik_Vorbereitung_WS_21_22.pdf Gleichung 4.5)
      Bei dem Algorithmus wird dieser Faktor außer acht gelassen, das kommt u.a. davon, dass man eine diskrete Fouriertrafo verwendet.

Schreibe einen Kommentar

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

Diese Seite verwendet Cookies. Mit der weiteren Verwendung stimmst du dem zu. This website uses cookies. If further used, you agree to the policy.