Python numpy.radians() 使用实例

The following are code examples for showing how to use . They are extracted from open source Python projects. You can vote up the examples you like or vote down the exmaples you don’t like. You can also save this page to your account.

Example 1

def angle_wrap(angle,radians=False):
    '''
    Wraps the input angle to 360.0 degrees.

    if radians is True: input is assumed to be in radians, output is also in
    radians

    '''

    if radians:
        wrapped = angle % (2.0*PI)
        if wrapped < 0.0:
            wrapped = 2.0*PI + wrapped

    else:

        wrapped = angle % 360.0
        if wrapped < 0.0:
            wrapped = 360.0 + wrapped

    return wrapped 

Example 2

def great_circle_dist(p1, p2):
    """Return the distance (in km) between two points in
    geographical coordinates.
    """
    lon0, lat0 = p1
    lon1, lat1 = p2
    EARTH_R = 6372.8
    lat0 = np.radians(float(lat0))
    lon0 = np.radians(float(lon0))
    lat1 = np.radians(float(lat1))
    lon1 = np.radians(float(lon1))
    dlon = lon0 - lon1
    y = np.sqrt(
        (np.cos(lat1) * np.sin(dlon)) ** 2
        + (np.cos(lat0) * np.sin(lat1)
           - np.sin(lat0) * np.cos(lat1) * np.cos(dlon)) ** 2)
    x = np.sin(lat0) * np.sin(lat1) + \
        np.cos(lat0) * np.cos(lat1) * np.cos(dlon)
    c = np.arctan2(y, x)
    return EARTH_R * c 

Example 3

def get_data(filename,headers,ph_units):
    # Importation des données .DAT
    dat_file = np.loadtxt("%s"%(filename),skiprows=headers,delimiter=',')
    labels = ["freq", "amp", "pha", "amp_err", "pha_err"]
    data = {l:dat_file[:,i] for (i,l) in enumerate(labels)}
    if ph_units == "mrad":
        data["pha"] = data["pha"]/1000                    # mrad to rad
        data["pha_err"] = data["pha_err"]/1000              # mrad to rad
    if ph_units == "deg":
        data["pha"] = np.radians(data["pha"])               # deg to rad
        data["pha_err"] = np.radians(data["pha_err"])       # deg to rad
    data["phase_range"] = abs(max(data["pha"])-min(data["pha"])) # Range of phase measurements (used in NRMS error calculation)
    data["Z"]  = data["amp"]*(np.cos(data["pha"]) + 1j*np.sin(data["pha"]))
    EI = np.sqrt(((data["amp"]*np.cos(data["pha"])*data["pha_err"])**2)+(np.sin(data["pha"])*data["amp_err"])**2)
    ER = np.sqrt(((data["amp"]*np.sin(data["pha"])*data["pha_err"])**2)+(np.cos(data["pha"])*data["amp_err"])**2)
    data["Z_err"] = ER + 1j*EI
    # Normalization of amplitude
    data["Z_max"] = max(abs(data["Z"]))  # Maximum amplitude
    zn, zn_e = data["Z"]/data["Z_max"], data["Z_err"]/data["Z_max"] # Normalization of impedance by max amplitude
    data["zn"] = np.array([zn.real, zn.imag]) # 2D array with first column = real values, second column = imag values
    data["zn_err"] = np.array([zn_e.real, zn_e.imag]) # 2D array with first column = real values, second column = imag values
    return data 

Example 4

def get_data(filename,headers,ph_units):
    # Importation des données .DAT
    dat_file = np.loadtxt("%s"%(filename),skiprows=headers,delimiter=',')
    labels = ["freq", "amp", "pha", "amp_err", "pha_err"]
    data = {l:dat_file[:,i] for (i,l) in enumerate(labels)}
    if ph_units == "mrad":
        data["pha"] = data["pha"]/1000                    # mrad to rad
        data["pha_err"] = data["pha_err"]/1000              # mrad to rad
    if ph_units == "deg":
        data["pha"] = np.radians(data["pha"])               # deg to rad
        data["pha_err"] = np.radians(data["pha_err"])       # deg to rad
    data["phase_range"] = abs(max(data["pha"])-min(data["pha"])) # Range of phase measurements (used in NRMS error calculation)
    data["Z"]  = data["amp"]*(np.cos(data["pha"]) + 1j*np.sin(data["pha"]))
    EI = np.sqrt(((data["amp"]*np.cos(data["pha"])*data["pha_err"])**2)+(np.sin(data["pha"])*data["amp_err"])**2)
    ER = np.sqrt(((data["amp"]*np.sin(data["pha"])*data["pha_err"])**2)+(np.cos(data["pha"])*data["amp_err"])**2)
    data["Z_err"] = ER + 1j*EI
    # Normalization of amplitude
    data["Z_max"] = max(abs(data["Z"]))  # Maximum amplitude
    zn, zn_e = data["Z"]/data["Z_max"], data["Z_err"]/data["Z_max"] # Normalization of impedance by max amplitude
    data["zn"] = np.array([zn.real, zn.imag]) # 2D array with first column = real values, second column = imag values
    data["zn_err"] = np.array([zn_e.real, zn_e.imag]) # 2D array with first column = real values, second column = imag values
    return data 

Example 5

def plotArc(start_angle, stop_angle, radius, width, **kwargs):
    """ write a docstring for this function"""
    numsegments = 100
    theta = np.radians(np.linspace(start_angle+90, stop_angle+90, numsegments))
    centerx = 0
    centery = 0
    x1 = -np.cos(theta) * (radius)
    y1 = np.sin(theta) * (radius)
    stack1 = np.column_stack([x1, y1])
    x2 = -np.cos(theta) * (radius + width)
    y2 = np.sin(theta) *  (radius + width)
    stack2 = np.column_stack([np.flip(x2, axis=0), np.flip(y2,axis=0)])
    #add the first values from the first set to close the polygon
    np.append(stack2, [[x1[0],y1[0]]], axis=0)
    arcArray = np.concatenate((stack1,stack2), axis=0)
    return patches.Polygon(arcArray, True, **kwargs), ((x1, y1), (x2, y2)) 

