Вращающиеся элементы в многомерном массиве MATLAB?

Question

У меня есть форма цилиндрической линии, которую я использую для создания трехмерного ядра; Я делаю это, сначала создавая простое двумерное ядро:

A = 625; 
B = 25; 

%2D kernal 
grid2d = zeros(101,101); 
c = 51; %centre 

for m = 1:101 
    for n = 1:101 

    r(m,n) = sqrt((m - c).^2 + (n - c).^2); %distance of point on grid to centre 

    %Populating the grid as a kernal 

    if r(m,n) <= 5 

     grid2d(m,n) = 100; 

    elseif r(m,n) >= 25 

     grid2d(m,n) = 0; 

    else 

     grid2d(m,n) = A./r(m,n) - B; 

    end 
    end 
end 

Это дает мне двумерное ядро. Теперь, если я определяю 3D-версию как 9 элементов, расположенных вдоль оси z внутри большей сетки, я могу создать 3D-ядро следующим:

gz = 147:155; %9 elements in the z axis 
H = length(gz); 
kernel3D = zeros(301,301,301); % 

for n = 1:H 

    kernel3D(151-50:1:151+50,151-50:1:151+50,gz(n)) = grid2d; 
end 

Это отлично работает, если я с помощью вертикальных линий источников, но мне интересно, можно ли вращать элементы этого массива в любом желаемом направлении, так что я мог бы генерировать ядра для наклонных источников; например, предположим, что я хотел повернуть этот массив на 45 градусов относительно плоскости XY и на 60 градусов по плоскости XZ через центр источника линии (151,151,151)?

Есть ли элегантный способ сделать это, возможно, используя матрицы вращения?

Answer 1

Я думаю, что нашел рабочее, но неэлегантное решение - в этой версии я сгенерировал 2D-ядро и вложил его в трехмерный массив. Затем я использую ротационные операторы для сопоставления соответствующих значений в позиции с преобразованием вращения. Это уродливо, но, похоже, работает.

%generate a test 2d kernel as before... 

A = 625; B = 25; c = 51; cen = 151; grid2d = zeros(101,101); 

for m = 1:101 
for n = 1:101 


     r(m,n) = sqrt((m - c).^2 + (n - c).^2); %distance of point on grid to centre 

%Populating the grid as a kernal 

if r(m,n) <= 5 

    grid2d(m,n) = 100; 

elseif r(m,n) >= 25 

    grid2d(m,n) = 0; 

else 

    grid2d(m,n) = A./r(m,n) - B; 

    end 
    end 
    end 

%make a bigger array for this to nest in 
flat3d = zeros(301,301,301); 
range = 101:201; %centre range of new grid 

flat3d(range,range,cen) = grid2d; % inserts grid2d at centre of new flat3d grid! 

%now the incline angles 
theta = input('Enter incline in XZ (vertical) plane in degrees : '); 
phi = input('Enter incline in XY (horizontal) plane in degrees : '); 

%Now a loop which inclines each point relative to theta around (x,y,z) = (cen,cen,cen); 

g = range(1); 

inc3d = zeros(301,301,301); %new inclined grid initialised 

for m = 1:101 
for n = 1:101 


    %the g - cen value is added to the index to map arrays sync up correctly and 
    %keep loop reasonable/avoid overlap outside index ! 

    xpart(m,n) = cen + cosd(theta).*(g + m -cen) - sind(theta).*(g + n - cen); 
    xn(m,n) = cen + cosd(phi).*(xpart(m,n) - cen) - sind(phi).*(g + n - cen); %this maps a value of x to new position theta and phi 
    yn(m,n) = cen + sind(phi).*(xn(m,n) - cen) + cosd(phi).*(g + n - cen); %maps a value of y to new position relative to phi 
    zn(m,n) = cen + sind(theta).*(g + m - cen) + cosd(theta).*(g + n - cen); %maps z to new position relative to theta (doesnt change relative to phi) 

    xn = round(xn); 
    yn = round(yn); %round to avoid integer problems with indexing 
    zn = round(zn); 


    inc3d(xn(m,n),yn(m,n),zn(m,n)) = flat3d(g + m,g + n,cen); 

end 
end 

Это, как представляется, работает; если я создаю изображение plot3 тэта = 30 и phi = 15, например, я получаю что-то вроде этого;

С этой 2D плоского изображения, с уравнением линии, я думаю, можно было бы просто суммировать и порождают ядро, для которых требуется когда-либо регионы. Это грязно, но, похоже, работает. Если кто-нибудь знает более чистый способ, обязательно добавьте его :)