create_cycdeg_image.m 2.7 KB

12345678910111213141516171819202122232425262728293031323334353637383940414243444546474849505152535455
  1. function D = create_cycdeg_image( im_size, pix_per_deg )
  2. % CREATE_CYCDEG_IMAGE (internal) create matrix that contains frequencies,
  3. % given in cycles per degree.
  4. %
  5. % D = create_cycdeg_image( im_size, pix_per_deg )
  6. % im_size - [height width] vector with image size
  7. % pix_per_deg - pixels per degree for both horizontal and vertical axis
  8. % (assumes square pixels)
  9. %
  10. % Useful for constructing Fourier-domain filters based on OTF or CSF data.
  11. %
  12. % (C) Rafal Mantiuk <mantiuk@gmail.com>
  13. % This is an experimental code for internal use. Do not redistribute.
  14. nyquist_freq = 0.5 * pix_per_deg;
  15. half_size = floor(im_size/2);
  16. odd = mod( im_size, 2 );
  17. freq_step = nyquist_freq./half_size;
  18. if( odd(2) )
  19. xx = [ linspace( 0, nyquist_freq, half_size(2)+1 ) linspace( -nyquist_freq, -freq_step(2), half_size(2) ) ];
  20. else
  21. xx = [ linspace( 0, nyquist_freq-freq_step(2), half_size(2) ) linspace( -nyquist_freq, -freq_step(2), half_size(2) ) ];
  22. end
  23. if( odd(1) )
  24. yy = [ linspace( 0, nyquist_freq, half_size(1)+1 ) linspace( -nyquist_freq, -freq_step(1), half_size(1) ) ];
  25. else
  26. yy = [ linspace( 0, nyquist_freq-freq_step(1), half_size(1) ) linspace( -nyquist_freq, -freq_step(1), half_size(1) ) ];
  27. end
  28. [XX YY] = meshgrid( xx, yy );
  29. D = sqrt( XX.^2 + YY.^2 );
  30. %[XX YY] = meshgrid( linspace( 0, nyquist_freq, half_size(2)+even(2) ), linspace( 0, nyquist_freq, half_size(1)+even(1) ) );
  31. %D1 = sqrt( XX.^2 + YY.^2 );
  32. %[XX YY] = meshgrid( linspace( nyquist_freq-frec_step(2), frec_step(2), half_size(2) ), linspace( 0, nyquist_freq, half_size(1)+even(1) ) );
  33. %D2 = sqrt( XX.^2 + YY.^2 );
  34. %[XX YY] = meshgrid( linspace( 0, nyquist_freq, half_size(2)+even(2) ), linspace( nyquist_freq-frec_step(1), frec_step(1), half_size(1) ) );
  35. %D3 = sqrt( XX.^2 + YY.^2 );
  36. %[XX YY] = meshgrid( linspace( nyquist_freq-frec_step(2), frec_step(2), half_size(2) ), linspace( nyquist_freq-frec_step(1), frec_step(1), half_size(1) ) );
  37. %D4 = sqrt( XX.^2 + YY.^2 );
  38. %D = [ D1 D2; D3 D4 ];
  39. %[XX YY] = meshgrid( linspace( 0, nyquist_freq, half_size(2)+even(2) ), linspace( 0, nyquist_freq, half_size(1)+even(1) ) );
  40. %D1 = sqrt( XX.^2 + YY.^2 );
  41. %[XX YY] = meshgrid( linspace( nyquist_freq-frec_step(2), frec_step(2), half_size(2) ), linspace( 0, nyquist_freq, half_size(1)+even(1) ) );
  42. %D2 = sqrt( XX.^2 + YY.^2 );
  43. %[XX YY] = meshgrid( linspace( 0, nyquist_freq, half_size(2)+even(2) ), linspace( nyquist_freq-frec_step(1), frec_step(1), half_size(1) ) );
  44. %D3 = sqrt( XX.^2 + YY.^2 );
  45. %[XX YY] = meshgrid( linspace( nyquist_freq-frec_step(2), frec_step(2), half_size(2) ), linspace( nyquist_freq-frec_step(1), frec_step(1), half_size(1) ) );
  46. %D4 = sqrt( XX.^2 + YY.^2 );
  47. end