From 27752838308df13ac174047e0079e38d55f990cb Mon Sep 17 00:00:00 2001 From: jakobsj Date: Tue, 24 Apr 2018 14:44:06 +0100 Subject: Ip removal (#109) * Removing IP data and files. Still in IP_parking branch. #102 * Also remove IP pyc file --- Wrappers/Python/wip/DemoRecIP.py | 110 ---------------------------- data/IP_data70channels.mat | Bin 8710147 -> 0 bytes data/__pycache__/read_IPdata.cpython-35.pyc | Bin 1674 -> 0 bytes data/read_IPdata.py | 58 --------------- 4 files changed, 168 deletions(-) delete mode 100644 Wrappers/Python/wip/DemoRecIP.py delete mode 100644 data/IP_data70channels.mat delete mode 100644 data/__pycache__/read_IPdata.cpython-35.pyc delete mode 100644 data/read_IPdata.py diff --git a/Wrappers/Python/wip/DemoRecIP.py b/Wrappers/Python/wip/DemoRecIP.py deleted file mode 100644 index 442e40e..0000000 --- a/Wrappers/Python/wip/DemoRecIP.py +++ /dev/null @@ -1,110 +0,0 @@ -#!/usr/bin/env python3 -# -*- coding: utf-8 -*- -""" -Reading multi-channel data and reconstruction using FISTA modular -""" - -import numpy as np -import matplotlib.pyplot as plt - -#import sys -#sys.path.append('../../../data/') -from read_IPdata import read_IPdata - -from ccpi.astra.astra_ops import AstraProjectorSimple, AstraProjectorMC -from ccpi.reconstruction.funcs import Norm2sq, Norm1, BaseFunction -from ccpi.reconstruction.algs import FISTA -#from ccpi.reconstruction.funcs import BaseFunction - -from ccpi.framework import ImageData, AcquisitionData, AcquisitionGeometry, ImageGeometry - -# read IP paper data into a dictionary -dataDICT = read_IPdata('..\..\..\data\IP_data70channels.mat') - -# Set ASTRA Projection-backprojection class (fan-beam geometry) -DetWidth = dataDICT.get('im_size')[0] * dataDICT.get('det_width')[0] / \ - dataDICT.get('detectors_numb')[0] -SourceOrig = dataDICT.get('im_size')[0] * dataDICT.get('src_to_rotc')[0] / \ - dataDICT.get('dom_width')[0] -OrigDetec = dataDICT.get('im_size')[0] * \ - (dataDICT.get('src_to_det')[0] - dataDICT.get('src_to_rotc')[0]) /\ - dataDICT.get('dom_width')[0] - -N = dataDICT.get('im_size')[0] - -vg = ImageGeometry(voxel_num_x=dataDICT.get('im_size')[0], - voxel_num_y=dataDICT.get('im_size')[0], - channels=1) - -pg = AcquisitionGeometry('cone', - '2D', - angles=(np.pi/180)*dataDICT.get('theta')[0], - pixel_num_h=dataDICT.get('detectors_numb')[0], - pixel_size_h=DetWidth, - dist_source_center=SourceOrig, - dist_center_detector=OrigDetec, - channels=1) - - -sino = dataDICT.get('data_norm')[0][:,:,34] # select mid-channel -b = AcquisitionData(sino,geometry=pg) - -# Initial guess -x_init = ImageData(np.zeros((N, N)),geometry=vg) - - - - - -Aop = AstraProjectorSimple(vg,pg,'gpu') -f = Norm2sq(Aop,b,c=0.5) - -# Run FISTA for least squares without regularization -opt = {'tol': 1e-4, 'iter': 10} -x_fista0, it0, timing0, criter0 = FISTA(x_init, f, None, opt) - -plt.imshow(x_fista0.array) -plt.show() - -# Now least squares plus 1-norm regularization -g1 = Norm1(10) - -# Run FISTA for least squares plus 1-norm function. -x_fista1, it1, timing1, criter1 = FISTA(x_init, f, g1, opt) - -plt.imshow(x_fista1.array) -plt.show() - -# Multiple channels -sino_mc = dataDICT.get('data_norm')[0][:,:,32:37] # select mid-channel - -vg_mc = ImageGeometry(voxel_num_x=dataDICT.get('im_size')[0], - voxel_num_y=dataDICT.get('im_size')[0], - channels=5) - -pg_mc = AcquisitionGeometry('cone', - '2D', - angles=(np.pi/180)*dataDICT.get('theta')[0], - pixel_num_h=dataDICT.get('detectors_numb')[0], - pixel_size_h=DetWidth, - dist_source_center=SourceOrig, - dist_center_detector=OrigDetec, - channels=5) - -b_mc = AcquisitionData(np.transpose(sino_mc,(2,0,1)), - geometry=pg_mc, - dimension_labels=("channel","angle","horizontal")) - -# ASTRA operator using volume and sinogram geometries -Aop_mc = AstraProjectorMC(vg_mc, pg_mc, 'gpu') - -f_mc = Norm2sq(Aop_mc,b_mc,c=0.5) - -# Initial guess -x_init_mc = ImageData(np.zeros((5, N, N)),geometry=vg_mc) - - -x_fista0_mc, it0_mc, timing0_mc, criter0_mc = FISTA(x_init_mc, f_mc, None, opt) - -plt.imshow(x_fista0_mc.as_array()[4,:,:]) -plt.show() \ No newline at end of file diff --git a/data/IP_data70channels.mat b/data/IP_data70channels.mat deleted file mode 100644 index 19b9421..0000000 Binary files a/data/IP_data70channels.mat and /dev/null differ diff --git a/data/__pycache__/read_IPdata.cpython-35.pyc b/data/__pycache__/read_IPdata.cpython-35.pyc deleted file mode 100644 index b062cf9..0000000 Binary files a/data/__pycache__/read_IPdata.cpython-35.pyc and /dev/null differ diff --git a/data/read_IPdata.py b/data/read_IPdata.py deleted file mode 100644 index a7565d7..0000000 --- a/data/read_IPdata.py +++ /dev/null @@ -1,58 +0,0 @@ -#!/usr/bin/env python3 -# -*- coding: utf-8 -*- -""" -function to read IP data and provide a dictionary with data and parameters as an output -""" -from scipy import io -import numpy as np -from collections import defaultdict - -def read_IPdata(datafile): - # read data from mat file (specify the location) - alldata = io.loadmat(datafile) - data_raw = alldata.get('Data_raw') # here is raw projection data - Phantom_ideal = alldata.get('Phantom_ideal') # here is 70 channels ideal phantom - Photon_flux = alldata.get('Photon_flux') # photon flux for normalization - del alldata - - # extract geometry-related parameters - proj_numb,detectors_numb,channels = np.shape(data_raw) - im_size = np.size(Phantom_ideal,1) - - theta = np.linspace(0,proj_numb-1,proj_numb)*360/proj_numb # projection angles - dom_width = 1.0 # width of domain in cm - src_to_rotc = 3.0 # dist. from source to rotation center - src_to_det = 5.0 # dist. from source to detector - det_width = 2.0 # detector width - - # negative log normalisation of the raw data (avoiding of log(0)) - data_norm = np.zeros(np.shape(data_raw)) - for i in range(0,channels): - slice1 = data_raw[:,:,i] - indx = np.nonzero(slice1>0) - slice2 = np.zeros((proj_numb,detectors_numb), 'float32') - slice2[indx] = -np.log(slice1[indx]/Photon_flux[i]) - indx2 = np.nonzero(slice1==0) - slice3 = np.zeros((proj_numb,detectors_numb), 'float32') - slice3[indx2] = np.log(slice2[indx2]+Photon_flux[i]) - data_norm[:,:,i] = slice2 + slice3 - del indx, indx2, slice1, slice2, slice3 - data_norm = np.float32(data_norm*(im_size/dom_width)) - - #build a dictionary for data and related parameters - dataDICT = defaultdict(list) - dataDICT['data_norm'].append(data_norm) - dataDICT['data_raw'].append(data_raw) - dataDICT['Photon_flux'].append(Photon_flux) - dataDICT['Phantom_ideal'].append(Phantom_ideal) - dataDICT['theta'].append(theta) - dataDICT['proj_numb'].append(proj_numb) - dataDICT['detectors_numb'].append(detectors_numb) - dataDICT['channels'].append(channels) - dataDICT['im_size'].append(im_size) - dataDICT['dom_width'].append(dom_width) - dataDICT['src_to_rotc'].append(src_to_rotc) - dataDICT['src_to_det'].append(src_to_det) - dataDICT['det_width'].append(det_width) - - return (dataDICT) \ No newline at end of file -- cgit v1.2.3