Python numpy.radians() Examples

The following are 30 code examples of numpy.radians(). You can vote up the ones you like or vote down the ones you don't like, and go to the original project or source file by following the links above each example. You may also want to check out all available functions/classes of the module numpy , or try the search function .
Example #1
Source File: rotated_mapping_tools.py    From LSDMappingTools with MIT License 7 votes vote down vote up
def HaversineDistance(lon1,lat1,lon2,lat2):
    """
    Function to calculate the great circle distance between two points
    using the Haversine formula
    """
    R = 6371. #Mean radius of the Earth

    # convert decimal degrees to radians
    lon1, lat1, lon2, lat2 = map(np.radians, [lon1, lat1, lon2, lat2])

    # haversine formula
    dlon = lon2 - lon1
    dlat = lat2 - lat1
    a = np.sin(dlat/2.)**2. + np.cos(lat1) * np.cos(lat2) * np.sin(dlon/2.)**2.
    c = 2.*np.arcsin(np.sqrt(a))
    distance = R * c

    return distance 
Example #2
Source File: omni_hro.py    From pysat with BSD 3-Clause "New" or "Revised" License 6 votes vote down vote up
def calculate_dayside_reconnection(inst):
    """ Calculate the dayside reconnection rate (Milan et al. 2014)

    Parameters
    -----------
    inst : pysat.Instrument
        Instrument with OMNI HRO data, requires BYZ_GSM and clock_angle

    Notes
    --------
    recon_day = 3.8 Re (Vx / 4e5 m/s)^1/3 Vx B_yz (sin(theta/2))^9/2
    """
    rearth = 6371008.8
    sin_htheta = np.power(np.sin(np.radians(0.5 * inst['clock_angle'])), 4.5)
    byz = inst['BYZ_GSM'] * 1.0e-9
    vx = inst['flow_speed'] * 1000.0

    recon_day = 3.8 * rearth * vx * byz * sin_htheta * np.power((vx / 4.0e5),
                                                                1.0/3.0)
    inst['recon_day'] = pds.Series(recon_day, index=inst.data.index)
    return 
Example #3
Source File: make_hpx_files.py    From gamma-astro-data-formats with Creative Commons Attribution 4.0 International 6 votes vote down vote up
def coords_to_vec(lon, lat):
    """ Converts longitute and latitude coordinates to a unit 3-vector

    return array(3,n) with v_x[i],v_y[i],v_z[i] = directional cosines
    """
    phi = np.radians(lon)
    theta = (np.pi / 2) - np.radians(lat)
    sin_t = np.sin(theta)
    cos_t = np.cos(theta)

    xVals = sin_t * np.cos(phi)
    yVals = sin_t * np.sin(phi)
    zVals = cos_t

    # Stack them into the output array
    out = np.vstack((xVals, yVals, zVals)).swapaxes(0, 1)
    return out 
Example #4
Source File: intensity_measures.py    From gmpe-smtk with GNU Affero General Public License v3.0 6 votes vote down vote up
def rotipp(acceleration_x, time_step_x, acceleration_y, time_step_y, periods,
        percentile, damping=0.05, units="cm/s/s", method="Nigam-Jennings"):
    """
    Returns the rotationally independent spectrum RotIpp as defined by
    Boore (2010)
    """
    if np.fabs(time_step_x - time_step_y) > 1E-10:
        raise ValueError("Record pair must have the same time-step!")
    acceleration_x, acceleration_y = equalise_series(acceleration_x,
                                                     acceleration_y)
    target, rota, rotv, rotd, angles = rotdpp(acceleration_x, time_step_x,
                                              acceleration_y, time_step_y,
                                              periods, percentile, damping,
                                              units, method)
    locn, penalty = _get_gmrotd_penalty(
        np.hstack([target["PGA"],target["Pseudo-Acceleration"]]),
        rota)
    target_theta = np.radians(angles[locn])
    arotpp = acceleration_x * np.cos(target_theta) +\
        acceleration_y * np.sin(target_theta)
    spec = get_response_spectrum(arotpp, time_step_x, periods, damping, units,
        method)[0]
    spec["GMRot{:2.0f}".format(percentile)] = target
    return spec 
