Double integral over array

Hi,

Evaluating a double integral over a large number of pixels takes an hour, even with the parallel toolbox. I can't figure out how to vectoirze this calculation to speed it up...suggestions welcome!

parfor iRow = 1:nRows

for iCol = 1:nCols

X = pitch_mm * (iRow - nRows/2 - 0.5);

Y = pitch_mm * (iCol - nCols/2 - 0.5);

% projected solid angle differential

dOmega = @(Xap,Yap) (apDist_mm ./ ((apDist_mm)^2 + (X-Xap).^2 + (Y-Yap).^2)).^2;

% convert to polar & integrate

dOmegaPolar = @(theta,r) dOmega(r.*cos(theta),r.*sin(theta)).*r;

cos4(iRow,iCol) = integral2(dOmegaPolar,0,2*pi,0,apDiam_mm/2);

end

end

Answer 1

[David Hill]

X=pitch_mm*([1:nRows]-nRows/2-.5);
Y=pitch_mm*([1:nCols]-nCols/2-.5);
dOmegaPolar = @(theta,r) dOmega(r.*cos(theta),r.*sin(theta)).*r;
dOmega = @(Xap,Yap) (apDist_mm ./ ((apDist_mm)^2 + (x-Xap).^2 + (y-Yap).^2)).^2;
cos4=zeros(nRows,nCols);%preallocate
parfor iRow = 1:nRows
    x=X(iRow);
    for iCol = 1:nCols
        y=Y(iCol);
        cos4(iRow,iCol) = integral2(dOmegaPolar,0,2*pi,0,apDiam_mm/2);
    end
end

Above should speed things up some.

[David Hill]

X=pitch_mm*([1:nRows]-nRows/2-.5);
Y=pitch_mm*([1:nCols]-nCols/2-.5);
dOmegaPolar = @(theta,r) dOmega(r.*cos(theta),r.*sin(theta)).*r;
cos4=zeros(nRows,nCols);%preallocate
parfor iRow = 1:nRows
    for iCol = 1:nCols
        dOmega = @(Xap,Yap) (apDist_mm ./ ((apDist_mm)^2 + (X(iRow)-Xap).^2 + (Y(iCol)-Yap).^2)).^2;%this needs to move inside loop
        cos4(iRow,iCol) = integral2(dOmegaPolar,0,2*pi,0,apDiam_mm/2);
    end
end

[steve solomon]

thanks David. That looks like it should be faster but it's not, about an hour of run time.