Open Access
Add to BibSonomyAbstract
We propose a novel simulation method based on rendering to evaluate the modulation transfer function (MTF) of optical imaging systems. The new simulation method corresponds to an experimental measurement of the MTF using imaging resolution test charts, and therefore allows an analysis of the resolving power of shift-variant optical systems that are difficult to evaluate with conventional methods based on the point spread function (PSF). Furthermore, the effects of stray light, such as from reflection or scattering, on the imaging performance can be analyzed. In contrast to methods based on illumination optics using Monte Carlo methods, the proposed method calculates the intensities on an image surface with rendering techniques used in three-dimensional computer graphics (3D CG), which results in calculations that are faster and have a higher precision. The proposed method is highly effective in analyzing the MTF of optical imaging systems through simulations.
© 2013 Optical Society of America
1. Introduction
Measurement of the modulation transfer function (MTF), which is the absolute value of the optical transfer function (OTF), is a method to quantitatively evaluate optical imaging systems. The MTF is an index of the resolving power for each spatial frequency and is used to evaluate optical systems together with aberration charts. The MTF of an optical imaging system is typically evaluated by imaging a test chart that has a contrast variation with a period equal to that of the designated spatial frequency and then measuring to what extent the contrast of the original test chart is reproduced. The MTF is often derived from the point spread function (PSF) in simulations. This assumes that the PSF is independent of the light source, and it also requires that the PSF is shift invariant; hence it does not change between points typically less than 1 mm apart when compared to the results from test chart measurements. Correct evaluation is impossible when the PSF changes significantly in the test chart interval because large fluctuations appear in the PSF-derived MTF. The obtained MTF is also affected by significant changes in the radiance distribution peaks when test charts are imaged in such optical systems.
Yoshida et al. (including the author) [1] have proposed a simulation technique based on illumination methods to address these problems. This method generates light rays from bright areas in the test chart using Monte Carlo methods and obtains the irradiance from the number of rays arriving at the image surface after passing through the optical imaging system. As the method is based on measurements that image test charts, interactions in shift-variant systems are reflected in the radiance distribution. As a consequence, it has become possible to analyze systems in which the simulated and experimental measurements correspond. However, one drawback of this method is the variance from the use of Monte Carlo methods that appears as high-frequency noise. Because a large number of light rays is necessary in the simulations, a satisfactory accuracy could not be obtained with the conventional method despite long computation times. Therefore, in this paper, rendering techniques used in three-dimensional computer graphics (3D CG) are applied to improve both the speed and accuracy. Rendering is a method to calculate the radiance at each pixel in an image, and ray tracing is typically carried out from the viewpoint (image surface), which is possible because of radiance invariance. While two-dimensional calculations are necessary to obtain the irradiance when a test chart is imaged, the radiance distribution necessary to derive the MTF can be obtained from a one-dimensional calculation when the radiance is calculated from the viewpoint. Therefore, highly efficient simulations as compared to methods based on illumination optics are possible. This paper outlines the formalism, simulation procedure, and algorithm of the new method.
2. Conventional methods
As examples of conventional simulation methods to analyze the MTF, we describe here the method using the PSF [2] and a method based on illumination optics [1].
2.1. Obtaining the transfer function using the PSF
The OTF is a quantitative measure of the reproduction capability of an optical system, and the MTF of an optical imaging system can be derived from the OTF. The OTF is the Fourier transform of the intensity distribution on the image surface normalized by the total amount of light. The OTF under incoherent illumination is defined as follows:
Here, I is the intensity distribution in the u–v plane, and fs and ft are the spatial frequencies in the sagittal and tangential directions, respectively. The absolute value of the OTF is the MTF, which can be written in terms of the PSF as In simulations of geometrical optics, rays of light from a point source are traced to the image plane, and the spot diagram in the image plane is used to calculate the MTF.2.2. Method based on illumination optics
Methods based on illumination optics correspond to experiments that measure the MTF from imaging test charts. These methods can analyze shift-variant optical systems and effects from stray light that cannot be analyzed with PSFs, which are delta functions. Figure 1 shows a schematic of the test chart and optical system arrangement.
Fig. 1 Simulation method based on illumination optics. The number of light rays emitted from a test chart illuminated by a light source is determined by Monte Carlo methods. The irradiance is obtained from the number of light rays that arrive at the image surface via the optical imaging system, and the MTF is derived from the irradiance distribution.
Incoherent illumination from a light source is reproduced based on illumination optics theory, and the irradiance distribution on the image surface is obtained by ray tracing using Monte Carlo methods. The MTF of an optical imaging system can ultimately be derived from the irradiance distribution using the equation
To be precise, a square-wave response function (SWRF) is obtained instead of the MTF when a square-wave test chart is used. Details are given in section 3.2.Illumination optics methods have been devised for an MTF analysis of the rod lens, which is an array of gradient index (GRIN) lenses. The rod lenses used in image scanners have small diameters and a periodic arrangement, and thus have significant changes in the PSF within 1 mm, which is the size of the test chart. Therefore, the MTF cannot be analyzed from the PSF because there are large differences in the intensity peaks at the image surface. Simulations based on illumination optics have been related to experimental systems, and information about the image of the test chart is useful in an evaluation of the radiance distribution. However, variance is a problem in Monte Carlo integration because the ray tracing uses Monte Carlo methods. Using too few rays increases the variance and affects the results as high-frequency noise. The number of rays needs to be increased four-fold to halve the variance, meaning that a vast number of rays are necessary in the simulations. In summary, such methods consume significant computational time and the accuracy is lowered because of noise.
3. New method using rendering
We propose a new method that uses rendering to increase the speed and accuracy compared with methods based on illumination optics. Rendering is a method to draw 3D CG images based on 3D objects and light sources given as numerical data. Physics-based rendering uses ray tracing and obtains the radiance by beginning the ray tracing from the viewpoint instead of deriving the irradiance from the light source. Figure 2 shows the new simulation method for test charts in optical systems.
Fig. 2 Simulation method using rendering. The radiance of pixels is obtained by ray tracing from the light receiving surface instead of the light source. One-dimensional ray tracing for each pixel on the light receiving surface is used to obtain the one-dimensional radiance distribution necessary to calculate the MTF.
The biggest benefit of ray tracing from the light receiving surface is that the radiance distribution can be obtained by one-dimensional ray tracing. Methods using illumination optics need to take into account all rays that arrive at the image surface, meaning that rays from all possible positions and directions need to be randomly emitted from the test-chart surface, and hence rays have to be generated two-dimensionally. Methods using rendering carry out ray tracing in the reverse direction, i.e., from the image surface, and so the MTF can be obtained using a one-dimensional calculation, which results in a significant reduction in the computational order. In other words, accuracy is increased by increasing the number of rays and using a slightly longer computational time. Two-dimensional ray tracing can obtain the image on the imaging surface and is more efficient than using illumination optics. This is because light rays that have been vignetted by the optical system are integrated as the radiance.
3.1. Formalism using the rendering equation
Problems in 3D CG rendering are reduced to the issue of how to solve the rendering equation. The rendering equation as proposed by Kajiya [3] shows the equilibrium state of light transmission in a model without participating media and is given as follows.
Here, Lo, Le, and Lr are the outgoing, emitted, and reflected radiances, respectively, p is the spatial coordinate, and ωo is a vector in the outgoing direction. The reflected radiance can be described as the integral of the incident radiance on a hemisphere Ω at a point p: Here, Li is the incident radiance, fr is the bidirectional reflectance distribution function (BRDF), ωi is a vector in the incident direction, and n is the normal vector. Many methods for solving the rendering equation have been proposed. Representative methods are the path tracing [3], radiosity [4], and photon mapping [5] methods. The difference in these methods lies in how to solve the integral in Eq. (5), for instance, by solving the equation probabilistically as a path integral using Monte Carlo methods or by using finite element methods. Global illumination can be achieved, and indirect illumination effects can be reproduced in all these methods, although there are advantages and disadvantages to each method. The optical system applied to the test chart in here illuminates the test chart with a light source and then detects the intensity on the image surface of the optical imaging system. Only the primary light from the light source that illuminates the test chart is relevant, and therefore indirect light does not need to be taken into account. Rendering considering complicated global illumination is thus unnecessary, and a simple rendering of local illumination by primary light is sufficient.Rendering of the primary light only was proposed by Whitted [6] in 1980 as a recursive ray tracing. Recursive ray tracing tracks light rays that recursively reflect and refract, and the integral term in Eq. (5) becomes
Here, V is the visibility function, which is 1 if the light source is visible from the observation point and 0 otherwise. Whitted’s recursive ray tracing method in Eq. (6) can be simplified further when applied to optical systems with test charts. If the test chart is uniformly irradiated, then the incident radiance Li is constant and there is no irregularity in the radiance. If the visibility function V is defined as 1 on the test chart and 0 (vignetted) elsewhere, then Li is always constant, and therefore If the illumination from the light source is considered as Lambertian diffuse lighting model, the BRDF is 1 regardless of the direction vector, and hence Direct light from the light source does not exist in the optical system for test charts, and so the emitted radiance Le does not need to be considered. This reduces the rendering equation in optical systems with test charts to Eq. (8). The incident radiance Li can change with the radiant flux and distance to the test chart; however Li can be set to 1 when only the contrast ratio is relevant for obtaining the MTF. Furthermore, taking the direction normal to the test chart as the z direction (n = {0,0,1}) means that the necessary reflected radiance is determined only by the z component of the incident light rays arriving at the test chart. In summary, the rendering equation becomes very simple in optical systems for test charts.3.2. Test chart
Rendering in 3D CG uses a 3D space for objects to be drawn called a scene. The scene for test-chart optical systems considered in this research simply becomes a two-dimensional plane perpendicular to the optical axis. Test charts with contrasts corresponding to different spatial frequencies are used. Two test charts are considered here. Sine-wave charts, where the contrast intensity distribution changes continuously in a sinusoidal manner, are used to obtain the MTF. Square-wave charts with a discontinuous contrast distribution, however, are often used because sine-wave charts are difficult to make experimentally. The SWRF is obtained in measurements using square-wave charts. Application of the Coltman correction [7] to a number of SWRFs at different spatial frequencies is necessary to obtain the MTF from SWRFs. While simulations can freely choose between square waves and sine waves, a sine-wave chart can simply be chosen to obtain the MTF, but the simulation of square-wave charts is necessary for a comparison with experimentally measured data. The contrast ratio from a square-wave chart is referred to as the SWRF to avoid confusion with the MTF, which is the contrast ratio from a sine-wave chart.
The contrast intensity distribution in the x direction of a sine-wave chart is given as follows.
Here, f is the spatial frequency and ϕ is the phase difference (brightness at each position). The square-wave chart can be denoted as follows in piecewise form. Here, T is the period (T = 1/f), and a is the position. The periodicity condition O(x + T) = O(x) is satisfied.3.3. Sampling and simulation procedure
Our model determines the origin of the light ray and the direction of emission based on sampling the pixels on the film surface and at the rear-most element of the optical imaging system. The sampling point for general lenses such as camera lenses becomes the exit pupil of the lens. This corresponds to a simulation of the depth of field in 3D CG rendering. There are many sampling methods such as stratified sampling [8] and the quasi-Monte Carlo method using a low-discrepancy sequence (LDS) [9]. These sampling methods are known to converge faster in general compared to using pseudo-random numbers, and hence simulation using an LDS is selected for this paper. If the optical imaging system contains a camera lens, then the exit pupil has a circular aperture, and it is necessary to convert the distribution for efficient sampling. Shirley [10] proposed a mapping from concentric squares to concentric circles. Kolb et al.[11] covered 3D CG rendering through a camera lens, its description in a simplified model, and discussed sampling in detail. Appropriate filters must be used when considering the properties of the charge-coupled device (CCD) or film used in experimental optical systems. Various filters have been proposed to model the film properties [8, 12]; however, the effects of the filter are not thought to affect the results significantly in the test-chart optical systems in this paper.
Figure 3 shows a block diagram of the entire simulation procedure including the sampling. The optical coordinate system and the algorithm are given in Figs. 4 and 5, respectively. The algorithm is much simpler than standard 3D CG rendering.
Fig. 3 Block diagram of the MTF analysis simulation using rendering. The test chart (position, spatial frequency, type) and film (size, position, number of pixels) are set and then the pixels and rear-most element of the optical imaging systems are sampled. The light rays generated by the sampling are traced, and the radiance distribution is obtained from the rays that arrive at the test chart.
4. Results
The performance analysis was carried out for an optical imaging system with a Tessar camera lens and another systems with a rod lens to demonstrate the MTF analysis using rendering. The simulation results for the Tessar lens are presented first. Figure 6 shows the specifications of the Tessar lens with data taken from reference [13]. The back focal length (distance from the exit pupil to the film surface) is 191.79 mm, and the working distance (distance from the entrance pupil to the test chart) is 191.29 mm. The test chart and film are aligned such that the center of both lies on the optical axis, and the size of the film is changed according to the spatial frequency of the test chart. The number of pixels on the film is fixed at 256 × 256, and the image is composed of a 256-step grayscale with the highest step set at the maximum radiance. The wavelength is set to the d-line value of 587.6 nm. A total of 2,000 samples per pixel (spp) are used for the two-dimensional image generation and 10,000 spp for the one-dimensional radiance distribution generation to increase the accuracy. Figure 7 shows the simulation results. The PSF does not change over the test-chart area for camera lenses such as the Tessar lens, and thus the radiance distribution peaks match. The MTF obtained from the radiance distribution agrees with that obtained from the PSF. However, lens-specific properties may affect the radiance distribution even in the case of camera lenses if the reflection of light by the lens material or stray light from absorption is taken into account.
Fig. 6 Specifications of the Tessar lens. Quantities nd and νd are the refractive index and Abbe number at the d-line (587.6 nm), respectively, and ap is the aperture size. Units are in mm.
Fig. 7 Simulated images and radiance distributions of square-wave and sine-wave charts at various spatial frequencies. The spatial frequencies are (a) 6 lp/mm, (b) 30 lp/mm, and (c) 50 lp/mm. For the radiance distribution, the solid and dashed lines are those of the square-wave and sine-wave charts, respectively.
We now present the simulation results for a rod lens in the optical imaging system. Rod lenses can be used as line sensors such as the read elements in image scanners. An MTF evaluation of rod lenses is effective and widely carried out. The rod lens for the simulation consists of 11 GRIN lenses with the refractive index distribution shown in Fig. 8; lens data from reference [14] are used. The settings for the test chart and film are the same as those for the Tessar lens simulations, and the sampling was 10,000 spp for the two-dimensional imaging and 150,000 spp for the one-dimensional radiance distribution. Figure 9 shows the simulation results.
Fig. 9 Simulated images and radiance distributions of the square-wave and sine-wave charts using a rod lens. The spatial frequency is 6 lp/mm. For the radiance distribution, the solid and dashed lines are those of the square-wave and sine-wave charts, respectively.
Because the rod lens in Fig. 8 has a smaller aperture than the test chart, the refractive index distribution is distorted and the resolution is low. The peaks in the radiance distributions in Fig. 9 are therefore different in some positions. Such optical systems cannot be evaluated correctly because the value of the MTF differs significantly in some positions. I recommend that radiance distributions with different peaks should be evaluated in terms of radiance distribution properties instead of using a single parameter, that is, the MTF.
Figure 10 shows comparative results between the conventional method based on illumination optics and the proposed method using rendering. The number of sampling is identical in both methods, and the one-dimensional radiance distribution is assumed as 150,000 spp × 256 pixels. Light rays were generated from bright areas in the square-wave chart, and the slit height of the chart along the y axis is 0.6 mm.
Fig. 10 Radiance distributions of the square-wave charts through a rod lens. The spatial frequency is 6 lp/mm. Solid and dashed lines are results of the conventional method based on illumination optics and the proposed method using rendering, respectively.
In the conventional method, high-frequency noise appears on the radiance distribution due to the variance of the Monte-Carlo method. These results indicate accuracy of the proposed method is higher than the conventional method in the same sampling number.
5. Conclusion
This paper has proposed a new method to simulate the MTF using rendering. Calculating the radiance in optical systems with test charts allows for an analysis of the MTF and radiance distribution for optical systems that are difficult to analyze using the PSF. The field of 3D CG rendering is making significant advances; however, many of the core technologies are from conventional optics, and photorealistic CG are made though physics-based simulations such as ray tracing. This research takes the opposite approach where rendering, which was developed in the field of CG, is applied to the analysis of imaging optics. This method is pioneering work in optical analysis using rendering.
References and links
1. S. Yoshida, S. Horiuchi, Z. Ushiyama, and M. Yamamoto, “A numerical analysis method for evaluating rod lenses using the Monte Carlo method,” Opt. Express 18, 27016–27027 (2010) [CrossRef] .
2. J. W. Goodman, Introduction to Fourier Optics, 3rd ed. (Roberts and Company Publishers, 2005).
3. J. T. Kajiya, “The rendering equation,” SIGGRAPH 1986 , 143–150 (1986) [CrossRef] .
4. T. Nishita, I. Okamura, and E. Nakamae, “Shading models for point and linear sources,” ACM Trans. Graph. 4, 124–146 (1985) [CrossRef] .
5. H. W. Jensen, Realistic image synthesis using photon mapping (AK Peters, Ltd., 2001).
6. T. Whitted, “An improved illumination model for shaded display,” Communications of the ACM 23, 343–349 (1980) [CrossRef] .
7. J. W. Coltman, “The Specification of Imaging Properties by Response to a Sine Wave Input,” J. Opt. Soc. Am. 44, 468–469 (1954) [CrossRef] .
8. M. Pharr and G. Humphreys, Physically based rendering: From theory to implementation (Morgan Kaufmann, 2010).
9. C. Lemieux, Monte Carlo and Quasi-Monte Carlo Sampling (Springer, 2009).
10. P. Shirley, “Discrepancy as a quality measure for sample distributions,” Proceedings of Eurographics 91, 183–193 (1991).
11. C. Kolb, D. Mitchell, and P. Hanrahan, “A realistic camera model for computer graphics,” Proceedings of the 22nd annual conference on Computer graphics and interactive techniques, 317–324 (1995).
12. P. Shirley and R. K. Morley, Realistic ray tracing (AK Peters, Ltd., 2003).
13. J. S. Warren, Modern Lens Design (Washington: SPIE Press, 2005).
14. S. Horiuchi, S. Yoshida, Z. Ushiyama, and M. Yamamoto, “Performance evaluation of GRIN lenses by ray tracing method,” Opt. Quant. Electron. 42, 81–88 (2010) [CrossRef] .
