TITLE,ā€¯ThomX_10/11/17 from Irynaā€™s ATā€¯;

option, warn, -info, echo;
Assign, echo=ThomX.echo;

E0:= 0.05;
beam, particle= ELECTRON, EX= 5.0E-08, EY = 5.0E-08, ENERGY= E0; !, RADIATE; !for EMIT

//Length of the ring
//L0 =  17.986715999999994;


//Dipoles, fringe as AT

BEND: SBEND, L=0.296351, ANGLE= 0.785398, E1=0.0332, E2=0.0332, FINT=0.5, hgap=0.01757532; 


!sept1: SBEND, L=0.251072, ANGLE= -0.16, E1=-0.16, E2=0, FINT= 0, HGAP= 0;
!kick1: SBEND, L=0.200004, ANGLE= -0.011, E1=-0.011, E2=0, FINT= 0, HGAP= 0;
!kick2: SBEND, L=0.200004, ANGLE= -0.011, E1=0, E2=-0.011, FINT= 0, HGAP= 0;
!sept2: SBEND, L=0.251072, ANGLE= -0.16, E1=0, E2=-0.16, FINT= 0, HGAP= 0;

//Drifts

dL := (0.296233-0.27646)/2;
dL_track := (0.296351 - 0.296233)/2;
dL_quad := (0.15731-0.15)/2;

SD0: drift, l:= 0.1000000E+00-dL - dL_track; 
SD1S1B : drift, l:= 0.600000E-01 - dL_quad; 
SD1S1: drift, l:=  0.190000E+00-LSX/2;
SD1S: drift, l:=  0.7500000E-01-LSX/2 - dL_quad;
SD3S2: drift, l:= 0.1250000E+00-LSX/2 - dL - dL_track;
SD2: drift, l:= 0.2100000E+00 - 2*dL_quad;
SD31: drift, l:= 0.7537500E+00 - dL - dL_quad - dL_track;
SD5B: drift, l:= 0.74;  !!contains RF/2, old= 0.6400000E+00; 
SD5: drift, l:=  0.0600000E+00 - dL_quad;
SD3S: drift, l:=  0.7500000E-01-LSX/2 - dL_quad;
SD3S1: drift, l:= 0.1350000E+00-LSX/2 - dL_quad;
SD2SB: drift, l:= 0.4000000E-01 - dL_quad;
SD2S: drift, l:=  0.1950000E+00 - dL - dL_track;
SD5Brf: drift, l:= 0.64;  !!contains RF/2, old= 0.7400000E+00; 



//Quadrupoles

QP1: quadrupole, l=0.15731, k1= -.4861690E+01;
QP2: quadrupole, l=0.15731, k1= 0.9441459E+01;
QP3: quadrupole, l=0.15731, k1= -.1854283E+02;
QP4: quadrupole, l=0.15731, k1= 0.1646009E+02;
QP31: quadrupole, l=0.15731, k1= -.8900146E+01;
QP41: quadrupole, l=0.15731, k1= 0.5107448E+01;

//Sextupoles (flag ON or OFF)

SX_ON = 1.;
lse:=1e-06;

!ā€¯thickā€¯ sextupoles
!k2SX1 :=  -325568.10*2*SX_ON;
!k2SX2 :=  2896110.38*2*SX_ON;
!k2SX3 := -8674118.16*2*SX_ON;

!SX1: sextupole, l=lse, k2= k2SX1; 
!SX2: sextupole, l=lse, k2= k2SX2;
!SX3: sextupole, l=lse, k2= k2SX3;

!Multipole/sextupoles
!AT sextupole values result in Cx=0.1, Cy=0.14
k2SX1 := -325568.10*lse*2*SX_ON;
k2SX2 :=  2896110.38*lse*2*SX_ON;
k2SX3 := -8674118.16*lse*2*SX_ON;

!Zero chromaticity in MADX: (k2*L)
!k2sx1              =       -0.265912254 ;
!k2sx2              =        5.743932173 ;
!k2sx3              =       -17.34823632 ;

SX1: multipole, lrad=0., knl:={0,0,k2SX1,0,0};
SX2: multipole, lrad=0., knl:={0,0,k2SX2,0,0};
SX3: multipole, lrad=0., knl:={0,0,k2SX3,0,0};

value,k2SX1,K2SX2,K2SX3;

!RF Cavity (flag ON or OFF)
RF_ON := 0;
RFVOLT:= 0.3*RF_ON;
RFC:   RFCAVITY, L=0.20, VOLT:=RFVOLT, LAG=0.5, HARMON=30;

