0

I am writing a function to draw an approximate circle on a square array (in Matlab, but the problem is mainly algorithmic).
The goal is to produce a mask for integrating light that falls on a portion of a CCD sensor from a diffraction-limited point source (whose diameter corresponds to a few pixels on the CCD array). In summary, the CCD sensor sees a pattern with revolution-symmetry, that has of course no obligation to be centered on one particular pixel of the CCD (see example image below).

Here is the algorithm that I currently use to produce my discretized circular mask, and which works partially (Matlab/Octave code):

xt = linspace(-xmax, xmax, npixels_cam); % in physical coordinates (meters)
[X Y] = meshgrid(xt-center(1), xt-center(2)); % shifted coordinate matrices
[Theta R] = cart2pol(X,Y);
R = R'; % cart2pol uses a different convention for lines/columns
mask = (R<=radius);

As you can see, my algorithm selects (sets to 1) all the pixels whose physical distance (in meters) is smaller or equal to a radius, which doesn't need to be an integer.

I feel like my algorithm may not be the best solution to this problem. In particular, I would like it to include the pixel in which the center is present, even when the radius is very small.
Any ideas ?

(See https://i.stack.imgur.com/3mZ5X.png for an example image of a diffraction-limited spot on a CCD camera).

calvin tiger
  • 391
  • 2
  • 7
  • 18
  • uh, do you really mean a pseudo-circle http://en.wikipedia.org/wiki/Pseudo-circle ? In that case my answer may be void. – Johan Lundberg Jan 26 '12 at 11:43
  • No no no, I just meant a "circle-like" thing :) You did understand my question correctly.. – calvin tiger Jan 26 '12 at 13:08
  • Your answer is very interesting and helpful, but I am also waiting for input from people who could suggest "standard" techniques used in such situations (I believe that my problem is quite common in digital image processing, and I am therefore surprised that no one else has commented/voted on this page). In any case, I will try your solution, and I am sure that it will solve my problem even if no straightforward method exists. Thank you again for your help ! – calvin tiger Jan 27 '12 at 14:04

2 Answers2

1

if you like to select pixels if and only if they contain any part of the circle C:

in each pixel place a small circle A with the radius = halv size of the pixel, and another one around it with R=sqrt(2)*half size of the circle (a circumscribed circle)

To test if two circles touch each other you just calculate the center to center distances and subtract the sum of the two radii.

If the test circle C is within A then you select the pixel. If it's within B but not C you need to test all four pixel sides for overlap like this Circle line-segment collision detection algorithm?

A brute force approximate method is to make a much finer grid within each pixel and test each center point in that grid.

Community
  • 1
  • 1
Johan Lundberg
  • 26,184
  • 12
  • 71
  • 97
  • Also, unless you are really studying sub-pixel effects at the edge of the mask, just doing what you have done should be fine. If you just want to make sure to cover all pixels which touch the circle and in addition you are not worried by having some 1-pixel overlaps, just do mask = (R<=radius+fudge) where fudge is 0.5*sqrt(pixelsize_x^2+pixelsize_y^2) – Johan Lundberg Jan 27 '12 at 18:17
  • Also calculating and comparing R2=sqrt(X.^2+Y.^2) is almost 5 times faster than using car2pol for a 5000x5000 matrix. – Johan Lundberg Jan 27 '12 at 18:34
1

This is a well-studied problem. Several levels of optimization are possible:

  • You can brute-force check if every pixel is inside the circle. (r^2 >= (x-x0)^2 + (y-y0)^2)

  • You can brute-force check if every pixel in a square bounding the circle is inside the circle. (r^2 >= (x-x0)^2 + (y-y0)^2 where |x-x0| < r and |y-y0| < r)

  • You can go line-by-line (where |y-y0| < r) and calculate the starting x ending x and fill all the lines in between. (Although square roots aren't cheap.)

  • There's an infinite possibility of more sophisticated algorithms. Here's a common one: http://en.wikipedia.org/wiki/Midpoint_circle_algorithm (filling in the circle is left as an exercise)

It really depends on how sophisticated you want to be based on how imperative good performance is.

Kaganar
  • 6,540
  • 2
  • 26
  • 59
  • A simple to understand method that performs reasonably well is to check every pixel in a square bounding area but with an optimized check -- instead of checking r^2 >= (x-x0)^2 + (y-y0)^2, check this once per each line on the left side of the bounding square and then add the difference between (x-x0)^2 and (x+1-x0)^2 to a value as you march across the line until it's over r^2 -- (x+1-x0)^2 - (x-x0)^2 = 2x + 1 (if I did my math right), and that's a lot faster to compute per pixel... but that's only true of low level optimization. Not sure if MATLAB extends that far. – Kaganar Jan 27 '12 at 16:08
  • Additionally, you may want to look into MATLAB's built-in functions. For example, if there's a function for filling a horizontal line in your square array with one call, it likely trumps other methods of filling pixels -- so you're better off computing the left and right bounds of the circle on each line and filling it via MATLAB's call. With this method you can likely even use the aforementioned Midpoint algorithm to make it even more efficient. – Kaganar Jan 27 '12 at 16:09