pro mkgrid, n, x, y, c=c, d=d, s=s ;+ ; NAME: ; mkgrid ;; ; PURPOSE: ; Make arrays containing the values of the x and y coordinates of a grid. ; ; INPUTS: ; n1,n2 - Integers. Array dimensions (number of pixels). ; ; OPTIONAL INPUTS: ; c0 - Floating point array [c1,c2] or c1 in both dimensions. Pixel coordinate of the ; [0.0,0.0] origin. Default: the centre of the array. ; dist - make distance array. each pixel contains distance to origin of the array. ; ; KEYWORD PARAMETERS: ; None. ; ; OUTPUTS: ; x,y - Floating point [n1*n2] array containing x or y values (all of the values ; in the other dimension are identical) ;- COMPILE_OPT idl2 ; Set reference point to the centre if not specified if n_elements(n) eq 1 then dim=[n,n] else dim=n>1 if n_elements(c) eq 0 then c=.5*float(dim)-0.5 else c=[0,0]+c ; Make x and y arrays ; but only if the actual keyword is set if keyword_set(routine_names(x,arg_name=-1)) then x = rebin(findgen(dim[0])-c[0],dim[0],dim[1],/sample) if keyword_set(routine_names(y,arg_name=-1)) then y = rebin(findgen(1,dim[1])-c[1],dim[0],dim[1],/sample) ; for distance array it is faster to recalculate than to take the sqrt of the x and y arrays if keyword_set(routine_names(d,arg_name=-1)) then \$ d=sqrt(rebin((findgen(dim[0])-c[0])^2,dim[0],dim[1],/sample) + rebin((findgen(1,dim[1])-c[1])^2,dim[0],dim[1],/sample)) ; even faster return the d^2 arrays if keyword_set(routine_names(s,arg_name=-1)) then s=rebin((findgen(dim[0])-c[0])^2,dim[0],dim[1],/sample) + rebin((findgen(1,dim[1])-c[1])^2,dim[0],dim[1],/sample) end pro getpsf, fpsf, ra, dec, psf, pa=pa, display=display, pixscl=pixscl, tsz=tsz, totcor=totcor ;+ ; NAME: ; getpsf ;; ; PURPOSE: ; Reconstruct the IRAC PSF at the specific location using the kinformation stored in the GREATS psf file. ; ; INPUTS: ; fpsf - string. Name of the psf file ; ra,dec - Floating point array. List of positions, in degrees, where the PSF need to be reconstruted. ; ; OPTIONAL INPUTS: ; pa - Floating point array. List of additional position angles requested for the final PSF. Default: 0.0 ; pixscl - Floating point value. Pixel scale of the output PSF. Default: adopt the same pixels scale as in the GREATS psf file ; tsz - Floating point value. Size of the output PSF, in pixels. Default: adopt the same size as the master PSF in the GREATS psf file (81 pixels). ; ; KEYWORD PARAMETERS: ; display - Keyword. If set, the PSF is displayed at the end of the process. Default: False. ; ; OUTPUTS: ; psf - Floating point array containing the PSFs. ;- common getpsf_aor_com, mpsf, hpsf, map, hmap, current_psf if not keyword_set(pa) then pa = 0.0 if n_elements(current_psf) le 0 then current_psf='' if not strcmp(current_psf,fpsf) then begin print,'getpsf --> restoring ',fpsf mpsf = readfits(fpsf, hpsf, exten=1,/silent) map = readfits(fpsf, hmap, exten=2, /silent) current_psf = fpsf end szpsf = size(mpsf,/dim) sz = size(map,/dim) len = sz[3] xy0 = (szpsf-1)/2. mkgrid, szpsf, d=d if keyword_set(pixscl) then begin if abs(pixscale(hpsf)/pixscl-1) gt 1e-3 then begin rscl = pixscale(hpsf)/pixscl if keyword_set(tsz) then sznew=tsz else sznew = oddsize(round(rscl*szpsf[0]),div=3) _gpsf=frebin(mpsf,rscl*szpsf[0],rscl*szpsf[0]) szz=size(_gpsf,/dim) sznew=oddsize(szz[0],div=3) if sznew ne szz[0] then begin tmp=fltarr([sznew, sznew]) tmp[0,0]=_gpsf _gpsf=tmp endif gcntrd,_gpsf,sznew/2.,sznew/2.,xc,yc,4 gpsf=fshift(_gpsf,sznew/2.-xc-0.5,sznew/2.-yc-0.5) end else begin gpsf = mpsf rscl=1.00 sznew = oddsize(round(rscl*szpsf[0]),div=3) ;endif endelse endif else begin gpsf = mpsf rscl=1.00 sznew = oddsize(round(rscl*szpsf[0]),div=3) end adxy, hmap, ra[*], dec[*], x, y ix = round(x) iy = round(y) nx=n_elements(x) psf = fltarr([size(gpsf,/dim),nx]) tmp = fltarr(size(gpsf,/dim)) for j=0,nx-1 do begin if total(map[ix,iy,0,1:*]) eq 0 then return mkgrid, sznew, d=d for i=0,len-1 do begin aor_w = map[ix,iy,0,i] aor_pa = map[ix,iy,1,i] if aor_w[0] le 0 then continue tmp += rot(gpsf,-aor_pa[0]-pa,rscl[0],xy0[0],xy0[1],cubic=-0.5,missing=0.0)*aor_w[0]/rscl^2 if keyword_set(display) then begin print, i+1, aor_w, aor_pa, format='(i,3f)' tvs, rot(gpsf,-aor_pa[0]-pa,rscl,xy0[0],xy0[1],cubic=-0.5,missing=0.0), fac=10, pos=0 tvs, tmp, fac=10, pos=1 end end psf[*,*,j] = tmp/total(tmp*(d lt (sznew-1.0)/2.)) end if keyword_set(display) then print, ra,dec, rscl, total(psf) end