/usr/share/gnudatalanguage/astrolib/substar.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 | pro substar,image,x,y,mag,id,psfname,VERBOSE = verbose ;Subtract scaled PSF stars
;+
; NAME:
; SUBSTAR
; PURPOSE:
; Subtract a scaled point spread function at specified star position(s).
; EXPLANATION:
; Part of the IDL-DAOPHOT photometry sequence
;
; CALLING SEQUENCE:
; SUBSTAR, image, x, y, mag, [ id, psfname, /VERBOSE]
;
; INPUT-OUTPUT:
; IMAGE - On input, IMAGE is the original image array. A scaled
; PSF will be subtracted from IMAGE at specified star positions.
; Make a copy of IMAGE before calling SUBSTAR, if you want to
; keep a copy of the unsubtracted image array
;
; INPUTS:
; X - REAL Vector of X positions found by NSTAR (or FIND)
; Y - REAL Vector of Y positions found by NSTAR (or FIND)
; MAG - REAL Vector of stellar magnitudes found by NSTAR (or APER)
; Used to scale the PSF to match intensity at star position.
; Stars with magnitude values of 0.0 are assumed missing and
; ignored in the subtraction.
;
; OPTIONAL INPUTS:
; ID - Index vector indicating which stars are to be subtracted. If
; omitted, (or set equal to -1), then stars will be subtracted
; at all positions specified by the X and Y vectors.
;
; PSFNAME - Name of the FITS file containing the PSF residuals, as
; generated by GETPSF. SUBSTAR will prompt for this parameter
; if not supplied.
;
; OPTIONAL INPUT KEYWORD:
; VERBOSE - If this keyword is set and nonzero, then SUBSTAR will
; display the star that it is currently processing
;
; COMMON BLOCKS:
; The RINTER common block is used (see RINTER.PRO) to save time in the
; PSF calculations
;
; PROCEDURES CALLED:
; DAO_VALUE(), READFITS(), REMOVE, SXOPEN, SXPAR(), SXREAD()
; REVISION HISTORY:
; Written, W. Landsman August, 1988
; Added VERBOSE keyword January, 1992
; Fix star subtraction near edges, W. Landsman May, 1996
; Assume the PSF file is in FITS format W. Landsman July, 1997
; Converted to IDL V5.0 W. Landsman September 1997
;-
common rinter,c1,c2,c3,init ;Save time in RINTER
if N_params() LT 4 then begin
print,'Syntax - SUBSTAR, image, x, y, mag,[ id, psfname, /VERBOSE]'
return
endif
s = size(image)
if s[0] NE 2 then $
message, 'ERROR - Input array (first parameter) must be 2 dimensions'
npts = N_elements(image)
if N_elements(psfname) NE 1 then begin
psfname = ''
read, 'Enter name of the FITS file containing PSF residuals: ', psfname
endif
if N_params() LT 5 then id = indgen( N_elements(x) ) else begin
if min(id) LT 0 then id = indgen( N_elements(x) ) ;Subtract all stars?
endelse
psf = readfits(psfname, hpsf)
nstar = N_elements(id) ;Number of stars to subtract
gauss = sxpar( hpsf, 'GAUSS*' )
psfmag = sxpar( hpsf, 'PSFMAG' )
psfrad = sxpar( hpsf, 'PSFRAD' )
fitrad = sxpar( hpsf, 'FITRAD' )
npsf = sxpar( hpsf, 'NAXIS1' )
nbox = ( 2*fix( psfrad + 0.5 ) + 1) > ((npsf-7)/2)
nhalf = (nbox-1)/2
psfrsq = psfrad^2
lx = fix( x[id] + 0.5 ) - nhalf
ly = fix( y[id] + 0.5 ) - nhalf
smag = mag[id]
scale = 10^(-0.4*(smag- psfmag))
xx = x[id] - lx
yy = y[id] - ly
bad = where( (smag EQ 0.0), Nbad) ;Any stars with missing magnitudes?
if Nbad GT 0 then begin
nstar = nstar - Nbad
remove,bad,lx,ly,xx,yy,scale
endif
rsq = fltarr( nbox, nbox)
boxgen = indgen(nbox)
; Compute RINTER common block arrays
p_1 = shift(psf,1,0) & p1 = shift(psf,-1,0) & p2 = shift(psf,-2,0)
c1 = 0.5*(p1-p_1)
c2 = 2.*p1 + p_1 - 0.5*(5.*psf + p2)
c3 = 0.5 *(3.*(psf-p1) + p2 - p_1)
init = 1
verbose = keyword_set(VERBOSE)
cr = string("15b)
for i = 0L,nstar-1 do begin
dx = boxgen - xx[i]
dy = boxgen - yy[i]
dx2 = dx^2 & dy2 = dy^2
for j = 0,nbox-1 do rsq[0,j] = dx2 + dy2[j]
good = where( rsq LT psfrsq)
xgood = good mod nbox & ygood = good/nbox
dx = dx[xgood] & dy = dy[ygood]
goodbig = ( xgood + lx[i] ) + ( ygood + ly[i] )*s[1]
bad = where( (goodbig LT 0) or (goodbig GE npts), Nbad)
if nbad GT 0 then remove,bad,goodbig,dx,dy
image[goodbig] = image[goodbig] - scale[i] * dao_value( dx,dy,gauss,psf )
if VERBOSE then $
print,f="($,'SUBSTAR: Processing Star',I5,A)",id[i],cr
endfor
return
end
|