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

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

项目:prysm    作者:brandondube    | 项目源码 | 文件源码
def generate_mtf(pixel_aperture=1, azimuth=0, num_samples=128):
    ''' generates the 1D diffraction-limited MTF for a given pixel size and azimuth.

    Args:
        pixel_aperture (`float`): aperture of the pixel, in microns.  Pixel is
            assumed to be square.

        azimuth (`float`): azimuth to retrieve the MTF at, in degrees.

        num_samples (`int`): number of samples in the output array.

    Returns:
        `tuple` containing:
            `numpy.ndarray`: array of units, in cy/mm.

            `numpy.ndarray`: array of MTF values (rel. 1.0).

    Notes:
        Azimuth is not actually implemented yet.
    '''
    pitch_unit = pixel_aperture / 1e3
    normalized_frequencies = np.linspace(0, 2, num_samples)
    otf = np.sinc(normalized_frequencies)
    mtf = np.abs(otf)
    return normalized_frequencies / pitch_unit, mtf
项目:nexrad_basemap    作者:rskschrom    | 项目源码 | 文件源码
def smPhiDP(phiDP, ran):
    # smooth phiDP field and take derivative
    # calculate lanczos filter weights
    numRan = ran.shape[0]
    numK = 31
    fc = 0.015
    kt = np.linspace(-(numK-1)/2, (numK-1)/2, numK)
    w = np.sinc(2.*kt*fc)*(2.*fc)*np.sinc(kt/(numK/2))

    #smoothPhiDP = convolve1d(phiDP, w, axis=1, mode='constant', cval=-999.)
    smoothPhiDP = conv(phiDP, w)
    #smoothPhiDP = np.ma.masked_where(smoothPhiDP==-999., smoothPhiDP)

    return smoothPhiDP

# function for estimating kdp
#----------------------------------
项目:tensorpac    作者:EtienneCmb    | 项目源码 | 文件源码
def n_even_fcn(f, o, w, l):
    """Even case."""
    # Variables :
    k = np.array(range(0, int(l) + 1, 1)) + 0.5
    b = np.zeros(k.shape)

    # # Run Loop :
    for s in range(0, len(f), 2):
        m = (o[s + 1] - o[s]) / (f[s + 1] - f[s])
        b1 = o[s] - m * f[s]
        b = b + (m / (4 * np.pi * np.pi) * (np.cos(2 * np.pi * k * f[
            s + 1]) - np.cos(2 * np.pi * k * f[s])) / (
            k * k)) * abs(np.square(w[round((s + 1) / 2)]))
        b = b + (f[s + 1] * (m * f[s + 1] + b1) * np.sinc(2 * k * f[
            s + 1]) - f[s] * (m * f[s] + b1) * np.sinc(2 * k * f[s])) * abs(
            np.square(w[round((s + 1) / 2)]))

    a = (np.square(w[0])) * 4 * b
    h = 0.5 * np.concatenate((np.flipud(a), a))

    return h
项目:brainpipe    作者:EtienneCmb    | 项目源码 | 文件源码
def NevenFcn(F, M, W, L):  # N is even
    # Variables :
    k = np.array(range(0, int(L) + 1, 1)) + 0.5
    b = np.zeros(k.shape)

    # # Run Loop :
    for s in range(0, len(F), 2):
        m = (M[s + 1] - M[s]) / (F[s + 1] - F[s])
        b1 = M[s] - m * F[s]
        b = b + (m / (4 * np.pi * np.pi) * (np.cos(2 * np.pi * k * F[
            s + 1]) - np.cos(2 * np.pi * k * F[s])) / (
            k * k)) * abs(np.square(W[round((s + 1) / 2)]))
        b = b + (F[s + 1] * (m * F[s + 1] + b1) * np.sinc(2 * k * F[
          s + 1]) - F[s] * (m * F[s] + b1) * np.sinc(2 * k * F[s])) * abs(
            np.square(W[round((s + 1) / 2)]))

    a = (np.square(W[0])) * 4 * b
    h = 0.5 * np.concatenate((np.flipud(a), a))

    return h


####################################################################
# - Filt the signal :
####################################################################
项目:Least-Squared-Error-Based-FIR-Filters    作者:fourier-being    | 项目源码 | 文件源码
def lpfls(N,wp,ws,W):
    M = (N-1)/2
    nq = np.arange(0,2*M+1)
    nb = np.arange(0,M+1)
    q = (wp/np.pi)*np.sinc((wp/np.pi)*nq) - W*(ws/np.pi)*np.sinc((ws/np.pi)*nq)
    b = (wp/np.pi)*np.sinc((wp/np.pi)*nb)
    b[0] = wp/np.pi
    q[0] = wp/np.pi + W*(1-ws/np.pi) # since sin(pi*n)/pi*n = 1, not 0
    b = b.transpose()

    Q1 = ln.toeplitz(q[0:M+1])
    Q2 = ln.hankel(q[0:M+1],q[M:])
    Q = Q1+Q2

    a = ln.solve(Q,b)
    h = list(nq)
    for i in nb:
        h[i] = 0.5*a[M-i]
        h[N-1-i] = h[i]
    h[M] = 2*h[M]
    hmax = max(np.absolute(h))
    for i in nq:
        h[i] = (8191/hmax)*h[i]
    return h
项目:Least-Squared-Error-Based-FIR-Filters    作者:fourier-being    | 项目源码 | 文件源码
def bpfls(N,ws1,wp1,wp2,ws2,W):
    M = (N-1)/2
    nq = np.arange(0,2*M+1)
    nb = np.arange(0,M+1)
    q = W*np.sinc(nq) - (W*ws2/np.pi) * np.sinc(nq* (ws2/np.pi)) + (wp2/np.pi) * np.sinc(nq*(wp2/np.pi)) - (wp1/np.pi) * np.sinc(nq*(wp1/np.pi)) + (W*ws1/np.pi) * np.sinc(nq*(ws1/np.pi))
    b = (wp2/np.pi)*np.sinc((wp2/np.pi)*nb) - (wp1/np.pi)*np.sinc((wp1/np.pi)*nb)
    b[0] = wp2/np.pi - wp1/np.pi
    q[0] = W - W*ws2/np.pi + wp2/np.pi - wp1/np.pi + W*ws1/np.pi # since sin(pi*n)/pi*n = 1, not 0
    b = b.transpose()

    Q1 = ln.toeplitz(q[0:M+1])
    Q2 = ln.hankel(q[0:M+1],q[M:])
    Q = Q1+Q2

    a = ln.solve(Q,b)
    h = list(nq)
    for i in nb:
        h[i] = 0.5*a[M-i]
        h[N-1-i] = h[i]
    h[M] = 2*h[M]
    hmax = max(np.absolute(h))
    for i in nq:
        h[i] = (8191/hmax)*h[i]
    return h