Example 6

def orthogonalization_matrix(lengths, angles):
    """Return orthogonalization matrix for crystallographic cell coordinates.

    Angles are expected in degrees.

    The de-orthogonalization matrix is the inverse.

    >>> O = orthogonalization_matrix([10, 10, 10], [90, 90, 90])
    >>> numpy.allclose(O[:3, :3], numpy.identity(3, float) * 10)
    True
    >>> O = orthogonalization_matrix([9.8, 12.0, 15.5], [87.2, 80.7, 69.7])
    >>> numpy.allclose(numpy.sum(O), 43.063229)
    True

    """
    a, b, c = lengths
    angles = numpy.radians(angles)
    sina, sinb, _ = numpy.sin(angles)
    cosa, cosb, cosg = numpy.cos(angles)
    co = (cosa * cosb - cosg) / (sina * sinb)
    return numpy.array([
        [ a*sinb*math.sqrt(1.0-co*co),  0.0,    0.0, 0.0],
        [-a*sinb*co,                    b*sina, 0.0, 0.0],
        [ a*cosb,                       b*cosa, c,   0.0],
        [ 0.0,                          0.0,    0.0, 1.0]]) 

Example 7

def computeRotMatrix(self,Phi=False):
    
        #######################################
        # COMPUTE ROTATION MATRIX SUCH THAT m(t) = A*L(t)*A'*Hp
        # Default set such that phi1,phi2 = 0 is UXO pointed towards North
        
        if Phi is False:
            phi1 = np.radians(self.phi[0])
            phi2 = np.radians(self.phi[1])
            phi3 = np.radians(self.phi[2])
        else:
            phi1 = np.radians(Phi[0])           # Roll (CCW)
            phi2 = np.radians(Phi[1])           # Inclination (+ve is nose pointing down)
            phi3 = np.radians(Phi[2])           # Declination (degrees CW from North)
        
        # A1 = np.r_[np.c_[np.cos(phi1),-np.sin(phi1),0.],np.c_[np.sin(phi1),np.cos(phi1),0.],np.c_[0.,0.,1.]] # CCW Rotation about z-axis
        # A2 = np.r_[np.c_[1.,0.,0.],np.c_[0.,np.cos(phi2),np.sin(phi2)],np.c_[0.,-np.sin(phi2),np.cos(phi2)]] # CW Rotation about x-axis (rotates towards North)
        # A3 = np.r_[np.c_[np.cos(phi3),-np.sin(phi3),0.],np.c_[np.sin(phi3),np.cos(phi3),0.],np.c_[0.,0.,1.]] # CCW Rotation about z-axis (direction of head of object)

        A1 = np.r_[np.c_[np.cos(phi1),np.sin(phi1),0.],np.c_[-np.sin(phi1),np.cos(phi1),0.],np.c_[0.,0.,1.]] # CW Rotation about z-axis
        A2 = np.r_[np.c_[1.,0.,0.],np.c_[0.,np.cos(phi2),np.sin(phi2)],np.c_[0.,-np.sin(phi2),np.cos(phi2)]] # CW Rotation about x-axis (rotates towards North)
        A3 = np.r_[np.c_[np.cos(phi3),np.sin(phi3),0.],np.c_[-np.sin(phi3),np.cos(phi3),0.],np.c_[0.,0.,1.]] # CW Rotation about z-axis (direction of head of object)
        
        return np.dot(A3,np.dot(A2,A1)) 

Example 8

def orthogonalization_matrix(lengths, angles):
    """Return orthogonalization matrix for crystallographic cell coordinates.

    Angles are expected in degrees.

    The de-orthogonalization matrix is the inverse.

    >>> O = orthogonalization_matrix([10, 10, 10], [90, 90, 90])
    >>> numpy.allclose(O[:3, :3], numpy.identity(3, float) * 10)
    True
    >>> O = orthogonalization_matrix([9.8, 12.0, 15.5], [87.2, 80.7, 69.7])
    >>> numpy.allclose(numpy.sum(O), 43.063229)
    True

    """
    a, b, c = lengths
    angles = numpy.radians(angles)
    sina, sinb, _ = numpy.sin(angles)
    cosa, cosb, cosg = numpy.cos(angles)
    co = (cosa * cosb - cosg) / (sina * sinb)
    return numpy.array([
        [ a*sinb*math.sqrt(1.0-co*co),  0.0,    0.0, 0.0],
        [-a*sinb*co,                    b*sina, 0.0, 0.0],
        [ a*cosb,                       b*cosa, c,   0.0],
        [ 0.0,                          0.0,    0.0, 1.0]]) 

Example 9

def test_projection(self):
        """Tests the electric field projection."""

        projection = self.field.projection

        # Top-right quadrant
        a = radians(45)
        self.assertTrue(isclose(projection([0, 0], a), -4*cos(a)))
        self.assertTrue(isclose(projection([3, 0], a), 0.375*cos(a)))
        self.assertTrue(isclose(projection([0, 1], a), -sqrt(2)*cos(a)))
        self.assertTrue(isclose(projection([[0, 0], [3, 0], [0, 1]], a),
                                array([-4, 0.375, -sqrt(2)])*cos(a)).all())


        # Bottom-left quadrant
        a1 = radians(-135)
        a2 = radians(45)
        self.assertTrue(isclose(projection([0, 0], a1), 4*cos(a2)))
        self.assertTrue(isclose(projection([3, 0], a1), -0.375*cos(a2)))
        self.assertTrue(isclose(projection([0, 1], a1),
                                sqrt(2)*cos(a2)))
        self.assertTrue(isclose(projection([[0, 0], [3, 0], [0, 1]], a1),
                                array([4, -0.375, sqrt(2)])*cos(a2)).all()) 

Example 10

