Python fits.getheader函数代码示例

    def _get_fits_header(self, result, options):
        """Populate FITS Header keywords"""
        super(WFIRSTInstrument, self)._get_fits_header(result, options)
        pupil_hdr = fits.getheader(self.pupil)
        apodizer_hdr = fits.getheader(self._apodizer_fname)
        fpm_hdr = fits.getheader(self._fpm_fname)
        lyotstop_hdr = fits.getheader(self._lyotstop_fname)

        result[0].header.set('MODE', self.mode, comment='Observing mode')
        result[0].header.set('CAMERA', self.camera, comment='Imager or IFS')
        result[0].header.set('APODIZER', self.apodizer, comment='Apodizer')
        result[0].header.set('APODTRAN', os.path.basename(self._apodizer_fname),
                             comment='Apodizer transmission')
        result[0].header.set('PUPLSCAL', apodizer_hdr['PUPLSCAL'],
                             comment='Apodizer pixel scale in m/pixel')
        result[0].header.set('PUPLDIAM', apodizer_hdr['PUPLDIAM'],
                             comment='Full apodizer array size, incl padding.')
        result[0].header.set('FPM', self.fpm, comment='Focal plane mask')
        result[0].header.set('FPMTRAN', os.path.basename(self._fpm_fname),
                             comment='FPM transmission')
        result[0].header.set('FPMSCAL', fpm_hdr['PIXSCALE'], comment='FPM spatial sampling, arcsec/pix')
        result[0].header.set('LYOTSTOP', self.lyotstop, comment='Lyot stop')
        result[0].header.set('LSTRAN', os.path.basename(self._lyotstop_fname),
                             comment='Lyot stop transmission')
        result[0].header.set('PUPLSCAL', lyotstop_hdr['PUPLSCAL'],
                             comment='Lyot stop pixel scale in m/pixel')
        result[0].header.set('PUPLDIAM', lyotstop_hdr['PUPLDIAM'],
                             comment='Lyot stop array size, incl padding.')

def SExtractorDualImage(drPath,tract,bands,patches,prefix,sexdir,dotsex,zps):
    '''
    Run SExtractor in dual image mode for each patch,band using chi2 image for detection
    '''
    os.chdir(sexdir)
    for band in bands:
        for patch in patches:
            # Move Variance extension to new file for SExtractor
            imname = drPath+'deepCoadd/'+band+'/'+tract+'/'+patch+'/'+\
                        'calexp-'+band+'-'+tract+'-'+patch[0]+patch[-1]
            varHead = fits.getheader(imname+'.fits','VARIANCE')
            maskHead = fits.getheader(imname+'.fits','MASK')
            im = fits.open(imname+'.fits')
            fits.writeto(imname+'_var.fits',im['VARIANCE'].data,varHead,clobber=True)
            fits.writeto(imname+'_mask.fits',im['MASK'].data.astype(float),maskHead,clobber=True)
            inImage = imname+'.fits[1]'
            outCat = drPath+'deepCoadd/'+band+'/'+tract+'/'+patch+'/'+band+'-'+tract+'-'+patch[0]+patch[-1]+'-chi2.cat'
            chi2Image = drPath+'deepCoadd/'+prefix+'chi2/'+tract+'/'+patch+'/'+prefix+'chi2-'+tract+'-'+patch+'.fits'
            # Run SExtractor
            os.system('sex '+chi2Image+','+inImage+' -c '+dotsex+' -CATALOG_NAME '+outCat+' -WEIGHT_IMAGE None,'+\
                        imname+'_var.fits -WEIGHT_TYPE NONE,MAP_VAR -MAG_ZEROPOINT '+str(zps[bands.index(band)]))
            # Add flags to catalog
            os.system('./venice -m '+imname+'_mask.fits -cat '+outCat+ ' -f all -xcol 2 -ycol 3 -o '+\
                    drPath+'deepCoadd/'+band+'/'+tract+'/'+patch+'/'+band+'-'+tract+'-'+patch[0]+patch[-1]+'-chi2-flags.cat')
            
            os.system('rm '+imname+'_var.fits')
            os.system('rm '+imname+'_mask.fits')

