diff options
-rw-r--r-- | build/linux/Makefile.in | 1 | ||||
-rw-r--r-- | build/msvc/gen.py | 2 | ||||
-rw-r--r-- | include/astra/ParallelBeamDistanceDrivenProjector2D.h | 205 | ||||
-rw-r--r-- | include/astra/ParallelBeamDistanceDrivenProjector2D.inl | 227 | ||||
-rw-r--r-- | include/astra/Projector2DImpl.inl | 1 | ||||
-rw-r--r-- | include/astra/ProjectorTypelist.h | 7 | ||||
-rw-r--r-- | src/ParallelBeamDistanceDrivenProjector2D.cpp | 174 |
7 files changed, 615 insertions, 2 deletions
diff --git a/build/linux/Makefile.in b/build/linux/Makefile.in index 7c4089a..078a1a2 100644 --- a/build/linux/Makefile.in +++ b/build/linux/Makefile.in @@ -168,6 +168,7 @@ BASE_OBJECTS=\ src/Globals.lo \ src/Logging.lo \ src/ParallelBeamBlobKernelProjector2D.lo \ + src/ParallelBeamDistanceDrivenProjector2D.lo \ src/ParallelBeamLinearKernelProjector2D.lo \ src/ParallelBeamLineKernelProjector2D.lo \ src/ParallelBeamStripKernelProjector2D.lo \ diff --git a/build/msvc/gen.py b/build/msvc/gen.py index 42d0e0e..2200ecc 100644 --- a/build/msvc/gen.py +++ b/build/msvc/gen.py @@ -249,6 +249,7 @@ P_astra["filters"]["Projectors\\source"] = [ "src\\FanFlatBeamLineKernelProjector2D.cpp", "src\\FanFlatBeamStripKernelProjector2D.cpp", "src\\ParallelBeamBlobKernelProjector2D.cpp", +"src\\ParallelBeamDistanceDrivenProjector2D.cpp", "src\\ParallelBeamLinearKernelProjector2D.cpp", "src\\ParallelBeamLineKernelProjector2D.cpp", "src\\ParallelBeamStripKernelProjector2D.cpp", @@ -393,6 +394,7 @@ P_astra["filters"]["Projectors\\headers"] = [ "include\\astra\\FanFlatBeamLineKernelProjector2D.h", "include\\astra\\FanFlatBeamStripKernelProjector2D.h", "include\\astra\\ParallelBeamBlobKernelProjector2D.h", +"include\\astra\\ParallelBeamDistanceDrivenProjector2D.h", "include\\astra\\ParallelBeamLinearKernelProjector2D.h", "include\\astra\\ParallelBeamLineKernelProjector2D.h", "include\\astra\\ParallelBeamStripKernelProjector2D.h", diff --git a/include/astra/ParallelBeamDistanceDrivenProjector2D.h b/include/astra/ParallelBeamDistanceDrivenProjector2D.h new file mode 100644 index 0000000..67dd21b --- /dev/null +++ b/include/astra/ParallelBeamDistanceDrivenProjector2D.h @@ -0,0 +1,205 @@ +/* +----------------------------------------------------------------------- +Copyright: 2010-2018, imec Vision Lab, University of Antwerp + 2014-2018, CWI, Amsterdam + +Contact: astra@astra-toolbox.com +Website: http://www.astra-toolbox.com/ + +This file is part of the 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/>. + +----------------------------------------------------------------------- +*/ + +#ifndef _INC_ASTRA_PARALLELDISTANCEDRIVENPROJECTOR +#define _INC_ASTRA_PARALLELDISTANCEDRIVENPROJECTOR + +#include "ParallelProjectionGeometry2D.h" +#include "ParallelVecProjectionGeometry2D.h" +#include "Float32Data2D.h" +#include "Projector2D.h" + +#include "Float32ProjectionData2D.h" +#include "Float32VolumeData2D.h" + +namespace astra +{ + +/** This class implements a "distance driven" two-dimensional projector. + * + * Reference: + * De Man, Bruno, and Samit Basu. “Distance-Driven Projection and Backprojection in Three Dimensions.” Physics in Medicine and Biology 49, no. 11 (2004): 2463. + * + * \par XML Configuration + * \astra_xml_item{ProjectionGeometry, xml node, The geometry of the projection.} + * \astra_xml_item{VolumeGeometry, xml node, The geometry of the volume.} + * + * \par MATLAB example + * \astra_code{ + * cfg = astra_struct('distance_driven');\n + * cfg.ProjectionGeometry = proj_geom;\n + * cfg.VolumeGeometry = vol_geom;\n + * proj_id = astra_mex_projector('create'\, cfg);\n + * } + */ +class _AstraExport CParallelBeamDistanceDrivenProjector2D : public CProjector2D { + +protected: + + /** Initial clearing. Only to be used by constructors. + */ + virtual void _clear(); + + /** Check the values of this object. If everything is ok, the object can be set to the initialized state. + * The following statements are then guaranteed to hold: + * - no NULL pointers + * - all sub-objects are initialized properly + * - blobvalues are ok + */ + virtual bool _check(); + +public: + + // type of the projector, needed to register with CProjectorFactory + static std::string type; + + /** Default constructor. + */ + CParallelBeamDistanceDrivenProjector2D(); + + /** Constructor. + * + * @param _pProjectionGeometry Information class about the geometry of the projection. Will be HARDCOPIED. + * @param _pReconstructionGeometry Information class about the geometry of the reconstruction volume. Will be HARDCOPIED. + */ + CParallelBeamDistanceDrivenProjector2D(CParallelProjectionGeometry2D* _pProjectionGeometry, + CVolumeGeometry2D* _pReconstructionGeometry); + + /** Destructor, is virtual to show that we are aware subclass destructor are called. + */ + ~CParallelBeamDistanceDrivenProjector2D(); + + /** Initialize the projector with a config object. + * + * @param _cfg Configuration Object + * @return initialization successful? + */ + virtual bool initialize(const Config& _cfg); + + /** Initialize the projector. + * + * @param _pProjectionGeometry Information class about the geometry of the projection. Will be HARDCOPIED. + * @param _pReconstructionGeometry Information class about the geometry of the reconstruction volume. Will be HARDCOPIED. + * @return initialization successful? + */ + virtual bool initialize(CParallelProjectionGeometry2D* _pProjectionGeometry, + CVolumeGeometry2D* _pVolumeGeometry); + + /** Clear this class. + */ + virtual void clear(); + + /** Returns the number of weights required for storage of all weights of one projection. + * + * @param _iProjectionIndex Index of the projection (zero-based). + * @return Size of buffer (given in SPixelWeight elements) needed to store weighted pixels. + */ + virtual int getProjectionWeightsCount(int _iProjectionIndex); + + /** Compute the pixel weights for a single ray, from the source to a detector pixel. + * + * @param _iProjectionIndex Index of the projection + * @param _iDetectorIndex Index of the detector pixel + * @param _pWeightedPixels Pointer to a pre-allocated array, consisting of _iMaxPixelCount elements + * of type SPixelWeight. On return, this array contains a list of the index + * and weight for all pixels on the ray. + * @param _iMaxPixelCount Maximum number of pixels (and corresponding weights) that can be stored in _pWeightedPixels. + * This number MUST be greater than the total number of pixels on the ray. + * @param _iStoredPixelCount On return, this variable contains the total number of pixels on the + * ray (that have been stored in the list _pWeightedPixels). + */ + virtual void computeSingleRayWeights(int _iProjectionIndex, + int _iDetectorIndex, + SPixelWeight* _pWeightedPixels, + int _iMaxPixelCount, + int& _iStoredPixelCount); + + /** Create a list of detectors that are influenced by point [_iRow, _iCol]. + * + * @param _iRow row of the point + * @param _iCol column of the point + * @return list of SDetector2D structs + */ + virtual std::vector<SDetector2D> projectPoint(int _iRow, int _iCol); + + + /** Policy-based projection of all rays. This function will calculate each non-zero projection + * weight and use this value for a task provided by the policy object. + * + * @param _policy Policy object. Should contain prior, addWeight and posterior function. + */ + template <typename Policy> + void project(Policy& _policy); + + /** Policy-based projection of all rays of a single projection. This function will calculate + * each non-zero projection weight and use this value for a task provided by the policy object. + * + * @param _iProjection Which projection should be projected? + * @param _policy Policy object. Should contain prior, addWeight and posterior function. + */ + template <typename Policy> + void projectSingleProjection(int _iProjection, Policy& _policy); + + /** Policy-based projection of a single ray. This function will calculate each non-zero + * projection weight and use this value for a task provided by the policy object. + * + * @param _iProjection Which projection should be projected? + * @param _iDetector Which detector should be projected? + * @param _policy Policy object. Should contain prior, addWeight and posterior function. + */ + template <typename Policy> + void projectSingleRay(int _iProjection, int _iDetector, Policy& _policy); + + /** Return the type of this projector. + * + * @return identification type of this projector + */ + virtual std::string getType(); + + +protected: + /** Internal policy-based projection of a range of angles and range. + * (_i*From is inclusive, _i*To exclusive) */ + template <typename Policy> + void projectBlock_internal(int _iProjFrom, int _iProjTo, + int _iDetFrom, int _iDetTo, Policy& _policy); + +}; + +//---------------------------------------------------------------------------------------- + +inline std::string CParallelBeamDistanceDrivenProjector2D::getType() +{ + return type; +} + + +} // namespace astra + + +#endif + diff --git a/include/astra/ParallelBeamDistanceDrivenProjector2D.inl b/include/astra/ParallelBeamDistanceDrivenProjector2D.inl new file mode 100644 index 0000000..6bf3b56 --- /dev/null +++ b/include/astra/ParallelBeamDistanceDrivenProjector2D.inl @@ -0,0 +1,227 @@ +/* +----------------------------------------------------------------------- +Copyright: 2010-2018, imec Vision Lab, University of Antwerp + 2014-2018, CWI, Amsterdam + +Contact: astra@astra-toolbox.com +Website: http://www.astra-toolbox.com/ + +This file is part of the 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/>. + +----------------------------------------------------------------------- +*/ + +#define policy_weight(p,rayindex,volindex,weight) do { if (p.pixelPrior(volindex)) { p.addWeight(rayindex, volindex, weight); p.pixelPosterior(volindex); } } while (false) + +template <typename Policy> +void CParallelBeamDistanceDrivenProjector2D::project(Policy& p) +{ + projectBlock_internal(0, m_pProjectionGeometry->getProjectionAngleCount(), + 0, m_pProjectionGeometry->getDetectorCount(), p); +} + +template <typename Policy> +void CParallelBeamDistanceDrivenProjector2D::projectSingleProjection(int _iProjection, Policy& p) +{ + projectBlock_internal(_iProjection, _iProjection + 1, + 0, m_pProjectionGeometry->getDetectorCount(), p); +} + +template <typename Policy> +void CParallelBeamDistanceDrivenProjector2D::projectSingleRay(int _iProjection, int _iDetector, Policy& p) +{ + projectBlock_internal(_iProjection, _iProjection + 1, + _iDetector, _iDetector + 1, p); +} + + + + + +template <typename Policy> +void CParallelBeamDistanceDrivenProjector2D::projectBlock_internal(int _iProjFrom, int _iProjTo, int _iDetFrom, int _iDetTo, Policy& p) +{ + // get vector geometry + const CParallelVecProjectionGeometry2D* pVecProjectionGeometry; + if (dynamic_cast<CParallelProjectionGeometry2D*>(m_pProjectionGeometry)) { + pVecProjectionGeometry = dynamic_cast<CParallelProjectionGeometry2D*>(m_pProjectionGeometry)->toVectorGeometry(); + } else { + pVecProjectionGeometry = dynamic_cast<CParallelVecProjectionGeometry2D*>(m_pProjectionGeometry); + } + + // precomputations + const float32 pixelLengthX = m_pVolumeGeometry->getPixelLengthX(); + const float32 pixelLengthY = m_pVolumeGeometry->getPixelLengthY(); + const float32 inv_pixelLengthX = 1.0f / pixelLengthX; + const float32 inv_pixelLengthY = 1.0f / pixelLengthY; + const int colCount = m_pVolumeGeometry->getGridColCount(); + const int rowCount = m_pVolumeGeometry->getGridRowCount(); + + // Performance note: + // This is not a very well optimizated version of the distance driven + // projector. The CPU projector model in ASTRA requires ray-driven iteration, + // which limits re-use of intermediate computations. + + // loop angles + for (int iAngle = _iProjFrom; iAngle < _iProjTo; ++iAngle) { + + const SParProjection * proj = &pVecProjectionGeometry->getProjectionVectors()[iAngle]; + + float32 detSize = sqrt(proj->fDetUX * proj->fDetUX + proj->fDetUY * proj->fDetUY); + + const bool vertical = fabs(proj->fRayX) < fabs(proj->fRayY); + + const float32 Ex = m_pVolumeGeometry->getWindowMinX() + pixelLengthX*0.5f; + const float32 Ey = m_pVolumeGeometry->getWindowMaxY() - pixelLengthY*0.5f; + + // loop detectors + for (int iDetector = _iDetFrom; iDetector < _iDetTo; ++iDetector) { + + const int iRayIndex = iAngle * m_pProjectionGeometry->getDetectorCount() + iDetector; + + // POLICY: RAY PRIOR + if (!p.rayPrior(iRayIndex)) continue; + + const float32 Dx = proj->fDetSX + (iDetector+0.5f) * proj->fDetUX; + const float32 Dy = proj->fDetSY + (iDetector+0.5f) * proj->fDetUY; + + if (vertical && true) { + + const float32 RxOverRy = proj->fRayX/proj->fRayY; + // TODO: Determine det/pixel scaling factors + const float32 lengthPerRow = m_pVolumeGeometry->getPixelLengthX() * m_pVolumeGeometry->getPixelLengthY(); + const float32 deltac = -pixelLengthY * RxOverRy * inv_pixelLengthX; + const float32 deltad = fabs((proj->fDetUX - proj->fDetUY * RxOverRy) * inv_pixelLengthX); + + // calculate c for row 0 + float32 c = (Dx + (Ey - Dy)*RxOverRy - Ex) * inv_pixelLengthX + 0.5f; + + // loop rows + for (int row = 0; row < rowCount; ++row, c+= deltac) { + + // horizontal extent of ray in center of this row: + // [ c - deltad/2 , c + deltad/2 ] + + // |-gapBegin-*---|------|----*-gapEnd-| + // * = ray extent intercepts; c - deltad/2 and c + deltad/2 + // | = pixel column edges + + const int colBegin = (int)floor(c - deltad/2.0f); + const int colEnd = (int)ceil(c + deltad/2.0f); + + // TODO: Optimize volume edge checks + + int iVolumeIndex = row * colCount + colBegin; + if (colBegin + 1 == colEnd) { + + if (colBegin >= 0 && colBegin < colCount) + policy_weight(p, iRayIndex, iVolumeIndex, + deltad * lengthPerRow); + } else { + const float gapBegin = (c - deltad/2.0f) - (float32)colBegin; + const float gapEnd = (float32)colEnd - (c + deltad/2.0f); + float tot = 1.0f - gapBegin; + if (colBegin >= 0 && colBegin < colCount) { + policy_weight(p, iRayIndex, iVolumeIndex, + (1.0f - gapBegin) * lengthPerRow); + } + iVolumeIndex++; + + for (int col = colBegin + 1; col + 1 < colEnd; ++col, ++iVolumeIndex) { + tot += 1.0f; + if (col >= 0 && col < colCount) { + policy_weight(p, iRayIndex, iVolumeIndex, lengthPerRow); + } + } + assert(iVolumeIndex == row * colCount + colEnd - 1); + tot += 1.0f - gapEnd; + if (colEnd > 0 && colEnd <= colCount) { + policy_weight(p, iRayIndex, iVolumeIndex, + (1.0f - gapEnd) * lengthPerRow); + } + assert(fabs(tot - deltad) < 0.0001); + } + + } + + } else if (!vertical && true) { + + const float32 RyOverRx = proj->fRayY/proj->fRayX; + // TODO: Determine det/pixel scaling factors + const float32 lengthPerCol = m_pVolumeGeometry->getPixelLengthX() * m_pVolumeGeometry->getPixelLengthY(); + const float32 deltar = -pixelLengthX * RyOverRx * inv_pixelLengthY; + const float32 deltad = fabs((proj->fDetUY - proj->fDetUX * RyOverRx) * inv_pixelLengthY); + + // calculate r for col 0 + float32 r = -(Dy + (Ex - Dx)*RyOverRx - Ey) * inv_pixelLengthY + 0.5f; + + // loop columns + for (int col = 0; col < colCount; ++col, r+= deltar) { + + // vertical extent of ray in center of this column: + // [ r - deltad/2 , r + deltad/2 ] + + const int rowBegin = (int)floor(r - deltad/2.0f); + const int rowEnd = (int)ceil(r + deltad/2.0f); + + // TODO: Optimize volume edge checks + + int iVolumeIndex = rowBegin * colCount + col; + if (rowBegin + 1 == rowEnd) { + + if (rowBegin >= 0 && rowBegin < rowCount) + policy_weight(p, iRayIndex, iVolumeIndex, + deltad * lengthPerCol); + } else { + const float gapBegin = (r - deltad/2.0f) - (float32)rowBegin; + const float gapEnd = (float32)rowEnd - (r + deltad/2.0f); + float tot = 1.0f - gapBegin; + + if (rowBegin >= 0 && rowBegin < rowCount) { + policy_weight(p, iRayIndex, iVolumeIndex, + (1.0f - gapBegin) * lengthPerCol); + } + iVolumeIndex += colCount; + + for (int row = rowBegin + 1; row + 1 < rowEnd; ++row, iVolumeIndex += colCount) { + tot += 1.0f; + if (row >= 0 && row < rowCount) { + policy_weight(p, iRayIndex, iVolumeIndex, lengthPerCol); + } + } + assert(iVolumeIndex == (rowEnd - 1) * colCount + col); + tot += 1.0f - gapEnd; + if (rowEnd > 0 && rowEnd <= rowCount) { + policy_weight(p, iRayIndex, iVolumeIndex, + (1.0f - gapEnd) * lengthPerCol); + } + assert(fabs(tot - deltad) < 0.0001); + } + + } + + } + + // POLICY: RAY POSTERIOR + p.rayPosterior(iRayIndex); + + } + } + + if (dynamic_cast<CParallelProjectionGeometry2D*>(m_pProjectionGeometry)) + delete pVecProjectionGeometry; +} diff --git a/include/astra/Projector2DImpl.inl b/include/astra/Projector2DImpl.inl index b7edb1f..1308141 100644 --- a/include/astra/Projector2DImpl.inl +++ b/include/astra/Projector2DImpl.inl @@ -26,6 +26,7 @@ along with the ASTRA Toolbox. If not, see <http://www.gnu.org/licenses/>. */ +#include "ParallelBeamDistanceDrivenProjector2D.inl" #include "ParallelBeamLinearKernelProjector2D.inl" #include "ParallelBeamLineKernelProjector2D.inl" #include "ParallelBeamStripKernelProjector2D.inl" diff --git a/include/astra/ProjectorTypelist.h b/include/astra/ProjectorTypelist.h index 1a65aad..063693d 100644 --- a/include/astra/ProjectorTypelist.h +++ b/include/astra/ProjectorTypelist.h @@ -35,6 +35,7 @@ along with the ASTRA Toolbox. If not, see <http://www.gnu.org/licenses/>. #include "Projector2D.h" #include "ParallelBeamLineKernelProjector2D.h" #include "ParallelBeamLinearKernelProjector2D.h" +#include "ParallelBeamDistanceDrivenProjector2D.h" #include "ParallelBeamBlobKernelProjector2D.h" #include "ParallelBeamStripKernelProjector2D.h" #include "SparseMatrixProjector2D.h" @@ -46,9 +47,10 @@ namespace astra{ #ifdef ASTRA_CUDA - typedef TYPELIST_8( + typedef TYPELIST_9( CFanFlatBeamLineKernelProjector2D, CFanFlatBeamStripKernelProjector2D, + CParallelBeamDistanceDrivenProjector2D, CParallelBeamLinearKernelProjector2D, CParallelBeamLineKernelProjector2D, CParallelBeamBlobKernelProjector2D, @@ -59,9 +61,10 @@ namespace astra{ #else - typedef TYPELIST_7( + typedef TYPELIST_8( CFanFlatBeamLineKernelProjector2D, CFanFlatBeamStripKernelProjector2D, + CParallelBeamDistanceDrivenProjector2D, CParallelBeamLinearKernelProjector2D, CParallelBeamLineKernelProjector2D, CParallelBeamBlobKernelProjector2D, diff --git a/src/ParallelBeamDistanceDrivenProjector2D.cpp b/src/ParallelBeamDistanceDrivenProjector2D.cpp new file mode 100644 index 0000000..a4ac450 --- /dev/null +++ b/src/ParallelBeamDistanceDrivenProjector2D.cpp @@ -0,0 +1,174 @@ +/* +----------------------------------------------------------------------- +Copyright: 2010-2018, imec Vision Lab, University of Antwerp + 2014-2018, CWI, Amsterdam + +Contact: astra@astra-toolbox.com +Website: http://www.astra-toolbox.com/ + +This file is part of the 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/>. + +----------------------------------------------------------------------- +*/ + +#include "astra/ParallelBeamDistanceDrivenProjector2D.h" + +#include <cmath> +#include <algorithm> + +#include "astra/DataProjectorPolicies.h" + +using namespace std; +using namespace astra; + +#include "astra/ParallelBeamDistanceDrivenProjector2D.inl" + +// type of the projector, needed to register with CProjectorFactory +std::string CParallelBeamDistanceDrivenProjector2D::type = "distance_driven"; + +//---------------------------------------------------------------------------------------- +// default constructor +CParallelBeamDistanceDrivenProjector2D::CParallelBeamDistanceDrivenProjector2D() +{ + _clear(); +} + + +//---------------------------------------------------------------------------------------- +// constructor +CParallelBeamDistanceDrivenProjector2D::CParallelBeamDistanceDrivenProjector2D(CParallelProjectionGeometry2D* _pProjectionGeometry, + CVolumeGeometry2D* _pReconstructionGeometry) + +{ + _clear(); + initialize(_pProjectionGeometry, _pReconstructionGeometry); +} + +//---------------------------------------------------------------------------------------- +// destructor +CParallelBeamDistanceDrivenProjector2D::~CParallelBeamDistanceDrivenProjector2D() +{ + clear(); +} + +//--------------------------------------------------------------------------------------- +// Clear - CParallelBeamDistanceDrivenProjector2D +void CParallelBeamDistanceDrivenProjector2D::_clear() +{ + CProjector2D::_clear(); + m_bIsInitialized = false; +} + +//--------------------------------------------------------------------------------------- +// Clear - Public +void CParallelBeamDistanceDrivenProjector2D::clear() +{ + CProjector2D::clear(); + m_bIsInitialized = false; +} + +//--------------------------------------------------------------------------------------- +// Check +bool CParallelBeamDistanceDrivenProjector2D::_check() +{ + // check base class + ASTRA_CONFIG_CHECK(CProjector2D::_check(), "ParallelBeamDistanceDrivenProjector2D", "Error in Projector2D initialization"); + + ASTRA_CONFIG_CHECK(dynamic_cast<CParallelProjectionGeometry2D*>(m_pProjectionGeometry) || dynamic_cast<CParallelVecProjectionGeometry2D*>(m_pProjectionGeometry), "ParallelBeamDistanceDrivenProjector2D", "Unsupported projection geometry"); + + ASTRA_CONFIG_CHECK(abs(m_pVolumeGeometry->getPixelLengthX() / m_pVolumeGeometry->getPixelLengthY()) - 1 < eps, "ParallelBeamDistanceDrivenProjector2D", "Pixel height must equal pixel width."); + + // success + return true; +} + +//--------------------------------------------------------------------------------------- +// Initialize, use a Config object +bool CParallelBeamDistanceDrivenProjector2D::initialize(const Config& _cfg) +{ + ASTRA_ASSERT(_cfg.self); + + // if already initialized, clear first + if (m_bIsInitialized) { + clear(); + } + + // initialization of parent class + if (!CProjector2D::initialize(_cfg)) { + return false; + } + + // success + m_bIsInitialized = _check(); + return m_bIsInitialized; +} + +//--------------------------------------------------------------------------------------- +// Initialize +bool CParallelBeamDistanceDrivenProjector2D::initialize(CParallelProjectionGeometry2D* _pProjectionGeometry, + CVolumeGeometry2D* _pVolumeGeometry) +{ + // if already initialized, clear first + if (m_bIsInitialized) { + clear(); + } + + // hardcopy geometries + m_pProjectionGeometry = _pProjectionGeometry->clone(); + m_pVolumeGeometry = _pVolumeGeometry->clone(); + + // success + m_bIsInitialized = _check(); + return m_bIsInitialized; +} + +//---------------------------------------------------------------------------------------- +// Get maximum amount of weights on a single ray +int CParallelBeamDistanceDrivenProjector2D::getProjectionWeightsCount(int _iProjectionIndex) +{ + int maxDim = max(m_pVolumeGeometry->getGridRowCount(), m_pVolumeGeometry->getGridColCount()); + int scale = m_pProjectionGeometry->getDetectorWidth() / min(m_pVolumeGeometry->getPixelLengthX(), m_pVolumeGeometry->getPixelLengthY()); + return maxDim * scale * 10 + 1; +} + +//---------------------------------------------------------------------------------------- +// Single Ray Weights +void CParallelBeamDistanceDrivenProjector2D::computeSingleRayWeights(int _iProjectionIndex, + int _iDetectorIndex, + SPixelWeight* _pWeightedPixels, + int _iMaxPixelCount, + int& _iStoredPixelCount) +{ + ASTRA_ASSERT(m_bIsInitialized); + StorePixelWeightsPolicy p(_pWeightedPixels, _iMaxPixelCount); + projectSingleRay(_iProjectionIndex, _iDetectorIndex, p); + _iStoredPixelCount = p.getStoredPixelCount(); +} + +//---------------------------------------------------------------------------------------- +// Splat a single point +std::vector<SDetector2D> CParallelBeamDistanceDrivenProjector2D::projectPoint(int _iRow, int _iCol) +{ + // TODO + assert(false); + + std::vector<SDetector2D> res; + + // return result vector + return res; + +} |