Python numpy 模块,sin() 实例源码

我们从Python开源项目中,提取了以下50个代码示例,用于说明如何使用numpy.sin()

项目:MKLMM    作者:omerwe    | 项目源码 | 文件源码
def deriveKernel(self, params, i):
        self.checkParamsI(params, i)
        ell = np.exp(params[0])
        p = np.exp(params[1])

        #compute d2
        if (self.K_sq is None): d2 = sq_dist(self.X_scaled.T / ell) #precompute squared distances
        else: d2 = self.K_sq / ell**2

        #compute dp
        dp = self.dp/p

        K = np.exp(-d2 / 2.0)
        if (i==0): return d2*K*np.cos(2*np.pi*dp)
        elif (i==1): return 2*np.pi*dp*np.sin(2*np.pi*dp)*K
        else: raise Exception('invalid parameter index:' + str(i))
项目:Modeling_Preparation    作者:Yangruipis    | 项目源码 | 文件源码
def _generate_data():
    """
    ?????
    ????u(k-1) ? y(k-1)?????y(k)
    """
    # u = np.random.uniform(-1,1,200)
    # y=[]
    # former_y_value = 0
    # for i in np.arange(0,200):
    #     y.append(former_y_value)
    #     next_y_value = (29.0 / 40) * np.sin(
    #         (16.0 * u[i] + 8 * former_y_value) / (3.0 + 4.0 * (u[i] ** 2) + 4 * (former_y_value ** 2))) \
    #                    + (2.0 / 10) * u[i] + (2.0 / 10) * former_y_value
    #     former_y_value = next_y_value
    # return u,y
    u1 = np.random.uniform(-np.pi,np.pi,200)
    u2 = np.random.uniform(-1,1,200)
    y = np.zeros(200)
    for i in range(200):
        value = np.sin(u1[i]) + u2[i]
        y[i] =  value
    return u1, u2, y
项目:pointnet    作者:charlesq34    | 项目源码 | 文件源码
def rotate_point_cloud(batch_data):
    """ Randomly rotate the point clouds to augument the dataset
        rotation is per shape based along up direction
        Input:
          BxNx3 array, original batch of point clouds
        Return:
          BxNx3 array, rotated batch of point clouds
    """
    rotated_data = np.zeros(batch_data.shape, dtype=np.float32)
    for k in range(batch_data.shape[0]):
        rotation_angle = np.random.uniform() * 2 * np.pi
        cosval = np.cos(rotation_angle)
        sinval = np.sin(rotation_angle)
        rotation_matrix = np.array([[cosval, 0, sinval],
                                    [0, 1, 0],
                                    [-sinval, 0, cosval]])
        shape_pc = batch_data[k, ...]
        rotated_data[k, ...] = np.dot(shape_pc.reshape((-1, 3)), rotation_matrix)
    return rotated_data
项目:pointnet    作者:charlesq34    | 项目源码 | 文件源码
def rotate_point_cloud_by_angle(batch_data, rotation_angle):
    """ Rotate the point cloud along up direction with certain angle.
        Input:
          BxNx3 array, original batch of point clouds
        Return:
          BxNx3 array, rotated batch of point clouds
    """
    rotated_data = np.zeros(batch_data.shape, dtype=np.float32)
    for k in range(batch_data.shape[0]):
        #rotation_angle = np.random.uniform() * 2 * np.pi
        cosval = np.cos(rotation_angle)
        sinval = np.sin(rotation_angle)
        rotation_matrix = np.array([[cosval, 0, sinval],
                                    [0, 1, 0],
                                    [-sinval, 0, cosval]])
        shape_pc = batch_data[k, ...]
        rotated_data[k, ...] = np.dot(shape_pc.reshape((-1, 3)), rotation_matrix)
    return rotated_data
项目:pycma    作者:CMA-ES    | 项目源码 | 文件源码
def monotoneTFosc(f):
    """Maps [-inf,inf] to [-inf,inf] with different constants
    for positive and negative part.

    """
    if np.isscalar(f):
        if f > 0.:
            f = np.log(f) / 0.1
            f = np.exp(f + 0.49 * (np.sin(f) + np.sin(0.79 * f))) ** 0.1
        elif f < 0.:
            f = np.log(-f) / 0.1
            f = -np.exp(f + 0.49 * (np.sin(0.55 * f) + np.sin(0.31 * f))) ** 0.1
        return f
    else:
        f = np.asarray(f)
        g = f.copy()
        idx = (f > 0)
        g[idx] = np.log(f[idx]) / 0.1
        g[idx] = np.exp(g[idx] + 0.49 * (np.sin(g[idx]) + np.sin(0.79 * g[idx])))**0.1
        idx = (f < 0)
        g[idx] = np.log(-f[idx]) / 0.1
        g[idx] = -np.exp(g[idx] + 0.49 * (np.sin(0.55 * g[idx]) + np.sin(0.31 * g[idx])))**0.1
        return g
项目:psola    作者:jcreinhold    | 项目源码 | 文件源码
def test_pitch_estimation(self):
        """
        test pitch estimation algo with contrived small example
        if pitch is within 5 Hz, then say its good (for this small example,
        since the algorithm wasn't made for this type of synthesized signal)
        """
        cfg = ExperimentConfig(pitch_strength_thresh=-np.inf)
        # the next 3 variables are in Hz
        tolerance = 5
        fs = 48000
        f = 150
        # create a sine wave of f Hz freq sampled at fs Hz
        x = np.sin(2*np.pi * f/fs * np.arange(2**10))
        # estimate the pitch, it should be close to f
        p, t, s = pest.pitch_estimation(x, fs, cfg)
        self.assertTrue(np.all(np.abs(p - f) < tolerance))
