Getting image data values with coordinates in circular region

Question

I have a number of gray level images (600, with size: 2400x2400 pixel, small instance see above). In this image there is a circular region of interest. The aim is to fit a function to the gray value data within region of interest (fit to a Zernike polynomial). For the fit with FindFit I need a list with points {{x1, y1,"image intensity value"},...} (x, y-being pixel positions) only from pixels within the region of interest. I only have a brute force approach: 1) Get the gray level data by ImageData[image], 2) Go through each pixel by Table (ParallelTable makes no difference), checking for each pixel if it is in the region of interest, if yes I store the value, if not I set the value to zero, 3) I then Select those data points that are not zero. Step two takes a loooong time (5sec with (345x345) image on windows 7 64bit, i7 laptop, dual core). After having the values of the image I use it in FindFit to fit Zernike polynomial. I so far only had the idea using procedural programming. What I am looking for is a significant increase in speed. My first idea was using functional programing or using matrix operation ... but I have not found something useful yet.

Does someone have an idea to increase the speed?

Note: I am aware of alternative techniques of retrieving the Zernike coefficients by integration technique. But this will be second step. Last step would be using CUDA for it.

This is the brute force code I have:

subDat=ParallelTable[{N[Sqrt[((i-centKoor[[1]])/roiRad)^2+((j-centKoor[[2]])/roiRad)^2]],N[ArcTanB[(j-centKoor[[2]]),(i-centKoor[[1]])]],If[Sqrt[(i-centKoor[[1]])^2+(j-centKoor[[2]])^2]<=roiRad,N[bildDats[[j,i]]],0]},{j,centKoor[[2]]-roiRad,centKoor[[2]]+roiRad},{i,centKoor[[1]]-roiRad,centKoor[[1]]+roiRad}];
subFitDat = Select[Flatten[subDat, 1], #[[3]] > 0 &];
ergSubuNorm=FindFit[subFitDat,term,par,{r,\[Phi]}];

centKoor->center coordinate of region of interest; roiRad->radius in pixel of radius of interest; bildDats->image data r, phi->unit radius and azimutal angle (polar coordinates)

Answer 1

Try to avoid Table where you can and use vectorized expressions instead. For example, if you create 2d arrays for x and y-coordinates:

dims = {500, 500};
center = {100, 100};
ys = Array[N[#1] &, dims] - center[[1]];
xs = Array[N[#2] &, dims] - center[[2]];

you can then use xs and ys in expressions, and Sqrt, ArcTan and friends will automatically thread over all values. This is usually very fast:

radius = Sqrt[xs^2 + ys^2];
angle = ArcTan[xs, ys];

If you need a boolean comparison, use UnitStep instead of Table/If where possible:

mask = UnitStep[100 - radius];

You can later transpose these arrays to the format FindFit needs, e.g.:

data = RandomReal[1, dims];
Pick[Transpose[Flatten /@ {radius, angle, data}], Flatten[mask], 1]

Answer 2

If your original image is image1 and your mask is maskimage then you can do the following:

newImage = ImageMultiply[image1, maskimage];

newimage contains only the interesting gray level pixels where your maskimage is white.