Fitting code contribution: SR_LUMC Matlab code for IVIM, DTI, and IVIM-DTI fitting

[s-rauh]

Fitting code contribution in src/original/fitting. The code contains Matlab fitting routines for IVIM, IVIM-DTI, and DTI.

[Devguru-codes]

I went through the code and spotted a few things:

1. remove_zeros.m - zero mask applied too late
   The zeros get removed from data on line 13, but then the options.bval filtering on line 16 checks data==0 after the removal - so there are no zeros left and it's a no-op. The mask should be saved before modifying data
   code -
   % something like:
   idx = data==0;
   data(idx) = [];
   if ~isempty(options.bval)
   options.bval(idx) = [];
   varargout{1} = options.bval;
   end
   Doesn't affect the current fitting routines since none of them pass options.bval, but it'll bite anyone who tries to use that feature.

2. fit_ivim.m - segmented - v == 1 guard on data selection
   code -
   case 'segmented'
   if options.seg_data && v == 1
   tmpdat = [tmpdat(tmpb==0); tmpdat(tmpb>=options.bcut)];
   tmpb = [tmpb(tmpb==0); tmpb(tmpb>=options.bcut)];
   end

tmpdat and tmpb are re-initialized from scratch at the start of every loop iteration (tmpdat = datafit(:,v); tmpb = bval;), so the v == 1 condition means only the first voxel gets the b-value range selection. All other voxels fit the full range. Compare with how v == 1 is correctly used elsewhere in the same file - for one-time initialization of bounds that persist:
code -
case 'two_step'
if v == 1 % Correct: bounds persist across iterations
lb_twostep = lb([1 3 4]); ...
end

Looks like the v == 1 condition was carried over from the two_step case where it's used for one-time bounds init, but here the data selection needs to happen every iteration.

3. norm_diffdata.m - wrong normalization after b-value scaling
   if min(bval) >= 1
   mS0 = mean(data((bval == min(bval)),:), 1); % normalize to min bval
   else
   mS0 = mean(data(bval<1, :), 1); % normalize to "b=0"
   end

This is called after bval_scaling() divides b-values by 1000 (when max(bval) > 100). So typical clinical b-values [0, 50, 200, 800] become [0, 0.05, 0.2, 0.8]. Since min(bval) = 0 < 1, it uses bval < 1, which selects ALL b-values (they're all < 1 after scaling). The normalization factor becomes the mean of ALL signals across ALL b-values, not just b=0.

This effectively destroys the normalization - it should be bval == min(bval) or bval < 0.001 in the else branch.
Also, I noticed the fitting routines call @ivimfun, @dtifun, etc. but those model function files don't seem to be included - are they in a models/ directory that didn't get added?