//Markers

RING.END: Marker;
SEPT: MARKER;
Inj_KICKER: MARKER;
RING.START: MARKER;
Mqp31: MARKER;
Msd1s1: MARKER;
BPMx: marker;
BPMz: marker;
IP: marker;
HCOR: marker;
VCOR: marker;
//RFC: marker;

str1: line=(SD5B, BPMx, BPMz, SD5, QP1, SD2, QP2, SD31);
str1rf: line=(SD5Brf, BPMx, BPMz, SD5, QP1, SD2, QP2, SD31);
str2: line=(SD31, QP2, SD2, QP1, SD5, BPMx, BPMz, SD5B);
str2rf: line=(SD31, QP2, SD2, QP1, SD5, BPMx, BPMz, SD5Brf);

cell1: line=(BEND, SD3S2, SX1, HCOR, VCOR, SD3S, QP31, SD3S, SX3, HCOR, VCOR, SD3S1, QP41, SD1S1B, BPMx, BPMz, SD1S1, SX2, HCOR, VCOR, SD1S, QP4, SD2, QP3, SD2SB, BPMx, BPMz, SD2S, BEND, SD0); 
cell2: line=( SD0, BEND, SD2S, BPMx, BPMz, SD2SB, QP3, SD2, QP4, SD1S, SX2, HCOR, VCOR, SD1S1, BPMx, BPMz, SD1S1B, QP41, SD3S1, SX3, HCOR, VCOR, SD3S, QP31, SD3S, SX1, HCOR, VCOR, SD3S2, BEND);

 
RING: line=(RING.START, SEPT, str1, cell1, IP, cell2, str2rf, RFC, str1rf, cell1, cell2, str2,RING.END);


use, sequence=RING;

//AT values: betx=3.455893, alfx=0.3167609E-15 , bety=2.662101,alfy=0.2883261E-15,
//dx=0.0290735, dpx=-0.1604619E-15, 

select,flag = twiss,column = name,s,betx,alfx,dx,bety,alfy;
twiss,  save, file = ThomX.twiss;

plot, table=twiss, haxis= s, vaxis1= betx,bety, vaxis2 = dx, colour= 100, file= ThomX.plots;



!-------- MADX CLASSIC MATCHING (FASTER, NOT ROBUST, CHROMA ONLY) -------

match, sequence= RING;
       vary, name=K2SX1, step=0.01;
       vary, name=K2SX2, step=0.01;
!        vary, name=K2SX3, step=0.01;
        global, sequence=RING, DQ1=0.0;
        global, sequence=RING, DQ2=0.0;
        lmdif, calls = 1000, tolerance = 1.e-21;
        migrad, calls = 1000, tolerance = 1.e-21;
ENDMATCH;
value,K2SX1,K2SX2,K2SX3;

!TWISS, SAVE;

//Adds calculation of Hx invariant and sigmas to Twiss table

select, flag = twiss,full, column = sigx,sigpx,sigy, sigpy,Hinvx;

sigx:= sqrt(beam->ex*table(twiss,betx));
sigpx:= sqrt(beam->ex/table(twiss,betx));
sigy:= sqrt(beam->ey*table(twiss,bety));
sigpy:= sqrt(beam->ey/table(twiss,bety));

Gammax:= (1+table(twiss,alfx)*table(twiss,alfx)/table(twiss,betx));

Hinvx:= Gammax*table(twiss,dx)*table(twiss,dx)+2*table(twiss,alfx)*table(twiss,dx)*table(twiss,dpx)+table(twiss,betx)*table(twiss,dpx)*table(twiss,dpx);

twiss,  save, file = ThomX.sigma;
plot, table=twiss, haxis= s, vaxis= sigx, sigy, colour= 100;
plot, table=twiss, haxis= s, vaxis= Hinvx, colour= 100;

//Note: to have the emittance calculation parameter LAG in the RF needs to be NOT zero!
EMIT;


