Python units.m方法代码示例

# 需要导入模块: from astropy import units [as 别名]
# 或者: from astropy.units import m [as 别名]
def test_log_occulterResults(self):
        """
        Test that log_occulter_results returns proper dictionary with keys
        """

        atts_list = ['slew_time', 'slew_angle', 'slew_dV', 'slew_mass_used', 'scMass']
        for mod in self.allmods:
            if 'log_occulterResults' in mod.__dict__:
                with RedirectStreams(stdout=self.dev_null):
                    if 'SotoStarshade' in mod.__name__:
                        obj = mod(f_nStars=4, **copy.deepcopy(self.spec))
                    else:
                        obj = mod(**copy.deepcopy(self.spec))
                DRM = {}
                slewTimes = np.ones((5,))*u.day
                sInds = np.arange(5)
                sd = np.ones((5,))*u.rad
                dV = np.ones((5,))*u.m/u.s

                DRM = obj.log_occulterResults(DRM, slewTimes, sInds, sd, dV)

                for att in atts_list:
                    self.assertTrue(att in DRM, 'Missing key in log_occulterResults for %s' % mod.__name__) 

# 需要导入模块: from astropy import units [as 别名]
# 或者: from astropy.units import m [as 别名]
def s_a(self, lat, lon):
        """ Standard deviation of terrain heights (m) within a 110 km × 110 km
        area with a 30 s resolution (e.g. the Globe “gtopo30” data).
        The value for the mid-path may be obtained from an area roughness
        with 0.5 × 0.5 degree resolution of geographical coordinates
        using bi-linear interpolation.
        """
        if not self._s_a:
            vals = load_data(os.path.join(dataset_dir, '530/v16_gtopo_30.txt'))
            lats = load_data(os.path.join(dataset_dir, '530/v16_lat.txt'))
            lons = load_data(os.path.join(dataset_dir, '530/v16_lon.txt'))
            self._Pr6 = bilinear_2D_interpolator(lats, lons, vals)

        return self._Pr6(
            np.array([lat.ravel(), lon.ravel()]).T).reshape(lat.shape)

    ###########################################################################
    #                               Section 2.2                               #
    ########################################################################### 

# 需要导入模块: from astropy import units [as 别名]
# 或者: from astropy.units import m [as 别名]
def XPD_outage_precipitation(self, lat, lon, d, f, el, C0_I, tau=45,
                                 U0=15, XPIF=0):
        """ Implementation of 'XPD_outage_precipitation' method for recommendation
        ITU-P R.530-16. See documentation for function
        'ITUR530.XPD_outage_precipitation'
        """
        # Step 1: Determine the path attenuation, A0.01 (dB), exceeded
        # for 0.01% of the time
        A001 = self.rain_attenuation(lat, lon, d, f, el, 0.01)

        # Step 2: Determine the equivalent path attenuation, Ap
        U = U0 + 30 * np.log10(f)
        V = np.where(f < 20, 12.8 * f**0.19, 22.6)
        Ap = 10 ** ((U - C0_I + XPIF) / V)                      # Eq. 112

        # Step 3: Determine parameters m and n
        m = min(23.26 * np.log10(Ap / (0.12 * A001)), 40)      # Eq. 113
        n = (-12.7 + np.sqrt(161.23 - 4 * m)) / 2              # Eq. 114

        # Step 4 : Determine the outage probability
        P_XPR = 10**(n - 2)                                     # Eq. 115 [%]
        return P_XPR 

# 需要导入模块: from astropy import units [as 别名]
# 或者: from astropy.units import m [as 别名]
def update(frame):
        x, y = np.random.uniform(-fov, fov, size=2)
        width = np.random.uniform(0.01, maxwid)
        length = np.random.uniform(width, maxlen)
        angle = np.random.uniform(0, 180)
        intens = width * length * (5e4 + 1e5 * np.random.exponential(2))

        model = toymodel.Gaussian(
            x=x * u.m,
            y=y * u.m,
            width=width * u.m,
            length=length * u.m,
            psi=angle * u.deg,
        )
        image, _, _ = model.generate_image(geom, intensity=intens, nsb_level_pe=5,)
        disp.image = image 