def cie_relative_luminance(sky_elevation, sky_azimuth=None, sun_elevation=None,
                           sun_azimuth=None, type='soc'):
    """ cie relative luminance of a sky element relative to the luminance
    at zenith
    
    angle in radians
    type is one of 'soc' (standard overcast sky), 'uoc' (uniform radiance)
    or 'clear_sky' (standard clear sky low turbidity)
    """

    if type == 'clear_sky' and (
                sun_elevation is None or sun_azimuth is None or sky_azimuth is None):
        raise ValueError, 'Clear sky requires sun position'

    if type == 'soc':
        return cie_luminance_gradation(sky_elevation, 4, -0.7)
    elif type == 'uoc':
        return cie_luminance_gradation(sky_elevation, 0, -1)
    elif type == 'clear_sky':
        return cie_luminance_gradation(sky_elevation, -1,
                                       -0.32) * cie_scattering_indicatrix(
            sun_azimuth, sun_elevation, sky_azimuth, sky_elevation, 10, -3,
            0.45)
    else:
        raise ValueError, 'Unknown sky type' 

Example 11

def ecliptic_longitude(hUTC, dayofyear, year):
    """ Ecliptic longitude

    Args:
        hUTC: fractional hour (UTC time)
        dayofyear (int):
        year (int):

    Returns:
        (float) the ecliptic longitude (degrees)

    Details:
        World Meteorological Organization (2006).Guide to meteorological
        instruments and methods of observation. Geneva, Switzerland.
    """

    jd = julian_date(hUTC, dayofyear, year)
    n = jd - 2451545
    # mean longitude (deg)
    L = numpy.mod(280.46 + 0.9856474 * n, 360)
    # mean anomaly (deg)
    g = numpy.mod(357.528 + 0.9856003 * n, 360)

    return L + 1.915 * numpy.sin(numpy.radians(g)) + 0.02 * numpy.sin(
        numpy.radians(2 * g)) 

Example 12

def actual_sky_irradiances(dates, ghi, dhi=None, Tdew=None, longitude=_longitude, latitude=_latitude, altitude=_altitude, method='dirint'):
    """ derive a sky irradiance dataframe from actual weather data"""

    df = sun_position(dates=dates, latitude=latitude, longitude=longitude, altitude=altitude, filter_night=False)
    df['am'] = air_mass(df['zenith'], altitude)
    df['dni_extra'] = sun_extraradiation(df.index)
    if dhi is None:
        pressure = pvlib.atmosphere.alt2pres(altitude)
        dhi = diffuse_horizontal_irradiance(ghi, df['elevation'], dates, pressure=pressure, temp_dew=Tdew, method=method)
    df['ghi'] = ghi
    df['dhi'] = dhi
    el = numpy.radians(df['elevation'])
    df['dni'] = (df['ghi'] - df['dhi']) / numpy.sin(el)

    df['brightness'] = brightness(df['am'], df['dhi'], df['dni_extra'])
    df['clearness'] = clearness(df['dni'], df['dhi'], df['zenith'])

    return df.loc[(df['elevation'] > 0) & (df['ghi'] > 0) , ['azimuth', 'zenith', 'elevation', 'am', 'dni_extra', 'clearness', 'brightness', 'ghi', 'dni', 'dhi' ]] 

Example 13

def rotation_matrix(alpha, beta, gamma):
    """
    Return the rotation matrix used to rotate a set of cartesian
    coordinates by alpha radians about the z-axis, then beta radians
    about the y'-axis and then gamma radians about the z''-axis.
    """
    ALPHA = np.array([[np.cos(alpha), -np.sin(alpha), 0],
                        [np.sin(alpha), np.cos(alpha), 0],
                        [0, 0, 1]])
    BETA = np.array([[np.cos(beta), 0, np.sin(beta)],
                        [0, 1, 0],
                        [-np.sin(beta), 0, np.cos(beta)]])
    GAMMA = np.array([[np.cos(gamma), -np.sin(gamma), 0],
                        [np.sin(gamma), np.cos(gamma), 0],
                        [0, 0, 1]])
    R = ALPHA.dot(BETA).dot(GAMMA)
    return(R) 

Example 14

def __init__(self, lat0, lon0, depth0, nlat, nlon, ndepth, dlat, dlon, ddepth):
# NOTE: Origin of spherical coordinate system and geographic coordinate
# system is not the same!
# Geographic coordinate system
        self.lat0, self.lon0, self.depth0 =\
                seispy.coords.as_geographic([lat0, lon0, depth0])
        self.nlat, self.nlon, self.ndepth = nlat, nlon, ndepth
        self.dlat, self.dlon, self.ddepth = dlat, dlon, ddepth
# Spherical/Pseudospherical coordinate systems
        self.nrho = self.ndepth
        self.ntheta = self.nlambda = self.nlat
        self.nphi = self.nlon
        self.drho = self.ddepth
        self.dtheta = self.dlambda = np.radians(self.dlat)
        self.dphi = np.radians(self.dlon)
        self.rho0 = seispy.constants.EARTH_RADIUS\
                - (self.depth0 + (self.ndepth - 1) * self.ddepth)
        self.lambda0 = np.radians(self.lat0)
        self.theta0 = ?/2 - (self.lambda0 + (self.nlambda - 1) * self.dlambda)
        self.phi0 = np.radians(self.lon0) 

Example 15

def _add_triangular_sides(self, xy_mask, angle, y_top_right, y_bot_left,
                              x_top_right, x_bot_left, n_material):
        angle = np.radians(angle)
        trap_len = (y_top_right - y_bot_left) / np.tan(angle)
        num_x_iterations = round(trap_len/self.x_step)
        y_per_iteration = round(self.y_pts / num_x_iterations)

        lhs_x_start_index = int(x_bot_left/ self.x_step + 0.5)
        rhs_x_stop_index = int(x_top_right/ self.x_step + 1 + 0.5)

        for i, _ in enumerate(xy_mask):
            xy_mask[i][:lhs_x_start_index] = False
            xy_mask[i][lhs_x_start_index:rhs_x_stop_index] = True

            if i % y_per_iteration == 0:
                lhs_x_start_index -= 1
                rhs_x_stop_index += 1

        self.n[xy_mask] = n_material
        return self.n 

Example 16