def load_uvot_images(galaxy, band, imdir=''):
    """ugh:
    :returns hdr:

    :returns images:
        List of images  [cps, exposure, sensitivity]

    :returns masks:
        Dictionary of mask images, keyed by strings:
        ['bkg', 'reg', 'reg1', 'reg5', 'reg20']
    """
    imroot = os.path.join(imdir, galaxy, galaxy)
    info = {'root': imroot, 'band': band}
    
    names = ['{root}_t{band}.fits', '{root}_t{band}_ex.fits', '{root}_t{band}_lss.fits']
    #print(names[0].format(**info))
    try:
        hdr = pyfits.getheader(names[0].format(**info), 1)
    except:
        hdr = pyfits.getheader(names[0].format(**info), 0)
    images = [pyfits.getdata(n.format(**info)) for n in names]

    masknames = ['bkg', 'reg', 'reg1', 'reg5', 'reg20']
    masks = [pyfits.getdata('{root}_{mask}.fits'.format(root=info['root'], mask=n))
            for n in masknames]

    return hdr, images, dict(zip(masknames, masks))

def hdr_get(gal,data_mode='12m'):
    if isinstance(gal,Galaxy):
        name = gal.name.lower()
    elif isinstance(gal,str):
        name = gal.lower()
    else:
        raise ValueError("'gal' must be a str or galaxy!")
    if data_mode == '7m':
        data_mode = '7m'
        conbeam=None
        print( 'WARNING: SFR maps come in 12m sizes only.') #(!!!) What about for all the new 15" maps?
        print( 'WARNING: Convolution forcibly disabled.')
    elif data_mode in ['12m','12m+7m']:
        data_mode = '12m+7m'  
    
    hdr = None
    hdr_found = False
    if name=='m33':
        for filename in [\
        'notphangsdata/M33_14B-088_HI.clean.image.GBT_feathered.pbcov_gt_0.5_masked.peakvels.fits']:
            if os.path.isfile(filename):
                hdr = fits.getheader(filename)
                hdr_found = True
    else:
        for filename in [\
        'phangsdata/'+name+'_co21_'+data_mode+'+tp_mom0.fits',\
        'phangsdata/'+name+'_co21_'+data_mode+'+tp_mom1.fits',\
        'phangsdata/'+name+'_co21_'+data_mode+'+tp_tpeak.fits']:
            if os.path.isfile(filename):
                hdr = fits.getheader(filename)
                hdr_found = True
    if hdr_found == False:
        print('WARNING: No header was found!')
        hdr = None
    return hdr

def test_adding_overscan_apogee_u9(make_overscan_test_files):
    original_dir = getcwd()

    apogee = ApogeeAltaU9()
    print(getcwd())
    oscan_dir, has_oscan, has_no_oscan = make_overscan_test_files
    print(getcwd())

    chdir(path.join(_test_dir, oscan_dir))
    # first, does requesting *not* adding overscan actually leave it alone?
    ph.patch_headers(dir='.', new_file_ext='', overwrite=True, purge_bad=False,
                     add_time=False, add_apparent_pos=False,
                     add_overscan=False, fix_imagetype=False)
    header_no_oscan = fits.getheader(has_no_oscan)
    assert 'biassec' not in header_no_oscan
    assert 'trimsec' not in header_no_oscan

    # Now add overscan
    ph.patch_headers(dir='.', new_file_ext='', overwrite=True, purge_bad=False,
                     add_time=False, add_apparent_pos=False,
                     add_overscan=True, fix_imagetype=False)
    print(_test_dir)
    header_no_oscan = fits.getheader(has_no_oscan)
    # This image had no overscan, so should be missing the relevant keywords.
    assert 'biassec' not in header_no_oscan
    assert 'trimsec' not in header_no_oscan
    header_yes_oscan = fits.getheader(has_oscan)
    # This one as overscan, so should include both of the overscan keywords.
    assert header_yes_oscan['biassec'] == apogee.useful_overscan
    assert header_yes_oscan['trimsec'] == apogee.trim_region
    print(getcwd())
    chdir(original_dir)
    print(getcwd())

