diff --git a/AT/StabilityDiagram/main.m b/AT/StabilityDiagram/main.m
index b80cf576a93f3741483f948696a1c54fe180a992..61d58aacd2e11b5d07e3aadf850d2370555365c2 100644
--- a/AT/StabilityDiagram/main.m
+++ b/AT/StabilityDiagram/main.m
@@ -6,8 +6,8 @@ delta = 0;
%% Lattice to be analysed
-ring = atreduce(sr.model().ring);
-% ring = rmin;
+load('../Lattices/ESRF_standard_cell.mat')
+ring = ARCA;% ring = rmin;
%% Sextupolar resonances
if sext
diff --git a/AT/StabilityDiagram/sext_resonances.m b/AT/StabilityDiagram/sext_resonances.m
index ea8705d48c53dd3ae3b9e385b0a97b95fdba3d34..8c81cb4bec0b6a05c262aac1d2c6f76eaf9a03d4 100644
--- a/AT/StabilityDiagram/sext_resonances.m
+++ b/AT/StabilityDiagram/sext_resonances.m
@@ -4,20 +4,24 @@ function sext_resonances(ring, Denergy, varargin)
% of sextupolar resonances on the dynamic aperture at a given
% energy deviation.
%
-% Follwoing the paper "Nonlinear dynamics with sextupoles in
+% Following the paper "Nonlinear dynamics with sextupoles in
% low-emittance light source storage rings", R. Nagaoka, K. Yoshida and
% M. Hara, 1991
%
-% Aumented with the inclusion of the energy deviation component.
-
+%
+% - Resonances of the type Qx = M to be confirmed
+% - To be augmented with the inclusion of the energy deviation component.
+% - Higher-order to be included
+%
+%
% Inputs : - ring, lattice to be analysed, can full or parts
% - Denergy, energy deviation
%
% Optional outputs : - orders, stating the order of sextupole resonance
% required. DEFAULT = linear resonances 1,3.
-
-
-
+%
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%
% if nargin <3
% order = [1;3];
% end
@@ -52,6 +56,7 @@ betay=beta(:,2);
betax0 = betax(1);
betay0 = betay(1);
+% Averaged twiss functions to be closer to the thin lens approximation
[lindata,avebeta,avemu,avedisp,nu,xsi]=atavedata(RING1,0,1:length(RING1)+1);
beta = avebeta;
disp = avedisp;
@@ -60,6 +65,7 @@ betax=beta(:,1);
betay=beta(:,2);
alpha=cat(1,lindata.alpha);
alphax=alpha(:,1);
+alphay=alpha(:,2);
nux=muxy(length(RING1)+1,1)/2/pi;
nuy=muxy(length(RING1)+1,2)/2/pi;
phix=muxy(:,1);
@@ -127,8 +133,6 @@ for k=1:N
delta = Qx-k/3;
xa_l_m= delta/A3M(k)/6*sqrt(betax0); %left boundary
xa_r_m= -2*xa_l_m; %right boundary
-
-
if(xa_l_m > -b/10^3 && xa_l_m < b/10^3 )
x(1,1) = xa_l_m*10^3; y(1,1) = 0;
x(1,2) = xa_l_m*10^3; y(1,2) = a;
@@ -144,14 +148,13 @@ xa_r_m= -2*xa_l_m; %right boundary
plot(x,y, 'g')
str = ['(3,', num2str(k),')'];
text(xa_r_m*10^3, a*rand, str, 'FontSize',10,'FontWeight','bold', 'Color', [0, 0, 0]);
- end
-
+ end
end
-%% Resonances of type Qx = M
+%% Resonances of type Qx = M // extracted from CATS
x1st=[];
y1st=[];
for k=1:N
- ystep = 0.05;
+ ystep = 0.001;
m3 = 3*k;
if(m3 >= 0 && m3 <= N)
delta = Qx-k;
@@ -206,108 +209,25 @@ for k=1:N
plot(x1stN(1,1:ic-1)*1000,y1stN(1,1:ic-1)*1000,'Color',[0.8,0.1,0.1],'LineWidth',2.5);
ResLab=['[1,',num2str(k),']'];
ic1pos = round((ic1-1+1)/2) + round((ic1-1-1)/3.0*(2*rand-1.0)/3.0);
- text(x1stN(1,ic1pos),y1stN(1,ic1pos),...
- texlabel(ResLab,'literal'),'FontSize',20,'FontWeight','bold','Color',[0.8,0.1,0.1])
+ text(x1stN(1,ic1pos)*1000,y1stN(1,ic1pos)*1000,...
+ texlabel(ResLab,'literal'),'FontSize',14,'FontWeight','bold','Color',[0.8,0.1,0.1])
end
end
end
-Iy = y.^2/2/betay0;
-FF = alphaM.*(B1M./(delta.^2))*12^2/2/betay0;
-GF = lambda.*(B1M./(delta.^2))*6^2/betay0;
-%F(m,y) = @(m,y) sqrt(1+ FF(m)*y^2);
-%xa_r_1(m,y)=@(m,y) -(delta(m)/alphaM(m))'/6*(1+F(m,y))*sqrt(betax0); %left boundary
-%xa_l_1(m,y)=@(m,y) +(delta(m)/alphaM(m))/6*(2*F(m,y)-1)*sqrt(betax0); %right boundary
-
-% for k = 1:N+1,
-% C = (alphaM(k+N)/2/lambda(k+N)+1)/(GF(k+N));
-% test = sqrt(C)-a/1000;
-% if (test<0)
-% ymin = 0;
-% ymax = a/1000;
-% else
-% ymin = 0;
-% ymax = +sqrt(C);
-% end
-%
-% y = linspace(ymin, ymax, 100)*sqrt(betay0);
-% F = sqrt(1*ones(1,100)+ FF(k+N)*(y.^2));
-% xa_r_1=-(delta(k+N)/alphaM(k+N))/6*(1*ones(1,100)+F)*sqrt(betax0); %left boundary
-% xa_l_1= (delta(k+N)/alphaM(k+N))/6*(2*F-1*ones(1,100))*sqrt(betax0); %right boundary
-% plot(xa_l_1*1000, y*1000, 'r')
-% plot(xa_r_1*1000, y*1000, 'r')
-% str = ['[1,',num2str(M(k+N)),']'];
-%
-% % if (xa_r_1(1)*1000<b)&&(-b<xa_r_1(1)*1000)
-% % text(xa_r_1(1)*1000 + 0.15, a/3, str, 'FontSize',10,'FontWeight','bold', 'Color', [1, 0, 0]);
-% % end
-%
-% if (xa_l_1(1)*1000<b)&&(-b<xa_l_1(1)*1000)
-% text(xa_l_1(1)*1000 - 0.15, a/3, str, 'FontSize',10,'FontWeight','bold', 'Color', [1, 0, 0]);
-% end
-% end
-%
-%
-% % for k = 1:N+1,
-% % if (FF(k+N)>0)
-% % ymin = 0;
-% % ymax = a/1000;
-% y = linspace(ymin, ymax, 100)*sqrt(betay0);
-% F = sqrt(1*ones(1,100)+ FF(k+N)*(y.^2));
-% xa_r_1=-(delta(k+N)/alphaM(k+N))/6*(1*ones(1,100)+F)*sqrt(betax0); %left boundary
-% xa_l_1= (delta(k+N)/alphaM(k+N))/6*(2*F-1*ones(1,100))*sqrt(betax0); %right boundary
-% plot(xa_l_1*1000, y*1000, 'r')
-% plot(xa_r_1*1000, y*1000, 'r')
-% str = ['[1,',num2str(M(k+N)),']'];
-% if (xa_r_1(1)*1000<b)&&(-b<xa_r_1(1)*1000)
-% text(xa_r_1(1)*1000 + 0.15, a/3, str, 'FontSize',10,'FontWeight','bold', 'Color', [1, 0, 0]);
-% end
-%
-% if (xa_l_1(1)*1000<b)&&(-b<xa_l_1(1)*1000)
-% text(xa_l_1(1)*1000 - 0.15, a/3, str, 'FontSize',10,'FontWeight','bold', 'Color', [1, 0, 0]);
-% end
-% elseif (FF(k+N)<0)
-% ymin = 0;
-% ymax = +sqrt(-1/FF(k+N));
-% y = linspace(ymin, ymax, 100)*sqrt(betay0);
-% F = sqrt(1*ones(1,100)+ FF(k+N)*(y.^2));
-% xa_r_1=-(delta(k+N)/alphaM(k+N))/6*(1*ones(1,100)+F)*sqrt(betax0); %left boundary
-% xa_l_1= (delta(k+N)/alphaM(k+N))/6*(2*F-1*ones(1,100))*sqrt(betax0); %right boundary
-% plot(xa_l_1*1000, y*1000, 'r');
-% plot(xa_r_1*1000, y*1000, 'r');
-% str = ['[1,',num2str(M(k+N)),']'];
-%
-% if (xa_r_1(1)*1000<b)&&(-b<xa_r_1(1)*1000)
-% text(xa_r_1(1)*1000 + 0.15, a/3, str, 'FontSize',10,'FontWeight','bold', 'Color', [1, 0, 0]);
-% end
-%
-% if (xa_l_1(1)*1000<b)&&(-b<xa_l_1(1)*1000)
-% text(xa_l_1(1)*1000 - 0.15, a/3, str, 'FontSize',10,'FontWeight','bold', 'Color', [1, 0, 0]);
-% end
-% end
-% end
-%
-
-
-
-
-
-%% Fichiers de sortie
+%% Output files
fid = fopen('resonances_sextupoles_out.txt','w');
-%ouvre un fichier ou le creer
-%�crit dans ce fichier, fid est sa reference pour matlab
fprintf(fid,'%s\n','Betatron tune (Qx, Qy)');
fprintf(fid,'%i\t %i%i\n',Qx);
fprintf(fid,'%i\n\n',Qy);
fprintf(fid,'%s\n','Initial Twiss and dispersion alphax betax alphay betay disp0 disp0p');
-fprintf(fid,'%i\t %i\t %i\t %i\n',alphax(1), betax0, alpha(2,1), betay0);
+fprintf(fid,'%i\t %i\t %i\t %i\n',alphax(1), betax0, alphay(1), betay0);
fprintf(fid,'%i\t %i\n\n',disp(1,1), disp(2,1));
fprintf(fid,'%s\n','Nonlinear magnets list');
-
fprintf(fid,'%s\n','Sextupoles');
fprintf(fid,'%s\n','Index Phix/2pi Phiy/2pi SK SKB s/R');
for k=1:n_sext
@@ -315,17 +235,6 @@ SK = betax(indice(k))^(3/2)*k2l(k);
SKB = betax(indice(k))^(1/2)* betay(indice(k))*k2l(k);
fprintf(fid,'%i\t %i\t %i\t %i\t %i\t %i\n\n', 3, phix(indice(k)), phiy(indice(k)), SK, SKB, s(indice(k))/R);
end
-%
-% fprintf(fid,'%s\n','Octupoles');
-% fprintf(fid,'%s\n','Index Phix/2pi Phiy/2pi TX TY TB s/R');
-% for k=1:n_oct,
-% TX = betax(indice_oct(k))^(2)*k3l(indice_oct(k));
-% TY = betay(indice_oct(k))^(2)*k3l(indice_oct(k));
-% TB = betax(indice_oct(k))*betay(indice_oct(k))*k3l(indice_oct(k));
-% fprintf(fid,'%i\t %i\t %i\t %i\t %i\t %i\t %i\n\n', 4, phix(indice_oct(k)), phiy(indice_oct(k)), TX, TY, TB, s(indice_oct(k))/R);
-% end
-%n'oublie pas de fermer le fichier sinon tu ne peux pas le lire
fclose(fid)
-save 'test.txt'