项目:stcad    作者:feschmidt    | 项目源码 | 文件源码
def make_wafer(self,wafer_r,frame,label,labelloc,labelwidth):
        """
        Generate wafer with primary flat on the left. From https://coresix.com/products/wafers/ I estimated that the angle defining the wafer flat to arctan(flat/2 / radius)
        """
        angled = 18
        angle = angled*np.pi/180
        circ = cad.shapes.Circle((0,0), wafer_r, width=self.boxwidth, initial_angle=180+angled, final_angle=360+180-angled, layer=self.layer_box)
        flat = cad.core.Path([(-wafer_r*np.cos(angle),wafer_r*np.sin(angle)),(-wafer_r*np.cos(angle),-wafer_r*np.sin(angle))], width=self.boxwidth, layer=self.layer_box)

        date = time.strftime("%d/%m/%Y")
        if labelloc==(0,0):
                    labelloc=(-2e3,wafer_r-1e3)
        # The label is added 100 um on top of the main cell
        label_grid_chip = cad.shapes.LineLabel( self.name + "  " +\
                                         date,500,position=labelloc,
                                         line_width=labelwidth,
                                         layer=self.layer_label)


        if frame==True:
            self.add(circ)
            self.add(flat)
        if label==True:
            self.add(label_grid_chip)
项目:demos    作者:jnez71    | 项目源码 | 文件源码
def ani_update(arg, ii=[0]):

    i = ii[0]  # don't ask...

    if np.isclose(t_arr[i], np.around(t_arr[i], 1)):
        fig2.suptitle('Evolution (Time: {})'.format(t_arr[i]), fontsize=24)

    graphic_floor[0].set_data([-floor_lim*np.cos(incline_history[i]) + radius*np.sin(incline_history[i]), floor_lim*np.cos(incline_history[i]) + radius*np.sin(incline_history[i])], [-floor_lim*np.sin(incline_history[i])-radius*np.cos(incline_history[i]), floor_lim*np.sin(incline_history[i])-radius*np.cos(incline_history[i])])
    graphic_wheel.center = (x_history[i], y_history[i])
    graphic_ind[0].set_data([x_history[i], x_history[i] + radius*np.sin(w_history[i])],
                            [y_history[i], y_history[i] + radius*np.cos(w_history[i])])
    graphic_pend[0].set_data([x_history[i], x_history[i] - cw_to_cm[1]*np.sin(q_history[i, 2])],
                             [y_history[i], y_history[i] + cw_to_cm[1]*np.cos(q_history[i, 2])])
    graphic_dist[0].set_data([x_history[i] - cw_to_cm[1]*np.sin(q_history[i, 2]), x_history[i] - cw_to_cm[1]*np.sin(q_history[i, 2]) - pscale*p_history[i]*np.cos(q_history[i, 2])],
                             [y_history[i] + cw_to_cm[1]*np.cos(q_history[i, 2]), y_history[i] + cw_to_cm[1]*np.cos(q_history[i, 2]) - pscale*p_history[i]*np.sin(q_history[i, 2])])

    ii[0] += int(1 / (timestep * framerate))
    if ii[0] >= len(t_arr):
        print("Resetting animation!")
        ii[0] = 0

    return [graphic_floor, graphic_wheel, graphic_ind, graphic_pend, graphic_dist]

# Run animation
项目:PyGPS    作者:gregstarr    | 项目源码 | 文件源码
def solveIter(mu,e):
    """__solvIter returns an iterative solution for Ek
    Mk = Ek - e sin(Ek)
    """
    thisStart = np.asarray(mu-1.01*e)
    thisEnd = np.asarray(mu + 1.01*e)
    bestGuess = np.zeros(mu.shape)

    for i in range(5): 
        minErr = 10000*np.ones(mu.shape)
        for j in range(5):
            thisGuess = thisStart + j*(thisEnd-thisStart)/10.0
            thisErr = np.asarray(abs(mu - thisGuess + e*np.sin(thisGuess)))
            mask = thisErr<minErr
            minErr[mask] = thisErr[mask]
            bestGuess[mask] = thisGuess[mask]

        # reset for next loop
        thisRange = thisEnd - thisStart
        thisStart = bestGuess - thisRange/10.0
        thisEnd = bestGuess + thisRange/10.0

    return(bestGuess)
项目:s2g    作者:caesar0301    | 项目源码 | 文件源码
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
项目:autolab_core    作者:BerkeleyAutomation    | 项目源码 | 文件源码
def x_axis_rotation(theta):
        """Generates a 3x3 rotation matrix for a rotation of angle
        theta about the x axis.

        Parameters
        ----------
        theta : float
            amount to rotate, in radians

        Returns
        -------
        :obj:`numpy.ndarray` of float
            A random 3x3 rotation matrix.
        """
        R = np.array([[1, 0, 0,],
                      [0, np.cos(theta), -np.sin(theta)],
                      [0, np.sin(theta), np.cos(theta)]])
        return R
项目:autolab_core    作者:BerkeleyAutomation    | 项目源码 | 文件源码
def y_axis_rotation(theta):
        """Generates a 3x3 rotation matrix for a rotation of angle
        theta about the y axis.

        Parameters
        ----------
        theta : float
            amount to rotate, in radians

        Returns
        -------
        :obj:`numpy.ndarray` of float
            A random 3x3 rotation matrix.
        """
        R = np.array([[np.cos(theta), 0, np.sin(theta)],
                      [0, 1, 0],
                      [-np.sin(theta), 0, np.cos(theta)]])
        return R