def logfile_sci(logfile, imlist, SaveSteps=False):
    '''Add science images and info to logfile
        
       logfile  - Log file name
       imlist   - List of science files reduced 
     SaveSteps  - Set to True if each step in the reduction is saved
                  as a new file
    '''
    imlist.sort()
    nfiles = len(imlist)
    log = open(logfile, 'a')
    log.write('# Science Images\n')
    if SaveSteps:   
        # if user is saving interim reduction steps to new files, 
        # include the extensions for each caliration type to header
        log.write('#           *.fits - raw image \n')
        log.write('#        *.bs.fits - bias-subtracted image \n')
        log.write('#     *.bs.ff.fits - bias-subtracted and flat-fielded \n')
    log.write('# Image  Object  Exptime  Filter\n')
    # for each file, print the image name, object type, exptime, and filter
    for i in range(nfiles):
        imname = os.path.splitext(os.path.basename(imlist[i]))[0]
        obj = fits.getheader(imlist[i])['OBJECT']
        exptime = fits.getheader(imlist[i])['EXPTIME']
        filt = fits.getheader(imlist[i])['FILTER']
        log.write('%s    %s   %.3f   %s \n' % (imname, obj, exptime, filt))
    log.write('\n')
    log.write('\n')
    log.write('# Comments: \n')
    log.write('\n')
    log.close()

def find_best_flat(flt_fits, verbose=True):
    """
    Find the most recent PFL file in $IREF for the filter used for the 
    provided FLT image.  Doesn't do any special check on USEAFTER date, just
    looks for the most-recently modified file. 
    """
    import glob
    import os.path
    import time
    
    IREF = os.environ["iref"]

    the_filter = fits.getheader(flt_fits,0).get('FILTER')
    
    pfls = glob.glob(IREF+'*pfl.fits')
    latest = 0
    best_pfl = None
    
    for pfl in pfls:
        head = fits.getheader(pfl)
        if head.get('FILTER') != the_filter:
            continue    
        
        this_created = os.path.getmtime(pfl)
        if this_created > latest:
            best_pfl = pfl
            latest = this_created
            
        if verbose:
            print '%s %s %s' %(pfl, the_filter, time.ctime(latest))
    
    return best_pfl #, the_filter, time.ctime(latest)

    def test_GetHeaderConvienceFunction(self, filename, ext, naxis1, naxis2):
        """Test the getheader convience function in both the fits and
        stpyfits namespace."""

        if ext is None:
            hd = stpyfits.getheader(self.data(filename))
            hd1 = fits.getheader(self.data(filename))
        else:
            hd = stpyfits.getheader(self.data(filename), ext)
            hd1 = fits.getheader(self.data(filename), ext)

        assert hd['NAXIS'] == 2
        assert hd1['NAXIS'] == 0
        assert hd['NAXIS1'] == naxis1
        assert hd['NAXIS2'] == naxis2

        for k in ('NAXIS1', 'NAXIS2'):
            with pytest.raises(KeyError):
                hd1[k]
        for k in ('NPIX1', 'NPIX2'):
            with pytest.raises(KeyError):
                hd[k]

        assert hd1['NPIX1'] == naxis1
        assert hd1['NPIX2'] == naxis2