# 需要导入模块: from astropy import units [as 别名]
# 或者: from astropy.units import m [as 别名]
def main():
    fig, axs = plt.subplots(1, 2, constrained_layout=True, figsize=(6, 3))

    model = Gaussian(0 * u.m, 0.1 * u.m, 0.3 * u.m, 0.05 * u.m, 25 * u.deg)
    cam = CameraGeometry.from_name("FlashCam")
    image, *_ = model.generate_image(cam, 2500)

    CameraDisplay(cam, ax=axs[0], image=image)
    CameraDisplay(
        cam.transform_to(EngineeringCameraFrame()), ax=axs[1], image=image,
    )

    axs[0].set_title("CameraFrame")
    axs[1].set_title("EngineeringCameraFrame")

    plt.show() 

# 需要导入模块: from astropy import units [as 别名]
# 或者: from astropy.units import m [as 别名]
def cam_to_nom():
    pix_x = np.ones(2048) * u.m
    pix_y = np.ones(2048) * u.m

    pointing_direction = SkyCoord(alt=70 * u.deg, az=180 * u.deg, frame=AltAz())
    camera_frame = CameraFrame(
        focal_length=15 * u.m, telescope_pointing=pointing_direction
    )
    camera_coord = SkyCoord(pix_x, pix_y, frame=camera_frame)

    # In this case we bypass the telescope system
    nominal_frame = NominalFrame(origin=AltAz(alt=75 * u.deg, az=180 * u.deg))
    nom_coord = camera_coord.transform_to(nominal_frame)

    horizon = camera_coord.transform_to(AltAz())

    print("Nominal Coordinate", nom_coord)
    print("Horizon coordinate", horizon)


# Once we are at the nominal system where most reconstruction will be done we
# can then convert to AltAz (currently we cannot transform directly from camera) 

# 需要导入模块: from astropy import units [as 别名]
# 或者: from astropy.units import m [as 别名]
def test_construct():
    """ Check we can make a CameraGeometry from scratch """
    x = np.linspace(-10, 10, 100)
    y = np.linspace(-10, 10, 100)
    geom = CameraGeometry(
        camera_name="Unknown",
        pix_id=np.arange(100),
        pix_x=x * u.m,
        pix_y=y * u.m,
        pix_area=x * u.m ** 2,
        pix_type="rectangular",
        pix_rotation="10d",
        cam_rotation="12d",
    )

    assert geom.camera_name == "Unknown"
    assert geom.pix_area is not None
    assert (geom.pix_rotation.deg - 10) < 1e-5
    assert (geom.cam_rotation.deg - 10) < 1e-5 

# 需要导入模块: from astropy import units [as 别名]
# 或者: from astropy.units import m [as 别名]
def test_rectangle_patch_neighbors():
    """" test that a simple rectangular camera has the expected neighbors """
    pix_x = np.array([-1.1, 0.1, 0.9, -1, 0, 1, -0.9, -0.1, 1.1]) * u.m
    pix_y = np.array([1.1, 1, 0.9, -0.1, 0, 0.1, -0.9, -1, -1.1]) * u.m
    cam = CameraGeometry(
        camera_name="testcam",
        pix_id=np.arange(pix_x.size),
        pix_x=pix_x,
        pix_y=pix_y,
        pix_area=None,
        pix_type="rectangular",
    )

    assert np.all(cam.neighbor_matrix.T == cam.neighbor_matrix)
    assert cam.neighbor_matrix.sum(axis=0).max() == 4
    assert cam.neighbor_matrix.sum(axis=0).min() == 2 

# 需要导入模块: from astropy import units [as 别名]
# 或者: from astropy.units import m [as 别名]
def test_pixel_width():
    geom = CameraGeometry(
        "test",
        pix_id=[1],
        pix_area=[2] * u.cm ** 2,
        pix_x=[0] * u.m,
        pix_y=[0] * u.m,
        pix_type="hex",
    )

    assert np.isclose(geom.pixel_width.to_value(u.cm), [2 * np.sqrt(1 / np.sqrt(3))])

    geom = CameraGeometry(
        "test",
        pix_id=[1],
        pix_area=[2] * u.cm ** 2,
        pix_x=[0] * u.m,
        pix_y=[0] * u.m,
        pix_type="rect",
    )

    assert np.isclose(geom.pixel_width.to_value(u.cm), [np.sqrt(2)]) 