项目:Least-Squared-Error-Based-FIR-Filters    作者:fourier-being    | 项目源码 | 文件源码
def hpfls(N,ws,wp,W):
    M = (N-1)/2
    nq = np.arange(0,2*M+1)
    nb = np.arange(0,M+1)
    b = 1 - (wp/np.pi)* np.sinc(nb * wp/np.pi)
    b[0] = 1- wp/np.pi
    q = 1 - (wp/np.pi)* np.sinc(nq * wp/np.pi) + W * (ws/np.pi) * np.sinc(nq * ws/np.pi) # since sin(pi*n)/pi*n = 1, not 0
    q[0] = b[0] + W* ws/np.pi
    b = b.transpose()

    Q1 = ln.toeplitz(q[0:M+1])
    Q2 = ln.hankel(q[0:M+1],q[M:])
    Q = Q1+Q2

    a = ln.solve(Q,b)
    h = list(nq)
    for i in nb:
        h[i] = 0.5*a[M-i]
        h[N-1-i] = h[i]
    h[M] = 2*h[M]
    hmax = max(np.absolute(h))
    for i in nq:
        h[i] = (8191/hmax)*h[i]
    return h
项目:zorro    作者:C-CINA    | 项目源码 | 文件源码
def lanczosSubPixKernel( subPixShift, kernelShape=3, lobes=None  ):
    """
    Generate a kernel suitable for ni.convolve to subpixally shift an image.
    """
    kernelShape = np.array( [kernelShape], dtype='int' )
    if kernelShape.ndim == 1: # make it 2-D
        kernelShape = np.array( [kernelShape[0], kernelShape[0]], dtype='int' )

    if lobes is None:
        lobes = (kernelShape[0]+1)/2

    x_range = np.arange(-kernelShape[1]/2,kernelShape[1]/2)+1.0-subPixShift[1]
    x_range = ( 2.0 / kernelShape[1] ) * x_range 
    y_range = np.arange(-kernelShape[1]/2,kernelShape[0]/2)+1.0-subPixShift[0]
    y_range = ( 2.0 /kernelShape[0] ) * y_range
    [xmesh,ymesh] = np.meshgrid( x_range, y_range )


    lanczos_filt = np.sinc(xmesh * lobes) * np.sinc(xmesh) * np.sinc(ymesh * lobes) * np.sinc(ymesh)

    lanczos_filt = lanczos_filt / np.sum(lanczos_filt) # Normalize filter output
    return lanczos_filt
项目:ASP    作者:TUIlmenauAMS    | 项目源码 | 文件源码
def sincinterp(x):
        """
        Sinc interpolation for computation of fractional transformations.
        As appears in :
        -https://github.com/audiolabs/frft/
        ----------
        Args:
            f       : (array) Complex valued input array
            a       : (float) Alpha factor
        Returns:
            ret     : (array) Real valued synthesised data
        """
        N = len(x)
        y = np.zeros(2 * N - 1, dtype=x.dtype)
        y[:2 * N:2] = x
        xint = fftconvolve( y[:2 * N], np.sinc(np.arange(-(2 * N - 3), (2 * N - 2)).T / 2),)
        return xint[2 * N - 3: -2 * N + 3]
项目:ScopeFoundry    作者:ScopeFoundry    | 项目源码 | 文件源码
def collect_pixel(self, pixel_i, frame_i,  k,j,i):
        #print pixel_i, k,j,i
        t0 = time.time()
        #px_data = np.random.rand()
        #px_data = t0 - self.prev_px
        x_hw = self.stage.settings.x_position.read_from_hardware(send_signal=False)
        y_hw = self.stage.settings.y_position.read_from_hardware(send_signal=False)

        theta = np.pi/10. * frame_i
        x = x_hw*np.cos(theta) - y_hw*np.sin(theta)
        y = x_hw*np.sin(theta) + y_hw*np.cos(theta)

        px_data = np.sinc((x-50)*0.05)**2 * np.sinc(0.05*(y-50))**2 #+ 0.05*np.random.random()
        #px_data = (x-xhw)**2 + ( y-yhw)**2
        #if px_data > 1:
        #    print('hw', x, xhw, y, yhw)
        self.display_image_map[k,j,i] = px_data
        if self.settings['save_h5']:
            self.test_data[frame_i, k,j,i] = px_data 
        time.sleep(self.settings['pixel_time'])
        #self.prev_px = t0
项目:urh    作者:jopohl    | 项目源码 | 文件源码
def design_windowed_sinc_lpf(fc, bw):
        N = Filter.get_filter_length_from_bandwidth(bw)

        # Compute sinc filter impulse response
        h = np.sinc(2 * fc * (np.arange(N) - (N - 1) / 2.))

        # We use blackman window function
        w = np.blackman(N)

        # Multiply sinc filter with window function
        h = h * w

        # Normalize to get unity gain
        h_unity = h / np.sum(h)

        return h_unity
