Scene-Based Nonuniformity Correction Method Using the Inverse Covariance Form of the Kalman

Scene-Based Nonuniformity Correction Method Using the Inverse Covariance Formof the Kalman FilterSergio Torres I.(e-mail:storres@die.udec.cl),Jorge Pezoa N.(e-mail:jpezoa@die.udec.cl) Department of Electrical Engineering,Universidad de Concepci´o n,CHILE Abstract A scene-based algorithm for nonuniformity correction(NUC)in focal-plane arrays(FPA)detectors has been developed.The NUC technique is based in theinverse covariance form(ICF)of the Kalmanfilter.The gain and the offset of eachdetector of the FPA are modeled by discrete-time Gauss-Markov processes.These pa-rameters are taken as constant within a given sequence of frames,corresponding toa certain time and operational conditions,but they randomly drift from one sequenceto another in response to new operational conditions.For each detector and each se-quence of frames,the ICFfilter input is an observation vector consisting of detector’sread-out values.The output of ICFfilter for any sequence of infrared frames is thedetectors’gain and offset.The efficacy of the ICF of the Kalmanfilter to compensatefor nonuniformity noise in infrared imagery is demonstrated using sequences of in-frared imagery with both artificial nonuniformity and artificial drift in the detectors’parameters.It is shown that the ICFfilter and the Kalmanfilter generate similar re-ductions of nonuniformity.However,the ICFfilter compensates the noisy images withless number of operations per pixel and per frame than the Kalmanfilter.1IntroductionFocal-plane array(FPA)detectors are used in the visible and infrared spectral bands imaging systems for a variety of applications.A FPA consist of a two-dimensional array of photode-tectors placed in the focal plane of the imaging lens.The performance of a FPA is strongly affected by the nonuniformity response of the detectors.The problem of nonuniformity in array detectors is attributed to the fact that each detec-tor in the array has a different responsivity to light.Thus,for an imaging sensor viewing a spatially constant scene,the sensor output may not be constant from pixel to pixel(as one would hope).This nonuniformity is also referred to asfixed-pattern noise and it may severely corrupt the images.In severe cases,the result is like looking at the world through a“dirty window”.As one might expect,thefixed-pattern noise reduces the resolving capability of an imaging system.For example,it may not be possible to capture minute pixel-to-pixel tem-perature variations.In addition,thefixed-pattern noise causes the performance degradation of many frequently used post-processing operations which are typically employed to identify object motion.What makes matters worse is that thefixed-pattern noise slowly changes over time(in the course of minutes to hours in some cases)preventing a simple one-time laboratory calibration.Of course,one can calibrate frequently;however,this requires sophisticated and expensive calibration targets(i.e.,black-body radiation sources)and more important,requireshalting normal imaging operations during the calibration process.In contrast to calibration-based nonuniformity correction(NUC),scene-based NUC is the process of compensation for the nonuniformity by using only the scenes frames being imaged.Hayat’s group[1–9]has been actively pursuing the development of novel algorithms for NUC based on statistical estimation theory.In particular,they developed a Gauss-Markov model for the slow variation in thefixed-pattern noise and utilized the model to estimate the nonuniformity parameters using a Kalmanfilter.Generally speaking,Kalmanfilters are com-putationally efficient estimation techniques that facilitate the calculation of a current estimate from the past estimate and current data.In the context of array detectors,blocks of frames (each containing hundreds of frames)separated by long time intervals(several minutes,in practice)are used to estimate the nonuniformity parameters.The nonuniformity parameters slightly vary from block to block.As new blocks arrive,the nonuniformity parameters corre-sponding to each new block are estimated by updating the old parameters(corresponding to the old block of frames).In this way,the valuable information contained in the old estimates are preserved and efficiently used in forming the current state of nonuniformity.In this paper,we propose a scene-based NUC method which optimally estimate the gain and the offset of each detector in the FPA by means of the inverse covariance form(ICF) of the Kalmanfilter.The ICF of the Kalmanfilter provides an algorithm that is both math-ematically equivalent to the conventional Kalmanfilter algorithm and computationally more efficient,when the dimension of the measurements is much greater than the dimension of the process noise.The matrix inverse operation is a very time consuming operation,and in the ICFfilter we compute the inverse of a matrix with the dimensions of the noise process, while in the Kalmanfilter we compute the inverse of a matrix of the dimension of the obser-vations.Further,the inverse covariancefilter is better suited for problems where there exists no previous knowledge of the initial system state.This paper is organized as follows.In section2the model of the system and the origi-nal Kalmanfilter are presented,and the ICF is developed following the standard procedures given in[12].In section3the NUC technique is tested with sequences of infrared data with simulated nonuniformity and drift in time,and three performance parameters are computed to measure the performance of the method.In Section4,we present the main conclusions. 2Estimation of the Gain and Offset Using the Inverse Covariance Form of the Kalman FilterWhen the FPA is operating in its linear input-output range,every detector may be modeled using a irradiance-linear model for the detector response,presented in[10]by Holst:Y ij(t)=X(1)ij (t)T ij(t)+X(2)ij(t)+V ij(t)(1)where the subscripts ij mean the position of the pixel in the array and t represents thetemporal variation of the parameters in the model,X(1)ij (t)and X(2)ij(t)are ij-th the detector’sgain and offset,respectively,at any time t.T ij(t)represents the ij-th input irradiance,at any time t,over the detector during integration time,V ij(t)is the ij-th additive readout noise and Y ij(t)is the ij-th readout signal for the detector response,both at any time t.From now on, the pixel subscripts ij will be omitted with the understanding that all operations are performed on a pixel-by-pixel basis.2.1The Kalman FilterTorres et al,[5][7],developed a Kalmanfilter that estimates the gain and the offset of each detector in the FPA.To obtain thefilter,the system states variables(gain and offset)were modeled using a Gauss-Markov random process.This model effectively captures the slow temporal variation inherent in the gain and the offset observed in many practical applications.The state equations for the Gauss-Markov process at the k-th block time are given by:X(1)ij(k+1)X(2)ij (k+1)=αk00βkX(1)ij(k)X(2)ij(k)+1001W(1)ij(k)W(2)ij(k)(2)X k+1=Φk X k+G k W k(3)again,X(1)(k)and X(2)(k)are detectors’gain and offset,W(1)(k)and W(2)(k)are the driving noise for the gain and the offset at the k th block.αk andβk parameters represent the amount of drift,in gain and offset,between k an k+1block time.Values ofαk andβk near to 1generate slow drift and values near to0produce high drift in gain and offset.Based on the state equation(3),on the extension of the observation equation(1)and assuming that for each frame in a block,the collected irradiance is a uniformly distributed random variable in some range[T min,T max],common to all detectors and all time,the fol-lowing Kalman Filter was derived[5,7]:ˆX−k =Φk−1ˆX k−1+M Tk−1(4)P−k =Φk−1P k−1ΦTk−1+G k−1Q k−1G Tk−1(5)K k=P−k H Tk[H k P−kH Tk+S k]−1(6)S k=R k+σ2T (σ2X(1)+X(1))I lk,l k(7)ˆX k =ˆX−k+K k(Y k−H kˆX−k)(8)P k=(I2,2−K k H k)P−k(9) with initail conditions:ˆX 0=E(X0)=X(1)X(2),P0=Λ=σ2X(1)0σ2X(2)(10)Equations(4)y(5)are the Kalmanfilter time updates,(6)is the Kalman gain and equa-tions(8)and(9)are the measurement updates.ˆX−k is the a priori state estimate,ˆX k is thecurrent state estimate.P−k is the a priori error covariance matrix,P k is the current error co-variance matrix.K k is the Kalman gain,Φk−1is the state transition matrix,Q k−1and R k−1 are the auto and cross covariance functions of the driving noise and additive noise,respec-tively.G k−1is the system associated matrix and Y k is the vector of observations.H k is the mean of the observation matrix.M k is a vector containing the mean of the driver noise,σ2Tis the variance of the inputinfrared irradiance,X(1)0andσ2X(1)0are the mean and variance of the initial condition of thegain,respectively.I lk,l k is an identity matrix with dimensions of the block length.The ma-trixes R k−1,S k−1are a diagonal square matrix of dimension of the observation vector.Table1:Number of operations,per pixel and per block of frames,required for one iteration.n represents the dimension of the system state variables,m the dimension of the observations and p the dimension of the process noiseKalman Filter Inverse Covariance FormAdds3n3+n2(2m+p−2)+4n3+n2(4p−m−1)++n(m−1+3m2+p2+2p)+m3+n(2p2−3p−m+2m2−1)+m+m2+p3 Multiplies3n3+n2(5m+p+2)+4n3+n2(4p+m+3)++n(4m+2m2+p2+p)+m3+n(2p2+m+2m2)+2m+p32.2Derivation of the Inverse Covariance Form of the Kalman FilterThe ICF or information matrix form of the Kalmanfilter is an alternative method of thefilter that is particularly well suited to problems where the measurement dimension is large or where there exists no a priori knowledge of the initial system state.In NUC methods for FPAs,the number of observations(l k),per detector,normally are between300to2000frames[5,7],so the dimension of the vector of observations is much greater than dimension of the noise and system states,the gain and offset in our case.Table 1shows the number of operations required for the k th iteration for the Kalman Filter and for the ICF of thefiing a block length of500frames the Kalmanfilter needs126754038 adds and126014052multiplies while inverse covariancefilter calculates1247570adds and 1004100multiplies.We can see that the inverse covariancefilter provides a computationally more efficient algorithm than the conventional Kalmanfilter.Further,for estimation problems where there exists no a priori knowledge of the initial system state,this algorithm is less susceptible to saturation problems within the equations and in this sense may provide a numerically more stable approach[12].The ICF is an alternative form of the Kalmanfilter in which the inverse of the errorcovariance matrix,P−1kand the estimate:ˆa k P−1kˆX k(11) are propagated in each iteration of thefilter,rather than the error covariance matrix and the estimateˆX k.The estimate of the system state,ˆX k,can be recovered fromˆa k by multiplying ˆa k by the matrix P k.The derivation of thefilter was made following the standard procedure given in[12],con-sisting in the application of the matrix inversion lemmas:(L+MN T)−1=L−1−L−1M(I+ N T L−1M)−1N T L−1and(L+MN T)−1M=L−1M(I+N T L−1M)and the combination of the Kalmanfilter equations using the definition11.The equations for the ICF of the Kalmanfilter are:ˆa−k ={I−B k−1G Tk−1}{Φ−1k−1ˆa k−1+A k−1M Tk−1}(12)(P−k )−1=(I−B k−1G Tk−1)A k−1(13)A k−1 Φ−T k−1P−1k−1Φ−1k−1(14)B k−1 A k−1G k−1{Q−1k−1+G T k−1A k−1G k−1}−1(15)P−1 k =(P−k)−1+H TkS−1kH k(16)ˆa k=ˆa−k +H TkS−1kY k(17)The ICF of the Kalman filter is given by the the time updates (equations (12)and (13)),the measurement updates (equations (17)and (16))and the initial conditions:P −1k =Λ−1and ˆa 0=Λ−1X 0In cases when there is no a priori knowledge of the initial system state,the inverse covari-ance filter can be utilized with the initial condition P −1k =0,corresponding to a P 0of infinity,situation that leads the traditional filter to saturation problems and consequently a significant and possibly catastrophic loss in numerical accuracy.Furthermore,with the Kalman filter,we must invert the matrix [H k P −k H T k +S k ]−1,how-ever with in the ICF the matrix to invert is [S k ]−1,both have a dimension of the observation vector,but S k −1is a diagonal matrix,so the computing time is considerably reduced.3Applications to Simulated Infrared Data 00.050.10.150.20.250.30.3500.10.20.30.40.50.60.7∆σGain ρRoughness Parameter 02040608010012014016018020000.20.40.60.811.21.41.61.8∆σOffset ρRoughness ParameterFigure 1:Roughness parameter for Kalman and ICF filter versus variations in the gain (left)and the offset (right)standard deviations,respectivelly.Solid (dotted)lines represents low (high)drift in parameters.On both pictures,the upper line is the roughness parameter of the uncorrected block.For the right graph,the middle lines are roughness obtained for ICF filter and lower lines roughness for Kalman filter.The performance of the inverse covariance filter is studied and compared with the perfor-mance of the traditional Kalman filter using infrared image sequences that are corrupted by simulated nonuniformity.NUC is performed by simple subtracting the estimated offset from the data and dividing the outcome by the estimated gain.To study the performance,we use the mean-square error (MSE)for the estimated gain and offset,averaged over all detectors.NUC capability is examined in terms of the root-mean-square error (RMSE)and roughness parameter (ρ)for a corrected image,metrics commonly used in NUC [2,10,11].All the assumptions for implementing the simulations are common for testing NUC methods and are given in [1–7,11].The inverse covariance filter assumed an unknown initial condition,P −1o =0.3.1Performance on a Block of Infrared FramesFigure (1)shows the computed roughness parameter,the upper line in both graphics is the roughness of the uncorrected image.In can be seen that,under variations in gain,even under00.050.10.150.20.250.30.3500.050.10.150.20.25M S E G a i n MSE of Gain∆σGain M S E O f f s e t 02040608010012014016018020000.10.20.30.40.50.60.7M S E G a i n MSE of Gain∆σOffset M S E Figure 2:Computed MSE of the estimated gain (left)and offset (right)for the Kalman and ICF filter versus variations in the gain and offset standard deviations.Solid (dotted)lines represents low (high)drift in parameters,’o’sign means Kalman filter and ’+’sign means ICF filterhigh drift in the parameters,there is a great reduction in the nonuniformity and it’s similar for both filters.However,under offset’s variations,the inverse covariance filter obtain a lower performance,situation that can be enhanced taking a larger number of frames.Figure (2)shows the mean-square error of the estimated gain and offset.Again,similar behavior for the estimation are obtained under variations in gain and offset,but we must note that,under off-set’s variations,the MSE of the gain is greater for the inverse covariance filter.However the offset’s MSE is approximated the same for both filters.In Figure (3)we can see the computed root mean-square error for the corrected block of frames.Note that again the inverse covari-ance performance is less than traditional filter for the case with offset nonuniformity.As an example,figures (4)and (5)show the true image,the true image with simulated nonunifor-mity and the ICF corrected image with mainly gain nonuniformity and with mainly offset nonuniformity,respectively.00.050.10.150.20.250.30.350.050.10.150.20.250.30.350.4∆σGain R M S E RMSE 0204060801001201401601802000.10.150.20.250.30.350.40.450.50.550.6∆σOffsetR M S E RMSE Figure 3:RMSE of the corrected image for traditional and inverse covariance filters under variations in the gain (left)and the offset (right).Solid (dotted)lines represents low (high)drift in parameters.Figure4:The True image(left),the True image with artificial nonuniformity(middle)and the ICF Corrected (right)frame.The nonuniformity is mainly generated for the gain.Figure5:The True image(left),the True image with artificial nonuniformity(middle)and the ICF Corrected (right)frame.The nonuniformity is mainly generated for the offset.Figure(6)shows the CPU time consumed by the Kalmanfilter and ICFfilter,versus the block length.By reducing the number of operations,the CPU time for the ICFfilter has been reduced several times compared with the Kalmanfilter for any frame size and block lengths greater than200frames.For example,reductions of55%in time is obtained for block lengths of1000frames.The tests were made using a Pentium IV1.6GHz processor and a 768MBytes RAM and using the Matlab’s function cputime.3.2Performance on a Set of Blocks of Infrared FramesThe ICF of the Kalmanfilter is evaluated using several blocks with simulated nonuniformity patterns at different levels of drift and each block with afixed length of500frames.The results of the tests are presented in Table2.The RMSE parameter is approximated the same for each k-th block,independent of the values ofαandβ(i.e.,for high,moderated and low levels of drift).Further,the MSE of the estimated gain and offset prove that the ICF of the Kalmanfilter can effectively estimatimates the system parameters,independent of the values ofαandβ.This demonstrate that the ICFfilter is capable of taking advantage of information contained in a previous block of frames and effectively updating this information using the current block.100200300400500600700800900100002468101214161820BlockLenght l k C P U T i m e CPU Time EvaluationFigure 6:The CPU time consumed by the traditional Kalman and ICF filter for frames of 32x32(dotted line),64x64(solid line)and 128x128(dashdot line)pixels,’o’sign means Kalman filter and ’+’sign means ICF filterTable2:Performance of the inverse covariancefilter estimating states on several blocks.k=1k=2k=3k=4k=5k=6k=7k=8k=9k=10α=β0.950.100.990.500.300.960.800.400.200.30 MSEˆX(1)0.01430.01520.01390.01820.01510.01450.01700.01630.01360.0178 MSEˆX(2)86.916499.243898.0473102.3034104.955097.555095.118492.403999.995298.6594 RMSE0.83120.84040.87710.82480.82440.83130.82270.82500.82020.85004ConclusionsWe have shown that the ICF of the Kalmanfilter is an alternative algorithm for the Kalman filter that is both mathematically equivalent to the Kalmanfilter and computationally more efficient for NUC on FPA.This mainly is because the number of observations needed to com-pensate corrupted block of frames is much greater than the number of system state variables. Further,the inverse covariancefilter is better suited for problems where there exists no pre-vious knowledge of the initial system states,fact well suited for practical situations in NUC applications.It was also shown that the NUC performance of the inverse covariancefilter is similar to the Kalmanfilter.The roughness parameter,the MSE for estimated states variables and the RMSE for the corrected frames are very similar to those obtained using Kalmanfilter.In offset estimation,the proposed ICF results in low performance than traditionalfilter,but with minor increments in the length of the block or using sampled frames,the performance of the filter can be improved.Thefilter,as traditional Kalmanfilter,is capable to capture the drift in nonuniformity.The improve of computational efficiency using the ICF of thefilter was demonstrated comparing bothfilters and significative reductions in time are obtained for block lengths greater than200frames suiting the algorithm as an excellent option to be calculated on line. 5AcknowledgmentsThis work was supported by the‘Fondo Nacional de Ciencia y Tecnolog´ıa’FONDECYT of the Chilean government.The authors thanks Ernest E.Armstrong at the United States Air Force Research Laboratory,Wright-Patterson Air Force Base and Majeed M.Hayat at University of New Mexico,USA for many valuable comments.References[1]M.Hayat,S.Torres,E.Armstrong,Model Base Real-Time Nonuniformity Correction in Focal Plane ArrayDetectors,Proc.SPIE,vol.3377,1998.[2]M.Hayat,S.Torres,E.Armstrong,B.Yasuda.Statistical Algorithm for Non-Uniformity Correction inFocal Plane Arrays,OSA Applied Optics,vol.39,2000.[3] E.Armstrong,M.Hayat,R.Hardie,S.Torres,B.Yasuda,Non-Uniformity Correction for Improved Regis-tration and High Resolution Image Reconstruction in IR Imagenary,Proc.SPIE,vol.3808,1999.[4]S.Torres,M.Hayat,E.Armstrong,B.Yasuda,On the Performance Analysis of a Recent Statistical Algo-rithm for Non-Uniformity Correction in Focal Plane Arrays,Proc.CISST,vol.3703,2000.[5]S.Torres,M.Hayat,E.Armstrong,B.Yasuda,A Kalman-Filtering Approach for Non-Uniformity Correc-tion in Focal-Plane Array Sensors,SPIE vol.4030,2000.[6]R.Hardie,M.Hayat,E.Armstrong,B.Yasuda,Scene Based Non-Uniformity Correction Using VideoSequences and Registration,OSA Applied Optics,vol.39,2000.[7]S.Torres,M.Hayat,Compensation for Gain and Bias Nonuniformity and Drift in Array Detectors:AKalman-Filtering Approach,Submitted to IEEE Trans.Image Processing,Feb.2001.[8]S.Cain,M.Hayat,E.Armstrong,Projection-Based Image Registration in the Presence of Fixed-PatternNoise,IEEE Trans.on Image Proc.,vol.10,num.12,2001.[9] B.Ratliff,M.Hayat,R.Hardie,An Algebraic Algorithm for Nonuniformity Correction in Focal-PlaneArrays,Proc.SPIE,vol.4372,2001.[10]G.Holst,CCD Arrays,Cameras and Displays,SPIE Optical Engineering Press,1996[11]M.Schultz,L.Caldwell,Nonuniformity Correction and Correctability of Infrared Focal Plane Arrays,Proc.SPIE,vol.2470,1995.[12]G.Minkler,J.Minkler,Theory and Application of the Kalman Filtering,Magellan Book Company,1993.[13] D.Scribner,M.Kruer,J.Killiany,Infrared Focal Plane Array Technology,Proc.of the IEEE,vol.79,1991.[14] ton,F.Barone,M.Kruer,Influence of Nonuniformity on Infrared Focal Plane Array Performance,SPIE Opt.Eng.,vol.24,1985.。