# 需要导入模块: from astropy import units [as 别名]
# 或者: from astropy.units import m [as 别名]
def test_construct_optics():
    """ create an OpticsDescription and make sure it
    fails if units are missing """
    OpticsDescription(
        name="test",
        num_mirrors=1,
        num_mirror_tiles=100,
        mirror_area=u.Quantity(550, u.m ** 2),
        equivalent_focal_length=u.Quantity(10, u.m),
    )

    with pytest.raises(TypeError):
        OpticsDescription(
            name="test",
            num_mirrors=1,
            num_mirror_tiles=100,
            mirror_area=550,
            equivalent_focal_length=10,
        ) 

# 需要导入模块: from astropy import units [as 别名]
# 或者: from astropy.units import m [as 别名]
def test_roundtrip_camera_horizon():
    from ctapipe.coordinates import CameraFrame, TelescopeFrame

    telescope_pointing = SkyCoord(alt=70 * u.deg, az=0 * u.deg, frame=AltAz())
    camera_frame = CameraFrame(
        focal_length=28 * u.m, telescope_pointing=telescope_pointing
    )

    cam_coord = SkyCoord(x=0.5 * u.m, y=0.1 * u.m, frame=camera_frame)
    telescope_coord = cam_coord.transform_to(TelescopeFrame())
    horizon_coord = telescope_coord.transform_to(AltAz())

    back_telescope_coord = horizon_coord.transform_to(TelescopeFrame())
    back_cam_coord = back_telescope_coord.transform_to(camera_frame)

    fov_lon = back_telescope_coord.fov_lon.to_value(u.deg)
    fov_lat = back_telescope_coord.fov_lat.to_value(u.deg)
    assert fov_lon == approx(telescope_coord.fov_lon.to_value(u.deg))
    assert fov_lat == approx(telescope_coord.fov_lat.to_value(u.deg))

    assert back_cam_coord.x.to_value(u.m) == approx(cam_coord.x.to_value(u.m))
    assert back_cam_coord.y.to_value(u.m) == approx(cam_coord.y.to_value(u.m)) 

# 需要导入模块: from astropy import units [as 别名]
# 或者: from astropy.units import m [as 别名]
def test_cam_to_hor():
    from ctapipe.coordinates import CameraFrame

    # Coordinates in any frame can be given as a numpy array of the xyz positions
    # e.g. in this case the position on pixels in the camera
    pix_x = [1] * u.m
    pix_y = [1] * u.m

    focal_length = 15000 * u.mm

    # first define the camera frame
    pointing = SkyCoord(alt=70 * u.deg, az=0 * u.deg, frame=AltAz())
    camera_frame = CameraFrame(focal_length=focal_length, telescope_pointing=pointing)

    # transform
    camera_coord = SkyCoord(pix_x, pix_y, frame=camera_frame)
    altaz_coord = camera_coord.transform_to(AltAz())

    # transform back
    altaz_coord2 = SkyCoord(az=altaz_coord.az, alt=altaz_coord.alt, frame=AltAz())
    camera_coord2 = altaz_coord2.transform_to(camera_frame)

    # check transform
    assert np.isclose(camera_coord.x.to_value(u.m), camera_coord2.y.to_value(u.m)) 

# 需要导入模块: from astropy import units [as 别名]
# 或者: from astropy.units import m [as 别名]
def test_all_to_value():
    """test all_to_value"""
    x_m = np.arange(5) * u.m
    y_mm = np.arange(5) * 1000 * u.mm
    z_km = np.arange(5) * 1e-3 * u.km
    nono_deg = np.arange(5) * 1000 * u.deg

    # one argument
    x = all_to_value(x_m, unit=u.m)
    assert (x == np.arange(5)).all()

    # two arguments
    x, y = all_to_value(x_m, y_mm, unit=u.m)
    assert (x == np.arange(5)).all()
    assert (y == np.arange(5)).all()

    # three
    x, y, z = all_to_value(x_m, y_mm, z_km, unit=u.m)
    assert (x == np.arange(5)).all()
    assert (y == np.arange(5)).all()
    assert (z == np.arange(5)).all()

    # cannot be converted
    with pytest.raises(u.UnitConversionError):
        all_to_value(x_m, nono_deg, unit=x_m.unit) 