项目:capriqorn    作者:bio-phys    | 项目源码 | 文件源码
def intensitiesFFIntraAtom(nq, dq, partNrs, nameList, ffDict):
    """ uses atomistic form factors """
    dIntensity = np.zeros((nq, 2), dtype=float)
    partInt = np.zeros((nq, len(partNrs) + 1))
    qList = np.zeros(nq)
    qList[:] = [float(i * dq) for i in range(nq)]
    dIntensity[:, 0] = qList[:]
    partInt[:, 0] = qList[:]
    # for j in range(0,len(dIntegrand)):
    #    r=dIntegrand[j,0]
    #    sinc=j0(q*r)
    k = 0
    formFacProd = np.zeros((nq, len(partNrs) + 1))
    # partNrsProd=getPartNrsProd(partNrs)
    # print "partNrsProd ", partNrsProd
    for i in range(nq):
        for k in range(1, len(partNrs) + 1):
            # print k
            formFacProd[i, k] = ff.fiveGaussian(
                ffDict[nameList[k - 1]], qList[i]) ** 2
            partInt[i, k] += partNrs[k - 1] * formFacProd[i, k]
    for i in range(nq):
        dIntensity[i, 1] = partInt[i, 1:].sum()
    return partInt, dIntensity
项目:pynufft    作者:jyhmiinlin    | 项目源码 | 文件源码
def iterate_l1(L, alpha, arg, beta, K, N, rr):
    oversample_ratio = (1.0 * K / N)

    for l1 in range(-L, L + 1):
        alf = alpha[abs(l1)] * 1.0
        if l1 < 0:
            alf = numpy.conj(alf)
    #             r1 = numpy.sinc(1.0*(arg+1.0*l1*beta)/(1.0*K/N))
        input_array = (arg + 1.0 * l1 * beta) / oversample_ratio
        r1 = dirichlet(input_array.astype(numpy.float32))
        rr = iterate_sum(rr, alf, r1)
    return rr
项目:pynufft    作者:jyhmiinlin    | 项目源码 | 文件源码
def nufft_r(om, N, J, K, alpha, beta):
    '''
    equation (30) of Fessler's paper

    '''
    def iterate_sum(rr, alf, r1):
        rr = rr + alf * r1
        return rr
    def iterate_l1(L, alpha, arg, beta, K, N, rr):
        oversample_ratio = (1.0 * K / N)
        import time
        t0=time.time()
        for l1 in range(-L, L + 1):
            alf = alpha[abs(l1)] * 1.0
    #         if l1 < 0:
    #             alf = numpy.conj(alf)
        #             r1 = numpy.sinc(1.0*(arg+1.0*l1*beta)/(1.0*K/N))
            input_array = (arg + 1.0 * l1 * beta) / oversample_ratio
            r1 = dirichlet(input_array)
            rr = iterate_sum(rr, alf, r1)
        return rr

    M = numpy.size(om)  # 1D size
    gam = 2.0 * numpy.pi / (K * 1.0)
    nufft_offset0 = nufft_offset(om, J, K)  # om/gam -  nufft_offset , [M,1]
    dk = 1.0 * om / gam - nufft_offset0  # om/gam -  nufft_offset , [M,1]
    arg = outer_sum(-numpy.arange(1, J + 1) * 1.0, dk)
    L = numpy.size(alpha) - 1
#     print('alpha',alpha)
    rr = numpy.zeros((J, M), dtype=numpy.float32)
    rr = iterate_l1(L, alpha, arg, beta, K, N, rr)
    return (rr, arg)
项目:prysm    作者:brandondube    | 项目源码 | 文件源码
def analytic_ft(self, unit_x, unit_y):
        ''' Analytic fourier transform of a pixel aperture.

        Args:
            unit_x (`numpy.ndarray`): sample points in x axis.

            unit_y (`numpy.ndarray`): sample points in y axis.

        Returns:
            `numpy.ndarray`: 2D numpy array containing the analytic fourier transform.

        '''
        xq, yq = np.meshgrid(unit_x, unit_y)
        return (sinc(xq * self.size_x / 1e3) *
                sinc(yq * self.size_y / 1e3)).astype(config.precision)
项目:untwist    作者:IoSR-Surrey    | 项目源码 | 文件源码
def make_bin_weights(self):
        erb_max = hz2erb(self.sr/2.0)
        ngrid = 1000
        erb_grid = np.arange(ngrid) * erb_max / (ngrid - 1)
        hz_grid = (np.exp(erb_grid / 9.26) - 1) / 0.00437
        resp = np.zeros((ngrid, self.n_bins))
        for b in range(self.n_bins):
            w = self.widths[b]
            r =  (2.0 * w + 1.0) / self.sr * (hz_grid - self.hz_freqs[b])
            resp[:,b] = np.square(np.sinc(r)+ 0.5 * np.sinc(r + 1.0) + 0.5 * np.sinc(r  - 1.0));
        self.weights = np.dot(linalg.pinv(resp), np.ones((ngrid,1)))
项目:smrt    作者:smrt-model    | 项目源码 | 文件源码
def autocorrelation_function(self, r):
        """compute the real space autocorrelation function for the Gaussian random field model

        """
        # compute the cut-level parameter beta
        beta = np.sqrt(2) * erfinv(2 * (1-self.frac_volume) - 1)

        # the covariance of the GRF
        acf_psi = (np.exp(-r/self.corr_length) * (1 + r / self.corr_length)
                   * np.sinc(2*r/self.repeat_distance))

        # integral discretization. henning says: the resolution 1e-2 is ad hoc, test required,
        # the integrand has a (integrable) singularity for t=1 and acf_psi = 1, so an adaptive
        # discretization seems preferable -> TODO
        dt = 1e-2
        t = np.arange(0, 1, dt)

        # the gridded integrand, via change of integration variable
        # compared to the wp-2 docu, to enable array-based computation
        t_gridded, acf_psi_gridded = np.meshgrid(t, acf_psi)
        integrand_gridded = (acf_psi_gridded / np.sqrt(1 - (t_gridded * acf_psi_gridded)**2)
                             * np.exp( - beta**2 / (1 + t_gridded * acf_psi_gridded)))

        acf = 1.0 / (2 * np.pi) * np.trapz(integrand_gridded, x=t_gridded)

        return acf

    # ft not known analytically: deligate
    # def ft_autocorrelation_function(self, k):
项目:smrt    作者:smrt-model    | 项目源码 | 文件源码
def autocorrelation_function(self, r):
        """compute the real space autocorrelation function for the Teubner Strey model

        """

        acf = self.corr_func_at_origin * np.exp(-r/self.corr_length) * np.sinc(2*r/self.repeat_distance)
        return acf
