Python OCLArray.empty方法代码示例

# 需要导入模块: from gputools import OCLArray [as 别名]
# 或者: from gputools.OCLArray import empty [as 别名]
def _deconv_rl_gpu_conv(data_g, h_g, Niter = 10):
    """ 
    using convolve

    """
        
    #set up some gpu buffers
    u_g = OCLArray.empty(data_g.shape,np.float32)

    u_g.copy_buffer(data_g)
    
    tmp_g = OCLArray.empty(data_g.shape,np.float32)
    tmp2_g = OCLArray.empty(data_g.shape,np.float32)

    #fix this
    hflip_g = OCLArray.from_array((h_g.get()[::-1,::-1]).copy())

    for i in range(Niter):
        convolve(u_g, h_g,
                 res_g = tmp_g)


        _divide_inplace(data_g,tmp_g)

        # return data_g, tmp_g
        
        convolve(tmp_g, hflip_g,
                 res_g = tmp2_g)
        _multiply_inplace(u_g,tmp2_g)

    return u_g

# 需要导入模块: from gputools import OCLArray [as 别名]
# 或者: from gputools.OCLArray import empty [as 别名]
    def _setup_gpu(self):
        dev = get_device()
        self._queue = dev.queue
        self._ctx = dev.context
        prog = OCLProgram(absPath("kernels/bpm_3d_kernels.cl"))

        # the buffers/ images
        Nx, Ny = self.simul_xy
        Nx0, Ny0 = self.shape[:2]

        self._plan = fft_plan((Ny, Nx), **self.fftplan_kwargs)
        self._buf_plane = OCLArray.empty((Ny, Nx), np.complex64)
        self._buf_H = OCLArray.empty((Ny, Nx), np.complex64)
        self._img_xy = OCLImage.empty((Ny, Nx), dtype=np.float32, num_channels=2)

        # buffer for the weighted dn average
        self.intens_g = OCLArray.empty((1, Ny, Nx), dtype=Bpm3d._real_type)
        self.intens_dn_g = OCLArray.empty((1, Ny, Nx), dtype=Bpm3d._real_type)
        self.intens_sum_g = OCLArray.zeros((), dtype=Bpm3d._real_type)
        self.intens_dn_sum_g = OCLArray.zeros((), dtype=Bpm3d._real_type)

        # the kernels
        self._kernel_compute_propagator = prog.compute_propagator
        self._kernel_compute_propagator.set_scalar_arg_dtypes((None,)+(np.float32,)*5)
        self._kernel_compute_propagator_buf = prog.compute_propagator_buf
        self._kernel_compute_propagator_buf.set_scalar_arg_dtypes((None,)+(np.float32,)*5+(None,)*2)

        self._kernel_mult_complex = prog.mult

        self._kernel_im_to_buf_field = prog.img_to_buf_field
        self._kernel_im_to_buf_intensity = prog.img_to_buf_intensity
        self._kernel_im_to_im_intensity = prog.img_to_img_intensity
        self._kernel_buf_to_buf_field = prog.buf_to_buf_field
        self._kernel_buf_to_buf_intensity = prog.buf_to_buf_intensity

        self._kernel_mult_dn_img_float = prog.mult_dn_image
        self._kernel_mult_dn_buf_float = prog.mult_dn
        self._kernel_mult_dn_img_complex = prog.mult_dn_image_complex
        self._kernel_mult_dn_buf_complex = prog.mult_dn_complex

        self._kernel_mult_dn_img_float_local = prog.mult_dn_image_local
        self._kernel_mult_dn_buf_float_local = prog.mult_dn_local
        self._kernel_mult_dn_img_complex_local = prog.mult_dn_image_complex_local
        self._kernel_mult_dn_buf_complex_local = prog.mult_dn_complex_local

        self._kernel_reduction = OCLMultiReductionKernel(np.float32,
                                                         neutral="0", reduce_expr="a+b",
                                                         map_exprs=["a[i]", "b[i]"],
                                                         arguments="__global float *a, __global float *b")

        self._fill_propagator(self.n0)