# 需要导入模块: from astropy import units [as 别名]
# 或者: from astropy.units import m [as 别名]
def make_interpolant(self, param, value, unit):
        """Make an interpolating function."""

        # the generic quadratic function we have used to test
        # always positive, reaches a maximum of 1.0 at x=0.5
        quadratic = lambda x: 0.25 + 3.0 * (1-x) * x

        if isinstance(value, (numbers.Number, np.ndarray)):
            if param == 'QE':
                return lambda lam: value*unit
            elif param in ('core_thruput', 'core_contrast', 'PSF'):
                # for PSF, will be an ndarray
                return lambda lam, WA: value*unit
        elif isinstance(value, basestring):
            # for most ty
            if param == 'QE':
                return lambda lam: quadratic(lam)*unit
            elif param in ('core_thruput', 'core_contrast'):
                return lambda lam, WA: quadratic(WA)*unit
            elif param == 'PSF':
                # this rather messy construct uses a value like
                # "psf_5x5.fits" to recover a PSF matrix to use, else,
                # it uses a fixed value.  The pattern matches
                # psf_NNNxMMM.fits where NNN and MMM are digit sequences.
                m = re.search("psf_(\d+)x(\d+)\.fits", value)
                if m is None:
                    # use a fixed value, we won't need it anyway
                    a_value = np.array([1])
                else:
                    # this is the size, like [5,5]
                    s = [int(n) for n in m.groups()]
                    # this is the value, which is always a progression
                    # from 0...prod(s)-1, reshaped to be the size asked for
                    a_value = np.arange(np.prod(s)).reshape(s)
                return lambda lam, WA: a_value*unit
        else:
            assert False, "unknown interpolant needed" 

# 需要导入模块: from astropy import units [as 别名]
# 或者: from astropy.units import m [as 别名]
def find_known_plans(self):
        """
        Find and return list of known RV stars and list of stars with earthlike planets
        """
        TL = self.TargetList
        SU = self.SimulatedUniverse

        c = 28.4 *u.m/u.s
        Mj = 317.8 * u.earthMass
        Mpj = SU.Mp/Mj                     # planet masses in jupiter mass units
        Ms = TL.MsTrue[SU.plan2star]
        Teff = TL.stellarTeff(SU.plan2star)
        mu = const.G*(SU.Mp + Ms)
        T = (2.*np.pi*np.sqrt(SU.a**3/mu)).to(u.yr)
        e = SU.e

        t_filt = np.where((Teff.value > 3000) & (Teff.value < 6800))[0]    # planets in correct temp range

        K = (c / np.sqrt(1 - e[t_filt])) * Mpj[t_filt] * np.sin(SU.I[t_filt]) * Ms[t_filt]**(-2/3) * T[t_filt]**(-1/3)

        K_filter = (T[t_filt].to(u.d)/10**4).value
        K_filter[np.where(K_filter < 0.03)[0]] = 0.03
        k_filt = t_filt[np.where(K.value > K_filter)[0]]               # planets in the correct K range

        a_filt = k_filt[np.where((SU.a[k_filt] > .95*u.AU) & (SU.a[k_filt] < 1.67*u.AU))[0]]   # planets in habitable zone
        r_filt = a_filt[np.where(SU.Rp.value[a_filt] < 1.75)[0]]                               # rocky planets
        self.earth_candidates = np.union1d(self.earth_candidates, r_filt).astype(int)

        known_stars = np.unique(SU.plan2star[k_filt])
        known_rocky = np.unique(SU.plan2star[r_filt])
        return known_stars.astype(int), known_rocky.astype(int)