Example #5
Source File: topography.py    From typhon with MIT License 6 votes vote down vote up
def _latlon_to_cart(lat, lon, R = typhon.constants.earth_radius):
    """
    Simple conversion of latitude and longitude to Cartesian coordinates.
    Approximates the Earth as sphere with radius :code:`R` and computes
    cartesian x, y, z coordinates with the center of the Earth as origin.

    Args:
        lat: Array of latitude coordinates.
        lon: Array of longitude coordinates.
        R: The radius to assume.
    Returns:
        Tuple :code:`(x, y, z)` of arrays :code:`x, y, z` containing the
        resulting x-, y- and z-coordinates.
    """
    lat = np.radians(lat)
    lon = np.radians(lon)
    x = R * np.cos(lat) * np.cos(lon)
    y = R * np.cos(lat) * np.sin(lon)
    z = R * np.sin(lat)
    return x, y, z 
Example #6
Source File: transform.py    From sixd_toolkit with MIT License 6 votes vote down vote up
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
Source File: lib.py    From sea_ice_drift with GNU General Public License v3.0 6 votes vote down vote up
def get_displacement_km(n1, x1, y1, n2, x2, y2):
    ''' Find displacement in kilometers using Haversine
        http://www.movable-type.co.uk/scripts/latlong.html
    Parameters
    ----------
        n1 : First Nansat object
        x1 : 1D vector - X coordinates of keypoints on image 1
        y1 : 1D vector - Y coordinates of keypoints on image 1
        n2 : Second Nansat object
        x1 : 1D vector - X coordinates of keypoints on image 2
        y1 : 1D vector - Y coordinates of keypoints on image 2
    Returns
    -------
        h : 1D vector - total displacement, km
    '''
    lon1, lat1 = n1.transform_points(x1, y1)
    lon2, lat2 = n2.transform_points(x2, y2)

    lt1, ln1, lt2, ln2 = map(np.radians, (lat1, lon1, lat2, lon2))
    dlat = lt2 - lt1
    dlon = ln2 - ln1
    d = (np.sin(dlat * 0.5) ** 2 +
         np.cos(lt1) * np.cos(lt2) * np.sin(dlon * 0.5) ** 2)
    return 2 * AVG_EARTH_RADIUS * np.arcsin(np.sqrt(d)) 
Example #8
Source File: rototranslation.py    From differentiable-renderer with MIT License 6 votes vote down vote up
def __init__(self, rotation: Vector, translation: Vector, angle_unit: str, notation: str='XYZ'):
        self.rotation    = rotation
        self.translation = translation

        self.angle_unit = angle_unit
        if self.angle_unit == 'degrees':
            self.rotation = Vector(*[radians(alpha) for alpha in rotation])

        self.R_x = None
        self.R_y = None
        self.R_z = None
        self.T   = None

        self.matrix = None
        self.notation = notation

        self._update_matrix() 
Example #9
Source File: transform.py    From patch_linemod with BSD 2-Clause "Simplified" License 6 votes vote down vote up
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 #10
Source File: test_angles.py    From orbit-predictor with MIT License 6 votes vote down vote up
def test_true_to_mean(self):
        # Data from Schlesinger & Udick, 1912
        data = [
            # ecc, M (deg), ta (deg)
            (0.0, 0.0, 0.0),
            (0.05, 10.0, 11.06),
            (0.06, 30.0, 33.67),
            (0.04, 120.0, 123.87),
            (0.14, 65.0, 80.50),
            (0.19, 21.0, 30.94),
            (0.35, 65.0, 105.71),
            (0.48, 180.0, 180.0),
            (0.75, 125.0, 167.57)
        ]
        for row in data:
            ecc, expected_M, ta = row

            M = angles.ta_to_M(radians(ta), ecc)

            self.assertAlmostEqual(degrees(M), expected_M, places=1) 
