VIEWS: 94 PAGES: 28 POSTED ON: 12/14/2011 Public Domain
Image Deconvolution Applied to CCD astronomical images with astrometric purposes Maria Teresa Merino Espasa Universidad de Barcelona (Spain) Summary Introduction to Image Restoration Adaptative deconvolution with wavelets Applications of our deconvolution algorithm with astrometric purposes: Increase of faint detections in Baker-Nunn images for faint asteroid discovering Quest aplication to increase astrometrical resolution to detect gravitational lenses candidates Work in progress to increse faint detection + increase astrometrical resolution Conclusions 1. Introduction to image Restoration Objective: Presentation of our algorithm for MLE deconvolution in linear space-invariant PSF systems with complete basic CCD images statistics (Poisson + Gaussian) Image restoration Our detected data images are not exactly the radiation from sources of interest. It’s important to consider that may be affected by: atmosphere effects optical distorsions and aberrations different pixel sensitivity additional noise (in emision, travelling and/or detection) … Image Restoration is the removal or reduction of those degradations that were incurred while the digital image was being obtained to bring the image toward what it would have been if it had been recorded without degradation. In astronomy, the most usual objectives are: deblurring removal atmospheric seeing radiation image sharpening fusion of data from diferent instruments increase S/N of sources astrometric resolution improvement … Linear image formation system with additive noise p(x,y): projection (measured) data p( x, y) g ( x, y) n( x, y) g(x,y): blurred image n(x,y): additive noise g ( x, y ) h( x, y; x1 , y1 )·d ( x1 , y1 )dx1dy1 h(x,y): blurring effects d(x,y): source+background radiation When space-invariant PSF: g ( x, y ) h( x x1 , y y1 )·d ( x1 , y1 )dx1dy1 h( x, y ) d ( x, y ) p( x, y) h( x, y) d ( x, y) n( x, y) In Fourier Space: P(u, v) H (u, v)·D(u, v) N (u, v) Deconvolution ill-posed problem Deconvolution: Image restoration based on some knowledge of the degradation processes effects and/or the expected source distribution using models of the statistical and blurring effects of the image formation system. ……...? P(u, v) N (u, v) Why not….. D(u, v) H (u, v) H (u, v) 1. at those points the division Because: operation is undefined or results in H(u,v) has zero o near zero meaningless values !!!! values in most of (u,v) range 2. for points having very small |H(u,v)|, although the division can be done, the noise will be amplified Conclusion: to an intorelable extent !!!! Image deconvolution is an inverse solution ill-posed problem !!! All the algorithms which try to solve it should be concerned about obtaining a solution with the three following properties: • existence • uniqueness • stability Deconvolution algorithms 1. Image formation, degradation and recovering model basic hypotesis Classical approach: • Direct inverse filter in Fourier space • Linear Regularized Filters (Wiener) • Minimum root mean squares • CLEAN Bayessian approach: • Fixed Model: MLE, MAP • Variable Model: Pixons 2. Regularization hypotesis (positiveness, energy conservation, frequency bounders, …) 3. Convergence iterative algorithm (gradient methods, EM, successive substitutions,…) 4. Iteration Stopping protocol (a posteriori, feasibility tests, crossvalidation maximum,…) The Bayesian Approach (fixed models) P(a|p): “a posterior probability”, radiation a given the recorded data p ( p | a ) ( a ) P(a): “a prior probability” of the source of ( a | p ) radiation ( p ) P(p|a): likelihood P(p): Normalizing constant Bayes Hypotesis: The most probable image a(x,y), given measured data p(x,y) is obtained by maximizing P(a|p) with a suposed P(p|a) image formation model (and maybe also a P(a) “a priory” knowledge). The most usual deconvolution algorithms division regarding how they maximize P(a|p) are: • MAP (Maximum A Posteriori), maximizes P(p|a)·P(a) • MLE (Maximum Likelihood Estimator), maximizes P(p|a) with P(a)=cte • MEM (Maximum Entropy Method) is MAP with p(a) given by maximum entropy hypotesis (=minimum data uncertainty) G+P_MLE deconvolution • sucessive sustitutions convergence method Algorithm B B D f a D D • energy conservation: qi ai p j b j ji i basics: i 1 i 1 j 1 Cj j 1 j 1 • Poisson photon gathering+Gauss readout noise: hj hj k ( p j k )2 B f 1 h j ij ·ai b j P k hj e P( p j | k ) e 2 2 k! 2 i 1 C j m D f ji p ' j (s) 1 ( s 1) ai ai qi f jl al C j b j B j 1 (s) l 1 D ( k p j ) 2 (h j ) k Iterative and nonlinear algorithm ( p bj ) ke j 2 2 j 1 (s) k! m k 0 p' j ( k p j ) 2 NOTE: This resulting algorithm, B D f ji p' j (h j ) k without flat and background, is a i B (s) e 2 2 k! similar to MLE Poisson (Richardson Lucy) computation i 1 j 1 l 1 f jl al C j b j (s) k 0 structure if p’=p and K=1 2. Adaptative Deconvolution with Wavelets Objective: Presentation of our multiresolution adaptative deconvolution algorithm based on wavelets Wavelets descomposition The wavelet transform replaces the Fourier transform’s sinusoidal waves by a family generated by translations and dilations of a function called Mother Wavelet. - b space and a scale - a ,b ( x ) 1 WT ( f (a, b)) a 2 f ( x) ( x b ) dx a functions base Discrete wavelet transform (DWT) algorithm specific properties: DWT is computed with à trous algorithm. Mother wavelet function is derived from a B3 cubic spline scaling function. This is a dyadic decomposition (scale 2i , i=1,...,n). No subsampling is applied wavelet planes have same size as original image. General properties of wavelets descomposition: Spatial and frequential content is almost decoupled. Better noise vs.signal discrimination than Fourier transform. Good processing flexibility. Example of Wavelets descomposition c0 w1 w2 w3 w4 c4 Original Residual plane image Wavelet planes at scale 4 co = w1 + w2 + w3 + w4 + c4 High frequency Low frequency Adaptative deconvolution with wavelets To deconvolve the image, this is decomposed in wavelet planes at every iteration of the method. The signal detection masks only applies deconvolution to those pixels with significant signal. Objectives of using wavelets: Multiescale deconvolution Objectives of adaptative deconvolution: - Do not amplify noise. - To detect and deconvolute only those features with signal. Signal detection mask: 3 2 w nf v, j v, j w vh,(js ) r p' 2 1 exp 2 v, j if vr, j 0 v, j t mv , j 2( vr, j ) 2 v, j nf with noise std. desv. at 0 if v, j r v, j 0 v, j r wavelet plane n sorrounding pixel j m Nw f ji w vh,(js ) mv , j w v , jj w vh,(js ) p' D (s) 1 G+P_AWMLE : ai ( s 1) ai qi v B j 1 f jl al C j b j (s) l 1 3. Applications of our deconvolution algorithm with astrometric purposes: Objective: Show astrometric objectives, aplications and deconvolution artifacts and effects Deconvolution of Baker-Nunn (RAO) images Data overview -50cm - f/1 Baker-Nunn instrument!!! -very large FOV 4.2ºx4.2º, -moderately undersampled data: 3.9”/pixel, -significant object blending, -bright stars cause overblooming, -PSF is low correlated with seeing and dominated by systematic distortions across FOV, which can be modelled -due to large pixel size, background level is dominated by Poisson noise. 0.8degx0.8deg subframe: note object blending and Deconvolution objective: blooming around bright stars Increase limit magnitude to detect faint RAO: Rothney Astrophysical Observatory NEOs with minimum false detections University of Calgary (Canada) Protocol for Baker-Nunn images 1) CCD image calibration: bias, dark and flatfield correction iraf.ccdproc 2) Compute astrometric plate transformation with reference catalogue wcstools.sua2, sextractor,iraf.ccmap 3) Perform aperture photometry of selected stars (bright but not saturated) iraf.phot, iraf.pstselect 4) PSF fitting with list created in 3) iraf.psf, iraf.seepsf 5) Image deconvolution 6) sextractor object detection and matching with USNOA2.0 catalogue 7) Compute limiting magnitude gain and rest of analysis Deconvolution of Baker-Nunn images Benefits of image deconvolution 1) Increase of SNR Limiting magnitude gain 2) Increase of resolution Object deblending Original image: in blue matched Deconvoluted image 60 iterations: detections yellow and red are new matched detections Deconvolution of Baker-Nunn images Matched vs. non-matched detections Evolution of raw and matched Adaptative wavelet deconvolution detections with the number of iterations. introduces small number of non-matched Asymptotic converge translated into detections. stable number of detections for large In contrast to classical Richardson-Lucy number of iterations. algorithms, false detections increase linearly, not exponentially. Deconvolution of Baker-Nunn images Where do the non-matched detections come from? Left: original image with explanation of the origin of non-matched detections displayed in red in the right frame, which corresponds to a 60 iteration deconvoluted image. Deconvolution of Baker-Nunn images Limiting magnitude gain Image # matched detections Original 1559 Deconvoluted 2365 180 iterations But deriving magnitude gain directly from ratio of the number of detections is not allways reliable, due to intrinsic magnitude distribution of studied stellar field. From estimating the completeness of magnitude histograms with respect to USNOA2.0 distribution: USNOA2.0 R magnitude distribution for the catalogue itself (top), deconvoluted image 180 iterations (middle) and original mR 0.6 image (bottom) Deconvolution of QUEST images QUEST= QUasar Ecuatorial Survey Team Data overview QUEST images to deconvolve WIYN images to evaluate results V filter, frames from 2 CCDs 2 MiniMosaic WIYN 1”/pixel, 0.141”/pixel 2 diferent fields in 2 single CCDs Covers aprox. the same field Limiting magnitude (S/N>10) 6+38 QSO candidates detected by aprox. 19.2 varability criteria in those fields Deconvolution objective: Detect a subsample of QSO candidates likely to be lensed by: 1. improving the resolution (deblend close companions) 2. improving the limiting magnitude Deconvolution of QUEST images Aplication of image deconvolution to 6 QSO candidates in the field. Each panel includes original QUEST image, deconvolved QUEST (400 iterations) and high resolution WIYN image, repectively. Those in green correspond to objecte present in all three images, in black those only resolved by WIYN and in red those resolved in both deconvolved QUEST and WIYN. Algorithm properties Powerful technique for increasing the number of useful science objects from the faint part of magnitude distribution. Enhances spatial resolution which translates into better object deblending discrimination. Generalizable for whatever observing facility employing a CCD. Cheap: only requires fast computers. The algorithm is highly parallelizable to be run in distributed computers. Low speed convergence and heavy CPU and RAM usage. Systematic usage best suited for astrometric surveys with moderate data throughput. Selective usage is suited for all kind of images, either high or low resolution. Work in progress Use multi frame reduction for removing cosmics and hot pixels. Improve PSF fitting to decrease false detections around bright stars. Compare internal astrometric accuracy of original and deconvoluted images. Accurate determination of CCD flats, gain and readout noise for better noise distribution modelling when deconvoluting. Define strategy to assess space variant PSF across the FOV Further evaluation of astrometric improvements and photometric effects. To sum up: In this presentation we have presented our selected algorithm: Gausian and Poison Adaptative Wavelets Maximum Likehood Estimator This algorithm has been proved with two diferent astrometrical type of data and diferent scientific goals Those examples have shown how our algorithm can produce in the particular case of astrometry of point sources: • recover faint objects • improve resolution Final Conclusions Deconvolution can be a powerful technique in astrometry projects. But it’s crucial to select the proper algorithm carefully considering: 1. your data 2. your scientifical objectives