Python numpy.fmod函数代码示例

def linearinterpolation(nstep,ndata,dt):
    """used for lienar interpolation between transportation matrices.
    returns weights alpha, beta and indices of matrices.
    Parameters
    -------
    nstep   : int Number of timesteps
    ndata   : int Number of matrices
    Returns
    -------
    alpha,beta : array 
                    coefficients for interpolation
    jalpha,jbeta : array
                    indices for interpolation
    """

    t = np.zeros(nstep,dtype=np.float_)
    for i in range(nstep):
        t[i] = np.fmod(0 + i*dt, 1.0)


    beta    = np.array(nstep,dtype=np.float_)
    alpha   = np.array(nstep,dtype=np.float_)

    w       = t * ndata+0.5
    beta    = np.float_(np.fmod(w, 1.0))
    alpha   = np.float_(1.0-beta)
    jalpha  = np.fmod(np.floor(w)+ndata-1.0,ndata).astype(int)
    jbeta   = np.fmod(np.floor(w),ndata).astype(int)

    return alpha,beta,jalpha,jbeta

def fundamental_arguments(t):

    """Compute the fundamental arguments (mean elements) of Sun and Moon.

    `t` - TDB time in Julian centuries since J2000.0, as float or NumPy array

    Outputs fundamental arguments, in radians:
          a[0] = l (mean anomaly of the Moon)
          a[1] = l' (mean anomaly of the Sun)
          a[2] = F (mean argument of the latitude of the Moon)
          a[3] = D (mean elongation of the Moon from the Sun)
          a[4] = Omega (mean longitude of the Moon's ascending node);
                 from Simon section 3.4(b.3),
                 precession = 5028.8200 arcsec/cy)

    """
    a = fa4 * t
    a += fa3
    a *= t
    a += fa2
    a *= t
    a += fa1
    a *= t
    a += fa0
    fmod(a, ASEC360, out=a)
    a *= ASEC2RAD
    if getattr(t, 'shape', ()):
        return a
    return a[:,0]

 def _visiting(self, step):
     """
     Assignement of components values based on visiting distribution.
     The way of exploring space depends on the Markov chain stepping
     """
     # It it is the first part of the markov chain
     # Changing all components at the same time
     if step < self._x.size:
         visits = np.array([self._visita() for _ in range(self._x.size)])
         visits[visits > 1.e8] = 1.e8 * self._random_state.random_sample()
         visits[visits < -1e8] = -1.e8 * self._random_state.random_sample()
         self._x = visits + self._xbackup
         a = self._x - self._lower
         b = np.fmod(a, self._xrange) + self._xrange
         self._x = np.fmod(b, self._xrange) + self._lower
         self._x[np.fabs(self._x - self._lower) < 1.e-10] += 1.e-10
     else:
         # Second part of the markov chain
         # Now change only one component at a time
         visit = self._visita()
         if visit > 1.e8:
             visit = 1.e8 * self._random_state.random_sample()
         elif visit < -1e8:
             visit = -1.e8 * self._random_state.random_sample()
         index = step - self._x.size
         self._x[index] = visit + self._xbackup[index]
         a = self._x[index] - self._lower[index]
         b = np.fmod(
             a, self._xrange[index]) + self._xrange[index]
         self._x[index] = np.fmod(
             b, self._xrange[index]) + self._lower[index]
         if np.fabs(self._x[index] - self._lower[
                 index]) < 1.e-10:
             self._x[index] += 1.e-10

def mod180deg(x):
    if x >= 0:
        retval = fmod(x+180.0, 360.0)-180.0
    else:
        retval = -(fmod(-x+180.0, 360.0)-180.0)
    assert -180.0 <= retval <= 180.0
    return retval