Example #11
Source File: test_angles.py    From orbit-predictor with MIT License 6 votes vote down vote up
def test_mean_to_true(self):
        # Data from Schlesinger & Udick, 1912
        data = [
            # ecc, M (deg), ta (deg)
            (0.0, 0.0, 0.0),
            (0.05, 10.0, 11.06),
            (0.06, 30.0, 33.67),
            (0.04, 120.0, 123.87),
            (0.14, 65.0, 80.50),
            (0.19, 21.0, 30.94),
            (0.35, 65.0, 105.71),
            (0.48, 180.0, 180.0),
            (0.75, 125.0, 167.57)
        ]
        for row in data:
            ecc, M, expected_ta = row

            ta = angles.M_to_ta(radians(M), ecc)

            self.assertAlmostEqual(degrees(ta), expected_ta, places=2) 
Example #12
Source File: transformations.py    From rai-python with MIT License 6 votes vote down vote up
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) 
Example #13
Source File: test_angles.py    From orbit-predictor with MIT License 6 votes vote down vote up
def test_true_to_eccentric(self):
        # Data from NASA-TR-R-158
        data = [
            # ecc, E (deg), ta(deg)
            (0.0, 0.0, 0.0),
            (0.05, 10.52321, 11.05994),
            (0.10, 54.67466, 59.49810),
            (0.35, 142.27123, 153.32411),
            (0.61, 161.87359, 171.02189)
        ]
        for row in data:
            ecc, expected_E, ta = row

            E = angles.ta_to_E(radians(ta), ecc)

            self.assertAlmostEqual(degrees(E), expected_E, places=4) 
Example #14
Source File: transform.py    From AugmentedAutoencoder with MIT License 6 votes vote down vote up
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 #15
Source File: transform.py    From AugmentedAutoencoder with MIT License 6 votes vote down vote up
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 #16
Source File: transformations.py    From visual_dynamics with MIT License 6 votes vote down vote up
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
Source File: transformations.py    From 6-PACK with MIT License 6 votes vote down vote up
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 #18
Source File: transform.py    From bop_toolkit with MIT License 6 votes vote down vote up
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 #19
Source File: test_angles.py    From orbit-predictor with MIT License 5 votes vote down vote up
def test_rotate_simple(self):
        vec = np.array([1, 0, 0])

        assert_allclose(rotate(vec, 0, np.radians(90)), np.array([1, 0, 0]), atol=1e-16)
        assert_allclose(rotate(vec, 1, np.radians(90)), np.array([0, 0, -1]), atol=1e-16)
        assert_allclose(rotate(vec, 2, np.radians(90)), np.array([0, 1, 0]), atol=1e-16) 
Example #20
Source File: view.py    From pysheds with GNU General Public License v3.0 5 votes vote down vote up
def _view_rectspherebivariate(cls, data, data_view, target_view, coords_in_radians=False,
                                  kx=3, ky=3, s=0, x_tolerance=1e-3, y_tolerance=1e-3):
        t_xmin, t_ymin, t_xmax, t_ymax = target_view.bbox
        d_xmin, d_ymin, d_xmax, d_ymax = data_view.bbox
        nodata = target_view.nodata
        yx_tolerance = np.sqrt(x_tolerance**2 + y_tolerance**2)
        target_dx, target_dy = target_view.affine.a, target_view.affine.e
        data_dx, data_dy = data_view.affine.a, data_view.affine.e
        viewrows, viewcols = target_view.grid_indices(col_ascending=True,
                                                      row_ascending=True)
        rows, cols = data_view.grid_indices(col_ascending=True,
                                            row_ascending=True)
        viewrows += target_dy
        viewcols += target_dx
        rows += data_dy
        cols += data_dx
        row_bool = (rows <= t_ymax + y_tolerance) & (rows >= t_ymin - y_tolerance)
        col_bool = (cols <= t_xmax + x_tolerance) & (cols >= t_xmin - x_tolerance)
        if not coords_in_radians:
            rows = np.radians(rows) + np.pi/2
            cols = np.radians(cols) + np.pi
            viewrows = np.radians(viewrows) + np.pi/2
            viewcols = np.radians(viewcols) + np.pi
        rsbs_interpolator = (interpolate.
                            RectBivariateSpline(rows[row_bool],
                                                cols[col_bool],
                                                data[np.ix_(row_bool[::-1], col_bool)],
                                                kx=kx, ky=ky, s=s))
        xy_query = np.vstack(np.dstack(np.meshgrid(viewrows, viewcols, indexing='ij')))
        view = rsbs_interpolator.ev(xy_query[:,0], xy_query[:,1]).reshape(target_view.shape)
        return view 