def orthogonalization_matrix(lengths, angles):
    """Return orthogonalization matrix for crystallographic cell coordinates.

    Angles are expected in degrees.

    The de-orthogonalization matrix is the inverse.

    >>> O = orthogonalization_matrix([10, 10, 10], [90, 90, 90])
    >>> numpy.allclose(O[:3, :3], numpy.identity(3, float) * 10)
    True
    >>> O = orthogonalization_matrix([9.8, 12.0, 15.5], [87.2, 80.7, 69.7])
    >>> numpy.allclose(numpy.sum(O), 43.063229)
    True

    """
    a, b, c = lengths
    angles = numpy.radians(angles)
    sina, sinb, _ = numpy.sin(angles)
    cosa, cosb, cosg = numpy.cos(angles)
    co = (cosa * cosb - cosg) / (sina * sinb)
    return numpy.array([
        [ a*sinb*math.sqrt(1.0-co*co),  0.0,    0.0, 0.0],
        [-a*sinb*co,                    b*sina, 0.0, 0.0],
        [ a*cosb,                       b*cosa, c,   0.0],
        [ 0.0,                          0.0,    0.0, 1.0]]) 

Example 17

def semiMajAmp(m1,m2,inc,ecc,p):
    """
    K = [(2*pi*G)/p]^(1/3) [m2*sin(i)/m2^(2/3)*sqrt(1-e^2)]
    units:
    K [m/s]
    m1 [Msun]
    m2 [Mj]
    p [yrs]
    inc [deg]
    """
    pSecs = p*sec_per_year
    m1KG = m1*const.M_sun.value
    A = ((2.0*np.pi*const.G.value)/pSecs)**(1.0/3.0)
    B = (m2*const.M_jup.value*np.sin(np.radians(inc)))
    C = m1KG**(2.0/3.0)*np.sqrt(1.0-ecc**2.0)
    print('Resulting K is '+repr(A*(B/C)))
    #return A*(B/C) 

Example 18

def inc_prior_fn(self, inc):
        ret = 1.0
        if (self.inc_prior==False) or (self.inc_prior=="uniform"):
            ret = self.uniform_fn(inc,7)
        else:
            if inc not in [0.0,90.0,180.0]:
                mn = self.mins_ary[7]
                mx = self.maxs_ary[7]
                # Account for case where min = -(max), as it causes error in denominator
                if mn == -1*mx:
                    mn = mn-0.1
                mn_rad = np.radians(mn)
                mx_rad = np.radians(mx)
                inc_rad = np.radians(inc)
                if (self.inc_prior == True) or (self.inc_prior == 'sin'):
                    ret = np.abs(np.sin(inc_rad)) / np.abs(np.cos(mn_rad)-np.cos(mx_rad))
                elif self.inc_prior == 'cos':
                    ret =  np.abs(np.cos(inc_rad)) / np.abs(np.cos(mn_rad)-np.cos(mx_rad))
        #if ret==0: ret=-np.inf
        return ret 

Example 19

def inc_prior_fn(self, inc):
        ret = 1.0
        if (self.inc_prior==False) or (self.inc_prior=="uniform"):
            ret = self.uniform_fn(inc,7)
        else:
            if inc not in [0.0,90.0,180.0]:
                mn = self.mins_ary[7]
                mx = self.maxs_ary[7]
                # Account for case where min = -(max), as it causes error in denominator
                if mn == -1*mx:
                    mn = mn-0.1
                mn_rad = np.radians(mn)
                mx_rad = np.radians(mx)
                inc_rad = np.radians(inc)
                if (self.inc_prior == True) or (self.inc_prior == 'sin'):
                    ret = np.abs(np.sin(inc_rad)) / np.abs(np.cos(mn_rad)-np.cos(mx_rad))
                elif self.inc_prior == 'cos':
                    ret =  np.abs(np.cos(inc_rad)) / np.abs(np.cos(mn_rad)-np.cos(mx_rad))
        #if ret==0: ret=-np.inf
        return ret 

Example 20

def hav_dist(locs1, locs2):
    """
    Return a distance matrix between two set of coordinates.
    Use geometric distance (default) or haversine distance (if longlat=True).

    Parameters
    ----------
    locs1 : numpy.array
        The first set of coordinates as [(long, lat), (long, lat)].
    locs2 : numpy.array
        The second set of coordinates as [(long, lat), (long, lat)].

    Returns
    -------
    mat_dist : numpy.array
        The distance matrix between locs1 and locs2
    """
    locs1 = np.radians(locs1)
    locs2 = np.radians(locs2)
    cos_lat1 = np.cos(locs1[..., 0])
    cos_lat2 = np.cos(locs2[..., 0])
    cos_lat_d = np.cos(locs1[..., 0] - locs2[..., 0])
    cos_lon_d = np.cos(locs1[..., 1] - locs2[..., 1])
    return 6367000 * np.arccos(
        cos_lat_d - cos_lat1 * cos_lat2 * (1 - cos_lon_d)) 

Example 21

def _get_rotation_matrix(axis, angle):
    """
    Helper function to generate a rotation matrix for an x, y, or z axis
    Used in get_major_angles
    """
    cos = np.cos
    sin = np.sin
    angle = np.radians(angle)
    if axis == 2:
        # z axis
        return np.array([[cos(angle), -sin(angle), 0], [sin(angle), cos(angle), 0], [0, 0, 1]])
    if axis == 1:
        # y axis
        return np.array([[cos(angle), 0, sin(angle)], [0, 1, 0], [-sin(angle), 0, cos(angle)]])
    else:
        # x axis
        return np.array([[1, 0, 0], [0, cos(angle), -sin(angle)], [0, sin(angle), cos(angle)]]) 

Example 22

def euler(xyz, order='xyz', units='deg'):
    if not hasattr(xyz, '__iter__'):
        xyz = [xyz]
    if units == 'deg':
        xyz = np.radians(xyz)
    r = np.eye(3)
    for theta, axis in zip(xyz, order):
        c = np.cos(theta)
        s = np.sin(theta)
        if axis == 'x':
            r = np.dot(np.array([[1, 0, 0], [0, c, -s], [0, s, c]]), r)
        if axis == 'y':
            r = np.dot(np.array([[c, 0, s], [0, 1, 0], [-s, 0, c]]), r)
        if axis == 'z':
            r = np.dot(np.array([[c, -s, 0], [s, c, 0], [0, 0, 1]]), r)
    return r 

