% Dirac Talbot with octave


%

NT = 512; % Number of nodes in time and angle

% Momentum list
momL = [[1,2]; [2,5]; [2,11]]
n1t = momL(:,1);
n2t = momL(:,2);
momLt = [[n1t,n2t];[n1t,-n2t];[-n1t,n2t];[-n1t,-n2t]];
momL = [[momLt(:,1),momLt(:,2)];[momLt(:,2),momLt(:,1)]];

% q^2
Q = 3

% time list
tL = [0.0,0.1,0.2,0.3]
tL = [0.5-0.05,0.5]
tL = [0.0,0.05]
%tL = [0.5,0.5-0.05]


trev = 2*pi/sqrt(momL(1,1)^2+momL(1,2)^2+Q)


th1 = linspace(-pi,pi,NT);
th2 = linspace(-pi,pi,NT);
[TH1,TH2] = meshgrid(th1,th2);

Nt = length(tL) 


for k = 1:Nt
	t = tL(k)
	S1 = 0;
	S2 = 0;
	for j = 1:size(momL,1)
		n1 = momL(j,1);
		n2 = momL(j,2);
		wm = sqrt(n1^2+n2^2+Q);
		% normalization coefficient
		normc = sqrt(1+ sqrt(Q)/wm)/sqrt(2);
		normc2 = normc*(n1+i*n2)/(sqrt(Q)+wm);
		S1 = S1 + normc*exp(i*n1*TH1+i*n2*TH2-i*wm*t*trev);
		S2 = S2 + normc2*exp(i*n1*TH1+i*n2*TH2-i*wm*t*trev);
	end

	R = sqrt(abs(S1).^2 + abs(S2).^2);


	figure(k)
	surf(TH1,TH2,R)
	colormap(jet(256));
	%colormap(gray);
	shading interp;
	view(0,90)
	%colorbar
	axis([-pi pi -pi pi])
	axis('square','off')
%	set(gca,'xtick',[])
%	set(gca,'xticklabel',[])
%	set(gca,'ytick',[])
%	set(gca,'yticklabel',[])
	grid off
	
	% https://www.mathworks.com/matlabcentral/answers/402194-how-to-save-a-plot-without-margin-of-figure
	set(gca, 'units', 'normalized'); %Just making sure it's normalized
	Tight = get(gca, 'TightInset');  %Gives you the bording spacing between plot box and any axis labels, [Left Bottom Right Top] spacing
	NewPos = [Tight(1)*0 Tight(2) 1-Tight(1)-Tight(3) 1-Tight(2)-Tight(4)]; %New plot position [X Y W H]
	set(gca, 'Position', NewPos);
	
	namefig = './fig9_%d.png';
	namefig = sprintf(namefig,k);
%	saveas(gca,namefig)

end

%pos = get(gcf,'Position');
%f = 1;
%f = 3;
%set(gcf, 'Position', [-2,-1,pos(3),pos(3)]); 
saveas(gcf,'borri.png')