Example #21
Source File: test_uvbeam.py    From pyuvdata with BSD 2-Clause "Simplified" License 5 votes vote down vote up
def test_interp_healpix_nside(cst_efield_2freq_cut, cst_efield_2freq_cut_healpix):
    efield_beam = cst_efield_2freq_cut

    efield_beam.interpolation_function = "az_za_simple"

    # test calling interp with healpix parameters directly gives same result
    min_res = np.min(
        np.array(
            [np.diff(efield_beam.axis1_array)[0], np.diff(efield_beam.axis2_array)[0]]
        )
    )
    nside_min_res = np.sqrt(3 / np.pi) * np.radians(60.0) / min_res
    nside = int(2 ** np.ceil(np.log2(nside_min_res)))

    new_efield_beam = cst_efield_2freq_cut_healpix
    assert new_efield_beam.nside == nside

    new_efield_beam.interpolation_function = "healpix_simple"

    # check error with cut sky
    with pytest.raises(
        ValueError, match="simple healpix interpolation requires full sky healpix maps."
    ):
        new_efield_beam.interp(
            az_array=efield_beam.axis1_array,
            za_array=efield_beam.axis2_array,
            az_za_grid=True,
            new_object=True,
        ) 
Example #22
Source File: double_difference.py    From seisflows with BSD 2-Clause "Simplified" License 5 votes vote down vote up
def distance(self, x1, y1, x2, y2):
        if PAR.UNITS in ['lonlat']:
            dlat = np.radians(y2-y1)
            dlon = np.radians(x2-x1)
            a = np.sin(dlat/2) * np.sin(dlat/2) + np.cos(np.radians(y1)) \
                * np.cos(np.radians(y2)) * np.sin(dlon/2) * np.sin(dlon/2)
            D = 2 * np.arctan2(np.sqrt(a), np.sqrt(1-a))
            D *= 180/np.pi
            return D

        else:
            return ((x1-x2)**2 + (y1-y2)**2)**0.5 
Example #23
Source File: occurrence_blackbox.py    From ibeis with Apache License 2.0 5 votes vote down vote up
def haversine(latlon1, latlon2):
    r"""
    Calculate the great circle distance between two points
    on the earth (specified in decimal degrees)

    Args:
        latlon1 (tuple): (lat, lon)
        latlon2 (tuple): (lat, lon)

    Returns:
        float : distance in kilometers

    References:
        en.wikipedia.org/wiki/Haversine_formula
        gis.stackexchange.com/questions/81551/matching-gps-tracks
        stackoverflow.com/questions/4913349/haversine-distance-gps-points

    Doctest:
        >>> from ibeis.algo.preproc.occurrence_blackbox import *  # NOQA
        >>> import scipy.spatial.distance as spdist
        >>> import functools
        >>> latlon1 = [-80.21895315, -158.81099213]
        >>> latlon2 = [  9.77816711,  -17.27471498]
        >>> kilometers = haversine(latlon1, latlon2)
        >>> result = ('kilometers = %s' % (kilometers,))
        >>> print(result)
        kilometers = 11930.909364189827
    """
    # convert decimal degrees to radians
    lat1, lon1 = np.radians(latlon1)
    lat2, lon2 = np.radians(latlon2)
    # haversine formula
    dlon = lon2 - lon1
    dlat = lat2 - lat1
    a = (np.sin(dlat / 2) ** 2) + np.cos(lat1) * np.cos(lat2) * (np.sin(dlon / 2) ** 2)
    c = 2 * np.arcsin(np.sqrt(a))
    # convert to kilometers
    EARTH_RADIUS_KM = 6367
    kilometers = EARTH_RADIUS_KM * c
    return kilometers 
