From 495903529d473a9968c1333d5a515e3b94732f0b Mon Sep 17 00:00:00 2001
From: Willem Jan Palenstijn <Willem.Jan.Palenstijn@cwi.nl>
Date: Wed, 23 Mar 2016 15:29:43 +0100
Subject: Move CUDA algorithm initialization to its own function

---
 src/CudaReconstructionAlgorithm2D.cpp | 22 +++++++++++++++-------
 src/CudaSirtAlgorithm.cpp             | 35 +++++++++++------------------------
 2 files changed, 26 insertions(+), 31 deletions(-)

(limited to 'src')

diff --git a/src/CudaReconstructionAlgorithm2D.cpp b/src/CudaReconstructionAlgorithm2D.cpp
index 5a1910c..2798434 100644
--- a/src/CudaReconstructionAlgorithm2D.cpp
+++ b/src/CudaReconstructionAlgorithm2D.cpp
@@ -328,6 +328,20 @@ bool CCudaReconstructionAlgorithm2D::setupGeometry()
 	return true;
 }
 
+//----------------------------------------------------------------------------------------
+
+void CCudaReconstructionAlgorithm2D::initCUDAAlgorithm()
+{
+	bool ok;
+
+	ok = setupGeometry();
+	ASTRA_ASSERT(ok);
+
+	ok = m_pAlgo->allocateBuffers();
+	ASTRA_ASSERT(ok);
+}
+
+
 //----------------------------------------------------------------------------------------
 // Iterate
 void CCudaReconstructionAlgorithm2D::run(int _iNrIterations)
@@ -339,13 +353,7 @@ void CCudaReconstructionAlgorithm2D::run(int _iNrIterations)
 	const CVolumeGeometry2D& volgeom = *m_pReconstruction->getGeometry();
 
 	if (!m_bAlgoInit) {
-
-		ok = setupGeometry();
-		ASTRA_ASSERT(ok);
-
-		ok = m_pAlgo->allocateBuffers();
-		ASTRA_ASSERT(ok);
-
+		initCUDAAlgorithm();
 		m_bAlgoInit = true;
 	}
 
diff --git a/src/CudaSirtAlgorithm.cpp b/src/CudaSirtAlgorithm.cpp
index 33e381a..7beb30e 100644
--- a/src/CudaSirtAlgorithm.cpp
+++ b/src/CudaSirtAlgorithm.cpp
@@ -113,36 +113,23 @@ bool CCudaSirtAlgorithm::initialize(CProjector2D* _pProjector,
 }
 
 //----------------------------------------------------------------------------------------
