/usr/share/gnudatalanguage/astrolib/gsssadxy.pro is in gdl-astrolib 2018.02.16+dfsg-1.
This file is owned by root:root, with mode 0o644.
The actual contents of the file can be viewed below.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 | pro GSSSadxy,gsa,ra,dec,x,y, PRINT = print
;+
; NAME:
; GSSSADXY
; PURPOSE:
; Converts RA and DEC (J2000) to (X,Y) for an STScI GuideStar image.
; EXPLANATION:
; The sky coordinates may be printed and/or returned in variables.
;
; CALLING SEQUENCE:
; GSSSADXY, GSA, Ra,Dec, [ X, Y, /Print ]
; INPUT:
; GSA - the GSSS Astrometry structure created by GSSSEXTAST
; RA - the RA coordinate(s) in *degrees*, scalar or vector
; DEC - the DEC coordinate(s) in *degrees*, scalar or vector
;
; OPTIONAL KEYWORD INPUT:
; /PRINT - If this keyword is set and non-zero, then coordinates will be
; displayed at the terminal
; OUTPUT:
; X - the corresponding X pixel coordinate(s), double precision
; Y - the corresponding Y pixel coordinate(s), double precision
;
; X and Y will be in IDL convention (first pixel 0,0)
; EXAMPLE:
; Given a FITS header, hdr, from the STScI Guidestar Survey, determine
; the X,Y coordinates of 3C 273 (RA = 12 29 6.7 +02 03 08)
;
; IDL> GSSSEXTAST, hdr, gsa ;Extract astrometry structure
; IDL> GSSSADXY, gsa, ten(12,29,6.7)*15,ten(2,3,8),/print
;
; NOTES:
; For most purpose users can simply use ADXY, which will call GSSSADXY
; if it is passed a GSSS header.
;
; PROCEDURES CALLED:
; ASTDISP - Print RA, Dec in standard format
; HISTORY:
; 10-JUL-90 Version 1 written by Eric W. Deutsch
; Derived from procedures written by Brian McLean
; Vectorized code W. Landsman March, 1991
; 14-AUG-91 Fixed error which caused returned X and Y to be .5 pixels too
; large. Now X,Y follows same protocol as ADXY.
; June 1994 - Dropped PRFLAG parameter, added /PRINT W. Landsman (HSTX)
; Converted to IDL V5.0 W. Landsman September 1997
; 29-JUN-99 Added support for AMD[X,Y]1[2-3] for DSS images by E. Deutsch
; Reduce memory requirements for large arrays D. Finkbeiner April 2004
; Remove
;-
On_error,2
arg = N_params()
if (arg lt 5) then begin
print,'Syntax - GSSSADXY, GSSS_Astrom_struct, ra, dec, x, y, print_flag'
print,'e.g.: IDL> GSSSADXY, gsa, ra, dec, x, y, 1'
return
endif
; Set Constants
iters = 0 & maxiters=50 & tolerance=0.0000005
radeg = 180.0d/!DPI & arcsec_per_radian= 3600.0d*radeg
pltdec = gsa.crval[1]/radeg
dec_rad = dec/radeg
cosd = cos(dec_rad)
sind = sin(temporary(dec_rad))
ra_dif = ra/radeg - gsa.crval[0]/radeg
div = ( sind*sin(pltdec) + cosd*cos(pltdec)*cos(ra_dif))
xi = cosd*sin(ra_dif)*arcsec_per_radian/div
eta = ( sind*cos(pltdec)-cosd*sin(pltdec)*cos(ra_dif))* $
(arcsec_per_radian/temporary(div))
ra_dif = 0
cosd = 0 & sind = 0
obx = xi/gsa.pltscl
oby = eta/gsa.pltscl
repeat begin
iters++
f= gsa.amdx[0]*obx+ $
gsa.amdx[1]*oby+ $
gsa.amdx[2]+ $
gsa.amdx[3]*obx*obx+ $
gsa.amdx[4]*obx*oby+ $
gsa.amdx[5]*oby*oby+ $
gsa.amdx[6]*(obx*obx+oby*oby)+ $
gsa.amdx[7]*obx*obx*obx+ $
gsa.amdx[8]*obx*obx*oby+ $
gsa.amdx[9]*obx*oby*oby+ $
gsa.amdx[10]*oby*oby*oby+ $
gsa.amdx[11]*obx*(obx*obx+oby*oby)+ $
gsa.amdx[12]*obx*(obx*obx+oby*oby)^2
fx=gsa.amdx[0]+ $
gsa.amdx[3]*2.0*obx+ $
gsa.amdx[4]*oby+ $
gsa.amdx[6]*2.0*obx+ $
gsa.amdx[7]*3.0*obx*obx+ $
gsa.amdx[8]*2.0*obx*oby+ $
gsa.amdx[9]*oby*oby+ $
gsa.amdx[11]*(3.0*obx*obx+oby*oby)+ $
gsa.amdx[12]*(5.0*obx^4 + 6.0*obx^2*oby^2 + oby^4)
fy=gsa.amdx[1]+ $
gsa.amdx[4]*obx+ $
gsa.amdx[5]*2.0*oby+ $
gsa.amdx[6]*2.0*oby+ $
gsa.amdx[8]*obx*obx+ $
gsa.amdx[9]*obx*2.0*oby+ $
gsa.amdx[10]*3.0*oby*oby+ $
gsa.amdx[11]*2.0*obx*oby+ $
gsa.amdx[12]*(4.0*obx^3*oby + 4.0*obx*oby^3)
g= gsa.amdy[0]*oby+ $
gsa.amdy[1]*obx+ $
gsa.amdy[2]+ $
gsa.amdy[3]*oby*oby+ $
gsa.amdy[4]*oby*obx+ $
gsa.amdy[5]*obx*obx+ $
gsa.amdy[6]*(obx*obx+oby*oby)+ $
gsa.amdy[7]*oby*oby*oby+ $
gsa.amdy[8]*oby*oby*obx+ $
gsa.amdy[9]*oby*obx*obx+ $
gsa.amdy[10]*obx*obx*obx+ $
gsa.amdy[11]*oby*(obx*obx+oby*oby)+ $
gsa.amdy[12]*oby*(obx*obx+oby*oby)^2
gx=gsa.amdy[1]+ $
gsa.amdy[4]*oby+ $
gsa.amdy[5]*2.0*obx+ $
gsa.amdy[6]*2.0*obx+ $
gsa.amdy[8]*oby*oby+ $
gsa.amdy[9]*oby*2.0*obx+ $
gsa.amdy[10]*3.0*obx*obx+ $
gsa.amdy[11]*2.0*obx*oby+ $
gsa.amdy[12]*(4.0*obx^3*oby + 4.0*obx*oby^3)
gy=gsa.amdy[0]+ $
gsa.amdy[3]*2.0*oby+ $
gsa.amdy[4]*obx+ $
gsa.amdy[6]*2.0*oby+ $
gsa.amdy[7]*3.0*oby*oby+ $
gsa.amdy[8]*2.0*oby*obx+ $
gsa.amdy[9]*obx*obx+ $
gsa.amdy[11]*(3.0*oby*oby+obx*obx)+ $
gsa.amdy[12]*(5.0*oby^4 + 6.0*obx^2*oby^2 + obx^4)
f -= xi
g -= eta
deltx = (-f*gy+g*fy) / (fx*gy-fy*gx)
delty = (-g*fx+f*gx) / (fx*gy-fy*gx)
obx += deltx
oby += delty
;print,deltx,delty,tolerance,iters,maxiters
endrep until (min(abs([deltx,delty])) lt tolerance) || (iters gt maxiters)
eta = 0 & xi = 0 & deltx = 0 & delty = 0
x = (gsa.ppo3-obx*1000.0)/gsa.xsz-gsa.xll - 0.5
y = (gsa.ppo6+oby*1000.0)/gsa.ysz-gsa.yll - 0.5
if keyword_set(PRINT) then AstDisp, x, y, ra, dec
return
end
|