Example #24
Source File: hypath.py    From pysplit with BSD 3-Clause "New" or "Revised" License 5 votes vote down vote up
def distance_between2pts(self, coord0, coord1, in_xy=False):
        """
        Calculate distance between two sets of coordinates.

        Parameters
        ----------
        coord0 : tuple of floats
            Coordinate pair in degrees
        coord1 : tuple of floats
            Coordinate pair in degrees
        in_xy : Boolean
            Default False.  If True, the pair is in (lon, lat)
            rather than (lat, lon)

        Returns
        -------
        distance : float
            Great circle distance in meters.

        """
        coord0 = np.radians(coord0)
        coord1 = np.radians(coord1)

        coord_order = {False: [0, 1],
                       True: [1, 0]}

        a, b = coord_order[in_xy]

        distance = (np.arccos(np.sin(coord1[a]) * np.sin(coord0[a]) +
                              np.cos(coord1[a]) * np.cos(coord0[a]) *
                              np.cos(coord0[b] - coord1[b])) * 6371) * 1000

        return distance 
Example #25
Source File: tests.py    From pycircstat with MIT License 5 votes vote down vote up
def cmtest(*args, **kwargs):
    """
    Non parametric multi-sample test for equal medians. Similar to a
    Kruskal-Wallis test for linear data.

    H0: the s populations have equal medians
    HA: the s populations have unequal medians

    :param alpha1: angles in radians
    :param alpha2: angles in radians
    :returns: p-value and test statistic of the common median test


    References: [Fisher1995]_

    """
    axis = kwargs.get('axis', None)
    if axis is None:
        alpha = list(map(np.ravel, args))
    else:
        alpha = args

    s = len(alpha)
    n = [(0 * a + 1).sum(axis=axis) for a in alpha]
    N = sum(n)

    med = descriptive.median(np.concatenate(alpha, axis=axis), axis=axis)
    if axis is not None:
        med = np.expand_dims(med, axis=axis)

    m = [np.sum(descriptive.cdiff(a, med) < 0, axis=axis) for a in alpha]
    if np.any([nn < 10 for nn in n]):
        warnings.warn('Test not applicable. Sample size in at least one group to small.')
    M = sum(m)
    P = (N ** 2. / (M * (N - M))) * sum([mm ** 2. / nn for mm, nn in zip(m, n)]) - N * M / (N - M)
    pval = stats.chi2.sf(P, df=s - 1)
    return pval, P 
Example #26
Source File: tests.py    From pycircstat with MIT License 5 votes vote down vote up
def _sample_cdf(alpha, resolution=100., axis=None):
    """

    Helper function for circ_kuipertest.
    Evaluates CDF of sample in thetas.

    :param alpha: sample (in radians)
    :param resolution: resolution at which the cdf is evaluated (default 100)
    :param axis: axis along which the cdf is computed
    :returns: points at which cdf is evaluated, cdf values

    """

    if axis is None:
        alpha = alpha.ravel()
        axis = 0
    bins = np.linspace(0, 2 * np.pi, resolution + 1)
    old_shape = alpha.shape
    alpha = alpha % (2 * np.pi)

    alpha = alpha.reshape((alpha.shape[0], int(np.prod(alpha.shape[1:])))).T
    cdf = np.array([np.histogram(a, bins=bins)[0]
                    for a in alpha]).cumsum(axis=1) / float(alpha.shape[1])
    cdf = cdf.T.reshape((len(bins) - 1,) + old_shape[1:])

    return bins[:-1], cdf 