# 需要导入模块: from gputools import OCLArray [as 别名]
# 或者: from gputools.OCLArray import empty [as 别名]
def focus_field_cylindrical(shape,units,lam = .5,NA = .3, n0=1.,
                            n_integration_steps = 100):
    """computes focus field of cylindrical lerns with given NA

    see:
    Colin J. R. Sheppard,
    Cylindrical lenses—focusing and imaging: a review

    Appl. Opt. 52, 538-545 (2013)

    return u,ex,ey,ez   with u being the intensity
    """

    p = OCLProgram(absPath("kernels/psf_cylindrical.cl"),build_options = str("-I %s -D INT_STEPS=%s"%(absPath("."),n_integration_steps)))

    
    Nx, Ny, Nz = shape
    dx, dy, dz = units

    alpha = np.arcsin(NA/n0)
    
    u_g = OCLArray.empty((Nz,Ny),np.float32)
    ex_g = OCLArray.empty((Nz,Ny),np.complex64)
    ey_g = OCLArray.empty((Nz,Ny),np.complex64)
    ez_g = OCLArray.empty((Nz,Ny),np.complex64)

    t = time.time()
    
    p.run_kernel("psf_cylindrical",u_g.shape[::-1],None,
                 ex_g.data,
                 ey_g.data,
                 ez_g.data,
                 u_g.data,
                 np.float32(-dy*(Ny-1)/2.),np.float32(dy*(Ny-1)/2.),
                 np.float32(-dz*(Nz-1)/2.),np.float32(dz*(Nz-1)/2.),
                 np.float32(lam/n0),
                 np.float32(alpha))

    u = np.array(np.repeat(u_g.get()[...,np.newaxis],Nx,axis=-1))
    ex = np.array(np.repeat(ex_g.get()[...,np.newaxis],Nx,axis=-1))
    ey = np.array(np.repeat(ey_g.get()[...,np.newaxis],Nx,axis=-1))
    ez = np.array(np.repeat(ez_g.get()[...,np.newaxis],Nx,axis=-1))

    
    print "time in secs:" , time.time()-t
    

    return u, ex, ey, ez

# 需要导入模块: from gputools import OCLArray [as 别名]
# 或者: from gputools.OCLArray import empty [as 别名]
def gpu_kuwahara(data, N=5):
    """Function to convolve an imgage with the Kuwahara filter on GPU."""
    # create numpy arrays


    if (N%2==0):       
        raise ValueError("Data has to be a (2n+1)x(2n+1) array.")

    
    data_g = OCLArray.from_array(data.astype(float32)) 
       
    res_g = OCLArray.empty((data.shape[0],data.shape[1]),float32) 
    
    prog = OCLProgram("./OpenCL/gpu_kernels/gpu_kuwahara.cl")
    
    # start kernel on gput
    prog.run_kernel("kuwahara",   # the name of the kernel in the cl file
                   data_g.shape[::-1], # global size, the number of threads e.g. (128,128,) 
                    None,   # local size, just leave it to None
                    data_g.data,res_g.data,
                    int32(N)) 
                    
    
#                    
    
    return res_g.get()

# 需要导入模块: from gputools import OCLArray [as 别名]
# 或者: from gputools.OCLArray import empty [as 别名]
def create_dn_buffer(size, units,points,
                     dn_inner = .0, rad_inner = 0,
                     dn_outer = .1, rad_outer = .4):

    Nx, Ny, Nz = size
    dx, dy, dz = units

    program = OCLProgram(absPath("kernels/bpm_3d_spheres.cl"))


    dn_g = OCLArray.empty((Nz,Ny,Nx),dtype=np.float32)

    # sort by z
    ps = np.array(points)
    ps = ps[np.argsort(ps[:,2]),:]

    Np = ps.shape[0]

    pointsBuf = OCLArray.from_array(ps.flatten().astype(np.float32))

    program.run_kernel("fill_dn",(Nx,Ny,Nz),None,dn_g.data,
                       pointsBuf.data,np.int32(Np),
                       np.float32(dx),np.float32(dy),np.float32(dz),
                       np.float32(dn_inner),np.float32(rad_inner),
                       np.float32(dn_outer),np.float32(rad_outer))


    return dn_g

# 需要导入模块: from gputools import OCLArray [as 别名]
# 或者: from gputools.OCLArray import empty [as 别名]
def _fft_convolve_gpu(data_g, h_g, res_g = None,
                      plan = None, inplace = False,
                      kernel_is_fft = False):
    """ fft convolve for gpu buffer
    """

    _complex_multiply_kernel = OCLElementwiseKernel(
        "cfloat_t *a, cfloat_t * b",
        "a[i] = cfloat_mul(b[i],a[i])","mult")


    dev = get_device()

    assert_bufs_type(np.complex64,data_g,h_g)

    if data_g.shape != h_g.shape:
        raise ValueError("data and kernel must have same size! %s vs %s "%(str(data_g.shape),str(h_g.shape)))


    if plan is None:
        plan = fft_plan(data_g.shape)

    if inplace:
        res_g = data_g
    else:
        if res_g is None:
            res_g = OCLArray.empty(data_g.shape,data_g.dtype)
            
        res_g.copy_buffer(data_g)
        
    if not kernel_is_fft:
        kern_g = OCLArray.empty(h_g.shape,h_g.dtype)
        kern_g.copy_buffer(h_g)
        fft(kern_g,inplace=True, plan = plan)
    else:
        kern_g = h_g


    fft(res_g,inplace=True, plan = plan)


    #multiply in fourier domain
    _complex_multiply_kernel(res_g,kern_g)

    fft(res_g,inplace = True, inverse = True, plan = plan)

    return res_g

