pro slitshft_sbsp,singleframe,cc1,cc2,sp1,sp2,sl1,sl2,ssum, $ ssl1,ssl2,nsearch,shftsl,delsh ; routine to compute the spatial shifts of the image along the slit ; direction as a function of spectral pixel number. When applied, ; these shifts rectify the images so that spectral variation is strictly ; along the spectral pixel direction. This process is known as removal ; of the spectral skew. The relative shift of the second CCDSIDE is found ; and folded into the output array. These are computed in advance of ; flat-fielding from the raw flat-field images. ; For HINODE SP, the continuum contrast of granulation plus the slit ; width variations of the raw data are sufficient for one frame of data ; to derive these shifts. No dark correction is necessary, since the ; dark levels are so uniform. ; INPUTS: singleframe - path to single flat field QS observation ; cc1,cc2 - pixel range for averaging spectral continuum ; sp1,sp2 - active spectral pixel range, each CCSDIDE ; sl1,sl2 - active pixel range along slit in CCD units ; nsearch - correlation pixel search range ; OUTPUTS: ; ssum - pixel summing along slit ; ssl1,ssl2 - active pixel range along slit, data units ; shftsl - shift to be applied for each spectral pixel ; delsh - shift between two CCDSIDES, to be applied to ; CCDSIDE2 ; read in the single frame data dat = float(readfits(singleframe,hdr)) ; EXTRACT PARAMETERS FROM HEADER ; Extract the number of pixels binned along the slit (along serial direction) ssum = sxpar(hdr,'CAMSSUM') ; Range of pixels along the slit ; Note that the pixel indices along the slit direction are reversed from ; those during pre-launch testing. The images are now reversed in ; the vertical direction during reformatting. Also, pixel indices in ; the header like SPCCDIX0, SPCCDIX1 refer to the pre-launch CCD coordinates, ; not the coordinates of LEVEL0 data. spccdix0 = sxpar(hdr,'SPCCDIX0') spccdix1 = sxpar(hdr,'SPCCDIX1') ; define these pixel positions in terms of the LEVEL0 image coordinates ; reverse order so that low limit of active range is in spccd0 spccd0 = 1023-spccdix1 spccd1 = 1023-spccdix0 ; bit shifting for SP data bitshft = sxpar(hdr,'SPBSHFT') ; Derive active range along slit for shifting in recorded pixels, not CCD pixels ssl1 = max([sl1(1),spccd0]) - spccd0 ; always truncate the limit in the case of summing along slit ssl2 = (min([sl2(1),spccd1]) - spccd0)/ssum if ssum eq 2 then begin if (ssl1 mod 2) eq 1 then begin ssl1 = (ssl1+1)/ssum endif else begin ssl1 = ssl1/ssum endelse endif ; correct for digital wrap-around of Stokes I, if wrap-around is present fltav = dat(*,*,*,0) whr = where(fltav lt 0., countwrap) if countwrap ne 0 then fltav(whr) = fltav(whr) + 65536. ; account for bit shifting of the data by multiplying the appropriate ; Stokes images by 2 if bitshft ne 0 then fltav = 2.*fltav ; working arrays extended next power of two larger in length nyy = sizeof_sbsp(fltav,2) for i = 0,20 do begin np2 = 2^i if np2 gt nyy then begin ntrp = np2 goto, kilroy endif endfor kilroy: sltref = fltarr(ntrp) sltvar = fltarr(ntrp) ; slit variation extended to next power of 2 nxx = sizeof_sbsp(fltav,1) ; full spectral dimension of data xx = findgen(nxx) ; spectral index shftsl = fltarr(nxx,2) ; final fitted shifts for each spectral pixel sdel = fltarr(nxx) & sdel(*) = 0. ; computed shifts for each spectral slop = fltarr(2) ; output array for slope of spectrum for two sides cctvar = fltarr(ntrp) ; continuum intensity variation along slit im1 = fltav(*,*,0) im2 = fltav(*,*,1) ; loop over two CCDSIDES window,xs=700,ys=500,/free for ks = 0,1 do begin if ks eq 0 then img = im1 else img = im2 ; Step 1: determine the variations of the spectrum skew in the ; spectral direction. Take a continuum fluctuation along the slit ; dimension and cross-correlate as a function of spectral pixel. ; Compute the relative shifts between the two sides also so that ; the resulting shifts applied to the data will put the two sides ; at the same horizontal pixel ny = ssl2 - ssl1 + 1 & nw = sp2(ks) - sp1(ks) + 1 ; get reference variation along slit for j = ssl1,ssl2 do cctvar(j) = mean_sbsp(img(cc1:cc2,j)) ; extend for smooth transition to average, remove mean sltref = extend_sbsp(cctvar,ssl1+1,ssl2-1) sltref = sltref - mean_sbsp(sltref) ; Fourier transform the reference variation along slit fftera = fft(sltref,-1) if ks eq 0 then fftera1 = fftera ; do cross-correlation with each slit variation temp = fltarr(ntrp) for i = 0,nw-1 do begin ; extend data, remove mean sltvar(ssl1:ssl2) = img(i,ssl1:ssl2) sltvar = extend_sbsp(sltvar,ssl1+1,ssl2-1) sltvar = sltvar - mean_sbsp(sltvar) ffterb = fft(sltvar,-1) ;; Find the maximum of cross-correlation. cormax_sbsp,fftera,ffterb,nsearch,sh sdel(i+sp1(ks)) = sh endfor ; collapse spectral variation of image to an slit-average spectrum avspec = fltarr(nw) for i = 0,nw-1 do avspec(i) = mean_sbsp(img(i,*)) whr = where(avspec gt 0.99*median(avspec)) xxp = findgen(nw) xxp = xxp(whr) + sp1(ks) yyp = sdel(sp1(ks):sp2(ks)) yyp = yyp(whr) ; fit the variation with a linear function clin = linfit(xxp,yyp) slop(ks) = clin(1) fitlin = clin(0) + clin(1)*xx if ks eq 0 then begin plot,xx(sp1(ks):sp2(ks)),sdel(sp1(ks):sp2(ks)), $ xrange=[0,nxx],yrange=[-.5,.5], $ title = 'Spectral Skewness and linear fits to two CCDSIDES' oplot,xx(sp1(ks):sp2(ks)),fitlin(sp1(ks):sp2(ks)),linestyle=2 endif else begin oplot,xx(sp1(ks):sp2(ks)),sdel(sp1(ks):sp2(ks)) oplot,xx(sp1(ks):sp2(ks)),fitlin(sp1(ks):sp2(ks)),linestyle=2 endelse ; Step 2: now compute the relative shifts for each of the original ; wavelength pixels, relative to the reference for the first CCDSIDE ; positive means shift the slit variation towards larger number shftsl(*,ks) = -fitlin if ks eq 1 then begin ; cross-correlate the slit variations for the two sides cormax_sbsp,fftera1,fftera,nsearch,sh ; adjust for shift between two sides shftsl(*,ks) = -fitlin - sh endif ; end loop over CCDSIDE endfor ; print offsets between ccdsides print,' offsets between ccdsides: ',sh delsh = -sh print,' pixel offset to be applied to CCDSIDE2 to align with CCDSIDE1:',delsh return end