项目:tensorpac    作者:EtienneCmb    | 项目源码 | 文件源码
def n_odd_fcn(f, o, w, l):
    """Odd case."""
    # Variables :
    b0 = 0
    m = np.array(range(int(l + 1)))
    k = m[1:len(m)]
    b = np.zeros(k.shape)

    # Run Loop :
    for s in range(0, len(f), 2):
        m = (o[s + 1] - o[s]) / (f[s + 1] - f[s])
        b1 = o[s] - m * f[s]
        b0 = b0 + (b1 * (f[s + 1] - f[s]) + m / 2 * (
            f[s + 1] * f[s + 1] - f[s] * f[s])) * abs(
            np.square(w[round((s + 1) / 2)]))
        b = b + (m / (4 * np.pi * np.pi) * (
            np.cos(2 * np.pi * k * f[s + 1]) - np.cos(2 * np.pi * k * f[s])
        ) / (k * k)) * abs(np.square(w[round((s + 1) / 2)]))
        b = b + (f[s + 1] * (m * f[s + 1] + b1) * np.sinc(2 * k * f[
            s + 1]) - f[s] * (m * f[s] + b1) * np.sinc(2 * k * f[s])) * abs(
            np.square(w[round((s + 1) / 2)]))

    b = np.insert(b, 0, b0)
    a = (np.square(w[0])) * 4 * b
    a[0] = a[0] / 2
    aud = np.flipud(a[1:len(a)]) / 2
    a2 = np.insert(aud, len(aud), a[0])
    h = np.concatenate((a2, a[1:] / 2))

    return h
项目:brainpipe    作者:EtienneCmb    | 项目源码 | 文件源码
def NoddFcn(F, M, W, L):  # N is odd
    # Variables :
    b0 = 0
    m = np.array(range(int(L + 1)))
    k = m[1:len(m)]
    b = np.zeros(k.shape)

    # Run Loop :
    for s in range(0, len(F), 2):
        m = (M[s + 1] - M[s]) / (F[s + 1] - F[s])
        b1 = M[s] - m * F[s]
        b0 = b0 + (b1 * (F[s + 1] - F[s]) + m / 2 * (
            F[s + 1] * F[s + 1] - F[s] * F[s])) * abs(
            np.square(W[round((s + 1) / 2)]))
        b = b + (m / (4 * np.pi * np.pi) * (
            np.cos(2 * np.pi * k * F[s + 1]) - np.cos(2 * np.pi * k * F[s])
        ) / (k * k)) * abs(np.square(W[round((s + 1) / 2)]))
        b = b + (F[s + 1] * (m * F[s + 1] + b1) * np.sinc(2 * k * F[
          s + 1]) - F[s] * (m * F[s] + b1) * np.sinc(2 * k * F[s])) * abs(
            np.square(W[round((s + 1) / 2)]))

    b = np.insert(b, 0, b0)
    a = (np.square(W[0])) * 4 * b
    a[0] = a[0] / 2
    aud = np.flipud(a[1:len(a)]) / 2
    a2 = np.insert(aud, len(aud), a[0])
    h = np.concatenate((a2, a[1:] / 2))

    return h


# Even case
项目:Least-Squared-Error-Based-FIR-Filters    作者:fourier-being    | 项目源码 | 文件源码
def lpfls2notch(N,wp,ws,wn1,wn2,W):
    M = (N-1)/2
    nq = np.arange(0,2*M+1)
    nb = np.arange(0,M+1)
    q = (wp/np.pi)*np.sinc((wp/np.pi)*nq) - W*(ws/np.pi)*np.sinc((ws/np.pi)*nq)
    b = (wp/np.pi)*np.sinc((wp/np.pi)*nb)
    q[0] = wp/np.pi + W*(1-ws/np.pi) # since sin(pi*n)/pi*n = 1, not 0
    b = np.asmatrix(b)
    b = b.transpose()

    Q1 = ln.toeplitz(q[0:M+1])
    Q2 = ln.hankel(q[0:M+1],q[M:])
    Q = Q1+Q2

    G1 = np.cos(wn1*nb)
    G2 = np.cos(wn2*nb)
    G = np.matrix([G1,G2])

    d = np.array([0,0])
    d = np.asmatrix(d)
    d = d.transpose()

    c = np.asmatrix(ln.solve(Q,b))

    mu = ln.solve(G*ln.inv(Q)*G.transpose(),G*c - d)

    a = c - ln.solve(Q,G.transpose()*mu)
    h = np.zeros(N)
    for i in nb:
        h[i] = 0.5*a[M-i]
        h[N-1-i] = h[i]
    h[M] = 2*h[M]
    hmax = max(np.absolute(h))
    for i in nq:
        h[i] = (8191/hmax)*h[i]
    return h
项目:Least-Squared-Error-Based-FIR-Filters    作者:fourier-being    | 项目源码 | 文件源码
def lpfls1notch(N,wp,ws,wn1,W):
    M = (N-1)/2
    nq = np.arange(0,2*M+1)
    nb = np.arange(0,M+1)
    q = (wp/np.pi)*np.sinc((wp/np.pi)*nq) - W*(ws/np.pi)*np.sinc((ws/np.pi)*nq)
    b = (wp/np.pi)*np.sinc((wp/np.pi)*nb)
    q[0] = wp/np.pi + W*(1-ws/np.pi) # since sin(pi*n)/pi*n = 1, not 0
    b = np.asmatrix(b)
    b = b.transpose()

    Q1 = ln.toeplitz(q[0:M+1])
    Q2 = ln.hankel(q[0:M+1],q[M:])
    Q = Q1+Q2

    G1 = np.cos(wn1*nb)
    G = np.matrix([G1])

    d = np.array([0])
    d = np.asmatrix(d)

    c = np.asmatrix(ln.solve(Q,b))

    mu = ln.solve(G*ln.inv(Q)*G.transpose(),G*c - d)

    a = c - ln.solve(Q,G.transpose()*mu)
    h = np.zeros(N)
    for i in nb:
        h[i] = 0.5*a[M-i]
        h[N-1-i] = h[i]
    h[M] = 2*h[M]
    hmax = max(np.absolute(h))
    for i in nq:
        h[i] = (8191/hmax)*h[i]
    return h
