From 0759d7e2d1151ffcb03ba33fa0d0b50361b74295 Mon Sep 17 00:00:00 2001 From: Martin Uecker Date: Fri, 20 Sep 2024 13:15:53 +0200 Subject: [PATCH] change scaling and command line for T2 moba to be consistent with mobafit --- src/moba.c | 16 ++++++++++++---- src/moba/T2fun.c | 20 ++++++++------------ src/moba/recon.c | 2 +- tests/moba.mk | 10 +++++----- 4 files changed, 26 insertions(+), 22 deletions(-) diff --git a/src/moba.c b/src/moba.c index 915a70ab7..0ad1eaa5f 100644 --- a/src/moba.c +++ b/src/moba.c @@ -195,13 +195,16 @@ int main_moba(int argc, char* argv[argc]) opt_reg_init(&ropts); + bool t2_old_flag = false; + const struct opt_s opts[] = { //FIXME: Sort options into optimization and others interface { 'r', NULL, true, OPT_SPECIAL, opt_reg_moba, &ropts, ":A:B:C", "generalized regularization options (-rh for help)" }, OPT_SELECT('L', enum mdb_t, &conf.mode, MDB_T1, "T1 mapping using model-based look-locker"), OPT_SELECT('P', enum mdb_t, &conf.mode, MDB_T1_PHY, "T1 mapping using reparameterized (M0, R1, alpha) model-based look-locker (TR required!)"), - OPT_SELECT('F', enum mdb_t, &conf.mode, MDB_T2, "T2 mapping using model-based Fast Spin Echo"), + OPT_SET('F', &t2_old_flag, "(T2 mapping using model-based Fast Spin Echo)"), + OPT_SELECT('T', enum mdb_t, &conf.mode, MDB_T2, "T2 mapping using model-based Fast Spin Echo"), OPT_SELECT('G', enum mdb_t, &conf.mode, MDB_MGRE, "T2* mapping using model-based multiple gradient echo"), OPT_SELECT('D', enum mdb_t, &conf.mode, MDB_IR_MGRE, "Joint T1 and T2* mapping using model-based IR multiple gradient echo"), OPTL_SELECT(0, "bloch", enum mdb_t, &conf.mode, MDB_BLOCH, "Bloch model-based reconstruction"), @@ -209,7 +212,6 @@ int main_moba(int argc, char* argv[argc]) OPT_PINT('l', &conf.opt_reg, "\b1/-l2", " toggle l1-wavelet or l2 regularization."), // extra spaces needed because of backsapce \b earlier OPT_PINT('i', &conf.iter, "iter", "Number of Newton steps"), OPTL_FLOAT('R', "reduction", &conf.redu, "redu", "reduction factor"), - OPT_FLOAT('T', &conf.damping, "damp", "damping on temporal frames"), OPT_FLOAT('j', &conf.alpha_min, "minreg", "Minimum regularization parameter"), OPT_FLOAT('u', &conf.rho, "rho", "ADMM rho [default: 0.01]"), OPT_PINT('C', &conf.inner_iter, "iter", "inner iterations"), @@ -230,6 +232,7 @@ int main_moba(int argc, char* argv[argc]) OPTL_PINT(0, "pusteps", &conf.pusteps, "ud", "Number of partial update steps for IRGNM"), OPTL_FLOAT(0, "ratio", &conf.ratio, "f:[0;1]", "Ratio of partial updates: ratio* + (1-ratio)*"), OPTL_FLOAT(0, "l1val", &conf.l1val, "f", "Regularization scaling of l1 wavelet (default: 1.)"), + OPTL_FLOAT(0, "temporal_damping", &conf.damping, "f", "Temporal damping factor."), OPTL_INT(0, "multi-gpu", &(conf.num_gpu), "num", "(number of gpus to use)"), OPT_INFILE('I', &init_file, "init", "File for initialization"), OPT_INFILE('t', &traj_file, "traj", "K-space trajectory"), @@ -269,6 +272,8 @@ int main_moba(int argc, char* argv[argc]) if (use_compat_to_version("v0.6.00")) conf.scaling_M0 = 2.; + if (t2_old_flag) + conf.mode = MDB_T2; if (conf.ropts->r > 0) conf.algo = ALGO_ADMM; @@ -287,7 +292,10 @@ int main_moba(int argc, char* argv[argc]) complex float* kspace_data = load_cfl(ksp_file, DIMS, ksp_dims); long TI_dims[DIMS]; - const complex float* TI = load_cfl(TI_file, DIMS, TI_dims); + complex float* TI = load_cfl(TI_file, DIMS, TI_dims); + + if (t2_old_flag) + md_zsmul(DIMS, TI_dims, TI, TI, 10.); assert(TI_dims[TE_DIM] == ksp_dims[TE_DIM]); assert(1 == ksp_dims[MAPS_DIM]); @@ -405,7 +413,7 @@ int main_moba(int argc, char* argv[argc]) md_clear(DIMS, coil_dims, sens, CFL_SIZE); } - md_zfill(DIMS, img_dims, img, 1.0); + md_zfill(DIMS, img_dims, img, 1.); complex float* k_grid_data = anon_cfl("", DIMS, grid_dims); diff --git a/src/moba/T2fun.c b/src/moba/T2fun.c index 815062f7a..16c5c1ea3 100644 --- a/src/moba/T2fun.c +++ b/src/moba/T2fun.c @@ -49,13 +49,11 @@ struct T2_s { complex float* dz; struct multiplace_array_s* TE; - - float scaling_z; }; DEF_TYPEID(T2_s); -// Calculate Model: rho .*exp(-scaling_z.*z.*TE) +// Calculate Model: rho .*exp(-z.*TE) static void T2_fun(const nlop_data_t* _data, complex float* dst, const complex float* src) { struct T2_s* data = CAST_DOWN(T2_s, _data); @@ -83,11 +81,11 @@ static void T2_fun(const nlop_data_t* _data, complex float* dst, const complex f complex float* tmp_map = md_alloc_sameplace(data->N, data->map_dims, CFL_SIZE, dst); - // -1*scaling_z.*z - md_zsmul(data->N, data->map_dims, tmp_map, data->z, -1. * data->scaling_z); + // -1*z + md_zsmul(data->N, data->map_dims, tmp_map, data->z, -1.); complex float* tmp_exp = md_alloc_sameplace(data->N, data->out_dims, CFL_SIZE, dst); - // exp(-TE.*scaling_z.*z) + // exp(-TE.*z) md_zmul2(data->N, data->out_dims, data->out_strs, tmp_exp, data->map_strs, tmp_map, data->TE_strs, multiplace_read(data->TE, dst)); md_free(tmp_map); @@ -99,13 +97,13 @@ static void T2_fun(const nlop_data_t* _data, complex float* dst, const complex f md_zsmul(data->N, data->out_dims, data->drho, tmp_exp, 1.0); // model: - // rho.*exp(-TE.*scaling_z.*z) + // rho.*exp(-TE.*z) md_zmul2(data->N, data->out_dims, data->out_strs, dst, data->map_strs, data->rho, data->out_strs, tmp_exp); - // dz: z' = -rho.*scaling_z.*TE.*exp(-TE.*scaling_z.*z) - // TE.*exp(-TE.*scaling_z.*z), + // dz: z' = -rho.*TE.*exp(-TE.*z) + // TE.*exp(-TE.*z), md_zmul2(data->N, data->out_dims, data->out_strs, tmp_exp, data->out_strs, tmp_exp, data->TE_strs, multiplace_read(data->TE, dst)); - md_zsmul(data->N, data->out_dims, tmp_exp, tmp_exp, -1. * data->scaling_z); + md_zsmul(data->N, data->out_dims, tmp_exp, tmp_exp, -1.); md_zmul2(data->N, data->out_dims, data->out_strs, data->dz, data->map_strs, data->rho, data->out_strs, tmp_exp); md_free(tmp_exp); @@ -241,7 +239,5 @@ struct nlop_s* nlop_T2_create(int N, const long map_dims[N], const long out_dims data->dz = NULL; data->TE = multiplace_move(N, TE_dims, CFL_SIZE, TE); - data->scaling_z = 10.; - return nlop_create(N, out_dims, N, in_dims, CAST_UP(PTR_PASS(data)), T2_fun, T2_der, T2_adj, NULL, NULL, T2_del); } diff --git a/src/moba/recon.c b/src/moba/recon.c index 198cddd3c..012d1acaf 100644 --- a/src/moba/recon.c +++ b/src/moba/recon.c @@ -242,7 +242,7 @@ static struct mobamod exp_create(const long dims[DIMS], const complex float* mas complex float* TE2 = md_alloc(DIMS, edims, CFL_SIZE); - md_zsmul(DIMS, edims, TE2, TE, -10.); + md_zsmul(DIMS, edims, TE2, TE, -1.); // chain T2 model diff --git a/tests/moba.mk b/tests/moba.mk index 14d5c9285..319d6b125 100644 --- a/tests/moba.mk +++ b/tests/moba.mk @@ -167,13 +167,13 @@ tests/test-moba-t2: phantom signal fmac fft ones index scale moba slice invert n $(TOOLDIR)/ones 6 16 16 1 1 1 16 psf.ra ;\ $(TOOLDIR)/index 5 16 tmp1.ra ;\ $(TOOLDIR)/scale 0.01 tmp1.ra TE.ra ;\ - $(TOOLDIR)/moba -F -i10 -f1 -C30 --scale_data=5000. --scale_psf=1000. --normalize_scaling -d4 -p psf.ra k_space.ra TE.ra reco.ra ;\ + $(TOOLDIR)/moba -T -i10 -f1 -C50 --scale_data=5000. --scale_psf=1000. --normalize_scaling -d4 -p psf.ra k_space.ra TE.ra reco.ra ;\ $(TOOLDIR)/slice 6 1 reco.ra R2.ra ;\ $(TOOLDIR)/invert R2.ra T2.ra ;\ $(TOOLDIR)/phantom -x16 -c circ.ra ;\ $(TOOLDIR)/fmac T2.ra circ.ra masked.ra ;\ - $(TOOLDIR)/scale -- 0.9 circ.ra ref.ra ;\ - $(TOOLDIR)/nrmse -t 0.0008 masked.ra ref.ra ;\ + $(TOOLDIR)/scale -- 0.09 circ.ra ref.ra ;\ + $(TOOLDIR)/nrmse -t 0.002 ref.ra masked.ra ;\ rm *.ra ; cd .. ; rmdir $(TESTS_TMP) touch $@ @@ -189,7 +189,7 @@ tests/test-moba-meco-noncart-r2s: traj scale phantom signal fmac index extract m $(TOOLDIR)/index 5 8 tmp1.ra ;\ $(TOOLDIR)/scale 1.6 tmp1.ra tmp2.ra ;\ $(TOOLDIR)/extract 5 1 8 tmp2.ra TE.ra ;\ - $(TOOLDIR)/moba -G -m3 -rQ:1 -rS:0 -rW:3:64:1 -i10 -C100 -u0.0001 -R3 -T 0.9 -o1.5 -k --kfilter-2 -t traj.ra data.ra TE.ra reco.ra ;\ + $(TOOLDIR)/moba -G -m3 -rQ:1 -rS:0 -rW:3:64:1 -i10 -C100 -u0.0001 -R3 --temporal_damping 0.9 -o1.5 -k --kfilter-2 -t traj.ra data.ra TE.ra reco.ra ;\ $(TOOLDIR)/slice 6 1 reco.ra R2S.ra ;\ $(TOOLDIR)/resize -c 0 8 1 8 R2S.ra R2S_crop.ra ;\ $(TOOLDIR)/phantom -x8 -c circ.ra ;\ @@ -211,7 +211,7 @@ tests/test-moba-meco-noncart-wfr2s: traj scale phantom signal fmac index extract $(TOOLDIR)/index 5 8 tmp1.ra ;\ $(TOOLDIR)/scale 1.6 tmp1.ra tmp2.ra ;\ $(TOOLDIR)/extract 5 1 8 tmp2.ra TE.ra ;\ - $(TOOLDIR)/moba -G -m1 -rQ:1 -rS:0 -rW:3:64:1 -i10 -C100 -u0.0001 -R3 -T0.9 -o1.5 -k --kfilter-2 -t traj.ra data.ra TE.ra reco.ra ;\ + $(TOOLDIR)/moba -G -m1 -rQ:1 -rS:0 -rW:3:64:1 -i10 -C100 -u0.0001 -R3 --temporal_damping 0.9 -o1.5 -k --kfilter-2 -t traj.ra data.ra TE.ra reco.ra ;\ $(TOOLDIR)/resize -c 0 8 1 8 reco.ra reco_crop.ra ;\ $(TOOLDIR)/slice 6 0 reco_crop.ra W.ra ;\ $(TOOLDIR)/slice 6 1 reco_crop.ra F.ra ;\