Secure your code as it's written. Use Snyk Code to scan source code in minutes - no build needed - and fix issues immediately.
pixel variances.
ds9: `pyds9.DS9`
Show intermediate steps of the drizzling
Returns
-------
hdu : `~astropy.io.fits.HDUList`
FITS HDUList with the drizzled 2D spectrum and weight arrays
"""
from drizzlepac.astrodrizzle import adrizzle
from astropy import log
log.setLevel('ERROR')
#log.disable_warnings_logging()
adrizzle.log.setLevel('ERROR')
NX = int(np.round(np.diff(wlimit)[0]*1.e4/dlam)) // 2
center = np.mean(wlimit[:2])*1.e4
out_header, output_wcs = utils.make_spectrum_wcsheader(center_wave=center,
dlam=dlam, NX=NX,
spatial_scale=spatial_scale, NY=NY)
sh = (out_header['NAXIS2'], out_header['NAXIS1'])
outsci = np.zeros(sh, dtype=np.float32)
outwht = np.zeros(sh, dtype=np.float32)
outctx = np.zeros(sh, dtype=np.int32)
outvar = np.zeros(sh, dtype=np.float32)
outwv = np.zeros(sh, dtype=np.float32)
outcv = np.zeros(sh, dtype=np.int32)
def apply_tweak_shifts(wcs_ref, shift_dict, grism_matches={}, verbose=True):
"""
"""
from drizzlepac import updatehdr
hdu = wcs_ref.to_fits(relax=True)
file0 = list(shift_dict.keys())[0].split('.fits')[0]
tweak_file = '{0}_tweak_wcs.fits'.format(file0)
hdu.writeto(tweak_file, clobber=True)
for file in shift_dict:
updatehdr.updatewcs_with_shift(file, tweak_file,
xsh=shift_dict[file][0],
ysh=shift_dict[file][1],
rot=0., scale=1.,
wcsname='SHIFT', force=True,
reusename=True, verbose=verbose,
sciext='SCI')
### Bug in astrodrizzle? Dies if the FLT files don't have MJD-OBS
### keywords
im = pyfits.open(file, mode='update')
im[0].header['MJD-OBS'] = im[0].header['EXPSTART']
im.flush()
# Update paired grism exposures
if file in grism_matches:
for grism_file in grism_matches[file]:
for drz_file in drz_files[1:]:
root = drz_file.split('_drz')[0]
result = align_drizzled_image(root=root, mag_limits=mag_limits,
radec=rd_ref,
NITER=5, clip=20)
orig_wcs, drz_wcs, out_shift, out_rot, out_scale = result
im = pyfits.open(drz_file)
files = []
for i in range(im[0].header['NDRIZIM']):
files.append(im[0].header['D{0:03d}DATA'.format(i+1)].split('[')[0])
for file in files:
updatehdr.updatewcs_with_shift(file, drz_files[0],
xsh=out_shift[0], ysh=out_shift[1],
rot=out_rot, scale=out_scale,
wcsname=ref_catalog, force=True,
reusename=True, verbose=True,
sciext='SCI')
im = pyfits.open(file, mode='update')
im[0].header['MJD-OBS'] = im[0].header['EXPSTART']
im.flush()
### Second drizzle
if len(files) > 1:
AstroDrizzle(files, output=root, clean=True, context=False, preserve=False, skysub=True, driz_separate=False, driz_sep_wcs=False, median=False, blot=False, driz_cr=False, driz_cr_corr=False, driz_combine=True, final_bits=576, coeffs=True, resetbits=0, build=False, final_wht_type='IVM')
else:
AstroDrizzle(files, output=root, clean=True, final_scale=None, final_pixfrac=1, context=False, final_bits=576, preserve=False, driz_separate=False, driz_sep_wcs=False, median=False, blot=False, driz_cr=False, driz_cr_corr=False, driz_combine=True, build=False, final_wht_type='IVM')
### Contamination-cleaned
adrizzle.do_driz(beam_data, beam_wcs, wht, output_wcs,
outsci, outwht, outctx, 1., 'cps', 1,
wcslin_pscale=beam.grism.wcs.pscale, uniqid=1,
pixfrac=pixfrac, kernel=kernel, fillval=0,
stepsize=10, wcsmap=None)
### Continuum
adrizzle.do_driz(beam_continuum, beam_wcs, wht, output_wcs,
coutsci, coutwht, coutctx, 1., 'cps', 1,
wcslin_pscale=beam.grism.wcs.pscale, uniqid=1,
pixfrac=pixfrac, kernel=kernel, fillval=0,
stepsize=10, wcsmap=None)
### Contamination
adrizzle.do_driz(beam.contam, beam_wcs, wht, output_wcs,
xoutsci, xoutwht, xoutctx, 1., 'cps', 1,
wcslin_pscale=beam.grism.wcs.pscale, uniqid=1,
pixfrac=pixfrac, kernel=kernel, fillval=0,
stepsize=10, wcsmap=None)
### Direct thumbnail
if direct_extension == 'REF':
if beam.direct['REF'] is None:
direct_extension = 'SCI'
if direct_extension == 'REF':
thumb = beam.direct['REF']
thumb_wht = np.cast[np.float32]((thumb != 0)*1)
else:
thumb = beam.direct[direct_extension]#/beam.direct.photflam
thumb_wht = 1./(beam.direct.data['ERR']/beam.direct.photflam)**2
wht = 1./im['ERR'].data**2
wht[(im['ERR'].data == 0) | (dq > 0) | (~np.isfinite(wht))] = 0
wht[im['SCI'].data < -3*im['ERR'].data] = 0
wht = np.cast[np.float32](wht)
exp_wcs = pywcs.WCS(im[1].header, relax=True)
exp_wcs.pscale = utils.get_wcs_pscale(exp_wcs)
#pf = 0.5
# import drizzlepac.wcs_functions as dwcs
# xx = out_wcs.deepcopy()
# #xx.all_pix2world = xx.wcs_world2pix
# map = dwcs.WCSMap(exp_wcs, xx)
astrodrizzle.adrizzle.do_driz(im['SCI'].data, exp_wcs, wht,
inter_wcs, outsci, outwht, outctx,
1., 'cps', 1,
wcslin_pscale=exp_wcs.pscale,
uniqid=1,
pixfrac=pixfrac, kernel=kernel,
fillval=0, stepsize=10,
wcsmap=SIP_WCSMap)
if ds9 is not None:
ds9.view(outsci, header=outh)
#outsci /= out_wcs.pscale**2
rms = 1/np.sqrt(outwht)
mask = (outwht == 0) | (rms > 100)
rms[mask] = 0
outsci[mask] = 0.
line_wcs = pywcs.WCS(h, relax=True)
line_wcs.pscale = utils.get_wcs_pscale(line_wcs)
# Science and wht arrays
sci = flt.grism['SCI'] - flt.model
wht = 1/(flt.grism['ERR']**2)
scl = np.exp(-(fcontam*np.abs(flt.model)/flt.grism['ERR']))
wht *= scl
wht[~np.isfinite(wht)] = 0
# Drizzle it
if verbose:
print('Drizzle {0} to wavelength {1:.2f}'.format(flt.grism.parent_file, wave))
adrizzle.do_driz(sci, line_wcs, wht, out_wcs,
outsci, outwht, outctx, 1., 'cps', 1,
wcslin_pscale=line_wcs.pscale, uniqid=1,
pixfrac=pixfrac, kernel=kernel, fillval=0,
stepsize=10, wcsmap=None)
# Done!
return outsci, outwht
wht[~np.isfinite(data_i+scl)] = 0
contam_weight[~np.isfinite(data_i+scl)] = 0
data_i[~np.isfinite(data_i+scl)] = 0
###### Go drizzle
### Contamination-cleaned
adrizzle.do_driz(data_i, beam_wcs, wht, output_wcs,
outsci, outwht, outctx, 1., 'cps', 1,
wcslin_pscale=1.0, uniqid=1,
pixfrac=pixfrac, kernel=kernel, fillval=0,
stepsize=10, wcsmap=None)
# For variance
adrizzle.do_driz(contam_weight, beam_wcs, wht, output_wcs,
outvar, outwv, outcv, 1., 'cps', 1,
wcslin_pscale=1.0, uniqid=1,
pixfrac=pixfrac, kernel=kernel, fillval=0,
stepsize=10, wcsmap=None)
if ds9 is not None:
ds9.view(outsci/output_wcs.pscale**2, header=out_header)
### Correct for drizzle scaling
outsci /= output_wcs.pscale**2
# variance
outvar /= outwv*output_wcs.pscale**2
outwht = 1/outvar
outwht[(outvar == 0) | (~np.isfinite(outwht))] = 0
else:
wht *= sens**2
beam_data /= sens
beam_continuum /= sens
###### Go drizzle
### Contamination-cleaned
adrizzle.do_driz(beam_data, beam_wcs, wht, output_wcs,
outsci, outwht, outctx, 1., 'cps', 1,
wcslin_pscale=beam.grism.wcs.pscale, uniqid=1,
pixfrac=pixfrac, kernel=kernel, fillval=0,
stepsize=10, wcsmap=None)
### Continuum
adrizzle.do_driz(beam_continuum, beam_wcs, wht, output_wcs,
coutsci, coutwht, coutctx, 1., 'cps', 1,
wcslin_pscale=beam.grism.wcs.pscale, uniqid=1,
pixfrac=pixfrac, kernel=kernel, fillval=0,
stepsize=10, wcsmap=None)
### Contamination
adrizzle.do_driz(beam.contam, beam_wcs, wht, output_wcs,
xoutsci, xoutwht, xoutctx, 1., 'cps', 1,
wcslin_pscale=beam.grism.wcs.pscale, uniqid=1,
pixfrac=pixfrac, kernel=kernel, fillval=0,
stepsize=10, wcsmap=None)
### Direct thumbnail
if direct_extension == 'REF':
if beam.direct['REF'] is None:
direct_extension = 'SCI'
dlam = np.interp(wave, beam.beam.lam[1:], np.diff(beam.beam.lam))
# 1e-17 erg/s/cm2 #, scaling closer to e-/s
sens *= 1./dlam
if sens == 0:
continue
else:
wht *= sens**2
beam_data /= sens
beam_continuum /= sens
###### Go drizzle
### Contamination-cleaned
adrizzle.do_driz(beam_data, beam_wcs, wht, output_wcs,
outsci, outwht, outctx, 1., 'cps', 1,
wcslin_pscale=beam.grism.wcs.pscale, uniqid=1,
pixfrac=pixfrac, kernel=kernel, fillval=0,
stepsize=10, wcsmap=None)
### Continuum
adrizzle.do_driz(beam_continuum, beam_wcs, wht, output_wcs,
coutsci, coutwht, coutctx, 1., 'cps', 1,
wcslin_pscale=beam.grism.wcs.pscale, uniqid=1,
pixfrac=pixfrac, kernel=kernel, fillval=0,
stepsize=10, wcsmap=None)
### Contamination
adrizzle.do_driz(beam.contam, beam_wcs, wht, output_wcs,
xoutsci, xoutwht, xoutctx, 1., 'cps', 1,
wcslin_pscale=beam.grism.wcs.pscale, uniqid=1,
scl = convert_to_flambda*beam.beam.total_flux/1.e-17
scl *= 1./beam.flat_flam.reshape(beam.beam.sh_beam).sum(axis=0)
#scl = convert_to_flambda/beam.beam.sensitivity
data_i *= scl
wht *= (1/scl)**2
#contam_weight *= scl
wht[~np.isfinite(data_i+scl)] = 0
contam_weight[~np.isfinite(data_i+scl)] = 0
data_i[~np.isfinite(data_i+scl)] = 0
###### Go drizzle
### Contamination-cleaned
adrizzle.do_driz(data_i, beam_wcs, wht, output_wcs,
outsci, outwht, outctx, 1., 'cps', 1,
wcslin_pscale=1.0, uniqid=1,
pixfrac=pixfrac, kernel=kernel, fillval=0,
stepsize=10, wcsmap=None)
# For variance
adrizzle.do_driz(contam_weight, beam_wcs, wht, output_wcs,
outvar, outwv, outcv, 1., 'cps', 1,
wcslin_pscale=1.0, uniqid=1,
pixfrac=pixfrac, kernel=kernel, fillval=0,
stepsize=10, wcsmap=None)
if ds9 is not None:
ds9.view(outsci/output_wcs.pscale**2, header=out_header)
### Correct for drizzle scaling