From 5edb35edc2c721b458334a65512b534912c2c542 Mon Sep 17 00:00:00 2001 From: Willem Jan Palenstijn Date: Wed, 23 Mar 2016 15:30:56 +0100 Subject: Add relaxation parameters to SIRT, SART --- src/SartAlgorithm.cpp | 8 +++++++- 1 file changed, 7 insertions(+), 1 deletion(-) (limited to 'src/SartAlgorithm.cpp') diff --git a/src/SartAlgorithm.cpp b/src/SartAlgorithm.cpp index 9346160..403f851 100644 --- a/src/SartAlgorithm.cpp +++ b/src/SartAlgorithm.cpp @@ -151,6 +151,9 @@ bool CSartAlgorithm::initialize(const Config& _cfg) CC.markOptionParsed("ProjectionOrderList"); } + m_fLambda = _cfg.self.getOptionNumerical("Relaxation", 1.0f); + CC.markOptionParsed("Relaxation"); + // create data objects m_pTotalRayLength = new CFloat32ProjectionData2D(m_pProjector->getProjectionGeometry()); m_pTotalPixelWeight = new CFloat32VolumeData2D(m_pProjector->getVolumeGeometry()); @@ -246,6 +249,7 @@ map CSartAlgorithm::getInformation() { map res; res["ProjectionOrder"] = getInformation("ProjectionOrder"); + res["Relaxation"] = getInformation("Relaxation"); return mergeMap(CReconstructionAlgorithm2D::getInformation(), res); }; @@ -253,6 +257,8 @@ map CSartAlgorithm::getInformation() // Information - Specific boost::any CSartAlgorithm::getInformation(std::string _sIdentifier) { + if (_sIdentifier == "Relaxation") + return m_fLambda; if (_sIdentifier == "ProjectionOrder") { vector res; for (int i = 0; i < m_iProjectionCount; i++) { @@ -286,7 +292,7 @@ void CSartAlgorithm::run(int _iNrIterations) m_pProjector, SinogramMaskPolicy(m_pSinogramMask), // sinogram mask ReconstructionMaskPolicy(m_pReconstructionMask), // reconstruction mask - SIRTBPPolicy(m_pReconstruction, m_pDiffSinogram, m_pTotalPixelWeight, m_pTotalRayLength), // SIRT backprojection + SIRTBPPolicy(m_pReconstruction, m_pDiffSinogram, m_pTotalPixelWeight, m_pTotalRayLength, m_fLambda), // SIRT backprojection m_bUseSinogramMask, m_bUseReconstructionMask, true // options on/off ); -- cgit v1.2.3 From f9cc36d3507f7cde4d20165836d65a584ced720f Mon Sep 17 00:00:00 2001 From: Willem Jan Palenstijn Date: Wed, 23 Mar 2016 18:21:42 +0100 Subject: Fix accumulating multiple raylengths in SART Thanks to @mohamedadaly for noticing. --- src/SartAlgorithm.cpp | 24 +++++++++++++++++++----- 1 file changed, 19 insertions(+), 5 deletions(-) (limited to 'src/SartAlgorithm.cpp') diff --git a/src/SartAlgorithm.cpp b/src/SartAlgorithm.cpp index 9346160..f80df61 100644 --- a/src/SartAlgorithm.cpp +++ b/src/SartAlgorithm.cpp @@ -272,9 +272,8 @@ void CSartAlgorithm::run(int _iNrIterations) m_bShouldAbort = false; - int iIteration = 0; - // data projectors + CDataProjectorInterface* pFirstForwardProjector; CDataProjectorInterface* pForwardProjector; CDataProjectorInterface* pBackProjector; @@ -292,7 +291,7 @@ void CSartAlgorithm::run(int _iNrIterations) // first time forward projection data projector, // also computes total pixel weight and total ray length - pForwardProjector = dispatchDataProjector( + pFirstForwardProjector = dispatchDataProjector( m_pProjector, SinogramMaskPolicy(m_pSinogramMask), // sinogram mask ReconstructionMaskPolicy(m_pReconstructionMask), // reconstruction mask @@ -303,16 +302,30 @@ void CSartAlgorithm::run(int _iNrIterations) m_bUseSinogramMask, m_bUseReconstructionMask, true // options on/off ); + // forward projection data projector + pForwardProjector = dispatchDataProjector( + m_pProjector, + SinogramMaskPolicy(m_pSinogramMask), // sinogram mask + ReconstructionMaskPolicy(m_pReconstructionMask), // reconstruction mask + CombinePolicy( // 2 basic operations + DiffFPPolicy(m_pReconstruction, m_pDiffSinogram, m_pSinogram), // forward projection with difference calculation + TotalPixelWeightPolicy(m_pTotalPixelWeight)), // calculate the total pixel weights + m_bUseSinogramMask, m_bUseReconstructionMask, true // options on/off + ); + // iteration loop - for (; iIteration < _iNrIterations && !m_bShouldAbort; ++iIteration) { + for (int iIteration = 0; iIteration < _iNrIterations && !m_bShouldAbort; ++iIteration) { int iProjection = m_piProjectionOrder[m_iIterationCount % m_iProjectionCount]; // forward projection and difference calculation m_pTotalPixelWeight->setData(0.0f); - pForwardProjector->projectSingleProjection(iProjection); + if (iIteration < m_iProjectionCount) + pFirstForwardProjector->projectSingleProjection(iProjection); + else + pForwardProjector->projectSingleProjection(iProjection); // backprojection pBackProjector->projectSingleProjection(iProjection); // update iteration count @@ -325,6 +338,7 @@ void CSartAlgorithm::run(int _iNrIterations) } + ASTRA_DELETE(pFirstForwardProjector); ASTRA_DELETE(pForwardProjector); ASTRA_DELETE(pBackProjector); -- cgit v1.2.3