# 需要导入模块: from gputools import OCLArray [as 别名]
# 或者: from gputools.OCLArray import empty [as 别名]
def _ocl_fft_gpu(plan, ocl_arr,res_arr = None, inverse = False, batch = 1):

    assert_bufs_type(np.complex64,ocl_arr)
    if res_arr is None:
        res_arr = OCLArray.empty(ocl_arr.shape,np.complex64)
    plan.execute(ocl_arr.data,res_arr.data, inverse = inverse, batch = batch)

    return res_arr

# 需要导入模块: from gputools import OCLArray [as 别名]
# 或者: from gputools.OCLArray import empty [as 别名]
def _deconv_rl_gpu_fft(data_g, h_g, Niter = 10):
    """ 
    using fft_convolve

    """


    if data_g.shape != h_g.shape:
        raise ValueError("data and h have to be same shape")

        
    #set up some gpu buffers
    u_g = OCLArray.empty(data_g.shape,np.complex64)

    u_g.copy_buffer(data_g)
    
    tmp_g = OCLArray.empty(data_g.shape,np.complex64)

    #fix this
    hflip_g = OCLArray.from_array((h_g.get()[::-1,::-1]).copy())

    plan = fft_plan(data_g.shape)

    #transform psf
    fft(h_g,inplace = True)
    fft(hflip_g,inplace = True)

    for i in range(Niter):
        print i
        fft_convolve(u_g, h_g,
                     res_g = tmp_g,
                     kernel_is_fft = True)


        _complex_divide_inplace(data_g,tmp_g)

        
        fft_convolve(tmp_g,hflip_g,
                     inplace = True,
                     kernel_is_fft = True)

        _complex_multiply_inplace(u_g,tmp_g)

    return u_g

# 需要导入模块: from gputools import OCLArray [as 别名]
# 或者: from gputools.OCLArray import empty [as 别名]
def focus_field_debye_at(x,y,z,lam, NA, n0 = 1., n_integration_steps = 200):
    """ the same as focus_field_debye but for the coordinates given in x, y, z (arrays of same shape)

        slower than focus_field_debye as it doesnt assume the coordinates to be on a grid
    """

    print absPath("kernels/psf_debye.cl")
    p = OCLProgram(absPath("kernels/psf_debye.cl"),
                   build_options = str("-I %s -D INT_STEPS=%s"%(absPath("."),n_integration_steps)))

    if np.isscalar(NA):
        NA = [0.,NA]

    alphas = np.arcsin(np.array(NA)/n0)
    assert len(alphas)%2 ==0

    assert x.shape == y.shape == z.shape
    dshape =x.shape
    N = np.prod(dshape)

    x_g = OCLArray.from_array(x.flatten().astype(np.float32))
    y_g = OCLArray.from_array(y.flatten().astype(np.float32))
    z_g = OCLArray.from_array(z.flatten().astype(np.float32))

    u_g = OCLArray.empty(N,np.float32)
    ex_g = OCLArray.empty(N,np.complex64)
    ey_g = OCLArray.empty(N,np.complex64)
    ez_g = OCLArray.empty(N,np.complex64)

    alpha_g = OCLArray.from_array(alphas.astype(np.float32))

    p.run_kernel("debye_wolf_at",(N,),None,
                 x_g.data,y_g.data,z_g.data,
                 ex_g.data,ey_g.data,ez_g.data, u_g.data,
                 np.float32(1.),np.float32(0.),
                 np.float32(lam/n0),
                 alpha_g.data, np.int32(len(alphas)))

    u = u_g.get().reshape(dshape)
    ex = ex_g.get().reshape(dshape)
    ey = ey_g.get().reshape(dshape)
    ez = ez_g.get().reshape(dshape)

    return u, ex, ey, ez

# 需要导入模块: from gputools import OCLArray [as 别名]
# 或者: from gputools.OCLArray import empty [as 别名]
def time_gpu(dshape, niter=100, fast_math=False):
    d_g = OCLArray.empty(dshape, np.complex64)
    get_device().queue.finish()
    plan = fft_plan(dshape, fast_math=fast_math)
    t = time()
    for _ in xrange(niter):
        fft(d_g, inplace=True, plan=plan)
    get_device().queue.finish()
    t = (time()-t)/niter
    print "GPU (fast_math = %s)\t%s\t\t%.2f ms"%(fast_math, dshape, 1000.*t)

