spm_dot.m 1.9 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566
  1. function [X] = spm_dot(X,x,i)
  2. % Multidimensional dot (inner) product
  3. % FORMAT [Y] = spm_dot(X,x,[DIM])
  4. %
  5. % X - numeric array
  6. % x - cell array of numeric vectors
  7. % DIM - dimensions to omit (asumes ndims(X) = numel(x))
  8. %
  9. % Y - inner product obtained by summing the products of X and x along DIM
  10. %
  11. % If DIM is not specified the leading dimensions of X are omitted.
  12. % If x is a vector the inner product is over the leading dimension of X
  13. %
  14. % See also: spm_cross
  15. %__________________________________________________________________________
  16. % Copyright (C) 2015 Wellcome Trust Centre for Neuroimaging
  17. % Karl Friston
  18. % $Id: spm_dot.m 7314 2018-05-19 10:13:25Z karl $
  19. % initialise dimensions
  20. %--------------------------------------------------------------------------
  21. if iscell(x)
  22. DIM = (1:numel(x)) + ndims(X) - numel(x);
  23. else
  24. DIM = 1;
  25. x = {x};
  26. end
  27. % omit dimensions specified
  28. %--------------------------------------------------------------------------
  29. if nargin > 2
  30. DIM(i) = [];
  31. x(i) = [];
  32. end
  33. % inner product using recursive summation (and bsxfun)
  34. %--------------------------------------------------------------------------
  35. for d = 1:numel(x)
  36. s = ones(1,ndims(X));
  37. s(DIM(d)) = numel(x{d});
  38. X = bsxfun(@times,X,reshape(full(x{d}),s));
  39. X = sum(X,DIM(d));
  40. end
  41. % eliminate singleton dimensions
  42. %--------------------------------------------------------------------------
  43. X = squeeze(X);
  44. return
  45. % NB: alternative scheme using outer product
  46. %==========================================================================
  47. % outer product and sum
  48. %--------------------------------------------------------------------------
  49. x = spm_cross(x);
  50. s = ones(1,ndims(X));
  51. S = size(X);
  52. s(DIM) = S(DIM);
  53. x = reshape(full(x),s);
  54. X = bsxfun(@times,X,x);
  55. for d = 1:numel(DIM)
  56. X = sum(X,DIM(d));
  57. end
  58. X = squeeze(X);