[RELEASE] v2.0.2

CHANGELOG.md

@@ -1,6 +1,14 @@
 # Change Log
 ### All notable changes to `COMMIT` will be documented in this file.
 
+## `v2.0.2`<br>_2023-12-06_
+### 🐛Fixed
+- Merged PR [#127](https://github.com/daducci/COMMIT/pull/127): [FIX] Operator is not recompiled when some parameters are changed
+- Merged PR [#129](https://github.com/daducci/COMMIT/pull/129): [FIX] Fista update
+
+---
+---
+
 ## `v2.0.1`<br>_2023-10-31_
 ### 🛠️Changed
 - Bump dependencies

commit/core.pyx

@@ -1,6 +1,5 @@
 #!python
 #cython: language_level=3, boundscheck=False, wraparound=False, nonecheck=False, cdivision=True, initializedcheck=False, binding=False
-from __future__ import print_function
 cimport cython
 import numpy as np
 cimport numpy as np
@@ -16,6 +15,7 @@ import commit.models
 import commit.solvers
 import amico.scheme
 import amico.lut
+from importlib import reload, invalidate_caches
 import pyximport
 from pkg_resources import get_distribution
 
@@ -170,7 +170,7 @@ cdef class Evaluation :
             ERROR( 'Scheme does not match with input data' )
         if self.scheme.dwi_count == 0 :
             ERROR( 'There are no DWI volumes in the data' )
-        
+
         # Check for Nan or Inf values in raw data
         if np.isnan(self.niiDWI_img).any() or np.isinf(self.niiDWI_img).any():
             if replace_bad_voxels is not None:
@@ -728,14 +728,14 @@ cdef class Evaluation :
             config.nISO       = self.KERNELS['iso'].shape[0]
             config.build_dir  = build_dir
 
+            sys.dont_write_bytecode = True
             pyximport.install( reload_support=True, language_level=3, build_dir=build_dir, build_in_temp=True, inplace=False )
 
-            if not 'commit.operator.operator' in sys.modules :
-                import commit.operator.operator
-            else :
-                reload( sys.modules['commit.operator.operator'] )
+        if 'commit.operator.operator' in sys.modules :
+            del sys.modules['commit.operator.operator']
+        import commit.operator.operator
 
-        self.A = sys.modules['commit.operator.operator'].LinearOperator( self.DICTIONARY, self.KERNELS, self.THREADS )
+        self.A = commit.operator.operator.LinearOperator( self.DICTIONARY, self.KERNELS, self.THREADS )
 
         LOG( '   [ %.1f seconds ]' % ( time.time() - tic ) )
 
@@ -1127,7 +1127,7 @@ cdef class Evaluation :
                 temp_weights[temp_weights>0] = xic[self.DICTIONARY['TRK']['kept']>0]
                 unravel_weights[ordered_idx] = temp_weights
                 xic = unravel_weights
-                
+
         else:
             if dictionary_info['blur_gauss_extent'] > 0 or dictionary_info['blur_core_extent'] > 0:
                 xic[ self.DICTIONARY['TRK']['kept']==1 ] *= self.DICTIONARY['TRK']['lenTot'] / self.DICTIONARY['TRK']['len']

commit/info.py

@@ -3,7 +3,7 @@
 # Format version as expected by setup.py (string of form "X.Y.Z")
 _version_major = 2
 _version_minor = 0
-_version_micro = 1
+_version_micro = 2
 _version_extra = '' #'.dev'
 __version__    = "%s.%s.%s%s" % (_version_major,_version_minor,_version_micro,_version_extra)
 

commit/proximals.pyx

@@ -1,11 +1,5 @@
 #!python
 #cython: language_level=3, boundscheck=False, wraparound=False, profile=False
-"""
-Author: Matteo Frigo - lts5 @ EPFL and Dep. of CS @ Univ. of Verona
-
-This structure is based on the previous work of Rafael Carrillo and was
-supported by the LTS5 laboratory at EPFL, Lausanne.
-"""
 cimport cython
 import numpy as np
 cimport numpy as np
@@ -32,10 +26,16 @@ cpdef soft_thresholding(double [::1] x, double lam, int compartment_start, int c
     cdef:
         int i
     for i in xrange(compartment_start, compartment_start+compartment_size):
-        if x[i] <= lam:
-            x[i] = 0.0
-        else:
+        # if x[i] <= lam:
+        #     x[i] = 0.0
+        # else:
+        #     x[i] = x[i] - lam
+        if x[i] > lam:
             x[i] = x[i] - lam
+        elif x[i] < -lam:
+            x[i] = x[i] + lam
+        else:
+            x[i] = 0.0
     return np.asarray( x )
 
 
@@ -60,7 +60,7 @@ cpdef projection_onto_l2_ball(double [::1] x, double lam, int compartment_start,
     return np.asarray( x )
 
 
-cpdef omega_group_sparsity(double [::1] x, int [::1] group_idx, int [::1] group_size, double [::1] group_weight, double lam, double n) :
+cpdef omega_group_sparsity(double [::1] x, int [::1] group_idx, int [::1] group_size, double [::1] group_weight, double lam) :
     """
     References:
         [1] Jenatton et al. - `Proximal Methods for Hierarchical Sparse Coding`
@@ -71,29 +71,18 @@ cpdef omega_group_sparsity(double [::1] x, int [::1] group_idx, int [::1] group_
         double omega = 0.0, gNorm, x_i
 
     if lam != 0:
-        if n == 2:
-            for k in xrange(nG):
-                N = group_size[k]
-                gNorm = 0.0
-                for i in xrange(j,j+N) :
-                    x_i = x[group_idx[i]]
-                    gNorm += x_i*x_i
-                omega += group_weight[k] * sqrt( gNorm )
-                j += N
-        elif n == np.inf:
-            for k in xrange(nG):
-                N = group_size[k]
-                gNorm = x[group_idx[j]]
-                for i in xrange(j+1,j+N) :
-                    x_i = x[group_idx[i]]
-                    if x_i > gNorm :
-                        gNorm = x_i
-                omega += group_weight[k] * gNorm
-                j += N
+        for k in xrange(nG):
+            N = group_size[k]
+            gNorm = 0.0
+            for i in xrange(j,j+N) :
+                x_i = x[group_idx[i]]
+                gNorm += x_i*x_i
+            omega += group_weight[k] * sqrt( gNorm )
+            j += N
     return lam*omega
 
 
-cpdef prox_group_sparsity( double [::1] x, int [::1] group_idx, int [::1] group_size, double [::1] group_weight, double lam, double n ) :
+cpdef prox_group_sparsity( double [::1] x, int [::1] group_idx, int [::1] group_size, double [::1] group_weight, double lam) :
     """
     References:
         [1] Jenatton et al. - `Proximal Methods for Hierarchical Sparse Coding`
@@ -109,33 +98,21 @@ cpdef prox_group_sparsity( double [::1] x, int [::1] group_idx, int [::1] group_
             x[i] = 0.0
 
     if lam != 0:
-        if n == 2 :
-            for k in xrange(nG) :
-                N = group_size[k]
-                gNorm = 0.0
-                for i in xrange(j,j+N) :
-                    x_i = x[group_idx[i]]
-                    gNorm += x_i*x_i
-                gNorm = sqrt( gNorm )
+        for k in xrange(nG) :
+            N = group_size[k]
+            gNorm = 0.0
+            for i in xrange(j,j+N) :
+                x_i = x[group_idx[i]]
+                gNorm += x_i*x_i
+            gNorm = sqrt( gNorm )
 
-                wl = group_weight[k] * lam
-                if gNorm <= wl :
-                    for i in xrange(j,j+N) :
-                        x[ group_idx[i] ] = 0.0
-                else :
-                    wl = (gNorm-wl)/gNorm
-                    for i in xrange(j,j+N) :
-                        x[ group_idx[i] ] *= wl
-                j += N
-        # elif n == np.inf :
-        # [TODO] TO be correctly implemented
-        #     for k in range(nG) :
-        #         idx = subtree[k]
-        #         # xn = max( v[idx] )
-        #         r = weight[k] * lam
-        #         for i in idx :
-        #             if v[i] <= r:
-        #                 v[i] = 0.0
-        #             else :
-        #                 v[i] -= r
+            wl = group_weight[k] * lam
+            if gNorm <= wl :
+                for i in xrange(j,j+N) :
+                    x[ group_idx[i] ] = 0.0
+            else :
+                wl = (gNorm-wl)/gNorm
+                for i in xrange(j,j+N) :
+                    x[ group_idx[i] ] *= wl
+            j += N
     return np.asarray( x )