1. Introduction

In a basic scheme of a line-scan hyperspectral imaging system (HIS) the observed sample is illuminated by a broadband light source. The diffused or transmitted light, containing information on the observed sample, is collected by the input optics, then passed through a narrow slit, and finally dispersed across the detector array forming one spatial slice of the hyperspectral image. For proper interpretation of the recorded information either an absolute relation between the physical properties of the observed sample and the recorded image or at least instrument-level uniformity with respect to the spectral and spatial position is required. Due to the optical properties and aberrations of the integrated lens (e.g. astigmatism, distortion, chromatic aberration, etc.), and aberrations arising from the dispersing element, the imaging system response is strongly spectrally and spatially dependent [1]. In an ideal image, free of geometric coregistration errors, the spectral and the spatial information would be aligned with the image columns (spatial channels) and rows (spectral channels), respectively. However, due to the higher order effects the recorded image is geometrically distorted. Speaking in relative terms, the recorded image appears wrapped in the spectral direction due to change of the dispersion angle with the spatial position (smile). Similar effect, arising from the wavelength-dependent magnification, is observed in the spatial direction (keystone). Images can be additionally distorted in both directions by rotational and angular misalignment between the spectrograph and the camera. On top of the geometric distortions the uniformity of the recorded image is further degraded by a spectrally and spatially dependent point spread function (PSF). Therefore, a simple global linear model (scaling and shifting of the image) is insufficient for the proper spectral and spatial calibration of the HIS.

Considering that the artifacts arising from the geometric misalignments and position variant PSF are well known issues in the field of remote sensing, the problem is surprisingly rarely addressed in the field of chemical imaging. To the best of our knowledge, there is only one publication [2] studying and dealing with the geometric coregistration problem, which was conducted in the field of Raman spectroscopy. The study has shown that even minor geometric misalignments lead to significant nonlinear intensity errors in the acquired spectra, which degrade the subsequent hyperspectral image analysis. The findings were in agreement with the studies made in the field of remote sensing, where numerous methods were proposed to both analyze and alleviate the impact of the problem [3–6]. However, these methods are specific to the recorded data, employed equipment and challenges in the field of remote sensing, and usually discuss only the spectral calibration.

Polder et al. [7] proposed a set of practical procedures suitable for characterization of several laboratory HIS setups. The spectral calibration was performed by finding the polynomial relation between measured peak locations of several known spectra and corresponding real spectral line wavelengths. Spectral resolution was reported as the full width at half maximum (FWHM) of a laser spectral line at 670 nm. The spatial resolution was assessed in terms of FWHM at three wavelengths by several methods. Most reliable results were obtained by calculating the derivative of the b-spline interpolated grey levels of a black-to-white step edge. The image position variant spectral response function (SRF) is usually estimated by stepping trough the spectral range by means of a spectrally tunable light source ([8, 9]). Unfortunately, this approach requires a calibrated monochromator and a relatively complex measurement setup. Analogous approach was used for the spatial calibration by observing a distant, physically thin, light source at predetermined positions [10]. The proposed practical method [2] for mitigation of the geometric coregistration problem in a Raman imaging spectrograph, utilized a feature-based image registration. The required control points were extracted from a pair of reference images corresponding to a neon light source and patterned white-illumination. Unfortunately, the SRF and the spatial response function (SiRF) were not assessed in this work. Most popular methods for the measurement of the optical system line spread function (LSF), closely related to the SiRF, are based on precise sub-pixel shifting of the calibration object [11, 12]. While this method is straightforward and accurate, it requires a precise mechanical actuator. Latter can pose a problem either when it is necessary to characterize an optical assembly prior to integration into the final setup or when scanning is performed on a system lacking the ability of precise positioning, e.g. conveyor belt.

In summary, the previously published methods either do not consider the complete spatial-spectral dependence of the HIS or require complex measurement setups, which render them impractical. To alleviate this problem, we propose a novel method based on modeling of the HIS response to a scene exhibiting well-defined spectral and spatial features. For a reliable and unambiguous modeling, spatial and spectral features were retrieved by acquisition and analysis of images corresponding to different spectral reference lamps and a spatial reference object. The proposed method requires a relatively simple measurement setup and does not depend on any special equipment.

