From 6f00d2f9f3c7b58b26ba4973ca36ac5f21723a1d Mon Sep 17 00:00:00 2001
From: Michiel Cottaar <MichielCottaar@protonmail.com>
Date: Mon, 16 Aug 2021 13:16:19 +0100
Subject: [PATCH] Made the dtifit tests deterministic

The random element of running GPS is taken out by using precomputed gradient orientations.
The old implementation set a random subset of bvecs to b=0 (or to b=2000).
This is now fixed ensuring that the full bvecs are used for each shell.
---
 unit_tests/fdt/dtifit/bvecs    | 100 +++++++++++++++++++++++++++++++++
 unit_tests/fdt/dtifit/feedsRun |  31 ++++++----
 2 files changed, 120 insertions(+), 11 deletions(-)
 create mode 100644 unit_tests/fdt/dtifit/bvecs

diff --git a/unit_tests/fdt/dtifit/bvecs b/unit_tests/fdt/dtifit/bvecs
new file mode 100644
index 0000000..01bd597
--- /dev/null
+++ b/unit_tests/fdt/dtifit/bvecs
@@ -0,0 +1,100 @@
+1  0  0  
+-0.5759015626  -0.7557308713  -0.3117823605  
+0.365980027  -0.7618208505  0.5344975318  
+0.8266615166  -0.1676571833  0.5371422586  
+-0.05315328551  0.9970703725  0.05500364094  
+-0.892798637  -0.417852473  -0.1682554742  
+0.3087651681  -0.9424873375  -0.127990975  
+0.7332427031  0.4098806082  -0.5425431092  
+-0.1825675796  -0.5151226918  -0.8374471275  
+0.09477098885  0.2595515517  0.9610678705  
+0.4475465206  -0.1412286865  0.8830382608  
+0.2530429794  -0.9593392114  0.1250501023  
+0.7354607419  0.5792546654  -0.3515132  
+-0.5011088611  0.8646571561  -0.03546705109  
+-0.3494413307  -0.8410040687  -0.4130410546  
+0.5847739281  -0.673936579  0.451496335  
+-0.02897174538  -0.6915161468  -0.721779784  
+-0.4770470052  0.5391935114  0.6940435952  
+0.4255848778  0.5245407073  0.7373835896  
+0.009732249087  0.9778708631  -0.2089829141  
+-0.744303127  -0.5609345609  -0.3624434763  
+0.9602209556  -0.2697422875  0.07221367413  
+0.1393337386  0.9597085901  0.2440195309  
+0.7522079351  0.1906474061  0.6307430451  
+-0.3567824679  0.728888815  0.5843178638  
+-0.7286370555  0.6829724737  0.05134823829  
+-0.5349212614  -0.4088516965  0.7393913269  
+0.529750921  0.7762409725  -0.3417805061  
+-0.569654883  -0.3145360974  -0.7593157167  
+0.6217432613  -0.06371685415  -0.7806250569  
+0.3453124252  0.7834611175  -0.5166701137  
+-0.8535268378  -0.5129379489  0.09157836941  
+-0.5587028258  -0.6151797057  0.5562419277  
+0.8519057956  0.3209637223  0.4138101066  
+0.7595775019  -0.2538503571  -0.5988338791  
+-0.07088030006  -0.2600860223  0.9629803965  
+0.1887815711  -0.9086574383  0.3724287559  
+0.3395785924  0.3083633785  0.8885934989  
+0.7939510312  -0.004256524139  -0.60796681  
+0.9655440027  0.09666497466  -0.2416209045  
+-0.4870237311  -0.1400514288  0.8620867025  
+-0.6162092547  0.2793353638  -0.7363816326  
+-0.8894556659  0.3184317103  -0.3278259663  
+0.07781834695  0.8865281032  -0.4560835748  
+0.5172584768  0.09327887008  0.850730698  
+-0.2525751629  -0.07009201098  0.9650351792  
+-0.1623217266  0.2030117623  0.9656282315  
+-0.9618908747  0.08459748658  -0.2600177119  
+0.4390180431  -0.8494682801  0.2926889115  
+-0.5639861908  0.696073956  0.4442979005  
+0.5206750927  0.6613037007  0.5399767247  
+-0.3735930616  -0.9135138536  -0.160998955  
+0.9657303691  -0.1908806804  -0.1758676206  
+0.2130723856  -0.6134846857  -0.7604187655  
+-0.05575164709  0.6279947277  -0.776217995  
+0.2012065496  0.9790943217  -0.02983678162  
+0.6831070125  -0.488981839  -0.5424588192  
+0.5649523258  -0.3142452495  -0.7629408842  
+0.7622387032  0.6365113381  0.1176667999  
+0.7470292555  -0.5960070066  -0.2944875881  
+-0.1310894929  0.4209333665  -0.8975692986  
+0.9622650553  0.2686836542  -0.04307037652  
+0.1061646581  -0.9424636145  -0.3170037867  
+-0.1291652717  0.8215637452  0.555292126  
+-0.5834669993  -0.8093404121  -0.06733764226  
+-0.1192284364  0.7942416415  -0.5957892201  
+0.8518057006  -0.5121976914  0.1099116614  
+-0.6837200427  -0.2064174643  0.6999419502  
+0.6851120126  0.7141287529  -0.1436720377  
+0.6585541815  -0.7279225646  0.1908798837  
+0.1957313245  -0.1849897966  0.9630514128  
+-0.06212628167  0.4530442329  0.8893206667  
+-0.2695387498  -0.06120002388  -0.9610428812  
+-0.3833288505  0.383118637  -0.8404041304  
+0.9090382753  -0.09865802712  -0.4048654194  
+-0.4058148335  0.1340505374  0.9040712219  
+-0.8626199411  -0.1899531769  0.4688332622  
+0.1080390675  0.863494933  0.4926500387  
+0.5327300523  -0.5300131134  0.659761162  
+0.3265434473  0.3324700658  -0.8847785216  
+-0.8871676959  -0.06216033205  -0.4572412628  
+0.3023205804  -0.6078512357  0.7342473302  
+-0.8859213785  0.4451050057  0.1304792892  
+0.7313288468  -0.4074607226  0.5469313278  
+0.3546148915  0.5889792358  -0.7261898777  
+-0.1470996957  -0.7322905532  0.6649152015  
+0.3241330104  -0.3892888966  -0.8622018015  
+-4.417301429e-10  4.004270839e-06  1  
+0.6761752544  0.4542177365  0.5800631631  
+-0.1433995357  -0.4987553831  0.8547980118  
+0.884545539  0.363778833  -0.2919660085  
+0.8552783312  -0.3584224946  -0.3742088875  
+0.9626357675  0.1700108002  0.2107811826  
+0.4405758073  0.88966168  -0.1199793859  
+0.7666600818  -0.5445902196  0.3400791258  
+0.2683263442  0.9167243162  -0.2960025356  
+0.2733316952  0.7152873417  0.6431592363  
+0.5569962534  -0.8043659644  -0.2067621074  
+-0.3508023035  0.8620202068  0.3658673351  
+0.6952687261  -0.03853288858  0.7177162496  
diff --git a/unit_tests/fdt/dtifit/feedsRun b/unit_tests/fdt/dtifit/feedsRun
index bea1b60..a0b3c18 100755
--- a/unit_tests/fdt/dtifit/feedsRun
+++ b/unit_tests/fdt/dtifit/feedsRun
@@ -56,24 +56,33 @@ def gen_data():
 
     diffusion_tensor = (eigen_values[..., None] * eigen_vectors[..., None, :] * eigen_vectors[..., :, None]).sum(-3)
 