项目:autolab_core    作者:BerkeleyAutomation    | 项目源码 | 文件源码
def z_axis_rotation(theta):
        """Generates a 3x3 rotation matrix for a rotation of angle
        theta about the z axis.

        Parameters
        ----------
        theta : float
            amount to rotate, in radians

        Returns
        -------
        :obj:`numpy.ndarray` of float
            A random 3x3 rotation matrix.
        """
        R = np.array([[np.cos(theta), -np.sin(theta), 0],
                      [np.sin(theta), np.cos(theta), 0],
                      [0, 0, 1]])
        return R
项目:autolab_core    作者:BerkeleyAutomation    | 项目源码 | 文件源码
def sph2cart(r, az, elev):
    """ Convert spherical to cartesian coordinates.

    Attributes
    ----------
    r : float
        radius
    az : float
        aziumth (angle about z axis)
    elev : float
        elevation from xy plane

    Returns
    -------
    float
        x-coordinate
    float
        y-coordinate
    float
        z-coordinate
    """
    x = r * np.cos(az) * np.sin(elev)
    y = r * np.sin(az) * np.sin(elev)
    z = r * np.cos(elev)
    return x, y, z
项目:mdct    作者:nils-werner    | 项目源码 | 文件源码
def mdst(x, odd=True):
    """ Calculate modified discrete sine transform of input signal in an
    inefficient pure-Python method.

    Use only for testing.

    Parameters
    ----------
    X : array_like
        The input signal
    odd : boolean, optional
        Switch to oddly stacked transform. Defaults to :code:`True`.

    Returns
    -------
    out : array_like
        The output signal

    """
    return trans(x, func=numpy.sin, odd=odd) * numpy.sqrt(2)
项目:mdct    作者:nils-werner    | 项目源码 | 文件源码
def imdst(X, odd=True):
    """ Calculate inverse modified discrete sine transform of input
    signal in an inefficient pure-Python method.

    Use only for testing.

    Parameters
    ----------
    X : array_like
        The input signal
    odd : boolean, optional
        Switch to oddly stacked transform. Defaults to :code:`True`.

    Returns
    -------
    out : array_like
        The output signal

    """
    return itrans(X, func=numpy.sin, odd=odd) * numpy.sqrt(2)
项目:mdct    作者:nils-werner    | 项目源码 | 文件源码
def cmdct(x, odd=True):
    """ Calculate complex modified discrete cosine transform of input
    inefficient pure-Python method.

    Use only for testing.

    Parameters
    ----------
    X : array_like
        The input signal
    odd : boolean, optional
        Switch to oddly stacked transform. Defaults to :code:`True`.

    Returns
    -------
    out : array_like
        The output signal

    """
    return trans(x, func=lambda x: numpy.cos(x) - 1j * numpy.sin(x), odd=odd)
项目:mdct    作者:nils-werner    | 项目源码 | 文件源码
def icmdct(X, odd=True):
    """ Calculate inverse complex modified discrete cosine transform of input
    signal in an inefficient pure-Python method.

    Use only for testing.

    Parameters
    ----------
    X : array_like
        The input signal
    odd : boolean, optional
        Switch to oddly stacked transform. Defaults to :code:`True`.

    Returns
    -------
    out : array_like
        The output signal

    """
    return itrans(X, func=lambda x: numpy.cos(x) + 1j * numpy.sin(x), odd=odd)
项目:kernel_goodness_of_fit    作者:karlnapf    | 项目源码 | 文件源码
def test_with_fake_log_prob(self):
        np.random.seed(42)


        def grad_log_prob(x):
            return -(x/2.0 + np.sin(x))*(1.0/2.0 + np.cos(x))

        def fake_log_prob(x):
            return -(x/5.0 + np.sin(x) )**2.0/2.0

        generator = mh_generator(log_density=fake_log_prob,x_start=1.0)
        tester = GaussianSteinTest(grad_log_prob,41)

        selector = SampleSelector(generator, sample_size=1000,thinning=20,tester=tester, max_iterations=5)

        data,converged = selector.points_from_stationary()

        assert converged is False
项目:kernel_goodness_of_fit    作者:karlnapf    | 项目源码 | 文件源码
def test_with_ugly(self):
        np.random.seed(42)


        def grad_log_prob(x):
            return -(x/5.0 + np.sin(x))*(1.0/5.0 + np.cos(x))

        def log_prob(x):
            return -(x/5.0 + np.sin(x) )**2.0/2.0

        generator = mh_generator(log_density=log_prob,x_start=1.0)
        tester = GaussianSteinTest(grad_log_prob,41)

        selector = SampleSelector(generator, sample_size=1000,thinning=20,tester=tester, max_iterations=5)

        data,converged = selector.points_from_stationary()

        tester = GaussianSteinTest(grad_log_prob,41)
        assert tester.compute_pvalue(data)>0.05

        assert converged is True
项目:DenoiseAverage    作者:Pella86    | 项目源码 | 文件源码
def bandpass(self, rin, sin, rout, sout):
        ''' To create a band pass two circle images are created, one inverted
        and pasted into dthe other'''

        # if radius zero dont create the inner circle
        if rin != 0:
            self.create_circle_mask(rin, sin)
        else:
            self.data = np.zeros(self.data.shape)

        # create the outer circle
        bigcircle = deepcopy(self)
        bigcircle.create_circle_mask(rout, sout)
        bigcircle.invert() 

        # sum the two pictures
        m = (self + bigcircle)

        # limit fro 0 to 1 and invert 
        m.limit(1)
        m.invert()  

        self.data = m.data