2. Methodology

Let $\{g(\theta ,\lambda )\in {ℝ}^2\}$ represent the true object defined in object coordinates (θ, λ), where θ is a spatial position and λ a wavelength. Firstly, the geometric properties of the HIS are captured by mapping the true object into the image space by a forward tuning function$\{T:(\theta ,\lambda )\to (\upsilon ,\omega )\in {ℝ}^2\}$. Secondly, presence of the blur is modeled by convolving the transformed object and the PSF of the HIS. Assuming the 2D PSF is separable it can be broken into a pair of 1D line spread functions, namely the spectral response function (SRF) and the spatial response function (SiRF). SRF and SiRF functions take into account the properties of the entire HIS, i.e. the lens, spectrograph and camera. Therefore, the acquired image I(u,w) of size U × W pixels can be estimated within error ε as:

 

Fig. 1 Example of a spectral (Ne) (a) and a spatial (b) reference image.

Download Full Size | PDF

The image (S_A) of the spectral reference object (g_S) provides a good definition of the properties in the spectral direction and is vastly insensitive to the distortions in the spatial direction, while the opposite holds true for the image (Si_A) of the spatial reference object (g_Si). In general, the SRF and SiRF are aligned with the spatial and the spectral axis of the spectrograph and not necessarily with the columns and rows of the acquired image. However, to lower the computational complexity of the proposed method, the alignment of the spatial and spectral axes between the spectrograph and the camera is assumed. The assumption holds if the initial geometric coregistration errors are small, which is realistic, since the proposed method is used after the proper adjustments are made by other, sufficiently accurate procedures [1]. Furthermore, the SRF and SiRF are position variant. Therefore, the Eq. (2) and Eq. (3) should be directly computed only for small regions, where the local invariance can be adopted.

The first goal of the characterization is the assessment of the forward (υ,ω) = (T_S(θ,λ),T_Si(θ,λ)) and of the corresponding inverse tuning function (θ,λ) = (Z_S(υ,ω),Z_Si(υ,ω)). The first is used to estimate the image of a given object, while the second is required for a physically meaningful interpretation of the acquired image. The tuning functions are determined by establishing a correspondence between the features in the object space and the corresponding positions in the acquired image. For this purpose, the feature locations have to be extracted from the acquired image. Furthermore, the quality of HIS is also defined by its resolution, therefore the position dependent response functions SRF(υ,ω) and SiRF(υ,ω) are also determined. With respect to the HIS properties and technological limitations of the reference object design, solving the two equations, namely Eq. (2) and Eq. (3), requires a different approach.

2.1 Spectral feature extraction and SRF characterization

Since the atomic emission spectral lines are significantly narrower than the width of the instrument’s SRF, these can be considered as shifted Dirac delta functions. Assuming locally invariant SRF and considering the delayed delta function shifting role with respect to the convolution, the measured response to a single spectral line is a good approximation of the SRF. In the majority of publications SRF is modeled as a Gaussian function. However, for the studied HIS, such model exhibited poor fit of the acquired data. Significantly better results were produced using the symmetric normalized Voigt profile V(ω-w₀;γ_L,γ_G), defined as a convolution of the centered Gaussian and Lorentzian broadening functions. Therefore, spectral response of the HIS to a single spectral line, centered at a position w₀ = T_S(θ,λ₀), is completely defined by two spectrally and spatially variant SRF parameters (γ_L,γ_G). Using a fast approximation formula proposed in [13] this model retains global parameterization and reasonable computational efficiency.