+    partial_bvecs = np.genfromtxt('bvecs')
+
     for multi_shell in (False, True):  # we should get the same results for single or multi-shell data
         directory[1] = 'multi' if multi_shell else 'single'
+
+        if multi_shell:
+            bvals = np.zeros(210)
+            bvals[10:110] = 1000.
+            bvals[110:] = 2000.
+            bvecs = np.concatenate((
+                np.zeros((10, 3)),
+                partial_bvecs,
+                partial_bvecs,
+            ), axis=0)
+        else:
+            bvals = np.zeros(110)
+            bvals[10:] = 1000.
+            bvecs = np.concatenate((
+                np.zeros((10, 3)),
+                partial_bvecs,
+            ), axis=0)
+
         for kurt in (0, 1, 2):  # 0: no kurtosis; 1: single kurtosis; 2: parallel and perpendicular kurtosis
             if kurt != 0 and not multi_shell:
                 continue
             directory[2] = ['', 'kurt', 'kurtdir'][kurt]
 
-            bvals = np.full(50, 1000.)
-            if multi_shell:
-                bvals[25:] = 2000.
-            bvals[::10] = 0.
-
-            bvecs_fn = os.path.join(OUTDIR, "bvecs")
-            run(['gps', '--ndir=50', '--out={}'.format(bvecs_fn)], check=True)
-            bvecs = np.genfromtxt(bvecs_fn)
-            assert bvecs.shape == (50, 3), f"GPS produced bvecs-file with shape {bvecs.shape} rather than the expected (50, 3)"
-            testing.assert_allclose((bvecs ** 2).sum(-1), 1.)
-
             if kurt == 2:
                 kvec = np.array([0.1, 0.05, 0.2])
                 beig = (
-- 
GitLab