项目:sand-glyphs    作者:inconvergent    | 项目源码 | 文件源码
def random_points_in_circle(n,xx,yy,rr):
  """
  get n random points in a circle.
  """

  rnd = random(size=(n,3))
  t = TWOPI*rnd[:,0]
  u = rnd[:,1:].sum(axis=1)
  r = zeros(n,'float')
  mask = u>1.
  xmask = logical_not(mask)
  r[mask] = 2.-u[mask]
  r[xmask] = u[xmask]
  xyp = reshape(rr*r,(n,1))*column_stack( (cos(t),sin(t)) )
  dartsxy  = xyp + array([xx,yy])
  return dartsxy
项目:BISIP    作者:clberube    | 项目源码 | 文件源码
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
项目:BISIP    作者:clberube    | 项目源码 | 文件源码
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
项目:sound_field_analysis-py    作者:QULab    | 项目源码 | 文件源码
def genSphCoords():
    """ Generates cartesian (x,y,z) and spherical (theta, phi) coordinates of a sphere
    Returns
    -------
    coords : named tuple
        holds cartesian (x,y,z) and spherical (theta, phi) coordinates
    """
    coords = namedtuple('coords', ['x', 'y', 'z', 'az', 'el'])
    az = _np.linspace(0, 2 * _np.pi, 360)
    el = _np.linspace(0, _np.pi, 181)
    coords.x = _np.outer(_np.cos(az), _np.sin(el))
    coords.y = _np.outer(_np.sin(az), _np.sin(el))
    coords.z = _np.outer(_np.ones(360), _np.cos(el))

    coords.el, coords.az = _np.meshgrid(_np.linspace(0, _np.pi, 181),
                                        _np.linspace(0, 2 * _np.pi, 360))
    return coords
项目:sound_field_analysis-py    作者:QULab    | 项目源码 | 文件源码
def sph2cartMTX(vizMTX):
    """ Converts the spherical vizMTX data to named tuple contaibubg .xs/.ys/.zs

    Parameters
    ----------
    vizMTX : array_like
       [180 x 360] matrix that hold amplitude information over phi and theta

    Returns
    -------
    V : named_tuple
       Contains .xs, .ys, .zs cartesian coordinates
    """
    rs = _np.abs(vizMTX.reshape((181, -1)).T)

    coords = genSphCoords()
    V = namedtuple('V', ['xs', 'ys', 'zs'])
    V.xs = rs * _np.sin(coords.el) * _np.cos(coords.az)
    V.ys = rs * _np.sin(coords.el) * _np.sin(coords.az)
    V.zs = rs * _np.cos(coords.el)
    return V
项目:MulensModel    作者:rpoleski    | 项目源码 | 文件源码
def convert_cof_mag2mass(t0, te, u0, alpha, s, q):
    """
    function to convert from center of magnification to center of mass
    coordinates. Note that this function is for illustration only. It has
    not been tested and may have sign errors.
    """
    if s <= 1.0:
        return t0, u0
    else:
        delta = q / (1. + q) / s
        delta_u0 = delta * np.sin(alpha * np.pi / 180.)
        delta_tau = delta * np.cos(alpha * np.pi / 180.)
        t0_prime = t0 + delta_tau * te
        u0_prime = u0 + delta_u0
        return t0_prime, u0_prime

#Define model parameters in CoMAGN system
项目:MulensModel    作者:rpoleski    | 项目源码 | 文件源码
def _B_0_function(self, z):
        """
        calculate B_0(z) function defined in:

        Gould A. 1994 ApJ 421L, 71 "Proper motions of MACHOs
        http://adsabs.harvard.edu/abs/1994ApJ...421L..71G

        Yoo J. et al. 2004 ApJ 603, 139 "OGLE-2003-BLG-262: Finite-Source
        Effects from a Point-Mass Lens"
        http://adsabs.harvard.edu/abs/2004ApJ...603..139Y

        """
        out = 4. * z / np.pi
        function = lambda x: (1.-value**2*np.sin(x)**2)**.5

        for (i, value) in enumerate(z):
            if value < 1.:
                out[i] *= ellipe(value*value)
            else:
                out[i] *= integrate.quad(function, 0., np.arcsin(1./value))[0]
        return out
项目:Image-Processing-Algorithms    作者:machinelearninggod    | 项目源码 | 文件源码
def get_orientation_sector(dx,dy):
    # rotate (dx,dy) by pi/8
    rotation = np.array([[np.cos(np.pi/8), -np.sin(np.pi/8)],
                          [np.sin(np.pi/8), np.cos(np.pi/8)]])
    rotated = np.dot(rotation, np.array([[dx], [dy]]))

    if rotated[1] < 0:
        rotated[0] = -rotated[0]
        rotated[1] = -rotated[1]

    s_theta = -1
    if rotated[0] >= 0 and rotated[0] >= rotated[1]:
        s_theta = 0
    elif rotated[0] >= 0 and rotated[0] < rotated[1]:
        s_theta = 1
    elif rotated[0] < 0 and -rotated[0] < rotated[1]:
        s_theta = 2
    elif rotated[0] < 0 and -rotated[0] >= rotated[1]:
        s_theta = 3

    return s_theta