Coverage of the utilized spectral range by a single spectral reference source is limited. Therefore, separate spectra of L different sources, each contributing J_l unique spectral lines, were simultaneously modeled. A given j-th spectral line in the spectrum of l-th reference lamp is defined by its wavelength λ_0,l,j and intensity. By modeling the spectral line with delta function, transforming it by forward spectral tuning function, and modeling intensity α_l,j, the T_s(g_S(θ,λ)) in Eq. (2) is substituted by α_l,j∙δ_u(ω-w_0,l,j,u). Using the Voigt profile as a SRF model in Eq. (2), and taking into account the convolution rule of SRF and shifted delta function, the modeled spectrum S_M,l,u of the l-th reference lamp corresponding to a spatial channel u of S_A,l can be written as a sum of responses to the individual spectral lines, i.e. shifted Voigt profiles α_l,j∙V(ω-w_0,l,j,u,γ_G(u,w_0,l,j,u),γ_L(u,w_0,l,j,u)). Approximating the integral in Eq. (2) over a rectangular area of size δu × δw by the rectangle rule [14], and adding a polynomial model of the baseline shift, the final spectrum model can be written as:

The three utilized spectral reference lamps (L = 3) shown in Fig. 2 provided 55 useful spectral lines. The most informative spectral range and the corresponding spectral lines were identified for each lamp. Three mercury lines of the Hg(Ar) lamp in the range from 527 nm to 588 nm were used and the Ne-lamp provided further 28 spectral lines in the range from 565 nm to 772 nm. Since the Hg(Ar) lamp exhibited a large intensity contrast between the argon and the mercury spectral lines, the additional Ar lamp was utilized in the range from 732 nm to 990 nm, providing 24 spectral lines. The accurate wavelengths of all the spectral lines were taken from the NIST spectral database [16].

 

Fig. 2 Weighted acquired spectra of the three utilized reference sources. The problem of low spectral resolution is clearly demonstrated by the two recorded spectral peaks of the Hg(Ar) lamp. The spectral peak Hg-I corresponds to the single Hg spectral line at 546.07 nm while the Hg-II peak is a response to two Hg spectral lines at 576.96 nm and 579.07 nm, unresolvable due to the limited instrument resolution and discrete sampling. This is the reason for the difference between the number of observed peaks and the number of tabulated spectral lines. Furthermore, at shorter wavelengths peaks are approximately Gaussian. However, towards the NIR end of the observed spectral range, the spectral peaks (e.g. Ar-I at 965.78 nm) become wider and their tails longer than the Gaussian function would model.

Download Full Size | PDF

Computation of the optimal parameters was performed by a nested approach. Parameters defining the properties of the HIS (b_G, b_L, b_S) were optimized by outer Broyden-Fletcher-Goldfarb-Shanno (BFGS) [17] nonlinear optimization method. This gradient-based iterative quasi-Newton method exhibits quadratic convergence and employs the observed behavior of the dissimilarity χ_S and its gradient for construction of the Hessian matrix approximation and step size selection. In each BFGS iteration the optimal auxiliary parameters α_l,j and β_l,k,u, capturing the unknown experimental variations (HIS independent), were computed by the constrained linear least-squares imposing the non-negativity of the spectral line intensities. Initial approximation of b_S were computed by a linear fit between the real and the acquired positions of the distinct spectral peaks (Hg:546.10 nm, Ne:671.70 nm, Ar:763.51 nm, 965.78 nm). While the initial SRF parameters b_G and b_L were set to values fitting a single Ar peak at the central part of the acquired image (Ar:763.51 nm). Prior to optimization the parameters were scaled to improve convergence and aid the proper selection of stopping criterion. Furthermore, convergence to local minima was prevented by iteratively increasing the number of breaks employed by the subsequent runs of the optimization. The number of breaks K_B and K_G act as regularization parameters. For values less than four the model exhibited underfitting. Between five and seven the optimization problem was properly regularized and best results were achieved. Further increasing the number of breaks led to overfitting. Since the computation time increases with the number of breaks, both K_B and K_G were set to five. Optimal parameters for a single spectrum were computed in approximately 11 s on a current generation personal computer.

The elements of the weight arrays µ_l(w), corresponding to the lamp specific informative spectral range were set to unit and left zero otherwise. In the HIS under study, the baseline shift was mainly contributed by the thermal drift of the camera sensor. The approximate shape of the baseline was assessed by recording an image of white reference material illuminated by a halogen broadband light and a series of images corresponding to the camera dark current. The frequency at which the dark frames were captured was low enough that the camera drift was measurable. For each frame a difference between the corresponding dark frame and the first dark frame in the series was computed. This difference was further divided by the difference between the white frame and the first dark frame in the series. The resulting series of images represented the time-lapse of baseline drift for each spatial channel. By fitting each baseline in the series it was determined that a baseline can be approximated by second order polynomial. Since the weights were nonzero only within narrow spectral bands, baseline was locally modeled by a first order polynomial (K = 1).