def fold(spike_trains, period):
    """Fold `spike_trains` by `period`."""

    # data = {key:[] for key in spike_trains.dtype.names}

    rows = []
    for i,row in spike_trains.iterrows():
        period_num = int( np.ceil(row['duration'] / period) )
        last_period = np.fmod(row['duration'], period)

        spikes = row['spikes']
        for idx in range(period_num):
            lo = idx * period
            hi = (idx+1) * period
            sec = spikes[(spikes>=lo) & (spikes<hi)]
            sec = np.fmod(sec, period)

            r = row.copy()
            r['spikes'] = sec
            r['duration'] = period

            rows.append(r)

        if last_period > 0:
            rows[-1]['duration'] = last_period


    folded_trains = pd.DataFrame(rows)
    folded_trains = folded_trains.reset_index(drop=True)

    return folded_trains

    def render(self, model, params, frame):

        # Scalar animation parameter, based on height and distance
        d = model.edgeCenters[:,2] + 0.5 * model.edgeDistances
        numpy.multiply(d, 1/self.height, d)

        # Add global offset for Z scrolling over time
        numpy.add(d, params.time * self.speed, d)

        # Add an offset that depends on which tree we're in
        numpy.add(d, numpy.choose(model.edgeTree, self.offsets), d)

        # Periodic animation, stored in our color table. Linearly interpolate.
        numpy.fmod(d, self.period, d)
        color = numpy.empty((model.numLEDs, 3))
        color[:,0] = numpy.interp(d, self.colorX, self.colorY[:,0])
        color[:,1] = numpy.interp(d, self.colorX, self.colorY[:,1])
        color[:,2] = numpy.interp(d, self.colorX, self.colorY[:,2])

        # Random flickering noise
        noise = numpy.random.rand(model.numLEDs).reshape(-1, 1)
        numpy.multiply(noise, 0.25, noise)
        numpy.add(noise, 0.75, noise)

        numpy.multiply(color, noise, color)
        numpy.add(frame, color, frame)

def _space_grid(Sm, spacing):
    '''
    Zero everywhere, except for the points at a certain grid interval.
    '''
    Ix, Iy, Iz = np.indices(Sm.shape)
    return Sm * (np.fmod(Ix, spacing) == 0) * (np.fmod(Iy, spacing) == 0) * \
           (np.fmod(Iz, spacing) == 0)

def rebin_counts_2D_indexing_new(x, y, I, xo, yo, Io):
    # use this one for 2D.
    csx = np.empty((len(x), len(y)-1))
    csx[0, :] = 0.0
    csx[1:, :] = np.cumsum(I, axis=0)

    xindices = np.interp(xo, x, np.arange(len(x), dtype='float'))
    xindices_whole = np.floor(xindices).astype(int)
    xindices_frac = np.fmod(xindices, 1.0)
    # the only way an index will match the highest bin edge is for lookups at or outside the range of the
    # source.  In this case the fractional portion will always == 0.0, so clipping this to the bin edge below
    # is safe and allows us to use the source weights unmodified
    xintegral = csx[xindices_whole, :] + I[xindices_whole.clip(max=len(x)-2), :]*xindices_frac[:, None]
    # rebinned over x
    ix = np.diff(xintegral, axis=0)

    csy = np.empty((len(xo)-1, len(y)))
    csy[:, 0] = 0.0
    csy[:, 1:] = np.cumsum(ix, axis=1)

    yindices = np.interp(yo, y, np.arange(len(y), dtype='float'))
    yindices_whole = np.floor(yindices).astype(int)
    yindices_frac = np.fmod(yindices, 1.0)
    yintegral = csy[:, yindices_whole] + ix[:, yindices_whole.clip(max=len(y)-2)]*yindices_frac[None, :]
    # rebinned over x and y
    ixy = np.diff(yintegral, axis=1)

    Io += ixy

def historyids_to_indices(historyids, historyScores): 
	hids=np.array(historyids, dtype=int)
	itr=np.fmod(hids, Global_BINWIDTH)
	sim=np.round(hids/Global_BINWIDTH)
	runlen=max(np.fmod(historyScores[:,0], Global_BINWIDTH))+1
	newi=itr+sim*runlen
	return newi.astype('int')