# 需要导入模块: from gputools import OCLArray [as 别名]
# 或者: from gputools.OCLArray import empty [as 别名]
def _deconv_rl_np_fft(data, h, Niter = 10, 
                h_is_fftshifted = False):
    """ deconvolves data with given psf (kernel) h

    data and h have to be same shape

    
    via lucy richardson deconvolution
    """

    if data.shape != h.shape:
        raise ValueError("data and h have to be same shape")

    if not h_is_fftshifted:
        h = np.fft.fftshift(h)


    hflip = h[::-1,::-1]
        
    #set up some gpu buffers
    y_g = OCLArray.from_array(data.astype(np.complex64))
    u_g = OCLArray.from_array(data.astype(np.complex64))
    
    tmp_g = OCLArray.empty(data.shape,np.complex64)

    hf_g = OCLArray.from_array(h.astype(np.complex64))
    hflip_f_g = OCLArray.from_array(hflip.astype(np.complex64))

    # hflipped_g = OCLArray.from_array(h.astype(np.complex64))
    
    plan = fft_plan(data.shape)

    #transform psf
    fft(hf_g,inplace = True)
    fft(hflip_f_g,inplace = True)

    for i in range(Niter):
        print i
        fft_convolve(u_g, hf_g,
                     res_g = tmp_g,
                     kernel_is_fft = True)

        _complex_divide_inplace(y_g,tmp_g)

        fft_convolve(tmp_g,hflip_f_g,
                     inplace = True,
                     kernel_is_fft = True)

        _complex_multiply_inplace(u_g,tmp_g)
        

    return np.abs(u_g.get())

# 需要导入模块: from gputools import OCLArray [as 别名]
# 或者: from gputools.OCLArray import empty [as 别名]
def _ocl_fft_gpu(ocl_arr,res_arr = None,inverse = False, plan = None):

    assert_bufs_type(np.complex64,ocl_arr)

    if plan is None:
        plan = Plan(ocl_arr.shape, queue = get_device().queue)

    if res_arr is None:
        res_arr = OCLArray.empty(ocl_arr.shape,np.complex64)
        
    plan.execute(ocl_arr.data,res_arr.data, inverse = inverse)

    return res_arr

# 需要导入模块: from gputools import OCLArray [as 别名]
# 或者: from gputools.OCLArray import empty [as 别名]
def perlin2(size, units, repeat = (10.,)*2):
    wx, wy = repeat
    dx, dy = units

    prog = OCLProgram(abspath("perlin.cl"))

    d = OCLArray.empty(size[::-1],np.float32)
    prog.run_kernel("perlin2d",d.shape[::-1],None,
                    d.data,
                    np.float32(dx),np.float32(dy),
                    np.float32(wx),np.float32(wy))

    return d.get()

# 需要导入模块: from gputools import OCLArray [as 别名]
# 或者: from gputools.OCLArray import empty [as 别名]
def focus_field_cylindrical_plane(shape = (128,128),
                            units = (.1,.1),
                            z = 0.,
                            lam = .5, NA = .6, n0 = 1.,
                            ex_g = None,
                            n_integration_steps = 200):
    """
    calculates the x component of the electric field  at a given z position z for a perfect, aberration free optical system
    via the vectorial debye diffraction integral for a cylindrical lens

    see
    Colin J. R. Sheppard,
    Cylindrical lenses—focusing and imaging: a review

    Appl. Opt. 52, 538-545 (2013)


    if ex_g is a valid OCLArray it fills it and returns None
    otherwise returns ex as a numpy array


    """

    p = OCLProgram(absPath("kernels/psf_cylindrical.cl"),build_options = str("-I %s -D INT_STEPS=%s"%(absPath("."),n_integration_steps)))


    Nx, Ny = shape
    dx, dy = units

    alpha = np.arcsin(NA/n0)

    if ex_g is None:
        use_buffer = False
        ex_g = OCLArray.empty((Ny,Nx),np.complex64)
    else:
        use_buffer = True

    assert ex_g.shape[::-1] == shape


    p.run_kernel("psf_cylindrical_plane",(Nx,Ny),None,
                 ex_g.data,
                 np.float32(-dy*(Ny-1)/2.),np.float32(dy*(Ny-1)/2.),
                 np.float32(z),
                 np.float32(lam/n0),
                 np.float32(alpha))

    if not use_buffer:
        return ex_g.get()

# 需要导入模块: from gputools import OCLArray [as 别名]
# 或者: from gputools.OCLArray import empty [as 别名]
def resample_buf(data, new_shape):
    """resamples d"""

    d1_g = OCLArray.from_array(data)
    d2_g = OCLArray.empty(new_shape,data.dtype)

    if data.dtype.type == np.float32:
        im = OCLImage.empty(data.shape[::1],dtype = np.float32)
    elif data.dtype.type == np.complex64:
        im = OCLImage.empty(data.shape[::1],dtype = np.float32, num_channels=2)

    im.copy_buffer(d1_g)
    d2_g.copy_image_resampled(im)

    return d2_g.get()