2.2 Spatial feature extraction and nonparametric SiRF characterization

The spatial reference object can be designed to suit the task at hand and, unlike the spectral lines, the modeled features can be well separated. Due to the manufacturing limitations the delta function approximation of these features is hard to achieve in an affordable reference object. Furthermore, in contrast to the SRF characterization, modeling the spatial resolution is significantly more demanding due to the SiRF complexity. However, reliable and computationally straightforward results were obtained by using a series of black-to-white step edges. Each spatial profile of the acquired spatial reference image at a given spectral channel w, may be interpreted as a sequence of blurred edges of alternating orientations. The model of the single edge with transition located at the image position u₀ is written as a Heaviside step function H(υ-u₀) modulated by a multiplicative c_M and an additive term c_A, arising from uneven illumination, temperature drift of the sensor and the reflected light. Therefore, in Eq. (3), the term T_Si(g_Si(θ,λ)) is substituted by the edge model at the wavelength λ

Since the SiRF model is nonparametric it was, under the approximation of local invariance, fitted to a number N_R of small (2∙∆u_R + 1) × (2∙∆w_R + 1) sized image regions R_r centered at (u_R,r,w_R,r). A region and corresponding variables are shown in Fig. 3. The two keys for selecting the region denoted by index r were capturing the information relevant to the convolution of the SiRF and keeping the edge e_r, approximately in the center of the region. The spatial region size ∆u_R was chosen to be the rounded half of the average distance between two adjacent edge approximations${\tilde{u}}_{0,e,w}{}_{}$:

 

Fig. 3 Illustration of the selected spatial region.

Download Full Size | PDF

The optimization of the main parameters (u_0,r, ρ_r, σ_r, µ_{G r}) was performed by a nonlinear constrained sequential quadratic programming (SQP) [18] algorithm, while the auxiliary parameters c_M,r and c_A,r were computed by a pseudo-inverse. At each iteration of the SQP method an approximation of the Hessian of the Lagrangian function is made using the BFGS update method. Step direction is determined by a quadratic programing solution of a subproblem, which is formed using the approximated Hessian of the Lagrangian. The step size is determined by a line search procedure. The smoothness of SiRF was implicitly enforced by setting M = 3 and imposing the lower limit of the parameters σ_m,r>0.5. The optimization problem was further constrained by solving the Eq. (12) and limiting the weights to nonnegative values (ρ_m,r≥0). Initial values of the edge positions were set to match the corresponding region centers u_0,r = u_R,r. while the parameters σ_m,r, ρ_m,r and υ_m,r were determined by fitting the given edge with simplified Gaussian SiRF model. According to the earlier described region size selection criterion, ∆w_R and ∆u_R were set to 16 and 2, respectively. Optimal parameters for a single region were computed in approximately 0.9 s.

2.3 Global modeling of the computed parameters

The parameters were actually computed only on a mesh of points formed by spectral lines and spatial edges. Therefore, the final step of the characterization is regularization and computation of global tuning function and interpolation functions for the SRF and the SiRF. Starting with the spectral direction, the first step is the definition of features in the image space. This is accomplished by inserting values of the optimized parameters into the constraint equation for each spatial channel, resulting in a 2D grid of points (u,w_0,l,j,u). By assigning the values of the target properties to each point of the grid, three sets of data are formed: one for the spectral tuning function λ(u, w_0,l,j,u), and two for the SRF shape parameters γ_G(u,w_0,l,j,u) and γ_L(u,w_0,l,j,u). By fitting 2D polynomial surfaces to each of the data sets, the interpolation functions (Γ_G, Γ_L) and the inverse spectral tuning function Z_S were computed:

Similar grid was assembled for the image of the spatial reference object by combining the edge positions and their corresponding region centers. Using the true edge positions (θ_0,e) the correspondence pairs θ_0,e(u_0,r, w_R,r) were formed and the inverse spatial tuning function Z_Si was obtained by a procedure analogous to the computation of the Z_S:

2.4 Hyperspectral imaging system and experimental setup

Calibration and validation images were acquired by a hyperspectral imaging system, shown in Fig. 4(a), consisting of a near infrared enhanced MV1-D1312(l)-160-CL-12 camera (Photonfocus) with 1312 × 1082 light sensitive pixels, an ImSpector V10E (from 400 nm to 1000 nm) spectrograph (Specim Ltd, Finland), and S23-f/2.4 lens (Specim Ltd, Finland). The broadband illumination source comprised a research QTH lamp housing (Newport, Oriel Instruments – model: 66881), and a radiometric power supply (Newport, Oriel Instruments – model: 69931). The spectral reference and validation images were recordings of the Hg(Ar), Ar, Ne, and Kr pen spectral reference lamps (Newport, Oriel instruments) models: 6035, 6030, 6032, 6031 respectively. Since a significant portion of the energy emitted by these lamps is in a form of ultraviolet light, optical 525 nm long-pass filter (Edmund Optics, model: NT64-634) was used to attenuate higher energy radiation, preventing hazards and unwanted effects of the light induced fluorescence of organic contaminants (mainly cellulose dust particles). The modulated light was dispersed by a tilted solid thermoplastic diffuse reflectance standard (T, Spectralon, Labsphere) and a thin diffuser (D) placed on the imaging plane. Despite the modeling of the multiplicative and additive effects in Eq. (4), care was taken to ensure a reasonably homogeneous illumination to satisfy the premise of minimal spatial variation of the spectral reference images. By replacing filtered gas emission lamps with broadband light source, and placing the spatial reference object on the imaging plane the system was prepared for the acquisition of the spatial reference image. The employed spatial reference object was a custom laser-machined 150 µm thick stainless-steel mask, shown in Fig. 4(b), with overall manufacturing accuracy of ± 8 µm. Prior to modeling or analysis all the acquired images were corrected by point-to-point subtraction of the dark image and division by the difference between the bright and the dark image. With regards to the position dependent signal-to-noise ratio, the spectral axis was cropped to the spectral range from 527 nm to 980 nm (W = 581 pixel). Likewise, due to the design of the studied HIS, the spatial image dimension was cropped to cover approximately 155 mm wide range (U = 1000 pixel).

 

Fig. 4 Schematic of the experimental setup (a). Spectralon tile (T), long-pass filter (F), diffuser (D). The spatial reference object (b).

Download Full Size | PDF

3. Results and validation of the proposed characterization method

In a perfect HIS we expect that spectra at different spatial channels are mutually aligned and relation between the image spectral positions ω and their corresponding wavelengths λ can be described by a first order polynomial model λ = z_S,1∙ω + z_S,0 + E_S. In such case, the only contribution to the fit error would be the spatially and spectrally independent noise E_S. However, due to the coregistration errors, the first order polynomial model has to be expanded by a term that includes spectral and spatial dependence of the inverse tuning function: λ = z_S,1∙ω + z_S,0 + E_S + δZ_S(υ,ω). The parameters z_S,1 and z_S,0 were assessed by a first order polynomial fit between the true wavelengths of the spectral lines λ_0,l,j and their corresponding locations w_0,l,j,u in all spatial channels of the acquired image (z_S,1 = 0.8004, z_S,0 = 524.3097). Example of the estimated spectral line positions is shown in Figs. 5(a) and 5(b). The position dependent residual δZ_S(υ,ω) is shown in Fig. 5(c). Likewise, the relation between true spatial position θ and the corresponding spatial position in the acquired image υ can be written as θ = z_Si,1∙υ + z_Si,0 + E₀ + δZ_Si(υ,ω). The parameters z_Si,1 and z_Si,0 were assessed by a first order polynomial fit between the true edge positions θ_0,e and their corresponding locations u_0,r in all regions R_r of the acquired image (z_Si,1 = 0.1525, z_Si,0 = 0.6108). Distribution of the extracted edge positions and the residual δZ_Si, are presented in Figs. 5(d)-5(f).

 

