Python astra.create_projector函数代码示例

def recon_sirt_fbp(im, angles, iter, temppath ):
	"""Reconstruct a sinogram with the SIRT-FBP algorithm (Pelt, 2015).

	Parameters
	----------
	im : array_like
		Sinogram image data as numpy array.

    iter : int
		Number of iterations to be used for the computation of SIRT filter.

	angles : double
		Value in radians representing the number of angles of the input sinogram.
	
	"""
	# Create ASTRA geometries:
	vol_geom = astra.create_vol_geom(im.shape[1] , im.shape[1])
	proj_geom = astra.create_proj_geom('parallel', 1.0, im.shape[1], linspace(0,angles,im.shape[0],False))
	proj = astra.create_projector('cuda', proj_geom, vol_geom)
	p = astra.OpTomo(proj)

	# Register plugin with ASTRA
	astra.plugin.register(sirtfbp.plugin)

	# Create the ASTRA projector
	im_rec = p.reconstruct('SIRT-FBP',im, iter, extraOptions={'filter_dir':temppath})
	
	return im_rec.astype(float32)

def astra_rec_cpu(tomo, center, recon, theta, vol_geom, niter, proj_type, opts):
    # Lazy import ASTRA
    import astra as astra_mod
    nslices, nang, ndet = tomo.shape
    cfg = astra_mod.astra_dict(opts['method'])
    if 'extra_options' in opts:
        cfg['option'] = opts['extra_options']
    proj_geom = astra_mod.create_proj_geom('parallel', 1.0, ndet, theta.astype(np.float64))
    pid = astra_mod.create_projector(proj_type, proj_geom, vol_geom)
    sino = np.zeros((nang, ndet), dtype=np.float32)
    sid = astra_mod.data2d.link('-sino', proj_geom, sino)
    cfg['ProjectorId'] = pid
    cfg['ProjectionDataId'] = sid
    for i in range(nslices):
        shft = int(np.round(ndet / 2. - center[i]))
        if not shft == 0:
            sino[:] = np.roll(tomo[i], shft)
            l = shft
            r = ndet + shft
            if l < 0:
                l = 0
            if r > ndet:
                r = ndet
            sino[:, :l] = 0
            sino[:,  r:] = 0
        else:
            sino[:] = tomo[i]
        vid = astra_mod.data2d.link('-vol', vol_geom, recon[i])
        cfg['ReconstructionDataId'] = vid
        alg_id = astra_mod.algorithm.create(cfg)
        astra_mod.algorithm.run(alg_id, niter)
        astra_mod.algorithm.delete(alg_id)
        astra_mod.data2d.delete(vid)
    astra_mod.data2d.delete(sid)
    astra_mod.projector.delete(pid)

def bp0(sino, det_col=111, num_angles=222, voxel_size_mm=1):
    """Wrapper for astra forward projector

    :param sino:
    :param projector_id:
    :return: backprojected sinogram
    """

    vol_geom = astra.create_vol_geom(sino.shape)
    proj_id = astra.create_projector(
        'cuda',
        astra.create_proj_geom('parallel', 1.0, det_col,
                               np.linspace(0, np.pi, num_angles, False)),
        vol_geom)


    rec_id, backprojection = astra.create_backprojection(sino * voxel_size_mm,
                                                   proj_id)

    # rec_id, backprojection = astra.create_backprojection(sino, projector_id)

    # backprojection /= sino.shape[0]
    # backprojection *= np.pi
    astra.data2d.delete(rec_id)
    astra.projector.delete(proj_id)

    return backprojection

def astra_rec_cuda(tomo, center, recon, theta, vol_geom, niter, proj_type, gpu_index, opts):
    # Lazy import ASTRA
    import astra as astra_mod
    nslices, nang, ndet = tomo.shape
    cfg = astra_mod.astra_dict(opts['method'])
    if 'extra_options' in opts:
        # NOTE: we are modifying 'extra_options' and so need to make a copy
        cfg['option'] = copy.deepcopy(opts['extra_options'])
    else:
        cfg['option'] = {}
    if gpu_index is not None:
        cfg['option']['GPUindex'] = gpu_index
    oc = None
    const_theta = np.ones(nang)
    proj_geom = astra_mod.create_proj_geom(
        'parallel', 1.0, ndet, theta.astype(np.float64))
    for i in range(nslices):
        if center[i] != oc:
            oc = center[i]
            proj_geom['option'] = {
                'ExtraDetectorOffset':
                (center[i] - ndet / 2.) * const_theta}
        pid = astra_mod.create_projector(proj_type, proj_geom, vol_geom)
        cfg['ProjectorId'] = pid
        sid = astra_mod.data2d.link('-sino', proj_geom, tomo[i])
        cfg['ProjectionDataId'] = sid
        vid = astra_mod.data2d.link('-vol', vol_geom, recon[i])
        cfg['ReconstructionDataId'] = vid
        alg_id = astra_mod.algorithm.create(cfg)
        astra_mod.algorithm.run(alg_id, niter)
        astra_mod.algorithm.delete(alg_id)
        astra_mod.data2d.delete(vid)
        astra_mod.data2d.delete(sid)
        astra_mod.projector.delete(pid)

