1.

Introduction

The retinal nerve fiber layer (RNFL) exhibits birefringence that is due to the cylindrical orientation of ganglion cell axons consisting of microtubules. A decrease in the RNFL phase retardation is caused by a decrease in the number of microtubules of the RNFL.¹ The changes in the RNFL birefringence correlate with glaucomatous damage.^{2, 3, 4} The retardation distribution of the optic nerve fiber layer is important information for carrying out glaucoma diagnosis. Because the optic nerve fiber layer has birefringence, the retardation has been measured by polarization analysis using circularly polarized light.^{5, 6} GDx-VCC (Carl Zeiss Meditec), which uses the measurement technique of scanning laser polarimetry (SLP), is commercially available. Corneal birefringence is an important source of variance in SLP.^{3, 7} GDx-VCC compensates for the corneal birefringence by monitoring a bow-tie pattern in the macular region while using a Babinet-Soleil compensator.^{2, 3} Then it statistically diagnoses glaucoma using a normative database. It has been reported that the ratio between the retardation of the RNFL around the optic disk and the thickness is not constant in all positions.⁸ Nevertheless, the retardation map is very useful information for diagnostic purposes.

Recently, a new method using polarization-sensitive spectral-domain optical coherence tomography (PS-SD-OCT) has been proposed to measure the retardation of the RNFL.^{9, 10} In this method, all of the elements of the Jones matrix of the sample arm, including the cornea, are obtained. By using the matrix diagonalization method developed by Park, the effect of the birefringence of the single mode fiber and the cornea is compensated.¹¹ Finally, the phase retardation and relative orientation of the RNFL are calculated. Although PS-SD-OCT systems have been reportedly used in retinal imaging,^{9, 10, 11} they are not yet commercially available.

In this paper, we present the measurement results of phase retardation in normal eyes using a newly developed polarization-sensitive fundus camera method. Our imaging system, which is based on a polarization analysis method,¹² is capable of detecting the Stokes vector onto a fundus camera instead of the usual CCD camera. If the Stokes vector expressing the whole state of polarization can be detected by a one-shot imaging system in which the interference method is not used, we can numerically eliminate the influence of the background scattering and of the retardation caused by the cornea. Because the retardation process of the eye optics can be represented by the numerical equation using the retardation matrix of each component and also the nonpolarized background scattering light, it can be calculated using the Stokes vector. The present method enables the detection of the Stokes vector for the reflection from the macula region. Therefore, it is possible to numerically eliminate the influence of the retardation caused by the cornea.

2.

Method

2.1.

Measuring Principle

The principle of the present method has been previously described.¹³ In this method, the circularly polarized light enters into the eye and comes out as elliptically polarized light. The circularly polarized light passes through the cornea and RNFL and is then partially reflected by the retinal pigment epithelium. The reflected light again passes through the RNFL and the cornea and exits the eye. Because the optic nerve fiber layer and the cornea have birefringence, the circularly polarized light returned from the eye becomes elliptically polarized, while the reflection and backscattering also carry the nonpolarized light component. The polarization of the reflected light can be presented with a Stokes vector that has four parameters, S = (s ₀, s ₁, s ₂, s ₃)^T, where T is the transposed matrix. s ₀ is the averaged intensity, s ₁ is the linearly polarized light component on 0 and 90 deg axes, s ₂ is the linear polarized light component on 45 and 135 deg axes, and s ₃ is the circularly polarized light component. Generally, the nonpolarized light is dominant, and [TeX:] $s_0 > \sqrt {s_1^2 + s_2^2 + s_3^2 }$ $s_0>\sqrt{s_1^2+s_2^2+s_3^2}$ . In this study, we took into account s ₁, s ₂, and s ₃, while the scattered light was not included in the analysis.

The Stokes vector is transposed with the Mueller matrix defined by Eq. 1. The vector on the right side is the Stokes vector for an incident light ray that enters into the eye. The vector on the left side is the Stokes vector transposed through the eye optics