def crossmatchtwofiles(img1, img2, radius=3):
    ''' This module is crossmatch two images:
        It run sextractor transform the pixels position of the the sources in coordinates and crossmatch them  
        The output is a dictionary with the objects in common
    '''

    from numpy import array, argmin, min, sqrt
    import agnkey

    hd1 = pyfits.getheader(img1)
    hd2 = pyfits.getheader(img2)
    wcs1 = pywcs.WCS(hd1)
    wcs2 = pywcs.WCS(hd2)

    xpix1, ypix1, fw1, cl1, cm1, ell1, bkg1, fl1 = agnkey.agnastrodef.sextractor(img1)
    xpix2, ypix2, fw2, cl2, cm2, ell2, bkg2, fl2 = agnkey.agnastrodef.sextractor(img2)
    xpix1, ypix1, xpix2, ypix2 = array(xpix1, float), array(ypix1, float), array(xpix2, float), array(ypix2, float)

    bb = wcs1.wcs_pix2sky(zip(xpix1, ypix1), 1)  # transform pixel in coordinate
    xra1, xdec1 = zip(*bb)
    bb = wcs2.wcs_pix2sky(zip(xpix2, ypix2), 1)  # transform pixel in coordinate
    xra2, xdec2 = zip(*bb)

    xra1, xdec1, xra2, xdec2 = array(xra1, float), array(xdec1, float), array(xra2, float), array(xdec2, float)
    distvec, pos1, pos2 = agnkey.agnastrodef.crossmatch(xra1, xdec1, xra2, xdec2, radius)
    # dict={}
    dict = {'ra1': xra1[pos1], 'dec1': xdec1[pos1], 'ra2': xra2[pos2], 'dec2': xdec2[pos2],
            'xpix1': xpix1[pos1], 'ypix1': ypix1[pos1], 'xpix2': xpix2[pos2], 'ypix2': ypix2[pos2]}

    out = next(tempfile._get_candidate_names())+'.list'
    np.savetxt(out, zip(xpix1[pos1], ypix1[pos1]), fmt='%10.10s\t%10.10s')
    return out, dict

def get_header_info(pipeline, cat_name, frame):
    """
    Get the exposure time, airmass, and gain from a DECam fits header
    
    Parameters
    ----------
    pipeline: :class:`.Pipeline`
        Object that contains parameters about a decam pipeline
    cat_name: str
        Name of a catalog that contains a list of sources found using decamtoyz
    frame: int
        CCD number of the header to load
    
    Returns
    -------
        exptime: float
            Exposure time of the image
        airmass: float
            Airmass (sec(zenith distance))
        gain: float
            Gain in electrons/adu
    """
    from astropyp import index
    from astropy.io import fits
    import numpy as np
    
    expnum = os.path.basename(cat_name).split('-')[0]
    sql = "select * from decam_files where EXPNUM={0} and PRODTYPE='image'".format(int(expnum))
    files = index.query(sql, pipeline.idx_connect_str)
    header = fits.getheader(os.path.join(pipeline.paths['decam'], files.iloc[0].filename), ext=0)
    exptime = header['exptime']
    airmass = 1/np.cos(np.radians(header['ZD']))
    header = fits.getheader(os.path.join(pipeline.paths['decam'], files.iloc[0].filename), ext=frame)
    gain = header['arawgain']
    return exptime, airmass, gain

