Smooth projection in SpectralIndices Index compute module

[Epta13]

When selecting the option of Smooth projection, this error appers:

---

LinAlgError Traceback (most recent call last) File ~.conda\envs\sediment\lib\site-packages\napari_sediment\spectral_indices_widget.py:1136, in SpectralIndexWidget._on_add_index_map_to_viewer(self=, event=False, force_recompute=True) 1133 pbr.set_description("Computing index") 1135 index_series = [x for key, x in self.index_collection.items() if self.index_pick_boxes[key].isChecked()] -> 1136 self.compute_selected_indices_map_and_proj([x.index_name for x in index_series], force_recompute=force_recompute) force_recompute = True index_series = [SpectralIndex(index_name='RABD620', index_type='RABD', left_band=580, right_band=640, middle_band=619, left_band_default=580, right_band_default=640, middle_band_default=619, index_map=None, index_proj=None, index_map_range=None, colormap='viridis')] self = <napari_sediment.spectral_indices_widget.SpectralIndexWidget object at 0x000002591B8B4700> 1137 for i in index_series: 1138 computed_index = self.index_collection[i.index_name].index_map

File ~.conda\envs\sediment\lib\site-packages\napari_sediment\spectral_indices_widget.py:870, in SpectralIndexWidget.compute_selected_indices_map_and_proj(self=, index_names=['RABD620'], force_recompute=True) 868 if (self.index_collection[index_name].index_map is None) or force_recompute: 869 computed_index = self.compute_and_clean_index(index_name) --> 870 proj = compute_index_projection( self = <napari_sediment.spectral_indices_widget.SpectralIndexWidget object at 0x000002591B8B4700> computed_index = <class 'numpy.ndarray'> (12155, 798) float32 self.viewer = Viewer(camera=Camera(center=(0.0, 2679.9685052061136, 661.14302400438555), zoom=0.33107879560242037, angles=(0.0, 0.0, 90.0), perspective=0.0, mouse_pan=True, mouse_zoom=True), cursor=Cursor(position=(1.0, 4055.7739011973995, 2938.546500136239), scaled=True, size=1, style=<CursorStyle.STANDARD: 'standard'>), dims=Dims(ndim=3, ndisplay=2, last_used=0, range=((0.0, 3.0, 1.0), (0.0, 12156.0, 1.0), (0.0, 798.0, 1.0)), current_step=(1, 6077, 398), order=(0, 1, 2), axis_labels=('0', '1', '2')), grid=GridCanvas(stride=1, shape=(-1, -1), enabled=False), layers=[<Image layer 'imcube' at 0x25924e55460>, <Image layer 'red' at 0x25924e1cb50>, <Image layer 'green' at 0x259280bd7f0>, <Image layer 'blue' at 0x259282895e0>, <Labels layer 'mask' at 0x25924e1ca60>, <Shapes layer 'rois' at 0x25912f6e6a0>], help='use <2> for transform, use for add rectangles, use for add ellipses, use for add lines, use for add path, use

for add polygons, use for add polygons lasso, use <4> for select vertices, use <5> for select shapes, use <2> for insert vertex, use <1> for remove vertex', status={'layer_base': 'rois', 'source_type': '', 'plugin': '', 'coordinates': ' [4056 2939]: [, ]'}, tooltip=Tooltip(visible=True, text=''), theme='dark', title='napari', mouse_over_canvas=True, mouse_move_callbacks=[], mouse_drag_callbacks=[], mouse_double_click_callbacks=[<bound method SpectralIndexWidget._add_analysis_roi of <napari_sediment.spectral_indices_widget.SpectralIndexWidget object at 0x000002591B8B4700>>], mouse_wheel_callbacks=[<function dims_scroll at 0x000002590DD88B80>], _persisted_mouse_event={}, _mouse_drag_gen={}, _mouse_wheel_gen={}, keymap={}) colmin = 398 colmax = 418 871 computed_index, self.viewer.layers['mask'].data, 872 colmin=colmin, colmax=colmax, 873 smooth_window=self.get_smoothing_window()) 874 self.index_collection[index_name].index_map = computed_index 875 self.index_collection[index_name].index_proj = proj

File ~.conda\envs\sediment\lib\site-packages\napari_sediment\spectralindex.py:251, in compute_index_projection(index_image=<class 'numpy.ndarray'> (12155, 798) float32, mask=<class 'numpy.ndarray'> (12155, 798) uint8, colmin=398, colmax=418, smooth_window=16) 248 proj = np.nanmean(index_image[:,colmin:colmax],axis=1) 250 if smooth_window is not None: --> 251 proj = savgol_filter(proj, window_length=smooth_window, polyorder=3) proj = <class 'numpy.ndarray'> (12155,) float32 smooth_window = 16 254 return proj

File ~.conda\envs\sediment\lib\site-packages\scipy\signal_savitzky_golay.py:352, in savgol_filter(x=<class 'numpy.ndarray'> (12155,) float32, window_length=16, polyorder=3, deriv=0, delta=1.0, axis=-1, mode='interp', cval=0.0) 348 # Do not pad. Instead, for the elements within window_length // 2 349 # of the ends of the sequence, use the polynomial that is fitted to 350 # the last window_length elements. 351 y = convolve1d(x, coeffs, axis=axis, mode="constant") --> 352 _fit_edges_polyfit(x, window_length, polyorder, deriv, delta, axis, y) x = <class 'numpy.ndarray'> (12155,) float32 y = <class 'numpy.ndarray'> (12155,) float32 window_length = 16 polyorder = 3 deriv = 0 delta = 1.0 axis = -1 353 else: 354 # Any mode other than 'interp' is passed on to ndimage.convolve1d. 355 y = convolve1d(x, coeffs, axis=axis, mode=mode, cval=cval)

File ~.conda\envs\sediment\lib\site-packages\scipy\signal_savitzky_golay.py:223, in _fit_edges_polyfit(x=<class 'numpy.ndarray'> (12155,) float32, window_length=16, polyorder=3, deriv=0, delta=1.0, axis=-1, y=<class 'numpy.ndarray'> (12155,) float32) 216 """ 217 Use polynomial interpolation of x at the low and high ends of the axis 218 to fill in the halflen values in y. 219 220 This function just calls _fit_edge twice, once for each end of the axis. 221 """ 222 halflen = window_length // 2 --> 223 _fit_edge(x, 0, window_length, 0, halflen, axis, halflen = 8 window_length = 16 x = <class 'numpy.ndarray'> (12155,) float32 axis = -1 polyorder = 3 deriv = 0 delta = 1.0 y = <class 'numpy.ndarray'> (12155,) float32 224 polyorder, deriv, delta, y) 225 n = x.shape[axis] 226 _fit_edge(x, n - window_length, n, n - halflen, n, axis, 227 polyorder, deriv, delta, y)

File ~.conda\envs\sediment\lib\site-packages\scipy\signal_savitzky_golay.py:193, in _fit_edge(x=<class 'numpy.ndarray'> (12155,) float32, window_start=0, window_stop=16, interp_start=0, interp_stop=8, axis=-1, polyorder=3, deriv=0, delta=1.0, y=<class 'numpy.ndarray'> (12155,) float32) 189 xx_edge = xx_edge.reshape(xx_edge.shape[0], -1) 191 # Fit the edges. poly_coeffs has shape (polyorder + 1, -1), 192 # where '-1' is the same as in xx_edge. --> 193 poly_coeffs = np.polyfit(np.arange(0, window_stop - window_start), np.polyfit = <function polyfit at 0x000002590596BEB0> window_start = 0 window_stop = 16 np = <module 'numpy' from 'C:\Users\pzahajska\.conda\envs\sediment\lib\site-packages\numpy\init.py'> window_stop - window_start = 16 xx_edge = <class 'numpy.ndarray'> (16, 1) float32 polyorder = 3 194 xx_edge, polyorder) 196 if deriv > 0: 197 poly_coeffs = _polyder(poly_coeffs, deriv)

File ~.conda\envs\sediment\lib\site-packages\numpy\lib\polynomial.py:669, in polyfit(x=<class 'numpy.ndarray'> (16,) float64, y=<class 'numpy.ndarray'> (16, 1) float32, deg=3, rcond=3.5527136788005009e-15, full=False, w=None, cov=False) 667 scale = NX.sqrt((lhs*lhs).sum(axis=0)) 668 lhs /= scale --> 669 c, resids, rank, s = lstsq(lhs, rhs, rcond) lhs = <class 'numpy.ndarray'> (16, 4) float64 rhs = <class 'numpy.ndarray'> (16, 1) float32 rcond = 3.5527136788005009e-15 lstsq = <function lstsq at 0x000002590591EBF0> 670 c = (c.T/scale).T # broadcast scale coefficients 672 # warn on rank reduction, which indicates an ill conditioned matrix

File ~.conda\envs\sediment\lib\site-packages\numpy\linalg\linalg.py:2326, in lstsq(a=<class 'numpy.ndarray'> (16, 4) float64, b=<class 'numpy.ndarray'> (16, 1) float32, rcond=3.5527136788005009e-15) 2323 if n_rhs == 0: 2324 # lapack can't handle n_rhs = 0 - so allocate the array one larger in that axis 2325 b = zeros(b.shape[:-2] + (m, n_rhs + 1), dtype=b.dtype) -> 2326 x, resids, rank, s = gufunc(a, b, rcond, signature=signature, extobj=extobj) signature = 'ddd->ddid' extobj = [8192, 1536, <function _raise_linalgerror_lstsq at 0x000002590591B790>] a = <class 'numpy.ndarray'> (16, 4) float64 b = <class 'numpy.ndarray'> (16, 1) float32 rcond = 3.5527136788005009e-15 gufunc = <ufunc 'lstsq_n'> 2327 if m == 0: 2328 x[...] = 0

File ~.conda\envs\sediment\lib\site-packages\numpy\linalg\linalg.py:124, in _raise_linalgerror_lstsq(err='invalid value', flag=8) 123 def _raise_linalgerror_lstsq(err, flag): --> 124 raise LinAlgError("SVD did not converge in Linear Least Squares")

LinAlgError: SVD did not converge in Linear Least Squares

[guiwitz]

Does this happen all the time or on a specific dataset? It's an algorithmic issue with computing the smoothing and it might be a very specific problem for a specific data point.

[Epta13]

So far it happend on all the cores which I tried, n=3

[guiwitz]

Can you print the output of conda list in your terminal?

[Epta13]

Bug seems to be only on Windows, pixed by: conda install "libblas==openblas" in the activated environment