项目:arlpy    作者:org-arl    | 项目源码 | 文件源码
def rcosfir(beta, sps, span=None):
    """Generates a raised cosine FIR filter.

    :param beta: shape of the raised cosine filter (0-1)
    :param sps: number of samples per symbol
    :param span: length of the filter in symbols (None => automatic selection)

    >>> import arlpy
    >>> rc = arlpy.comms.rcosfir(0.25, 6)
    >>> bb = arlpy.comms.modulate(arlpy.comms.random_data(100), arlpy.comms.psk())
    >>> pb = arlpy.comms.upconvert(bb, 6, 27000, 18000, rc)
    """
    if beta < 0 or beta > 1:
        raise ValueError('Beta must be between 0 and 1')
    if span is None:
        # from http://www.commsys.isy.liu.se/TSKS04/lectures/3/MichaelZoltowski_SquareRootRaisedCosine.pdf
        # since this recommendation is for root raised cosine filter, it is conservative for a raised cosine filter
        span = 33-int(44*beta) if beta < 0.68 else 4
    delay = int(span*sps/2)
    t = _np.arange(-delay, delay+1, dtype=_np.float)/sps
    denom = 1 - (2*beta*t)**2
    eps = _np.finfo(float).eps
    idx1 = _np.nonzero(_np.abs(denom) > _sqrt(eps))
    b = _np.full_like(t, beta*_sin(_pi/(2*beta))/(2*sps))
    b[idx1] = _np.sinc(t[idx1]) * _cos(_pi*beta*t[idx1])/denom[idx1] / sps
    b /= _sqrt(_np.sum(b**2))
    return b
项目:cebl    作者:idfah    | 项目源码 | 文件源码
def lanczos(x, freq, radius=3):
    l = 0.5*freq * np.sinc(0.5*freq*x) * np.sinc(x/radius)
    l[(x < -radius) | (x > radius)] = 0.0
    return l
项目:cebl    作者:idfah    | 项目源码 | 文件源码
def pollData(self):
        """Poll for new data.  This method sleeps in order to ensure
        that self.pollSize observations are generated at a realistic rate.
        """
        # figure time between polls using exponentially weighted moving average
        curTime = time.time()
        if self.pollDelay >= 0.0:
            self.pollDelay = 0.8*self.pollDelay + 0.2*(curTime - self.lastPollTime)
        else:
            self.pollDelay = 0.0
            self.lastPollTime = curTime - self.shift
        sleepTime = np.max((0.0, self.shift - self.pollDelay))
        self.lastPollTime = curTime + sleepTime
        time.sleep(sleepTime)

        with self.lock:
            # generate some random data
            data = np.random.uniform(-self.scale.value, self.scale.value,
                                     size=(self.pollSize,self.nChan))

            erpEnd = self.erpStart.value +\
                self.erpSpeed.value *\
                data.shape[0]/float(self.sampRate)

            erp = np.linspace(self.erpStart.value, erpEnd, data.shape[0])
            erp = np.repeat(erp, data.shape[1]).reshape((-1,data.shape[1]))
            erp = erp * 0.5*(np.arange(data.shape[1])+1.0)
            erp = np.sinc(erp)

            data += erp

            self.erpStart.value = erpEnd

        return data
项目:cebl    作者:idfah    | 项目源码 | 文件源码
def lanczos(n, radius=3):
    taps = np.linspace(-radius,radius, n)
    return np.sinc(taps/radius)
项目:cebl    作者:idfah    | 项目源码 | 文件源码
def sinc(n, radius=3, freq=1.0):
    taps = np.linspace(-radius,radius, n)
    return freq * np.sinc(freq*taps)
项目:cebl    作者:idfah    | 项目源码 | 文件源码
def initImpulseResponse(self, window):
        if self.bandType == 'allpass':
            self.impulseResponse = windows.kroneckerDelta(self.order+1)

        elif self.bandType == 'allstop':
            self.impulseResponse = np.zeros_like(window)

        elif self.bandType == 'lowpass':
            hightaps = self.high*self.taps
            self.impulseResponse = self.high*np.sinc(hightaps) * window

        elif self.bandType == 'highpass':
            lowtaps = self.low*self.taps
            self.impulseResponse = (-self.low*np.sinc(lowtaps) * window +
                windows.kroneckerDelta(self.order+1))

        elif self.bandType == 'bandpass':
            lowtaps = self.low*self.taps
            hightaps = self.high*self.taps
            self.impulseResponse = (self.high*np.sinc(hightaps) -
                self.low*np.sinc(lowtaps)) * window

        elif self.bandType == 'bandstop':
            lowtaps = self.low*self.taps
            hightaps = self.high*self.taps
            self.impulseResponse = ((self.high*np.sinc(hightaps) -
                self.low*np.sinc(lowtaps)) * window +
                windows.kroneckerDelta(self.order+1))

        else:
            raise Exception('Invalid bandType: ' + str(self.bandType))

        self.impulseResponse = self.impulseResponse.astype(self.dtype, copy=False)
项目:cebl    作者:idfah    | 项目源码 | 文件源码
def BandpassFilter(lowFreq, highFreq, sampRate=1.0, order=None, filtType='butter', **kwargs):
    filtType = filtType.lower()
    if filtType in ('butter', 'cheby1', 'cheby2', 'ellip', 'bessel'):
        if order is None: order = 3
        return BandpassFilterIIR(lowFreq=lowFreq, highFreq=highFreq,
                    sampRate=sampRate, order=order, filtType=filtType, **kwargs)
    elif filtType in ('lanczos', 'sinc-blackman', 'sinc-hamming', 'sinc-hann'):
        if order is None: order = 20
        return BandpassFilterFIR(lowFreq=lowFreq, highFreq=highFreq,
                    sampRate=sampRate, order=order, filtType=filtType, **kwargs)
    else:
        raise Exception('Invalid filter type: ' + str(filtType) + '.')
