summaryrefslogtreecommitdiffstats
diff options
context:
space:
mode:
Diffstat
-rw-r--r--build/linux/Makefile.in6
-rw-r--r--matlab/mex/astra_mex_data3d_c.cpp543
-rw-r--r--matlab/mex/mexCopyDataHelpFunctions.cpp257
-rw-r--r--matlab/mex/mexCopyDataHelpFunctions.h53
-rw-r--r--matlab/mex/mexDataManagerHelpFunctions.cpp371
-rw-r--r--matlab/mex/mexDataManagerHelpFunctions.h90
6 files changed, 827 insertions, 493 deletions
diff --git a/build/linux/Makefile.in b/build/linux/Makefile.in
index 49bc981..a7ac1f2 100644
--- a/build/linux/Makefile.in
+++ b/build/linux/Makefile.in
@@ -35,6 +35,8 @@ endif
ifeq ($(matlab),yes)
CPPFLAGS+=-I$(MATLAB_ROOT)/extern/include -DMATLAB_MEX_FILE
+CXXFLAGS+=-fopenmp
+LDFLAGS+=-fopenmp
endif
BOOST_CPPFLAGS=
@@ -196,7 +198,9 @@ TEST_OBJECTS=\
tests/test_XMLDocument.o
MATLAB_CXX_OBJECTS=\
- matlab/mex/mexHelpFunctions.o
+ matlab/mex/mexHelpFunctions.o \
+ matlab/mex/mexCopyDataHelpFunctions.o \
+ matlab/mex/mexDataManagerHelpFunctions.o
MATLAB_MEX=\
matlab/mex/astra_mex_algorithm_c.$(MEXSUFFIX) \
diff --git a/matlab/mex/astra_mex_data3d_c.cpp b/matlab/mex/astra_mex_data3d_c.cpp
index 35f3480..7c2af8f 100644
--- a/matlab/mex/astra_mex_data3d_c.cpp
+++ b/matlab/mex/astra_mex_data3d_c.cpp
@@ -32,6 +32,8 @@ $Id$
*/
#include <mex.h>
#include "mexHelpFunctions.h"
+#include "mexCopyDataHelpFunctions.h"
+#include "mexDataManagerHelpFunctions.h"
#include <list>
@@ -68,186 +70,38 @@ void astra_mex_data3d_create(int& nlhs, mxArray* plhs[], int& nrhs, const mxArra
return;
}
- string sDataType = mex_util_get_string(prhs[1]);
- CFloat32Data3DMemory* pDataObject3D = NULL;
+ const mxArray * const geometry = prhs[2];
+ const mxArray * const data = nrhs > 3 ? prhs[3] : NULL;
- if (nrhs >= 4 && !(mex_is_scalar(prhs[3]) || mxIsDouble(prhs[3]) || mxIsSingle(prhs[3]))) {
- mexErrMsgTxt("Data must be single or double.");
+ if (!checkStructs(geometry)) {
+ mexErrMsgTxt("Argument 3 is not a valid MATLAB struct.\n");
return;
}
- mwSize dims[3];
-
- // SWITCH DataType
- if (sDataType == "-vol") {
-
- // Read geometry
- if (!mxIsStruct(prhs[2])) {
- mexErrMsgTxt("Argument 3 is not a valid MATLAB struct.\n");
- }
- Config cfg;
- XMLDocument* xml = struct2XML("VolumeGeometry", prhs[2]);
- if (!xml)
- return;
- cfg.self = xml->getRootNode();
- CVolumeGeometry3D* pGeometry = new CVolumeGeometry3D();
- if (!pGeometry->initialize(cfg)) {
- mexErrMsgTxt("Geometry class not initialized. \n");
- delete pGeometry;
- delete xml;
- return;
- }
- delete xml;
-
- // If data is specified, check dimensions
- if (nrhs >= 4 && !mex_is_scalar(prhs[3])) {
- get3DMatrixDims(prhs[3], dims);
- if (pGeometry->getGridColCount() != dims[0] || pGeometry->getGridRowCount() != dims[1] || pGeometry->getGridSliceCount() != dims[2]) {
- mexErrMsgTxt("The dimensions of the data do not match those specified in the geometry. \n");
- delete pGeometry;
- return;
- }
- }
-
- // Initialize data object
- pDataObject3D = new CFloat32VolumeData3DMemory(pGeometry);
- delete pGeometry;
- }
-
- else if (sDataType == "-sino" || sDataType == "-proj3d") {
-
- // Read geometry
- if (!mxIsStruct(prhs[2])) {
- mexErrMsgTxt("Argument 3 is not a valid MATLAB struct.\n");
- }
- XMLDocument* xml = struct2XML("ProjectionGeometry", prhs[2]);
- if (!xml)
- return;
- Config cfg;
- cfg.self = xml->getRootNode();
-
- // FIXME: Change how the base class is created. (This is duplicated
- // in Projector2D.cpp.)
- std::string type = cfg.self->getAttribute("type");
- CProjectionGeometry3D* pGeometry = 0;
- if (type == "parallel3d") {
- pGeometry = new CParallelProjectionGeometry3D();
- } else if (type == "parallel3d_vec") {
- pGeometry = new CParallelVecProjectionGeometry3D();
- } else if (type == "cone") {
- pGeometry = new CConeProjectionGeometry3D();
- } else if (type == "cone_vec") {
- pGeometry = new CConeVecProjectionGeometry3D();
- } else {
- mexErrMsgTxt("Invalid geometry type.\n");
- return;
- }
-
- if (!pGeometry->initialize(cfg)) {
- mexErrMsgTxt("Geometry class not initialized. \n");
- delete pGeometry;
- delete xml;
- return;
- }
- delete xml;
-
- // If data is specified, check dimensions
- if (nrhs >= 4 && !mex_is_scalar(prhs[3])) {
- get3DMatrixDims(prhs[3], dims);
- if (pGeometry->getDetectorColCount() != dims[0] || pGeometry->getProjectionCount() != dims[1] || pGeometry->getDetectorRowCount() != dims[2]) {
- mexErrMsgTxt("The dimensions of the data do not match those specified in the geometry. \n");
- delete pGeometry;
- return;
- }
- }
-
- // Initialize data object
- pDataObject3D = new CFloat32ProjectionData3DMemory(pGeometry);
- delete pGeometry;
- }
-
- else if (sDataType == "-sinocone") {
- // Read geometry
- if (!mxIsStruct(prhs[2])) {
- mexErrMsgTxt("Argument 3 is not a valid MATLAB struct.\n");
- }
- XMLDocument* xml = struct2XML("ProjectionGeometry", prhs[2]);
- if (!xml)
- return;
- Config cfg;
- cfg.self = xml->getRootNode();
- CConeProjectionGeometry3D* pGeometry = new CConeProjectionGeometry3D();
- if (!pGeometry->initialize(cfg)) {
- mexErrMsgTxt("Geometry class not initialized. \n");
- delete xml;
- delete pGeometry;
- return;
- }
- delete xml;
- // If data is specified, check dimensions
- if (nrhs >= 4 && !mex_is_scalar(prhs[3])) {
- get3DMatrixDims(prhs[3], dims);
- if (pGeometry->getDetectorRowCount() != dims[2] || pGeometry->getProjectionCount() != dims[1] || pGeometry->getDetectorColCount() != dims[0]) {
- mexErrMsgTxt("The dimensions of the data do not match those specified in the geometry. \n");
- delete pGeometry;
- return;
- }
- }
- // Initialize data object
- pDataObject3D = new CFloat32ProjectionData3DMemory(pGeometry);
- delete pGeometry;
- }
- else {
- mexErrMsgTxt("Invalid datatype. Please specify '-vol' or '-proj3d'. \n");
+ if (data && !checkDataType(data)) {
+ mexErrMsgTxt("Data must be single or double.");
return;
}
- // Check initialization
- if (!pDataObject3D->isInitialized()) {
- mexErrMsgTxt("Couldn't initialize data object.\n");
- delete pDataObject3D;
+ const string sDataType = mex_util_get_string(prhs[1]);
+
+ // step2: Allocate data
+ CFloat32Data3DMemory* pDataObject3D =
+ allocateDataObject(sDataType, geometry, data);
+ if (!pDataObject3D) {
+ // Error message was already set by the function
return;
}
- // Store data
-
- // fill with scalar value
- if (nrhs < 4 || mex_is_scalar(prhs[3])) {
- float32 fValue = 0.0f;
- if (nrhs >= 4)
- fValue = (float32)mxGetScalar(prhs[3]);
- for (int i = 0; i < pDataObject3D->getSize(); ++i) {
- pDataObject3D->getData()[i] = fValue;
- }
- }
- // fill with array value
- else if (mxIsDouble(prhs[3])) {
- double* pdMatlabData = mxGetPr(prhs[3]);
- int i = 0;
- int col, row, slice;
- for (slice = 0; slice < dims[2]; ++slice) {
- for (row = 0; row < dims[1]; ++row) {
- for (col = 0; col < dims[0]; ++col) {
- // TODO: Benchmark and remove triple indexing?
- pDataObject3D->getData3D()[slice][row][col] = pdMatlabData[i];
- ++i;
- }
- }
- }
- }
- else if (mxIsSingle(prhs[3])) {
- const float* pfMatlabData = (const float*)mxGetData(prhs[3]);
- int i = 0;
- int col, row, slice;
- for (slice = 0; slice < dims[2]; ++slice) {
- for (row = 0; row < dims[1]; ++row) {
- for (col = 0; col < dims[0]; ++col) {
- // TODO: Benchmark and remove triple indexing?
- pDataObject3D->getData3D()[slice][row][col] = pfMatlabData[i];
- ++i;
- }
- }
- }
+ // step3: Initialize data
+ if (!data) {
+ mxArray * emptyArray = mxCreateDoubleMatrix(0, 0, mxREAL);
+ copyMexToCFloat32Array(emptyArray, pDataObject3D->getData(),
+ pDataObject3D->getSize());
+ mxDestroyArray(emptyArray);
+ } else {
+ copyMexToCFloat32Array(data, pDataObject3D->getData(),
+ pDataObject3D->getSize());
}
pDataObject3D->updateStatistics();
@@ -258,7 +112,6 @@ void astra_mex_data3d_create(int& nlhs, mxArray* plhs[], int& nrhs, const mxArra
if (1 <= nlhs) {
plhs[0] = mxCreateDoubleScalar(iIndex);
}
-
}
//-----------------------------------------------------------------------------------------
@@ -274,209 +127,38 @@ void astra_mex_data3d_create(int& nlhs, mxArray* plhs[], int& nrhs, const mxArra
*/
#ifdef USE_MATLAB_UNDOCUMENTED
-extern "C" {
-mxArray *mxCreateSharedDataCopy(const mxArray *pr);
-bool mxUnshareArray(const mxArray *pr, const bool noDeepCopy);
-mxArray *mxUnreference(const mxArray *pr);
-#if 0
-// Unsupported in Matlab R2014b
-bool mxIsSharedArray(const mxArray *pr);
-#endif
-}
-
-class CFloat32CustomMemoryMatlab3D : public CFloat32CustomMemory {
-public:
- // offset allows linking the data object to a sub-volume (in the z direction)
- // offset is measured in floats.
- CFloat32CustomMemoryMatlab3D(const mxArray* _pArray, bool bUnshare, size_t iOffset) {
- //fprintf(stderr, "Passed:\narray: %p\tdata: %p\n", (void*)_pArray, (void*)mxGetData(_pArray));
- // First unshare the input array, so that we may modify it.
- if (bUnshare) {
-#if 0
- // Unsupported in Matlab R2014b
- if (mxIsSharedArray(_pArray)) {
- fprintf(stderr, "Performance note: unsharing shared array in link\n");
- }
-#endif
- mxUnshareArray(_pArray, false);
- //fprintf(stderr, "Unshared:\narray: %p\tdata: %p\n", (void*)_pArray, (void*)mxGetData(_pArray));
- }
- // Then create a (persistent) copy so the data won't be deleted
- // or changed.
- m_pLink = mxCreateSharedDataCopy(_pArray);
- //fprintf(stderr, "SharedDataCopy:\narray: %p\tdata: %p\n", (void*)m_pLink, (void*)mxGetData(m_pLink));
- mexMakeArrayPersistent(m_pLink);
- m_fPtr = (float *)mxGetData(m_pLink);
- m_fPtr += iOffset;
- }
- virtual ~CFloat32CustomMemoryMatlab3D() {
- // destroy the shared array
- //fprintf(stderr, "Destroy:\narray: %p\tdata: %p\n", (void*)m_pLink, (void*)mxGetData(m_pLink));
- mxDestroyArray(m_pLink);
- }
-private:
- mxArray* m_pLink;
-};
void astra_mex_data3d_link(int& nlhs, mxArray* plhs[], int& nrhs, const mxArray* prhs[])
{
- // TODO: Reduce code duplication with _create!
- // TODO: Also allow empty argument to let this function create its own mxArray
+ // TODO: Allow empty argument to let this function create its own mxArray
// step1: get datatype
if (nrhs < 4) {
mexErrMsgTxt("Not enough arguments. See the help document for a detailed argument list. \n");
return;
}
- if (!mxIsStruct(prhs[2])) {
- mexErrMsgTxt("Argument 3 is not a valid MATLAB struct.\n");
- }
-
- const mxArray* _pArray = prhs[3];
-
- string sDataType = mex_util_get_string(prhs[1]);
- CFloat32Data3DMemory* pDataObject3D = NULL;
+ const mxArray * const geometry = prhs[2];
+ const mxArray * const data = nrhs > 3 ? prhs[3] : NULL;
+ const mxArray * const unshare = nrhs > 4 ? prhs[4] : NULL;
+ const mxArray * const zIndex = nrhs > 5 ? prhs[5] : NULL;
- if (mex_is_scalar(_pArray) || !mxIsSingle(_pArray)) {
- mexErrMsgTxt("Data must be a single array.");
+ if (!checkStructs(geometry)) {
+ mexErrMsgTxt("Argument 3 is not a valid MATLAB struct.\n");
return;
}
- unsigned int iZ = 0;
- bool bUnshare = true;
- if (nrhs > 4) {
- if (!mex_is_scalar(prhs[4])) {
- mexErrMsgTxt("Argument 5 (read-only) must be scalar");
- return;
- }
- // unshare the array if we're not linking read-only
- bUnshare = !(bool)mxGetScalar(prhs[4]);
- }
- if (nrhs > 5) {
- if (!mex_is_scalar(prhs[5])) {
- mexErrMsgTxt("Argument 6 (Z) must be scalar");
- return;
- }
- iZ = (unsigned int)mxGetScalar(prhs[5]);
- }
-
- mwSize dims[3];
-
- // SWITCH DataType
- if (sDataType == "-vol") {
- // Read geometry
- Config cfg;
- XMLDocument* xml = struct2XML("VolumeGeometry", prhs[2]);
- if (!xml)
- return;
- cfg.self = xml->getRootNode();
- CVolumeGeometry3D* pGeometry = new CVolumeGeometry3D();
- if (!pGeometry->initialize(cfg)) {
- mexErrMsgTxt("Geometry class not initialized. \n");
- delete pGeometry;
- delete xml;
- return;
- }
- delete xml;
-
- // If data is specified, check dimensions
- if (nrhs >= 4) {
- get3DMatrixDims(prhs[3], dims);
- bool ok = pGeometry->getGridColCount() == dims[0] && pGeometry->getGridRowCount() == dims[1];
- if (nrhs == 4) {
- ok &= pGeometry->getGridSliceCount() == dims[2];
- } else if (nrhs > 4) {
- // sub-volume in Z direction
- ok &= iZ + pGeometry->getGridSliceCount() <= dims[2];
- }
- if (!ok) {
- mexErrMsgTxt("The dimensions of the data do not match those specified in the geometry. \n");
- delete pGeometry;
- return;
- }
- }
-
- size_t iOffset = iZ;
- iOffset *= dims[0];
- iOffset *= dims[1];
- CFloat32CustomMemoryMatlab3D* pHandle = new CFloat32CustomMemoryMatlab3D(_pArray, bUnshare, iOffset);
-
- // Initialize data object
- pDataObject3D = new CFloat32VolumeData3DMemory(pGeometry, pHandle);
- delete pGeometry;
- }
- else if (sDataType == "-sino" || sDataType == "-proj3d" || sDataType == "-sinocone") {
-
- // Read geometry
- if (!mxIsStruct(prhs[2])) {
- mexErrMsgTxt("Argument 3 is not a valid MATLAB struct.\n");
- }
- XMLDocument* xml = struct2XML("ProjectionGeometry", prhs[2]);
- if (!xml)
- return;
- Config cfg;
- cfg.self = xml->getRootNode();
-
- // FIXME: Change how the base class is created. (This is duplicated
- // in Projector2D.cpp.)
- std::string type = cfg.self->getAttribute("type");
- CProjectionGeometry3D* pGeometry = 0;
- if (type == "parallel3d") {
- pGeometry = new CParallelProjectionGeometry3D();
- } else if (type == "parallel3d_vec") {
- pGeometry = new CParallelVecProjectionGeometry3D();
- } else if (type == "cone") {
- pGeometry = new CConeProjectionGeometry3D();
- } else if (type == "cone_vec") {
- pGeometry = new CConeVecProjectionGeometry3D();
- } else {
- mexErrMsgTxt("Invalid geometry type.\n");
- return;
- }
-
- if (!pGeometry->initialize(cfg)) {
- mexErrMsgTxt("Geometry class not initialized. \n");
- delete pGeometry;
- delete xml;
- return;
- }
- delete xml;
-
- // If data is specified, check dimensions
- if (nrhs >= 4) {
- get3DMatrixDims(prhs[3], dims);
- bool ok = pGeometry->getDetectorColCount() == dims[0] && pGeometry->getProjectionCount() == dims[1];
- if (nrhs == 4) {
- ok &= pGeometry->getDetectorRowCount() == dims[2];
- } else if (nrhs > 4) {
- // sub-volume in Z direction
- ok &= iZ + pGeometry->getDetectorRowCount() <= dims[2];
- }
- if (!ok) {
- mexErrMsgTxt("The dimensions of the data do not match those specified in the geometry. \n");
- delete pGeometry;
- return;
- }
- }
-
- size_t iOffset = iZ;
- iOffset *= dims[0];
- iOffset *= dims[1];
- CFloat32CustomMemoryMatlab3D* pHandle = new CFloat32CustomMemoryMatlab3D(_pArray, bUnshare, iOffset);
-
- // Initialize data object
- pDataObject3D = new CFloat32ProjectionData3DMemory(pGeometry, pHandle);
- delete pGeometry;
- }
- else {
- mexErrMsgTxt("Invalid datatype. Please specify '-vol' or '-proj3d'. \n");
+ if (data && !checkDataType(data)) {
+ mexErrMsgTxt("Data must be single or double.");
return;
}
- // Check initialization
- if (!pDataObject3D->isInitialized()) {
- mexErrMsgTxt("Couldn't initialize data object.\n");
- delete pDataObject3D;
+ string sDataType = mex_util_get_string(prhs[1]);
+
+ // step2: Allocate data
+ CFloat32Data3DMemory* pDataObject3D =
+ allocateDataObject(sDataType, geometry, data, unshare, zIndex);
+ if (!pDataObject3D) {
+ // Error message was already set by the function
return;
}
@@ -489,10 +171,7 @@ void astra_mex_data3d_link(int& nlhs, mxArray* plhs[], int& nrhs, const mxArray*
if (1 <= nlhs) {
plhs[0] = mxCreateDoubleScalar(iIndex);
}
-
}
-
-
#endif
@@ -547,43 +226,8 @@ void astra_mex_data3d_create_cache(int nlhs, mxArray* plhs[], int nrhs, const mx
*/
void astra_mex_data3d_get(int nlhs, mxArray* plhs[], int nrhs, const mxArray* prhs[])
-{
- // step1: input
- if (nrhs < 2) {
- mexErrMsgTxt("Not enough arguments. See the help document for a detailed argument list. \n");
- return;
- }
- int iDataID = (int)(mxGetScalar(prhs[1]));
-
- // step2: get data object
- CFloat32Data3DMemory* pDataObject = dynamic_cast<CFloat32Data3DMemory*>(astra::CData3DManager::getSingleton().get(iDataID));
- if (!pDataObject || !pDataObject->isInitialized()) {
- mexErrMsgTxt("Data object not found or not initialized properly.\n");
- return;
- }
-
- // create output
- if (1 <= nlhs) {
- mwSize dims[3];
- dims[0] = pDataObject->getWidth();
- dims[1] = pDataObject->getHeight();
- dims[2] = pDataObject->getDepth();
-
- plhs[0] = mxCreateNumericArray(3, dims, mxDOUBLE_CLASS, mxREAL);
- double* out = mxGetPr(plhs[0]);
-
- int i = 0;
- for (int slice = 0; slice < pDataObject->getDepth(); slice++) {
- for (int row = 0; row < pDataObject->getHeight(); row++) {
- for (int col = 0; col < pDataObject->getWidth(); col++) {
- // TODO: Benchmark and remove triple indexing?
- out[i] = pDataObject->getData3D()[slice][row][col];
- ++i;
- }
- }
- }
- }
-
+{
+ generic_astra_mex_data3d_get<mxDOUBLE_CLASS>(nlhs, plhs, nrhs, prhs);
}
//-----------------------------------------------------------------------------------------
@@ -596,43 +240,8 @@ void astra_mex_data3d_get(int nlhs, mxArray* plhs[], int nrhs, const mxArray* pr
*/
void astra_mex_data3d_get_single(int nlhs, mxArray* plhs[], int nrhs, const mxArray* prhs[])
-{
- // step1: input
- if (nrhs < 2) {
- mexErrMsgTxt("Not enough arguments. See the help document for a detailed argument list. \n");
- return;
- }
- int iDataID = (int)(mxGetScalar(prhs[1]));
-
- // step2: get data object
- CFloat32Data3DMemory* pDataObject = dynamic_cast<CFloat32Data3DMemory*>(astra::CData3DManager::getSingleton().get(iDataID));
- if (!pDataObject || !pDataObject->isInitialized()) {
- mexErrMsgTxt("Data object not found or not initialized properly.\n");
- return;
- }
-
- // create output
- if (1 <= nlhs) {
- mwSize dims[3];
- dims[0] = pDataObject->getWidth();
- dims[1] = pDataObject->getHeight();
- dims[2] = pDataObject->getDepth();
-
- plhs[0] = mxCreateNumericArray(3, dims, mxSINGLE_CLASS, mxREAL);
- float* out = (float *)mxGetData(plhs[0]);
-
- int i = 0;
- for (int slice = 0; slice < pDataObject->getDepth(); slice++) {
- for (int row = 0; row < pDataObject->getHeight(); row++) {
- for (int col = 0; col < pDataObject->getWidth(); col++) {
- // TODO: Benchmark and remove triple indexing?
- out[i] = pDataObject->getData3D()[slice][row][col];
- ++i;
- }
- }
- }
- }
-
+{
+ generic_astra_mex_data3d_get<mxSINGLE_CLASS>(nlhs, plhs, nrhs, prhs);
}
@@ -652,70 +261,20 @@ void astra_mex_data3d_store(int nlhs, mxArray* plhs[], int nrhs, const mxArray*
mexErrMsgTxt("Not enough arguments. See the help document for a detailed argument list. \n");
return;
}
- int iDataID = (int)(mxGetScalar(prhs[1]));
- // step2: get data object
- CFloat32Data3DMemory* pDataObject = dynamic_cast<CFloat32Data3DMemory*>(astra::CData3DManager::getSingleton().get(iDataID));
- if (!pDataObject || !pDataObject->isInitialized()) {
- mexErrMsgTxt("Data object not found or not initialized properly.\n");
+ if (!checkDataType(prhs[2])) {
+ mexErrMsgTxt("Data must be single or double.");
return;
}
- if (!(mex_is_scalar(prhs[2]) || mxIsDouble(prhs[2]) || mxIsSingle(prhs[2]))) {
- mexErrMsgTxt("Data must be single or double.");
+ // step2: get data object
+ CFloat32Data3DMemory* pDataObject = NULL;
+ if (!checkID(mxGetScalar(prhs[1]), pDataObject)) {
+ mexErrMsgTxt("Data object not found or not initialized properly.\n");
return;
}
- // fill with scalar value
- if (mex_is_scalar(prhs[2])) {
- float32 fValue = (float32)mxGetScalar(prhs[2]);
- for (int i = 0; i < pDataObject->getSize(); ++i) {
- pDataObject->getData()[i] = fValue;
- }
- }
- // fill with array value
- else if (mxIsDouble(prhs[2])) {
- mwSize dims[3];
- get3DMatrixDims(prhs[2], dims);
- if (dims[0] != pDataObject->getWidth() || dims[1] != pDataObject->getHeight() || dims[2] != pDataObject->getDepth()) {
- mexErrMsgTxt("Data object dimensions don't match.\n");
- return;
-
- }
- double* pdMatlabData = mxGetPr(prhs[2]);
- int i = 0;
- int col, row, slice;
- for (slice = 0; slice < dims[2]; ++slice) {
- for (row = 0; row < dims[1]; ++row) {
- for (col = 0; col < dims[0]; ++col) {
- // TODO: Benchmark and remove triple indexing?
- pDataObject->getData3D()[slice][row][col] = pdMatlabData[i];
- ++i;
- }
- }
- }
- }
- else if (mxIsSingle(prhs[2])) {
- mwSize dims[3];
- get3DMatrixDims(prhs[2], dims);
- if (dims[0] != pDataObject->getWidth() || dims[1] != pDataObject->getHeight() || dims[2] != pDataObject->getDepth()) {
- mexErrMsgTxt("Data object dimensions don't match.\n");
- return;
-
- }
- const float* pfMatlabData = (const float *)mxGetData(prhs[2]);
- int i = 0;
- int col, row, slice;
- for (slice = 0; slice < dims[2]; ++slice) {
- for (row = 0; row < dims[1]; ++row) {
- for (col = 0; col < dims[0]; ++col) {
- // TODO: Benchmark and remove triple indexing?
- pDataObject->getData3D()[slice][row][col] = pfMatlabData[i];
- ++i;
- }
- }
- }
- }
+ copyMexToCFloat32Array(prhs[2], pDataObject->getData(), pDataObject->getSize());
pDataObject->updateStatistics();
}
diff --git a/matlab/mex/mexCopyDataHelpFunctions.cpp b/matlab/mex/mexCopyDataHelpFunctions.cpp
new file mode 100644
index 0000000..bbcad30
--- /dev/null
+++ b/matlab/mex/mexCopyDataHelpFunctions.cpp
@@ -0,0 +1,257 @@
+/*
+-----------------------------------------------------------------------
+Copyright 2013 iMinds-Vision Lab, University of Antwerp
+
+Contact: astra@ua.ac.be
+Website: http://astra.ua.ac.be
+
+
+This file is part of the
+All Scale Tomographic Reconstruction Antwerp Toolbox ("ASTRA Toolbox").
+
+The ASTRA Toolbox is free software: you can redistribute it and/or modify
+it under the terms of the GNU General Public License as published by
+the Free Software Foundation, either version 3 of the License, or
+(at your option) any later version.
+
+The ASTRA Toolbox is distributed in the hope that it will be useful,
+but WITHOUT ANY WARRANTY; without even the implied warranty of
+MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+GNU General Public License for more details.
+
+You should have received a copy of the GNU General Public License
+along with the ASTRA Toolbox. If not, see <http://www.gnu.org/licenses/>.
+
+-----------------------------------------------------------------------
+$Id$
+*/
+
+#include "mexCopyDataHelpFunctions.h"
+
+#include "mexHelpFunctions.h"
+
+#define HAVE_OMP
+
+#define ROUND_DOWN(x, s) ((x) & ~((s)-1))
+
+#ifdef HAVE_OMP
+# include <omp.h>
+#endif
+
+#if defined(__SSE2__)
+# include <emmintrin.h>
+
+# define STORE_32F_64F_CORE8(in, out, count, base) \
+ {\
+ const __m128 inV0 = *((const __m128 *) &in[count + 0 + base]);\
+ const __m128 inV1 = *((const __m128 *) &in[count + 4 + base]);\
+\
+ *((__m128d *)&out[count + 0 + base]) = _mm_cvtps_pd(inV0);\
+ *((__m128d *)&out[count + 2 + base]) = _mm_cvtps_pd(_mm_movehl_ps(inV0, inV0));\
+\
+ *((__m128d *)&out[count + 4 + base]) = _mm_cvtps_pd(inV1);\
+ *((__m128d *)&out[count + 6 + base]) = _mm_cvtps_pd(_mm_movehl_ps(inV1, inV1));\
+ }
+# define STORE_64F_32F_CORE8(in, out, count, base) \
+ {\
+ const __m128d inV0 = *((const __m128d *) &in[count + 0 + base]);\
+ const __m128d inV1 = *((const __m128d *) &in[count + 2 + base]);\
+\
+ const __m128d inV2 = *((const __m128d *) &in[count + 4 + base]);\
+ const __m128d inV3 = *((const __m128d *) &in[count + 6 + base]);\
+\
+ *((__m128 *)&out[count + 0 + base]) = _mm_movelh_ps(_mm_cvtpd_ps(inV0), _mm_cvtpd_ps(inV1));\
+ *((__m128 *)&out[count + 4 + base]) = _mm_movelh_ps(_mm_cvtpd_ps(inV2), _mm_cvtpd_ps(inV3));\
+ }
+# define STORE_32F_32F_CORE8(in, out, count, base) \
+ {\
+ *((__m128 *)&out[count + 0 + base]) = *((const __m128 *)&in[count + 0 + base]);\
+ *((__m128 *)&out[count + 4 + base]) = *((const __m128 *)&in[count + 4 + base]);\
+ }
+# define STORE_CONST_32F_CORE8(val, out, count, base) \
+ {\
+ *((__m128 *)&out[count + 0 + base]) = val;\
+ *((__m128 *)&out[count + 4 + base]) = val;\
+ }
+#else
+# define STORE_32F_64F_CORE8(in, out, count, base) \
+ {\
+ out[count + 0 + base] = (double)in[count + 0 + base];\
+ out[count + 1 + base] = (double)in[count + 1 + base];\
+ out[count + 2 + base] = (double)in[count + 2 + base];\
+ out[count + 3 + base] = (double)in[count + 3 + base];\
+ out[count + 4 + base] = (double)in[count + 4 + base];\
+ out[count + 5 + base] = (double)in[count + 5 + base];\
+ out[count + 6 + base] = (double)in[count + 6 + base];\
+ out[count + 7 + base] = (double)in[count + 7 + base];\
+ }
+# define STORE_64F_32F_CORE8(in, out, count, base) \
+ {\
+ out[count + 0 + base] = (float)in[count + 0 + base];\
+ out[count + 1 + base] = (float)in[count + 1 + base];\
+ out[count + 2 + base] = (float)in[count + 2 + base];\
+ out[count + 3 + base] = (float)in[count + 3 + base];\
+ out[count + 4 + base] = (float)in[count + 4 + base];\
+ out[count + 5 + base] = (float)in[count + 5 + base];\
+ out[count + 6 + base] = (float)in[count + 6 + base];\
+ out[count + 7 + base] = (float)in[count + 7 + base];\
+ }
+# define STORE_32F_32F_CORE8(in, out, count, base) \
+ {\
+ out[count + 0 + base] = in[count + 0 + base];\
+ out[count + 1 + base] = in[count + 1 + base];\
+ out[count + 2 + base] = in[count + 2 + base];\
+ out[count + 3 + base] = in[count + 3 + base];\
+ out[count + 4 + base] = in[count + 4 + base];\
+ out[count + 5 + base] = in[count + 5 + base];\
+ out[count + 6 + base] = in[count + 6 + base];\
+ out[count + 7 + base] = in[count + 7 + base];\
+ }
+#endif
+
+const char * warnDataTypeNotSupported = "Data type not supported: nothing was copied";
+
+void
+_copyMexToCFloat32Array(const mxArray * const inArray, astra::float32 * const out)
+{
+ const long long tot_size = mxGetNumberOfElements(inArray);
+ const long long block = 32;
+ const long long totRoundedPixels = ROUND_DOWN(tot_size, block);
+
+ if (mxIsDouble(inArray)) {
+ const double * const pdMatlabData = mxGetPr(inArray);
+
+#pragma omp for nowait
+ for (long long count = 0; count < totRoundedPixels; count += block) {
+ STORE_64F_32F_CORE8(pdMatlabData, out, count, 0);
+ STORE_64F_32F_CORE8(pdMatlabData, out, count, 8);
+ STORE_64F_32F_CORE8(pdMatlabData, out, count, 16);
+ STORE_64F_32F_CORE8(pdMatlabData, out, count, 24);
+ }
+#pragma omp for nowait
+ for (long long count = totRoundedPixels; count < tot_size; count++) {
+ out[count] = pdMatlabData[count];
+ }
+ }
+ else if (mxIsSingle(inArray)) {
+ const float * const pfMatlabData = (const float *)mxGetData(inArray);
+
+#pragma omp for nowait
+ for (long long count = 0; count < totRoundedPixels; count += block) {
+ STORE_32F_32F_CORE8(pfMatlabData, out, count, 0);
+ STORE_32F_32F_CORE8(pfMatlabData, out, count, 8);
+ STORE_32F_32F_CORE8(pfMatlabData, out, count, 16);
+ STORE_32F_32F_CORE8(pfMatlabData, out, count, 24);
+ }
+#pragma omp for nowait
+ for (long long count = totRoundedPixels; count < tot_size; count++) {
+ out[count] = pfMatlabData[count];
+ }
+ }
+ else {
+#pragma omp single nowait
+ mexWarnMsgIdAndTxt("ASTRA_MEX:wrong_datatype", warnDataTypeNotSupported);
+ }
+}
+
+void
+_copyCFloat32ArrayToMex(const astra::float32 * const in, mxArray * const outArray)
+{
+ const long long tot_size = mxGetNumberOfElements(outArray);
+ const long long block = 32;
+ const long long tot_rounded_size = ROUND_DOWN(tot_size, block);
+
+ if (mxIsDouble(outArray)) {
+ double * const pdMatlabData = mxGetPr(outArray);
+
+#pragma omp for nowait
+ for (long long count = 0; count < tot_rounded_size; count += block) {
+ STORE_32F_64F_CORE8(in, pdMatlabData, count, 0);
+ STORE_32F_64F_CORE8(in, pdMatlabData, count, 8);
+ STORE_32F_64F_CORE8(in, pdMatlabData, count, 16);
+ STORE_32F_64F_CORE8(in, pdMatlabData, count, 24);
+ }
+#pragma omp for nowait
+ for (long long count = tot_rounded_size; count < tot_size; count++) {
+ pdMatlabData[count] = in[count];
+ }
+ }
+ else if (mxIsSingle(outArray)) {
+ float * const pfMatlabData = (float *) mxGetData(outArray);
+
+#pragma omp for nowait
+ for (long long count = 0; count < tot_rounded_size; count += block) {
+ STORE_32F_32F_CORE8(in, pfMatlabData, count, 0);
+ STORE_32F_32F_CORE8(in, pfMatlabData, count, 8);
+ STORE_32F_32F_CORE8(in, pfMatlabData, count, 16);
+ STORE_32F_32F_CORE8(in, pfMatlabData, count, 24);
+ }
+#pragma omp for nowait
+ for (long long count = tot_rounded_size; count < tot_size; count++) {
+ pfMatlabData[count] = in[count];
+ }
+ }
+ else {
+#pragma omp single nowait
+ mexWarnMsgIdAndTxt("ASTRA_MEX:wrong_datatype", warnDataTypeNotSupported);
+ }
+}
+
+void
+_initializeCFloat32Array(const astra::float32 & val, astra::float32 * const out,
+ const size_t & tot_size)
+{
+#ifdef __SSE2__
+ const long long block = 32;
+ const long long tot_rounded_size = ROUND_DOWN(tot_size, block);
+
+ const __m128 vecVal = _mm_set1_ps(val);
+
+ {
+#pragma omp for nowait
+ for (long long count = 0; count < tot_rounded_size; count += block) {
+ STORE_CONST_32F_CORE8(vecVal, out, count, 0);
+ STORE_CONST_32F_CORE8(vecVal, out, count, 8);
+ STORE_CONST_32F_CORE8(vecVal, out, count, 16);
+ STORE_CONST_32F_CORE8(vecVal, out, count, 24);
+ }
+#else
+ const long long tot_rounded_size = 0;
+ {
+#endif
+#pragma omp for nowait
+ for (long long count = tot_rounded_size; count < (long long)tot_size; count++) {
+ out[count] = val;
+ }
+ }
+}
+
+void
+copyMexToCFloat32Array(const mxArray * const in,
+ astra::float32 * const out, const size_t &tot_size)
+{
+#pragma omp parallel
+ {
+ // fill with scalar value
+ if (mex_is_scalar(in)) {
+ astra::float32 fValue = 0.f;
+ if (!mxIsEmpty(in)) {
+ fValue = (astra::float32)mxGetScalar(in);
+ }
+ _initializeCFloat32Array(fValue, out, tot_size);
+ }
+ // fill with array value
+ else {
+ _copyMexToCFloat32Array(in, out);
+ }
+ }
+}
+
+void
+copyCFloat32ArrayToMex(const float * const in, mxArray * const outArray)
+{
+#pragma omp parallel
+ {
+ _copyCFloat32ArrayToMex(in, outArray);
+ }
+}
diff --git a/matlab/mex/mexCopyDataHelpFunctions.h b/matlab/mex/mexCopyDataHelpFunctions.h
new file mode 100644
index 0000000..d0cf3c6
--- /dev/null
+++ b/matlab/mex/mexCopyDataHelpFunctions.h
@@ -0,0 +1,53 @@
+/*
+-----------------------------------------------------------------------
+Copyright 2013 iMinds-Vision Lab, University of Antwerp
+
+Contact: astra@ua.ac.be
+Website: http://astra.ua.ac.be
+
+
+This file is part of the
+All Scale Tomographic Reconstruction Antwerp Toolbox ("ASTRA Toolbox").
+
+The ASTRA Toolbox is free software: you can redistribute it and/or modify
+it under the terms of the GNU General Public License as published by
+the Free Software Foundation, either version 3 of the License, or
+(at your option) any later version.
+
+The ASTRA Toolbox is distributed in the hope that it will be useful,
+but WITHOUT ANY WARRANTY; without even the implied warranty of
+MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+GNU General Public License for more details.
+
+You should have received a copy of the GNU General Public License
+along with the ASTRA Toolbox. If not, see <http://www.gnu.org/licenses/>.
+
+-----------------------------------------------------------------------
+$Id$
+*/
+
+#ifndef MEXCOPYDATAHELPFUNCTIONS_H_
+#define MEXCOPYDATAHELPFUNCTIONS_H_
+
+#include <mex.h>
+
+#include "astra/Globals.h"
+
+#include <vector>
+
+void copyMexToCFloat32Array(const mxArray * const, astra::float32 * const,
+ const size_t &);
+void copyCFloat32ArrayToMex(const astra::float32 * const, mxArray * const);
+
+template<mxClassID MType, class AType>
+mxArray * createEquivMexArray(const AType * const pDataObj)
+{
+ mwSize dims[3];
+ dims[0] = pDataObj->getWidth();
+ dims[1] = pDataObj->getHeight();
+ dims[2] = pDataObj->getDepth();
+
+ return mxCreateNumericArray(3, dims, MType, mxREAL);
+}
+
+#endif /* MEXCOPYDATAHELPFUNCTIONS_H_ */
diff --git a/matlab/mex/mexDataManagerHelpFunctions.cpp b/matlab/mex/mexDataManagerHelpFunctions.cpp
new file mode 100644
index 0000000..2985a9d
--- /dev/null
+++ b/matlab/mex/mexDataManagerHelpFunctions.cpp
@@ -0,0 +1,371 @@
+/*
+-----------------------------------------------------------------------
+Copyright 2013 iMinds-Vision Lab, University of Antwerp
+
+Contact: astra@ua.ac.be
+Website: http://astra.ua.ac.be
+
+
+This file is part of the
+All Scale Tomographic Reconstruction Antwerp Toolbox ("ASTRA Toolbox").
+
+The ASTRA Toolbox is free software: you can redistribute it and/or modify
+it under the terms of the GNU General Public License as published by
+the Free Software Foundation, either version 3 of the License, or
+(at your option) any later version.
+
+The ASTRA Toolbox is distributed in the hope that it will be useful,
+but WITHOUT ANY WARRANTY; without even the implied warranty of
+MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+GNU General Public License for more details.
+
+You should have received a copy of the GNU General Public License
+along with the ASTRA Toolbox. If not, see <http://www.gnu.org/licenses/>.
+
+-----------------------------------------------------------------------
+$Id$
+*/
+
+#include "mexDataManagerHelpFunctions.h"
+
+#include "mexHelpFunctions.h"
+
+#include "astra/ParallelProjectionGeometry3D.h"
+#include "astra/ParallelVecProjectionGeometry3D.h"
+#include "astra/ConeProjectionGeometry3D.h"
+#include "astra/ConeVecProjectionGeometry3D.h"
+#include "astra/Float32VolumeData3DMemory.h"
+#include "astra/Float32ProjectionData3DMemory.h"
+
+#define USE_MATLAB_UNDOCUMENTED
+
+#ifdef USE_MATLAB_UNDOCUMENTED
+extern "C" {
+mxArray *mxCreateSharedDataCopy(const mxArray *pr);
+bool mxUnshareArray(const mxArray *pr, const bool noDeepCopy);
+mxArray *mxUnreference(const mxArray *pr);
+#if 0
+// Unsupported in Matlab R2014b
+bool mxIsSharedArray(const mxArray *pr);
+#endif
+}
+
+class CFloat32CustomMemoryMatlab3D : public astra::CFloat32CustomMemory {
+public:
+ // offset allows linking the data object to a sub-volume (in the z direction)
+ // offset is measured in floats.
+ CFloat32CustomMemoryMatlab3D(const mxArray* _pArray, bool bUnshare, size_t iOffset)
+ {
+ // Convert from slice to offset
+ mwSize dims[3];
+ get3DMatrixDims(_pArray, dims);
+ iOffset *= dims[0];
+ iOffset *= dims[1];
+
+ //fprintf(stderr, "Passed:\narray: %p\tdata: %p\n", (void*)_pArray, (void*)mxGetData(_pArray));
+ // First unshare the input array, so that we may modify it.
+ if (bUnshare) {
+#if 0
+ // Unsupported in Matlab R2014b
+ if (mxIsSharedArray(_pArray)) {
+ fprintf(stderr, "Performance note: unsharing shared array in link\n");
+ }
+#endif
+ mxUnshareArray(_pArray, false);
+ //fprintf(stderr, "Unshared:\narray: %p\tdata: %p\n", (void*)_pArray, (void*)mxGetData(_pArray));
+ }
+ // Then create a (persistent) copy so the data won't be deleted
+ // or changed.
+ m_pLink = mxCreateSharedDataCopy(_pArray);
+ //fprintf(stderr, "SharedDataCopy:\narray: %p\tdata: %p\n", (void*)m_pLink, (void*)mxGetData(m_pLink));
+ mexMakeArrayPersistent(m_pLink);
+ m_fPtr = (float *)mxGetData(m_pLink);
+ m_fPtr += iOffset;
+ }
+ virtual ~CFloat32CustomMemoryMatlab3D() {
+ // destroy the shared array
+ //fprintf(stderr, "Destroy:\narray: %p\tdata: %p\n", (void*)m_pLink, (void*)mxGetData(m_pLink));
+ mxDestroyArray(m_pLink);
+ }
+private:
+ mxArray* m_pLink;
+};
+#endif
+
+//-----------------------------------------------------------------------------------------
+bool
+checkID(const astra::int32 & id, astra::CFloat32Data3DMemory *& pDataObj)
+{
+ pDataObj = dynamic_cast<astra::CFloat32Data3DMemory *>(
+ astra::CData3DManager::getSingleton().get(id) );
+ return (pDataObj && pDataObj->isInitialized());
+}
+
+//-----------------------------------------------------------------------------------------
+bool
+checkDataType(const mxArray * const in)
+{
+ return (mex_is_scalar(in) || mxIsDouble(in) || mxIsSingle(in));
+}
+
+//-----------------------------------------------------------------------------------------
+bool
+checkStructs(const mxArray * const in)
+{
+ return mxIsStruct(in);
+}
+
+//-----------------------------------------------------------------------------------------
+bool
+checkDataSize(const mxArray * const mArray,
+ const astra::CProjectionGeometry3D * const geom)
+{
+ mwSize dims[3];
+ get3DMatrixDims(mArray, dims);
+ return (geom->getDetectorColCount() == dims[0]
+ && geom->getProjectionCount() == dims[1]
+ && geom->getDetectorRowCount() == dims[2]);
+}
+
+//-----------------------------------------------------------------------------------------
+bool
+checkDataSize(const mxArray * const mArray,
+ const astra::CVolumeGeometry3D * const geom)
+{
+ mwSize dims[3];
+ get3DMatrixDims(mArray, dims);
+ return (geom->getGridColCount() == dims[0]
+ && geom->getGridRowCount() == dims[1]
+ && geom->getGridSliceCount() == dims[2]);
+}
+
+//-----------------------------------------------------------------------------------------
+bool
+checkDataSize(const mxArray * const mArray,
+ const astra::CProjectionGeometry3D * const geom,
+ const mwIndex & zOffset)
+{
+ mwSize dims[3];
+ get3DMatrixDims(mArray, dims);
+ return (geom->getDetectorColCount() == dims[0]
+ && geom->getProjectionCount() == dims[1]
+ && (zOffset + geom->getDetectorRowCount()) <= dims[2]);
+}
+
+//-----------------------------------------------------------------------------------------
+bool
+checkDataSize(const mxArray * const mArray,
+ const astra::CVolumeGeometry3D * const geom,
+ const mwIndex & zOffset)
+{
+ mwSize dims[3];
+ get3DMatrixDims(mArray, dims);
+ return (geom->getGridColCount() == dims[0]
+ && geom->getGridRowCount() == dims[1]
+ && (zOffset + geom->getGridSliceCount()) <= dims[2]);
+}
+
+//-----------------------------------------------------------------------------------------
+void
+updateStatistics(const std::vector<astra::CFloat32Data3DMemory *> & vecIn)
+{
+ const size_t tot_size = vecIn.size();
+ for (size_t count = 0; count < tot_size; count++)
+ {
+ vecIn[count]->updateStatistics();
+ }
+}
+
+//-----------------------------------------------------------------------------------------
+void
+getDataPointers(const std::vector<astra::CFloat32Data3DMemory *> & vecIn,
+ std::vector<astra::float32 *> & vecOut)
+{
+ const size_t tot_size = vecIn.size();
+ vecOut.resize(tot_size);
+ for (size_t count = 0; count < tot_size; count++)
+ {
+ vecOut[count] = vecIn[count]->getData();
+ }
+}
+
+//-----------------------------------------------------------------------------------------
+void
+getDataSizes(const std::vector<astra::CFloat32Data3DMemory *> & vecIn,
+ std::vector<size_t> & vecOut)
+{
+ const size_t tot_size = vecIn.size();
+ vecOut.resize(tot_size);
+ for (size_t count = 0; count < tot_size; count++)
+ {
+ vecOut[count] = vecIn[count]->getSize();
+ }
+}
+
+//-----------------------------------------------------------------------------------------
+astra::CFloat32Data3DMemory *
+allocateDataObject(const std::string & sDataType,
+ const mxArray * const geometry, const mxArray * const data,
+ const mxArray * const unshare, const mxArray * const zIndex)
+{
+ astra::CFloat32Data3DMemory* pDataObject3D = NULL;
+
+ bool bUnshare = true;
+ if (unshare)
+ {
+ if (!mex_is_scalar(unshare))
+ {
+ mexErrMsgTxt("Argument 5 (read-only) must be scalar");
+ return NULL;
+ }
+ // unshare the array if we're not linking read-only
+ bUnshare = !(bool)mxGetScalar(unshare);
+ }
+
+ mwIndex iZ = 0;
+ if (zIndex)
+ {
+ if (!mex_is_scalar(zIndex))
+ {
+ mexErrMsgTxt("Argument 6 (Z) must be scalar");
+ return NULL;
+ }
+ iZ = (mwSignedIndex)mxGetScalar(zIndex);
+ }
+
+ // SWITCH DataType
+ if (sDataType == "-vol")
+ {
+ // Read geometry
+ astra::XMLDocument* xml = struct2XML("VolumeGeometry", geometry);
+ if (!xml) {
+ return NULL;
+ }
+ astra::Config cfg;
+ cfg.self = xml->getRootNode();
+
+ astra::CVolumeGeometry3D* pGeometry = new astra::CVolumeGeometry3D();
+ if (!pGeometry->initialize(cfg))
+ {
+ mexErrMsgTxt("Geometry class not initialized. \n");
+ delete pGeometry;
+ delete xml;
+ return NULL;
+ }
+ delete xml;
+
+ // If data is specified, check dimensions
+ if (data && !mex_is_scalar(data))
+ {
+ if (! (zIndex
+ ? checkDataSize(data, pGeometry, iZ)
+ : checkDataSize(data, pGeometry)) )
+ {
+ mexErrMsgTxt("The dimensions of the data do not match those specified in the geometry. \n");
+ delete pGeometry;
+ return NULL;
+ }
+ }
+
+ // Initialize data object
+#ifdef USE_MATLAB_UNDOCUMENTED
+ if (unshare) {
+ CFloat32CustomMemoryMatlab3D* pHandle =
+ new CFloat32CustomMemoryMatlab3D(data, bUnshare, iZ);
+
+ // Initialize data object
+ pDataObject3D = new astra::CFloat32VolumeData3DMemory(pGeometry, pHandle);
+ }
+ else
+ {
+ pDataObject3D = new astra::CFloat32VolumeData3DMemory(pGeometry);
+ }
+#else
+ pDataObject3D = new astra::CFloat32VolumeData3DMemory(pGeometry);
+#endif
+ delete pGeometry;
+ }
+ else if (sDataType == "-sino" || sDataType == "-proj3d" || sDataType == "-sinocone")
+ {
+ // Read geometry
+ astra::XMLDocument* xml = struct2XML("ProjectionGeometry", geometry);
+ if (!xml) {
+ return NULL;
+ }
+ astra::Config cfg;
+ cfg.self = xml->getRootNode();
+
+ // FIXME: Change how the base class is created. (This is duplicated
+ // in Projector2D.cpp.)
+ std::string type = cfg.self->getAttribute("type");
+ astra::CProjectionGeometry3D* pGeometry = 0;
+ if (type == "parallel3d") {
+ pGeometry = new astra::CParallelProjectionGeometry3D();
+ } else if (type == "parallel3d_vec") {
+ pGeometry = new astra::CParallelVecProjectionGeometry3D();
+ } else if (type == "cone") {
+ pGeometry = new astra::CConeProjectionGeometry3D();
+ } else if (type == "cone_vec") {
+ pGeometry = new astra::CConeVecProjectionGeometry3D();
+ } else {
+ mexErrMsgTxt("Invalid geometry type.\n");
+ return NULL;
+ }
+
+ if (!pGeometry->initialize(cfg)) {
+ mexErrMsgTxt("Geometry class not initialized. \n");
+ delete pGeometry;
+ delete xml;
+ return NULL;
+ }
+ delete xml;
+
+ // If data is specified, check dimensions
+ if (data && !mex_is_scalar(data))
+ {
+ if (! (zIndex
+ ? checkDataSize(data, pGeometry, iZ)
+ : checkDataSize(data, pGeometry)) )
+ {
+ mexErrMsgTxt("The dimensions of the data do not match those specified in the geometry. \n");
+ delete pGeometry;
+ return NULL;
+ }
+ }
+
+ // Initialize data object
+#ifdef USE_MATLAB_UNDOCUMENTED
+ if (unshare)
+ {
+ CFloat32CustomMemoryMatlab3D* pHandle =
+ new CFloat32CustomMemoryMatlab3D(data, bUnshare, iZ);
+
+ // Initialize data object
+ pDataObject3D = new astra::CFloat32ProjectionData3DMemory(pGeometry, pHandle);
+ }
+ else
+ {
+ pDataObject3D = new astra::CFloat32ProjectionData3DMemory(pGeometry);
+ }
+#else
+ pDataObject3D = new astra::CFloat32ProjectionData3DMemory(pGeometry);
+#endif
+ delete pGeometry;
+ }
+ else
+ {
+ mexErrMsgTxt("Invalid datatype. Please specify '-vol' or '-proj3d'. \n");
+ return NULL;
+ }
+
+ // Check initialization
+ if (!pDataObject3D->isInitialized())
+ {
+ mexErrMsgTxt("Couldn't initialize data object.\n");
+ delete pDataObject3D;
+ return NULL;
+ }
+
+ return pDataObject3D;
+}
+
diff --git a/matlab/mex/mexDataManagerHelpFunctions.h b/matlab/mex/mexDataManagerHelpFunctions.h
new file mode 100644
index 0000000..36b74d8
--- /dev/null
+++ b/matlab/mex/mexDataManagerHelpFunctions.h
@@ -0,0 +1,90 @@
+/*
+-----------------------------------------------------------------------
+Copyright 2013 iMinds-Vision Lab, University of Antwerp
+
+Contact: astra@ua.ac.be
+Website: http://astra.ua.ac.be
+
+
+This file is part of the
+All Scale Tomographic Reconstruction Antwerp Toolbox ("ASTRA Toolbox").
+
+The ASTRA Toolbox is free software: you can redistribute it and/or modify
+it under the terms of the GNU General Public License as published by
+the Free Software Foundation, either version 3 of the License, or
+(at your option) any later version.
+
+The ASTRA Toolbox is distributed in the hope that it will be useful,
+but WITHOUT ANY WARRANTY; without even the implied warranty of
+MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+GNU General Public License for more details.
+
+You should have received a copy of the GNU General Public License
+along with the ASTRA Toolbox. If not, see <http://www.gnu.org/licenses/>.
+
+-----------------------------------------------------------------------
+$Id$
+*/
+
+#ifndef MEXDATAMANAGERHELPFUNCTIONS_H_
+#define MEXDATAMANAGERHELPFUNCTIONS_H_
+
+#include <mex.h>
+
+#include "astra/Globals.h"
+#include "astra/AstraObjectManager.h"
+#include "astra/Float32Data3DMemory.h"
+#include "astra/ProjectionGeometry3D.h"
+#include "astra/VolumeGeometry3D.h"
+
+#include "mexCopyDataHelpFunctions.h"
+
+bool checkID(const astra::int32 &, astra::CFloat32Data3DMemory *&);
+
+bool checkDataType(const mxArray * const);
+bool checkStructs(const mxArray * const);
+
+bool checkDataSize(const mxArray * const, const astra::CProjectionGeometry3D * const);
+bool checkDataSize(const mxArray * const, const astra::CVolumeGeometry3D * const);
+bool checkDataSize(const mxArray * const, const astra::CProjectionGeometry3D * const,
+ const mwIndex & zOffset);
+bool checkDataSize(const mxArray * const, const astra::CVolumeGeometry3D * const,
+ const mwIndex & zOffset);
+
+void updateStatistics(const std::vector<astra::CFloat32Data3DMemory *> &);
+
+void getDataPointers(const std::vector<astra::CFloat32Data3DMemory *> &,
+ std::vector<astra::float32 *> &);
+void getDataSizes(const std::vector<astra::CFloat32Data3DMemory *> &,
+ std::vector<size_t> &);
+
+astra::CFloat32Data3DMemory * allocateDataObject(const std::string & sDataType,
+ const mxArray * const geometry, const mxArray * const data,
+ const mxArray * const unshare = NULL, const mxArray * const zOffset = NULL);
+
+//-----------------------------------------------------------------------------------------
+template<mxClassID datatype>
+void generic_astra_mex_data3d_get(int nlhs, mxArray* plhs[], int nrhs,
+ const mxArray* prhs[])
+{
+ // step1: input
+ if (nrhs < 2) {
+ mexErrMsgTxt("Not enough arguments. See the help document for a detailed argument list. \n");
+ return;
+ }
+
+ // step2: get data object/s
+ astra::CFloat32Data3DMemory* pDataObject = NULL;
+ if (!checkID(mxGetScalar(prhs[1]), pDataObject)) {
+ mexErrMsgTxt("Data object not found or not initialized properly.\n");
+ return;
+ }
+
+ // create output
+ if (1 <= nlhs) {
+ plhs[0] = createEquivMexArray<datatype>(pDataObject);
+ copyCFloat32ArrayToMex(pDataObject->getData(), plhs[0]);
+ }
+}
+
+#endif /* MEXDATAMANAGERHELPFUNCTIONS_H_ */
generated by cgit v1.2.3 (git 2.25.1) at 2024-11-12 01:10:41 +0000