-// Iterate
-void CCudaSirtAlgorithm::run(int _iNrIterations)
-{
-	// check initialized
-	ASTRA_ASSERT(m_bIsInitialized);
 
-	if (!m_bAlgoInit) {
-		// We only override the initialisation step to copy the min/max masks
+void CCudaSirtAlgorithm::initCUDAAlgorithm()
+{
+	CCudaReconstructionAlgorithm2D::initCUDAAlgorithm();
 
-		bool ok = setupGeometry();
-		ASTRA_ASSERT(ok);
+	astraCUDA::SIRT* pSirt = dynamic_cast<astraCUDA::SIRT*>(m_pAlgo);
 
-		ok = m_pAlgo->allocateBuffers();
+	if (m_pMinMask || m_pMaxMask) {
+		const CVolumeGeometry2D& volgeom = *m_pReconstruction->getGeometry();
+		const float *pfMinMaskData = 0;
+		const float *pfMaxMaskData = 0;
+		if (m_pMinMask) pfMinMaskData = m_pMinMask->getDataConst();
+		if (m_pMaxMask) pfMaxMaskData = m_pMaxMask->getDataConst();
+		bool ok = pSirt->uploadMinMaxMasks(pfMinMaskData, pfMaxMaskData, volgeom.getGridColCount());
 		ASTRA_ASSERT(ok);
-
-		if (m_pMinMask || m_pMaxMask) {
-			const CVolumeGeometry2D& volgeom = *m_pReconstruction->getGeometry();
-			astraCUDA::SIRT* pSirt = dynamic_cast<astraCUDA::SIRT*>(m_pAlgo);
-			const float *pfMinMaskData = 0;
-			const float *pfMaxMaskData = 0;
-			if (m_pMinMask) pfMinMaskData = m_pMinMask->getDataConst();
-			if (m_pMaxMask) pfMaxMaskData = m_pMaxMask->getDataConst();
-			ok = pSirt->uploadMinMaxMasks(pfMinMaskData, pfMaxMaskData, volgeom.getGridColCount());
-			ASTRA_ASSERT(ok);
-		}
-
-		m_bAlgoInit = true;
 	}
 
-	CCudaReconstructionAlgorithm2D::run(_iNrIterations);
 }
 
 
-- 
cgit v1.2.3


From f03ceb16d2dbde0c43e8c90683c5feafe01e5356 Mon Sep 17 00:00:00 2001
From: Willem Jan Palenstijn <Willem.Jan.Palenstijn@cwi.nl>
Date: Wed, 23 Mar 2016 15:30:47 +0100
Subject: Rename ART lambda option to Relaxation

---
 src/ArtAlgorithm.cpp | 10 +++++++---
 1 file changed, 7 insertions(+), 3 deletions(-)

(limited to 'src')

diff --git a/src/ArtAlgorithm.cpp b/src/ArtAlgorithm.cpp
index b59bd93..526c263 100644
--- a/src/ArtAlgorithm.cpp
+++ b/src/ArtAlgorithm.cpp
@@ -156,8 +156,12 @@ bool CArtAlgorithm::initialize(const Config& _cfg)
 		return false;
 	}
 
+	// "Lambda" is replaced by the more descriptive "Relaxation"
 	m_fLambda = _cfg.self.getOptionNumerical("Lambda", 1.0f);
-	CC.markOptionParsed("Lambda");
+	m_fLambda = _cfg.self.getOptionNumerical("Relaxation", m_fLambda);
+	if (!_cfg.self.hasOption("Relaxation"))
+		CC.markOptionParsed("Lambda");
+	CC.markOptionParsed("Relaxation");
 
 	// success
 	m_bIsInitialized = _check();
@@ -232,7 +236,7 @@ map<string,boost::any> CArtAlgorithm::getInformation()
 {
 	map<string, boost::any> res;
 	res["RayOrder"] = getInformation("RayOrder");
-	res["Lambda"] = getInformation("Lambda");
+	res["Relaxation"] = getInformation("Relaxation");
 	return mergeMap<string,boost::any>(CReconstructionAlgorithm2D::getInformation(), res);
 };
 
@@ -240,7 +244,7 @@ map<string,boost::any> CArtAlgorithm::getInformation()
 // Information - Specific
 boost::any CArtAlgorithm::getInformation(std::string _sIdentifier) 
 {
-	if (_sIdentifier == "Lambda")	{ return m_fLambda; }
+	if (_sIdentifier == "Relaxation")	{ return m_fLambda; }
 	if (_sIdentifier == "RayOrder") {
 		vector<float32> res;
 		for (int i = 0; i < m_iRayCount; i++) {
-- 
cgit v1.2.3


From 5edb35edc2c721b458334a65512b534912c2c542 Mon Sep 17 00:00:00 2001
From: Willem Jan Palenstijn <Willem.Jan.Palenstijn@cwi.nl>
Date: Wed, 23 Mar 2016 15:30:56 +0100
Subject: Add relaxation parameters to SIRT, SART

---
 src/CudaSartAlgorithm.cpp | 17 ++++++++++++++++-
 src/CudaSirtAlgorithm.cpp |  6 ++++++
 src/SartAlgorithm.cpp     |  8 +++++++-
 src/SirtAlgorithm.cpp     | 11 +++++++++--
 4 files changed, 38 insertions(+), 4 deletions(-)

(limited to 'src')

diff --git a/src/CudaSartAlgorithm.cpp b/src/CudaSartAlgorithm.cpp
index d202847..bf97224 100644
--- a/src/CudaSartAlgorithm.cpp
+++ b/src/CudaSartAlgorithm.cpp
@@ -107,7 +107,8 @@ bool CCudaSartAlgorithm::initialize(const Config& _cfg)
 		CC.markOptionParsed("ProjectionOrderList");
 	}
 
-
+	m_fLambda = _cfg.self.getOptionNumerical("Relaxation", 1.0f);
+	CC.markOptionParsed("Relaxation");
 
 	return true;
 }
@@ -123,12 +124,26 @@ bool CCudaSartAlgorithm::initialize(CProjector2D* _pProjector,
 	if (!m_bIsInitialized)
 		return false;
 
+	m_fLambda = 1.0f;
+
 	m_pAlgo = new astraCUDA::SART();
 	m_bAlgoInit = false;
 
 	return true;
 }
 
+//----------------------------------------------------------------------------------------
+
+void CCudaSartAlgorithm::initCUDAAlgorithm()
+{
+	CCudaReconstructionAlgorithm2D::initCUDAAlgorithm();
+
+	astraCUDA::SART* pSart = dynamic_cast<astraCUDA::SART*>(m_pAlgo);
+
+	pSart->setRelaxation(m_fLambda);
+}
+
+
 
 } // namespace astra
 
diff --git a/src/CudaSirtAlgorithm.cpp b/src/CudaSirtAlgorithm.cpp
index 7beb30e..c8dc677 100644
--- a/src/CudaSirtAlgorithm.cpp
+++ b/src/CudaSirtAlgorithm.cpp
@@ -50,6 +50,8 @@ CCudaSirtAlgorithm::CCudaSirtAlgorithm()
 
 	m_pMinMask = 0;
 	m_pMaxMask = 0;
+
+	m_fLambda = 1.0f;
 }
 
 //----------------------------------------------------------------------------------------