Fig. 5 Extracted image positions of two spectral lines (a, b) and the corresponding positions estimated by the computed model, clearly demonstrate subpixel rotation of the image due to misalignment of the camera and the spectrograph slit, and the bending of the spectral lines which is in agreement with previously published work [1]. Furthermore, the spectral axis nonlinearity is evident from spectral dependency of the spectral tuning function residual δZ_S (c). Extracted positions of the two outermost edges at different spectral channels and the corresponding edge positions estimated by the computed characterization model (d, e) suggest rotation and the wavelength dependent magnification. The higher order spatial tuning function residual δZ_Si (f).

Download Full Size | PDF

The SRF, shown in Fig. 6, exhibits a significant variation with wavelength, starting as a Gaussian at short wavelengths and becoming increasingly influenced by the Lorentzian broadening at the near infrared (NIR). The maximum spectral half width at half maximum (HWHM) of the SRF was 2.17 nm assuming an average dispersion of:0.8 nm/pixel. Spatial resolution is the highest in the central part of the image. While the SiRF, depicted in Fig. 7, shows some variance with spatial position, the main variability is with respect to the wavelength. With increasing deviation towards shorter wavelengths, the complex structure and asymmetry of the SiRF becomes evident, while in the opposite direction the symmetry of the SiRF improves, but the resolution is lower. Both can be observed in Fig. 8.

 

Fig. 6 Plots of the Voigt profile parameter interpolation surfaces Γ_L(u,w) (a), and Γ_G(u,w) (b) show the increasing influence of Lorentzian broadening with longer wavelength. Spectral resolution in terms of Half Width at Half Maximum HWHM of the Voigt profile (c). At the NIR end of the utilized spectral range the resolution is significantly lower than at shorter wavelengths.

Download Full Size | PDF

 

Fig. 7 Normalized SiRFs at different image locations, plotted on 32 pixel wide intervals marked by black vertical lines. Best resolution is exhibited in the central part of the acquired image.

Download Full Size | PDF

 

Fig. 8 (a) Asymmetry and complexity of the SiRF are easily observed in plots of the SiRF for spatial channel 742 at two wavelengths corresponding to spectral channels 10 and 570. (b) Quality of the proposed model fit and wavelength dependent shifts of the edges due to coregistration errors are evident from overlay of the acquired (black) and the estimated (red) profiles of two edges in the image of the validation object. Plotted profiles correspond to the two depicted SiRFs.

Download Full Size | PDF

Accuracy of the characterization method was assessed by analyzing the validation image of Kr spectral reference source. Firstly, the forward tuning function and the SRF model were validated by predicting the acquired validation image, shown in Fig. 9(a), by the computed characterization model. Secondly, the positions of four well resolved spectral lines in the acquired validation image were extracted by fitting the SRF model to the acquired intensity profiles of the emission lines. The extracted positions w_V,j,u of j-th spectral line at spatial channel u were then transformed into corresponding wavelengths λ_V,j,u by the inverse tuning function. The wavelength prediction error was calculated as δλ_V,j,u = λ_V,j,u-λ_V,0,j, where λ_V,0,j are the true wavelengths of the spectral lines. For comparison, additional spectral calibration was performed according to the method proposed by Polder et.al [7]. The spectral line pixel positions were extracted from each spatial channel of the acquired Hg(Ar) reference lamp image. The relation between pixel locations and the true wavelengths was modeled by a third order polynomial fit. For the purpose of validation, the extracted pixel positions of the four well resolved Kr spectral lines were transformed by the calculated polynomial fit and compared to the corresponding true wavelengths. The prediction error statistics of the two methods for each of the four selected Kr spectral lines is presented in Table 1.

 

Fig. 9 Accuracy of the forward tuning function and SRF estimation are presented by comparison of the acquired image and the estimated image of the spectral and spatial validation objects. Estimated spectrum in the central part of the Kr validation image and the fit error expressed as a difference between the estimated and acquired spectrum (a). Estimated image profile in the central part of the spatial validation image and the fit error computed by subtracting the estimated and the acquired profile (b).

