pro euvi_point,ii,error=errs,gtfile=gfil,no_grh=no_grh,forceroll=forcr, $ quiet=quiet,version=make_ver ;+ ; $Id: euvi_point.pro,v 1.23 2024/06/28 15:55:50 nathan Exp $ ; ; Project : STEREO SECCHI ; ; Name : euvi_point ; ; Purpose : update CDELTi & PCj_i, and calculate CRPIXi from FPS and GT info ; ; Use : IDL> euvi_point,i,err=err ; ; Input/Output : ; i = FITS index structure. Attitude related keywords are updated, ; including .crpixi, .crvali, cdelti, .pci_j, .att_file (and others) ; ; Keywords input: ; gtfile : two options: ; 1. path/name of a single Level 0 secchi telemetry packet file ; (name must follow proper naming convention: *yyyy_doy_n_vv.???) ; 2. path/name of GT database file (created by scc_gtdbase_update.pro) ; (name must end in 'geny') ; If none supplied, looks in $SCC_DATA/gt/ then in current dir for ; 'secchi_gtdbase.geny' ; /no_grh : no ground receipt header (ignored if file is GT database file) ; /forceroll : Get roll info from SPICE, even if i.att_file not "NotFound" ; version : string. Use to apply an old version of the pointing correction. ; Any version prior to '1.16' will not correct for pointing drift. ; Default: latest version. ; ; Keywords output: ; error = byte. 0 if no error, otherwise nonzero error code: ; 7 : invalid observatory. Must be STEREO_A or STEREO_B. ; 3 : error reading GT calibration file. ; 1 : outdated GT calibration data ; Note: this keyword does not indicate whether a successful roll ; update was done. ; Instead, check the ATT_FILE keyword for a valid file name ; ; Common : ; ; Restrictions: ; 1. Assumes that the image has been rectified, i.e., ecliptic N up! ; ; Side effects: ; ; Category : Pipeline ; ; Prev. Hist.: ; Written : Jean-Pierre Wuelser, LMSAL, 9-Aug-2006 ; JPW, 6-Feb-2007 Major changes ; JPW, 16-Feb-2007 coefficient update, update crota with SPICE ; JPW, 10-Apr-2007 call get_stereo_att_file instead of ; get_stereo_spice kernel -> major speed gain ; JPW, 8-May-2007 coefficient updates: ; - GT-EUVI roll from offpoint analysis ; - EUVI-S/C roll from lunar transit ; checking history to see if roll recalc. needed ; added /quiet, fixed putting err into att_file ; JPW, 21-Jun-2007 changed logic on how to do the roll correction ; - checks whether spice/icy dlm present ; - adds EUVI roll offset even if no attit. file ; now allowing vector of header structure ; ; $Log: euvi_point.pro,v $ ; Revision 1.23 2024/06/28 15:55:50 nathan ; JPW - CRPIXi drift correction update post 2015 conjunction ; ; Revision 1.22 2023/08/15 22:15:50 secchia ; roll08/12 correction ; ; Revision 1.21 2015/07/16 18:49:38 nathan ; JPW - update for post-conjunction configuration ; ; Revision 1.20 2013/07/29 18:21:34 nathan ; JPW - moved the definitions of crpix0/cdelt0 after SECCHI header check ; ; Revision 1.19 2013/06/12 20:03:26 nathan ; add loud msg if recalculating roll ; ; Revision 1.18 2013/06/12 19:26:30 nathan ; add /loud output; use i.summed for determining binning ; ; Revision 1.17 2013/06/07 20:43:34 nathan ; print delta changes and put in history ; ; Revision 1.16 2013/06/07 18:26:46 nathan ; Completely new version of euvi_point.pro by J-P Wuelser, May 9, 2013. ; Main update is implementing pointing drift correction. Also has logic ; for outdates databases, SPICE not available, whether to do new roll ; correction. ; ; Revision 1.15 2012/04/05 21:29:06 nathan ; append VERSION keyword in input ; ; Revision 1.14 2011/02/02 12:54:32 mcnutt ; only applies erol correction to crota id att_file available to avoid the correction from being done twice ; ; Revision 1.13 2009/01/15 19:05:58 nathan ; only apply ins_r0 roll offset if recomputing SC roll angle ; ; Revision 1.12 2008/12/22 20:43:33 wuelser ; test icy before convert_stereo_lonlat ; ; Revision 1.11 2008-09-26 21:21:50 nathan ; moved computation of xycen to after PC matrix ; ; Revision 1.10 2007/08/30 18:58:24 secchia ; nr - make ins_?0 pixel independent ; ; Revision 1.9 2007/08/30 18:48:21 nathan ; added INS_[xyr]0 (Bug 91) ; ; Revision 1.8 2007/08/23 23:52:08 nathan ; obtain RA/Dec coordinates ; ; Revision 1.7 2007/06/22 19:47:38 nathan ; CDELT1A should be negative (Bug 191) ; ; Revision 1.6 2007/06/21 23:44:41 euvi ; new logic for roll update ; ; Revision 1.5 2007/05/08 22:50:38 euvi ; coefficient update ; ; Revision 1.4 2007/05/04 21:29:05 nathan ; left stop in ; ; Revision 1.3 2007/05/04 21:18:28 nathan ; add err to ATT_FILE and update_history ; ; Revision 1.2 2007/04/11 18:49:54 euvi ; use get_stereo_att_file, faster ; ; Revision 1.1 2007/02/19 23:09:03 euvi ; initial version ; ;- info="$Id: euvi_point.pro,v 1.23 2024/06/28 15:55:50 nathan Exp $" IF keyword_set(QUIET) THEN loud=0 ELSE loud=1 if (loud) then print,info icy = test_spice_icy_dlm() ; spice available? errs = bytarr(n_elements(ii) > 1) rollwarningflag = 0 for j=0L,n_elements(ii)-1 do begin i = ii[j] ; valid observatory? if tag_exist(i,'obsrvtry') then obss = strupcase(strtrim(i[0].obsrvtry,2)) $ else obss = '' if tag_exist(i,'detector') then dets = strupcase(strtrim(i[0].detector,2)) $ else dets = '' err = 0B if obss eq 'STEREO_A' and dets eq 'EUVI' then begin ; set various constants. For A: obs = 0 obs_s = 'a' flp1 = 1.0 ; E-W flipped during rectify, 1=true flp2 = 1.0 ; S-N flipped during rectify, 1=true ; initial values of off1,2 on 2007-02-07, prior to fit_limb ; p2a = 1.590 ; EUVI plate scale in arcsec ; off1 = 1027.0 ; GT axis location on EUVI CCD ; off2 = 1110.0 ; relative to p1col=0,p1row=0 ; erol = 0.0 / !radeg ; Roll of EUVI-N measured east from S/C N ; values updated on 2007-02-16: off1,2 based on go_scc_fit_limb1.pro, ; p2a_A/p2a_B & erol_A-erol_B based on go_scc_fit_rot1.pro with updated CROTA p2a = 1.590*0.9986 ; EUVI plate scale in arcsec (A=0.9986*B) ; off1 = 1028.0 ; GT axis location on EUVI CCD ; off2 = 1110.1 ; relative to p1col=0,p1row=0 ; erol = 1.245 / !radeg ; Roll of EUVI-N measured east from S/C N ; grol = 0.006 ; Roll of GT-N measured east from EUVI N ; value updates 2007-05-02: (assumes GT plate scale mod in scc_gt2sunvec) off1 = 1027.5 ; GT axis location on EUVI CCD off2 = 1110.6 ; relative to p1col=0,p1row=0 erol = (1.245-1.125) / !radeg ; Roll of EUVI-N measured east from S/C N grol = -0.005 ; Roll of GT-N measured east from EUVI N ; Post-conjunction changes: if tag_exist(i,'rectrota') then if i.rectrota eq 4 then begin flp1 = (1.0-flp1) flp2 = (1.0-flp2) grol = grol + !pi erol = erol + !pi endif ; define the pointing drift records for A: npdrec = 8 ; total number of pd records pdrec = replicate({obs:'',ver:0,ts:0d0,te:0d0,t0:0d0,torb:0d0, $ c1:fltarr(9),c2:fltarr(9)},npdrec) pdrec.obs = 'a' pdrec.torb = 344.6*86400.0 pdrec.t0 = anytim('2007-01-01T00:00') k = 0 ; first record: no drift correction pdrec[k].ver = 0 pdrec[k].ts = anytim('2000-01-01T00:00') pdrec[k].te = anytim('2099-01-01T00:00') pdrec[k].t0 = anytim('2007-01-01T00:00') pdrec[k].c1 = fltarr(9) ; no correction pdrec[k].c2 = fltarr(9) ; no correction k = 1 pdrec[k].ver = 1016 pdrec[k].ts = anytim('2007-04-01T00:00') pdrec[k].te = anytim('2009-04-01T00:00') pdrec[k].t0 = anytim('2007-01-01T00:00') pdrec[k].c1 = [ -0.19698, 0.40215, -0.02734, 0.02834, 0.02683, $ 0.01781, 0.00528, 0.00799, 0.70364] pdrec[k].c2 = [ 2.01189, -2.56016, 0.08225, 0.05025, -0.01835, $ 0.01585, 0.00151, 0.00655, 0.78864] k = 2 pdrec[k].ver = 1016 pdrec[k].ts = anytim('2009-04-01T00:00') pdrec[k].te = anytim('2011-04-01T00:00') pdrec[k].t0 = anytim('2009-01-01T00:00') pdrec[k].c1 = [ -1.21082, 1.17563, -0.03543, -0.01422, -0.00833, $ 0.01370, -0.00109, -0.01159, 0.33896] pdrec[k].c2 = [ 1.61051, 0.02998, 0.01900, 0.16401, -0.01370, $ 0.04178, 0.00161, -0.01293, -0.99862] k = 3 pdrec[k].ver = 1016 pdrec[k].ts = anytim('2011-04-01T00:00') pdrec[k].te = anytim('2099-01-01T00:00') pdrec[k].t0 = anytim('2011-01-01T00:00') pdrec[k].c1 = [ -1.59089, 0.96487, -0.00403, -0.00339, -0.02923, $ -0.04176, 0.01798, 0.01447, 0.26979] pdrec[k].c2 = [ 2.28952, -0.68690, -0.10883, 0.22979, -0.00715, $ 0.03479, 0.00254, 0.01736, 3.86424] k = 4 pdrec[k].ver = 1021 pdrec[k].ts = anytim('2015-05-19T00:00') ; post conjunction pdrec[k].te = anytim('2023-08-12T00:00') ;s/c flip pdrec[k].t0 = anytim('2011-01-01T00:00') pdrec[k].c1 = -pdrec[3].c1 ; polarity change pdrec[k].c1[8] = pdrec[3].c1[8] ; no polarity change for exponent term pdrec[k].c2 = -pdrec[3].c2 ; polarity change pdrec[k].c2[8] = pdrec[3].c2[8] ; no polarity change for exponent term k = 5 pdrec[k].ver = 1023 pdrec[k].ts = anytim('2015-11-18T00:00') pdrec[k].te = anytim('2021-08-01T00:00') pdrec[k].t0 = anytim('2015-11-18T00:00') pdrec[k].c1 = [ 1.16430, 0.11202, 0.05279, 0.05599, -0.01810, $ 0.01131, 0.00301, -0.00421, -0.40094] pdrec[k].c2 = [ -2.44524, -0.10207, 0.10195, 0.00222, -0.01788, $ 0.00797, -0.01294, 0.02516, -0.42572] k = 6 pdrec[k].ver = 1023 pdrec[k].ts = anytim('2021-08-01T00:00') pdrec[k].te = anytim('2099-01-01T00:00') pdrec[k].t0 = anytim('2021-08-01T00:00') pdrec[k].c1 = [ 2.94178, -0.19203, -0.00031, 0.07529, -0.02484, $ -0.00588, -0.00306, 0.00382, 1.35169] pdrec[k].c2 = [ -6.85677, 3.03264, 0.09955, 0.00847, 0.02142, $ 0.01430, -0.01344, 0.02884, 0.10000] k = 7 pdrec[k].ver = 1023 pdrec[k].ts = anytim('2023-08-12T00:00') ; post inferior conj., s/c flip pdrec[k].te = anytim('2031-10-01T00:00') ; 2nd superior conj.? pdrec[k].t0 = anytim('2021-08-01T00:00') pdrec[k].c1 = -pdrec[6].c1 ; polarity change pdrec[k].c1[8] = pdrec[6].c1[8] ; no polarity change for exponent term pdrec[k].c2 = -pdrec[6].c2 ; polarity change pdrec[k].c2[8] = pdrec[6].c2[8] ; no polarity change for exponent term endif else if obss eq 'STEREO_B' and dets eq 'EUVI' then begin ; constants for B: obs = 1 obs_s = 'b' flp1 = 0.0 ; E-W flipped during rectify, 1=true flp2 = 1.0 ; S-N flipped during rectify, 1=true ; initial values of off1,2 on 2007-02-07, prior to fit_limb ; p2a = 1.590 ; EUVI-B plate scale in arcsec ; off1 = 1035.0 ; GT axis location on EUVI CCD ; off2 = 1096.0 ; relative to p1col=0,p1row=0 ; erol = 0.0 / !radeg ; Roll of EUVI-N measured east from S/C N ; values updated on 2007-02-16: off1,2 based on go_scc_fit_limb1.pro, ; prelaunch p2a & erol=0 for Behind, awaiting lunar transit results p2a = 1.590 ; EUVI-B plate scale in arcsec ; off1 = 1034.6 ; GT axis location on EUVI CCD ; off2 = 1095.4 ; relative to p1col=0,p1row=0 ; erol = 0.000 / !radeg ; Roll of EUVI-N measured east from S/C N ; grol = 0.020 ; Roll of GT-N measured east from EUVI N ; value updates 2007-05-02: (assumes GT plate scale mod in scc_gt2sunvec) off1 = 1034.4 ; GT axis location on EUVI CCD off2 = 1095.7 ; relative to p1col=0,p1row=0 erol = -1.125 / !radeg ; Roll of EUVI-N measured east from S/C N grol = -0.015 ; Roll of GT-N measured east from EUVI N ; Post-conjunction changes: if tag_exist(i,'rectrota') then if i.rectrota eq 1 then begin flp1 = (1.0-flp1) flp2 = (1.0-flp2) grol = grol + !pi erol = erol + !pi endif ; define the pointing drift records for B: npdrec = 5 ; total number of pd records pdrec = replicate({obs:'',ver:0,ts:0d0,te:0d0,t0:0d0,torb:0d0, $ c1:fltarr(9),c2:fltarr(9)},npdrec) pdrec.obs = 'b' pdrec.torb = 389.0*86400.0 k = 0 ; first record: no drift correction pdrec[k].ver = 0 pdrec[k].ts = anytim('2000-01-01T00:00') pdrec[k].te = anytim('2099-01-01T00:00') pdrec[k].t0 = anytim('2007-01-01T00:00') pdrec[k].c1 = fltarr(9) ; no correction pdrec[k].c2 = fltarr(9) ; no correction k = 1 pdrec[k].ver = 1016 pdrec[k].ts = anytim('2007-05-01T00:00') pdrec[k].te = anytim('2009-06-01T00:00') pdrec[k].t0 = anytim('2007-01-31T00:00') pdrec[k].c1 = [ -3.61068, 3.71167, 0.04840, -0.34213, -0.00013, $ -0.02756, -0.00926, -0.00304, 0.10000] pdrec[k].c2 = [ -3.50374, 3.65205, 0.06367, -0.01608, 0.03177, $ -0.12241, -0.01370, 0.08606, 0.38241] k = 2 pdrec[k].ver = 1016 pdrec[k].ts = anytim('2009-06-01T00:00') pdrec[k].te = anytim('2011-03-01T00:00') pdrec[k].t0 = anytim('2009-03-03T00:00') pdrec[k].c1 = [ -7.40505, 6.93326, -0.06104, -0.46446, -0.04072, $ -0.09851, 0.00199, -0.02096, 0.10000] pdrec[k].c2 = [ -3.09684, 1.37612, 0.03446, 0.08630, -0.01247, $ -0.04378, 0.04979, 0.10301, 0.49868] k = 3 pdrec[k].ver = 1016 pdrec[k].ts = anytim('2011-03-01T00:00') pdrec[k].te = anytim('2099-01-01T00:00') pdrec[k].t0 = anytim('2010-12-01T00:00') pdrec[k].c1 = [ -8.56397, 7.15365, -0.49817, 0.52057, 0.02147, $ -0.05027, -0.00185, -0.02525, 0.10000] pdrec[k].c2 = [ -2.68742, 0.23191, 0.14165, -0.38149, 0.09431, $ 0.05550, 0.09927, 0.07492, 4.66914] k = 4 pdrec[k].ver = 1021 pdrec[k].ts = anytim('2015-05-19T00:00') ; post conjunction pdrec[k].te = anytim('2023-08-12T00:00') ;s/c flip pdrec[k].t0 = anytim('2011-01-01T00:00') pdrec[k].c1 = -pdrec[3].c1 ; polarity change pdrec[k].c1[8] = pdrec[3].c1[8] ; no polarity change for exponent term pdrec[k].c2 = -pdrec[3].c2 ; polarity change pdrec[k].c2[8] = pdrec[3].c2[8] ; no polarity change for exponent term endif else err = 7B if err ne 7B then begin crpix0=[i.crpix1,i.crpix2] cdelt0=[i.cdelt1,i.cdelt2] er1_1 = cos(erol) ; convert EUVI roll into PC matrix er1_2 = -sin(erol) er2_1 = sin(erol) er2_2 = cos(erol) ; determine status and version of previously applied euvi_point, if any ; relevant keywords: ATT_FILE, VERSION, HISTORY Id: euvi_point.pro ; variables set: ver, usever, oldver, olderr, att_file if n_elements(make_ver) ne 1 then make_ver=strtrim(strmid(info,22,5),2) usever = strsplit(make_ver,' .',/extract) ;help,usever if n_elements(usever) gt 1 then $ usever=fix(usever[0])*1000+fix(usever[1]) else usever=9999 if tag_exist(i,'version') then begin ver = strsplit(i.version,' p',/extract) if n_elements(ver) gt 1 then oldver = ver[1] if strmid(i.version,0,1) ne 'p' then ver=ver[0] else ver='' endif if n_elements(oldver) eq 0 && tag_exist(i,'history') then begin whis = where(strmid(i.history,0,20) eq 'Id: euvi_point.pro,v',nwhis) if nwhis gt 0 then $ oldver = strtrim(strmid(i.history[whis[nwhis-1]],21,5),2) endif if n_elements(oldver) gt 0 then begin oldver = strsplit(oldver,' .',/extract) if n_elements(oldver) gt 1 then $ oldver=fix(oldver[0])*1000+fix(oldver[1]) else oldver=fix(oldver) endif else oldver=0 if tag_exist(i,'att_file') then att_file=i.att_file else att_file='' att_file = strsplit(att_file,'+',/extract) olderr = att_file[n_elements(att_file)-1] if strmid(olderr,1,2,/reverse) eq 'GT' then begin olderr = strmid(olderr,2,1,/reverse) case olderr of '0': olderr=0 '1': olderr=1 else: olderr=3 endcase endif else olderr=3 ; bare att file name w(ithout +?GT) or NotFound: if strlen(att_file[0]) lt 8 then att_file='NotFound' else att_file=att_file[0] ;help,make_ver,usever,oldver,ver,att_file,olderr ; get the sun vector relative to the GT boresight (and the err flag!) time = anytim(i.date_obs)+i.exptime+2.0 ; end exp + 2 sec dsun = i.dsun_obs / 1.496e11 ; sun distance in AU fpsoff = [i.fpsoffy,i.fpsoffz] ; FPS pointing offset svecg = scc_gt2sunvec(time,dsun,fpsoff,obs=obs_s,file=gfil,error=err, $ no_grh=no_grh) ; svec in arcsec svec = [svecg[0]*cos(grol)-svecg[1]*sin(grol), $ ; rotate svec from GT svecg[0]*sin(grol)+svecg[1]*cos(grol)] ; to EUVI system ; Decide whether we have enough information to actually improve the pointing: ; 1. if (err le 1) && (err le olderr) -> we can do at least as well or better ; 2. if (olderr le 1) && (oldver le 1.15) -> adding drift improves old point. if (err le 1) && (err le olderr) then do_point=1 else $ if olderr le 1 && oldver le 1015 && oldver ge 1005 && usever gt 1015 $ then do_point=2 else do_point=0 ; Recalculate pointing and/or add pointing drift correction if do_point gt 0 then begin IF (loud) THEN message,'Recalculating pointing.',/info IF (loud) THEN help,do_point,err,olderr,oldver,usever,att_file ; find the appropriate pointing drift coefficient record wpd=where(pdrec.ver le usever and pdrec.ts lt time and pdrec.te ge time,npd) if npd lt 1 then wpd=0 else wpd=wpd[npd-1] ; use latest if multiple pd = pdrec[wpd] ; calculate the pointing correction ; consists of an offset c[0], an exponential c[1]*exp(-c[8]*t), ; and 3 sine an cosine terms at the 1st, 2nd, and 3rd harmonic of torb tpd = (time-pd.t0) / pd.torb ; orb. phase pd1 = pd.c1[0] + pd.c1[1]*exp(-pd.c1[8]*tpd) ; const+exp term pd2 = pd.c2[0] + pd.c2[1]*exp(-pd.c2[8]*tpd) ; const+exp term for k=1,3 do begin ppd = 2.0*!pi*k*tpd pd1 = pd1 + pd.c1[2*k]*sin(ppd) + pd.c1[2*k+1]*cos(ppd) ; sin,cos terms pd2 = pd2 + pd.c2[2*k]*sin(ppd) + pd.c2[2*k+1]*cos(ppd) ; sin,cos terms endfor ; get the summing factors, assuming image is rectified: sumfac1 = 2.^(i.summed - 1) ;Amount of image summing sumfac2 = sumfac1 IF (loud) THEN help,sumfac1 ; recalculate crpix from scratch if feasible (i.e., if do_point eq 1) ; otherwise (if do_point eq 2) use existing crpix (and add pointing drift) if do_point eq 1 then begin ; calculate the basic crpix (excl. pointing drift) ; The calculation has 3 steps, each using one line in the statements below: ; 1. Calculate (IDL) full resolution pixel coordinate of GT axis. ; This is essentially OFF1-P1ROW or OFF2-P1COL. ; If the axis had to be flipped during the rectification process, then ; this becomes, e.g., NAXIS1*SUMMED-1 - (OFF1-PROW). ; The actual equation uses FLP1 and FLP2 to incorporate both cases. ; 2. Apply the GT bias, and the FPS bias. The sign is opposite for A and B: ; +Y of Ahead points N, +Y of Behind points S ; +Z of Ahead points W, +Z of Behind points E ; a positive GT bias forces a positive Y/Z sun vector component: ; - on Ahead it moves suncenter N/W on the CCD, ; - on Behind it moves suncenter S/E on the CCD. ; 3. Apply a correction for summed pixels, divide by the summing, ; and add 1 (CRPIX uses FITS standard = IDL standard + 1). i.crpix1 = flp1*i.naxis1 + ((1.-2.*flp1)*(off1-i.p1row)-flp1 $ + (1.-2.*obs)*svec[1]/p2a $ - (sumfac1-1.)/2.) / sumfac1 + 1. i.crpix2 = flp2*i.naxis2 + ((1.-2.*flp2)*(off2-i.p1col)-flp2 $ + (1.-2.*obs)*svec[0]/p2a $ - (sumfac2-1.)/2.) / sumfac2 + 1. endif ; subtract pointing drift i.crpix1 = i.crpix1 - pd1/sumfac1 i.crpix2 = i.crpix2 - pd2/sumfac2 ; make sure crval1,2 are zero i.crval1 = 0.0 i.crval2 = 0.0 ; calculate cdelt. i.cdelt1 = p2a * sumfac1 i.cdelt2 = p2a * sumfac2 ; ; Set instrument offset keywords ; IF tag_exist(i,'ins_x0') THEN BEGIN ; i.ins_x0 = p2a*((2048+1)/2. - off1) ; pre 1.16 version ; i.ins_y0 = p2a*((2048+1)/2. - off2) ; pre 1.16 version ; polarity now post-rectify and includes pointing drift i.ins_x0 = p2a*((1.-2.*flp1)*((2048+1)/2.-off1)+pd1) ; arcsec i.ins_y0 = p2a*((1.-2.*flp2)*((2048+1)/2.-off2)+pd2) ; arcsec ENDIF if tag_exist(i,'crpix1a') then begin i.crpix1a = i.crpix1 i.crpix2a = i.crpix2 endif if tag_exist(i,'cdelt1a') then begin i.cdelt1a = -i.cdelt1/3.6e3 i.cdelt2a = i.cdelt2/3.6e3 endif ; Obtain RA/Dec coordinates if icy && tag_exist(i,'crval1a') && (oldver lt 1008 || keyword_set(forcr)) $ then begin yp=[0d,0d] convert_stereo_lonlat,i.date_obs,yp,'HPC','GEI',spacecraft=obs_s,/degrees i.crval1a = yp[0] i.crval2a = yp[1] endif ; Show changes. dcrpix=[i.crpix1,i.crpix2]-crpix0 dcdelt=[i.cdelt1,i.cdelt2]-cdelt0 crmsg='delta crpix=['+trim(dcrpix[0])+','+trim(dcrpix[1])+'] pixels' cdmsg='delta cdelt=['+trim(dcrpix[0])+','+trim(dcrpix[1])+'] arcsec' IF (loud) THEN BEGIN message,crmsg,/info IF (abs(dcdelt[0]+abs(dcdelt[1]) GT .0001)) THEN message,cdmsg,/info ENDIF ; update version and history keywords if tag_exist(i,'version') then if (size(i.version))[1] eq 7 then begin i.version = ver+'p'+strtrim(usever/1000,2)+'.'+strtrim(usever mod 1000,2) endif if tag_exist(i,'history') then begin len=strlen(info) histinfo=strmid(info,1,len-2) i=scc_update_history(i, histinfo) IF total(crpix0) NE 0 THEN i=scc_update_history(i,crmsg) IF total(cdelt0) NE 0 and (abs(dcdelt[0]+abs(dcdelt[1]) GT .0001)) THEN i=scc_update_history(i,cdmsg) endif endif ; end of conditional pointing update ; try to get SPICE roll info if: ; - att_file eq 'NotFound', or ; - olderr gt 1 ; - version lt 1.5 ; - keyword forceroll set, or ; - Sanity check between PC, CROTA and SC_ROLL fails ; if icy not available, and sanity check failed: ; - replace PC and CROTA with SC_ROLL derived number ; sanity check between existing PC, CROTA and SC_ROLL: roll_bad = 0 if tag_exist(i,'sc_roll') then begin scrota = i.sc_roll + erol*!radeg pcrota = !radeg*atan(i.pc2_1,i.pc1_1) if abs(scrota-pcrota) gt 0.01 then roll_bad=1 if tag_exist(i,'crota') then $ if abs(scrota-i.crota) gt 0.01 then roll_bad=1 endif ; figure out whether we want to and can do a roll update with SPICE if (att_file eq 'NotFound' || olderr gt 1 || oldver lt 1005 || $ keyword_set(forcr) || roll_bad) && icy then begin att_file_spice = get_stereo_att_file(i.date_obs,obs_s) if att_file_spice eq '' then begin att_file_spice = 'NotFound' ; do_roll = 0 ; do not use SPICE if att file not found do_roll = 1 ; use SPICE anyway: predicted roll better than no roll endif else do_roll=1 endif else do_roll=0 if (att_file eq 'NotFound' || olderr gt 1 || oldver lt 1005 || $ keyword_set(forcr) || roll_bad) && ~icy then rollwarningflag=1 ; do the roll update if do_roll || roll_bad then begin if do_roll then begin ; use SPICE IF (loud) THEN message,'Re-calculating roll from SPICE.',/info att_file = att_file_spice if tag_exist(i,'pc1_1a') then begin rollrada = get_stereo_roll(i.date_obs,obs_s,system='GEI') / !radeg i.pc1_1a = cos(rollrada) ; convert roll into PC matrix i.pc1_2a = -sin(rollrada) i.pc2_1a = sin(rollrada) i.pc2_2a = cos(rollrada) pc1_1 = i.pc1_1a pc1_2 = i.pc1_2a pc2_1 = i.pc2_1a pc2_2 = i.pc2_2a ; calculate updated pc matrix: equivalent to er # pc. i.pc1_1a = er1_1 * pc1_1 + er1_2 * pc2_1 i.pc1_2a = er1_1 * pc1_2 + er1_2 * pc2_2 i.pc2_1a = er2_1 * pc1_1 + er2_2 * pc2_1 i.pc2_2a = er2_1 * pc1_2 + er2_2 * pc2_2 endif hcv = get_stereo_hpc_point(i.date_obs,obs_s) rollrad = hcv[2]/!radeg if tag_exist(i,'crota') then $ i.crota = (hcv[2]+erol*!radeg+180.0) mod 360.0 - 180.0 endif else begin ; use SC_ROLL keyword instead rollrad = i.sc_roll/!radeg if tag_exist(i,'crota') then $ i.crota = (i.sc_roll+erol*!radeg+180.0) mod 360.0 - 180.0 endelse i.pc1_1 = cos(rollrad) ; convert roll into PC mtrx i.pc1_2 = -sin(rollrad) i.pc2_1 = sin(rollrad) i.pc2_2 = cos(rollrad) pc1_1 = i.pc1_1 pc1_2 = i.pc1_2 pc2_1 = i.pc2_1 pc2_2 = i.pc2_2 ; calculate updated pc matrix: equivalent to er # pc. i.pc1_1 = er1_1 * pc1_1 + er1_2 * pc2_1 i.pc1_2 = er1_1 * pc1_2 + er1_2 * pc2_2 i.pc2_1 = er2_1 * pc1_1 + er2_2 * pc2_1 i.pc2_2 = er2_1 * pc1_2 + er2_2 * pc2_2 if tag_exist(i,'ins_r0') then i.ins_r0 = erol*!radeg ; deg endif ; update att_file if tag_exist(i,'att_file') then begin if do_point ne 1 then err=olderr ; use old error code if pnt not recalc. i.att_file = att_file+'+'+strtrim(fix(err),2)+'GT' endif ; unfortunately, SECCHI supports a slew of secondary and obsolete keywords: if tag_exist(i,'xcen') then begin ceni = (i.naxis1+1)/2. - i.crpix1 cenj = (i.naxis2+1)/2. - i.crpix2 i.xcen = i.cdelt1*(i.pc1_1*ceni+i.pc1_2*cenj) i.ycen = i.cdelt2*(i.pc2_1*ceni+i.pc2_2*cenj) endif ii[j] = i endif errs[j] = err endfor if rollwarningflag && ~keyword_set(quiet) then $ print,'No roll update, spice/icy dlm not available' if n_elements(errs) eq 1 then errs = errs[0] end