def compute_lscore_over_time(events, historyScores, run, outputfn):
    runlen = max(np.fmod(historyScores[:, 0], Global_BINWIDTH)) + 1
    myrun = run
    myeventids = ["total"]
    myruncounts = [1]
    hidmin = myrun * Global_BINWIDTH
    hidmax = hidmin + runlen
    tmpi = np.where((historyScores[:, 0] >= hidmin) & (historyScores[:, 0] <= hidmax))[0]
    costsarray = np.mean(historyScores[tmpi, 1:3], axis=1)
    totalscores = np.cumsum(np.exp(-1 * Global_K * costsarray))
    mydata = totalscores / totalscores[runlen - 1]
    for event in events:
        event.histories = listout_ranges(event.histRanges)
        myeventids.append(event.id)
        sys.stderr.write("working on event: %s\n" % (event.id))
        myruncounts.append(event.numsims)
        lscores = np.zeros(runlen)
        hids = np.array(event.histories)
        myhids = hids[np.where((hids <= hidmax) & (hids >= hidmin))]
        myhis = np.fmod(hids, Global_BINWIDTH)
        for i in xrange(myhis.size):
            totalscore = totalscores[myhis[i]]
            lscore = compute_likelihood_histories(myhids[: (i + 1)], historyScores, totalscore)
            lscores[myhis[i]] = lscore
        for i in xrange(max(myhis), hidmax):
            totalscore = totalscores[i]
            lscore = compute_likelihood_histories(myhids, historyScores, totalscore)
            lscores[i] = lscore
        mydata = np.vstack((mydata, lscores))
    np.savetxt(
        outputfn, mydata.T, delimiter="\t", header="\t".join(myeventids) + "\n" + "\t".join(map(str, myruncounts))
    )

 def rotate(self, angle, mask=None):
     """Rotate the grids (arena centered)
     
     Grids to be rotated can be optionally specified by bool/index array
     *mask*, otherwise population is rotated. Specified *angle* can be a
     scalar value to be applied to the population or a population- or
     mask-sized array depending on whether *mask* is specified.
     """
     rot2D = lambda psi: [[cos(psi), sin(psi)], [-sin(psi),  cos(psi)]]
     if mask is not None and type(mask) is np.ndarray:
         if mask.dtype.kind == 'b':
             mask = mask.nonzero()[0]
         if type(angle) is np.ndarray and angle.size == mask.size:
             for i,ix in enumerate(mask):
                 self._phi[ix] = np.dot(self._phi[ix], rot2D(angle[i]))
         elif type(angle) in (int, float, np.float64):
             angle = float(angle)
             self._phi[mask] = np.dot(self._phi[mask], rot2D(angle))
         else:
             raise TypeError, 'angle must be mask-sized array or float'
         self._psi[mask] = np.fmod(self._psi[mask]+angle, 2*pi)
     elif mask is None:
         if type(angle) is np.ndarray and angle.size == self.num_maps:
             for i in xrange(self.num_maps):
                 self._phi[i] = np.dot(self._phi[i], rot2D(angle[i]))
         elif type(angle) in (int, float, np.float64):
             angle = float(angle)
             self._phi = np.dot(self._phi, rot2D(angle))
         else:
             raise TypeError, 'angle must be num_maps array or float'
         self._psi = np.fmod(self._psi+angle, 2*pi)
     else:
         raise TypeError, 'mask must be bool/index array'

def plot_fft_results(freq, pow, cum, angles, x, y, freqi, powi, angi):
	plt.figure(figsize=(10,15))
	plt.subplot(5,1,1)
	plt.plot(freq, pow, 'x')
	plt.gca().set_xscale('log')
	plt.xlabel('frequency')
	plt.ylabel('power')
	plt.hlines(powi, plt.xlim()[0], plt.xlim()[1], alpha=0.3)
	plt.subplot(5,1,2)
	plt.plot(freq, cum, 'x')
	plt.gca().set_xscale('log')
	plt.ylabel('cumulative probability')
	plt.xlabel('frequency')
	plt.vlines(freqi, 0, 1, alpha=0.3)
	plt.subplot(5,1,3)
	plt.plot(cum, pow**0.5, 'x')
	plt.xlabel('cumulative probability')
	plt.ylabel('amplitude')
	plt.hlines(powi**0.5, plt.xlim()[0], plt.xlim()[1], alpha=0.3)
	plt.subplot(5,1,4)
	plt.plot(cum, numpy.fmod(angles+10*numpy.pi, 2*numpy.pi), 'x')
	plt.xlabel('cumulative probability')
	plt.ylabel('phase')
	plt.hlines(numpy.fmod(angi+10*numpy.pi, 2*numpy.pi), plt.xlim()[0], plt.xlim()[1], alpha=0.3)
	plt.subplot(5,1,5)
	plt.plot(x, y, 'x')
	plt.xlabel('time')
	plt.ylabel('value')
	xnew = numpy.linspace(x.min(), x.max(), 200)
	ynew = powi**0.5 * numpy.cos(xnew * freqi * pi * 2 + angi)
	plt.plot(xnew, ynew, '-', alpha=0.5, lw=3, 
		label="freq=%.5f\npow=%.3f\nang=%.3f" % (freqi, powi, angi))
	plt.legend(loc='best', prop=dict(size=10))

 def test_impulse_negative(self):
     """ check the transform of a negative impulse at a random place """
     win_s = 256
     i = int(floor(random()*win_s))
     impulse = -.1
     f = fft(win_s)
     timegrain = fvec(win_s)
     timegrain[0] = 0
     timegrain[i] = impulse
     fftgrain = f ( timegrain )
     #self.plot_this ( fftgrain.phas )
     assert_almost_equal ( fftgrain.norm, abs(impulse), decimal = 5 )
     if impulse < 0:
         # phase can be pi or -pi, as it is not unwrapped
         #assert_almost_equal ( abs(fftgrain.phas[1:-1]) , pi, decimal = 6 )
         assert_almost_equal ( fftgrain.phas[0], pi, decimal = 6)
         assert_almost_equal ( np.fmod(fftgrain.phas[-1], pi), 0, decimal = 6)
     else:
         #assert_equal ( fftgrain.phas[1:-1] == 0, True)
         assert_equal ( fftgrain.phas[0], 0)
         assert_almost_equal ( np.fmod(fftgrain.phas[-1], pi), 0, decimal = 6)
     # now check the resynthesis
     synthgrain = f.rdo ( fftgrain )
     #self.plot_this ( fftgrain.phas.T )
     assert_equal ( fftgrain.phas <= pi, True)
     assert_equal ( fftgrain.phas >= -pi, True)
     #self.plot_this ( synthgrain - timegrain )
     assert_almost_equal ( synthgrain, timegrain, decimal = 6 )