Download Full Size | PDF

Table 1. The wavelength prediction error statistics of the four Kr spectral lines.

View Table

Likewise, the spatial validation image shown in Fig. 9(b), was produced by imaging the reference object displaced by a known shift. Quantitative assessment was made by finding the image-wise edge positions and calculating their expected locations in the object space. Analogically to the spectral direction, the 30 equally distributed edges covering the entire spatial region were located by locally fitting the modeled edge response to each of the 581 spectral channels resulting in a total of 17430 validation points. The residual misalignments, shown in Fig. 10, were computed as a difference between the expected and the true edge positions in the object space. The mean residual misalignment for all the validation points was 0.001 mm, with the standard deviation of 0.005 mm. Speaking in pixel-size terms, the mean error was 0.007 pixel and the standard deviation was 0.05 pixel. The validation shows that the expected error was in the range of 1/10 of a pixel, which is comparable to the spatial reference object uncertainty specified by the manufacturer.

 

Fig. 10 Residual misalignment between the measured and the true edge positions. Centers represent the mean error and the error bars correspond to the standard deviation. Dashed line denotes the reference object uncertainty. The edge locations correspond to the edges in the profile shown in Fig. 9(b).

Download Full Size | PDF

The computed forward tuning function can be used for calibration of the acquired image by means of geometric image transform [19]. First, a pixel grid (θ_q, λ_q) is formed in the object space. Each pixel q is then transformed by the forward tuning function into the acquired image space (υ_q,ω_q) = (T_S(θ_q, λ_q)_, T_Si(θ_q, λ_q)). Finally, the intensities of the object space image are computed by bicubic interpolation. Results of the validation image calibration are displayed in Fig. 11.

 

Fig. 11 The effects of the calibration observed in the spectral and spatial validation images. For better visualization of the small coregistration errors only spectra and spatial profiles limited to image subregions are shown. (a) Spectra of the 760.15 nm and 769.45 nm Kr spectral lines at spatial channels 10 and 990 corresponding to true spatial positions 2.07 mm and 152.16 mm. Due to coregistration errors, the spectra in the acquired image are shifted by approximately 0.25 pixel. (b) After transforming the acquired image into the object space the image is calibrated to real physical units and the coregistration improves. Similar observation holds for the spatial validation image. (c) Plot of the three spatial validation image profiles of a single edge at spectral channels 17, 282 and 555 corresponding to wavelengths 540 nm, 750 nm and 970 nm clearly shows wavelength dependent spatial distortions. (d) After the calibration the profiles are spatially aligned.

Download Full Size | PDF

4. Conclusion

Prior work has demonstrated the significant impact of the errors arising from poor uniformity of the acquired images on the accuracy of the subsequent hyperspectral image analysis. Various methods were proposed for the assessment and/or correction of the effects arising from the operation principle of the spectrograph based hyperspectral imaging system. However, the majority of these methods require equipment often not present on the application site or the characterization is based on a small number of image locations, effectively ignoring the spectral-spatial dependency. Moreover, many of the spectral lines emitted by the standard reference sources are unresolvable by the moderate spectral resolution of the imaging spectrographs. Consequently, the simple estimation of the spectral line position by its center of gravity is inaccurate. This significantly limits the number of well resolved spectral lines useful for the calibration, and thereby the degree of modeling, and consequently the accuracy of the fitting. We have alleviated this problem by using spectral matching method, which allowed taking into account all the significant spectral lines in the observed range and the spectral resolution of the device. Furthermore, a simple visual inspection suggested asymmetry of the spatial response function (SiRF). This was an important finding, as the edges extracted by the popular gradient based edge detection algorithms would be inevitably biased. Therefore, we proposed an extension of the SiRF model allowing greater flexibility. The results of the proposed characterization method can be used for calibration of the acquired images and/or as a tool for quantization of the expected errors in various hyperspectral imaging techniques. Furthermore, we expect, the proposed characterization of the SRF and SiRF can be used in future work for development of hyperspectral image restoration methods based on previously published work on hyperspectral image deconvolution [20, 21].