Example 23

def get_q_per_pixel(self):
        '''Gets the delta-q associated with a single pixel. This is computed in
        the small-angle limit, so it should only be considered approximate.
        For instance, wide-angle detectors will have different delta-q across
        the detector face.'''
        
        if self.q_per_pixel is not None:
            return self.q_per_pixel
        
        c = (self.pixel_size_um/1e6)/self.distance_m
        twotheta = np.arctan(c) # radians
        
        self.q_per_pixel = 2.0*self.get_k()*np.sin(twotheta/2.0)
        
        return self.q_per_pixel
    
    
    # Maps
    ######################################## 

Example 24

def add_motion_spikes(motion_mat,frequency,severity,TR):
	time = motion_mat[:,0]
	max_translation = 5 * severity / 1000 * np.sqrt(2*np.pi)#Max translations in m, factor of sqrt(2*pi) accounts for normalisation factor in norm.pdf later on
	max_rotation = np.radians(5 * severity) *np.sqrt(2*np.pi) #Max rotation in rad
	time_blocks = np.floor(time[-1]/TR).astype(np.int32)
	for i in range(time_blocks):
		if np.random.uniform(0,1) < frequency: #Decide whether to add spike
			for j in range(1,4):
				if np.random.uniform(0,1) < 1/6:
					motion_mat[:,j] = motion_mat[:,j] \
					+ max_translation * random.uniform(-1,1) \
					* norm.pdf(time,loc = (i+0.5)*TR,scale = TR/5)
			for j in range(4,7):
				if np.random.uniform(0,1) < 1/6:
					motion_mat[:,j] = motion_mat[:,j] \
					+ max_rotation * random.uniform(-1,1) \
					* norm.pdf(time,loc = (i+0.5 + np.random.uniform(-0.25,-.25))*TR,scale = TR/5)
	return motion_mat 

Example 25

def tiltFactor(self, midpointdepth=None,
                   printAvAngle=False):
        '''
        get tilt factor from inverse distance law
        https://en.wikipedia.org/wiki/Inverse-square_law
        '''
        # TODO: can also be only def. with FOV, rot, tilt
        beta2 = self.viewAngle(midpointdepth=midpointdepth)
        try:
            angles, vals = getattr(
                emissivity_vs_angle, self.opts['material'])()
        except AttributeError:
            raise AttributeError("material[%s] is not in list of know materials: %s" % (
                self.opts['material'], [o[0] for o in getmembers(emissivity_vs_angle)
                                        if isfunction(o[1])]))
        if printAvAngle:
            avg_angle = beta2[self.foreground()].mean()
            print('angle: %s DEG' % np.degrees(avg_angle))

        # use averaged angle instead of beta2 to not overemphasize correction
        normEmissivity = np.clip(
            InterpolatedUnivariateSpline(
                np.radians(angles), vals)(beta2), 0, 1)
        return normEmissivity 

Example 26

def setPose(self, obj_center=None, distance=None,
                rotation=None, tilt=None, pitch=None):
        tvec, rvec = self.pose()

        if distance is not None:
            tvec[2, 0] = distance
        if obj_center is not None:
            tvec[0, 0] = obj_center[0]
            tvec[1, 0] = obj_center[1]

        if rotation is None and tilt is None:
            return rvec
        r, t, p = rvec2euler(rvec)
        if rotation is not None:
            r = np.radians(rotation)
        if tilt is not None:
            t = np.radians(tilt)
        if pitch is not None:
            p = np.radians(pitch)
        rvec = euler2rvec(r, t, p)

        self._pose = tvec, rvec 

Example 27

def hav(alpha):
    """ Formula for haversine

    Parameters
    ----------
    alpha : (float)
        Angle in radians

    Returns
    --------
    hav_alpha : (float)
        Haversine of alpha, equal to the square of the sine of half-alpha
    """
    hav_alpha = np.sin(alpha * 0.5)**2

    return hav_alpha 

Example 28

def archav(hav):
    """ Formula for the inverse haversine

    Parameters
    -----------
    hav : (float)
        Haversine of an angle

    Returns
    ---------
    alpha : (float)
        Angle in radians
    """
    alpha = 2.0 * np.arcsin(np.sqrt(hav))

    return alpha 

Example 29

def spherical_to_cartesian(s,degrees=True,normalize=False):
   '''
   Takes a vector in spherical coordinates and converts it to cartesian.
   Assumes the input vector is given as [radius,colat,lon] 
   '''

   if degrees:
      s[1] = np.radians(s[1])  
      s[2] = np.radians(s[2])

   x1 = s[0]*np.sin(s[1])*np.cos(s[2])
   x2 = s[0]*np.sin(s[1])*np.sin(s[2])
   x3 = s[0]*np.cos(s[1])

   x = [x1,x2,x3]

   if normalize:
      x /= np.linalg.norm(x)
   return x 

Example 30

def rotate_delays(lat_r,lon_r,lon_0=0.0,lat_0=0.0,degrees=0):
   '''
   Rotates the source and receiver of a trace object around an
   arbitrary axis.
   '''

   alpha = np.radians(degrees)
   colat_r = 90.0-lat_r
   colat_0 = 90.0-lat_0

   x_r = lon_r - lon_0
   y_r = colat_0 - colat_r

   #rotate receivers
   lat_rotated = 90.0-colat_0+x_r*np.sin(alpha) + y_r*np.cos(alpha)
   lon_rotated = lon_0+x_r*np.cos(alpha) - y_r*np.sin(alpha)

   return lat_rotated, lon_rotated 

Example 31

