Python astropy.units 模块，Jy() 实例源码

def sky(self):
        """
        OSKAR sky model array converted from the input image.

        Columns
        -------
        ra : (J2000) right ascension (deg)
        dec : (J2000) declination (deg)
        flux : source (Stokes I) flux density (Jy)
        """
        idx = self.mask.flatten()
        ra, dec = self.ra_dec
        ra = ra.flatten()[idx]
        dec = dec.flatten()[idx]
        flux = self.image.flatten()[idx] * self.factor_K2JyPixel
        sky_ = np.column_stack([ra, dec, flux])
        return sky_

def write_sky_model(self, outfile, clobber=False):
        """
        Write the converted sky model for simulation.
        """
        if os.path.exists(outfile) and (not clobber):
            raise OSError("OSKAR sky model file already exists: %s" % outfile)
        sky = self.sky
        counts = sky.shape[0]
        percent = 100 * counts / self.image.size
        logger.info("Source counts: %d (%.1f%%)" % (counts, percent))
        header = ("Frequency = %.3f [MHz]\n" % self.freq +
                  "Pixel size = %.2f [arcsec]\n" % self.pixelsize +
                  "K2JyPixel = %.2f\n" % self.factor_K2JyPixel +
                  "RA0 = %.4f [deg]\n" % self.ra0 +
                  "Dec0 = %.4f [deg]\n" % self.dec0 +
                  "Minimum value = %.4e [K]\n" % self.minvalue +
                  "Maximum value = %.4e [K]\n" % self.maxvalue +
                  "Source counts = %d (%.1f%%)\n\n" % (counts, percent) +
                  "R.A.[deg]    Dec.[deg]    flux[Jy]")
        logger.info("Writing sky model ...")
        np.savetxt(outfile, sky, fmt='%.10e, %.10e, %.10e', header=header)
        logger.info("Wrote OSKAR sky model to file: %s" % outfile)

def fits_header(self):
        header = self.wcs.to_header()
        header["BUNIT"] = ("Jy/pixel", "Brightness unit")
        header["FREQ"] = (self.freq, "Frequency [MHz]")
        header["RA0"] = (self.ra0, "Center R.A. [deg]")
        header["DEC0"] = (self.dec0, "Center Dec. [deg]")
        header["PixSize"] = (self.pixelsize, "Pixel size [arcsec]")
        header["K2JyPix"] = (self.factor_K2JyPixel, "[K] -> [Jy/pixel]")
        header["MINVALUE"] = (self.minvalue, "[K] minimum threshold")
        if np.isfinite(self.maxvalue):
            header["MAXVALUE"] = (self.maxvalue, "[K] maximum threshold")
        return header

def factor_K2JyPixel(self):
        """
        Conversion factor from [K] to [Jy/pixel]
        """
        pixarea = (self.pixelsize * au.arcsec) ** 2
        freq = self.freq * au.MHz
        equiv = au.brightness_temperature(pixarea, freq)
        factor = au.K.to(au.Jy, equivalencies=equiv)
        return factor

def SubFluxOperation(y):
    print(galaxies[y])

    collection = ['F475W','F814W','F160W']

    flux = np.zeros([len(collection),len(radii)]) #*u.Jy
    subflux = np.zeros([len(collection),len(radii)])

    for i in range (0, len(collection)):

        # read in the images
        file = glob.glob(dir+galaxies[y]+'_final_'+collection[i]+'*sci.fits')
        hdu = fits.open(file[0])
        data[i], header[i] = hdu[0].data, hdu[0].header
        fnu[i] = header[i]['PHOTFNU']
        exp[i] = header[i]['EXPTIME']

        #define positions for photometry
        positions = [(xcen[y], ycen[y])]

        #do photometry on images
        #convert to proper units
        for j in range(0,len(radii)):
            aperture = CircularAperture(positions, radii[j])
            phot_table = aperture_photometry(data[i], aperture)
            #the next line changes the data from its original e-/s to something else that i don't remember.... because i suck.
            flux[i,j] = phot_table['aperture_sum'][0]*(fnu[i]/exp[i])
            if j == 0:
                subflux[i,j] = flux[i,j]
            else:
                subflux[i,j] = flux[i,j]-flux[i,j-1]
    return(subflux)
    # now let's put stuff in the data structure we created

def test_from_fnu(self):
        fluxd = 6.161081515869728 * u.mJy
        nu = 2.998e14 / 11.7 * u.Hz
        aper = 1 * u.arcsec
        eph = dict(rh=1.5 * u.au, delta=1.0 * u.au)
        S = 1.711e14 * u.Jy
        afrho = Afrho.from_fluxd(nu, fluxd, aper, eph, S=S)
        assert np.isclose(afrho.cm, 1000.0)

def test_from_fnu(self):
        fluxd = 6.0961896974549115 * u.mJy
        nu = 2.998e14 / 11.7 * u.Hz
        aper = 1 * u.arcsec
        eph = dict(rh=1.5 * u.au, delta=1.0 * u.au)
        S = 1.711e14 * u.Jy
        efrho = Efrho.from_fluxd(nu, fluxd, aper, eph)
        assert np.isclose(efrho.cm, 33.0)