Example #27
Source File: utils.py    From kaggle-taxi-ii with MIT License 5 votes vote down vote up
def heading(p1, p2):
    p1 = np.radians(p1)
    p2 = np.radians(p2)
    lat1, lon1 = p1[1], p1[0]
    lat2, lon2 = p2[1], p2[0]
    aa = np.sin(lon2 - lon1) * np.cos(lat2)
    bb = np.cos(lat1) * np.sin(lat2) - np.sin(lat1) * np.cos(lat2) * np.cos(lon2 - lon1)
    return np.arctan2(aa, bb) + np.pi 
Example #28
Source File: image.py    From PynPoint with GNU General Public License v3.0 5 votes vote down vote up
def rotate_coordinates(center: Tuple[float, float],
                       position: Union[Tuple[float, float], np.ndarray],
                       angle: float) -> Tuple[float, float]:
    """
    Function to rotate coordinates around the image center.

    Parameters
    ----------
    center : tuple(float, float)
        Image center (y, x).
    position : tuple(float, float)
        Position (y, x) in the image, or a 2D numpy array of positions.
    angle : float
        Angle (deg) to rotate in counterclockwise direction.

    Returns
    -------
    tuple(float, float)
        New position (y, x).
    """

    pos_y = (position[1] - center[1]) * math.sin(np.radians(angle)) + \
            (position[0] - center[0]) * math.cos(np.radians(angle))

    pos_x = (position[1] - center[1]) * math.cos(np.radians(angle)) - \
            (position[0] - center[0]) * math.sin(np.radians(angle))

    return center[0] + pos_y, center[1] + pos_x 
Example #29
Source File: image.py    From PynPoint with GNU General Public License v3.0 5 votes vote down vote up
def polar_to_cartesian(image: np.ndarray,
                       sep: float,
                       ang: float) -> Tuple[float, float]:
    """
    Function to convert polar coordinates to pixel coordinates.

    Parameters
    ----------
    image : numpy.ndarray
        Input image (2D or 3D).
    sep : float
        Separation (pix).
    ang : float
        Position angle (deg), measured counterclockwise with respect to the positive y-axis.

    Returns
    -------
    tuple(float, float)
        Cartesian coordinates (y, x). The bottom left corner of the image is (-0.5, -0.5).
    """

    center = center_subpixel(image)  # (y, x)

    x_pos = center[1] + sep * math.cos(math.radians(ang + 90))
    y_pos = center[0] + sep * math.sin(math.radians(ang + 90))

    return y_pos, x_pos 
Example #30
Source File: spore_context.py    From spore with MIT License 5 votes vote down vote up
def rotate_into(self, direction, rotation):
        """ slerp the given rotation values into the direction given
        by the brush_state
        @param direction MVector: the target direction
        @param rotation MVector: current euler rotation """

        vector_weight = self.brush_state.settings['strength']
        up_vector = om.MVector(0, 1, 0)
        local_up = up_vector.rotateBy(om.MEulerRotation(math.radians(rotation.x),
                                                        math.radians(rotation.y),
                                                        math.radians(rotation.z)))

        target_rotation = om.MQuaternion(local_up, direction, vector_weight)

        util = om.MScriptUtil()
        x_rot = np.radians(rotation.x)
        y_rot = np.radians(rotation.y)
        z_rot = np.radians(rotation.z)
        util.createFromDouble(x_rot, y_rot, z_rot)
        rotation_ptr = util.asDoublePtr()
        mat = om.MTransformationMatrix()
        mat.setRotation(rotation_ptr, om.MTransformationMatrix.kXYZ)

        mat = mat.asMatrix() * target_rotation.asMatrix()
        rotation = om.MTransformationMatrix(mat).rotation()

        return om.MVector(math.degrees(rotation.asEulerRotation().x),
                        math.degrees(rotation.asEulerRotation().y),
                        math.degrees(rotation.asEulerRotation().z))