@@ -86,6 +88,8 @@ bool CCudaSirtAlgorithm::initialize(const Config& _cfg)
 	}
 	CC.markOptionParsed("MaxMaskId");
 
+	m_fLambda = _cfg.self.getOptionNumerical("Relaxation", 1.0f);
+	CC.markOptionParsed("Relaxation");
 
 	m_pAlgo = new astraCUDA::SIRT();
 	m_bAlgoInit = false;
@@ -108,6 +112,7 @@ bool CCudaSirtAlgorithm::initialize(CProjector2D* _pProjector,
 
 	m_pAlgo = new astraCUDA::SIRT();
 	m_bAlgoInit = false;
+	m_fLambda = 1.0f;
 
 	return true;
 }
@@ -130,6 +135,7 @@ void CCudaSirtAlgorithm::initCUDAAlgorithm()
 		ASTRA_ASSERT(ok);
 	}
 
+	pSirt->setRelaxation(m_fLambda);
 }
 
 
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<string,boost::any> CSartAlgorithm::getInformation()
 {
 	map<string, boost::any> res;
 	res["ProjectionOrder"] = getInformation("ProjectionOrder");
+	res["Relaxation"] = getInformation("Relaxation");
 	return mergeMap<string,boost::any>(CReconstructionAlgorithm2D::getInformation(), res);
 };
 
@@ -253,6 +257,8 @@ map<string,boost::any> CSartAlgorithm::getInformation()
 // Information - Specific
 boost::any CSartAlgorithm::getInformation(std::string _sIdentifier) 
 {
+	if (_sIdentifier == "Relaxation")
+		return m_fLambda;
 	if (_sIdentifier == "ProjectionOrder") {
 		vector<float32> 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
 		); 
 
diff --git a/src/SirtAlgorithm.cpp b/src/SirtAlgorithm.cpp
index d9f3a65..ff25648 100644
--- a/src/SirtAlgorithm.cpp
+++ b/src/SirtAlgorithm.cpp
@@ -76,6 +76,7 @@ void CSirtAlgorithm::_clear()
 	m_pDiffSinogram = NULL;
 	m_pTmpVolume = NULL;
 
+	m_fLambda = 1.0f;
 	m_iIterationCount = 0;
 }
 
@@ -91,6 +92,7 @@ void CSirtAlgorithm::clear()
 	ASTRA_DELETE(m_pDiffSinogram);
 	ASTRA_DELETE(m_pTmpVolume);
 
+	m_fLambda = 1.0f;
 	m_iIterationCount = 0;
 }
 
@@ -128,6 +130,9 @@ bool CSirtAlgorithm::initialize(const Config& _cfg)
 		return false;
 	}
 
