summaryrefslogtreecommitdiffstats
diff options
context:
space:
mode:
authorEdoardo Pasca <edo.paskino@gmail.com>2017-10-19 12:47:02 +0100
committerEdoardo Pasca <edo.paskino@gmail.com>2017-10-19 12:47:02 +0100
commit99cfddc7764a292c438868a02c9d512dabefbdc4 (patch)
tree2604436119f3621906dee04ecd4d4071ab8a510e
parentcb8ef11f00e897b6f5b3049126dd32baa3c50cf9 (diff)
downloadregularization-99cfddc7764a292c438868a02c9d512dabefbdc4.tar.gz
regularization-99cfddc7764a292c438868a02c9d512dabefbdc4.tar.bz2
regularization-99cfddc7764a292c438868a02c9d512dabefbdc4.tar.xz
regularization-99cfddc7764a292c438868a02c9d512dabefbdc4.zip
moved code out of if-then-else
fixed #5
-rw-r--r--main_func/FISTA_REC.m54
1 files changed, 15 insertions, 39 deletions
diff --git a/main_func/FISTA_REC.m b/main_func/FISTA_REC.m
index bea1860..d177748 100644
--- a/main_func/FISTA_REC.m
+++ b/main_func/FISTA_REC.m
@@ -517,6 +517,18 @@ else
% subsets loop
counterInd = 1;
+ if (strcmp(proj_geom.type,'parallel') || strcmp(proj_geom.type,'fanflat') || strcmp(proj_geom.type,'fanflat_vec'))
+ % if geometry is 2D use slice-by-slice projection-backprojection routine
+ for kkk = 1:SlicesZ
+ [sino_id, sinoT] = astra_create_sino_cuda(X_t(:,:,kkk), proj_geomSUB, vol_geom);
+ sino_updt_Sub(:,:,kkk) = sinoT';
+ astra_mex_data2d('delete', sino_id);
+ end
+ else
+ % for 3D geometry (watch the GPU memory overflow in earlier ASTRA versions < 1.8)
+ [sino_id, sino_updt_Sub] = astra_create_sino3d_cuda(X_t, proj_geomSUB, vol_geom);
+ astra_mex_data3d('delete', sino_id);
+ end
for ss = 1:subsets
X_old = X;
t_old = t;
@@ -528,18 +540,7 @@ else
if (lambdaR_L1 > 0)
% Group-Huber fidelity (ring removal)
- if (strcmp(proj_geom.type,'parallel') || strcmp(proj_geom.type,'fanflat') || strcmp(proj_geom.type,'fanflat_vec'))
- % if geometry is 2D use slice-by-slice projection-backprojection routine
- for kkk = 1:SlicesZ
- [sino_id, sinoT] = astra_create_sino_cuda(X_t(:,:,kkk), proj_geomSUB, vol_geom);
- sino_updt_Sub(:,:,kkk) = sinoT';
- astra_mex_data2d('delete', sino_id);
- end
- else
- % for 3D geometry (watch the GPU memory overflow in earlier ASTRA versions < 1.8)
- [sino_id, sino_updt_Sub] = astra_create_sino3d_cuda(X_t, proj_geomSUB, vol_geom);
- astra_mex_data3d('delete', sino_id);
- end
+
residualSub = zeros(Detectors, numProjSub, SlicesZ,'single'); % residual for a chosen subset
for kkk = 1:numProjSub
@@ -550,18 +551,7 @@ else
elseif (studentt > 0)
% student t data fidelity
- if (strcmp(proj_geom.type,'parallel') || strcmp(proj_geom.type,'fanflat') || strcmp(proj_geom.type,'fanflat_vec'))
- % if geometry is 2D use slice-by-slice projection-backprojection routine
- for kkk = 1:SlicesZ
- [sino_id, sinoT] = astra_create_sino_cuda(X_t(:,:,kkk), proj_geomSUB, vol_geom);
- sino_updt_Sub(:,:,kkk) = sinoT';
- astra_mex_data2d('delete', sino_id);
- end
- else
- % for 3D geometry (watch the GPU memory overflow in earlier ASTRA versions < 1.8)
- [sino_id, sino_updt_Sub] = astra_create_sino3d_cuda(X_t, proj_geomSUB, vol_geom);
- astra_mex_data3d('delete', sino_id);
- end
+
% artifacts removal with Students t penalty
residualSub = squeeze(weights(:,CurrSubIndeces,:)).*(sino_updt_Sub - squeeze(sino(:,CurrSubIndeces,:)));
@@ -576,21 +566,7 @@ else
objective(i) = ff; % for the objective function output
else
% PWLS model
- if (strcmp(proj_geom.type,'parallel') || strcmp(proj_geom.type,'fanflat') || strcmp(proj_geom.type,'fanflat_vec'))
- % if geometry is 2D use slice-by-slice projection-backprojection routine
- for kkk = 1:SlicesZ
- [sino_id, sinoT] = astra_create_sino_cuda(X_t(:,:,kkk), proj_geomSUB, vol_geom);
- sino_updt_Sub(:,:,kkk) = sinoT';
- astra_mex_data2d('delete', sino_id);
- end
- else
- % for 3D geometry (watch the GPU memory overflow in earlier ASTRA versions < 1.8)
- [sino_id, sino_updt_Sub] = astra_create_sino3d_cuda(X_t, proj_geomSUB, vol_geom);
- astra_mex_data3d('delete', sino_id);
- end
- residualSub = squeeze(weights(:,CurrSubIndeces,:)).*(sino_updt_Sub - squeeze(sino(:,CurrSubIndeces,:)));
- objective(i) = 0.5*norm(residualSub(:)); % for the objective function output
- end
+
if (strcmp(proj_geom.type,'parallel') || strcmp(proj_geom.type,'fanflat') || strcmp(proj_geom.type,'fanflat_vec'))
% if geometry is 2D use slice-by-slice projection-backprojection routine