def rotate_delays(lat_r,lon_r,lon_0=0.0,lat_0=0.0,degrees=0):
   '''
   Rotates the source and receiver of a trace object around an
   arbitrary axis.
   '''

   alpha = np.radians(degrees)
   colat_r = 90.0-lat_r
   colat_0 = 90.0-lat_0

   x_r = lon_r - lon_0
   y_r = colat_0 - colat_r

   #rotate receivers
   lat_rotated = 90.0-colat_0+x_r*np.sin(alpha) + y_r*np.cos(alpha)
   lon_rotated = lon_0+x_r*np.cos(alpha) - y_r*np.sin(alpha)

   return lat_rotated, lon_rotated 

Example 32

def rotate_setup(self,lon_0=60.0,colat_0=90.0,degrees=0):

     alpha = np.radians(degrees)
     lon_s = self.sy
     lon_r = self.ry
     colat_s = self.sx
     colat_r = self.rx

     x_s = lon_s - lon_0
     y_s = colat_0 - colat_s
     x_r = lon_r - lon_0
     y_r = colat_0 - colat_r
 
     #rotate receiver
     self.rx = colat_0+x_r*np.sin(alpha) + y_r*np.cos(alpha)
     self.ry = lon_0+x_r*np.cos(alpha) - y_r*np.sin(alpha)

     #rotate source
     self.sx = colat_0+x_s*np.sin(alpha) + y_s*np.cos(alpha)
     self.sy = lon_0+x_s*np.cos(alpha) - y_s*np.sin(alpha)

  #########################################################################
  # Plot map of earthquake and station
  ######################################################################### 

Example 33

def _calcArea_(self, v1, v2):
        """
        Private method to calculate the area covered by a spherical
        quadrilateral with one corner defined by the normal vectors
        of the two intersecting great circles.

        INPUTS:
            v1, v2:  float array(3), the normal vectors

        RETURNS:
            area: float, the area given in square radians
        """
        angle = self.calcAngle(v1, v2)
        area = (4*angle - 2*np.math.pi)

        return area 

Example 34

def _rotVector_(self, v, angle, axis):
        """
        Rotate a vector by an angle around an axis
        INPUTS:
            v:    3-dim float array
            angle:    float, the rotation angle in radians
            axis: string, 'x', 'y', or 'z'

        RETURNS:
            float array(3):  the rotated vector
        """
        # axisd = {'x':[1,0,0], 'y':[0,1,0], 'z':[0,0,1]}

        # construct quaternion and rotate...
        rot = cgt.quat()
        rot.fromAngleAxis(angle, axis)
        return list(rot.rotateVec(v)) 

Example 35

def _geodetic_to_cartesian(lat, lon, alt):
    """Conversion from latitude, longitue and altitude coordinates to
    cartesian with respect to an ellipsoid

    Args:
        lat (float): Latitude in radians
        lon (float): Longitue in radians
        alt (float): Altitude to sea level in meters

    Return:
        numpy.array: 3D element (in meters)
    """
    C = Earth.r / np.sqrt(1 - (Earth.e * np.sin(lat)) ** 2)
    S = Earth.r * (1 - Earth.e ** 2) / np.sqrt(1 - (Earth.e * np.sin(lat)) ** 2)
    r_d = (C + alt) * np.cos(lat)
    r_k = (S + alt) * np.sin(lat)

    norm = np.sqrt(r_d ** 2 + r_k ** 2)
    return norm * np.array([
        np.cos(lat) * np.cos(lon),
        np.cos(lat) * np.sin(lon),
        np.sin(lat)
    ]) 

Example 36

def orthogonalization_matrix(lengths, angles):
    """Return orthogonalization matrix for crystallographic cell coordinates.

    Angles are expected in degrees.

    The de-orthogonalization matrix is the inverse.

    >>> O = orthogonalization_matrix([10, 10, 10], [90, 90, 90])
    >>> numpy.allclose(O[:3, :3], numpy.identity(3, float) * 10)
    True
    >>> O = orthogonalization_matrix([9.8, 12.0, 15.5], [87.2, 80.7, 69.7])
    >>> numpy.allclose(numpy.sum(O), 43.063229)
    True

    """
    a, b, c = lengths
    angles = numpy.radians(angles)
    sina, sinb, _ = numpy.sin(angles)
    cosa, cosb, cosg = numpy.cos(angles)
    co = (cosa * cosb - cosg) / (sina * sinb)
    return numpy.array([
        [ a*sinb*math.sqrt(1.0-co*co),  0.0,    0.0, 0.0],
        [-a*sinb*co,                    b*sina, 0.0, 0.0],
        [ a*cosb,                       b*cosa, c,   0.0],
        [ 0.0,                          0.0,    0.0, 1.0]]) 

Example 37

def orthogonalization_matrix(lengths, angles):
    """Return orthogonalization matrix for crystallographic cell coordinates.

    Angles are expected in degrees.

    The de-orthogonalization matrix is the inverse.

    >>> O = orthogonalization_matrix([10, 10, 10], [90, 90, 90])
    >>> numpy.allclose(O[:3, :3], numpy.identity(3, float) * 10)
    True
    >>> O = orthogonalization_matrix([9.8, 12.0, 15.5], [87.2, 80.7, 69.7])
    >>> numpy.allclose(numpy.sum(O), 43.063229)
    True

    """
    a, b, c = lengths
    angles = numpy.radians(angles)
    sina, sinb, _ = numpy.sin(angles)
    cosa, cosb, cosg = numpy.cos(angles)
    co = (cosa * cosb - cosg) / (sina * sinb)
    return numpy.array([
        [ a*sinb*math.sqrt(1.0-co*co),  0.0,    0.0, 0.0],
        [-a*sinb*co,                    b*sina, 0.0, 0.0],
        [ a*cosb,                       b*cosa, c,   0.0],
        [ 0.0,                          0.0,    0.0, 1.0]]) 

Example 38

