フィルターのクリア

How to vectorize this nested for loop?

3 ビュー (過去 30 日間)
Preetham Manjunatha
Preetham Manjunatha 2024 年 5 月 4 日
編集済み: Preetham Manjunatha 2024 年 5 月 5 日
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

採用された回答

Walter Roberson
Walter Roberson 2024 年 5 月 4 日
編集済み: Walter Roberson 2024 年 5 月 4 日
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 件のコメント
Preetham Manjunatha
Preetham Manjunatha 2024 年 5 月 5 日
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
Preetham Manjunatha 2024 年 5 月 5 日
編集済み: Preetham Manjunatha 2024 年 5 月 5 日
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!

サインインしてコメントする。

その他の回答 (0 件)

カテゴリ

Help Center および File ExchangeMedical Physics についてさらに検索

製品


リリース

R2023a

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by