use the proper configuration fista by default

mrirecon · Oct 22, 2024 · 7cded5c · 7cded5c
1 parent 11d81fb
commit 7cded5c
Show file tree

Hide file tree

Showing 4 changed files with 23 additions and 10 deletions.
diff --git a/src/iter/italgos.c b/src/iter/italgos.c
@@ -202,14 +202,13 @@ void ist(int maxiter, float epsilon, float tau, long N,
 		if (NULL != ist_continuation)
 			ist_continuation(&itrdata);
 
-		iter_op_p_call(thresh, itrdata.scale * itrdata.tau, x, x);
-
-
 		iter_op_call(op, r, x);		// r = A x
 		vops->xpay(N, -1., r, b);	// r = b - r = b - A x
 
 		itrdata.rsnew = vops->norm(N, r);
 
+		iter_op_p_call(thresh, itrdata.scale * itrdata.tau, x, x);
+
 		debug_printf(DP_DEBUG3, "#It %03d: %f \n", itrdata.iter, itrdata.rsnew / itrdata.rsnot);
 
 		if (itrdata.rsnew < epsilon)
@@ -292,7 +291,7 @@ void fista(int maxiter, float epsilon, float tau,
 		vops->axpy(N, x, itrdata.tau, r);
 	}
 
-	if (last) {
+	if (!last) {
 
 		iter_monitor(monitor, vops, x);
 		iter_op_p_call(thresh, itrdata.scale * itrdata.tau, x, x);

diff --git a/src/pics.c b/src/pics.c
@@ -327,7 +327,7 @@ int main_pics(int argc, char* argv[argc])
 		OPTL_STRING(0, "wavelet", &wtype_str, "name", "wavelet type (haar,dau2,cdf44)"),
 		OPTL_ULONG(0, "mpi", &mpi_flags, "flags", "distribute over this dimensions with use of MPI"),
 		OPTL_FLVEC3(0, "fista_pqr", &fista.params, "p:q:r", "parameters for FISTA acceleration"),
-		OPTL_SET(0, "fista_last", &fista.last, "end iteration with call to proximal op"),
+		OPTL_SET(0, "fista_last", &fista.last, "end iteration with call to data consistency"),
 		OPTL_INFILE(0, "motion-field", &motion_file, "file", "motion field"),
 	};
 

diff --git a/tests/pics.mk b/tests/pics.mk
@@ -598,6 +598,20 @@ tests/test-pics-eulermaruyama3: ones scale zeros pics var nrmse
 	touch $@
 
 
+tests/test-pics-fista: phantom upat squeeze fmac pics nrmse
+	set -e; mkdir $(TESTS_TMP) ; cd $(TESTS_TMP)						;\
+	$(TOOLDIR)/phantom -k -s8 k.ra								;\
+	$(TOOLDIR)/phantom -S8 s.ra								;\
+	$(TOOLDIR)/upat -y 2 p.ra								;\
+	$(TOOLDIR)/squeeze p.ra p2.ra								;\
+	$(TOOLDIR)/fmac k.ra p2.ra kp.ra							;\
+	$(TOOLDIR)/pics -w1. -i100 -l2 -r100000000. kp.ra s.ra x.ra				;\
+	$(TOOLDIR)/pics -w1. -i200 -S -e -l2 -r100000000. --fista kp.ra s.ra xI.ra		;\
+	$(TOOLDIR)/nrmse -t 0.001 x.ra xI.ra							;\
+	rm *.ra ; cd .. ; rmdir $(TESTS_TMP)
+	touch $@
+
+
 tests/test-pics-ist: phantom upat squeeze fmac pics nrmse
 	set -e; mkdir $(TESTS_TMP) ; cd $(TESTS_TMP)						;\
 	$(TOOLDIR)/phantom -k -s8 k.ra								;\
@@ -606,8 +620,8 @@ tests/test-pics-ist: phantom upat squeeze fmac pics nrmse
 	$(TOOLDIR)/squeeze p.ra p2.ra								;\
 	$(TOOLDIR)/fmac k.ra p2.ra kp.ra							;\
 	$(TOOLDIR)/pics -w1. -i100 -l2 -r100000000. kp.ra s.ra x.ra				;\
-	$(TOOLDIR)/pics -w1. -i500 -S -e -l2 -r100000000. -I kp.ra s.ra xI.ra			;\
-	$(TOOLDIR)/nrmse -t 0.006 x.ra xI.ra							;\
+	$(TOOLDIR)/pics -w1. -i400 -S -e -l2 -r100000000. --ist kp.ra s.ra xI.ra		;\
+	$(TOOLDIR)/nrmse -t 0.001 x.ra xI.ra							;\
 	rm *.ra ; cd .. ; rmdir $(TESTS_TMP)
 	touch $@
 
@@ -624,6 +638,6 @@ TESTS += tests/test-pics-wavl1-dau2 tests/test-pics-wavl1-cdf44 tests/test-pics-
 TESTS += tests/test-pics-noncart-lowmem tests/test-pics-noncart-lowmem-stack0 tests/test-pics-noncart-lowmem-stack1 tests/test-pics-noncart-lowmem-stack2 tests/test-pics-noncart-lowmem-no-toeplitz
 TESTS += tests/test-pics-phase
 TESTS += tests/test-pics-eulermaruyama tests/test-pics-eulermaruyama2 tests/test-pics-eulermaruyama3
-TESTS += tests/test-pics-ist
+TESTS += tests/test-pics-fista tests/test-pics-ist
 
 
diff --git a/utests/test_iter.c b/utests/test_iter.c
@@ -251,7 +251,7 @@ static bool test_iter_irgnm_l1(void)
 
 	nlop_apply(nlid, N, dims, dst1, N, dims, src1);
 
-	md_zfill(N, dims, src2, 0.975);
+	md_zfill(N, dims, src2, 0.5);
 	md_zfill(N, dims, src3, 0.);
 
 	const struct operator_p_s* lsqr = NULL;
@@ -299,7 +299,7 @@ static bool test_iter_irgnm_l1(void)
 	md_free(src2);
 	md_free(src3);
 
-	UT_RETURN_ASSERT(err < 1.E-10);
+	UT_RETURN_ASSERT(err < 1.E-7);
 }