!=======================================================================!
! (Chromatic) TWISS WITH MADX output in twiss* detuning* chroma*        !
!=======================================================================!
!/*
!**** TUNES, BETAS AND ALFAS @ SS CENTRE Vs ENERGY DEVIATION (+-2.0%)****
de_t   := 0.02;  ! scan range
DP     := -de_t; ! scan start
NMAX   =  20;    ! scan points

SAVEBETA, label=TWSSIP, place=#e;

option,-info,-warn;   ! otherwise you get redifiniton warnings
n = 1; 
create,table=twiss_end,column=DELTA_P,bx,by,ax,ay;
create,table=twiss_tune,column=DELTA_P,Qx,Qy;
while (n < NMAX+1)
{
  TWISS, DELTAP=DP;
  DELTA_P :=DP*1e2;
  bx:=table(twiss,RING.END,betx); 
  by:=table(twiss,RING.END,bety); 
  ax:=table(twiss,RING.END,alfx); 
  ay:=table(twiss,RING.END,alfy); 
  Qx:=table(twiss,RING.END,mux)-3; 
  Qy:=table(twiss,RING.END,muy)-1; 
  if(Bx > 1e-100 && By > 1e-100)
  {
   fill,table=twiss_end;
   fill,table=twiss_tune;
  }
  DP:=-de_t*(1-2/NMAX*n);
  n  = n + 1;
};

PLOT, TABLE=twiss_tune, HAXIS=DELTA_P, VAXIS1=QX, VAXIS2=QY, STYLE=100,
      symbol=4,!VMIN=0.10,0.60,VMAX=0.30,0.80,
      TITLE="Tune Vs Dp/p [%]",noversion,FILE="RING-plots";

PLOT, TABLE=twiss_end, HAXIS=DELTA_P, VAXIS1=BX, VAXIS2=BY, STYLE=100,
      symbol=4,!VMIN=3.0,1.4,VMAX=6.0,3.0,
      TITLE="Beta Vs Dp/p [%], SS centre",noversion,FILE="RING-plots";

PLOT, TABLE=twiss_end, HAXIS=DELTA_P, VAXIS1=AX, VAXIS2=AY, STYLE=100,
      symbol=4,!VMIN=3.0,1.4,VMAX=6.0,3.0,
      TITLE="Alfa Vs Dp/p [%], SS centre",noversion,FILE="RING-plots";

write,table=twiss_end;
write,table=twiss_tune; write,table=twiss_tune,file=chroma_tunes.dat;


!***************** TUNES Vs INITIAL CONDITIONS *************************

de_t   := 0.040;  ! scan range [m]
DP    =  -de_t;      ! scan start
NMAX   =  20;        ! scan points


TWISS, SAVE;
SAVEBETA, label=TWSSIP, place=#e;

n = 1; 
create,table=twiss_detuneAx,column=Ax,Ay,Qx,Qy;
while (n < NMAX+1)
{
  TWISS, DELTAP=0,X=DP;
  Ax :=DP*1e3;
  Ay :=0;
  Qx:=table(twiss,RING.END,mux)-3; 
  Qy:=table(twiss,RING.END,muy)-1; 
  if(Qx < 1 && Qy < 1)
   {
    value,n,Ax,Ay,Qx,Qy; 
    fill,table=twiss_detuneAx;
   }

//  DP:=de_t/NMAX*n;
 DP:=-de_t*(1-2/NMAX*n);
  n  = n + 1;
};

PLOT, TABLE=twiss_detuneAx, HAXIS=Ax, VAXIS1=QX, VAXIS2=QY, STYLE=100,
      symbol=4,!VMIN=0.10,0.60,VMAX=0.30,0.80,
      TITLE="Tunes Vs X position [mm]",noversion;


de_t   := 0.020;  ! scan range [m]
DP     := 0;      ! scan start

n = 1; 
create,table=twiss_detuneAy,column=Ax,Ay,Qx,Qy;
while (n < NMAX+1)
{
  TWISS, DELTAP=0,Y=DP;
  Ax :=0;
  Ay :=DP*1e3;
  Qx:=table(twiss,RING.END,mux)-3; 
  Qy:=table(twiss,RING.END,muy)-1; 
  if(Qx < 1 && Qy < 1)
   {
    value,n,Ax,Ay,Qx,Qy; 
    fill,table=twiss_detuneAy;
   }
  DP:=de_t/NMAX*n;
  n  = n + 1;
};

PLOT, TABLE=twiss_detuneAy, HAXIS=Ay, VAXIS1=QX, VAXIS2=QY, STYLE=100,
      symbol=4,!VMIN=0.10,0.60,VMAX=0.30,0.80,
      TITLE="Tunes Vs Y position [mm]",noversion;

write,table=twiss_detuneAx;write,table=twiss_detuneAx,file=detuning_x.dat;
write,table=twiss_detuneAy;write,table=twiss_detuneAy,file=detuning_y.dat;


*/
!=======================================================================!
! END OF (Chromatic) TWISS WITH MADX output in twiss* detuning* chroma* !
!=======================================================================!


return;

stop;