From 56809b0359af7e9108adeb1fd21823a225edf6fa Mon Sep 17 00:00:00 2001 From: Willem Jan Palenstijn Date: Wed, 20 Jan 2016 18:08:59 +0100 Subject: Remove dependency of libastra on libpython by refactoring PluginAlgorithm --- matlab/mex/astra_mex_plugin_c.cpp | 18 +++++++++++++++--- 1 file changed, 15 insertions(+), 3 deletions(-) (limited to 'matlab') diff --git a/matlab/mex/astra_mex_plugin_c.cpp b/matlab/mex/astra_mex_plugin_c.cpp index 177fcf4..0dc0edf 100644 --- a/matlab/mex/astra_mex_plugin_c.cpp +++ b/matlab/mex/astra_mex_plugin_c.cpp @@ -48,7 +48,11 @@ using namespace astra; */ void astra_mex_plugin_get_registered(int nlhs, mxArray* plhs[], int nrhs, const mxArray* prhs[]) { - astra::CPluginAlgorithmFactory *fact = astra::CPluginAlgorithmFactory::getSingletonPtr(); + astra::CPluginAlgorithmFactory *fact = astra::CPluginAlgorithmFactory::getFactory(); + if (!fact) { + mexPrintf("Plugin support not initialized."); + return; + } std::map mp = fact->getRegisteredMap(); for(std::map::iterator it=mp.begin();it!=mp.end();it++){ mexPrintf("%s: %s\n",it->first.c_str(), it->second.c_str()); @@ -62,9 +66,13 @@ void astra_mex_plugin_get_registered(int nlhs, mxArray* plhs[], int nrhs, const */ void astra_mex_plugin_register(int nlhs, mxArray* plhs[], int nrhs, const mxArray* prhs[]) { + astra::CPluginAlgorithmFactory *fact = astra::CPluginAlgorithmFactory::getFactory(); + if (!fact) { + mexPrintf("Plugin support not initialized."); + return; + } if (2 <= nrhs) { string class_name = mexToString(prhs[1]); - astra::CPluginAlgorithmFactory *fact = astra::CPluginAlgorithmFactory::getSingletonPtr(); fact->registerPlugin(class_name); }else{ mexPrintf("astra_mex_plugin('register', class_name);\n"); @@ -78,9 +86,13 @@ void astra_mex_plugin_register(int nlhs, mxArray* plhs[], int nrhs, const mxArra */ void astra_mex_plugin_get_help(int nlhs, mxArray* plhs[], int nrhs, const mxArray* prhs[]) { + astra::CPluginAlgorithmFactory *fact = astra::CPluginAlgorithmFactory::getFactory(); + if (!fact) { + mexPrintf("Plugin support not initialized."); + return; + } if (2 <= nrhs) { string name = mexToString(prhs[1]); - astra::CPluginAlgorithmFactory *fact = astra::CPluginAlgorithmFactory::getSingletonPtr(); mexPrintf((fact->getHelp(name)+"\n").c_str()); }else{ mexPrintf("astra_mex_plugin('get_help', name);\n"); -- cgit v1.2.3 From 7f11ed333cdc93e9e93d4677f857d97ffbbfec87 Mon Sep 17 00:00:00 2001 From: Willem Jan Palenstijn Date: Thu, 21 Jan 2016 12:14:44 +0100 Subject: Fix whitespace --- matlab/mex/astra_mex_plugin_c.cpp | 12 ++++++------ 1 file changed, 6 insertions(+), 6 deletions(-) (limited to 'matlab') diff --git a/matlab/mex/astra_mex_plugin_c.cpp b/matlab/mex/astra_mex_plugin_c.cpp index 0dc0edf..627e34a 100644 --- a/matlab/mex/astra_mex_plugin_c.cpp +++ b/matlab/mex/astra_mex_plugin_c.cpp @@ -128,12 +128,12 @@ void mexFunction(int nlhs, mxArray* plhs[], initASTRAMex(); // SWITCH (MODE) - if (sMode == std::string("get_registered")) { - astra_mex_plugin_get_registered(nlhs, plhs, nrhs, prhs); - }else if (sMode == std::string("get_help")) { - astra_mex_plugin_get_help(nlhs, plhs, nrhs, prhs); - }else if (sMode == std::string("register")) { - astra_mex_plugin_register(nlhs, plhs, nrhs, prhs); + if (sMode == std::string("get_registered")) { + astra_mex_plugin_get_registered(nlhs, plhs, nrhs, prhs); + }else if (sMode == std::string("get_help")) { + astra_mex_plugin_get_help(nlhs, plhs, nrhs, prhs); + }else if (sMode == std::string("register")) { + astra_mex_plugin_register(nlhs, plhs, nrhs, prhs); } else { printHelp(); } -- cgit v1.2.3 From 96fa285132bf88462c33c5b24cc1b241fb3d4d73 Mon Sep 17 00:00:00 2001 From: Willem Jan Palenstijn Date: Thu, 21 Jan 2016 14:30:32 +0100 Subject: Initialize Python plugins when available in Matlab --- matlab/mex/astra_mex_plugin_c.cpp | 25 ++++++++++++++++++++++++- matlab/mex/mexInitFunctions.cpp | 8 ++++++++ 2 files changed, 32 insertions(+), 1 deletion(-) (limited to 'matlab') diff --git a/matlab/mex/astra_mex_plugin_c.cpp b/matlab/mex/astra_mex_plugin_c.cpp index 627e34a..9bd4446 100644 --- a/matlab/mex/astra_mex_plugin_c.cpp +++ b/matlab/mex/astra_mex_plugin_c.cpp @@ -37,10 +37,31 @@ $Id$ #include "astra/PluginAlgorithm.h" +#include + using namespace std; using namespace astra; +//----------------------------------------------------------------------------------------- +/** astra_mex_plugin('init'); + * + * Initialize plugin support by initializing python and importing astra + */ +void astra_mex_plugin_init() +{ + if(!Py_IsInitialized()){ + Py_Initialize(); + PyEval_InitThreads(); + } + + // Importing astra may be overkill, since we only need to initialize + // PythonPluginAlgorithmFactory from astra.plugin_c. + PyObject *mod = PyImport_ImportModule("astra"); + Py_DECREF(mod); +} + + //----------------------------------------------------------------------------------------- /** astra_mex_plugin('get_registered'); * @@ -128,7 +149,9 @@ void mexFunction(int nlhs, mxArray* plhs[], initASTRAMex(); // SWITCH (MODE) - if (sMode == std::string("get_registered")) { + if (sMode == "init") { + astra_mex_plugin_init(); + } else if (sMode == std::string("get_registered")) { astra_mex_plugin_get_registered(nlhs, plhs, nrhs, prhs); }else if (sMode == std::string("get_help")) { astra_mex_plugin_get_help(nlhs, plhs, nrhs, prhs); diff --git a/matlab/mex/mexInitFunctions.cpp b/matlab/mex/mexInitFunctions.cpp index bd3df2c..7245af2 100644 --- a/matlab/mex/mexInitFunctions.cpp +++ b/matlab/mex/mexInitFunctions.cpp @@ -23,5 +23,13 @@ void initASTRAMex(){ if(!astra::CLogger::setCallbackScreen(&logCallBack)){ mexErrMsgTxt("Error initializing mex functions."); } + mexIsInitialized=true; + + + // If we have support for plugins, initialize them. + // (NB: Call this after setting mexIsInitialized, to avoid recursively + // calling initASTRAMex) + mexEvalString("if exist('astra_mex_plugin_c') == 3; astra_mex_plugin_c('init'); end"); + } -- cgit v1.2.3 From 5cb48f351fcc1a7c5b11180edcded083176a2431 Mon Sep 17 00:00:00 2001 From: Willem Jan Palenstijn Date: Wed, 27 Jan 2016 09:28:10 +0100 Subject: Use Py_XDECREF to guard against null --- matlab/mex/astra_mex_plugin_c.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) (limited to 'matlab') diff --git a/matlab/mex/astra_mex_plugin_c.cpp b/matlab/mex/astra_mex_plugin_c.cpp index 9bd4446..a279df8 100644 --- a/matlab/mex/astra_mex_plugin_c.cpp +++ b/matlab/mex/astra_mex_plugin_c.cpp @@ -58,7 +58,7 @@ void astra_mex_plugin_init() // Importing astra may be overkill, since we only need to initialize // PythonPluginAlgorithmFactory from astra.plugin_c. PyObject *mod = PyImport_ImportModule("astra"); - Py_DECREF(mod); + Py_XDECREF(mod); } -- cgit v1.2.3 From 4a75ab09962015d90acaa3e2be8357a1c917d5eb Mon Sep 17 00:00:00 2001 From: Willem Jan Palenstijn Date: Wed, 27 Jan 2016 10:05:57 +0100 Subject: Move fixLapackLoading to astra_mex_plugin --- matlab/mex/astra_mex_plugin_c.cpp | 33 +++++++++++++++++++++++++++++++++ 1 file changed, 33 insertions(+) (limited to 'matlab') diff --git a/matlab/mex/astra_mex_plugin_c.cpp b/matlab/mex/astra_mex_plugin_c.cpp index a279df8..4ed534e 100644 --- a/matlab/mex/astra_mex_plugin_c.cpp +++ b/matlab/mex/astra_mex_plugin_c.cpp @@ -42,6 +42,35 @@ $Id$ using namespace std; using namespace astra; +static void fixLapackLoading() +{ + // When running in Matlab, we need to force numpy + // to use its internal lapack library instead of + // Matlab's MKL library to avoid errors. To do this, + // we set Python's dlopen flags to RTLD_NOW|RTLD_DEEPBIND + // and import 'numpy.linalg.lapack_lite' here. We reset + // Python's dlopen flags afterwards. + PyObject *sys = PyImport_ImportModule("sys"); + if (sys != NULL) { + PyObject *curFlags = PyObject_CallMethod(sys, "getdlopenflags", NULL); + if (curFlags != NULL) { + PyObject *retVal = PyObject_CallMethod(sys, "setdlopenflags", "i", 10); // RTLD_NOW|RTLD_DEEPBIND + if (retVal != NULL) { + PyObject *lapack = PyImport_ImportModule("numpy.linalg.lapack_lite"); + if (lapack != NULL) { + Py_DECREF(lapack); + } + PyObject *retVal2 = PyObject_CallMethod(sys, "setdlopenflags", "O",curFlags); + if (retVal2 != NULL) { + Py_DECREF(retVal2); + } + Py_DECREF(retVal); + } + Py_DECREF(curFlags); + } + Py_DECREF(sys); + } +} //----------------------------------------------------------------------------------------- /** astra_mex_plugin('init'); @@ -55,6 +84,10 @@ void astra_mex_plugin_init() PyEval_InitThreads(); } +#ifndef _MSC_VER + fixLapackLoading(); +#endif + // Importing astra may be overkill, since we only need to initialize // PythonPluginAlgorithmFactory from astra.plugin_c. PyObject *mod = PyImport_ImportModule("astra"); -- cgit v1.2.3 From 8cfe13a410a051e5a6b9835e3a4cd6c808cf5d29 Mon Sep 17 00:00:00 2001 From: Willem Jan Palenstijn Date: Tue, 16 Feb 2016 15:34:08 +0100 Subject: Add astra_mex delete/info based on index manager --- matlab/mex/astra_mex_c.cpp | 49 ++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 49 insertions(+) (limited to 'matlab') diff --git a/matlab/mex/astra_mex_c.cpp b/matlab/mex/astra_mex_c.cpp index fdf4f33..a8623be 100644 --- a/matlab/mex/astra_mex_c.cpp +++ b/matlab/mex/astra_mex_c.cpp @@ -36,10 +36,14 @@ $Id$ #include "mexInitFunctions.h" #include "astra/Globals.h" +#include "astra/AstraObjectManager.h" + #ifdef ASTRA_CUDA #include "../cuda/2d/darthelper.h" #include "astra/CompositeGeometryManager.h" #endif + + using namespace std; using namespace astra; @@ -142,6 +146,47 @@ void astra_mex_version(int nlhs, mxArray* plhs[], int nrhs, const mxArray* prhs[ } } +//----------------------------------------------------------------------------------------- + +void astra_mex_info(int nlhs, mxArray* plhs[], int nrhs, const mxArray* prhs[]) +{ + if (nrhs < 2) { + mexErrMsgTxt("Usage: astra_mex('info', index/indices);\n"); + return; + } + + for (int i = 1; i < nrhs; i++) { + int iDataID = (int)(mxGetScalar(prhs[i])); + CAstraObjectManagerBase *ptr; + ptr = CAstraIndexManager::getSingleton().get(iDataID); + if (ptr) { + mexPrintf("%s\t%s\n", ptr->getType().c_str(), ptr->getInfo(iDataID).c_str()); + } + } + +} + +//----------------------------------------------------------------------------------------- + +void astra_mex_delete(int nlhs, mxArray* plhs[], int nrhs, const mxArray* prhs[]) +{ + if (nrhs < 2) { + mexErrMsgTxt("Usage: astra_mex('delete', index/indices);\n"); + return; + } + + for (int i = 1; i < nrhs; i++) { + int iDataID = (int)(mxGetScalar(prhs[i])); + CAstraObjectManagerBase *ptr; + ptr = CAstraIndexManager::getSingleton().get(iDataID); + if (ptr) + ptr->remove(iDataID); + } + +} + + + //----------------------------------------------------------------------------------------- static void printHelp() @@ -178,6 +223,10 @@ void mexFunction(int nlhs, mxArray* plhs[], astra_mex_credits(nlhs, plhs, nrhs, prhs); } else if (sMode == std::string("set_gpu_index")) { astra_mex_set_gpu_index(nlhs, plhs, nrhs, prhs); + } else if (sMode == std::string("info")) { + astra_mex_info(nlhs, plhs, nrhs, prhs); + } else if (sMode == std::string("delete")) { + astra_mex_delete(nlhs, plhs, nrhs, prhs); } else { printHelp(); } -- cgit v1.2.3 From da11f9e7c7559eb01cc422ac55c9f7a1cd4bb803 Mon Sep 17 00:00:00 2001 From: Willem Jan Palenstijn Date: Mon, 22 Feb 2016 12:17:53 +0100 Subject: Add missing help --- matlab/mex/astra_mex_c.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) (limited to 'matlab') diff --git a/matlab/mex/astra_mex_c.cpp b/matlab/mex/astra_mex_c.cpp index a8623be..f499528 100644 --- a/matlab/mex/astra_mex_c.cpp +++ b/matlab/mex/astra_mex_c.cpp @@ -192,7 +192,7 @@ void astra_mex_delete(int nlhs, mxArray* plhs[], int nrhs, const mxArray* prhs[] static void printHelp() { mexPrintf("Please specify a mode of operation.\n"); - mexPrintf(" Valid modes: version, use_cuda, credits\n"); + mexPrintf(" Valid modes: version, use_cuda, credits, set_gpu_index, info, delete\n"); } //----------------------------------------------------------------------------------------- -- cgit v1.2.3 From 5a81584cb37fe578d5985b8e78b511076fa75b1c Mon Sep 17 00:00:00 2001 From: Folkert Bleichrodt Date: Wed, 2 Mar 2016 16:47:25 +0100 Subject: opTomo: output type matches input type If opTomo is used to do a forward or backprojection, the precision of the output data now matches the precision of the input data. So the output will be single precision only if the input is stored in single precision. --- matlab/tools/opTomo.m | 86 +++++++++++++++++---------------------------------- 1 file changed, 29 insertions(+), 57 deletions(-) (limited to 'matlab') diff --git a/matlab/tools/opTomo.m b/matlab/tools/opTomo.m index 71dfb1e..33c8565 100644 --- a/matlab/tools/opTomo.m +++ b/matlab/tools/opTomo.m @@ -44,11 +44,9 @@ classdef opTomo < opSpot vol_id fp_alg_id bp_alg_id + proj_id % ASTRA IDs handle astra_handle - % geometries - proj_geom; - vol_geom; end % properties properties ( SetAccess = private, GetAccess = public ) @@ -139,6 +137,17 @@ classdef opTomo < opSpot error(['Only type ' 39 'cuda' 39 ' is supported ' ... 'for 3D geometries.']) end + + % setup projector + cfg = astra_struct('cuda3d'); + cfg.ProjectionGeometry = proj_geom; + cfg.VolumeGeometry = vol_geom; + cfg.option.GPUindex = gpu_index; + + % create projector + op.proj_id = astra_mex_projector3d('create', cfg); + % create handle to ASTRA object, for cleaning up + op.astra_handle = opTomo_helper_handle(op.proj_id); % create a function handle op.funHandle = @opTomo_intrnl3D; @@ -148,8 +157,6 @@ classdef opTomo < opSpot % pass object properties op.proj_size = proj_size; op.vol_size = vol_size; - op.proj_geom = proj_geom; - op.vol_geom = vol_geom; op.cflag = false; op.sweepflag = false; @@ -170,11 +177,6 @@ classdef opTomo < opSpot x = full(x); end - % convert input to single - if isa(x, 'single') == false - x = single(x); - end - % the multiplication y = op.funHandle(op, x, mode); @@ -194,7 +196,7 @@ classdef opTomo < opSpot function y = opTomo_intrnl2D(op,x,mode) if mode == 1 - % X is passed as a vector, reshape it into an image. + % x is passed as a vector, reshape it into an image. x = reshape(x, op.vol_size); % Matlab data copied to ASTRA data @@ -204,9 +206,13 @@ classdef opTomo < opSpot astra_mex_algorithm('iterate', op.fp_alg_id); % retrieve Matlab array - y = astra_mex_data2d('get_single', op.sino_id); + if isa(x, 'single') + y = astra_mex_data2d('get_single', op.sino_id); + else + y = astra_mex_data2d('get', op.sino_id); + end else - % X is passed as a vector, reshape it into a sinogram. + % x is passed as a vector, reshape it into a sinogram. x = reshape(x, op.proj_size); % Matlab data copied to ASTRA data @@ -216,8 +222,13 @@ classdef opTomo < opSpot astra_mex_algorithm('iterate', op.bp_alg_id); % retrieve Matlab array - y = astra_mex_data2d('get_single', op.vol_id); + if isa(x, 'single') + y = astra_mex_data2d('get_single', op.vol_id); + else + y = astra_mex_data2d('get', op.vol_id); + end end + end % opTomo_intrnl2D @@ -225,55 +236,16 @@ classdef opTomo < opSpot function y = opTomo_intrnl3D(op,x,mode) if mode == 1 - % X is passed as a vector, reshape it into an image + % x is passed as a vector, reshape it into an image x = reshape(x, op.vol_size); - % initialize output - y = zeros(op.proj_size, 'single'); - - % link matlab array to ASTRA - vol_id = astra_mex_data3d_c('link', '-vol', op.vol_geom, x, 0); - sino_id = astra_mex_data3d_c('link', '-sino', op.proj_geom, y, 1); - - % initialize fp algorithm - cfg = astra_struct('FP3D_CUDA'); - cfg.ProjectionDataId = sino_id; - cfg.VolumeDataId = vol_id; - - alg_id = astra_mex_algorithm('create', cfg); - % forward projection - astra_mex_algorithm('iterate', alg_id); - - % cleanup - astra_mex_data3d('delete', vol_id); - astra_mex_data3d('delete', sino_id); - astra_mex_algorithm('delete', alg_id); + y = astra_mex_direct('FP3D', op.proj_id, x); else - % X is passed as a vector, reshape it into projection data + % x is passed as a vector, reshape it into projection data x = reshape(x, op.proj_size); - - % initialize output - y = zeros(op.vol_size,'single'); - - % link matlab array to ASTRA - vol_id = astra_mex_data3d_c('link', '-vol', op.vol_geom, y, 1); - sino_id = astra_mex_data3d_c('link', '-sino', op.proj_geom, x, 0); - % initialize bp algorithm - cfg = astra_struct('BP3D_CUDA'); - cfg.ProjectionDataId = sino_id; - cfg.ReconstructionDataId = vol_id; - - alg_id = astra_mex_algorithm('create', cfg); - - % backprojection - astra_mex_algorithm('iterate', alg_id); - - % cleanup - astra_mex_data3d('delete', vol_id); - astra_mex_data3d('delete', sino_id); - astra_mex_algorithm('delete', alg_id); + y = astra_mex_direct('BP3D', op.proj_id, x); end end % opTomo_intrnl3D -- cgit v1.2.3 From 6205f3fcc68be3d49fa5b47a02a65c21ffc53930 Mon Sep 17 00:00:00 2001 From: Folkert Bleichrodt Date: Mon, 14 Mar 2016 13:18:02 +0100 Subject: Removed regression from opTomo --- matlab/tools/opTomo.m | 23 +++++++++++++---------- 1 file changed, 13 insertions(+), 10 deletions(-) (limited to 'matlab') diff --git a/matlab/tools/opTomo.m b/matlab/tools/opTomo.m index 33c8565..04b3634 100644 --- a/matlab/tools/opTomo.m +++ b/matlab/tools/opTomo.m @@ -176,12 +176,23 @@ classdef opTomo < opSpot if issparse(x) x = full(x); end + + if isa(x, 'double') + isdouble = true; + x = single(x); + else + isdouble = false; + end % the multiplication y = op.funHandle(op, x, mode); % make sure output is column vector y = y(:); + + if isdouble + y = double(y); + end end % multiply @@ -206,11 +217,7 @@ classdef opTomo < opSpot astra_mex_algorithm('iterate', op.fp_alg_id); % retrieve Matlab array - if isa(x, 'single') - y = astra_mex_data2d('get_single', op.sino_id); - else - y = astra_mex_data2d('get', op.sino_id); - end + y = astra_mex_data2d('get_single', op.sino_id); else % x is passed as a vector, reshape it into a sinogram. x = reshape(x, op.proj_size); @@ -222,11 +229,7 @@ classdef opTomo < opSpot astra_mex_algorithm('iterate', op.bp_alg_id); % retrieve Matlab array - if isa(x, 'single') - y = astra_mex_data2d('get_single', op.vol_id); - else - y = astra_mex_data2d('get', op.vol_id); - end + y = astra_mex_data2d('get_single', op.vol_id); end end % opTomo_intrnl2D -- cgit v1.2.3