So I have a profile that I am rotating to make a 3D surface like this
a=0.661/0.511;
r=2.8179;
vo=1;
m=0.511;
theta=linspace(0,pi);
vp=vo./(1+a.*(1-cos(theta)));
dOmega=r^2/2.*(vp./vo).^2.*(vo./vp+vp./vo-sin(theta).^2);
figure
plot(dOmega.*cos(theta),dOmega.*sin(theta),'r','LineWidth',2)
Which produces the following plot
I am trying to make a surface using the cylinder command like this
[X,Y,Z] = cylinder(dOmega.*sin(theta));
figure
surf(8.*Z,Y,X)
However, the profile negates the x-part of the 2D profile and instead of rounded ends it comes out pointy as shown below.
Does anyone know how to get around this? Thank you in advance for your help!
Richard

 Réponse acceptée

Richard McCulloch
Richard McCulloch le 25 Mar 2014

4 votes

OK, so I parametrized the functions and it worked out. Here is the completed code for those interested. It calculates the compton scattering cross section envelope at a specified gamma ray energy.
a=0.255/0.511;
r=2.8179;
vo=1;
m=0.511;
F = @(z)r^2/2.*((vo./(1+a.*(1-cos(z))))./vo).^2.*...
(vo./(vo./(1+a.*(1-cos(z))))+...
(vo./(1+a.*(1-cos(z))))./vo-sin(z).^2);
t = linspace(0,2*pi,55);
z = linspace(0,pi,55);
[T,U] = meshgrid(t,z);
X = F(U).*sin(U).*cos(T);
Y = F(U).*sin(U).*sin(T);
Z = F(U).*cos(U);
surf(Z,Y,X,X.^2+Y.^2+Z.^2)
shading interp
axis vis3d
axis equal
title('Differential Cross Section (10^{-26} cm^2/electron)')
This is the output
Richard

3 commentaires

Elvira Tanjung
Elvira Tanjung le 19 Jan 2016
Dear Mr. Richard, I have similar problem with ur case. I don't understand the parameter that you used in this codes, like what are the a, r, m, vo. vp? I would be very thankful if you could share with me. Thank you so much :: Vira ::
Ragavendiran R
Ragavendiran R le 2 Mai 2019
Could you please explain the parameters a,r,m,vo and vp.
Richard McCulloch
Richard McCulloch le 2 Mai 2019
This is the Klein–Nishina formula for Compton Scattering angle calculations. For the purpose of this code a, r, m, and vo are just function constants. “The Atomic Nucleus” by R. Evans has a good treatment for a more detailed understanding.

Connectez-vous pour commenter.

Plus de réponses (0)

Community Treasure Hunt

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

Start Hunting!

Translated by