def analytic_phase_reset_old(phi, dx=0., dy=0.,
        a=0., l_ccw=0.2, l_cw=-1., X=1., Y=1.):

    # this is an older expression derived with several assumptions
    assert(X == 1.)
    assert(Y == 1.)
    assert(l_cw == -1.)

    r0 = iris_fixedpoint(a, l_ccw, l_cw, X, Y)
    T = iris_period(a, l_ccw, l_cw, X, Y)
    if r0 == None:
        raise RuntimeError("No limit cycle found")
    else:
        quad1or2 = np.fmod(phi, 2 * math.pi) < math.pi
        quad1or3 = np.fmod(phi, math.pi) < math.pi/2
        du = (quad1or3 * dy + (1 - quad1or3) * dx) * (1 - 2*quad1or2)
        ds = (quad1or3 * dx + (1 - quad1or3) * dy) * (
                1 - 2*quad1or2*quad1or3 - 2*(1 - quad1or2)*(1 - quad1or3))
        t = np.fmod(phi, math.pi/2)/(math.pi/2) * T/4
        Q = 1/(l_ccw * r0) * np.exp(-l_ccw * t)
        dt0 = -Q * du
        dr = np.exp(-T/4) * (X * -dt0 + np.exp(t) * ds)
        return (
                (dt0
                + -1./(l_ccw * r0)
                * 1./(1 - 1./(l_ccw * r0) * (r0/Y)**(1/l_ccw)) * dr)
                / T * 2*math.pi
                )

def _level1_xfm_no_highpass(X, h0o, h1o, ext_mode):
    """Perform level 1 of the 3d transform discarding highpass subbands.

    """
    # Check shape of input according to ext_mode. Note that shape of X is
    # double original input in each direction.
    if ext_mode == 4 and np.any(np.fmod(X.shape, 2) != 0):
        raise ValueError('Input shape should be a multiple of 2 in each direction when self.ext_mode == 4')
    elif ext_mode == 8 and np.any(np.fmod(X.shape, 4) != 0):
        raise ValueError('Input shape should be a multiple of 4 in each direction when self.ext_mode == 8')

    out = np.zeros_like(X)

    # Loop over 2nd dimension extracting 2D slice from first and 3rd dimensions
    for f in xrange(X.shape[1]):
        # extract slice
        y = X[:, f, :].T
        out[:, f, :] = colfilter(y, h0o).T

  # Loop over 3rd dimension extracting 2D slice from first and 2nd dimensions
    for f in xrange(X.shape[2]):
        y = colfilter(out[:, :, f].T, h0o).T
        out[:, :, f] = colfilter(y, h0o)

    return out