项目:pdc-project    作者:ealain    | 项目源码 | 文件源码
def rootRaisedCosine(t):
    bit_period = 1/BIT_FREQUENCY

    if (t== bit_period/(4*BETA)):
        return (BETA/(np.pi*np.sqrt(2*bit_period)) * \
                ((np.pi + 2)*np.sin(np.pi/(4*BETA)) + (np.pi - 2)*np.cos(np.pi/(4*BETA))))
    else:
        return (4 * BETA / np.pi / np.sqrt(bit_period) * \
            (np.cos((1 + BETA) * np.pi * t / bit_period) + \
             (1 - BETA) * np.pi / (4 * BETA) * np.sinc((1-BETA)*t/bit_period)) / \
                (1 - (4*BETA*t/bit_period)**2))
项目:pdc-project    作者:ealain    | 项目源码 | 文件源码
def rootRaisedCosine(t):
    bit_period = 1/BIT_FREQUENCY

    if (t== bit_period/(4*BETA)):
        return (BETA/(np.pi*np.sqrt(2*bit_period)) * \
                ((np.pi + 2)*np.sin(np.pi/(4*BETA)) + (np.pi - 2)*np.cos(np.pi/(4*BETA))))
    else:
        return (4 * BETA / np.pi / np.sqrt(bit_period) * \
            (np.cos((1 + BETA) * np.pi * t / bit_period) + \
             (1 - BETA) * np.pi / (4 * BETA) * np.sinc((1-BETA)*t/bit_period)) / \
                (1 - (4*BETA*t/bit_period)**2))
项目:QuantumClassicalDynamics    作者:dibondar    | 项目源码 | 文件源码
def get_CML(self, q, t):
        """
        Calculate C, M, L forming the elements of T matrix
        :param q: a shifted coordinate grid
        :param t: time
        :return: tuple C, M, L
        """
        assert q is self.x_plus or q is self.x_minus, \
            "the shifted coordinate (q) must be either x_plus or x_minus"

        # get the difference of adiabatic potential curves
        Vg_minus_Ve = (self._Vg_plus_Ve_x_plus if q is self.x_plus else self._Vg_minus_Ve_x_minus)

        Veg = self.Veg(q, t)

        D = Veg**2 + 0.25*Vg_minus_Ve**2
        np.sqrt(D, out=D)

        S = np.sinc(D * self.dt / np.pi)
        S *= self.dt

        C = D * self.dt
        np.cos(C, out=C)

        M = S * Vg_minus_Ve
        M *= 0.5

        L = S * Veg

        return C, M, L
项目:VideoTemplateCreator    作者:markdtw    | 项目源码 | 文件源码
def ef_cascade(screenpos, i, nletters):
    v = np.array([0, -1])
    d = lambda t : 1 if t < 0 else abs(np.sinc(t) / (1 + t**4))
    return lambda t: screenpos + v * 400 * d(t - 0.15 * i)
项目:ScopeFoundry    作者:ScopeFoundry    | 项目源码 | 文件源码
def collect_pixel(self, pixel_i, k,j,i):
        #print pixel_i, k,j,i
        t0 = time.time()
        #px_data = np.random.rand()
        #px_data = t0 - self.prev_px
        x0,y0 = self.pos
        x_set = self.stage.settings['x_position']
        y_set = self.stage.settings['y_position']
        x_hw = self.stage.settings.x_position.read_from_hardware(send_signal=False)
        y_hw = self.stage.settings.y_position.read_from_hardware(send_signal=False)
        if np.abs(x_hw - x0) > 1:
            self.log.debug('='*60)
            self.log.debug('pos      {} {}'.format(x0, y0))
            self.log.debug('settings {} {}'.format(x_set, y_set))
            self.log.debug('hw       {} {}'.format(x_hw, y_hw))            
            self.log.debug('settings value delta {} {}'.format(x_set-x0, y_set-y0))
            self.log.debug('read_hw  value delta {} {}'.format(x_hw-x0, y_hw-y0))
            self.log.debug('='*60)

        x = x_hw
        y = y_hw

        px_data = np.sinc((x-50)*0.05)**2 * np.sinc(0.05*(y-50))**2 #+ 0.05*np.random.random()
        #px_data = (x-xhw)**2 + ( y-yhw)**2
        #if px_data > 1:
        #    print('hw', x, xhw, y, yhw)
        self.display_image_map[k,j,i] = px_data
        if self.settings['save_h5']:
            self.test_data[k,j,i] = px_data 
        time.sleep(self.settings['pixel_time'])
        #self.prev_px = t0
项目:ThinkX    作者:AllenDowney    | 项目源码 | 文件源码
def Sinc(freq=440, amp=1.0, offset=0):
    """Makes a Sinc function.

    freq: float frequency in Hz
    amp: float amplitude, 1.0 is nominal max
    offset: float phase offset in radians

    returns: Sinusoid object
    """
    return Sinusoid(freq, amp, offset, func=np.sinc)
项目:capriqorn    作者:bio-phys    | 项目源码 | 文件源码
def intensitiesFFFaster(nq, dq, dIntegrand, keys, ffDict, dr):
    """ uses atomistic form factors """
    dIntensity = np.zeros((nq, 2), dtype=float)
    partInt = np.zeros((nq, len(dIntegrand[0])))
    nameList = [k.split(",") for k in keys]
    # print "nameList=",nameList
    qList = np.zeros(nq)
    qList[:] = [float(i * dq) for i in range(nq)]
    dIntensity[:, 0] = qList[:]
    partInt[:, 0] = qList[:]
    rList = dIntegrand[:, 0]
    # for j in range(0,len(dIntegrand)):
    #    r=dIntegrand[j,0]
    #    sinc=j0(q*r)
    formFacProd = np.zeros((nq, len(dIntegrand[0])))
    for i in range(nq):
        sincList = np.sinc(rList * qList[i] / math.pi) * dr
        for k in range(1, len(dIntegrand[0])):
            # print k
            formFacProd[i, k] = ff.fiveGaussian(ffDict[nameList[k - 1][0]], qList[i])\
                * ff.fiveGaussian(ffDict[nameList[k - 1][1]], qList[i])
            partInt[i, k] += (sincList[:] * dIntegrand[:, k]
                              ).sum() * formFacProd[i, k]
    for i in range(nq):
        dIntensity[i, 1] = partInt[i, 1:].sum()
    return partInt, dIntensity
