function atrous, spec, filter = filter, n_scales = n_scales, check = check ;+ ; NAME: ; ATROUS ; PURPOSE: ; Perform a 1-D "a trous" wavelet decomposition on an image. ; ; CALLING SEQUENCE: ; decomposition = ATROUS( spectrum [, filter = filter, \$ ; n_scales = n_scales, /check] ; ; INPUTS: ; SPECTRUM -- A 1-D Image to Filter ; ; KEYWORD PARAMETERS: ; FILTER -- A 1D (!!) array representing the filter to be used. ; Defaults to a B_3 spline filter (see Stark & Murtaugh ; "Astronomical Image and Data Analysis", Spinger-Verlag, ; 2002, Appendix A) ; N_SCALES -- Set to name of variable to receive the number of ; scales performed by the decomposition. ; CHECK -- Check number of scales to be performed and return. ; OUTPUTS: ; DECOMPOSITION -- A 2-D array with scale running along the 2nd axis ; (large scales -> small scales). The first plane ; of the array is the smoothed image. To recover ; the image at any scale, just total the array ; along the 2nd dimension up to the scale desired. ; ; ; RESTRICTIONS: ; Uses FFT convolutions which edge-wrap instead of mirroring the ; edges as suggested by Stark & Mutaugh. Wait for it. ; ; MODIFICATION HISTORY: ; ; Mon Oct 6 11:49:50 2003, Erik Rosolowsky ; Written. ; 24 April 2005, JohnJohn: Rewrote as a function, converted to work ; with 1D spectra rather than 2D images. ; ;- on_error, 2 if n_elements(filter) eq 0 then filter = 1./[16, 4, 8/3., 4, 16] nel = n_elements(spec) if 1-keyword_set(n_scales) then \$ n_scales = floor((alog(nel)/n_elements(filter))/alog(2)) decomp = fltarr(nel, n_scales+1) sp = spec for k = 0, n_scales-1 do begin ; Smooth image with a convolution by a filter smooth = convol(sp, filter, /edge_wrap) ; num_conv, sp, filter, smooth decomp[0, n_scales-k] = sp-smooth sp = smooth ; Generate new filter newfilter = fltarr(2*n_elements(filter)-1) newfilter[2*findgen(n_elements(filter))] = filter ; note filter is padded with zeros between the images filter = newfilter endfor ; Stick last smoothed image into end of array decomp[*, 0] = smooth return, decomp end