Rotate Spherical Cap On Sphere

Dear All,

I am trying to produce a plot where I have a sphere with another spherical section representation a region of interest. I also want to be able to rotate this spherical section over any theta,phi. I will walk though what I have done and explain where I am stuck.

The first think I do is produce a sphere that is transparent. Here is the code for that.

%%%%%%%%%%% Generate a sphere consisting of 200 by 200 faces  %%%%%%%%%%%%%
[x,y,z]=sphere(200);
rad=1;
hSurface1=surf(rad*x,rad*y,rad*z);hold on
set(hSurface1,'FaceColor',[0 0 1],'FaceAlpha',0.25,'FaceLighting','gouraud','EdgeColor','none')
hSurface2=surf(rad*x*.7,rad*y*.7,rad*z*.7);hold on
set(hSurface2,'FaceColor',[0 1 0],'FaceAlpha',0.25,'FaceLighting','gouraud','EdgeColor','none')
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

This is one sphere with unit radius =1. The inner sphere represents the shell thickness. It is for visual purposes only. I next define a range for the spherical section. As expected the phi range is from 0-2pi. The theta range is determined by a region of interest for the problem. For example:

numPoints = 4;
thetaMin = 0
thetaMax = pi/2-0.68;
thetaRange = linspace(thetaMin, thetaMax, numPoints);

phiMin = 0;
phiMax = 2*pi;
phiRange = linspace(phiMin, phiMax, numPoints);

[theta,phi] = meshgrid(thetaRange,phiRange,numPoints);
[x1,y1,z1] = sph2cart(phi, theta, rad);
x1 = round(x1,4);
y1 = round(y1,4);
z1 = round(z1,4);
%surf(x1,y1,z1)
plot3(x1,y1,z1,'*w');

This will generate a spherical cap (or points) that sits on top of the sphere. My question is how do I rotate the cap (or indvidual points) using any theta and phi. Is there a unique way to transform the mesh points? I have tried to convert the individual points in the following approach:

theta_rot = 0;phi_rot = 0;

These are the new theta and phi I want to rotate mesh points to. I tried to find the correct difference in theta and phi to calculate a new vector b.

for i=1:length(x1)
for j=1:length(x1)
a=[x1(i,j) y1(i,j) z1(i,j)] 
[azimuth,elevation,r] = cart2sph(x1(i,j),y1(i,j),z1(i,j))
b=[sind(theta_rot+elevation)*cosd(phi_rot+azimuth) sind(theta_rot+elevation)*sind(phi_rot+azimuth) cosd(theta_rot+elevation)]

w=cross(a,b);
d=dot(a,b);
Id = [1 0 0;0 1 0;0 0 1];
vx = [0 -w(3) w(2);w(3) 0 w(1);-w(2) w(1) 0];
R = Id + vx + vx^2*((1-d)/norm(w)^2);

rot_vals = mtimes(R,[x1(i,j) y1(i,j) z1(i,j)]');
rot_x1(i,j) = rot_vals(1);
rot_y1(i,j) = rot_vals(2);
rot_z1(i,j) = rot_vals(3);
end
end

This approach does not work since I am unable to find a better way to construct vector (b) by finding the difference in the points. I hope this is not too confusing and any help is appreaciated.

Answer 1

[darova]

I rotated data about X axis first (theta angle). Then i am rotating about Z (phi angle)

clc,clear

p = linspace(0,2*pi,50);
t = linspace(pi/2,pi/6,10);
[phi,theta] = ndgrid(p,t);
[x1,y1,z1] = sph2cart(phi, theta, 1);


[X,Y,Z] = sphere(20);
h = surf(X,Y,Z);
set(h,'EdgeColor','none')
% set(h,'FaceCOlor', 'b');
alpha(h,0.5)

t = 45;
st1 = sind( t*sind(t) );
ct1 = cosd( t*sind(t) );
hold on
for t = linspace(0,180,100)
        % ----------    PLAY WITH THIS
%     st1 = sind( t*sind(t) );
%     ct1 = cosd( t*sind(t) );
    st3 = sind( t );
    ct3 = cosd( t );
        % ------------
    Rx = [1 0 0; 0 ct1 -st1; 0 st1 ct1];    % rotation matrix around X axis
    Rz = [ct3 -st3 0; st3 ct3 0; 0 0 1];    % rotation matrix around Z axis
    V = Rx*[x1(:) y1(:) z1(:)]';            % rotate data
    V = Rz*V;
    x2 = reshape(V(1,:),size(x1));
    y2 = reshape(V(2,:),size(x1));
    z2 = reshape(V(3,:),size(x1));
    h = plot3(x2,y2,z2,'g');
    pause(0.05)
    delete(h)
end

hold off
axis equal

[Chad]

I am unable to get this to work with theta below the equator or at theta of 135 degrees. Am I missing something?