项目:atoolbox    作者:liweitianux    | 项目源码 | 文件源码
def to_radec(coords, xc=0, yc=0):
    """
    Convert the generated coordinates to (ra, dec) (unit: degree).

    xc, yc: the center coordinate (ra, dec)
    """
    results = []
    for r, theta in coords:
        # FIXME: spherical algebra should be used!!!
        dx = r * np.cos(theta*np.pi/180)
        dy = r * np.sin(theta*np.pi/180)
        x = xc + dx
        y = yc + dy
        results.append((x, y))
    if len(results) == 1:
        return results[0]
    else:
        return results
项目:pauvre    作者:conchoecia    | 项目源码 | 文件源码
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))
项目:Parallel.GAMIT    作者:demiangomez    | 项目源码 | 文件源码
def ct2lg(dX, dY, dZ, lat, lon):

    n = dX.size
    R = np.zeros((3, 3, n))

    R[0, 0, :] = -np.multiply(np.sin(np.deg2rad(lat)), np.cos(np.deg2rad(lon)))
    R[0, 1, :] = -np.multiply(np.sin(np.deg2rad(lat)), np.sin(np.deg2rad(lon)))
    R[0, 2, :] = np.cos(np.deg2rad(lat))
    R[1, 0, :] = -np.sin(np.deg2rad(lon))
    R[1, 1, :] = np.cos(np.deg2rad(lon))
    R[1, 2, :] = np.zeros((1, n))
    R[2, 0, :] = np.multiply(np.cos(np.deg2rad(lat)), np.cos(np.deg2rad(lon)))
    R[2, 1, :] = np.multiply(np.cos(np.deg2rad(lat)), np.sin(np.deg2rad(lon)))
    R[2, 2, :] = np.sin(np.deg2rad(lat))

    dxdydz = np.column_stack((np.column_stack((dX, dY)), dZ))

    RR = np.reshape(R[0, :, :], (3, n))
    dx = np.sum(np.multiply(RR, dxdydz.transpose()), axis=0)
    RR = np.reshape(R[1, :, :], (3, n))
    dy = np.sum(np.multiply(RR, dxdydz.transpose()), axis=0)
    RR = np.reshape(R[2, :, :], (3, n))
    dz = np.sum(np.multiply(RR, dxdydz.transpose()), axis=0)

    return dx, dy, dz
项目:Parallel.GAMIT    作者:demiangomez    | 项目源码 | 文件源码
def distance(self, lon1, lat1, lon2, lat2):
        """
        Calculate the great circle distance between two points
        on the earth (specified in decimal degrees)
        """

        # convert decimal degrees to radians
        lon1 = lon1*pi/180
        lat1 = lat1*pi/180
        lon2 = lon2*pi/180
        lat2 = lat2*pi/180
        # 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))
        km = 6371 * c
        return km
项目:Parallel.GAMIT    作者:demiangomez    | 项目源码 | 文件源码
def distance(self, lon1, lat1, lon2, lat2):
        """
        Calculate the great circle distance between two points
        on the earth (specified in decimal degrees)
        """

        # convert decimal degrees to radians
        lon1 = lon1*pi/180
        lat1 = lat1*pi/180
        lon2 = lon2*pi/180
        lat2 = lat2*pi/180
        # haversine formula
        dlon = lon2 - lon1
        dlat = lat2 - lat1
        a = numpy.sin(dlat/2)**2 + numpy.cos(lat1) * numpy.cos(lat2) * numpy.sin(dlon/2)**2
        c = 2 * numpy.arcsin(numpy.sqrt(a))
        km = 6371 * c
        return km
项目:Parallel.GAMIT    作者:demiangomez    | 项目源码 | 文件源码
def ct2lg(dX, dY, dZ, lat, lon):

    n = dX.size
    R = np.zeros((3, 3, n))

    R[0, 0, :] = -np.multiply(np.sin(np.deg2rad(lat)), np.cos(np.deg2rad(lon)))
    R[0, 1, :] = -np.multiply(np.sin(np.deg2rad(lat)), np.sin(np.deg2rad(lon)))
    R[0, 2, :] = np.cos(np.deg2rad(lat))
    R[1, 0, :] = -np.sin(np.deg2rad(lon))
    R[1, 1, :] = np.cos(np.deg2rad(lon))
    R[1, 2, :] = np.zeros((1, n))
    R[2, 0, :] = np.multiply(np.cos(np.deg2rad(lat)), np.cos(np.deg2rad(lon)))
    R[2, 1, :] = np.multiply(np.cos(np.deg2rad(lat)), np.sin(np.deg2rad(lon)))
    R[2, 2, :] = np.sin(np.deg2rad(lat))

    dxdydz = np.column_stack((np.column_stack((dX, dY)), dZ))

    RR = np.reshape(R[0, :, :], (3, n))
    dx = np.sum(np.multiply(RR, dxdydz.transpose()), axis=0)
    RR = np.reshape(R[1, :, :], (3, n))
    dy = np.sum(np.multiply(RR, dxdydz.transpose()), axis=0)
    RR = np.reshape(R[2, :, :], (3, n))
    dz = np.sum(np.multiply(RR, dxdydz.transpose()), axis=0)

    return dx, dy, dz
项目:Parallel.GAMIT    作者:demiangomez    | 项目源码 | 文件源码
def distance(lon1, lat1, lon2, lat2):
    """
    Calculate the great circle distance between two points
    on the earth (specified in decimal degrees)
    """

    # convert decimal degrees to radians
    lon1 = lon1*pi/180
    lat1 = lat1*pi/180
    lon2 = lon2*pi/180
    lat2 = lat2*pi/180
    # 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))
    km = 6371 * c
    return km