def pid(number_of_angles):
    """Return projector ID.
    :param number_of_angles:
    :return:
    """
    return astra.create_projector('cuda', astra.create_proj_geom('parallel',
    1.0, vshape[0], np.linspace(0, np.pi, number_of_angles, False)), vol_geom)

 def run(self, iterations):
     self.fn = getFilterFile(self.fd, self.pg, iterations, self.rg, self.osz)
     if os.path.exists(self.fn):
         flt = np.load(self.fn)
         self.v[:] = self.customFBP(flt)
         return
     nd = self.nd
     if self.osz:
         nds = self.osz
     else:
         nds = nd
     na = len(self.ang)
     pgc = astra.create_proj_geom('parallel',1.0,nd,self.ang)
     vgc = astra.create_vol_geom((nds,nds))
     pidc = astra.create_projector('strip',pgc,vgc)
     x = np.zeros((nds,nds))
     xs = np.zeros((nds,nds))
     sf = np.zeros((na,nd))
     vid = astra.data2d.create('-vol',vgc)
     sid = astra.data2d.create('-sino',pgc)
     cfg = astra.astra_dict('FP')
     cfg['ProjectorId']=pidc
     cfg['ProjectionDataId']=sid
     cfg['VolumeDataId']=vid
     fpid = astra.algorithm.create(cfg)
     cfg = astra.astra_dict('BP')
     cfg['ProjectorId']=pidc
     cfg['ProjectionDataId']=sid
     cfg['ReconstructionDataId']=vid
     bpid = astra.algorithm.create(cfg)
     vc = astra.data2d.get_shared(vid)
     sc = astra.data2d.get_shared(sid)
     x[nds//2,nds//2]=1
     alp = 1./(na*nds)
     if self.rg:
         if self.rg*alp >=0.1:
             alp = 0.1/self.rg
     astra.log.info('Computing filter...')
     for i in range(iterations):
         if i%10==0: astra.log.info('{:.2f} % done'.format(100*float(i)/iterations))
         xs+=x
         vc[:] = x
         astra.algorithm.run(fpid)
         astra.algorithm.run(bpid)
         if self.rg:
             dx = x[:-1,:] - x[1:,:]
             dy = x[:,:-1] - x[:,1:]
             x[:-1,:] -= self.rg*dx*alp
             x[1:,:] += self.rg*dx*alp
             x[:,:-1] -= self.rg*dy*alp
             x[:,1:] += self.rg*dy*alp
         x -= vc*alp
     vc[:] = xs
     astra.algorithm.run(fpid)
     flt = sc.copy()*alp
     astra.algorithm.delete([fpid,bpid])
     astra.algorithm.delete([vid,sid])
     np.save(self.fn,flt)
     self.v[:] = self.customFBP(flt)
     astra.projector.delete(pidc)

def astrafp(im, ang, prj="cuda"):
    proj_geom = astra.create_proj_geom("parallel", 1.0, im.shape[0], np.array([ang, 0]))
    vol_geom = astra.create_vol_geom(im.shape)
    pid = astra.create_projector(prj, proj_geom, vol_geom)
    w = astra.OpTomo(pid)
    fpim = w * im
    astra.projector.delete(pid)
    return fpim[0 : im.shape[0]]

 def __init__(self, n_pixels, n_angles, rayperdetec=None):
     '''
     Initialize the ASTRA toolbox with a simple parallel configuration.
     The image is assumed to be square, and the detector count is equal to the number of rows/columns.
     '''
     self.vol_geom = astra.create_vol_geom(n_pixels, n_pixels)
     self.proj_geom = astra.create_proj_geom('parallel', 1.0, n_pixels, np.linspace(0,np.pi,n_angles,False))
     self.proj_id = astra.create_projector('cuda', self.proj_geom, self.vol_geom)

    def reconstruct(self, sinogram, centre_of_rotation, angles, shape, center):

        ctr = centre_of_rotation
        width = sinogram.shape[1]
        pad = 50

        sino = np.nan_to_num(sinogram)

        # pad the array so that the centre of rotation is in the middle
        alen = ctr
        blen = width - ctr
        mid = width / 2.0

        if (ctr > mid):
            plow = pad
            phigh = (alen - blen) + pad
        else:
            plow = (blen - alen) + pad
            phigh = pad

        logdata = np.log(sino+1)
        sinogram = np.pad(logdata, ((0, 0), (int(plow), int(phigh))),
                          mode='reflect')

        width = sinogram.shape[1]

        vol_geom = astra.create_vol_geom(shape[0], shape[1])
        proj_geom = astra.create_proj_geom('parallel', 1.0, width,
                                           np.deg2rad(angles))

        sinogram_id = astra.data2d.create("-sino", proj_geom, sinogram)

        # Create a data object for the reconstruction
        rec_id = astra.data2d.create('-vol', vol_geom)

        proj_id = astra.create_projector('strip', proj_geom, vol_geom)

        cfg = astra.astra_dict(self.parameters['reconstruction_type'])
        cfg['ReconstructionDataId'] = rec_id
        cfg['ProjectionDataId'] = sinogram_id
        cfg['ProjectorId'] = proj_id

        # Create the algorithm object from the configuration structure
        alg_id = astra.algorithm.create(cfg)
        # Run 20 iterations of the algorithm
        itterations = int(self.parameters['number_of_iterations'])
        # This will have a runtime in the order of 10 seconds.
        astra.algorithm.run(alg_id, itterations)
        # Get the result
        rec = astra.data2d.get(rec_id)

        # Clean up.
        astra.algorithm.delete(alg_id)
        astra.data2d.delete(rec_id)
        astra.data2d.delete(sinogram_id)

        return rec

 def astra_reconstruction(self, sino, vol_geom, proj_geom):
     # currently hard-coded - for 3D version only!
     #sino = np.transpose(sino, (1, 0, 2))
     self.sino_id = self.astra_function.create("-sino", proj_geom, sino)
     # Create a data object for the reconstruction
     self.rec_id = self.astra_function.create('-vol', vol_geom)
     cfg = self.cfg_setup()
     if self.alg_type is '2D' and not "CUDA" in self.name:
         proj_id = astra.create_projector('strip', proj_geom, vol_geom)
         cfg['ProjectorId'] = proj_id
     return self.run_astra(cfg)

def _basic_par2d_fp(type):
  import astra
  import numpy as np
  vg = astra.create_vol_geom(2, 32)
  pg = astra.create_proj_geom('parallel', 1, 32, [0])
  proj_id = astra.create_projector(type, pg, vg)
  vol = np.random.rand(2, 32)
  (sino_id, sino) = astra.create_sino(vol, proj_id)
  astra.data2d.delete(sino_id)
  astra.projector.delete(proj_id)
  err = np.max(np.abs(sino[0,:] - np.sum(vol,axis=0)))
  return err < 1e-6

 def set_config(self, rec_id, sino_id, proj_geom, vol_geom):
     cfg = astra.astra_dict(self.alg)
     cfg['ReconstructionDataId'] = rec_id
     cfg['ProjectionDataId'] = sino_id
     if 'FBP' in self.alg:
         fbp_filter = self.parameters['FBP_filter'] if 'FBP_filter' in \
             self.parameters.keys() else 'none'
         cfg['FilterType'] = fbp_filter
     if 'projector' in self.parameters.keys():
         proj_id = astra.create_projector(
             self.parameters['projector'], proj_geom, vol_geom)
         cfg['ProjectorId'] = proj_id
     cfg = self.set_options(cfg)
     return cfg

 def set_config(self, rec_id, sino_id, proj_geom, vol_geom):
     cfg = astra.astra_dict(self.alg)
     cfg['ReconstructionDataId'] = rec_id
     cfg['ProjectionDataId'] = sino_id
     if 'FBP' in self.alg:
         cfg['FilterType'] = self.parameters['FBP_filter']
     if 'projector' in self.parameters.keys():
         proj_id = astra.create_projector(
             self.parameters['projector'], proj_geom, vol_geom)
         cfg['ProjectorId'] = proj_id
     # mask not currently working correctly for SIRT or SART algorithms
     sirt_or_sart = [a for a in ['SIRT', 'SART'] if a in self.alg]
     if self.mask_id and not sirt_or_sart:
         cfg['option'] = {}
         cfg['option']['ReconstructionMaskId'] = self.mask_id
     cfg = self.set_options(cfg)
     return cfg

    def geom_setup_2D(self, dim_det_cols):
        in_pData = self.get_plugin_in_datasets()[0]
        cor = self.cor + self.pad_amount

        sino_shape = in_pData.get_shape()
        nCols = sino_shape[dim_det_cols]
        self.p_low, self.p_high = \
            self.array_pad(cor[0], nCols + 2*self.pad_amount)
        sino_width = \
            sino_shape[1] + self.p_low + self.p_high + 2*self.pad_amount

        vol_geom = astra.create_vol_geom(self.vol_shape[0], self.vol_shape[1])
        self.proj_geom = astra.create_proj_geom('parallel', 1.0, sino_width,
                                                np.deg2rad(self.angles))

        self.rec_id = self.astra_function.create('-vol', vol_geom)
        self.cfg = self.cfg_setup()
        if self.alg_type is '2D' and "CUDA" not in self.name:
            proj_id = astra.create_projector('strip', self.proj_geom, vol_geom)
            self.cfg['ProjectorId'] = proj_id

def fp0(image, det_col=111, num_angles=222, voxel_size_mm=1):
    """Wrapper for astra forward projector

    :param image:
    :return: sinogram
    """

    vol_geom = astra.create_vol_geom(image.shape)
    proj_id = astra.create_projector(
        'cuda',
        astra.create_proj_geom('parallel', 1.0, det_col,
                               np.linspace(0, np.pi, num_angles, False)),
        vol_geom)

    sino_id, sino = astra.create_sino(image, proj_id)
    sino *= voxel_size_mm

    astra.data2d.delete(sino_id)
    astra.projector.delete(proj_id)

    return sino