circ_samplecdf.m 1.7 KB

1234567891011121314151617181920212223242526272829303132333435363738394041424344454647484950515253545556575859606162636465666768697071727374
  1. function [phis, cdf, phiplot, cdfplot] = circ_samplecdf(thetas, resolution)
  2. % [phis, cdf, phiplot, cdfplot] = circ_samplecdf(thetas, resolution)
  3. %
  4. % Helper function for circ_kuipertest.
  5. % Evaluates CDF of sample in thetas.
  6. %
  7. % Input:
  8. % thetas sample (in radians)
  9. % resolution resolution at which the cdf is evaluated
  10. %
  11. % Output:
  12. % phis angles at which CDF is evaluated
  13. % cdf CDF values at these angles
  14. % phiplot as phi, for plotting
  15. % cdfplot as cdf, for plotting
  16. %
  17. %
  18. % Circular Statistics Toolbox for Matlab
  19. % By Marc J. Velasco, 2009
  20. % velasco@ccs.fau.edu
  21. if nargin < 2
  22. resolution = 100;
  23. end
  24. phis = 0;
  25. cdf = zeros(1, length(phis));
  26. phis = linspace(0,2*pi,resolution+1);
  27. phis = phis(1:end-1);
  28. % ensure all points in thetas are on interval [0, 2pi)
  29. x = thetas(thetas<0);
  30. thetas(thetas<0) = (2*pi-abs(x));
  31. % compute cdf
  32. thetas = sort(thetas);
  33. dprob = 1/length(thetas); %incremental change in probability
  34. cumprob = 0; %cumultive probability so far
  35. % for a little bit, we'll add on 2pi to the end of phis
  36. phis = [phis 2*pi];
  37. for j=1:resolution
  38. minang = phis(j);
  39. maxang = phis(j+1);
  40. currcount = sum(thetas >= minang & thetas < maxang);
  41. cdf(j) = cumprob + dprob*currcount;
  42. cumprob = cdf(j);
  43. end
  44. phis = phis(1:end-1);
  45. % for each point in x, duplicate it with the preceding value in y
  46. phis2 = phis;
  47. cdf2 = [0 cdf(1:end-1)];
  48. cdfplottable = [];
  49. phisplottable = [];
  50. for j=1:length(phis);
  51. phisplottable = [phisplottable phis(j) phis2(j)]; %#ok<AGROW>
  52. cdfplottable = [cdfplottable cdf2(j) cdf(j)]; %#ok<AGROW>
  53. end
  54. phiplot = [phisplottable 2*pi];
  55. cdfplot = [cdfplottable 1];