Eq. 1

[TeX:] \documentclass[12pt]{minimal}\begin{document}\begin{equation} \left(\begin{array}{l} s_0^{\prime} \\[5pt] s_1^{\prime} \\[5pt] s_2^{\prime} \\[5pt] s_3^1 \\[5pt] \end{array} \right) = M\left(\begin{array}{l} s_0 \\[5pt] s_1 \\[5pt] s_2 \\[5pt] s_3 \\[5pt] \end{array} \right),\quad M = \left[ T_{\theta} \right]\left[ C_{\alpha} \right]\left[ T_{ - \theta } \right]. \end{equation}\end{document} $$\left(\begin{array}{c}{s}_0^{\prime }\hfill \\ s_1^{\prime }\hfill \\ s_2^{\prime }\hfill \\ s_3^1\hfill \end{array}\right)=M\left(\begin{array}{c}{s}_0\hfill \\ s_1\hfill \\ s_2\hfill \\ s_3\hfill \end{array}\right),M=\left[T_{\theta }\right]\left[C_{\alpha }\right]\left[T_{-\theta }\right].$$

T _θ, and C _α are the matrixes given by

Eq. 2

[TeX:] \documentclass[12pt]{minimal}\begin{document}\begin{eqnarray} [ {T_\theta }] &=& \left[ {\begin{array}{c@\quad c@\quad c@\quad c} 1 & 0 & 0 & 0 \\ 0 & {\cos 2\theta } & { - \sin 2\theta } & 0 \\ 0 & {\sin 2\theta } & {\cos 2\theta } & 0 \\ 0 & 0 & 0 & 1 \\ \end{array}} \right]\nonumber\\ \left[ {C_\alpha } \right] &=& \left[ {\begin{array}{c@\quad c@\quad c@\quad c} 1 & 0 & 0 & 0 \\ 0 & 1 & 0 & 0 \\ 0 & 0 & {\cos \Gamma } & {\sin \Gamma } \\ 0 & 0 & { - \sin \Gamma } & {\cos \Gamma } \\ \end{array}} \right], \end{eqnarray}\end{document} $$\begin{array}{ccc}\hfill \left[T_{\theta }\right]& =& \left[\begin{array}{cccc}\hfill 1& \hfill 0& \hfill 0& 0\\ \hfill 0& \hfill \cos 2\theta & \hfill -\sin 2\theta & 0\\ \hfill 0& \hfill \sin 2\theta & \hfill \cos 2\theta & 0\\ \hfill 0& \hfill 0& \hfill 0& 1\end{array}\right]\hfill \\ \hfill \left[C_{\alpha }\right]& =& \left[\begin{array}{cccc}\hfill 1& \hfill 0& \hfill 0& 0\\ \hfill 0& \hfill 1& \hfill 0& 0\\ \hfill 0& \hfill 0& \hfill \cos \Gamma & \sin \Gamma \\ \hfill 0& \hfill 0& \hfill -\sin \Gamma & \cos \Gamma \end{array}\right],\hfill \end{array}$$

where θ is the angle between the fast axis and horizontal axis, Γ is the retardation that depends on the thickness (d) of the birefringence material, the index difference between the fast and slow axes, and the wave length λ

Eq. 3

[TeX:] \documentclass[12pt]{minimal}\begin{document}\begin{equation} \Gamma = \frac{{d(n_1 - n_2)}}{\lambda }2\pi. \end{equation}\end{document} $$\Gamma =\frac{d(n_1-n_2)}{\lambda }2\pi .$$

The optical model of the eye is shown in Fig. 1. M _c and M _r are the Mueller matrices of the cornea and the nerve fiber layer, respectively. If the polarized light reflected from the surface of the retinal pigment layer retains its polarization, the Stokes vector of the light from the eye optics may be shown as follows:

Eq. 4

