;+ ; NAME: ; iris_calc_wave() ; ; PURPOSE: ; A sub-routine of iris_spec_cal. ; ; Accurately determines the line centers by fitting Gaussians ; and fits a parabolic function to the wavelength solution in an ; image which has already be geometrically corrected, such that ; the spectral line centers all fall at the same pixel from the ; top to the bottom of the CCD. ; ; CALLING SEQUENCE: ; result = iris_calc_wave(img, wavelength, pixel, width, karray, wave='') ; ; INPUTS: ; img - The spectrum from which the wavelength will be fitted. ; wavelength - An array containing all the line wavelengths for ; the lines to be fitted. ; pixel - An array containing the first guess pixel positions of ; the lines. ; width - An array containing the width in pixels to fit ; around each line ; wave - A string containing the IRIS wavelenth ID for the ; spectrum, either 'nuv', 'fuvs', or 'fuvl'. ; ; OPTIONAL INPUTS: ; None ; ; KEYWORD PARAMETERS: ; None ; ; OUTPUTS: ; Returns a structure containing the tags: ; .fit - the array of fitted parabolic parameters which define ; the wavelength solution. ; e.g. wavelength = a[0] + a[1]*pixel + a[2]*pixel^2 ; .dl - The linear dispersion [Ang/pix] of the data ; after the wavelength and geometry correction ; (this is hard-coded) ; .dmapx - The x-component of the distortion map used ; to apply the wavelength+geometry solution ; .dmapy - The y-component of the distortion map used ; to apply the wavelength+geometry correction. ; .wave_coord - The array of wavelength coordinates for ; every pixel in the original image (in Angstroms) ; .spat_coord - The array of spatial coordinates for every ; pixel in the original image (in pixels) ; ; e.g. result = a where wavelength = a[0] + a[1]*pixel + a[2]*pixel^2 ; ; OPTIONAL OUTPUTS: ; None ; ; COMMON BLOCKS: ; None ; ; PROCEDURES USED: ; ; COMMENTS: ; ; EXAMPLES: ; IDL> result = iris_calc_wave(img, wavelength, pixel, wave='nuv') ; ; MODIFICATION HISTORY: ; Started 2013-June-3 by Sarah A. Jaeggli, Montana State University ; ; 2013-June-12 SAJ Removed dimension keyword in mean and ; stddev calls for IDL 7 compatibility ; ; 2013-Sept-16 SAJ changed philosophy on using ; polynomial solution, wavelength now produces ; distortion maps, as a result the spatial and wave ; coordinate is available for every pixel in Level0/1 ; ; 2014-Jun-13 SAJ wavelengths are now fitted from a ; spatially averaged spectrum, this greatly increases ; fitting accuracy for weak neutral lines in the FUV ; ; 2014-Jun-23 SAJ added width parameter for each line so ; fitted regions can be specified precisely ; ; 2014-July-28 SAJ standardized FUV1, FUV2 labeling to ; match image header standards ;- FUNCTION iris_calc_wave, img, wavelength, pixel, width, karray, wave=wave ;set up wavelength-specific parameters case wave of 'NUV' : begin nterms = 6 ;number of terms in the line fit tol = 5 ;number of pixels tolerance in result sign =-1. ;sign of line amplitude dl = 0.02546D ;desired dispersion of unwarped spectrum end 'FUV2' : begin nterms = 5 ;number of terms in the line fit tol = 5 ;number of pixels tolerance in result sign = 1. ;sign of line amplitude dl = 0.01272D end 'FUV1' : begin nterms = 5 ;number of terms in the line fit tol = 5 ;number of pixels tolerance in result sign = 1. ;sign of line amplitude dl = 0.01298D end endcase ;get image dimensions nx=(size(img))[1] ny=(size(img))[2] xarray=dindgen(nx) yarray=dindgen(ny) nlines=n_elements(pixel) ;open a new window for plotting results window, 0, xsize=500, ysize=300 ;set up fit results buffer xpos=fltarr(nlines) err_xpos=fltarr(nlines) ;average spectrum over entire spatial dimension linespec=mean(img, dimension=2, /nan) for i=0,nlines-1 do begin a=0 xx=xarray[pixel[i]-width[i]:pixel[i]+width[i]] ii=linespec[pixel[i]-width[i]:pixel[i]+width[i]] ; ifit=gaussfit(xx, ii, a, nterms=nterms) ifit=mpfitpeak(xx, ii, a, nterms=nterms, /gaussian, /nan) plot, xx, ii, xstyle=1, ystyle=1, $ title='Line Fit', xtitle='Pixel', ytitle='Intensity' oplot, xx, ifit, color=250 wait, 1 ;throw out absurd fits and populate fit array if sign*a[0] lt 0. or $ a[1] gt pixel[i]+tol or $ a[1] lt pixel[i]-tol then a[*]=-!values.f_nan xpos[i]=a[1] endfor good=where(finite(xpos) eq 1) xpos=xpos[good] wavelength=wavelength[good] fit=poly_fit(xpos, wavelength, 1, yfit=yfit1) fit=poly_fit(xpos, wavelength, 2, yfit=yfit2) plot, xpos, wavelength-yfit1, xstyle=2, ystyle=2, psym=4, $ title='Wavelength Fit', xtitle='Pixel', ytitle='Wavelength [Ang]' oplot, xpos, yfit2-yfit1, color=250 ;retrieve polywarp arrays from Karray variable kx=karray[*,*,0] & ky=karray[*,*,1] kxi=karray[*,*,2] & kyi=karray[*,*,3] nkx=(size(kx))[1] & nky=(size(kx))[2] ;set the dispersion to a constant value (defined in case statements) c1=fit[1]/dl c2=fit[2]/dl ;unwarped (original) coordinates x0=rebin(reform(xarray, nx, 1), nx, ny) y0=rebin(reform(yarray, 1, ny), nx, ny) xo=replicate(0.D,nx,ny) yo=replicate(0.D,nx,ny) xi=replicate(0.D,nx,ny) yi=replicate(0.D,nx,ny) for i=0,nkx-1 do begin for j=0,nky-1 do begin ;forward xo+=kx[i,j] * x0^j * y0^i yo+=ky[i,j] * x0^j * y0^i ;and backward transformations xi+=kxi[i,j] * x0^j * y0^i yi+=kyi[i,j] * x0^j * y0^i endfor endfor ;distortion array for geo-transformed image, linearized in wavelength delx=x0-(c1*x0 + c2*x0^2) dely=0.D*y0 delxp=interpolate(delx, xi, yi, cubic=-0.5) delyp=interpolate(dely, xi, yi, cubic=-0.5) dmapx=xo + delxp dmapy=yo + delyp ;also generate an array of spatial and wavelength coordinates wave_array=fit[0] + fit[1]*x0 + fit[2]*x0^2 spat_array=y0 ;transform these into the curved (original) coordinate frame wave_array_prime=interpolate(wave_array, xi, yi, cubic=-0.5) spat_array_prime=interpolate(spat_array, xi, yi, cubic=-0.5) RETURN, {fit:fit, dl:dl, dmapx:dmapx, dmapy:dmapy, $ wave_coord:wave_array_prime, spat_coord:spat_array_prime} END