Function iris_isp2solar,str,roll_angle=roll_angle,noisy=noisy,db=db,pzt2xy=p2xy ;+ ; Name : iris_isp2solar ; ; Purpose : return header structure with IRIS WCS keywords added ; ; Use : ; IDL> db = iris_mk_pointdb() ; IDL> xy = iris_isp2solar(index,db=db,roll_angle=roll_angle) ; ; Input: ; index = FITS header index structure. From level 0, 1, or 1.5 FITS header ; ; Required(!) input keywords: ; roll_angle = float. S/C roll angle. Required, even if 0. ; db = struct. Database of wcs constants. Required keyword. ; ; Optional input keywords: ; /pzt2xy : pzt pointing offset in pixels, primarily for Level 2 processing. ; y coord. along slit, xy does not include wobble and rot.tracking ; requires index to be a level 1.5 header structure ; /noisy : prints lots of diagnostics ; ; Output: ; xy = [solar_x,solar_y]. Location of slit center on sun in arcsec. ; in tang.coord: x: E-W, y:S-N ; if /pzt2xy set: offset in pixels relative to nominal pointing ; with y axis along slit, excl. wobble & rot.tracking ; ; 18-April-2013 TDT version 0.95 ; 20-May-2013 TDT version 0.96 fixes signs for wedges & PZTs (thanks to J-P), ; new regression for PZTs ; 20-May-2013 TDT version 0.97 added shift due to focus mechanism ; 31-May-2013 JPW version 0.98 added db keyword, e.g.,: ; crval=iris_isp2solar(str,db=iris_mk_pointdb()) ; Also added missing TOTALs in the focus correction, and ; moved open-loop GT bias from after to before S/C roll ; 23-July-2013 JPW made db a required keyword to avoid the use of outdated ; pointing parameters. ; 21-Sept-2013 JPW bug fixes in 3 rotations (overwriting old x) ; 05-Oct-2013 JPW added /pzt2xy and improved documentation ; 05-Dec-2013 JPW wedge motor encoder correction ; temperature correction of wm_yoff (temporary fix) ; 06-Dec-2013 JPW allow use of FITS keywords it01pmrf & it06telm if there ; 21-Mar-2013 JPW changed format of wm_yoff temp corr, includes it18sppz ;- ; str = structure made from Level 0 FITS header ; to make from raw header, str = fitshead2struct(header) ;if (n_elements(roll_angle) eq -1) then roll_angle = $ ; get_acs_roll_angle(str.t_obs) common iris_temp_pointdb, tempdb ; tempory implementation of wm_yoff fix version = 'Version 0.99' IF N_ELEMENTS(DB) EQ 0 THEN BEGIN print,'Error in iris_isp2solar.pro: pointing database not supplied return,-1 ENDIF ; IRIS pointing database parameters ; The purpose is to determine the solar pointing x,y, in arcseconds EW ; and NS from disk center, of the center of the slit (half way between ; the 2 fiducials) in the NUV continuum SJI CCD image ; The roll_angle passed to this routine or determined from ACS HK is ; assumed to be the ACS Z-axis angle wrt solar North, in degrees, ; with positive sense meaning positive rotation about S/C +Z axis ; (which is the line-of-sight to sun); ie, the top of the slit in ; the SJI image (large CCD Y coordinate) rotates into the ; NW quadrant of the Sun ; This roll angle = SAT_ROT in Level 1 keywords or ; = A_EULERBR_Z in S/C telemetry (APID 601 = 0x259, ID 122, 64 ; bits, starting byte/bit 116/00) ; ACS raw telemetry is in radians but engineering conversion using ; the HK database routines will convert to degrees ; assume it is degrees for SAT_ROT (& SAT_RTER) ; roll angle bias of slit when ACS roll angle = 0, with same ; sign convention as roll, ie positive means slit is in NW and SE ; quadrants when ACS roll angle = 0; ; in degrees slit_roll_bias = db.slit_roll_bias ; roll angle bias between SJI CCD x & y axes and SG slit, with same ; sign convention as roll, ie positive means CCD y is in NW and SE ; quadrants when slit is exactly NS; ; in degrees ; Preliminary value from ground test data is possible but not here yet sji_ccd_roll_bias = db.sji_ccd_roll_bias ; plate scale of NUV SJI continuum band images, which were used to ; calibrate PZT and WM offsets, units are arcseconds/pixel pzt_off_pixel = db.pzt_off_pixel ; home position of PZT offsets used in calibration of wedge motor ; offsets wm_x(y)off listed below, by 4-limb or other coalignment pzt_off_home = [db.pzt_off_homea,db.pzt_off_homeb,db.pzt_off_homec] ; rows of matrix converting PZT offset ABC to NUV SJI CCD x & y ; These give motion of a spot on a level 0 CCD image caused by PZT offsets ; Derived by regression from 7 Nov 2012 CPT results on raw CCD ; images and 8 Feb 2013 CPT results on transposed CCD images ; Nov 7 images have pixel 1,1 = SE corner, pixel 1096,1036 = NW ; corner, according to Russ Lindgren's diagram ; ; Nov 7 CCD X (1st) = Solar N/S with larger coordinates N = Solar Y ; Nov 7 CCD Y (2nd) = Solar E/W with larger coordinates W = Solar X ; Nov 7 CCD Y is slit scan direction in these raw images ; These Nov 7 images have NOT been rotated in IDL ; ; However, after November, 2012, level 0 CCD images became transposes ; of the earlier level 0 images, as shown by the hot pixels. ; Again, the analysis was done on images NOT rotated in IDL ; Feb 8 CCD X = Nov 7 CCD Y and Feb 8 CCD Y = Nov 7 CCD X ; Feb 8 CCD X (1st) = Solar E/W with larger coordinates W = Solar X ; Feb 8 CCD Y (2nd) = Solar N/S with larger coordinates N = Solar Y ; ; 1-sigma uncertainties from regression ~ 0.0003 ; see file ~tarbell/iris/GT/pzt_offs_regress.txt ; Note that solar coordinate changes of the center of the slit ; (which is fixed on the CCD) are -1 times these pzt_off_sjix = [db.pzt_off_sjixa,db.pzt_off_sjixb,db.pzt_off_sjixc] pzt_off_sjiy = [db.pzt_off_sjiya,db.pzt_off_sjiyb,db.pzt_off_sjiyc] ; For information, NUV SJI continuum pixels per DN of PZT offset for A, B, C ; pzt_off_scale = sqrt(pzt_off_sjix^2.+pzt_off_sjiy^2.) ; sign conversion from 2013 level 0 CCD X direction to solar x (EW, W = +) pzt_x_sign = db.pzt_x_sign ; sign conversion from 2013 level 0 CCD Y direction to solar y (NS, N = +) pzt_y_sign = db.pzt_y_sign ; radius of boresight deflection for WM1, from Akin/Tarbell ground ; cal, in arcsec wm1_rmax = db.wm1_rmax ; radius of boresight deflection for WM2, from Akin/Tarbell ground ; cal, in arcsec wm2_rmax = db.wm2_rmax ; step position of WM2 that cancels the image shift of M1 when WM1=0, ; from Akin/Tarbell ground calibration wm_offset2 = db.wm_offset2 ; step position of M1 that produces offset only in ground ; +V direction = ~IRIS +/-X direction = ~solar NS ; from Akin/Tarbell ground calibration wm_m1v = db.wm_m1v ; sign conversion from ground +H direction to solar x (EW, W = +) ; this has changed since the April 18 version (ver 0.95) ; J-P suggested and I verified that I determined the signs wrong, ; before the GT-ACS polarity test, relying on the incorrect statement ; in the ICD that GT axes were rotated 180 degrees from S/C axes wm_x_sign = db.wm_x_sign ; sign conversion from ground +V direction to solar y (NS, N = +) ; also changed for same reason as wm_x_sign wm_y_sign = db.wm_y_sign ; solar x pointing (at center of slit in SJI continuum image) when GT ; is at center of wedge polar coordinate system and roll=0, in ; arcsec EW (W is +), from 4-limb coalignment or equivalent wm_xoff = db.wm_xoff ; solar y pointing when GT is at center of wedge polar coordinate ; system and roll=0 , in arcsec NS (N is +), from 4-limb coalignment ; or equivalent wm_yoff = db.wm_yoff ; roll angle delta between Akin's ground horizontal & vertical ; coordinate system and SG slit, with same sign convention as roll, ; ie positive means WM vertical axis is in NW and SE quadrants; ; in degrees ; Preliminary value from ground test data is possible but not here yet wm_roll_bias = db.wm_roll_bias ; focus position at which wedge calibration to determine wm_x/yoff was ; done (ie, 4-limb coalignment), in focus mech steps wm_focus_cal = db.wm_focus_cal ; this is a guess, TBC ; Regression coefficients for NUV SJI X & Y motion of laser spot on ; image in pixels vs focus motor position in steps, ; linear & quadratic in (FMPOS+100) where FMPOS is in steps ; Note that solar coordinate changes of the center of the slit ; (which is fixed on the CCD) are -1 times these ; Preliminary values from ground test data of 7 Nov 2012, see ; ~tarbell/iris/GT/file focus_mech_regression.txt ; Note 2013 Level 0 data are transposes of Nov 2012 data--this has ; been taken care of in the coefficients below. ; This offset of +100. in the independent variable need not ; be equal to -wm_focus_cal--I just happened to use -100 as the ; "center" when doing the regression on Nov 7 data. focus_xregr = [db.focus_xregr1,db.focus_xregr2] focus_yregr = [db.focus_yregr1,db.focus_yregr2] foc_rgcl = db.focus_regcal ; Wedge motor encoder correction wm1_enco = db.wm1_enco wm2_enco = db.wm2_enco ; Coefficients for temperature correction of wm_yoff wm_ytfit0 = db.wm_ytfit0 wm_ytfit1 = db.wm_ytfit1 ; End of setting values from database ; Derive the X & Y shifts (px) due to focus, for the position wm_focus_cal focus_temp = [float(wm_focus_cal)-foc_rgcl,(float(wm_focus_cal)-foc_rgcl)^2] wm_focus_cal_x = total(focus_xregr*focus_temp) wm_focus_cal_y = total(focus_yregr*focus_temp) ; ******************************************* ; get the temperatures for the wm_yoff_t temperature correction ; - get them from the FITS header if available ; - otherwise use a ssw database read into the common block if total(abs(wm_ytfit1)) ne 0. then begin ; get the temps from the header, if valid if tag_exist(str,'it01pmrf') then if ~finite(str.it01pmrf,/nan) $ then it01pmrf=str.it01pmrf if tag_exist(str,'it06telm') then if ~finite(str.it06telm,/nan) $ then it06telm=str.it06telm if tag_exist(str,'it18sppz') then if ~finite(str.it18sppz,/nan) $ then it18sppz=str.it18sppz ; read the ssw database with the daily temperature averages if necessary if (n_elements(it01pmrf) eq 0 or n_elements(it06telm) eq 0 or $ n_elements(it18sppz) eq 0) and n_elements(tempdb) eq 0 then $ restgenx,file=concat_dir('$SSW_IRIS_DATA','iris_temp_pointdb.geny'),tempdb ; get temps from the database if n_elements(tempdb) gt 0 then begin if n_elements(it01pmrf) eq 0 then $ it01pmrf = interpol(tempdb.it01pmrf,tempdb.anytim,anytim(str.t_obs)) if n_elements(it06telm) eq 0 then $ it06telm = interpol(tempdb.it06telm,tempdb.anytim,anytim(str.t_obs)) if n_elements(it18sppz) eq 0 then $ it18sppz = interpol(tempdb.it18sppz,tempdb.anytim,anytim(str.t_obs)) ; write the temps into the header if they were NaN if tag_exist(str,'it01pmrf') then if finite(str.it01pmrf,/nan) $ then str.it01pmrf=it01pmrf if tag_exist(str,'it06telm') then if finite(str.it06telm,/nan) $ then str.it06telm=it06telm if tag_exist(str,'it18sppz') then if finite(str.it18sppz,/nan) $ then str.it18sppz=it18sppz endif endif ; get relevant pzt_offsets from ISP (header) ; isqpata/b/c includes both rotation tracking & orbital wobble ; orbital wobble does not change solar pointing pzts = [str.isqpzta, str.isqpztb,str.isqpztc] - $ [str.iwbpzta, str.iwbpztb,str.iwbpztc]- pzt_off_home ; if used for assembling level2 data, subtract rotation tracking offset: if keyword_set(p2xy) then pzts=pzts-[str.irtpzta,str.irtpztb,str.irtpztc] ; minus sign because solar coord change = -motion of image ccdx = -total(pzt_off_sjix*pzts)*pzt_off_pixel*pzt_x_sign ccdy = -total(pzt_off_sjiy*pzts)*pzt_off_pixel*pzt_y_sign if keyword_set(noisy) then print,'SQ PZTs - home ',pzts if keyword_set(noisy) then print,'SQ PZTs ccdx & ccdy: ',ccdx,ccdy ; add shift of image by focus mechanism using current position focus ; and subtract image shift by focus at the calibration pos focus_temp = [float(str.ifmpos)-foc_rgcl,(float(str.ifmpos)-foc_rgcl)^2] focusx = total(focus_temp*focus_xregr) - wm_focus_cal_x focusy = total(focus_temp*focus_yregr) - wm_focus_cal_y ; minus sign because solar coord change = -motion of image ccdx2 = ccdx - focusx*pzt_off_pixel ccdy2 = ccdy - focusy*pzt_off_pixel if keyword_set(noisy) then print,'PZTs & focus ccdx & ccdy: ',ccdx2,ccdy2 ; rotate these to solar X, Y at roll = 0 if keyword_set(p2xy) then begin ; relative to slit, so no slit_roll_bias sth = sin((sji_ccd_roll_bias)/!radeg) cth = cos((sji_ccd_roll_bias)/!radeg) endif else begin ; relative to S/C coord: w/slit_roll_bias sth = sin((slit_roll_bias+sji_ccd_roll_bias)/!radeg) cth = cos((slit_roll_bias+sji_ccd_roll_bias)/!radeg) endelse ccdx = ccdx2*cth + ccdy2*sth ccdy = -ccdx2*sth + ccdy2*cth if keyword_set(noisy) then print,'ccdx & ccdy at roll=0: ',ccdx,ccdy ; add wedge motor image shifts from formulae in iris_wedge2solar_tl ; unles p2xy set: if keyword_set(p2xy) then begin x2 = ccdx y2 = ccdy endif else begin fac = 2.*!pi/240. ;th1 = fac*(float(str.iwm1cpos)-wm_m1v) ;th2 = fac*(float(str.iwm2cpos)-wm_m1v+120.-wm_offset2) ; adding wedge motor encoder correction: mod1 = fix(str.iwm1cpos+0.5) mod 2 th1 = float(str.iwm1cpos) + mod1*wm1_enco th1 = fac*(th1-wm_m1v) mod2 = fix(str.iwm2cpos+0.5) mod 2 th2 = float(str.iwm2cpos) + mod2*wm2_enco th2 = fac*(th2-wm_m1v+120.-wm_offset2) xwm2 = wm_x_sign*(wm1_rmax*sin(th1) + wm2_rmax*sin(th2)) ywm2 = wm_y_sign*(wm1_rmax*cos(th1) + wm2_rmax*cos(th2)) if keyword_set(noisy) then print,'WM boresight shift x & y: ',xwm2,ywm2 ; rotate by roll bias wrt solar North at roll=0 degrees sth = sin((slit_roll_bias+wm_roll_bias)/!radeg) cth = cos((slit_roll_bias+wm_roll_bias)/!radeg) xwm = xwm2*cth + ywm2*sth ywm = -xwm2*sth + ywm2*cth if keyword_set(noisy) then print,'WM rotated shift x & y: ',xwm,ywm ; For the next step, calculate the temperature correction for wm_yoff: if n_elements(it01pmrf) eq 1 and n_elements(it06telm) eq 1 and $ n_elements(it18sppz) eq 1 and n_elements(wm_ytfit1) eq 3 then begin it_temps = [it01pmrf,it06telm,it18sppz] wm_yoff_t = wm_ytfit0 + total(it_temps*wm_ytfit1) endif else wm_yoff_t = 0. ; add the center of the WM coordinate system: note this is in solar ; X & Y, not in the wedge motor or GT diode or CCD X & Y xwm = xwm + wm_xoff ywm = ywm + wm_yoff + wm_yoff_t if keyword_set(noisy) then print,'WM x & y at roll=0: ',xwm,ywm ; co-add both contributions x2 = xwm + ccdx y2 = ywm + ccdy endelse ; add GT bias offsets in ACS packets, only if ISS loop is not closed ; if loop is closed, then PZTs will zero this (assuming they have ; enough range) ; this is very rare and may never be used during science observing, ; so ignore the small roll bias between the GT diodes coordinate ; system and the others if (str.iissloop ne 'CLOSED') then begin ; minus sign because solar coord change = -motion of image x2 = x2 - str.igtpoffx*0.01 y2 = y2 - str.igtpoffy*0.01 if keyword_set(noisy) then print,'solar X, Y after GT bias offsets: ',x2,y2 end ; Now rotate everything to non-zero roll angle ; unless p2xy set, in that case, scale from arcsec to pixels instead (no roll) if keyword_set(p2xy) then begin x = x2 / str.cdelt2 y = y2 / str.cdelt2 endif else begin sroll = sin(roll_angle/!radeg) croll = cos(roll_angle/!radeg) x = x2*croll + y2*sroll y = -x2*sroll + y2*croll endelse if keyword_set(noisy) then print,'solar X, Y, roll_angle: ',x,y,roll_angle return,[x,y] end