项目:Parallel.GAMIT    作者:demiangomez    | 项目源码 | 文件源码
def todictionary(self, time_series=False):
        # convert the ETM adjustment into a dirtionary
        # optionally, output the whole time series as well

        # start with the parameters
        etm = dict()
        etm['Linear'] = {'tref': self.Linear.tref, 'params': self.Linear.values.tolist()}
        etm['Jumps'] = [{'type':jump.type, 'year': jump.year, 'a': jump.a.tolist(), 'b': jump.b.tolist(), 'T': jump.T} for jump in self.Jumps.table]
        etm['Periodic'] = {'frequencies': self.Periodic.frequencies, 'sin': self.Periodic.sin.tolist(), 'cos': self.Periodic.cos.tolist()}

        if time_series:
            ts = dict()
            ts['t'] = self.ppp_soln.t.tolist()
            ts['x'] = self.ppp_soln.x.tolist()
            ts['y'] = self.ppp_soln.y.tolist()
            ts['z'] = self.ppp_soln.z.tolist()

            etm['time_series'] = ts

        return etm
项目:Parallel.GAMIT    作者:demiangomez    | 项目源码 | 文件源码
def ct2lg(self, dX, dY, dZ, lat, lon):

        n = dX.size
        R = numpy.zeros((3, 3, n))

        R[0, 0, :] = -numpy.multiply(numpy.sin(numpy.deg2rad(lat)), numpy.cos(numpy.deg2rad(lon)))
        R[0, 1, :] = -numpy.multiply(numpy.sin(numpy.deg2rad(lat)), numpy.sin(numpy.deg2rad(lon)))
        R[0, 2, :] = numpy.cos(numpy.deg2rad(lat))
        R[1, 0, :] = -numpy.sin(numpy.deg2rad(lon))
        R[1, 1, :] = numpy.cos(numpy.deg2rad(lon))
        R[1, 2, :] = numpy.zeros((1, n))
        R[2, 0, :] = numpy.multiply(numpy.cos(numpy.deg2rad(lat)), numpy.cos(numpy.deg2rad(lon)))
        R[2, 1, :] = numpy.multiply(numpy.cos(numpy.deg2rad(lat)), numpy.sin(numpy.deg2rad(lon)))
        R[2, 2, :] = numpy.sin(numpy.deg2rad(lat))

        dxdydz = numpy.column_stack((numpy.column_stack((dX, dY)), dZ))

        RR = numpy.reshape(R[0, :, :], (3, n))
        dx = numpy.sum(numpy.multiply(RR, dxdydz.transpose()), axis=0)
        RR = numpy.reshape(R[1, :, :], (3, n))
        dy = numpy.sum(numpy.multiply(RR, dxdydz.transpose()), axis=0)
        RR = numpy.reshape(R[2, :, :], (3, n))
        dz = numpy.sum(numpy.multiply(RR, dxdydz.transpose()), axis=0)

        return dx, dy, dz
项目:Parallel.GAMIT    作者:demiangomez    | 项目源码 | 文件源码
def distance(self, lon1, lat1, lon2, lat2):
        """
        Calculate the great circle distance between two points
        on the earth (specified in decimal degrees)
        """

        # convert decimal degrees to radians
        lon1 = lon1*pi/180
        lat1 = lat1*pi/180
        lon2 = lon2*pi/180
        lat2 = lat2*pi/180
        # haversine formula
        dlon = lon2 - lon1
        dlat = lat2 - lat1
        a = numpy.sin(dlat/2)**2 + numpy.cos(lat1) * numpy.cos(lat2) * numpy.sin(dlon/2)**2
        c = 2 * numpy.arcsin(numpy.sqrt(a))
        km = 6371 * c
        return km
项目:Parallel.GAMIT    作者:demiangomez    | 项目源码 | 文件源码
def distance(self, lon1, lat1, lon2, lat2):
        """
        Calculate the great circle distance between two points
        on the earth (specified in decimal degrees)
        """

        # convert decimal degrees to radians
        lon1 = lon1*pi/180
        lat1 = lat1*pi/180
        lon2 = lon2*pi/180
        lat2 = lat2*pi/180
        # haversine formula
        dlon = lon2 - lon1
        dlat = lat2 - lat1
        a = numpy.sin(dlat/2)**2 + numpy.cos(lat1) * numpy.cos(lat2) * numpy.sin(dlon/2)**2
        c = 2 * numpy.arcsin(numpy.sqrt(a))
        km = 6371 * c
        return km
项目:Parallel.GAMIT    作者:demiangomez    | 项目源码 | 文件源码
def ct2lg(self, dX, dY, dZ, lat, lon):

        n = dX.size
        R = numpy.zeros((3, 3, n))

        R[0, 0, :] = -numpy.multiply(numpy.sin(numpy.deg2rad(lat)), numpy.cos(numpy.deg2rad(lon)))
        R[0, 1, :] = -numpy.multiply(numpy.sin(numpy.deg2rad(lat)), numpy.sin(numpy.deg2rad(lon)))
        R[0, 2, :] = numpy.cos(numpy.deg2rad(lat))
        R[1, 0, :] = -numpy.sin(numpy.deg2rad(lon))
        R[1, 1, :] = numpy.cos(numpy.deg2rad(lon))
        R[1, 2, :] = numpy.zeros((1, n))
        R[2, 0, :] = numpy.multiply(numpy.cos(numpy.deg2rad(lat)), numpy.cos(numpy.deg2rad(lon)))
        R[2, 1, :] = numpy.multiply(numpy.cos(numpy.deg2rad(lat)), numpy.sin(numpy.deg2rad(lon)))
        R[2, 2, :] = numpy.sin(numpy.deg2rad(lat))

        dxdydz = numpy.column_stack((numpy.column_stack((dX, dY)), dZ))

        RR = numpy.reshape(R[0, :, :], (3, n))
        dx = numpy.sum(numpy.multiply(RR, dxdydz.transpose()), axis=0)
        RR = numpy.reshape(R[1, :, :], (3, n))
        dy = numpy.sum(numpy.multiply(RR, dxdydz.transpose()), axis=0)
        RR = numpy.reshape(R[2, :, :], (3, n))
        dz = numpy.sum(numpy.multiply(RR, dxdydz.transpose()), axis=0)

        return dx, dy, dz
