Merge branch 'master' into aniso

1 files changed, 26 insertions, 6 deletions

diff --git a/src/SartAlgorithm.cpp b/src/SartAlgorithm.cpp
index 9346160..8df3342 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++) {
@@ -272,9 +278,8 @@ void CSartAlgorithm::run(int _iNrIterations)
 
 	m_bShouldAbort = false;
 
-	int iIteration = 0;
-
 	// data projectors
+	CDataProjectorInterface* pFirstForwardProjector;
 	CDataProjectorInterface* pForwardProjector;
 	CDataProjectorInterface* pBackProjector;
 
@@ -286,13 +291,13 @@ 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
 		); 
 
 	// 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 +308,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<DiffFPPolicy, TotalPixelWeightPolicy>(					// 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 +344,7 @@ void CSartAlgorithm::run(int _iNrIterations)
 	}
 
 
+	ASTRA_DELETE(pFirstForwardProjector);
 	ASTRA_DELETE(pForwardProjector);
 	ASTRA_DELETE(pBackProjector);