[TeX:] \documentclass[12pt]{minimal}\begin{document}\begin{equation} S = M_c M_r M_r M_c A, \end{equation}\end{document} $$S=M_c{M}_r{M}_r{M}_{c}A,$$

where S is the detected Stokes vector and A is the Stokes vector of the circularly polarized light that enters the eye.

Fig. 1

Optical model of the eye. The circularly polarized light enters into the eye and comes out as elliptically polarized light. The circularly polarized light passes through the cornea and RNFL and is partially reflected by the retinal pigment epithelium. M _c and M _r are the Mueller matrices of the cornea and the nerve fiber layer. The polarized light reflected from the surface of the retinal pigment layer retains its polarization.

Provided that M _c can be measured in advance, we can multiply M _c ⁻¹ to both sides of Eq. 4 to yield

Eq. 5

[TeX:] \documentclass[12pt]{minimal}\begin{document}\begin{equation} S^{\prime} = M_R B, \end{equation}\end{document} $$S^{\prime }=M_{R}B,$$

where [TeX:] $S^\prime = M_c^{ - 1} S$ $S^{\prime }=M_c^{-1}S$ , B = M _c A, and M _R = M _r M _r. Since Stokes vectors S ^′ and B are known, Mueller matrix M _R can be derived using Eq. 5. There are two parameters, θ_R and Γ_R, in the matrix M _R. θ_R is the axis of the retardation caused by the nerve fiber layer, and Γ_R is the magnitude of the retardation. By performing an optimization to derive the two parameters, we may evaluate

Eq. 6

[TeX:] \documentclass[12pt]{minimal}\begin{document}\begin{equation} f(\theta,\Gamma) = \left\| {S^{\prime} - M(\theta,\Gamma)B} \right\|, \end{equation}\end{document} $$f(\theta ,\Gamma )=‖S^{\prime }-M(\theta ,\Gamma )B‖,$$

and then find the elements θ₀ and Γ₀ such thatf(θ₀, Γ₀) ⩽ f(θ, Γ). θ₀ and Γ₀ are regarded as the optimal solution of θ_R and Γ_R, and θ and Γ cover the space of the candidate solutions.

2.2.

Retardation Measurement

Our polarization analysis system is constructed with a CCD area image sensor, a linear polarizing plate, and a microphase plates array. As illustrated in Fig. 2, the retardations of microphase plates are the same, but the fast axis angles are different. There are four different axes microphase plates, and the polarization is analyzed with four pixels of data. In the present measurement, the angles and the retardation were determined in accordance with the study of Sabatke ¹⁴

Fig. 2

Construction of polarization analysis camera. The polarization analysis system is constructed with a CCD area image sensor, a linear polarizing plate, and a microphase plates array. The retardations of microphase plates are the same, while the fast axis angles are different. There are four different axes microphase plates, and the polarization is analyzed with four pixels of data.

The four fast axes angles of the plates are +51.7, −51.7, +15.7, and −15.7 deg, respectively. Meanwhile, the retardation for the microphase plates is 132 deg. If we define the retardation of wavelength plates as Δ and their fast axis angles as ϕ₁, ϕ₂, ϕ₃, and ϕ₄,theStokesvectorS = (S ₀, S ₁, S ₂, S ₃)^Tmay

be related to the intensities of the four pixels (I ₁, I ₂, I ₃, I ₄)^T by Eq. 7, and they can be derived with an inverse matrix calculation

Eq. 7