def orthogonalization_matrix(lengths, angles):
    """Return orthogonalization matrix for crystallographic cell coordinates.

    Angles are expected in degrees.

    The de-orthogonalization matrix is the inverse.

    >>> O = orthogonalization_matrix([10, 10, 10], [90, 90, 90])
    >>> numpy.allclose(O[:3, :3], numpy.identity(3, float) * 10)
    True
    >>> O = orthogonalization_matrix([9.8, 12.0, 15.5], [87.2, 80.7, 69.7])
    >>> numpy.allclose(numpy.sum(O), 43.063229)
    True

    """
    a, b, c = lengths
    angles = numpy.radians(angles)
    sina, sinb, _ = numpy.sin(angles)
    cosa, cosb, cosg = numpy.cos(angles)
    co = (cosa * cosb - cosg) / (sina * sinb)
    return numpy.array([
        [ a*sinb*math.sqrt(1.0-co*co),  0.0,    0.0, 0.0],
        [-a*sinb*co,                    b*sina, 0.0, 0.0],
        [ a*cosb,                       b*cosa, c,   0.0],
        [ 0.0,                          0.0,    0.0, 1.0]]) 

Example 39

def set_radmask(self, cluster, mpcscale):
        """
        Assign mask (0/1) values to maskgals for a given cluster

        parameters
        ----------
        cluster: Cluster object
        mpcscale: float
            scaling to go from mpc to degrees (check units) at cluster redshift

        results
        -------
        sets maskgals.mark

        """

        # note this probably can be in the superclass, no?
        ras = cluster.ra + self.maskgals.x/(mpcscale*SEC_PER_DEG)/np.cos(np.radians(cluster.dec))
        decs = cluster.dec + self.maskgals.y/(mpcscale*SEC_PER_DEG)
        self.maskgals.mark = self.compute_radmask(ras,decs) 

Example 40

def _calc_bkg_density(self, r, chisq, refmag):
        """
        Internal method to compute background filter

        parameters
        ----------
        bkg: Background object
           background
        cosmo: Cosmology object
           cosmology scaling info

        returns
        -------

        bcounts: float array
            b(x) for the cluster
        """
        mpc_scale = np.radians(1.) * self.cosmo.Dl(0, self.z) / (1 + self.z)**2
        sigma_g = self.bkg.sigma_g_lookup(self.z, chisq, refmag)
        return 2 * np.pi * r * (sigma_g/mpc_scale**2) 

Example 41

def orthogonalization_matrix(lengths, angles):
    """Return orthogonalization matrix for crystallographic cell coordinates.

    Angles are expected in degrees.

    The de-orthogonalization matrix is the inverse.

    >>> O = orthogonalization_matrix([10, 10, 10], [90, 90, 90])
    >>> numpy.allclose(O[:3, :3], numpy.identity(3, float) * 10)
    True
    >>> O = orthogonalization_matrix([9.8, 12.0, 15.5], [87.2, 80.7, 69.7])
    >>> numpy.allclose(numpy.sum(O), 43.063229)
    True

    """
    a, b, c = lengths
    angles = numpy.radians(angles)
    sina, sinb, _ = numpy.sin(angles)
    cosa, cosb, cosg = numpy.cos(angles)
    co = (cosa * cosb - cosg) / (sina * sinb)
    return numpy.array([
        [ a*sinb*math.sqrt(1.0-co*co),  0.0,    0.0, 0.0],
        [-a*sinb*co,                    b*sina, 0.0, 0.0],
        [ a*cosb,                       b*cosa, c,   0.0],
        [ 0.0,                          0.0,    0.0, 1.0]]) 

Example 42

def apply_grad_cartesian_tensor(grad_X, zmat_dist):
    """Apply the gradient for transformation to cartesian space onto zmat_dist.

    Args:
        grad_X (:class:`numpy.ndarray`): A ``(3, n, n, 3)`` array.
            The mathematical details of the index layout is explained in
            :meth:`~chemcoord.Cartesian.get_grad_zmat()`.
        zmat_dist (:class:`~chemcoord.Zmat`):
            Distortions in Zmatrix space.

    Returns:
        :class:`~chemcoord.Cartesian`: Distortions in cartesian space.
    """
    columns = ['bond', 'angle', 'dihedral']
    C_dist = zmat_dist.loc[:, columns].values.T
    try:
        C_dist = C_dist.astype('f8')
        C_dist[[1, 2], :] = np.radians(C_dist[[1, 2], :])
    except (TypeError, AttributeError):
        C_dist[[1, 2], :] = sympy.rad(C_dist[[1, 2], :])
    cart_dist = np.tensordot(grad_X, C_dist, axes=([3, 2], [0, 1])).T
    from chemcoord.cartesian_coordinates.cartesian_class_main import Cartesian
    return Cartesian(atoms=zmat_dist['atom'],
                     coords=cart_dist, index=zmat_dist.index) 

Example 43

def drawSkymapCatalog(ax,lon,lat,**kwargs):
    mapping = {
        'ait':'aitoff',
        'mol':'mollweide',
        'lam':'lambert',
        'ham':'hammer'
    }
    kwargs.setdefault('proj','aitoff')
    kwargs.setdefault('s',2)
    kwargs.setdefault('marker','.')
    kwargs.setdefault('c','k')

    proj = kwargs.pop('proj')
    projection = mapping.get(proj,proj)
    #ax.grid()
    # Convert from 
    # [0. < lon < 360] -> [-pi < lon < pi]
    # [-90 < lat < 90] -> [-pi/2 < lat < pi/2]
    lon,lat= np.radians([lon-360.*(lon>180),lat])
    ax.scatter(lon,lat,**kwargs) 

Example 44

def healpixMap(nside, lon, lat, fill_value=0., nest=False):
    """
    Input (lon, lat) in degrees instead of (theta, phi) in radians.
    Returns HEALPix map at the desired resolution 
    """

    lon_median, lat_median = np.median(lon), np.median(lat)
    max_angsep = np.max(ugali.utils.projector.angsep(lon, lat, lon_median, lat_median))
    
    pix = angToPix(nside, lon, lat, nest=nest)
    if max_angsep < 10:
        # More efficient histograming for small regions of sky
        m = np.tile(fill_value, healpy.nside2npix(nside))
        pix_subset = ugali.utils.healpix.angToDisc(nside, lon_median, lat_median, max_angsep, nest=nest)
        bins = np.arange(np.min(pix_subset), np.max(pix_subset) + 1)
        m_subset = np.histogram(pix, bins=bins - 0.5)[0].astype(float)
        m[bins[0:-1]] = m_subset
    else:
        m = np.histogram(pix, np.arange(hp.nside2npix(nside) + 1))[0].astype(float)
    if fill_value != 0.:
        m[m == 0.] = fill_value
    return m

