How to vectorize this nested for loop?

Question

Preetham Manjunatha 2024-5-4

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/2115151-how-to-vectorize-this-nested-for-loop

编辑： Preetham Manjunatha 2024-5-5

checker.jpg

I am finding the cylindrical projection of an image given the camera intrinsic matrix and the distortion coefficients. After finding the warpped coordinates, I need to find the image value for those cooridnates and assign them to an empty warpped image. I would like to vectorize the nested for loop written at the end. I tried with

sub2ind

had series of errors and gave up. Below is the code:

clc; close all; clear;
% Inputs
fileName = 'checker.jpg';
% Add grid line
add_grid = 0;
% Focal lengths
fx = 50;
fy = 50;
% Read image
image = (imread(fileName));
% Get image size
[ydim, xdim, bypixs] = size(image);
% Camera intrinsics
K = [fx, 0, xdim/2; 0, fy, ydim/2; 0, 0, 1];
% Distortion coefficients [k1, k2, k3, p1, p2]
DC = [0, 0, 0, 0, 0];
% Get distrotion coefficients
fx = K(1,1);
fy = K(2,2);
k1 = DC(1);
k2 = DC(2);
k3 = DC(3);
p1 = DC(4);
p2 = DC(5);
% Get image size
[ydim, xdim, bypixs] = size(image);
% Initialize an array
imageCylindrical = uint8(zeros(ydim, xdim, bypixs));
% Get the center of image
xc = xdim/2;
yc = ydim/2;
% Create X and Y coordinates grid
[X,Y] = meshgrid(1:xdim, 1:ydim);
% Perform the cylindrical projection
theta = (X - xc) / fx;
h = (Y - yc) / fy;
% Cylindrical coordinates to Cartesian
xcap = sin(theta);
ycap = h;
zcap = cos(theta);
xyz_cap = cat(3, xcap, ycap, zcap);
xyz_cap = reshape(xyz_cap,[],3);
% Normalized coords
xyz_cap_norm = (K * xyz_cap')';
xn = xyz_cap_norm(:,1) ./ xyz_cap_norm(:,3);
yn = xyz_cap_norm(:,2) ./ xyz_cap_norm(:,3);
% Radial and tangential distortion 
r = xn.^2 + yn.^2;
xd_r = xn .* (1 + k1 * r.^2 + k2 * r.^4 + k3 * r.^6);
yd_r = yn .* (1 + k1 * r.^2 + k2 * r.^4 + k3 * r.^6);
xd_t = 2 * p1 * xn .* yn + p2 * (r.^2 + 2 * xn.^2);
yd_t = p1 * (r.^2 + 2 * yn.^2) + 2 * p2 * xn .* yn;
xd = xd_r + xd_t;
yd = yd_r + yd_t;
% Reshape and clip coordinates
xd = reshape(ceil(xd),[ydim, xdim]);
yd = reshape(ceil(yd),[ydim, xdim]);
mask = xd > 0 & xd <= xdim & yd > 0 & yd <= ydim;
% Get masked coordinates
xd = xd .* mask;
yd = yd .* mask;
% Get projections
for i = 1:ydim
    for j = 1:xdim
        if yd(i,j) ~= 0 || xd(i,j)~=0
            imageCylindrical(i,j,:) = image(yd(i,j), xd(i,j),:);
        end
    end
end

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

Walter Roberson 2024-5-4

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/2115151-how-to-vectorize-this-nested-for-loop#answer_1452401

编辑：Walter Roberson 2024-5-4

在 MATLAB Online 中打开

ind = sub2ind(size(image), yd, xd, 1);
imageCylindrical(:,:,1) = image(ind + 0 * ydim*xdim);
imageCylindrical(:,:,2) = image(ind + 1 * ydim*xdim);
imageCylindrical(:,:,3) = image(ind + 2 * ydim*xdim);

4 个评论
显示 2更早的评论隐藏 2更早的评论

Walter Roberson 2024-5-4

在 MATLAB Online 中打开

checker.jpg

% Inputs

fileName = 'checker.jpg';

% Add grid line

add_grid = 0;

% Focal lengths

fx = 50;

fy = 50;

% Read image

image = (imread(fileName));

% Get image size

[ydim, xdim, bypixs] = size(image);

% Camera intrinsics

K = [fx, 0, xdim/2; 0, fy, ydim/2; 0, 0, 1];

% Distortion coefficients [k1, k2, k3, p1, p2]

DC = [0, 0, 0, 0, 0];

% Get distrotion coefficients

fx = K(1,1);

fy = K(2,2);

k1 = DC(1);

k2 = DC(2);

k3 = DC(3);

p1 = DC(4);

p2 = DC(5);

% Get image size

[ydim, xdim, bypixs] = size(image);

% Initialize an array

%imageCylindrical = uint8(zeros(ydim, xdim, bypixs));

% Get the center of image

xc = xdim/2;

yc = ydim/2;

% Create X and Y coordinates grid

[X,Y] = meshgrid(1:xdim, 1:ydim);

% Perform the cylindrical projection

theta = (X - xc) / fx;

h = (Y - yc) / fy;

% Cylindrical coordinates to Cartesian

xcap = sin(theta);

ycap = h;

zcap = cos(theta);

xyz_cap = cat(3, xcap, ycap, zcap);

xyz_cap = reshape(xyz_cap,[],3);

% Normalized coords

xyz_cap_norm = (K * xyz_cap')';

xn = xyz_cap_norm(:,1) ./ xyz_cap_norm(:,3);

yn = xyz_cap_norm(:,2) ./ xyz_cap_norm(:,3);

% Radial and tangential distortion

r = xn.^2 + yn.^2;

xd_r = xn .* (1 + k1 * r.^2 + k2 * r.^4 + k3 * r.^6);

yd_r = yn .* (1 + k1 * r.^2 + k2 * r.^4 + k3 * r.^6);

xd_t = 2 * p1 * xn .* yn + p2 * (r.^2 + 2 * xn.^2);

yd_t = p1 * (r.^2 + 2 * yn.^2) + 2 * p2 * xn .* yn;

xd = xd_r + xd_t;

yd = yd_r + yd_t;

% Reshape and clip coordinates

xd = reshape(ceil(xd),[ydim, xdim]);

yd = reshape(ceil(yd),[ydim, xdim]);

mask = xd > 0 & xd <= xdim & yd > 0 & yd <= ydim;

% Get masked coordinates

xd = xd .* mask;

yd = yd .* mask;

ind = sub2ind(size(image), yd(mask), xd(mask), 1);

IC1 = zeros(ydim, xdim, 'uint8');

IC1(mask) = image(ind + 0 * ydim * xdim);

if bypixs == 1

imageCylindrical = IC1;

else

IC2 = zeros(ydim, xdim, 'uint8');

IC3 = zeros(ydim, xdim, 'uint8');

IC2(mask) = image(ind + 1 * ydim * xdim);

IC3(mask) = image(ind + 2 * ydim * xdim);

imageCylindrical = cat(3, IC1, IC2, IC3);

end

imshow(imageCylindrical)

Preetham Manjunatha 2024-5-5

在 MATLAB Online 中打开

It is strange that when I run on my Windows 11 machine (MATLAB R2023a), it throws an error:

Error using sub2ind
The subscript vectors must all be of the same size.
Error in demo (line 123)
ind = sub2ind(size(image), yd(mask), xd(mask), 1);

However, when I was able to run on online MATLAB, it works and produces the image.

Preetham Manjunatha 2024-5-5

编辑：Preetham Manjunatha 2024-5-5

在 MATLAB Online 中打开

Adding and fixing the ind:

ind = sub2ind(size(image), yd(mask), xd(mask), ones(1,length(xd(mask)))');

did the trick. It seems R2024a can work with

ind = sub2ind(size(image), yd(mask), xd(mask), 1);

Thanks Walter!

请先登录，再进行评论。

How to vectorize this nested for loop?

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

采纳的回答

4 个评论
显示 2更早的评论隐藏 2更早的评论

更多回答（0 个）

另请参阅

类别

标签

产品

版本

Community Treasure Hunt

How to vectorize this nested for loop?

0 个评论 显示 -2更早的评论隐藏 -2更早的评论

采纳的回答

4 个评论 显示 2更早的评论隐藏 2更早的评论

更多回答（0 个）

另请参阅

类别

标签

产品

版本

Community Treasure Hunt

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

4 个评论
显示 2更早的评论隐藏 2更早的评论