项目:esys-pbi    作者:fsxfreak    | 项目源码 | 文件源码
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]])
项目:em_examples    作者:geoscixyz    | 项目源码 | 文件源码
def plotPlaceTxRxSphereXY(Ax,xtx,ytx,xrx,yrx,x0,y0,a):


    Xlim = Ax.get_xlim()
    Ylim = Ax.get_ylim()

    FS = 20

    Ax.scatter(xtx,ytx,s=100,color='k')
    Ax.text(xtx-0.75,ytx+1.5,'$\mathbf{Tx}$',fontsize=FS+6)
    Ax.scatter(xrx,yrx,s=100,color='k')
    Ax.text(xrx-0.75,yrx-4,'$\mathbf{Rx}$',fontsize=FS+6)

    xs = x0 + a*np.cos(np.linspace(0,2*np.pi,41))
    ys = y0 + a*np.sin(np.linspace(0,2*np.pi,41))

    Ax.plot(xs,ys,ls=':',color='k',linewidth=3)

    Ax.set_xbound(Xlim)
    Ax.set_ybound(Ylim)

    return Ax
项目:em_examples    作者:geoscixyz    | 项目源码 | 文件源码
def calc_IndCurrent_cos_range(self,f,t):
        """Induced current over a range of times"""

        Bpx = self.Bpx
        Bpz = self.Bpz
        a2  = self.a2
        azm = np.pi*self.azm/180.
        R   = self.R
        L   = self.L

        w = 2*np.pi*f

        Ax = np.pi*a2**2*np.sin(azm)
        Az = np.pi*a2**2*np.cos(azm)

        Phi = (Ax*Bpx + Az*Bpz)
        phi = np.arctan(R/(w*L))-np.pi  # This is the phase and not phase lag
        Is  = -(w*Phi/(R*np.sin(phi) + w*L*np.cos(phi)))*np.cos(w*t + phi)
        Ire = -(w*Phi/(R*np.sin(phi) + w*L*np.cos(phi)))*np.cos(w*t)*np.cos(phi)
        Iim =  (w*Phi/(R*np.sin(phi) + w*L*np.cos(phi)))*np.sin(w*t)*np.sin(phi)

        return Ire,Iim,Is,phi
项目:pyku    作者:dubvulture    | 项目源码 | 文件源码
def is_grid(self, grid, image):
        """
        Checks the "gridness" by analyzing the results of a hough transform.
        :param grid: binary image
        :return: wheter the object in the image might be a grid or not
        """
        #   - Distance resolution = 1 pixel
        #   - Angle resolution = 1° degree for high line density
        #   - Threshold = 144 hough intersections
        #        8px digit + 3*2px white + 2*1px border = 16px per cell
        #           => 144x144 grid
        #        144 - minimum number of points on the same line
        #       (but due to imperfections in the binarized image it's highly
        #        improbable to detect a 144x144 grid)
        lines = cv2.HoughLines(grid, 1, np.pi / 180, 144)

        if lines is not None and np.size(lines) >= 20:
            lines = lines.reshape((lines.size / 2), 2)
            # theta in [0, pi] (theta > pi => rho < 0)
            # normalise theta in [-pi, pi] and negatives rho
            lines[lines[:, 0] < 0, 1] -= np.pi
            lines[lines[:, 0] < 0, 0] *= -1

            criteria = (cv2.TERM_CRITERIA_EPS, 0, 0.01)
            # split lines into 2 groups to check whether they're perpendicular
            if cv2.__version__[0] == '2':
                density, clmap, centers = cv2.kmeans(
                    lines[:, 1], 2, criteria, 5, cv2.KMEANS_RANDOM_CENTERS)
            else:
                density, clmap, centers = cv2.kmeans(
                    lines[:, 1], 2, None, criteria,
                    5, cv2.KMEANS_RANDOM_CENTERS)

            if self.debug:
                self.save_hough(lines, clmap)

            # Overall variance from respective centers
            var = density / np.size(clmap)
            sin = abs(np.sin(centers[0] - centers[1]))
            # It is probably a grid only if:
            #   - centroids difference is almost a 90° angle (+-15° limit)
            #   - variance is less than 5° (keeping in mind surface distortions)
            return sin > 0.99 and var <= (5*np.pi / 180) ** 2
        else:
            return False
