circ_moment.m 1.6 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869
  1. function [mp rho_p mu_p] = circ_moment(alpha, w, p, cent, dim)
  2. % [mp cbar sbar] = circ_moment(alpha, w, p, cent, dim)
  3. % Calculates the complex p-th centred or non-centred moment
  4. % of the angular data in angle.
  5. %
  6. % Input:
  7. % alpha sample of angles
  8. % [w weightings in case of binned angle data]
  9. % [p p-th moment to be computed, default is p=1]
  10. % [cent if true, central moments are computed, default = false]
  11. % [dim compute along this dimension, default is 1]
  12. %
  13. % If dim argument is specified, all other optional arguments can be
  14. % left empty: circ_moment(alpha, [], [], [], dim)
  15. %
  16. % Output:
  17. % mp complex p-th moment
  18. % rho_p magnitude of the p-th moment
  19. % mu_p angle of th p-th moment
  20. %
  21. %
  22. % References:
  23. % Statistical analysis of circular data, Fisher, p. 33/34
  24. %
  25. % Circular Statistics Toolbox for Matlab
  26. % By Philipp Berens, 2009
  27. % berens@tuebingen.mpg.de
  28. if nargin < 5
  29. dim = 1;
  30. end
  31. if nargin < 4
  32. cent = false;
  33. end
  34. if nargin < 3 || isempty(p)
  35. p = 1;
  36. end
  37. if nargin < 2 || isempty(w)
  38. % if no specific weighting has been specified
  39. % assume no binning has taken place
  40. w = ones(size(alpha));
  41. else
  42. if size(w,2) ~= size(alpha,2) || size(w,1) ~= size(alpha,1)
  43. error('Input dimensions do not match');
  44. end
  45. end
  46. if cent
  47. theta = circ_mean(alpha,w,dim);
  48. v = size(alpha)./size(theta);
  49. alpha = circ_dist(alpha,repmat(theta,v));
  50. end
  51. n = size(alpha,dim);
  52. cbar = sum(cos(p*alpha).*w,dim)/n;
  53. sbar = sum(sin(p*alpha).*w,dim)/n;
  54. mp = cbar + i*sbar;
  55. rho_p = abs(mp);
  56. mu_p = angle(mp);