[TeX:] \documentclass[12pt]{minimal}\begin{document}\begin{equation} \left[ {\begin{array}{c} {I_1 } \\ {I_2 } \\ {I_3 } \\ {I_4 } \\ \end{array}} \right] = \frac{1}{2}\left[ {\begin{array}{c@\quad c@\quad c@\quad c} 1 & {1 - (1 - \cos \Delta)\sin ^2 2\phi _1 } & {(1 - \cos \Delta)\sin 2\phi _1 \cos 2\phi _1 } & { - \sin \Delta \sin 2\phi _1 } \\ 1 & {1 - (1 - \cos \Delta)\sin ^2 2\phi _2 } & {(1 - \cos \Delta)\sin 2\phi _2 \cos 2\phi _2 } & { - \sin \Delta \sin 2\phi _2 } \\ 1 & {1 - (1 - \cos \Delta)\sin ^2 2\phi _3 } & {(1 - \cos \Delta)\sin 2\phi _3 \cos 2\phi _3 } & { - \sin \Delta \sin 2\phi _3 } \\ 1 & {1 - (1 - \cos \Delta)\sin ^2 2\phi _4 } & {(1 - \cos \Delta)\sin 2\phi _4 \cos 2\phi _4 } & { - \sin \Delta \sin 2\phi _4 } \\ \end{array}} \right]\left[ {\begin{array}{c} {S_0 } \\ {S_1 } \\ {S_2 } \\ {S_3 } \\ \end{array}} \right]. \end{equation}\end{document} $$\left[\begin{array}{c}{I}_1\\ I_2\\ I_3\\ I_4\end{array}\right]=\frac{1}{2}\left[\begin{array}{cccc}\hfill 1& \hfill 1-(1-\cos \Delta ){\sin }^{2}2{\phi }_1& \hfill (1-\cos \Delta )\sin 2{\phi }_1\cos 2{\phi }_1& -\sin \Delta \sin 2{\phi }_1\\ \hfill 1& \hfill 1-(1-\cos \Delta ){\sin }^{2}2{\phi }_2& \hfill (1-\cos \Delta )\sin 2{\phi }_2\cos 2{\phi }_2& -\sin \Delta \sin 2{\phi }_2\\ \hfill 1& \hfill 1-(1-\cos \Delta ){\sin }^{2}2{\phi }_3& \hfill (1-\cos \Delta )\sin 2{\phi }_3\cos 2{\phi }_3& -\sin \Delta \sin 2{\phi }_3\\ \hfill 1& \hfill 1-(1-\cos \Delta ){\sin }^{2}2{\phi }_4& \hfill (1-\cos \Delta )\sin 2{\phi }_4\cos 2{\phi }_4& -\sin \Delta \sin 2{\phi }_4\end{array}\right]\left[\begin{array}{c}{S}_0\\ S_1\\ S_2\\ S_3\end{array}\right].$$

2.3.

Method to Negate the Influence Caused by the Birefringence of the Cornea

Our method is provided with the condition that there is very small polarization effect at the macular region. We presuppose that at the macular region, the optic nerve fiber layer is very thin, little retardation is caused by the NFL, the fiber layer of Henle that has a radial structure, and the thickness is the same on all meridians.

Henle's fiber layer causes retardation of the slow axis, which is in the radial direction. In the direction that is the same as the slow axis of the retardation caused by the cornea, the retardation is increased, and in the direction that is the same as the fast axis of cornea, the retardation is decreased.

From the difference between the decreased and increased retardation, we may calculate the retardation caused by the cornea Γ_c, from the maximum retardation Γ_max, and the minimum retardation Γ_min as follows:

Eq. 8

[TeX:] \documentclass[12pt]{minimal}\begin{document}\begin{equation} \Gamma _c = \textstyle{\frac{1}{2}}(\Gamma _{\max } \pm \Gamma _{\min }). \end{equation}\end{document} $${\Gamma }_c=\frac{1}{2}({\Gamma }_{\max }\pm {\Gamma }_{\min }).$$

As Γ_c is generally not significantly large, we choose to the minus sign in Eq. 8. As illustrated in Fig. 3, the illuminating light of the fundus camera passes through the ring aperture peripheral part of the cornea and the CCD camera captures the light that passes through the center of the cornea. The cornea exhibits different values of retardation at different positions.¹⁵ Γ_c is two times the magnitude of the mean retardation of the cornea Γ_cHole with pinhole aperture for observation and the cornea Γ_cRing with ring aperture for illumination. It has been reportedly estimated that the cornea ring aperture and cornea pinhole aperture have nearly the same axis of retardation.¹⁵ Here, we set Γ_cRing and Γ_cHole to half of Γ_c.

