/usr/share/gnudatalanguage/astrolib/hrebin.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 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208 209 210 211 212 213 214 215 216 217 218 219 220 221 222 223 224 225 226 227 228 229 230 231 232 233 234 235 236 237 238 239 240 241 242 243 244 245 246 247 248 249 250 251 252 253 254 255 256 257 258 259 260 261 262 263 264 265 266 267 268 269 270 271 272 273 274 275 276 277 | pro hrebin, oldim, oldhd, newim, newhd, newx, newy, TOTAL = total, $
SAMPLE=sample, OUTSIZE = outsize, ERRMSG = errmsg, ALT=alt
;+
; NAME:
; HREBIN
; PURPOSE:
; Expand or contract a FITS image using (F)REBIN and update the header
; EXPLANATION:
; If the output size is an exact multiple of the input size then REBIN is
; used, else FREBIN is used. User can either overwrite the input array,
; or write to new variables. By default, the counts/pixel is preserved,
; though one can preserve the total counts or surface flux by setting /TOTAL
;
; CALLING SEQUENCE:
; HREBIN, oldhd ;Special calling sequence to just update header
; HREBIN, oldim, oldhd, [ newim, newhd, newx, newy, OUTSIZE = ,/SAMPLE,
; ERRMSG = ]
;
; INPUTS:
; OLDIM - the original image array
; OLDHD - the original image FITS header, string array
;
; OPTIONAL INPUTS:
; NEWX - size of the new image in the X direction, integer scalar
; NEWY - size of the new image in the Y direction, integer scalar
; HREBIN will prompt for NEWX and NEWY if not supplied
;
; OPTIONAL OUTPUTS:
; NEWIM - the image after expansion or contraction with REBIN
; NEWHD - header for newim containing updated astrometry info
; If output parameters are not supplied, the program will modify
; the input parameters OLDIM and OLDHD to contain the new array and
; updated header.
;
; OPTIONAL INPUT KEYWORDS:
; ALT - Single character 'A' through 'Z' or ' ' specifying which astrometry
; system to modify in the FITS header. The default is to use the
; primary astrometry of ALT = ' '. See Greisen and Calabretta (2002)
; for information about alternate astrometry keywords.
;
; OUTSIZE - Two element integer vector which can be used instead of the
; NEWX and NEWY parameters to specify the output image dimensions
;
; /SAMPLE - Expansion or contraction is done using REBIN which uses
; bilinear interpolation when magnifying and boxaveraging when
; minifying. If the SAMPLE keyword is supplied and non-zero,
; then nearest neighbor sampling is used in both cases. Keyword
; has no effect when output size is not a multiple of input size.
;
; /TOTAL - If set then the output image will have the same total number of counts
; as the input image. Because HREBIN also updates the astrometry,
; use of the TOTAL keyword also preserves counts per surface area, e.g.
; counts/(arc sec)@
;
; OPTIONAL KEYWORD OUTPUT:
; ERRMSG - If this keyword is supplied, then any error mesasges will be
; returned to the user in this parameter rather than depending on
; on the MESSAGE routine in IDL. If no errors are encountered
; then a null string is returned.
; PROCEDURE:
; The parameters BSCALE, NAXIS1, NAXIS2, CRPIX1, and CRPIX2 and the CD
; (or CDELT) parameters are updated for the new FITS header.
;
; EXAMPLE:
; Compress a 2048 x 2048 image array IM, with FITS header HDR, to a
; 724 x 724 array. Overwrite the input variables with the compressed
; image and header.
;
; IDL> hrebin, im, hdr, OUT = [724, 724]
;
; PROCEDURES USED:
; CHECK_FITS, EXTAST, FREBIN, GSSS_STDAST, STRN(), SXPAR(), SXADDHIST,
; SXADDPAR, ZPARCHECK
;
; MODIFICATION HISTORY:
; Written, December 1990 W. Landsman, ST System Corp.
; Update CD1_1 keywords W. Landsman November 1992
; Check for a GSSS header W. Landsman June 1994
; Update BSCALE even if no astrometry present W. Landsman May 1997
; Converted to IDL V5.0 W. Landsman September 1997
; Use FREBIN to accept sizes that are not a integer multiple of the original
; size W. Landsman August 1998
; Correct for "edge" effects when expanding with REBIN W. Landsman Apr. 1999
; Fixed initialization of header only call broken in Apr 98 change May. 1999
; Remove reference to obsolete !ERR W. Landsman February 2000
; Use double precision formatting for CD matrix W. Landsman April 2000
; Recognize PC00n00m astrometry format W. Landsman December 2001
; Correct astrometry for integral contraction W. Landsman April 2002
; Fix output astrometry for non-equal plate scales for PC matrix or
; CROTA2 keyword, added ALT keyword. W. Landsman May 2005
; Update distortion parameters if present W. Landsman August 2007
; Don't update BSCALE/BZERO for unsigned integer W.Landsman Mar 2008
; Use post-V6.0 notation W. Landsman Nov 2011
; Write CRPIX values as double precision if necessary W. Landsman Oct. 2012
; Always call FREBIN, added TOTAL keyword W. Landsman Nov 2015
;-
On_error,2
compile_opt idl2
npar = N_params() ;Check # of parameters
if (npar EQ 3) || (npar EQ 5) || (npar EQ 0) then begin
print,'Syntax - HREBIN, oldim, oldhd,[ newim, newhd, OUTSIZE=, ' + $
'/SAMPLE, ERRMSG= ]'
return
endif
if ~keyword_set(SAMPLE) then sample = 0
save_err = arg_present(errmsg) ;Does user want to return error messages?
; If only 1 parameter is supplied, then assume it is a FITS header
if ( npar EQ 1 ) then begin
zparcheck, 'HREBIN', oldim, 1, 7, 1, 'Image header'
oldhd = oldim
xsize = sxpar( oldhd,'NAXIS1' )
ysize = sxpar( oldhd,'NAXIS2' )
endif else begin
check_FITS, oldim, oldhd, dimen, /NOTYPE, ERRMSG = errmsg
if errmsg NE '' then begin
if ~save_err then message,'ERROR - ' + errmsg,/CON
return
endif
if N_elements(dimen) NE 2 then begin
errmsg = 'Input image array must be 2-dimensional'
if ~save_err then message,'ERROR - ' + errmsg,/CON
return
endif
xsize = dimen[0] & ysize = dimen[1]
endelse
tname = size(oldim,/tname)
if ( npar LT 6 ) then begin
if ( N_elements(OUTSIZE) NE 2 ) then begin
tit = !MSG_PREFIX + 'HREBIN: '
print, tit, 'Original array size is '+ strn(xsize) + ' by ' + strn(ysize)
read, tit + 'Enter size of new image in the X direction: ',newx
read, tit + 'Enter size of new image in the Y direction: ',newy
endif else begin
newx = outsize[0]
newy = outsize[1]
endelse
endif
; Modified Nov 2015 to always call FREBIN. FREBIN() will call the IDL REBIN()
; function if we are changing dimensions by an exact multiple.
if npar GT 1 then begin
if npar GT 2 then newim = frebin( oldim, newx, newy,total=total) $
else oldim = frebin( oldim, newx, newy,total=total)
endif
if ( sample GT 0 ) then type = ' Nearest Neighbor Approximation' else begin
if ( newx LT xsize ) then type = ' Box Averaging' else $
type = ' Bilinear Interpolation'
endelse
newhd = oldhd
sxaddpar, newhd, 'NAXIS1', fix(newx)
sxaddpar, newhd, 'NAXIS2', fix(newy)
label = 'HREBIN: '+ strmid( systime(),4,20 )
sxaddpar,newhd,'history',label + ' Original Image Size Was '+ $
strn(xsize) +' by ' + strn(ysize)
if ( npar GT 1 ) then sxaddpar,newhd,'history',label+type
xratio = float(newx) / xsize ;Expansion or contraction in X
yratio = float(newy) / ysize ;Expansion or contraction in Y
lambda = yratio/xratio ;Measures change in aspect ratio.
pix_ratio = xratio*yratio ;Ratio of pixel areas
; Update astrometry info if it exists
extast, newhd, astr, noparams, ALT = alt
if noparams GE 0 then begin
if strmid(astr.ctype[0],5,3) EQ 'GSS' then begin
gsss_stdast, newhd
extast, newhd, astr, noparams
endif
; Correct the position of the reference pixel. Note that CRPIX values are
; given in FORTRAN (first pixel is (1,1)) convention
crpix = astr.crpix
; When expanding with REBIN with bilinear interpolation (SAMPLE = 0), edge
; effects are introduced, which require a different calculation of the updated
; CRPIX1 and CRPIX2 values.
exact = (~(xsize mod newx) || ~(newx mod xsize)) && $
(~(ysize mod newy) || ~(newy mod ysize))
if (exact) && (~keyword_set(SAMPLE)) && (xratio GT 1) then $
crpix1 = (crpix[0]-1.0)*xratio + 1.0 else $
crpix1 = (crpix[0]-0.5)*xratio + 0.5
if (exact) && (~keyword_set(SAMPLE)) && (yratio GT 1) then $
crpix2 = (crpix[1]-1.0)*yratio + 1.0 else $
crpix2 = (crpix[1]-0.5)*yratio + 0.5
if N_elements(alt) EQ 0 then alt = ''
sxaddpar, newhd, 'CRPIX1' + alt, crpix1
sxaddpar, newhd, 'CRPIX2' + alt, crpix2
if tag_exist(astr,'DISTORT') then begin
distort = astr.distort
message,'Updating SIP distortion parameters',/INF
update_distort,distort, [1./xratio,0],[1./yratio,0]
astr.distort= distort
add_distort, newhd, astr
endif
; Scale either the CDELT parameters or the CD1_1 parameters.
if (noparams NE 2) then begin
cdelt = astr.cdelt
sxaddpar, newhd, 'CDELT1' + alt, CDELT[0]/xratio
sxaddpar, newhd, 'CDELT2' + alt, CDELT[1]/yratio
; Adjust the PC matrix if aspect ratio has changed. See equation 187 in
; Calabretta & Greisen (2002)
if lambda NE 1.0 then begin
cd = astr.cd
if noparams EQ 1 then begin
;Can no longer use the simple CROTA2 convention, change to PC keywords
sxaddpar,newhd,'PC1_1'+alt, cd[0,0]
sxaddpar, newhd,'PC2_2'+alt, cd[1,1]
sxdelpar, newhd, ['CROTA2','CROTA1']
endif
sxaddpar, newhd, 'PC1_2'+alt, cd[0,1]/lambda
sxaddpar, newhd, 'PC2_1'+alt, cd[1,0]*lambda
endif
endif else begin ;CDn_m Matrix format
cd = astr.cd
sxaddpar, newhd, 'CD1_1'+alt, cd[0,0]/xratio
sxaddpar, newhd, 'CD1_2'+alt, cd[0,1]/yratio
sxaddpar, newhd, 'CD2_1'+alt, cd[1,0]/xratio
sxaddpar, newhd, 'CD2_2'+alt, cd[1,1]/yratio
endelse
endif
; Adjust BZERO and BSCALE for new pixel size, unless these values are used
; to define unsigned integer data types.
if ~keyword_set(TOTAL) then begin
bscale = sxpar( oldhd, 'BSCALE')
bzero = sxpar( oldhd, 'BZERO')
unsgn = (tname EQ 'UINT') || (tname EQ 'ULONG')
if ~unsgn then begin
if (bscale NE 0) && (bscale NE 1) then $
sxaddpar, newhd, 'BSCALE', bscale/pix_ratio, 'Calibration Factor'
if (bzero NE 0) then sxaddpar, newhd, 'BZERO', bzero/pix_ratio, $
' Additive Constant for Calibration'
endif
endif
pixelsiz = sxpar( oldhd,'PIXELSIZ' , Count = N_pixelsiz)
if N_pixelsiz GT 0 then sxaddpar, newhd, 'PIXELSIZ', pixelsiz/xratio
if npar EQ 2 then oldhd = newhd else $
if npar EQ 1 then oldim = newhd
return
end
|