项目:pyku    作者:dubvulture    | 项目源码 | 文件源码
def save_hough(self, lines, clmap):
        """
        :param lines: (rho, theta) pairs
        :param clmap: clusters assigned to lines
        :return: None
        """
        height, width = self.image.shape
        ratio = 600. * (self.step+1) / min(height, width)
        temp = cv2.resize(self.image, None, fx=ratio, fy=ratio,
                          interpolation=cv2.INTER_CUBIC)
        temp = cv2.cvtColor(temp, cv2.COLOR_GRAY2BGR)
        colors = [(0, 127, 255), (255, 0, 127)]

        for i in range(0, np.size(lines) / 2):
            rho = lines[i, 0]
            theta = lines[i, 1]
            color = colors[clmap[i, 0]]
            if theta < np.pi / 4 or theta > 3 * np.pi / 4:
                pt1 = (rho / np.cos(theta), 0)
                pt2 = (rho - height * np.sin(theta) / np.cos(theta), height)
            else:
                pt1 = (0, rho / np.sin(theta))
                pt2 = (width, (rho - width * np.cos(theta)) / np.sin(theta))
            pt1 = (int(pt1[0]), int(pt1[1]))
            pt2 = (int(pt2[0]), int(pt2[1]))
            cv2.line(temp, pt1, pt2, color, 5)

        self.save2image(temp)
项目:pyku    作者:dubvulture    | 项目源码 | 文件源码
def is_grid(self, grid, image):
        """
        Checks the "gridness" by analyzing the results of a hough transform.
        :param grid: binary image
        :return: wheter the object in the image might be a grid or not
        """
        #   - Distance resolution = 1 pixel
        #   - Angle resolution = 1° degree for high line density
        #   - Threshold = 144 hough intersections
        #        8px digit + 3*2px white + 2*1px border = 16px per cell
        #           => 144x144 grid
        #        144 - minimum number of points on the same line
        #       (but due to imperfections in the binarized image it's highly
        #        improbable to detect a 144x144 grid)

        lines = cv2.HoughLines(grid, 1, np.pi / 180, 144)

        if lines is not None and np.size(lines) >= 20:
            lines = lines.reshape((lines.size/2), 2)
            # theta in [0, pi] (theta > pi => rho < 0)
            # normalise theta in [-pi, pi] and negatives rho
            lines[lines[:, 0] < 0, 1] -= np.pi
            lines[lines[:, 0] < 0, 0] *= -1

            criteria = (cv2.TERM_CRITERIA_EPS, 0, 0.01)
            # split lines into 2 groups to check whether they're perpendicular
            if cv2.__version__[0] == '2':
                density, clmap, centers = cv2.kmeans(
                    lines[:, 1], 2, criteria,
                    5, cv2.KMEANS_RANDOM_CENTERS)
            else:
                density, clmap, centers = cv2.kmeans(
                    lines[:, 1], 2, None, criteria,
                    5, cv2.KMEANS_RANDOM_CENTERS)

            # Overall variance from respective centers
            var = density / np.size(clmap)
            sin = abs(np.sin(centers[0] - centers[1]))
            # It is probably a grid only if:
            #   - centroids difference is almost a 90° angle (+-15° limit)
            #   - variance is less than 5° (keeping in mind surface distortions)
            return sin > 0.99 and var <= (5*np.pi / 180) ** 2
        else:
            return False
项目:TDOSE    作者:kasperschmidt    | 项目源码 | 文件源码
def build_2D_cov_matrix(sigmax,sigmay,angle,verbose=True):
    """
    Build a covariance matrix for a 2D multivariate Gaussian

    --- INPUT ---
    sigmax          Standard deviation of the x-compoent of the multivariate Gaussian
    sigmay          Standard deviation of the y-compoent of the multivariate Gaussian
    angle           Angle to rotate matrix by in degrees (clockwise) to populate covariance cross terms
    verbose         Toggle verbosity
    --- EXAMPLE OF USE ---
    import tdose_utilities as tu
    covmatrix = tu.build_2D_cov_matrix(3,1,35)

    """
    if verbose: print ' - Build 2D covariance matrix with varinaces (x,y)=('+str(sigmax)+','+str(sigmay)+\
                      ') and then rotated '+str(angle)+' degrees'
    cov_orig      = np.zeros([2,2])
    cov_orig[0,0] = sigmay**2.0
    cov_orig[1,1] = sigmax**2.0

    angle_rad     = (180.0-angle) * np.pi/180.0 # The (90-angle) makes sure the same convention as DS9 is used
    c, s          = np.cos(angle_rad), np.sin(angle_rad)
    rotmatrix     = np.matrix([[c, -s], [s, c]])

    cov_rot       = np.dot(np.dot(rotmatrix,cov_orig),np.transpose(rotmatrix))  # performing rot * cov * rot^T

    return cov_rot
# = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = =
项目:pyballd    作者:Yurlungur    | 项目源码 | 文件源码
def residual(r,theta,u,d):
    out = np.empty_like(u)
    out[0] = (2*np.sin(theta)*r*d(u[0],1,0)
              + r*r*np.sin(theta)*d(u[0],2,0)
              + np.cos(theta)*d(u[0],0,1)
              + np.sin(theta)*d(u[1],0,2))
    out[1] = (2*np.sin(theta)*r*d(u[1],1,0)
              + r*r*np.sin(theta)*d(u[1],2,0)
              + np.cos(theta)*d(u[1],0,1)
              + np.sin(theta)*d(u[1],0,2))
    return out
项目:pyballd    作者:Yurlungur    | 项目源码 | 文件源码
def residual(r,theta,u,d):
    u = u[0]
    out = (2*np.sin(theta)*r*d(u,1,0)
           + r*r*np.sin(theta)*d(u,2,0)
           + np.cos(theta)*d(u,0,1)
           + np.sin(theta)*d(u,0,2))
    out = out.reshape(tuple([1]) + out.shape)
    return out