def image_limits(image, wispfield):
    '''Find the '''
    # read in combined image
    hdr = fits.getheader(image)
    hdr_wcs = pywcs.WCS(hdr)
    # list of images that were combined
    ims = hdr['IMCMB*']
    Nim = len(ims)
    xx = np.zeros((Nim,2), dtype=int)
    yy = np.zeros((Nim,2), dtype=int)
    for i in range(Nim):
        if os.path.splitext(ims[i])[1] != '.fits':
            ims[i] = ims[i] + '.fits'
        hd = fits.getheader(os.path.join(wispfield,ims[i]))
        nx1 = hd['NAXIS1']
        nx2 = hd['NAXIS2']
        hd_wcs = pywcs.WCS(hd)
        # find WCS coordinates of input image corners
        coo = hd_wcs.wcs_pix2sky([[1,1],[1,nx1-1],[nx2-1,1],[nx2-1,nx1-1]],1)
        # find corresponding pixel coords in combined image
        pix = hdr_wcs.wcs_sky2pix(coo,1)
        xx[i,0] = pix[0][1]
        xx[i,1] = pix[1][1]
        yy[i,0] = pix[0][0]
        yy[i,1] = pix[2][0]

    xmin,xmax = np.max(xx[:,0]), np.min(xx[:,1])
    ymin,ymax = np.max(yy[:,0]), np.min(yy[:,1])
    return xmin,xmax,ymin,ymax

    def testGetHeaderConvienceFunction(self):
        """Test the getheader convience function in both the fits and
           stpyfits namespace."""

        hd = stpyfits.getheader(self.data('cdva2.fits'))
        hd1 = fits.getheader(self.data('cdva2.fits'))

        assert_equal(hd['NAXIS'], 2)
        assert_equal(hd1['NAXIS'], 0)
        assert_equal(hd['NAXIS1'], 10)
        assert_equal(hd['NAXIS2'], 10)

        assert_raises(KeyError, lambda: hd1['NAXIS1'])
        assert_raises(KeyError, lambda: hd1['NAXIS2'])
        assert_raises(KeyError, lambda: hd['NPIX1'])
        assert_raises(KeyError, lambda: hd['NPIX2'])

        assert_equal(hd1['NPIX1'], 10)
        assert_equal(hd1['NPIX2'], 10)

        hd = stpyfits.getheader(self.data('o4sp040b0_raw.fits'), 2)
        hd1 = fits.getheader(self.data('o4sp040b0_raw.fits'), 2)

        assert_equal(hd['NAXIS'], 2)
        assert_equal(hd1['NAXIS'], 0)
        assert_equal(hd['NAXIS1'], 62)
        assert_equal(hd['NAXIS2'], 44)

        assert_raises(KeyError, lambda: hd1['NAXIS1'])
        assert_raises(KeyError, lambda: hd1['NAXIS2'])
        assert_raises(KeyError, lambda: hd['NPIX1'])
        assert_raises(KeyError, lambda: hd['NPIX2'])

        assert_equal(hd1['NPIX1'], 62)
        assert_equal(hd1['NPIX2'], 44)

def filter_combo(Palim, UVISims):
    '''Return the filter combination that should be used for 
       determination of the color term.
    '''
    # get Palomar SDSS filters
    hdr = fits.getheader(Palim[0])
    Pal_filt = hdr['FILTER'].rstrip("'")
    
    # get UVIS filters
    UVIS_filt = []
    for im in UVISims:
        hdr = fits.getheader(im)
        UVIS_filt.append(hdr['FILTER'])

    # are both UVIS filters available?
    if len(UVIS_filt) > 1:
        print '\nNOTE: Both UVIS filters -- %s -- are available for '+\
            'this field\n' % [x for x in UVIS_filt]
   
    # possible filter combinations for photometric calibration
    if Pal_filt == 'g':
        use_filt = [x for x in UVIS_filt if x == 'F600LP' or x == 'F814W']
        use_im = [x for x in UVISims if \
                  os.path.basename(x).split('_')[0] == 'F600LP']
    if Pal_filt == 'i':
        use_filt = [x for x in UVIS_filt if x == 'F475X' or x == 'F606W']
        use_im = [x for x in UVISims if \
                  os.path.basename(x).split('_')[0] == 'F475X']
    
    use_filt.append(Pal_filt)
    use_filt.sort()
    return use_filt, use_im