+	m_fLambda = _cfg.self.getOptionNumerical("Relaxation", 1.0f);
+	CC.markOptionParsed("Relaxation");
+
 	// init data objects and data projectors
 	_init();
 
@@ -152,6 +157,8 @@ bool CSirtAlgorithm::initialize(CProjector2D* _pProjector,
 	m_pSinogram = _pSinogram;
 	m_pReconstruction = _pReconstruction;
 
+	m_fLambda = 1.0f;
+
 	// init data objects and data projectors
 	_init();
 
@@ -248,7 +255,7 @@ void CSirtAlgorithm::run(int _iNrIterations)
 			x = 1.0f / x;
 		else
 			x = 0.0f;
-		pfT[i] = x;
+		pfT[i] = m_fLambda * x;
 	}
 	pfT = m_pTotalRayLength->getData();
 	for (int i = 0; i < m_pTotalRayLength->getSize(); ++i) {
@@ -296,7 +303,7 @@ void CSirtAlgorithm::run(int _iNrIterations)
 		m_pTmpVolume->setData(0.0f);
 		pBackProjector->project();
 
-		// divide by pixel weights
+		// multiply with relaxation factor divided by pixel weights
 		(*m_pTmpVolume) *= (*m_pTotalPixelWeight);
 		(*m_pReconstruction) += (*m_pTmpVolume);
 
-- 
cgit v1.2.3


From 16430239d04ff738a21146c410918c285552543f Mon Sep 17 00:00:00 2001
From: Willem Jan Palenstijn <Willem.Jan.Palenstijn@cwi.nl>
Date: Wed, 23 Mar 2016 15:50:24 +0100
Subject: Add relaxation parameters to SIRT3D

---
 src/CudaSirtAlgorithm3D.cpp | 8 ++++++++
 1 file changed, 8 insertions(+)

(limited to 'src')

diff --git a/src/CudaSirtAlgorithm3D.cpp b/src/CudaSirtAlgorithm3D.cpp
index 605c470..c819f8e 100644
--- a/src/CudaSirtAlgorithm3D.cpp
+++ b/src/CudaSirtAlgorithm3D.cpp
@@ -56,6 +56,7 @@ CCudaSirtAlgorithm3D::CCudaSirtAlgorithm3D()
 	m_iGPUIndex = -1;
 	m_iVoxelSuperSampling = 1;
 	m_iDetectorSuperSampling = 1;
+	m_fLambda = 1.0f;
 }
 
 //----------------------------------------------------------------------------------------
@@ -128,6 +129,8 @@ bool CCudaSirtAlgorithm3D::initialize(const Config& _cfg)
 		return false;
 	}
 
+	m_fLambda = _cfg.self.getOptionNumerical("Relaxation");
+
 	initializeFromProjector();
 
 	// Deprecated options
@@ -135,6 +138,7 @@ bool CCudaSirtAlgorithm3D::initialize(const Config& _cfg)
 	m_iDetectorSuperSampling = (int)_cfg.self.getOptionNumerical("DetectorSuperSampling", m_iDetectorSuperSampling);
 	m_iGPUIndex = (int)_cfg.self.getOptionNumerical("GPUindex", m_iGPUIndex);
 	m_iGPUIndex = (int)_cfg.self.getOptionNumerical("GPUIndex", m_iGPUIndex);
+
 	CC.markOptionParsed("VoxelSuperSampling");
 	CC.markOptionParsed("DetectorSuperSampling");
 	CC.markOptionParsed("GPUIndex");
@@ -164,6 +168,8 @@ bool CCudaSirtAlgorithm3D::initialize(CProjector3D* _pProjector,
 		clear();
 	}
 
+	m_fLambda = 1.0f;
+
 	// required classes
 	m_pProjector = _pProjector;
 	m_pSinogram = _pSinogram;
@@ -224,6 +230,8 @@ void CCudaSirtAlgorithm3D::run(int _iNrIterations)
 
 		ASTRA_ASSERT(ok);
 
+		m_pSirt->setRelaxation(m_fLambda);
+
 		m_bAstraSIRTInit = true;
 
 	}
-- 
cgit v1.2.3