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 onset(hSurface1,'FaceColor',[0 0 1],'FaceAlpha',0.25,'FaceLighting','gouraud','EdgeColor','none')hSurface2=surf(rad*x*.7,rad*y*.7,rad*z*.7);hold onset(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 = 0thetaMax = 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);endend
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.
Best Answer