############################################################ 

Example 45

def galToCel(ll, bb):
    """
    Converts Galactic (deg) to Celestial J2000 (deg) coordinates
    """
    bb = numpy.radians(bb)
    sin_bb = numpy.sin(bb)
    cos_bb = numpy.cos(bb)

    ll = numpy.radians(ll)
    ra_gp = numpy.radians(192.85948)
    de_gp = numpy.radians(27.12825)
    lcp = numpy.radians(122.932)

    sin_lcp_ll = numpy.sin(lcp - ll)
    cos_lcp_ll = numpy.cos(lcp - ll)

    sin_d = (numpy.sin(de_gp) * sin_bb) \
            + (numpy.cos(de_gp) * cos_bb * cos_lcp_ll)
    ramragp = numpy.arctan2(cos_bb * sin_lcp_ll,
                            (numpy.cos(de_gp) * sin_bb) \
                            - (numpy.sin(de_gp) * cos_bb * cos_lcp_ll))
    dec = numpy.arcsin(sin_d)
    ra = (ramragp + ra_gp + (2. * numpy.pi)) % (2. * numpy.pi)
    return numpy.degrees(ra), numpy.degrees(dec) 

Example 46

def ang2const(lon,lat,coord='gal'):
    import ephem

    scalar = np.isscalar(lon)
    lon = np.array(lon,copy=False,ndmin=1)
    lat = np.array(lat,copy=False,ndmin=1)

    if coord.lower() == 'cel':
        ra,dec = lon,lat
    elif coord.lower() == 'gal':
        ra,dec = gal2cel(lon,lat)
    else:
        msg = "Unrecognized coordinate"
        raise Exception(msg)

    x,y = np.radians([ra,dec])
    const = [ephem.constellation(coord) for coord in zip(x,y)]
    if scalar: return const[0]
    return const 

Example 47

def orthogonalization_matrix(lengths, angles):
    """Return orthogonalization matrix for crystallographic cell coordinates.

    Angles are expected in degrees.

    The de-orthogonalization matrix is the inverse.

    >>> O = orthogonalization_matrix([10, 10, 10], [90, 90, 90])
    >>> numpy.allclose(O[:3, :3], numpy.identity(3, float) * 10)
    True
    >>> O = orthogonalization_matrix([9.8, 12.0, 15.5], [87.2, 80.7, 69.7])
    >>> numpy.allclose(numpy.sum(O), 43.063229)
    True

    """
    a, b, c = lengths
    angles = numpy.radians(angles)
    sina, sinb, _ = numpy.sin(angles)
    cosa, cosb, cosg = numpy.cos(angles)
    co = (cosa * cosb - cosg) / (sina * sinb)
    return numpy.array([
        [ a*sinb*math.sqrt(1.0-co*co),  0.0,    0.0, 0.0],
        [-a*sinb*co,                    b*sina, 0.0, 0.0],
        [ a*cosb,                       b*cosa, c,   0.0],
        [ 0.0,                          0.0,    0.0, 1.0]]) 

Example 48

def __init__(self, survey_features=None, condition_features=None, time_lag=10.,
                 filtername='r', twi_change=-18.):
        """
        Paramters
        ---------
        time_lag : float (10.)
            If there is a gap between observations longer than this, let the filter change (minutes)
        twi_change : float (-18.)
            The sun altitude to consider twilight starting/ending
        """
        self.time_lag = time_lag/60./24.  # Convert to days
        self.twi_change = np.radians(twi_change)
        self.filtername = filtername
        if condition_features is None:
            self.condition_features = {}
            self.condition_features['Current_filter'] = features.Current_filter()
            self.condition_features['Sun_moon_alts'] = features.Sun_moon_alts()
            self.condition_features['Current_mjd'] = features.Current_mjd()
        if survey_features is None:
            self.survey_features = {}
            self.survey_features['Last_observation'] = features.Last_observation()

        super(Strict_filter_basis_function, self).__init__(survey_features=self.survey_features,
                                                           condition_features=self.condition_features) 

Example 49

def treexyz(ra, dec):
    """
    Utility to convert RA,dec postions in x,y,z space, useful for constructing KD-trees.

    Parameters
    ----------
    ra : float or array
        RA in radians
    dec : float or array
        Dec in radians

    Returns
    -------
    x,y,z : floats or arrays
        The position of the given points on the unit sphere.
    """
    # Note ra/dec can be arrays.
    x = np.cos(dec) * np.cos(ra)
    y = np.cos(dec) * np.sin(ra)
    z = np.sin(dec)
    return x, y, z 

Example 50

def orthogonalization_matrix(lengths, angles):
    """Return orthogonalization matrix for crystallographic cell coordinates.

    Angles are expected in degrees.

    The de-orthogonalization matrix is the inverse.

    >>> O = orthogonalization_matrix((10., 10., 10.), (90., 90., 90.))
    >>> numpy.allclose(O[:3, :3], numpy.identity(3, float) * 10)
    True
    >>> O = orthogonalization_matrix([9.8, 12.0, 15.5], [87.2, 80.7, 69.7])
    >>> numpy.allclose(numpy.sum(O), 43.063229)
    True

    """
    a, b, c = lengths
    angles = numpy.radians(angles)
    sina, sinb, _ = numpy.sin(angles)
    cosa, cosb, cosg = numpy.cos(angles)
    co = (cosa * cosb - cosg) / (sina * sinb)
    return numpy.array((
        ( a*sinb*math.sqrt(1.0-co*co),  0.0,    0.0, 0.0),
        (-a*sinb*co,                    b*sina, 0.0, 0.0),
        ( a*cosb,                       b*cosa, c,   0.0),
        ( 0.0,                          0.0,    0.0, 1.0)),
        dtype=numpy.float64) 
点赞