项目:pynufft    作者:jyhmiinlin    | 项目源码 | 文件源码
def dirichlet(x):
    return numpy.sinc(x)
项目:pynufft    作者:jyhmiinlin    | 项目源码 | 文件源码
def dirichlet(x):
    return numpy.sinc(x)
项目:cebl    作者:idfah    | 项目源码 | 文件源码
def initLanczos(self, filtOrder):
        self.filtOrder = filtOrder

        if self.filtOrder % 2 != 0:
             raise Exception('Invalid filtOrder: ' + str(self.filtOrder) +
                 ' Must be an even integer.')

        radius = self.filtOrder // 2
        win = np.sinc(np.linspace(-radius, radius, self.filtOrder+1) / float(radius)) # lanczos
        #win = spsig.hamming(self.filtOrder+1) # sinc-hamming

        # this should be automated somehow XXX - idfah
        if self.filtOrder <= 6:
            cutoff = 2*0.570

        elif self.filtOrder <= 8:
            cutoff = 2*0.676

        elif self.filtOrder <= 12:
            cutoff = 2*0.781

        elif self.filtOrder <= 16:
            cutoff = 2*0.836

        elif self.filtOrder <= 32:
            cutoff = 2*0.918

        elif self.filtOrder <= 64:
            cutoff = 2*0.959

        # need to fix for multiple pool sizes XXX - idfah
        cutoff /= float(self.poolSize)

        taps = cutoff * np.linspace(-radius, radius, self.filtOrder+1, dtype=self.dtype)
        impulseResponse = cutoff * np.sinc(taps) * win

        self.filters = []
        nReadoutLayers = 1 if self.nHidden is None else 2
        for ni, no in self.layerDims[:-nReadoutLayers]:
            noEmb = no*(self.filtOrder+1) # no outs after filter embedding

            filtMat = np.zeros(noEmb*2, dtype=self.dtype)
            filtMat[noEmb-1::-no] = impulseResponse

            # filters strided for embedding
            sz = filtMat.itemsize
            filtMat = npst.as_strided(filtMat, (no,noEmb), strides=(sz,sz))[::-1].T
            self.filters.append(filtMat.copy())
项目:pdc-project    作者:ealain    | 项目源码 | 文件源码
def rrc(t):
    '''
    Input: T, evaluation point (seconds)
    Output: value of root-raised-cosine at time T
    '''
    # Delay between two bits
    bit_period = 1/BIT_FREQUENCY
    # Total amount of bits to transmit
    nb_bits = len(LIST_OF_BITS)
    # To be returned (sum of contributions)
    s = 0.0
    # Max value of rrc
    m = 4*BETA/np.pi/np.sqrt(bit_period) + (1-BETA)/np.sqrt(bit_period) + sum(abs(2*rootRaisedCosine(i*bit_period)) for i in range(1, TRUNCATION))

    if(t < - TRUNCATION * bit_period  or t >= (nb_bits + TRUNCATION) * bit_period):
        # T out of support
        r = 0.0
    else:
        # Bits that will affect function at time T
        relevant_bits = np.zeros(2*TRUNCATION+1)
        for i in range(2*TRUNCATION+1):
            j = t/bit_period + i - TRUNCATION
            j = int(j) if int(j) <= j else int(j) - 1
            if(j >= 0 and j < nb_bits):
                relevant_bits[i] = -1 if LIST_OF_BITS[j] == '0' else 1

        for i in range(2*TRUNCATION+1):
            tt = t/bit_period
            tt = t - int(tt)*bit_period if int(tt) <= tt else t - (int(tt)-1)*bit_period
            if(t == bit_period * (1 / 4 / BETA + (i - TRUNCATION))):
                # L'Hospital's rule because of potential discontinuity
                s += relevant_bits[i] * BETA / np.pi / np.sqrt(2*bit_period) * 1 / m * \
                     ((np.pi + 2) * np.sin(np.pi/4/BETA) + \
                      (np.pi - 2) * np.cos(np.pi/4/BETA))
            else:
                # General case formula
                s += relevant_bits[i] * 4*BETA/np.pi/np.sqrt(bit_period) * 1 / m * \
                     (np.cos((1 + BETA) * np.pi * ((tt / bit_period - (i-TRUNCATION)))) + \
                      (1 - BETA) * np.pi / 4 / BETA * \
                      np.sinc((1 - BETA) * (tt / bit_period - (i-TRUNCATION))))/ \
                      (1 - (4*BETA*(tt / bit_period - (i-TRUNCATION)))**2)

    return s

###   ###   ###   ###   ###   ###   ###
项目:CVFisher    作者:radiohep    | 项目源码 | 文件源码
def noise_power(self, freq):
        """Returns a function to calculate the noise PS at the given freq.
        """

        z = freq_to_z(freq)

        beam_size = self.beam_size(freq)
        A_pix = beam_size**2
        A_survey = self.f_sky * 4 * np.pi

        tau = A_pix / A_survey * units.year * self.num_year * self.beam_num

        # Calculate the comoving size of a frequency bin (at the given freq)
        d = self.proper_distance
        dxf = (d(freq_to_z(freq - self.freq_width)) - d(z))

        # Define the window function in k-space for the parallel and perpendicular directions.
        # Use a sinc function for parallel as it is the FT of a top-hat bin. This is probably a bad choice.
        def window_par(kpar):
            y = kpar * dxf / (4 * np.pi)
            return np.sinc(y) * (np.abs(y) < 1.0)

        # Azimuthally average over the X and Y window functions to produce an
        # overall k_perp window function. Do this by averaging for a set number
        # of points k values and then generating an interpolating function to
        # appeoximate the full result.
        def _int(phi, k):
            # Integrand to average over
            x = (3e2 * k * d(z)) / (freq * 2 * np.pi)

            xx = x * np.cos(phi)
            xy = x * np.sin(phi)
            return (self.window_x(xx) * self.window_y(xy))**2

        def _w_xy_average(k):
            # Full averaged window function
            return scipy.integrate.fixed_quad(_int, 0, 2 * np.pi, args=(k,), n=1024)[0]

        # Generate a log interpolated approximation
        k_val = np.linspace(0, self.kmax, 256)
        int_val = np.array([_w_xy_average(k)**0.5 for k in k_val])
        _w_perp_interp = scipy.interpolate.interp1d(k_val, np.log(int_val))

        def window_perp(kperp):
            return np.exp(_w_perp_interp(kperp))

        # Calculate the comoving volume of a single pixel (beam)
        V_pix = A_pix * d(z)**2 * dxf

        # Receiver temperature contribution to instrumental Stokes I
        T_recv_I = self.T_recv / 2**0.5

        return inv_noise_ps_21cm(T_recv_I + self.T_sky(freq), tau, V_pix,
                                 self.freq_width, window_par, window_perp)