def crossmatchtwofiles(img1, img2, radius=3):
    ''' This module is crossmatch two images:
        It run sextractor transform the pixels position of the the sources in coordinates and crossmatch them  
        The output is a dictionary with the objects in common
    '''
    import lsc
    from astropy.wcs import WCS
    from numpy import array, argmin, min, sqrt

    hd1 = fits.getheader(img1)
    hd2 = fits.getheader(img2)
    wcs1 = WCS(hd1)
    wcs2 = WCS(hd2)

    xpix1, ypix1, fw1, cl1, cm1, ell1, bkg1, fl1 = lsc.lscastrodef.sextractor(img1)
    xpix2, ypix2, fw2, cl2, cm2, ell2, bkg2, fl2 = lsc.lscastrodef.sextractor(img2)
    xpix1, ypix1, xpix2, ypix2 = array(xpix1, float), array(ypix1, float), array(xpix2, float), array(ypix2, float)

    bb = wcs1.wcs_pix2world(zip(xpix1, ypix1), 1)  #   transform pixel in coordinate
    xra1, xdec1 = zip(*bb)
    bb = wcs2.wcs_pix2world(zip(xpix2, ypix2), 1)  #   transform pixel in coordinate
    xra2, xdec2 = zip(*bb)

    xra1, xdec1, xra2, xdec2 = array(xra1, float), array(xdec1, float), array(xra2, float), array(xdec2, float)
    distvec, pos1, pos2 = lsc.lscastrodef.crossmatch(xra1, xdec1, xra2, xdec2, radius)
    #dict={}
    dict = {'ra1': xra1[pos1], 'dec1': xdec1[pos1], 'ra2': xra2[pos2], 'dec2': xdec2[pos2], \
            'xpix1': xpix1[pos1], 'ypix1': ypix1[pos1], 'xpix2': xpix2[pos2], 'ypix2': ypix2[pos2]}
    np.savetxt('substamplist', zip(xpix1[pos1], ypix1[pos1]), fmt='%10.10s\t%10.10s')
    return 'substamplist', dict

    def setup(self, SE_only, no_conv):
        """Return lists of all images and filters used in this Par.
           We will use the unrotated images for use with a single psf
           Image filenames:
             ParXXX/DATA/UVIS/IRtoUVIS/FYYY_UVIS_sci.fits
             ParXXX/DATA/UVIS/IRtoUVIS/FYYY_UVIS_rms.fits
        """
        images = glob(os.path.join(self.imdir, 'F*_UVIS_sci.fits'))

        # dictionary of zero points
        zps = self.get_zp(images[0])

        # build table
        t = Table(data=None, 
                  names=['filt','image','convim','rms','wht','exptime','zp'],
                  dtype=['S10', 'S60', 'S60', 'S60', 'S60', float, float])
        for image in images:
            filt = fits.getheader(image)['FILTER']

            # weight map
            wht = image.split('_sci.fits')[0] + '_wht.fits'

            # clean image for convolution
            tmp = os.path.splitext(image)[0] + '_cln.fits'
            image_cln = os.path.join(self.outdir, os.path.basename(tmp))
            if SE_only is False:
                print 'Cleaning %s' % os.path.basename(image)
                if no_conv:
                    clean_image(image, image_cln, cln_by_wht=False, whtim=wht)
                else:
                    clean_image(image, image_cln, cln_by_wht=True, whtim=wht)
            
            # names of convolved images
            if filt == self.reddest_filt:
                convim = image_cln
            else:
                check = re.search('\d+', self.reddest_filt)
                rf = check.group(0)
                convim = os.path.join(self.outdir,'%s_convto%s.fits'%(filt,rf))

            # replace zeros with 1.e9 in rms analysis maps
            rms0 = image.split('_sci.fits')[0] + '_rms.fits'
            tmp = os.path.splitext(rms0)[0] + '_analysis.fits'
            rms_analysis = os.path.join(self.outdir, os.path.basename(tmp))
            self.fix_rms_map(rms0, rms_analysis, value=1.e10,rmstype='analysis')

            # for detection image, create detection RMS map as well
            if filt == self.detect_filt:
                tmp2 = os.path.splitext(rms0)[0] + '_detection.fits'
                rms_detect = os.path.join(self.outdir, os.path.basename(tmp2))
                self.fix_rms_map(rms0, rms_detect, value=0.01, 
                                 rmstype='detection', whtim=wht)
            
            exptime = fits.getheader(image)['EXPTIME']
            zp = zps[filt]

            t.add_row([filt, image_cln, convim, rms_analysis, wht, exptime, zp])
        # set detection RMS map
        self.detect_rms = rms_detect
        return t