Fig. 3

Illuminating light of the fundus camera passes through the ring aperture peripheral part of the cornea and the CCD camera captures the light that passes through the center of the cornea.

This method can calculate only the polarized portion.¹⁴ However, the total reflected light may contain those from subretinal layers, not just those from the fiber layer. Such effects may induce errors in the calculation and may be evaluated from the measured Stokes vector. We calculated with the assumption that such an effect does not exist. So, the result of retardation calculated from the Stokes vector may potentially include errors.

3.

Experiment

3.1.

Apparatus

The imaging system employed a polarization analysis camera designed with a central wavelength of 550 nm. To selectively use a portion of the spectrum of the broadband light source, which is a halogen lamp, a green pass filter centered at the wavelength of 550 nm and with a bandwidth of 30 nm was inserted in the illuminating path. A 1/4-wave plate and a polarization filter were also placed in the illuminating path for circularly polarized light illumination. Finally, we applied a polarization analysis camera instead of the usual CCD camera of the fundus camera. A photo of the prototype is shown in Fig. 4.

Fig. 4

A photo showing the outlook of the prototype. We applied a polarization analysis camera instead of the usual CCD camera of the fundus camera. A 1/4 wave plate and a polarization filter were also placed in the illuminating path for circularly polarized light illumination.

4.

Results

One of images taken with this system is shown in Fig. 5a. The striped pattern seen in the image is an interference pattern induced by the space between the phase, polarizing, and CCD plates. There is no influence that these stripes analyze with calibration beforehand by using the linear polarization and the circularly polarized light every four pixels.

The Stokes parameters S1, S2, and S3 were calculated from Fig. 5a, and their distribution maps are depicted in Figs. 5b, 5c, 5d, respectively, where the gray level corresponds to the value ranging from −1 to 1.

As we explained above, these three images do not show the interference pattern. This result provides evidence that the transformations from intensity to Stokes vector were performed well. In Fig. 5, the values of S1, S2, and S3 are displayed with the pixel intensity of values ranging from 0 to 255, where the original value of S1, S2, and S3 were between −1 and +1. The detected light was always partially polarized and the values of S0 were from 0.8 to 2.0 and the average of S0 was 1.2.

Fig. 5

Data sample of right eye of a 29-year old woman. (a) Image taken with polarization analysis camera. Distribution maps of Stokes parameter (b) S1, (c) S2, and (d) S3. The gray scale in the distribution maps corresponds to −1 to 1.

The distribution maps of the retardation taken with this system are shown in Fig. 6. The value of retardation is illustrated with color. The retardation of 0 to 90 nm is assigned to the color of dark blue to red as shown by the color bar in Fig. 6. The retardation in Fig. 6a is caused by the nerve fiber layer and the cornea. We are aware that there is almost no nerve fiber layer at the macula region marked as macular. Therefore, we may consider that the retardation shown at macular is caused by the cornea and the fiber layer of Henle with thin radial distribution and equal thickness. If we subtract the retardation of the cornea from the entire retardation [Fig. 6a] by numerical calculation, the result shown in Fig. 6b is the retardation caused only by the nerve fiber layer. The macula region of Fig. 6b is all blue, indicating that the compensation was probably correct.

Fig. 6

Distribution maps of retardation of the nerve fiber layer of a 29-year old woman. The color bar corresponds to the retardation in a range of 0 to 90 nm. (a) Distribution of retardation caused by the nerve fiber layer and the cornea. Macular: The retardation of the macula area is due to the cornea. (b) Distribution of retardation that was without the retardation caused by the cornea. The retardation of the cornea was subtracted from the entire retardation (a) by numerical calculation. The result shown in (b) is the retardation solely caused by the nerve fiber layer. Optic disk: The area between the two circles surrounding the optic disk (indicated by the arrow) is the area under measurement.