项目:SoftSAR    作者:eduardosufan    | 项目源码 | 文件源码
def focus_multiprocessing(self, row):
        """
        Focus SAR image with TDBP algorithm. NOTE: Image must be range compressed.
        It uses local squint angle (antenna coordinate system) and distances to target to focus.

        Parameters
        ----------
        row: int.
          image row to be focus.

        Returns
        -------
        list of numpy complex.
          List containing entire focused row (numpy complex data) calculated in parallel mode.
        """

        # Light speed.
        c = 300000000.0
        # SAR bandwidth, central frequency and lambda.
        sar_B = self.param.get_float_parameter("Radar/B")
        sar_f0 = self.param.get_float_parameter("Radar/f0")
        sar_lambda = c/sar_f0

        nt_fast_time = self.simulated_image.traj.nt
        nt_slow_time = self.simulated_image.Nt

        # Partial row calculated in parallel mode focusing.
        partial_row = np.empty(self.ny, dtype=np.complex128)

        x_foc_ind = row
        for y_foc_ind in range(self.ny):
            foc_lin_ind = x_foc_ind*self.ny + y_foc_ind

            # Synthetic range compressed data (matched 2D filter).
            # Antenna Enclosure (lobe).
            doppler_amplitude_lin = (np.sinc(self.local_squint_ref_traj[foc_lin_ind, :]/self.radar_beamwidth*0.886 ))**2
            doppler_amplitude = np.tile(doppler_amplitude_lin, [self.simulated_image.Nt, 1])

            # Range amplitude: range positions in raw data of backscattered signal. These are the sincs with range 
            # migration (range compressed image).
            range_amplitude = np.sinc( sar_B*( (np.tile(self.simulated_image.t_axis_fast_time, [nt_fast_time, 1])).transpose()
                                               - np.tile(2*self.distances_ref_traj[foc_lin_ind, :]/c, [nt_slow_time, 1]) ) )

            # Limit bandwidth to threshold given by a window. Use only 3dB of antenna lobe for azimuth, limited by squint threshold.
            doppler_threshold_win = np.absolute( np.tile(self.local_squint_ref_traj[foc_lin_ind, :], [nt_slow_time, 1]) ) < self.squint_threshold
            raw_amplitude = doppler_amplitude*range_amplitude*doppler_threshold_win

            # Phase of backscattered signal (2*pi*2*r/lambda).
            raw_phase = np.exp(-1j*4*np.pi/sar_lambda*np.tile(self.distances_ref_traj[foc_lin_ind, :], [nt_slow_time, 1]))

            # Get module of raw_amplitude (for every xn, yn).
            mod_raw_amplitude = np.sum(abs(raw_amplitude)**2)
            # Repeat over x,y (slow time and fast time) to normalize.
            mod_raw_amplitude = np.tile(mod_raw_amplitude, [nt_slow_time, nt_fast_time])

            # Get raw odographer with raw_amplitude and raw_phase, i.e. with amplitude and phase information, and normalize.
            raw_to_foc = (np.conjugate(raw_phase))*raw_amplitude/mod_raw_amplitude

            partial_row[y_foc_ind] = np.sum(self.new_raw*raw_to_foc)

        return list(partial_row)
项目:SoftSAR    作者:eduardosufan    | 项目源码 | 文件源码
def generate_img(self, param):
        """
        Generate range compressed image.

        Parameters
        ----------
        param: object (ConfigurationManager instance).
          ConfigurationManager instance to read parameters from file.

        Returns
        -------
        -.
        """

        # Light speed.
        c = 300000000.0

        # Load beamwidth, bandwidth and central frequency to use locally.
        sar_bmw = (param.get_float_parameter("Radar/beamwidth")*np.pi)/180.
        sar_B = param.get_float_parameter("Radar/B")
        sar_f0 = param.get_float_parameter("Radar/f0")

        # Get angles squint and look of view with respect to the antenna coordinate system.
        #self.get_angles_antenna()
        self.local_look, self.local_squint = Utils.get_angles_antenna(self.traj, self.nom_target)

        # Set fast time axis.
        start = 2*(min(self.distances))/c - self.fast_time_pixel_margin_mono*self.radar_dt
        end = 2*(max(self.distances))/c + self.fast_time_pixel_margin_mono*self.radar_dt
        step = self.radar_dt
        self.t_axis_fast_time = np.arange(start, end, step)

        # Number of elements in fast time axis.
        self.Nt = np.size(self.t_axis_fast_time)

        self.freq_axis_fftshift = Utils.freq_axis(self.radar_dt, self.Nt, False, True)

        sar_lambda = c/sar_f0

        # Doppler amplitude (envolvente de la antena).
        doppler_amplitude = (np.sinc( (np.tile(self.local_squint, [self.Nt, 1]))/sar_bmw*(2*0.443) ))**2

        # Range amplitude: range positions in raw data of backscattered signal.
        Nd = np.size(self.distances)
        range_amplitude = np.sinc( sar_B*( (np.tile(self.t_axis_fast_time, [Nd, 1])).transpose() - np.tile(2*self.distances/c, [self.Nt, 1]) ) )

        # Signal phase received: 2*pi*2*r/lambda.
        signal_phase = np.exp(-1j*4*np.pi/sar_lambda*np.tile(self.distances, [self.Nt, 1]))

        # Generate range compressed simulated image.
        self.image = doppler_amplitude*range_amplitude*signal_phase