Research Article
Vol. 54, No. 14 / May 10 2015 / Applied Optics
4561
Image reconstruction method based on CCD calibration in frequency domain SHENG-JUN XIONG,1 BIN XIANGLI,1,2 YANG HE,1
AND
ZE ZHANG1,*
1
Academy of Opto-electronics, China Academy of Science, Beijing 100094, China e-mail: [email protected] *Corresponding author: [email protected]
2
Received 28 January 2015; revised 14 April 2015; accepted 16 April 2015; posted 16 April 2015 (Doc. ID 233433); published 8 May 2015
We demonstrate a method to reconstruct CCD images by calculating out the pixel response function accurately with laser interference patterns. This method is proven theoretically to have the ability to improve image quality greatly and thus may find great application in high-quality imaging fields. © 2015 Optical Society of America OCIS codes: (040.1520) CCD, charge-coupled device; (100.2000) Digital image processing; (110.3010) Image reconstruction techniques. http://dx.doi.org/10.1364/AO.54.004561
1. INTRODUCTION It has been more than 40 years since charge-coupled devices (CCDs) were invented around 1970. Due to the advantages of high quality, large capacity, and low price, CCDs have replaced chemical-based photography in most image taking systems [1]. Unlike photography, the CCD mechanism divides an image into many parts and displays each part with a pixel [2]. In ideal conditions, one would expect that different parts inside a pixel, and different pixels could respond to light intensity linearly, and display varied intensities with different gray values accordingly. However, in reality, each pixel gives only one value no matter how complex the incident light field. Moreover, due to an imperfect industrial manufacturing level or other unknown reasons, different parts of a pixel and different pixels respond dissimilarly even under the same illumination. Therefore, aberrations are caused when a CCD is used to display optical images. For common usage, these aberrations can be neglected due to low accuracy requirements. However, in applications such as astro-observation, police evidence investigation, and military reconnaissance, one may need to retrieve very detailed information from images, and thus aberrations may need to be highly reduced. Fundamentally, the aberrations are caused by two factors: the intrapixel and interpixel response. In previous research, quantum efficiency (QE), response uniformity and linearity, and geometric properties of pixel arrays are measured with various methods [3–11]. Among these methods, the flat field correction method (FFCM) is commonly used [7–10]. FFCM is very rough and can only reduce interpixel aberration to a very limited extent. To reduce CCD aberration dramatically, the straightforward way is to measure the CCD response function experimentally, as demonstrated by Kavaldjiev and Ninkov 1559-128X/15/144561-05$15/0$15.00 © 2015 Optical Society of America
[12,13]. In the researchers’ experimental system, they use a microscope system to focus and scan a laser beam across target CCDs to observe and measure the response intensity at each point. However, due to limited measurement accuracy and the lengthy measurement procedure, this method did not get well developed. Zhai et al. demonstrated a method to calibrate CCDs with laser interference patterns and measure the centroid position of the point object [14]. It has been shown that this method can reach 10−6 pixels theoretically and 10−5 pixels experimentally with centroid positional accuracy [14,15]. Image reconstruction technology has become more important in many application areas [16,17]. In this paper, we focus on the image reconstruction method by calibrating the CCD response function up to the tenth coefficient in frequency domain with laser interference patterns (IRTWLIP), and thus super-high-quality images can be achieved. The second part of this paper gives the theory of this method, and the third and fourth parts simulate and analyze the working mechanism, respectively. 2. THEORY In an imaging procedure, a pixel is always assumed to be the smallest processing unit. Namely, it is assumed that the response inside a pixel is uniform [1,2]. However, due to complex microstructures, imperfect craft, and pixel–pixel crosstalk, the response varies over the physical detection area of a pixel, which is called intrapixel variation. In addition, the pixel size and pixel grid also vary slightly in different areas, which is called interpixel variation. These variations surely affect image quality and cause aberrations. The measurement of these variations is essential to achieve high-quality images.
4562
Research Article
Vol. 54, No. 14 / May 10 2015 / Applied Optics
The imaging process of a pixel can be written as Z Z I mn dx dySx; yR mn x; y;
R˜ mn 0; 0 qmn0 (1)
where, I mn is the output of pixel m; n; Rx; y is the pixel response function (PRF); and Sx; y is the incident light field. A CCD has aberrations because Rx; y is always taken with approximate expressions. For example, in FFCM, Rx; y is taken as a constant function; in some subpixel-accuracy measurement, Rx; y is taken as a slightly modulated Gaussian function [18]. If Rx; y is measured out more accurately, then less aberration would remain. We assume Rx; y is written as R mn x p ; y p exp−x 2p y 2p ∕r 2g × c 0 c 1 x p c 2 y p c 3 x 2p
dyR mn x; y:
(7)
× expik x x p k y y p dx p dy p ;
(8)
where x mn ; ymn is the theoretical center, and Δx mn ; Δy mn is the displacement between the actual and theoretical centers. To solve R˜ mn kx ; k y , we use the McLaughlin formula to simplify Eq. (8) as R˜ mn kx ; ky q expfikx x mn Δx mn k y y mn Δy mn g mn0
c 4 y2p
× 1 qmn1 k x qmn2 k y qmn3 k 2x
c 5 x p y p …:
Please notice that x p ; yp is the coordinate’s inside pixel m; n; r g denotes the size of the Gaussian spot; and c i is the modulation coefficient. If c i is measured out more precisely, then Rx; y can be more accurately achieved. To measure Rx; y in frequency domain, we suppose ˜ x ; ky is the Fourier transform of Sx; y, and then Sk ZZ ˜ Skx ; ky Sx; y exp−ikx x ky ydxdy; (3) Sx; y
Z dx
In the coordinate of pixel m; n, Eq. (6) can be written as ˜ R mn kx ; ky expfikx x mn Δx mn ky y mn Δy mn g ZZ × R mn x p x mn Δx mn ; y p y mn Δy mn
qmn4 k 2y qmn5 k x ky …:
(2)
ZZ
Z
˜ x ; ky expik x x ky ydkx dk y ; Sk
(4)
where kx 2π∕d x , ky 2π∕d y represents the spatial frequency; and d x , d y are the spatial period. Insert Eq. (4) into Eq. (1): Z Z Z Z ˜ x ; ky R mn x; y I mn dx dy dkx dk y Sk × expikx x ky y Z Z Z Z ˜ x ; k y dx dyR mn x; y dkx dky Sk × expikx x ky y Z Z ˜ x ; k y R˜ mn k x ; ky ; dkx dky Sk
(9)
We use laser interference patterns to measure the coefficients in Eq. (9). The interference pattern generated by two laser beams I 1 and I 2 can be written as pffiffiffiffiffiffiffiffiffi I x; y; φt I 1 I 2 2 I 1 I 2 coskx x ky y φt : (10) Subscribe Eq. (10) to Eq. (1): Z Z I mn φt dx dyIx; y; φt R mn x; y Z Z pffiffiffiffiffiffiffiffiffi dx dyI 1 I 2 2 I 1 I 2 × coskx x ky y φt R mn x; y I 1 I 2 R˜ mn 0; 0 pffiffiffiffiffiffiffiffiffi 2 I 1 I 2 ReR˜ mn k x ; ky expiφt :
(11)
Then, subscribe Eq. (9) into Eq. (11): I mn φt I mn kx ;k y ;φt qmn0 I 1 I 2 2qmn0 1 qmn1 k2x qmn2 k 2y qmn3 kx ky … pffiffiffiffiffiffiffiffiffi × I 1 I 2 cosk x x mn Δx mn ky y mn Δy mn φt : (12)
(5)
where R˜ mn kx ; k y represents the Fourier transform of Rx; y, that is, ZZ R mn x; y expik x x k y ydxdy: (6) R˜ mn k x ; k y Equation (5) shows the imaging process of a CCD in frequency domain. It is worthwhile to mention that if the incident light field is uniform, then namely Sx; y constant. Here, we suppose Sx; y 1, then I mn0 R R dx dyR mn x; y. This is the so-called flat field correction that indicates that the output is the response integration of each pixel. Since the spatial period of a flat field is infinite, the spatial frequency is zero, and then Eq. (6) can be written as
From Eq. (12), we can see that with k x ; ky ; φt and I mn k x ; ky ; φt , coefficients Δx mn , Δy mn , q mni can be calculated out using the least square method, namely R˜ mn k x ; ky can be calculated out. Suppose CCD has N × N pixels, then the corresponding discrete coordinates can be written as kg ; k h g; h2π∕N : Correspondingly, Eq. (5) is written as XX ˜ g ; kh ·R˜ mn k g ; kh · Δk g · Δk h ; I mn Sk g
(13)
(14)
h
where g; h −N ∕2; −N ∕2 1; …; N ∕2 − 1 and Δk g ˜ x ; ky can be solved by Δkh 2π∕N are discrete steps. Sk N × N equations of Eq. (14). Then, the incident light field of Sx; y can be reconstructed by doing the inverse Fourier transformation:
Research Article Sx; y
XX g
Vol. 54, No. 14 / May 10 2015 / Applied Optics
4563
˜ g ; kh expikg x kh y · Δkg · Δkh : Sk
h
(15)
3. SIMULATION In this part, we simulate the image reconstruction process of a 30 × 30 CCD. For this CCD, the response function of pixel m; n is given by setting a series of random values to c i and r g in Eq. (2). Our purpose is to measure out the corresponding coefficients Δx mn , Δy mn , qmni and reconstruct the incident light field. To begin, we generate 25 interference patterns with known spatial frequency k x ; ky 0, π∕24, π∕12, π∕8, π∕6, and initial phase φt , and display these patterns with the target CCD. Then, the output of every pixel I mn kx ; ky ; φt can be read out. As discussed in Section 2, the coefficients of Δx mn , Δy mn , qmni can be obtained if kx ; k y ; φt and I mn kx ; k y ; φt are fed into Eq. (12) and thus R˜ mn k x ; ky can be matched out. In the end, using Eqs. (14) and (15), the incident light field can be reconstructed. In simulation, the imaging process of an optical system can be viewed as a low-pass filtering process. To give examples, we generate images of Airy spot, line, and grids, and filter their high frequency parts, as shown in Figs. 1(a), 2(a)and 3(a), respectively. First, FFCM is used to optimize these images, as shown in Figs. 1(b), 2(b)and 3(b). Then, IRTWLIP is used, as shown in Figs. 1(c), 2(c)and 3(c). From Figs. 1(c), 2(c) and 3(c), we can see that images taken with IRTWLIP have almost the same shape as the input light fields. Meanwhile, the images taken by FFCM are discrete, which means that the aberration is still large. To quantize the imaging quality, we measure the conveyed energy of these images and compare them with the input energy as listed in
Fig. 1. Simulations for imaging results of an Airy spot of (a) an original image, (b) after FFCM and (c) after IRTWLIP.
Fig. 2. Simulations for imaging results of a line of (a) an original image, (b) after FFCM and (c) after IRTWLIP.
Fig. 3. Simulations for imaging results of grids of (a) an original image, (b) after FFCM and (c) after IRTWLIP.
Table 1. Measured Energy of Two Different Processes J∕s · m2 Input Light Field
Input Energy
FFCM Image
IRTWLIP Image
Airy disk Line Grids
10000 10000 10000
9996.97 9984.83 10008.12
10000.146 9999.994 10000.103
Table 1. It is proven that the IRTWLIP method works much better. 4. ANALYSIS In this part, we analyze two issues of our image-reconstruction method: resolution and aberration. The resolution of a common imaging system is dependent on two subsystems: the optical system and the electrical detection system. In frequency domain, the system resolution is represented by high cutoff spatial frequency. For an ideal optical system, the high cutoff frequency is mainly related to pupil size and focal length: 2πD : (16) K O max λf For a CCD, the high cutoff frequency is mainly related to pixel size, not pixel response function. Therefore, no matter how accurate the pixel response function is measured out, the high cutoff frequency doesn’t change. Actually, the high cutoff frequency depends on the Nyquist theorem, which states that “an object needs to be sampled at least twice so as to be recognized.” For this reason, the high cutoff frequency of a CCD system is written as K D max π∕pixel:
(17)
The detection system is usually arranged right on the focal plane of the optical system so K O max has priority over K D max . This means that if K O max ≤ K D max , the resolution of the whole imaging system depends on K O max ; otherwise, it depends on K D max . For this reason, one needs to match the detection system to the optical system. Generally speaking, the PSF radius of an optical system needs to be designed larger than a pixel, which means K O max ≤ K D max [19]. Smaller spots will be a waste of system performance. Although the resolution couldn’t be improved, the aberration of the detecting system could be highly reduced with IRTWLIP. Here, we use a merit transform function (MTF)
4564
Research Article
Vol. 54, No. 14 / May 10 2015 / Applied Optics
calibrated with FFCM. Apparently, due to these great working effects, this method has prospective application in both civilian and military areas such as astro-observation, police evidence investigation, and military reconnaissance. National Natural Science Foundation of China (NSFC) (61308101, 61475161). REFERENCES
Fig. 4. MTF of a CCD. The blue line shows the MTF of an ideal CCD; the green line shows the MTF with FFCM; and the red line shows the MTF with IRTWLIP.
to study the performance of the detecting system. For an ideal CCD, the MTF can be defined as MTFdetector sincpk [20], in which p is the pixel size, and k is the spatial frequency of the image in optical system. From the equation, we can see that the MTF is getting smaller and smaller when k increases, as shown with the blue line in Fig. 4. To calculate the MTF of the CCD after using IRTWLIP, we take pictures for given sinusoidal fringes [21,22], and the MTF of the detecting system at special frequency can be measured using MTFk V out ∕V in ;
(18)
where V out and V in represent the contrast of output and input images, respectively, of given sinusoidal fringes, which can be calculated using I − Im V out∕in M ; (19) IM Im where I M and I m denote the maximum and minimum intensity, respectively, of the image of given sinusoidal fringes. Using Eqs. (18) and (19), the MTF at different spatial frequencies can be measured out, as marked with the red circles in Fig. 4. We can see that the MTF of a calibrated CCD in frequency domain is nearly 1, even in a high frequency area. In Fig. 4, the MTF of a CCD calibrated with FFCM is also given with a green line. The MTF with FFCM is slightly lower than that in the ideal condition due to the aberration caused by the nonuniformity of intrapixel response and the slight location shifts of the pixels. 5. CONCLUSION In this paper, we demonstrated a method to reconstruct CCD images by calibrating pixel response function in frequency domain. For example, we simulated the imaging process for objects such as the Airy spot, line, and grids. It is observed that the reconstructed images have much better quality than those taken by FFCM. To quantize the imaging quality, we measured the MTF of a CCD calibrated with different methods by taking images of standard sinusoidal fringes at different spatial frequencies. It was proven that the MTFs with IRTWLIP at all distinguishable spatial frequencies are nearly 1, which is much higher than that of the ideal CCD and the CCD
1. G. E. Smith, “The invention and early history of the CCD,” J. Appl. Phys 109, 102421 (2011). 2. I. Haq and A. A. Mudassar, “Geometric superresolution of a CCD pixel,” Opt. Lett. 35, 2705–2707 (2010). 3. T. C. Williams and C. R. Shaddix, “Simultaneous correction of flat field and nonlinearity response of intensified charge-coupled devices,” Rev. Sci. Instrum. 78, 123702 (2007). 4. S. Hubert and V. Prévot, “Quantum efficiency measurements of an x-ray charge-coupled device in the 2–10 keV spectral region by means of a double crystal monochromator coupled to an x-ray tube,” Appl. Opt. 53, 8078–8085 (2014). 5. A. Meda, I. Ruo-Berchera, I. P. Degiovanni, G. Brida, M. L. Rastello, and M. Genovese, “Absolute calibration of a charge-coupled device camera with twin beams,” Appl. Phys. Lett. 105, 101113 (2014). 6. F. W. Leslie, S. S. Cha, and N. Ramachandran, “A simple approach for calibrating solid-state imaging systems,” Rev. Sci. Instrum. 71, 4149–4154 (2000). 7. S. L. Barna, M. W. Tate, S. M. Gruner, and E. F. Eikenberry, “Calibration procedures for charge-coupled device x-ray detectors,” Rev. Sci. Instrum. 70, 2927–2934 (1999). 8. K. P. Klaasen, M. F. A’Hearn, M. Baca, A. Delamere, M. Desnoyer, T. Farnham, O. Groussin, D. Hampton, S. Ipatov, J. Li, C. Lisse, N. Mastrodemos, S. McLaughlin, J. Sunshine, P. Thomas, and D. Wellnitz, “Invited Article: Deep impact instrument calibration,” Rev. Sci. Instrum. 79, 091301 (2008). 9. S. M. Gruner, M. W. Tate, and E. F. Eikenberry, “Charge-coupled device area x-ray detectors,” Rev. Sci. Instrum. 73, 2815–2841 (2002). 10. G. Filacchione, E. Ammannito, A. Coradini, F. Capaccioni, G. Piccioni, M. C. De Sanctis, M. Dami, and A. Barbis, “On-ground characterization of Rosetta/VIRTIS-M. II. Spatial and radiometric calibrations,” Rev. Sci. Instrum. 77, 103106 (2006). 11. A. P. Hammersley, S. O. Svensson, A. Thompson, H. Graafsma, A. Kvick, and J. P. Moy, “Calibration and correction of distortions in twodimensional detector systems,” Rev. Sci. Instrum. 66, 2729–2733 (1995). 12. D. Kavaldjiev and Z. Ninkov, “Subpixel sensitivity map for a chargecoupled device sensor,” Opt. Eng. 37, 948–954 (1998). 13. A. Piterman and Z. Ninkov, “Sub-pixel sensitivity maps for a backilluminated CCD and the effects of non-uniform response on measurement accuracy,” Opt. Eng. 41, 1192–1202 (2002). 14. C. Zhai, M. Shao, R. Goullioud, and B. Nemati, “Micro-pixel accuracy centroid displacement estimation and detector calibration,” Proc. R. Soc. A. 467, 3550–3569 (2011). 15. B. Nemati, M. Shao, C. Zhai, H. Erlig, X. Wang, and R. Goullioud, “Micro-pixel image position sensing testbed,” Proc. SPIE 8151, 81510W (2011). 16. D. A. Ralys and V. Dadurkevicius, “Image recovering from the modulated CCD imagers,” Astrophys. Space Sci. 139, 155–162 (1987). 17. Y. Ohneda, N. Baba, N. Miura, and T. Sakura, “Multiresolution approach to image reconstruction with phase-diversity technique,” Opt. Rev. 8, 32–36 (2001). 18. J. G. Ingalls, J. E. Krick, S. J. Garey, S. Laine, J. A. Surace, W. J. Glaccum, C. C. Grillmair, and P. J. Lowrance, “Intra-pixel gain variations and high-precision photometry with the infrared array camera (IRAC),” Proc. SPIE 8442, 84421Y (2012). 19. B. F. Alexander and K. C. Ng, “Elimination of systematic error in subpixel accuracy centroid estimation,” Opt. Eng. 30, 1320–1331 (1991).
Research Article 20. M. M. Blouke, J. R. Janesick, J. E. Hall, M. W. Cowens, and P. J. May, “800 × 800 charge-coupled device image sensor,” Opt. Eng. 22, 225607 (1983). 21. S. Najafi and K. Madanipour, “Measurement of the modulation transfer function of a charge-coupled device array by the combination
Vol. 54, No. 14 / May 10 2015 / Applied Optics
4565
of the self-imaging effect and slanted edge method,” Appl. Opt. 52, 4724–4727 (2013). 22. S. R. Shortes, W. W. Chan, W. C. Rhines, J. B. Barton, and D. R. Collins, “Characteristics of thinned backside-illuminated chargecoupled device imagers,” Appl. Phys. Lett. 24, 565–567 (1974).
Vol. 54, No. 14 / May 10 2015 / Applied Optics
4561
Image reconstruction method based on CCD calibration in frequency domain SHENG-JUN XIONG,1 BIN XIANGLI,1,2 YANG HE,1
AND
ZE ZHANG1,*
1
Academy of Opto-electronics, China Academy of Science, Beijing 100094, China e-mail: [email protected] *Corresponding author: [email protected]
2
Received 28 January 2015; revised 14 April 2015; accepted 16 April 2015; posted 16 April 2015 (Doc. ID 233433); published 8 May 2015
We demonstrate a method to reconstruct CCD images by calculating out the pixel response function accurately with laser interference patterns. This method is proven theoretically to have the ability to improve image quality greatly and thus may find great application in high-quality imaging fields. © 2015 Optical Society of America OCIS codes: (040.1520) CCD, charge-coupled device; (100.2000) Digital image processing; (110.3010) Image reconstruction techniques. http://dx.doi.org/10.1364/AO.54.004561
1. INTRODUCTION It has been more than 40 years since charge-coupled devices (CCDs) were invented around 1970. Due to the advantages of high quality, large capacity, and low price, CCDs have replaced chemical-based photography in most image taking systems [1]. Unlike photography, the CCD mechanism divides an image into many parts and displays each part with a pixel [2]. In ideal conditions, one would expect that different parts inside a pixel, and different pixels could respond to light intensity linearly, and display varied intensities with different gray values accordingly. However, in reality, each pixel gives only one value no matter how complex the incident light field. Moreover, due to an imperfect industrial manufacturing level or other unknown reasons, different parts of a pixel and different pixels respond dissimilarly even under the same illumination. Therefore, aberrations are caused when a CCD is used to display optical images. For common usage, these aberrations can be neglected due to low accuracy requirements. However, in applications such as astro-observation, police evidence investigation, and military reconnaissance, one may need to retrieve very detailed information from images, and thus aberrations may need to be highly reduced. Fundamentally, the aberrations are caused by two factors: the intrapixel and interpixel response. In previous research, quantum efficiency (QE), response uniformity and linearity, and geometric properties of pixel arrays are measured with various methods [3–11]. Among these methods, the flat field correction method (FFCM) is commonly used [7–10]. FFCM is very rough and can only reduce interpixel aberration to a very limited extent. To reduce CCD aberration dramatically, the straightforward way is to measure the CCD response function experimentally, as demonstrated by Kavaldjiev and Ninkov 1559-128X/15/144561-05$15/0$15.00 © 2015 Optical Society of America
[12,13]. In the researchers’ experimental system, they use a microscope system to focus and scan a laser beam across target CCDs to observe and measure the response intensity at each point. However, due to limited measurement accuracy and the lengthy measurement procedure, this method did not get well developed. Zhai et al. demonstrated a method to calibrate CCDs with laser interference patterns and measure the centroid position of the point object [14]. It has been shown that this method can reach 10−6 pixels theoretically and 10−5 pixels experimentally with centroid positional accuracy [14,15]. Image reconstruction technology has become more important in many application areas [16,17]. In this paper, we focus on the image reconstruction method by calibrating the CCD response function up to the tenth coefficient in frequency domain with laser interference patterns (IRTWLIP), and thus super-high-quality images can be achieved. The second part of this paper gives the theory of this method, and the third and fourth parts simulate and analyze the working mechanism, respectively. 2. THEORY In an imaging procedure, a pixel is always assumed to be the smallest processing unit. Namely, it is assumed that the response inside a pixel is uniform [1,2]. However, due to complex microstructures, imperfect craft, and pixel–pixel crosstalk, the response varies over the physical detection area of a pixel, which is called intrapixel variation. In addition, the pixel size and pixel grid also vary slightly in different areas, which is called interpixel variation. These variations surely affect image quality and cause aberrations. The measurement of these variations is essential to achieve high-quality images.
4562
Research Article
Vol. 54, No. 14 / May 10 2015 / Applied Optics
The imaging process of a pixel can be written as Z Z I mn dx dySx; yR mn x; y;
R˜ mn 0; 0 qmn0 (1)
where, I mn is the output of pixel m; n; Rx; y is the pixel response function (PRF); and Sx; y is the incident light field. A CCD has aberrations because Rx; y is always taken with approximate expressions. For example, in FFCM, Rx; y is taken as a constant function; in some subpixel-accuracy measurement, Rx; y is taken as a slightly modulated Gaussian function [18]. If Rx; y is measured out more accurately, then less aberration would remain. We assume Rx; y is written as R mn x p ; y p exp−x 2p y 2p ∕r 2g × c 0 c 1 x p c 2 y p c 3 x 2p
dyR mn x; y:
(7)
× expik x x p k y y p dx p dy p ;
(8)
where x mn ; ymn is the theoretical center, and Δx mn ; Δy mn is the displacement between the actual and theoretical centers. To solve R˜ mn kx ; k y , we use the McLaughlin formula to simplify Eq. (8) as R˜ mn kx ; ky q expfikx x mn Δx mn k y y mn Δy mn g mn0
c 4 y2p
× 1 qmn1 k x qmn2 k y qmn3 k 2x
c 5 x p y p …:
Please notice that x p ; yp is the coordinate’s inside pixel m; n; r g denotes the size of the Gaussian spot; and c i is the modulation coefficient. If c i is measured out more precisely, then Rx; y can be more accurately achieved. To measure Rx; y in frequency domain, we suppose ˜ x ; ky is the Fourier transform of Sx; y, and then Sk ZZ ˜ Skx ; ky Sx; y exp−ikx x ky ydxdy; (3) Sx; y
Z dx
In the coordinate of pixel m; n, Eq. (6) can be written as ˜ R mn kx ; ky expfikx x mn Δx mn ky y mn Δy mn g ZZ × R mn x p x mn Δx mn ; y p y mn Δy mn
qmn4 k 2y qmn5 k x ky …:
(2)
ZZ
Z
˜ x ; ky expik x x ky ydkx dk y ; Sk
(4)
where kx 2π∕d x , ky 2π∕d y represents the spatial frequency; and d x , d y are the spatial period. Insert Eq. (4) into Eq. (1): Z Z Z Z ˜ x ; ky R mn x; y I mn dx dy dkx dk y Sk × expikx x ky y Z Z Z Z ˜ x ; k y dx dyR mn x; y dkx dky Sk × expikx x ky y Z Z ˜ x ; k y R˜ mn k x ; ky ; dkx dky Sk
(9)
We use laser interference patterns to measure the coefficients in Eq. (9). The interference pattern generated by two laser beams I 1 and I 2 can be written as pffiffiffiffiffiffiffiffiffi I x; y; φt I 1 I 2 2 I 1 I 2 coskx x ky y φt : (10) Subscribe Eq. (10) to Eq. (1): Z Z I mn φt dx dyIx; y; φt R mn x; y Z Z pffiffiffiffiffiffiffiffiffi dx dyI 1 I 2 2 I 1 I 2 × coskx x ky y φt R mn x; y I 1 I 2 R˜ mn 0; 0 pffiffiffiffiffiffiffiffiffi 2 I 1 I 2 ReR˜ mn k x ; ky expiφt :
(11)
Then, subscribe Eq. (9) into Eq. (11): I mn φt I mn kx ;k y ;φt qmn0 I 1 I 2 2qmn0 1 qmn1 k2x qmn2 k 2y qmn3 kx ky … pffiffiffiffiffiffiffiffiffi × I 1 I 2 cosk x x mn Δx mn ky y mn Δy mn φt : (12)
(5)
where R˜ mn kx ; k y represents the Fourier transform of Rx; y, that is, ZZ R mn x; y expik x x k y ydxdy: (6) R˜ mn k x ; k y Equation (5) shows the imaging process of a CCD in frequency domain. It is worthwhile to mention that if the incident light field is uniform, then namely Sx; y constant. Here, we suppose Sx; y 1, then I mn0 R R dx dyR mn x; y. This is the so-called flat field correction that indicates that the output is the response integration of each pixel. Since the spatial period of a flat field is infinite, the spatial frequency is zero, and then Eq. (6) can be written as
From Eq. (12), we can see that with k x ; ky ; φt and I mn k x ; ky ; φt , coefficients Δx mn , Δy mn , q mni can be calculated out using the least square method, namely R˜ mn k x ; ky can be calculated out. Suppose CCD has N × N pixels, then the corresponding discrete coordinates can be written as kg ; k h g; h2π∕N : Correspondingly, Eq. (5) is written as XX ˜ g ; kh ·R˜ mn k g ; kh · Δk g · Δk h ; I mn Sk g
(13)
(14)
h
where g; h −N ∕2; −N ∕2 1; …; N ∕2 − 1 and Δk g ˜ x ; ky can be solved by Δkh 2π∕N are discrete steps. Sk N × N equations of Eq. (14). Then, the incident light field of Sx; y can be reconstructed by doing the inverse Fourier transformation:
Research Article Sx; y
XX g
Vol. 54, No. 14 / May 10 2015 / Applied Optics
4563
˜ g ; kh expikg x kh y · Δkg · Δkh : Sk
h
(15)
3. SIMULATION In this part, we simulate the image reconstruction process of a 30 × 30 CCD. For this CCD, the response function of pixel m; n is given by setting a series of random values to c i and r g in Eq. (2). Our purpose is to measure out the corresponding coefficients Δx mn , Δy mn , qmni and reconstruct the incident light field. To begin, we generate 25 interference patterns with known spatial frequency k x ; ky 0, π∕24, π∕12, π∕8, π∕6, and initial phase φt , and display these patterns with the target CCD. Then, the output of every pixel I mn kx ; ky ; φt can be read out. As discussed in Section 2, the coefficients of Δx mn , Δy mn , qmni can be obtained if kx ; k y ; φt and I mn kx ; k y ; φt are fed into Eq. (12) and thus R˜ mn k x ; ky can be matched out. In the end, using Eqs. (14) and (15), the incident light field can be reconstructed. In simulation, the imaging process of an optical system can be viewed as a low-pass filtering process. To give examples, we generate images of Airy spot, line, and grids, and filter their high frequency parts, as shown in Figs. 1(a), 2(a)and 3(a), respectively. First, FFCM is used to optimize these images, as shown in Figs. 1(b), 2(b)and 3(b). Then, IRTWLIP is used, as shown in Figs. 1(c), 2(c)and 3(c). From Figs. 1(c), 2(c) and 3(c), we can see that images taken with IRTWLIP have almost the same shape as the input light fields. Meanwhile, the images taken by FFCM are discrete, which means that the aberration is still large. To quantize the imaging quality, we measure the conveyed energy of these images and compare them with the input energy as listed in
Fig. 1. Simulations for imaging results of an Airy spot of (a) an original image, (b) after FFCM and (c) after IRTWLIP.
Fig. 2. Simulations for imaging results of a line of (a) an original image, (b) after FFCM and (c) after IRTWLIP.
Fig. 3. Simulations for imaging results of grids of (a) an original image, (b) after FFCM and (c) after IRTWLIP.
Table 1. Measured Energy of Two Different Processes J∕s · m2 Input Light Field
Input Energy
FFCM Image
IRTWLIP Image
Airy disk Line Grids
10000 10000 10000
9996.97 9984.83 10008.12
10000.146 9999.994 10000.103
Table 1. It is proven that the IRTWLIP method works much better. 4. ANALYSIS In this part, we analyze two issues of our image-reconstruction method: resolution and aberration. The resolution of a common imaging system is dependent on two subsystems: the optical system and the electrical detection system. In frequency domain, the system resolution is represented by high cutoff spatial frequency. For an ideal optical system, the high cutoff frequency is mainly related to pupil size and focal length: 2πD : (16) K O max λf For a CCD, the high cutoff frequency is mainly related to pixel size, not pixel response function. Therefore, no matter how accurate the pixel response function is measured out, the high cutoff frequency doesn’t change. Actually, the high cutoff frequency depends on the Nyquist theorem, which states that “an object needs to be sampled at least twice so as to be recognized.” For this reason, the high cutoff frequency of a CCD system is written as K D max π∕pixel:
(17)
The detection system is usually arranged right on the focal plane of the optical system so K O max has priority over K D max . This means that if K O max ≤ K D max , the resolution of the whole imaging system depends on K O max ; otherwise, it depends on K D max . For this reason, one needs to match the detection system to the optical system. Generally speaking, the PSF radius of an optical system needs to be designed larger than a pixel, which means K O max ≤ K D max [19]. Smaller spots will be a waste of system performance. Although the resolution couldn’t be improved, the aberration of the detecting system could be highly reduced with IRTWLIP. Here, we use a merit transform function (MTF)
4564
Research Article
Vol. 54, No. 14 / May 10 2015 / Applied Optics
calibrated with FFCM. Apparently, due to these great working effects, this method has prospective application in both civilian and military areas such as astro-observation, police evidence investigation, and military reconnaissance. National Natural Science Foundation of China (NSFC) (61308101, 61475161). REFERENCES
Fig. 4. MTF of a CCD. The blue line shows the MTF of an ideal CCD; the green line shows the MTF with FFCM; and the red line shows the MTF with IRTWLIP.
to study the performance of the detecting system. For an ideal CCD, the MTF can be defined as MTFdetector sincpk [20], in which p is the pixel size, and k is the spatial frequency of the image in optical system. From the equation, we can see that the MTF is getting smaller and smaller when k increases, as shown with the blue line in Fig. 4. To calculate the MTF of the CCD after using IRTWLIP, we take pictures for given sinusoidal fringes [21,22], and the MTF of the detecting system at special frequency can be measured using MTFk V out ∕V in ;
(18)
where V out and V in represent the contrast of output and input images, respectively, of given sinusoidal fringes, which can be calculated using I − Im V out∕in M ; (19) IM Im where I M and I m denote the maximum and minimum intensity, respectively, of the image of given sinusoidal fringes. Using Eqs. (18) and (19), the MTF at different spatial frequencies can be measured out, as marked with the red circles in Fig. 4. We can see that the MTF of a calibrated CCD in frequency domain is nearly 1, even in a high frequency area. In Fig. 4, the MTF of a CCD calibrated with FFCM is also given with a green line. The MTF with FFCM is slightly lower than that in the ideal condition due to the aberration caused by the nonuniformity of intrapixel response and the slight location shifts of the pixels. 5. CONCLUSION In this paper, we demonstrated a method to reconstruct CCD images by calibrating pixel response function in frequency domain. For example, we simulated the imaging process for objects such as the Airy spot, line, and grids. It is observed that the reconstructed images have much better quality than those taken by FFCM. To quantize the imaging quality, we measured the MTF of a CCD calibrated with different methods by taking images of standard sinusoidal fringes at different spatial frequencies. It was proven that the MTFs with IRTWLIP at all distinguishable spatial frequencies are nearly 1, which is much higher than that of the ideal CCD and the CCD
1. G. E. Smith, “The invention and early history of the CCD,” J. Appl. Phys 109, 102421 (2011). 2. I. Haq and A. A. Mudassar, “Geometric superresolution of a CCD pixel,” Opt. Lett. 35, 2705–2707 (2010). 3. T. C. Williams and C. R. Shaddix, “Simultaneous correction of flat field and nonlinearity response of intensified charge-coupled devices,” Rev. Sci. Instrum. 78, 123702 (2007). 4. S. Hubert and V. Prévot, “Quantum efficiency measurements of an x-ray charge-coupled device in the 2–10 keV spectral region by means of a double crystal monochromator coupled to an x-ray tube,” Appl. Opt. 53, 8078–8085 (2014). 5. A. Meda, I. Ruo-Berchera, I. P. Degiovanni, G. Brida, M. L. Rastello, and M. Genovese, “Absolute calibration of a charge-coupled device camera with twin beams,” Appl. Phys. Lett. 105, 101113 (2014). 6. F. W. Leslie, S. S. Cha, and N. Ramachandran, “A simple approach for calibrating solid-state imaging systems,” Rev. Sci. Instrum. 71, 4149–4154 (2000). 7. S. L. Barna, M. W. Tate, S. M. Gruner, and E. F. Eikenberry, “Calibration procedures for charge-coupled device x-ray detectors,” Rev. Sci. Instrum. 70, 2927–2934 (1999). 8. K. P. Klaasen, M. F. A’Hearn, M. Baca, A. Delamere, M. Desnoyer, T. Farnham, O. Groussin, D. Hampton, S. Ipatov, J. Li, C. Lisse, N. Mastrodemos, S. McLaughlin, J. Sunshine, P. Thomas, and D. Wellnitz, “Invited Article: Deep impact instrument calibration,” Rev. Sci. Instrum. 79, 091301 (2008). 9. S. M. Gruner, M. W. Tate, and E. F. Eikenberry, “Charge-coupled device area x-ray detectors,” Rev. Sci. Instrum. 73, 2815–2841 (2002). 10. G. Filacchione, E. Ammannito, A. Coradini, F. Capaccioni, G. Piccioni, M. C. De Sanctis, M. Dami, and A. Barbis, “On-ground characterization of Rosetta/VIRTIS-M. II. Spatial and radiometric calibrations,” Rev. Sci. Instrum. 77, 103106 (2006). 11. A. P. Hammersley, S. O. Svensson, A. Thompson, H. Graafsma, A. Kvick, and J. P. Moy, “Calibration and correction of distortions in twodimensional detector systems,” Rev. Sci. Instrum. 66, 2729–2733 (1995). 12. D. Kavaldjiev and Z. Ninkov, “Subpixel sensitivity map for a chargecoupled device sensor,” Opt. Eng. 37, 948–954 (1998). 13. A. Piterman and Z. Ninkov, “Sub-pixel sensitivity maps for a backilluminated CCD and the effects of non-uniform response on measurement accuracy,” Opt. Eng. 41, 1192–1202 (2002). 14. C. Zhai, M. Shao, R. Goullioud, and B. Nemati, “Micro-pixel accuracy centroid displacement estimation and detector calibration,” Proc. R. Soc. A. 467, 3550–3569 (2011). 15. B. Nemati, M. Shao, C. Zhai, H. Erlig, X. Wang, and R. Goullioud, “Micro-pixel image position sensing testbed,” Proc. SPIE 8151, 81510W (2011). 16. D. A. Ralys and V. Dadurkevicius, “Image recovering from the modulated CCD imagers,” Astrophys. Space Sci. 139, 155–162 (1987). 17. Y. Ohneda, N. Baba, N. Miura, and T. Sakura, “Multiresolution approach to image reconstruction with phase-diversity technique,” Opt. Rev. 8, 32–36 (2001). 18. J. G. Ingalls, J. E. Krick, S. J. Garey, S. Laine, J. A. Surace, W. J. Glaccum, C. C. Grillmair, and P. J. Lowrance, “Intra-pixel gain variations and high-precision photometry with the infrared array camera (IRAC),” Proc. SPIE 8442, 84421Y (2012). 19. B. F. Alexander and K. C. Ng, “Elimination of systematic error in subpixel accuracy centroid estimation,” Opt. Eng. 30, 1320–1331 (1991).
Research Article 20. M. M. Blouke, J. R. Janesick, J. E. Hall, M. W. Cowens, and P. J. May, “800 × 800 charge-coupled device image sensor,” Opt. Eng. 22, 225607 (1983). 21. S. Najafi and K. Madanipour, “Measurement of the modulation transfer function of a charge-coupled device array by the combination
Vol. 54, No. 14 / May 10 2015 / Applied Optics
4565
of the self-imaging effect and slanted edge method,” Appl. Opt. 52, 4724–4727 (2013). 22. S. R. Shortes, W. W. Chan, W. C. Rhines, J. B. Barton, and D. R. Collins, “Characteristics of thinned backside-illuminated chargecoupled device imagers,” Appl. Phys. Lett. 24, 565–567 (1974).
