MATLAB Plotting Beam Emittance with MATLAB Code

[1Keenan]

Hi,

I want to write a very simple code to plot the beam emittance knowing x,y,z and vx,vy,vz of each particle.
I have generated with a tracking code a very simple beam (cone distribution with 1° divergence) and I want to compare the emittance of this beam with the one calculated by a different code.
I'm using MATLAB and this is the code I have written. It's not working properly, but I'm not able to spot the error. COuld you please help me?

A= Vel ; %%%Vel is a matrix of n rows and columns: [x,y,z,vx,vy,vz] 
      
      
      
%%%%%   evaluate the momentum of the reference particle

c=3e8; %[m/s] light speed
e=1.602176e-19; %[C] electron charge 

E0=938.27; %[MeV] proton rest energy
M0=E0/c^2; %[MeV/c^2] proton rest mass
mkg=M0*10^6*e; %[kg] mass in kg e

T=5; %[MeV] Kinetic energy (it is a monochromatic beam)

Etot= E0+T; %Total Energy

gamma = Etot/E0; %gamma relativistic

beta= sqrt( 1 - (1/gamma^2)); %beta relativistic = v/c

v=beta*c;  

%%%%% momentum (not taking in account particle mass%%%%
V0=v*gamma;

angX=sin(A(:,4)./V0); %%%angle in radians
angX=angX.*10^3; %%% angle in mradians
figure()
plot(A(:,1),angX,'.')

angY=sin(A(:,5)./V0); %%%angle in radians
angY=angY.*10^3; %%% angle in mradians
figure()
plot(A(:,2),angY,'.')

 

[mfb]

It's not working properly

What is wrong?

Did you check variables like v? Your protons are non-relativistic, you can check their speed with a non-relativistic calculation.
What is A?
The sin of something will never give you an angle. Do you mean the arcsin there?

 

[1Keenan]

I check variables and their dimensions. I try with classic and relativistic formulas... what is wrong is the plot. Only the Vx-Vy plot seems to be ok.

The sin was wrong, you are right, but actually nothing change as I hav very small angles...

 

[mfb]

Looks like an odd distribution of your particles, unrelated to your code.
Your code just scales everything in A.

 

[1Keenan]

I don't know, the particles should be distributed around a circle.
I have the right prifile in xy ad vx-vy, but in the x-vx I have a straight line.
What I sent you is what I get in the focal point of a quadrupole

 

[mfb]

I would not expect lines at all - you should have a 2D-distribution, otherwise your phase space volume is zero.
There is nothing fundamentally wrong with a v_x,x-distribution that is very stretched. It just means you have to focus in that direction again.

 

[1Keenan]

The v_x,x distribution is not stretched, it is a line in my plot... that's why I'm concerned...