Clinicians commonly collect data for the region surrounding the optic disk. The place to be measured is marked optic disk in Fig. 6. The area between the two circles is the area to be measured. The average of the diameters of the two circles on the retina is approximately 3.4 mm. The result of the measured retardation of the above part is shown in Fig. 7. The total processing time with a personal computer, Intel Core™2 Duo (2.93GHz), was around 20 s per measurement, including the operation time on the software that provides retardation.

Fig. 7

Measured result of retardation for the area between the two circles surrounding the optic disk shown in Fig. 6b. A double-humped shape is observed.

For the second result, the CCD acquired an image of a 23-year old woman's left eye. The retardation images without and with corneal compensation are shown in the left and right sides of Fig. 8, respectively. The macula region of Fig. 8b is not all blue. The compensation on macula region is worse than that of the first result. The numerical data of retardation with corneal compensation is shown in Fig. 9.

Fig. 8

Analyzed result of a 23-year old woman's left eye. Result (a) without and (b) with correction. The cornea compensation by the proposed method was not perfect with this second eye. Note that the surroundings of the fovea will reach an almost constant low value shown in blue color if the influence from the cornea is perfectly compensated.

Fig. 9

Measured result of retardation for the area between the two circles surrounding the optic disk shown in Fig. 8b.

The images shown in Fig. 8 indicate that the cornea compensation by the proposed method was not perfect with the second eye. The surroundings of the fovea should reach an almost constant low value shown in blue color if the compensation was perfectly done. However, the actual result is different. From the inspection of the macula region after compensation and the comparison of color variation in the color bars in Fig. 8b, we estimate that this method has an error of 10 nm or less in retardation.

5.

Discussion

The calculated retardation with the proposed system shows a double hump distribution, which is similar to the result reported by Huang.¹⁶ Our system provided a similar result in the measurement of the retardation caused by the RNFL. This may offer some evidences that the present system is feasible for the measurement of the retardation primarily caused by the RNFL in the optic disk region.

Our system offers a unique advantage that the image, which covers the macula and optic disk regions, can be acquired in one shot. This feature enables us to calculate the Stokes vector of every point from the image. Meanwhile, we may compensate the retardation error caused by the cornea at the optic disk region by compensating with the retardation calculated at the macula region. The glaucoma community is looking for systems that can provide reliable numbers on glaucoma detection and glaucoma progression. The changes in the RNFL birefringence may correlate with damage in glaucoma.^{1, 2, 3} Our method has the possibility to provide nerve fiber layer retardation measurements with a simply constructed apparatus.

The fundus camera has a ring aperture for fundus illumination. We calculated the compensation with the assumption that there is no difference between the average of the birefringence characteristic of the ring aperture region of the cornea and that of the pinhole region of the cornea, where the pinhole is placed in the center of the ring aperture. However, the birefringence characteristic of the cornea is not perfectly uniform,¹⁵ and this may become a cause of error. A separate measurement of the birefringence axis and retardation of the central pinhole region of the cornea may enable a more accurate calculation of the retardation by the RNFL.

6.

Conclusion

A measurement technique using a polarization-sensitive fundus camera has been developed for the measurement of phase retardation in the RFNL. The present technique is capable of detecting the Stokes vector expressing the whole state of polarization in single-shot imaging, making it feasible to numerically eliminate the influence of the background scattering and of the retardation caused by the cornea. Our measurement results in normal eyes demonstrated that we could calculate the retardation mainly caused by the optic nerve fiber layer. This process is confirmed by the numerical difference between the measurements on the macula region, including the Henle's fiber layer and the region surrounding the optic disk. It is our hope that the current method can provide a simple yet efficient way to measure the retardation in the RFNL, which may be useful in Glaucoma detection and monitoring.