-
Notifications
You must be signed in to change notification settings - Fork 0
/
BeamOpticsSupportFunctions4D.html
280 lines (280 loc) · 73.6 KB
/
BeamOpticsSupportFunctions4D.html
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
<!DOCTYPE html PUBLIC "-//W3C//DTD HTML 4.01 Transitional//EN">
<html><head><meta http-equiv="Content-Type" content="text/html; charset=utf-8"><meta http-equiv="X-UA-Compatible" content="IE=edge,IE=9,chrome=1"><meta name="generator" content="MATLAB 2021a"><title>Beam optics support functions 4D (Section 3.5)</title><style type="text/css">.rtcContent { padding: 30px; } .S0 { margin: 2px 10px 9px 4px; padding: 0px; line-height: 21px; min-height: 0px; white-space: pre-wrap; color: rgb(0, 0, 0); font-family: Helvetica, Arial, sans-serif; font-style: normal; font-size: 14px; font-weight: normal; text-align: left; }
.S1 { margin: 15px 10px 5px 4px; padding: 0px; line-height: 28.8px; min-height: 0px; white-space: pre-wrap; color: rgb(213, 80, 0); font-family: Helvetica, Arial, sans-serif; font-style: normal; font-size: 24px; font-weight: normal; text-align: left; }
.S2 { margin: 15px 10px 5px 4px; padding: 0px; line-height: 18px; min-height: 0px; white-space: pre-wrap; color: rgb(60, 60, 60); font-family: Helvetica, Arial, sans-serif; font-style: normal; font-size: 17px; font-weight: bold; text-align: left; }
.S3 { margin: 10px 0px 20px; padding-left: 0px; font-family: Helvetica, Arial, sans-serif; font-size: 14px; }
.S4 { margin-left: 56px; line-height: 21px; min-height: 0px; text-align: left; white-space: pre-wrap; }
.CodeBlock { background-color: #F7F7F7; margin: 10px 0 10px 0;}
.S5 { border-left: 1px solid rgb(233, 233, 233); border-right: 1px solid rgb(233, 233, 233); border-top: 1px solid rgb(233, 233, 233); border-bottom: 0px none rgb(0, 0, 0); border-radius: 4px 4px 0px 0px; padding: 6px 45px 0px 13px; line-height: 17.234px; min-height: 18px; white-space: nowrap; color: rgb(0, 0, 0); font-family: Menlo, Monaco, Consolas, "Courier New", monospace; font-size: 14px; }
.S6 { border-left: 1px solid rgb(233, 233, 233); border-right: 1px solid rgb(233, 233, 233); border-top: 0px none rgb(0, 0, 0); border-bottom: 0px none rgb(0, 0, 0); border-radius: 0px; padding: 0px 45px 0px 13px; line-height: 17.234px; min-height: 18px; white-space: nowrap; color: rgb(0, 0, 0); font-family: Menlo, Monaco, Consolas, "Courier New", monospace; font-size: 14px; }
.S7 { border-left: 1px solid rgb(233, 233, 233); border-right: 1px solid rgb(233, 233, 233); border-top: 0px none rgb(0, 0, 0); border-bottom: 1px solid rgb(233, 233, 233); border-radius: 0px 0px 4px 4px; padding: 0px 45px 4px 13px; line-height: 17.234px; min-height: 18px; white-space: nowrap; color: rgb(0, 0, 0); font-family: Menlo, Monaco, Consolas, "Courier New", monospace; font-size: 14px; }
.S8 { margin: 3px 10px 5px 4px; padding: 0px; line-height: 18px; min-height: 0px; white-space: pre-wrap; color: rgb(60, 60, 60); font-family: Helvetica, Arial, sans-serif; font-style: normal; font-size: 17px; font-weight: bold; text-align: left; }
.S9 { margin: 10px 10px 9px 4px; padding: 0px; line-height: 21px; min-height: 0px; white-space: pre-wrap; color: rgb(0, 0, 0); font-family: Helvetica, Arial, sans-serif; font-style: normal; font-size: 14px; font-weight: normal; text-align: left; }</style></head><body><div class = rtcContent><div class = 'S0'><span>Companion software for "Volker Ziemann, </span><span style=' font-style: italic;'>Hands-on Accelerator physics using MATLAB, CRCPress, 2019</span><span>" (https://www.crcpress.com/9781138589940)</span></div><h1 class = 'S1'><span>Beam optics support functions 4D (Section 3.5)</span></h1><div class = 'S0'><span>Volker Ziemann, 211119</span></div><div class = 'S0'><span>In this live script we define the functions for the 4D beam optics calculations, such as </span><span style=' font-family: monospace;'>calcmat()</span><span> that a frequently used in other calculations. All described functions reside in the subdirectory </span><span style=' font-family: monospace;'>4D</span><span> that is contained in the archive </span><span style=' font-family: monospace;'>BeamOpticsSupportFile.zip</span><span>. Any scripts using these function need to include that subirectory with the command "</span><span style=' font-family: monospace;'>addpath ./4D</span><span>".</span></div><h3 class = 'S2'><span>The function </span><span style=' font-family: monospace;'>calcmat()</span><span> to calculate all transfer matrices</span></h3><div class = 'S0'><span>The following function receives the </span><span style=' font-family: monospace;'>beamline</span><span> description as input and returns</span></div><ul class = 'S3'><li class = 'S4'><span>Racc(4,4,nmat): transfer matrices from the start to the each of each segment, such that R(:,:,end) is the transfer matrix from the start to the end of the beamline.</span></li><li class = 'S4'><span>spos: position along the beamline after each segment, useful when plotting.</span></li><li class = 'S4'><span>nmat: number of segments </span></li><li class = 'S4'><span>nlines: number of lines in the beamline </span></li></ul><div class="CodeBlock"><div class="inlineWrapper"><div class = 'S5'><span style="white-space: pre"><span style="color: rgb(14, 0, 255);">function </span><span >[Racc,spos,nmat,nlines]=calcmat(beamline)</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >ndim=size(DD(1),1); </span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >nlines=size(beamline,1); </span><span style="color: rgb(2, 128, 9);">% number of lines in beamline</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >nmat=sum(beamline(:,2))+1; </span><span style="color: rgb(2, 128, 9);">% sum over repeat-count in column 2</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >Racc=zeros(ndim,ndim,nmat); </span><span style="color: rgb(2, 128, 9);">% matrices from start to element-end</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >Racc(:,:,1)=eye(ndim); </span><span style="color: rgb(2, 128, 9);">% initialize first with unit matrix</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >spos=zeros(nmat,1); </span><span style="color: rgb(2, 128, 9);">% longitudinal position</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >ic=1; </span><span style="color: rgb(2, 128, 9);">% element counter</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span style="color: rgb(14, 0, 255);">for </span><span >line=1:nlines </span><span style="color: rgb(2, 128, 9);">% loop over input elements</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > </span><span style="color: rgb(14, 0, 255);">for </span><span >seg=1:beamline(line,2) </span><span style="color: rgb(2, 128, 9);">% loop over repeat-count </span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > ic=ic+1; </span><span style="color: rgb(2, 128, 9);">% next element </span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > Rcurr=eye(4); </span><span style="color: rgb(2, 128, 9);">% matrix in next element</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > </span><span style="color: rgb(14, 0, 255);">switch </span><span >beamline(line,1) </span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > </span><span style="color: rgb(14, 0, 255);">case </span><span >1 </span><span style="color: rgb(2, 128, 9);">% drift</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > Rcurr=DD(beamline(line,3));</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > </span><span style="color: rgb(14, 0, 255);">case </span><span >2 </span><span style="color: rgb(2, 128, 9);">% thin quadrupole</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > Rcurr=Q(beamline(line,4)); </span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > </span><span style="color: rgb(14, 0, 255);">case </span><span >4 </span><span style="color: rgb(2, 128, 9);">% sector dipole</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > phi=beamline(line,4)*pi/180; </span><span style="color: rgb(2, 128, 9);">% convert to radians</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > rho=beamline(line,3)/phi;</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > Rcurr=SB(beamline(line,3),rho);</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > </span><span style="color: rgb(14, 0, 255);">case </span><span >5 </span><span style="color: rgb(2, 128, 9);">% thick quadrupole</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > Rcurr=QQ(beamline(line,3),beamline(line,4)); </span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > </span><span style="color: rgb(14, 0, 255);">case </span><span >20 </span><span style="color: rgb(2, 128, 9);">% coordinate roll</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > Rcurr=ROLL(beamline(line,4)); </span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > </span><span style="color: rgb(14, 0, 255);">otherwise</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > disp(</span><span style="color: rgb(170, 4, 249);">'unsupported code'</span><span >)</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > </span><span style="color: rgb(14, 0, 255);">end</span><span > </span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > Racc(:,:,ic)=Rcurr*Racc(:,:,ic-1); </span><span style="color: rgb(2, 128, 9);">% concatenate </span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > spos(ic)=spos(ic-1)+beamline(line,3); </span><span style="color: rgb(2, 128, 9);">% position of element </span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > </span><span style="color: rgb(14, 0, 255);">end</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span style="color: rgb(14, 0, 255);">end</span></span></div></div><div class="inlineWrapper"><div class = 'S7'><span style="white-space: pre"><span style="color: rgb(14, 0, 255);">end</span></span></div></div></div><h3 class = 'S8'><span>Transfer matrix for a drift space </span><span style=' font-family: monospace;'>DD(L)</span></h3><div class = 'S0'><span>The function </span><span style=' font-family: monospace;'>DD()</span><span> receives the length</span><span style=' font-family: monospace;'> L </span><span>of a drift space and resturns the 4x4 transfer matrix</span><span style=' font-family: monospace;'> out</span><span> for a drift space.</span></div><div class="CodeBlock"><div class="inlineWrapper"><div class = 'S5'><span style="white-space: pre"><span style="color: rgb(14, 0, 255);">function </span><span >out=DD(L)</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >out=eye(4);</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >out(1,2)=L;</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >out(3,4)=L;</span></span></div></div><div class="inlineWrapper"><div class = 'S7'><span style="white-space: pre"><span style="color: rgb(14, 0, 255);">end</span></span></div></div></div><h3 class = 'S2'><span>Transfer matrix for a thin-lens quadrupole </span><span style=' font-family: monospace;'>Q(F)</span></h3><div class = 'S0'><span>The function </span><span style=' font-family: monospace;'>Q()</span><span> receives the focal length</span><span style=' font-family: monospace;'> F </span><span>as input and returns the 4x4 transfer matrix</span><span style=' font-family: monospace;'> out</span><span> for a thin-lens quadrupole.</span></div><div class="CodeBlock"><div class="inlineWrapper"><div class = 'S5'><span style="white-space: pre"><span style="color: rgb(14, 0, 255);">function </span><span >out=Q(F)</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >out=eye(4);</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span style="color: rgb(14, 0, 255);">if </span><span >abs(F)<1e-8, </span><span style="color: rgb(14, 0, 255);">return</span><span >; </span><span style="color: rgb(14, 0, 255);">end</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >out(2,1)=-1/F;</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >out(4,3)=1/F;</span></span></div></div><div class="inlineWrapper"><div class = 'S7'><span style="white-space: pre"><span style="color: rgb(14, 0, 255);">end</span></span></div></div></div><h3 class = 'S2'><span>Transfer matrix for a thick quadrupole </span><span style=' font-family: monospace;'>Q(F)</span></h3><div class = 'S0'><span>The function </span><span style=' font-family: monospace;'>QQ()</span><span> receives the length</span><span style=' font-family: monospace;'> L </span><span>and </span><span style=' font-family: monospace;'>k1 </span><span>as input and returns the 4x4 transfer matrix</span><span style=' font-family: monospace;'> out</span><span> for a thick quadrupole.</span></div><div class="CodeBlock"><div class="inlineWrapper"><div class = 'S5'><span style="white-space: pre"><span style="color: rgb(14, 0, 255);">function </span><span >out=QQ(L,k)</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >ksq=sqrt(abs(k));</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >out=eye(4);</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span style="color: rgb(14, 0, 255);">if </span><span >abs(k) < 1e-6</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > out(1,2)=L;</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > out(3,4)=L;</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span style="color: rgb(14, 0, 255);">else</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > A=[cos(ksq*L),sin(ksq*L)/ksq;-ksq*sin(ksq*L),cos(ksq*L)]; </span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > B=[cosh(ksq*L),sinh(ksq*L)/ksq;ksq*sinh(ksq*L),cosh(ksq*L)];</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > </span><span style="color: rgb(14, 0, 255);">if </span><span >k>0</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > out(1:2,1:2)=A;</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > out(3:4,3:4)=B;</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > </span><span style="color: rgb(14, 0, 255);">else</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > out(1:2,1:2)=B;</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > out(3:4,3:4)=A; </span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > </span><span style="color: rgb(14, 0, 255);">end</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span style="color: rgb(14, 0, 255);">end</span></span></div></div><div class="inlineWrapper"><div class = 'S7'><span style="white-space: pre"><span style="color: rgb(14, 0, 255);">end</span></span></div></div></div><h3 class = 'S2'><span>Transfer matrix for a sector dipole </span><span style=' font-family: monospace;'>SB(L,rho)</span></h3><div class = 'S0'><span>The function </span><span style=' font-family: monospace;'>SB()</span><span> receives the length</span><span style=' font-family: monospace;'> L</span><span> and bending radius </span><span style=' font-family: monospace;'>rho</span><span> of a horizontally deflecting sector dipole magnet and returns its 4x4 transfer matrix</span><span style=' font-family: monospace;'> out</span><span>.</span></div><div class="CodeBlock"><div class="inlineWrapper"><div class = 'S5'><span style="white-space: pre"><span style="color: rgb(14, 0, 255);">function </span><span >out=SB(L,rho)</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >phi=L/rho;</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >out=eye(4);</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >out(3,4)=L;</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span style="color: rgb(14, 0, 255);">if </span><span >abs(phi)<1e-8</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > out(1,2)=L;</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span style="color: rgb(14, 0, 255);">else</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > out(1:2,1:2)=[cos(phi),rho*sin(phi); </span><span style="color: rgb(14, 0, 255);">...</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > -sin(phi)/rho,cos(phi)]; </span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span style="color: rgb(14, 0, 255);">end</span></span></div></div><div class="inlineWrapper"><div class = 'S7'><span style="white-space: pre"><span style="color: rgb(14, 0, 255);">end</span></span></div></div></div><h3 class = 'S2'><span>Transfer matrix for coordinate rotation </span><span style=' font-family: monospace;'>ROLL(phi)</span></h3><div class = 'S0'><span>The function </span><span style=' font-family: monospace;'>ROLL()</span><span> receives the roll angle</span><span style=' font-family: monospace;'> phi</span><span> (in degree) around the s-direction as input and returns the corresponding 4x4 transfer matrix </span><span style=' font-family: monospace;'>out</span><span>.</span></div><div class="CodeBlock"><div class="inlineWrapper"><div class = 'S5'><span style="white-space: pre"><span style="color: rgb(14, 0, 255);">function </span><span >out=ROLL(phi) </span><span style="color: rgb(2, 128, 9);">% phi in degree</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >c=cos(phi*pi/180); s=sin(phi*pi/180);</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >out=zeros(4);</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >out(1,1)=c; out(1,3)=s; out(2,2)=c; out(2,4)=s;</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >out(3,1)=-s; out(3,3)=c; out(4,2)=-s; out(4,4)=c;</span></span></div></div><div class="inlineWrapper"><div class = 'S7'><span style="white-space: pre"><span style="color: rgb(14, 0, 255);">end</span></span></div></div></div><h3 class = 'S2'><span>R2beta()</span></h3><div class = 'S0'><span>The function </span><span style=' font-family: monospace;'>R2beta()</span><span> receives a transfer matrix </span><span style=' font-family: monospace;'>R</span><span> as input and returns the "tune" </span><span texencoding="Q=\mu/2\pi" style="vertical-align:-5px"><img src="data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAIIAAAAmCAYAAAAIjkMFAAAICUlEQVR4Xu2aZaw0RRBFz4e7u7sGD+7u7hYIBHdNsODu7u4OQRPcIbi7u0OA4JJDepJhGOndt/t2B6b+vORtT3d1dXXVrVs9hEYaCwBDGis0FtACjSM0fvC3BRpHaByhcYRB9oGhgKeBzYAnB3ntyuWaiFBpoo4NWAq4ERgf+L5js3ZoosYROmTIiGkuBIYBNowYO+hD2nWExYDVgMWBSYCRgA+Be4DTQggc9M308YLa51NgbeD2HD09h62BFYFZgZGB54M9jwJ+7PbeWnWEJYEjgFmAC4DzgVeBUYBFgeOAccOYfYA/u72Bmsy/QbDNxMDvGZ1NFVcDCxfs5S1gFeDFbu411hGGA44EdgI+BlYFnshRbDLg5RAh9gcO6qbyNZr7NuAVYJccnW8C5geOB54BRgBmA7YHxgzjBZnzAr92a88xjmBeuwpYHfgGmBt4o0ShY4Fdg+dPD7zZLeVrMq833rQ5D/BURmdTxUnBph9kfjOyPgJMHf6/PnBFt/Zc5QhDA5cC6wYFVgD07jJxczqOchhgivg/y87AlsBMOUa4HjgbuLXAQGKGm8NvxwB7dMuQVY5geD8gLH4lsF6EIoawR8M4N7FyxDf/5SFyBtcAh2c2OWy4VJaVRSJo/C4Qf17IjbplqDJH0IPNTeKDn4AZgHcjFFk+5eGOnyLim34ZMhowKvAz8EWOUoZpU6UVgGmySmYMIE8bvJczWJLpj4pJvgTGioyu6j9pAO9l0z6eXbfMEe5PIVmrAxmxGHHcuWHgV8DYMR/1yZjLQ9Q7A9gmo5PA7XPAdGn1dHeEzqbGBQDL7XZk+FA6ek7OcV/BJBsDuwWQWbXO28BU2UFFjmAt+2xqsEBHL4qRE0J14VgriLzcmDfP0QEpx6xRNsbb541tR6yIJgjOYCpMy5ohxBsddQr/lom2fSdUTsnFaFWnJM1acRhdsmJEORnYtoWJT82zc5EjWMoIcpSXgJlbWMgSyPJHuQTQW2MkvWbM+KIxHmQ7jmCFo8GVCYFPMgucGUDfncDSEQp6gwXW6vNtxPi8IacA2wHiiLtyBhgF5HAE56YsL6+gcs8wVhbTM1gIeKjKa/N+N3wkud0DshyMEUseDyFxsE2Ai2I+7IMxWwGmhLwo5s3zdpt/9wqEWZXK5wDm7HWqBhb87lpewosLbrwOoIPcEL4/JFRoCwIPh//JYspRmJ5/a9URRLOGPTevWMIUlTfZuSVMZBcVa2dz0S9tGmKwP0vwwek5hk8DYHmUPDItra+kkBHFaChh1I74nfNYslcRSUYDgbn4ZaIABI1q74f15YBKJS81iIzThNHkBYg3b2Jp0AQT6BTihbpIgg/kTBIeJNHdEL9cqBS8XVVI3yhgLvZQqg4xzz6GdsP6IpFpRabXyCAnIWeh7A6Iu8QPOnfLjjAN8Hr4yl6ByDVmM8umGiqi6ymBH6oU6JPfLY1NCUoWY4gHzLPjAdcBgsYq8TabSnaoGpjz+1oBYIoxPov8Xr7GyJ0m/J4LPSHPs5LdzYsI6dRgPT1iRPPINCJItBml0ywD3Bu5iWRYL6sGO3/emiwwNuReFmjeOSNv1zjAR6H0fqxFG1iWWgXofKbWGBFL6HReOjGaZzZ74IByS8W8SYuqhteAacMHhjfDZpnYIHEDyubAeTE7yIzpZdUgh29KyOIDO63eLCOCtrKy0DZVttgRmK5FG8g3aDfxiAeYJ5atPmpJR+gTAdcTnG4RPrLZt18AvwkfYoQumrfwzaIhzWaIYqi6tmRT5jFLKiNJXXsLAjubQ+7bki3Zt7hAutzc6y23jWyK8P9F1ZBRwFB9cAuOMBcgWJWOt62fJ0Zm57UlnaTclQI2ECTOEaKy394RKgrLS8G7Ou+dogT+NX9RRBCtGib1ImlmkXK2j+5k5iQ5cJXUiBqsbiJR414VSSS7fOIAAZt5WsDrTXOfRjsf3six2APIilHUgzQv+44gRuRofNBj9EmTeMm3ntEYgek0GnnDTT/2L3wDojwQgGXyjfS4oNYIYGoTvAo+C0nBMopZ7xctm/etR82jlih+I/NoONo0eJ9OkADMmM330xgN6+Em4msg+wjmawGkNbm1uJjHkOybjARYZvdhSF4iEDgxe9RxpPIFqDEiTvFiWpYmEcnyXOY37UQSY6axRGx4GREKpar7aAPGRSWGXCzh2g2TevFZwAsxO+jjMUYBb4w3zEN8MFCw1uCKRNOhAbzp/EV8v2NF5z4tk4WsktHDczTBXox4mz0DRVAot+Mclum3ZCZwPzq3gNMzspQtlSpH8GPzj2FLrGCJ4m2IeYJmzm2H6q3SudO/q6Oh1lAa01EsWl9Gz0aUt/vrTivZ7fmqHEGQIaKeD/DdnSAkhhvwcYr0p7epnyXBB74cErANRKSntdcaA5mkV9+WOYKdL8Olr5QTMRLYVzcqmJMEKb62FXXLfZvDPHxf7cp2lfLbvdp0at2k7PV5nUxcu+KbDW0gmPTVUe2kyBE8TMsOw6Z5xtBpXZz0H4o2KpjSqPsWVBn9ZqCkzLIka7cn4J7k8m01mxbq0lv5x1nkOUJSqkhupHO8r2TkDCxZfHptT8KxMlmWSoIXAUps2dRrp9C5JcpsERv1qt4XlOkr9SwdbGVVS6nCCFWbkkSK6UNUzdOL330ImjByBw5AAdk+04IVR2nPfwBrdP3TgTpC1xWswQICael1SbfaSuMItT26zireOEJn7Vnb2RpHqO3RdVbxxhE6a8/azvYXt7yDNhm7FmgAAAAASUVORK5CYII=" width="65" height="19" /></span><span> for the transfer matrix R, as well as the periodic Twiss parameters </span><span style="font-family: STIXGeneral, STIXGeneral-webfont, serif; font-style: italic; font-weight: normal; color: rgb(0, 0, 0);">α</span><span>, </span><span style="font-family: STIXGeneral, STIXGeneral-webfont, serif; font-style: italic; font-weight: normal; color: rgb(0, 0, 0);">β</span><span>, and </span><span style="font-family: STIXGeneral, STIXGeneral-webfont, serif; font-style: italic; font-weight: normal; color: rgb(0, 0, 0);">γ</span><span> following Equation 3.60.</span></div><div class="CodeBlock"><div class="inlineWrapper"><div class = 'S5'><span style="white-space: pre"><span style="color: rgb(14, 0, 255);">function </span><span >[Q,alpha,beta,gamma]=R2beta(R)</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >mu=acos(0.5*(R(1,1)+R(2,2)));</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span style="color: rgb(14, 0, 255);">if </span><span >(R(1,2)<0), mu=2*pi-mu; </span><span style="color: rgb(14, 0, 255);">end</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >Q=mu/(2*pi);</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >beta=R(1,2)/sin(mu);</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >alpha=(0.5*(R(1,1)-R(2,2)))/sin(mu);</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >gamma=(1+alpha^2)/beta;</span></span></div></div><div class="inlineWrapper"><div class = 'S7'><span style="white-space: pre"><span style="color: rgb(14, 0, 255);">end</span></span></div></div></div><h3 class = 'S2'><span>plot_betas()</span></h3><div class = 'S0'><span>The function </span><span style=' font-family: monospace;'>plot_betas()</span><span> receives the </span><span style=' font-family: monospace;'>beamline</span><span> description and the initial 4x4 beam matrix </span><span style=' font-family: monospace;'>sigma0 </span><span>as input an produces a plot of the horizontal and the vertical beta function. This function assumes that the emittance of sigma0 is 1, or </span><span texencoding="\det\sigma_0=1" style="vertical-align:-6px"><img src="data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAH0AAAAoCAYAAADTw/VhAAAGCElEQVR4Xu2aZag1VRSGn8/ubsUWC+zADvCHIhZ2YIuJjS2KLRYiit2KKGKC8UOxUBEVsbuwFbuL52ONbIczd86cOR/fuTOz/tx7Zvbs2Xu9K9619kygk9ZpYELrdtxtmA70FhpBB3oHegs10MItd57egd63BuYENgP8e1HfT3UDR0IDVT39AGBnYC1gSuABYNOR2Em7F7EMcDLwHHBBmSqqgj4fsDpwN0xk/h3oZRqetPeXBk4CdgKmAI4Fzil7ZVXQs/neApYcMujTAdsDN5QturvPLMDRwBzACsC6oZNJCvoLwEpDBv2KiB77dqCWamAq4C/gn4i8z45H0PcErgGuAjrQSzH/34BFgPfHG+gSw2uBaTrQq6Edo4cO+pbAwcDywLfA28BxwE0l4X124AhgC2AJ4E/g5QD1+ghLrnk7wN/TxwYMWX/E/7cCew2gBnnKPMCCwNQ9nn8P+GKAeUf1kaGBbjlmqN0F2B8QgF+ADeO6LH7Ggpy+KnA/cAdwHfAdsDZwftT1XpdtagiZfAo4Z53wvg6wH7AVTCQ6RbIB8NioIjjAuoYG+pnh0XqaoTcVa/Sn4kK+ZJOBvw7cFxEifW6PZK5DgYuHBPrMwGVhoGU6+xWYDfitbGCf958J9tzn8J7D3gRWrDHBUEBfGHgD+ApYLOeR2do+Ahbq4eknAqdFOfdObiMq++uoJ71nyVfX0+cGHgRWjon0YJsTHwD3AjskBlpDr4WPZlVMnbl1kmVrTDAU0K+OXHp5hPZe6ykq2TJj+LFgEzPFdUuNWYEf4vcg4d0U9DCwUcwhzzAyZZzAtHQGIC95ooZSR/3RoYCeNV6OB84q2HEv0AVUEH+OfF9FWYOAfnjScrwncrnGlImNC6PVT8AawGtVFjSOxtYG3Zyskmzp6Sl6e7+eLoGz96vMFaG8X91VBV22/2G8xzxteMxq1fSdnweT1+NNPanIMXYEZPO2My8B7ux3wSM0rjboKu/V2NCpwCkVQLdDZwRQVkn+70c/VUHfOgHIaGRUyoulm2WmTP7KYPbZGMcfBXhQYflmefcScCBwWz8LHqExtUG3vv4mNmS5ZaesX09XuZZnisai0YwlGsmLMaAq6B7nWgEoyxWE7sWBjEweBFwa4+ePKHF6bo1Gg0OABYAiTpLfT2PYu5YvK/4M0Ip+74FcEZEzh3rGLku3CjC/9xLraT0xSx8Z6LZi9+7Di24HtgW+D0LY6xEPH4wC8oxFE2P2oOJcIF+vbwI8VJLW8u9pDHu/Edg1dterTveWZYZ58HFg/UQTFwKHJZFinzgYSJW1cTRh9FDzsWJeFZhHAO+XiV7r2b4dwqV6DJZUmuM1wDw3eTKaRdb3qUdL/DTW8XZcXDu8q7+sTpfU2YXbJhSR6Xa35PjTnGle/BL4G5gXeBeYIQZbQ+vNEjyV6tGpYVmjuisBKwPCSGFnznasXbxPCgja7tHt09Ot/1PWLgm1fLPrd3Y0mVK7sKx0nfb6U/E5u4SyfNvO40V0nldisWNVXP/tp+g8XXCs1/UYwTSMPR9eZbi3MePLFBshhsssZ3qA4pm4dXRerKGNBNnY7P4xAZC/jR7OaR0vYdMA8uLcGsqacSpn+1bR2OzO2ZKVlBmuU3G/rl9jNYXlRU4i+EaB8SJHAufFYjV2nXJMGesjCsO3eXG16J59DBj69R67XYbWW6LjlXqaL5Sk+ay9esEzfOvtkifbjnmxBDMieEBjyXhzfBEyVrvUZ04ANgdcm1FBImqLWOCz1JF/Vwa6hzJ5kcQ6b3YAVKa/yXnfL2Y0er9VTHEUG3WtDvNd0YnrHfTLmcm52brv1kD08mlzE6kLw7tRRubfWGkj6Pbn14tUIGfJxIhk2H80ae02Evg2gp7lQE8Ln05QtYQTcNu7jf6su42gW51YYZgT07MFG0pyBO/bN2istBF0wbRbaP1uBWJtbgln+9k+g2Sz0dJW0AXV5o4lod8O2Gvwi56iA6ZGGUGbQW8UkFU204FeRVsNGduB3hAgq2yjA72KthoytgO9IUBW2UYHehVtNWRsB3pDgKyyjX8BMGxTOIhstjAAAAAASUVORK5CYII=" width="62.5" height="20" /></span><span> in both planes, such that </span><span texencoding="\sigma_{11}=\beta_x" style="vertical-align:-6px"><img src="data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAGcAAAAoCAYAAADjT+RAAAAFOklEQVR4Xu2ad6gkRRCHv1NRFAMGTAcqZlQUFOOZEBTBnLP/mDAnFBOoYM6KmHPAjIqioogKJowI5qyY0ePEnOV71MLevNnd7n27zt7ZDQ8eszVV3fXrrvpV9UyijJH1wKSRnVmZGAWcEd4EBZwCzgh7YISnVk5OAWeEPTD4qc0P7AQsD8wGvAzcBfyaa6qcnFyPdZc/EtgLOAV4ElgbeAB4H9gI+C7HXAEnx1vdZS8CtgamAF+2id4G7AocC5ybY66Ak+OtzrJ7ALcA2wL3V8TODmAuBw7KMVfAyfFWveycwLvAN8AaNSJ3ADsD+wNX55gr4OR4q172BOB0wHxjaGsfcwGfAb8BKwLf55gr4OR4a7zs7MBXgAxtaeCjisjRwFnA9kEMsqzlguNOmAwsUGPlD+CVLOszvvA2wH3Ah8AysZyFgGnAlsBlwAH9AKOuFHDmBg4EdukQU1sufgLYZMb3d9YKWvnkGmA/4FZgd+Dv+DsiAPonS2sI9wJH9K+I09JLv9z+1F5CA/59K+DOAejcOwrFHFVzAFMBo8luwO2Ap2bxCG/nB2DXAfvkKG7JdgNHdiH9mwX4BTgzJnAosETEUXdIk6MVViY6B+sQT0HOMEo8Hi8sCnxdeXlz4OF4tjrwao7ybmFtU+CRAMbd4Q59NpSbBK14n4od8Xuu0ZlE3s16HPAGsErNmqx57o3n5p2rctddd3I8pjp/sYib6wAvVhRLGQ8HrgX2zTU6k8i/FDn4kvBFdVknAqfFw35OZi0hOAY4J5SabyQD1eEzmcifwLwR9qoyhpyTgDV7gJEqN0qYLhhFpyG/rivgXFvkwP9XiEI1aw11J+c1YFXgh+Du39ZolL+fF8+l1l+0yawHmGDNWRZfVtB1I1Uua0H/kbDMVQJgzhUoqXP7kCx8DJiL3gJWBrIZWxUciylzjEMW5CTqxk3RfbVl4QTaDQuGBMJj76Q6gZMq183fTbE1qbMMzHC/Vs0E3ZhXxvO+QprvVsHRma+H0qOAC2sMW/dYCUsbzwCMrXXDusd81Qmc1jupcnU2mmJrnoolo/o/vjIxQ907wLLAC8C6ccKyD3UVnIXbKOGeETerSmUoMpUPgNWAnxoEJ3vBA3jBSzSd79gMeKyi01JDkmBbx3xrb62vUZdzPBVLAeaVCypavTxyMjbwNg6AOhlOPRGpcn0tcAgvHQxcClhCmAZ+brOxYfjnR2AL4Pka+5IDayT12O4yPxtdTBUbRCoYu5SrA8c47p3E23GTJzFwrB/PLdY8yr06rKlOT5Ubgp/7UmkvzXAqUZIMCZJ+9AZUdvsmsGMQgm4GLEe83zE9HBJAWtzfDPzVCRyfGydPDiYi29C4d+AWUtWa5/90cvwmwF1t+WAEcQPPCiwXeeVG4J4ArBfyFvqPAjfESfRbg+lGr95aLwPdfk89EalyE5nLoN6V/j8TyiRPnpJ+xzxBwa/vVMgXcPJcazSxwftJ5OW8t6eXltV9GozO+55xo4CT596n4wOOTp2THG3enO4QIXGRAk6O68bLyswsus07MrGHJqDOLrXMzE62X+e02jvmrvdaeod5cmwDrQTYyuh2tZAqNwFfDOTVw4CLgc+jAB1jVBljvugq6Pzt4n+7K9ZBfqFjq8srGjfA2BgGOH48Z6faItYhg3kQuLuykFS5jPUPVdQd7tWACbxa/6UY9lRIJmS7fhHaqo/s7NtJ0WfPtSsaBjgpEy0yCR4o4CQ4qSmRAk5Tnk+wW8BJcFJTIgWcpjyfYLeAk+CkpkQKOE15PsHuv8FZ+ymCV/HXAAAAAElFTkSuQmCC" width="51.5" height="20" /></span><span> and </span><span texencoding="\sigma_{33}=\beta_y" style="vertical-align:-6px"><img src="data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAGgAAAAoCAYAAAASRL/NAAAF0UlEQVR4Xu2bB4glRRCGv1MxIWbEnCOKgjknUBQjYo5gDphRzAHMOZ05oJgDJsyICipiVswYMYsRc5bvqDnmZuftm5439/btOQ0Le7PVVd39d1f9Vd03hrYN9AqMGejRtYOjBWjAN0ELUAvQgK/AgA+vPUEtQAO+As0PbyZgK2BRYArgBeA24Lc6ptoTVGfVOvc5GNgJOAF4HFgJuBd4F1gL+CbVXAtQ6op1lj8P2BRYDfg8J3YTsC1wOHBmqrkWoNQVK5ffAbge2By4uyByeoBzCbBvqrkWoNQVGyo/DfAO8BWwXIm6W4CtgT2BK1LNtQClrthQ+aOAkwHjj24u36YFPgF+BxYHfkg11wKUumITyk8JfAHI3BYEPiioOxQ4DdgiyEKytVSA3BFzATOXWPoTeDF5BKO7w2bAXcD7wEIxlVmB74GNgYuBveqCo74qAE0H7ANs08HHZkv8GLDu6F7v5NFn8eVKYA/gBmB74J/4OShA+jdZc3ToBpC74NI4Nd1syP1P7CbU8N83AW5tQOfOkUymqJoK+BbQq2wH3Ax4euYMV3d2gHY1sFuK4rzscADJOqSGkwG/AqfGIPYH5g2/6k4ZyZa5mF7HYJ7iaUhpeotHo8PswJeFzhsAD8S3ZYGXUpRnsp0AWg94MMBxl7hTn45OBkYz4ydiZ/xRx/Ak0McNewTwOrBUyXzMie6M78ahy+vMuQwgj6wAzBF+dGXguYJy6eSBwFXA7nUMTwJ9no+YfEGsRXFKRwMnxcc6J3Rc1zKADgPOCMXGHwlCsflNhvIXMH24wEzG+tMpMfgfgdsBcwXdZL5VlRtELGeJxFT3X1Y9cMwZYfD3xSKZTZ5LGUCvAEsDLq7c/usSrfL7s+K7tPuz+H3FqDedD3wXzM/jfSFwQE5PVbnkCfWpg4xWUmAMFixpdb5JID4EjE1vAksCtZhcESATLmOOTXbkQMradVG1tbzhIDLj1wYQ+YzZ4DhDgJ3pqirXbb1HisVJq2Vmun43W7FJsC6Lj7Xdm/2LAIn0a6H4EODcEuPmRWbMUkpdmb7WNjkwN/BRoY9kQwA3TJTrBo5/HykW5+mYL6oERxYGqtt7G1gYeBZYJU5alfkMkSkCNFuOLu4YfrTYSeYig3kPWAb4eRjLHnUB8yTK+jq1qnK1JtlwJy/iBMC2PvBIQb9piMTBEtAKUYurPYSyGOTpmB8wzpxT0Gxgd0C6sLUDpE7Gpw43eQ+gS+hVrvYkG+64H3ARYHphSPglp3/NWJ+fgI2AZ0ps7wJ4PWHle/lc/W5vQK8lqONDRBlA+nXvNN6KG0HJgm31+G5C57EerjIruNLweWIix3S4rKoq1/Aa96TO2puuVfIkQRIo19GbVFnvG8CWQRI6GZLVvQqMDVCU0xuZ2K6TO6Eda3H6zeODochCHIB36iZbxZyobBCSAksenjLjlP82mJo75FtVuZ5WtMHOvjHw2trUQk/iJjb2LhJxRvJzR4DWzazvFFxn43bWvLbIYvq4b91qcd2MVPm7VYmHgWNziVtZv6pyVWxOLJlVgadCuYTK01K37RqJ/hIBtOAfFz/jdfYDII1Jx03cvNQarlWVq7sovfbTq1gUlvgYp3tpephPIy0xT5SOG7NkiH0H6ONwdRZfh2tV5XpZmF76PhmPQjpVWFJ1WxSQEUocTOR9vzBB68cJMl94OS60siS4bCJV5VIXoSl5GZsnXFckQ7u/AcXGdJ9jCbivf6TmEw2gGYH7gIciD7D8YVKrYe9EsspuVbkG5t+oCne4JSzdkpvp7wa0S61lcj7J8v5oSGvyBPm6xTgjMzEPkI5LQW+M6nhmvKpcA/NvVIUbzWuFa0ryw7qGXCtf+nhfVHpt0yRAdQf5f+5nIdlapSWh0tYCNHLbw9Nobqj779hagPoL0BrAApE/WRYqltImagzq71RHpzWLqF5TeCFa6YFNe4L6D7RF5Mr/FaUFqP8AJVlsAUparv4L/we+LCs4UyfUiwAAAABJRU5ErkJggg==" width="52" height="20" /></span><span> are the beta functions. It then uses Equation 3.43 to propagate </span><span style="font-family: STIXGeneral, STIXGeneral-webfont, serif; font-style: italic; font-weight: normal; color: rgb(0, 0, 0);">σ</span><span>.</span></div><div class="CodeBlock"><div class="inlineWrapper"><div class = 'S5'><span style="white-space: pre"><span style="color: rgb(14, 0, 255);">function </span><span >plot_betas(beamline,sigma0)</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >[Racc,spos,nmat,nlines]=calcmat(beamline);</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >betax=zeros(1,nmat); betay=betax;</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span style="color: rgb(14, 0, 255);">for </span><span >k=1:nmat</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > sigma=Racc(:,:,k)*sigma0*Racc(:,:,k)';</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > betax(k)=sigma(1,1); betay(k)=sigma(3,3);</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span style="color: rgb(14, 0, 255);">end</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >plot(spos,betax,</span><span style="color: rgb(170, 4, 249);">'k'</span><span >,spos,betay,</span><span style="color: rgb(170, 4, 249);">'r-.'</span><span >); </span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >xlabel(</span><span style="color: rgb(170, 4, 249);">' s[m]'</span><span >); ylabel(</span><span style="color: rgb(170, 4, 249);">'\beta_x,\beta_y [m]'</span><span >)</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >legend(</span><span style="color: rgb(170, 4, 249);">'\beta_x'</span><span >,</span><span style="color: rgb(170, 4, 249);">'\beta_y'</span><span >)</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >axis([0, max(spos), 0, 1.05*max([betax,betay])])</span></span></div></div><div class="inlineWrapper"><div class = 'S7'><span style="white-space: pre"><span style="color: rgb(14, 0, 255);">end</span></span></div></div></div><h3 class = 'S2'><span>plot_sigmas()</span></h3><div class = 'S0'><span>The function </span><span style=' font-family: monospace;'>plot_sigmas()</span><span> receives the </span><span style=' font-family: monospace;'>beamline</span><span> description and the initial 4x4 beam matrix </span><span style=' font-family: monospace;'>sigma0 </span><span>as input an produces a plot of the horizontal and the vertical beam sizes </span><span texencoding="\sigma_x" style="vertical-align:-6px"><img src="data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAB8AAAAoCAYAAAAG0SEsAAACGUlEQVRYR+3WS4iOYRQH8N9IkqRckktJLCTFwgIlt5ISFjYibEgspIhoyqVcSkKSW2xEYSeRspCNjRJSLMhObikkQujUmYlv3vlm3q93ZjbfWT7P+/z/z/9/znnO26IPo6UPuTXJ+8T9pu1N23vVgWbB9ardbWRlbR+EsRhWcNufeFhGRXfIB2MTVmB6HfC7WFAl+RKcSbVd4e7Fvq4++ne/nvINOI1++IZDuILNGIfl+F2GrPbbzsgX4nYSf8RS3M/DA/AC97AOPxq9QBF5FFWAj05lM/GghuA4tuAC1ldJvh2HEzDyHcVWG7F2Cr8wJNNS+g5Fyh9jKr5gAj4UoG7DkVyP1ntdmpkOv1FDETmOuJbtVYR7EWvwDqPwpwryKXiaQFtxrAA0+v4VRuAgWhshjjO1to/E2wRbjcsFwDuz7V5iGr5WRR44oWo8Iq9Ha4Bn4A4+YR7iArUxB3OzUHfjPCZmGj9jftuBooKLnr6O5wiyKLyI2bl+FbvyAvVEP8L3fHL34FZ2xo165LE3C3FgOJ5legLoXEHPd3aBA9iRdbMfofq/6M5gaTSli3Ezu+JSEUhPkk/K1MVLeKK3yWMoLco0xTjuED2lfCOeYC2WYQziAeuP910VXCN5jgcq2izGbAyneKBCcYzhlZic8759DFepfBVO4my2YggYmLa/Qey3q47NKslLu9UkL21ZFQeatlfhYmmMv5vIWykaNpXyAAAAAElFTkSuQmCC" width="15.5" height="20" /></span><span> and </span><span texencoding="\sigma_y" style="vertical-align:-6px"><img src="data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAB8AAAAoCAYAAAAG0SEsAAACNklEQVRYR+3WSajNYRjH8c81haQMyVCSbKRYWBiykJINCzYoShlKGYoUSdhQ5oU5CkVIyo4sxIKFUooQkjInEiKEHr2njnP/Z7r9z7mb8yz/w/N9xt/7tulEa+tEtha8U6rfKnur7E2tQGvgmlruAqzesvfGMPTPiPYn7tSTRS3wPliOuRhfwfk1TMsTPhOHU7bV/G7B1mofFb+vlPkyHEIXfMN2nMVKDMcc/K4HVvptOfh0XE7gD5iFm+nnHniC61iMHx0NIAseQxXOh6TMJuJ2CWAfVuM4luQJX4cdyWH0O4at1OLZQfxC39SWumPIyvwuxuIzRuJ9hte12JWex+q9qptMu2tUP0SPw86n9cryewoL8Q6D8ScP+BjcS47WYG+G09j7ZxiIbdjYEXD8U1r2QXibnC3A6QzH69PaPcU4fM0LHn4iqxGIvu4pcTwBV/EJUxEBFNtoLMJs3CjahG54hP3F1cwauNjpS3iIgMXghU1Jz89hQwqgXNInMC+JUcxF2IF0Jswv/FROZCZhMwbgQWrPdxzN2PmsAAqzsxTH0gejUgIR2D+r5WDpaEuf41aqQPgIWb6A182AH0n6H0PcNSni7uJMGpl5DN1FRAsmp7PiRbPgQ/EyyXOIV5yK/1kjMw/QG3xBHE7tZLrR8Cs4iTNZU9tIeC/sxIpy69JI+KakaB+bAe+OVWmqZyBKfr+SSOSZeZzrcQMKRYyb7uNq6pQnPFg9ETJck+UNrwnaDHmtGshfvbVhKT3gPaoAAAAASUVORK5CYII=" width="15.5" height="20" /></span><span>. It uses Equation 3.43 to propagate </span><span style="font-family: STIXGeneral, STIXGeneral-webfont, serif; font-style: italic; font-weight: normal; color: rgb(0, 0, 0);">σ</span><span>.</span></div><div class="CodeBlock"><div class="inlineWrapper"><div class = 'S5'><span style="white-space: pre"><span style="color: rgb(14, 0, 255);">function </span><span >plot_sigmas(beamline,sigma0)</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >[Racc,spos,nmat,nlines]=calcmat(beamline);</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >sigmax=zeros(nmat,1); sigmay=sigmax; </span><span style="color: rgb(2, 128, 9);">% allocate space</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span style="color: rgb(14, 0, 255);">for </span><span >k=1:nmat</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > sigma=Racc(:,:,k)*sigma0*Racc(:,:,k)';</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > sigmax(k)=sqrt(sigma(1,1)); sigmay(k)=sqrt(sigma(3,3));</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span style="color: rgb(14, 0, 255);">end</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >plot(spos,sigmax,</span><span style="color: rgb(170, 4, 249);">'k'</span><span >,spos,sigmay,</span><span style="color: rgb(170, 4, 249);">'k-.'</span><span >); </span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >xlabel(</span><span style="color: rgb(170, 4, 249);">' s[m]'</span><span >); ylabel(</span><span style="color: rgb(170, 4, 249);">'\sigma_x,\sigma_y [m]'</span><span >)</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >legend(</span><span style="color: rgb(170, 4, 249);">'\sigma_x'</span><span >,</span><span style="color: rgb(170, 4, 249);">'\sigma_y'</span><span >)</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >axis([0, max(spos), 0, 1.05*max([sigmax,sigmay])])</span></span></div></div><div class="inlineWrapper"><div class = 'S7'><span style="white-space: pre"><span style="color: rgb(14, 0, 255);">end</span></span></div></div></div><h3 class = 'S2'><span>periodic_beammatrix()</span></h3><div class = 'S0'><span>The function </span><span style=' font-family: monospace;'>periodic_beammatrix()</span><span> receives the 4x4 transfer matrix </span><span style=' font-family: monospace;'>Rend</span><span> and the emittances </span><span style=' font-family: monospace;'>epsx</span><span> and </span><span style=' font-family: monospace;'>epsy</span><span> as input and returns the 4x4 beam matrix </span><span style="font-family: STIXGeneral, STIXGeneral-webfont, serif; font-style: italic; font-weight: normal; color: rgb(0, 0, 0);">σ</span><span> that obeys </span><span texencoding="\sigma=R_{end}\sigma R_{end}^t" style="vertical-align:-8px"><img src="data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAK4AAAAsCAYAAAD8dIM9AAALNElEQVR4Xu3cA7AkWxIG4P+tbTvWtm3b3o21bb/dfW9t27atWNu2bXvjm80zUVNT3T19b9++3e9WRkzMdPGcPHky//wza/bLKKMG1lAD+63hmMchjxrIaLijEaylBkbDXctlGwc9Gu5oA2upgdFw13LZxkGPhjvawFpqYCca7hWS3DnJBZIcordqf0ry7yTfT/LtJM9O8rYk/13L1Z0+6FMleViSiyY5cu/SPyc5WJIfJvlBklcneWGSv66KHnai4TbdPzHJ7evH95JcqBbpGEluneT+SQ6Z5MNJLrFKi7Zg47lqktfUM/+Z5NJJ3pvk8EmumORxSY6d5Gelo28MvP8pSe6d5A8LHtvEx+1kw7UgdyrNMNKH9rTk933r2POS3GRZi7Lk91wpyevrne9JcvHe+/1+Vx37VpJTV1Ry6FBJHp3kvEnOvsxx72TDfV+SC5eyz53k4z3FW6Cv1rHfJjn6QRQyPDDJg2qe90ryiAED5G15XXLmJJ8vz/zkJCdL8s2KTHdJQldbLjvVcA+T5DdJDluKPmbHizSlH7Wuab+PleSXW74iy39BdwOfNclnB4bAUM9Yx6+Z5FX1b/iYl75Ikvcvc+g71XC74Q++u/qA0mHethiMHPY9qCVph6uNK+T/KonN2Z+jc38sWND1uP59YJI7JDlakn+ssuGa6PFroP1xAvafWebgN/Eu4fAedf8tkzxz4FmO3byOw3F338T7tvrWoyQ5YRLr05dfJ4FNh+QySd5aJ16W5LoDFzn2kjr+qSTn6FzziYpCl9vqCfafvy8e9wiVZV8rydmmDFDIETrWQWyws9RAT5rku71BM+an17EvJDlfElTZKsmJy9uJFieaMjAY9sETzj+2qEGnJZ+S0K4wUgwDG/h9kgsmoQ8CSvHScO0Tlq2YWYZ7+VpAXnaWAPj7z7qozkuEGmbax1v2ugwtc6YN3Czk/yLZ1WD088qI/5KEIeA2b5vknEn+k+S5Se7Zw7obeOXCb0E9YUJg9FkiAf3AhIu+mOT0de6SlXQdpxIujgqepSdOCXX45c5zGo12uiRfmTWIRZ+fZri3SPK0IqIRz8jql9cE7HADt7gbEQmA7HQz8rUkp9nAA66dRFgkf0ty6Fqc9qgPVfLx5iTf2cDzt/IWRQFe8Yb1EvwznImuEs7hdV50X4SB/rQuBPPYQrcgowDjmaBEn3Fx21OL5z1BGf+X9uWli7pmkuEi3N9eRisxUW36SL0UWIeZ7OKbLhuUL2Diz+lwspI0CdhdOzTQfWqTLuBVC3+EkM/Tks8lgVFRVcQm5kHvmERBYJbcoKphrsNZPyDJuZJ8tG58Z5JLTXmId6ELOQGsw1I3+ZDhAvgM87jlUXGcn+xN4PGlIEZws1kaWrHzyrkiBm8Lp/kbHfaT8jhfL5J9xYa9K7+wDtbsR0mE6H6lijGDTzzyi2ZMQAmX8ZIunGjUl2hKTz+e8BzvUirmvODgpcqQ4cqeH1mjkKAof/bFMaHiX0mOtEblUBgWxCDvrlJum9sbK7L4PVSQWObCMAZQ6pRJblcE/1uSXLYGAe68YmBAjsGlytTnnzFgG5Vzgu9t4EZnSbYeU/dOKkgsUxeD7xoy3LbjcHcybpljX4RWFBGRuFHCOggjeNKERenW7GH722zjhPRIvCHJeapihycFCQ5ekA3LMSRvSiKhVsli9JNEQibUE5AQ5GiiQsajw7ugwGm3UQ8TX9033G616JVJZJZD0sKM7BzIn5eY3y5WgTFoHCFq65/uTA52l6wwEmVL3ujv27hoEliw5sr153U1FhvKxhoSnVz4XJUtnneS6I5rSdzdOh62Xd82gN8Ylj5U3Ea1/P/VfcOFm1p2KGRoROkLTg/viVaS0bZGlHkmsx2sAg+CjAdtJJxwbZ8VUXtHhxH8aOuammdui7iWA1Fetga62IzJ2Miksix40xIr1FZrjBkaD6agednWe9C9ztxbWdd7WxfdIua2kGf0DVfJD7dJrt+pmHRfBvegxtAlEgG9m+sgwiuqi/AozfN2x45wVw2ads0y5spbwqsafSSL7bd3o5+GEiYhHwuAsrzOlEGKLDautkVFBZukHzFRhKKPczb78baJPcL62IRgD/ZiN64fwri86UmKIupzgugSO9mEZaKMd13k4VVMMF7Uz0MmDByZjlrSUK5IMolc58FhfZ5bFIIN4WShWvKqMse7t5Y/FBWju95AQ4qkUaFDrqAByDPpV1GE6MBiwDAudqFfWgcfPF92D+NOa/jWb6s5nrj+YhP0IPluifmksvhWrz37VGw5oD/vIcPF2cKCsm+GKkkjslTHWb2HMd51EYag6uNvgjlpyWV/DgyIkRNMg37VvjBAnu0dtQHOUKXQUxSVyNBQirJz2FT3FC+PQurDK3V+4ZhnkVTROc4c1agI1ORRSeBRNNeNyktaP0ZrYW1Gz5lWFGL4xtyMFVPB0IfEOD5WJ3jfkxcDsew1B1fNFzTdPbdJBQjZrBq3HlSZpevwnRpPVg6oz9Ak0h7XLNlqwiPxNubYTdCcFxYlORaZCMGaclqnGI8Io4NIrTTMO9nYrYytpq9Aw7uCUzaDvxluF4Lp7RDBLMyL633NI16jPpnpTs91t0qi0sXIj1jEP6ZkEt/a7gcfNBb1K5Y2CW/90gE9Wntwhch99CT4nGmZIrqoEIpmu2VWr8IyB7gu7+KpJU2M1mLydoxW9t9KqDaLhOr5Se5XoRs3a9FtINSWQk+rPnW9nq8yeGrRARbdyQJ+0SnI8ozRcDduChIb1SrcNkpKMssjN8zYniy759ng1Fbn59Fwos3jgVtgg6SxldPdr3ggOVrqpzAbV8muO3lD8INeQB4FEngdtIHJ4X1RSMGKd6cXc+/nUHICzVocqoKIPzb+Xh18o8edb8Ua5TStVVAZVCauSYcXJvTMy/LAYAPBXsCN3U+CGquhctmum2+Ey70abFL+Vxq+ShkaGg/zJLrYgKKSfAj9BoczVvAIe9UtpMgP5FB6MdCQoBo9qhXsVUwZDXe+hW60lNDV+nU9gafALPAyjFVCJqmT3BHYV0VS85JmFIuBvbF4ra3Q/RJfZV1eSzIJR69yEmwjXq1YGNU2+tEKKrLIExiuYo7EUqFK4UOJGQRinK1RvyXP9NV6X7Rsmr/8QcVzDxkNdz7DRU/pA26dWY0W5Dlw25ImoVASZTEaLdXoKrQXo5T9S/AYrTBqIYVEeA4+lhyDFN3NMd9It/5quBxT0ipvN06isw6W/2C9vn2ThiERReiuMS7dfgsMjf/ngsG3jWr+dMSTt6+Qd89qNNz5F5iHYaQ8CI4XPdVtI8S3Srq637FhCnjT15ZHwRH7po13gpN9cCiZk+Qx1hdUArjRfuf5ZzX/HZrLzcFGE3HQp8bdLZPTE8iD725V2Nbnbd6ggE0K46vAwr1NJL3uRYP9bvS48y/QeMfeGoA/eUaJ56QChrt8o6YJS3LVoo9SMkyMJyY8MV69W3oGE/wvOgpc7boRKoyWuGkNtNYAUUTIbyKCS7JAAs1KEjXcL/qQON8YGZVLFVpN7KITaNU6EVvDPA8sEXbdHh98jlBh02u4Yx8AGuBZYVPMAfzug0ssA6yugKLDEA3WMC8MK+l0jjGCWIwTd62wAnYpkkj0nlWVQ0krDldSt8cO2bGaHye+KQ0oc8PjQryQrjEfldWarhgbhqV9SeNlvCr8j1nB/aIIVSrBB0mqL4gle9gZMESVVtKHhRmhwqaWa7x5JTQwQoWVWIZxEPNqYDTceTU2Xr8SGhgNdyWWYRzEvBoYDXdejY3Xr4QGRsNdiWUYBzGvBkbDnVdj4/UroYH/Af4jR0sUHR2gAAAAAElFTkSuQmCC" width="87" height="22" /></span><span>. In other words, it is periodic.</span></div><div class="CodeBlock"><div class="inlineWrapper"><div class = 'S5'><span style="white-space: pre"><span style="color: rgb(14, 0, 255);">function </span><span >sigma=periodic_beammatrix(Rend,epsx,epsy)</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >[Qx,alphax,betax,gammax]=R2beta(Rend(1:2,1:2)); </span><span style="color: rgb(2, 128, 9);">% eq. 3.60</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >[Qy,alphay,betay,gammay]=R2beta(Rend(3:4,3:4));</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >sigma=zeros(4,4);</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >sigma(1:2,1:2)=epsx*[betax,-alphax;-alphax,gammax]; </span><span style="color: rgb(2, 128, 9);">% eq. 3.78</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >sigma(3:4,3:4)=epsy*[betay,-alphay;-alphay,gammay];</span></span></div></div><div class="inlineWrapper"><div class = 'S7'><span style="white-space: pre"><span style="color: rgb(14, 0, 255);">end</span></span></div></div></div><h3 class = 'S2'><span>tunes()</span></h3><div class = 'S0'><span>The function tunes() receives the 4x4 transfer matrix </span><span style=' font-family: monospace;'>Rend</span><span> and returns the horizontal and vertical tunes, </span><span texencoding="Q_x" style="vertical-align:-6px"><img src="data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAACYAAAAoCAYAAACSN4jeAAADW0lEQVRYR+3XaeilYxgG8N/YBlETg+z7EllSlpItoiGMsiUiPtjX8okiSii7EGFoFMIHCUWRLUX2JXtkTckXZY0u3Wc6c7znvO97JjXlPPWvf+e97vu+nnt/5lhOz5zllJcZsb6RmXnsf++x/bAQ+2MjrI6v8SxuxRt9PTSK75tjB+Aq7IhFuAcfYg3si+uwTmEuxl/TEuxKbBVcjfPwLY7Aaw1GN8EH5cFLcfl/SWwlPIQj8RN2wycTDF6LC/EntsWn05Br89iKuB/HlvJD8GSLoaPrIoFdiYS092kjlnBcVlofxHEdLOyBVwr3OA7rIPMvyCRi21d1Jb9+wXb4ooORBXiicMFv1kGmF7HnsXdJpPpO6WgguLsK+yPW7ii3FGycx3bCW0PI3fFqRwM3VPUGngqN53ufccSux/ml7X3s0EPzm9i58ItxYg/ZJdBxxD4fyo2QTPl3OWmu37NknToJ93URHMU0EVu5kn2FAh86lMxtNi6o7h9cRtQW+K1NqOl7E7EtRxropviyo/L3hnIqJJNvU50mYlvh49KWWTcXv3fQfjCeKtwP2Bw/d5BrhLSF8les1mEYJ+xJ+gz3XOIgPDctqciNS/6PsHUp3qAG9yQ7Z+PmApyKu5eF1CRi5+CmUn4UHplgaB88gxTN1LOxS1UGsyrSv5InWfqyUWRbGD0Z6hnyCXcuc2eLp45BdrrDkdGV8EdH5PI3mMsTn28b1iaRvElSn16zMuHPZDgXJ+PpIjUomLYozq/UyHp0b5GME9JeXh8It20Xa1bnTqPMWEq1ZRX6ptboO/BuG5OG7y9VVHLhxrWojVh0hkjW6eRamm3mX5eVeb2aAk28r8FF2BhfNQHaiK2LB7Anjq+wdelNWRYPxGljvJnfb6+Z+nZfYln4Hq5X0EA2nsoUiNeyfbyAd/BdPUh2LTJ5NeVd8EeD0bVwI07AmbitD7HcKC+eDOQkZRJ2Gwzm57i0SnNNUl8ypoojl6XgCrxcyZ5IpGcuVTxNoZyHM6pJhtjg5KbpWXmmZYHMTA020+GzGvR5U+b/0ZP1Oi1or4rCi+WptI+sResPLZf/yLblWFvBpal2maO3VPEkdI+W0l2QN8FjOGu0oJaVWBvxqb/PiPV13cxjM4/19UBf/N9iMpwp8XsoUgAAAABJRU5ErkJggg==" width="19" height="20" /></span><span> and </span><span texencoding="Q_y" style="vertical-align:-6px"><img src="data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAACYAAAAoCAYAAACSN4jeAAADjklEQVRYR+3XWcjtYxQG8N8xj6VjKmSeMlOGG/MUkuGYJWMhsxuJIsqQzJIhQ6IQScZChAtFIVPmKFNKLigzPVr7a9tn7/1/93fy+cpZV7v2s973+a93rfWsNccstTmzlJeFxCZ9mYUR+99HbBcciF2xBpbBF3geN+P1SSM0iJ80x3bHFdgcd+MuvI/lsDOuwcqFuQB/TpdgK7ElcCXOwlc4AK8NuXRNvFcRvAiX/JvEFsODOAjfY1t8NObCq3EufsdG+Hg65Loitijuw+F1+L54quOiQ+tDArsMedKJrYtYnuPiOvUBHNFww/Z4pXCPY/8Gn/kg44htUtWV/PoJG+Ozhkv2wZOFC37tBp+JiL2IHcsj1XdC4wXB3VHY77Bio98/YKMitgXe7ENuh1cbL7iuqjfwVGgiP7GNInYtzq7T3sWmE5z8BrYs/L04ZgLfKegoYp/25UZIpvxbLM31G6bGqWNxT4vjIGYYscUr2Rcp8H59ydx1xznV/YOLRK2LX7qchv0/jNh6Aw10LXzeePg7fTkVksm3adkwYuvjwzotWrckfm04fW88XbhvsQ5+bPAbCul6yp+xdIMY59mT9BH3fMReeGG6pOI3Kvk/wAZ18Gol3OPuOR03FuBE3LkgpMYROwM31OGH4OExF+2EZ5GimbY2tlRlMEsh/St5kqEvE0WmhUGLqEfk89z5mNvHfECGyqjCbji/r41sjWhqBoWXe/7jtHL1miSSN0nqU0or4xNlOBPH4Zki1SuYca8Y3X2rimKbAq6AJ4rorS3Eglm+OncaZWQp1ZZR6Msao2/D2xPm02m4qXK4N9edhOeQxv63dY09wYRIxunkWppt9K9lZF61VGCQd5puhsdTcUv9eXk97xS2i9gquB874Kh6tpbelGFxD5w8IprZExLpeUjf3AoP9WPHEcvAF3C2oJ4lUlGBRC3Tx0uVM1/XQpK8CZlsTdkLfhtB7HociXx4htFL8UcLsRyejSeCHM1bCRuip5+j0irNNTP/hSOquOeX6KeaDytdfqylXaRK8v5pkiHWs7lIz8qalgEymhps1OGTEvrslPndZZmGE/VHaz+dD9+VY10XpKm26OjgOYl8KjzpMnTjWlBiXcRH/b8Z9kRmvaH2XxBbFudV0o9sOzNFLBv6wZWHR+Mq/DAu3DNFLM05/fARHN9FKoRniljuymCQ/bTJZpJYE6EeaNYS+wsQRqQp/+7riwAAAABJRU5ErkJggg==" width="19" height="20" /></span><span>, respectively.</span></div><div class="CodeBlock"><div class="inlineWrapper"><div class = 'S5'><span style="white-space: pre"><span style="color: rgb(14, 0, 255);">function </span><span >Q=tunes(Rend);</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >[Qx,alphax,betax,gammax]=R2beta(Rend(1:2,1:2));</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >[Qy,alphay,betay,gammay]=R2beta(Rend(3:4,3:4));</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >Q=[Qx,Qy];</span></span></div></div><div class="inlineWrapper"><div class = 'S7'><span style="white-space: pre"><span style="color: rgb(14, 0, 255);">end</span></span></div></div></div><h3 class = 'S2'><span>drawmag()</span></h3><div class = 'S0'><span>The function</span><span style=' font-family: monospace;'> drawmag()</span><span> receives the beamline description and the vertical position </span><span style=' font-family: monospace;'>vpos</span><span> and </span><span style=' font-family: monospace;'>height</span><span> of the magnets on the plot as input and produces a graphical rendition of the quadrupoles and dipoles on a plot.</span></div><div class="CodeBlock"><div class="inlineWrapper"><div class = 'S5'><span style="white-space: pre"><span style="color: rgb(14, 0, 255);">function </span><span >drawmag(beamline,vpos,height)</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >hold </span><span style="color: rgb(170, 4, 249);">on</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span style="color: rgb(2, 128, 9);">% legend('AutoUpdate','off') % avoids an extra entry for the magnet drawings</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >nlines=size(beamline,1);</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >nmat=sum(beamline(:,2))+1;</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >spos=zeros(nmat,1);</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >ic=1;</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span style="color: rgb(14, 0, 255);">for </span><span >line=1:nlines </span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > </span><span style="color: rgb(14, 0, 255);">for </span><span >seg=1:beamline(line,2)</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > ic=ic+1;</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > </span><span style="color: rgb(14, 0, 255);">switch </span><span >beamline(line,1)</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > </span><span style="color: rgb(14, 0, 255);">case </span><span >2</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > dv=0.15*height*sign(beamline(line,4));</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > rectangle(</span><span style="color: rgb(170, 4, 249);">'Position'</span><span >,[spos(ic-1),vpos+dv,0.1,height])</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > </span><span style="color: rgb(14, 0, 255);">case </span><span >4</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > L=beamline(line,3);</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > rectangle(</span><span style="color: rgb(170, 4, 249);">'Position'</span><span >, </span><span style="color: rgb(14, 0, 255);">...</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > [spos(ic-1),vpos+0.25*height,L,0.5*height])</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > </span><span style="color: rgb(14, 0, 255);">case </span><span >5</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > L=beamline(line,3);</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > dv=0.15*height*sign(beamline(line,4));</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > rectangle(</span><span style="color: rgb(170, 4, 249);">'Position'</span><span >,[spos(ic-1),vpos+dv,L,height])</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > </span><span style="color: rgb(14, 0, 255);">end</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > spos(ic)=spos(ic-1)+beamline(line,3);</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > </span><span style="color: rgb(14, 0, 255);">end</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span style="color: rgb(14, 0, 255);">end</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >plot([spos(1),spos(end)],[vpos+0.5*height,vpos+0.5*height],</span><span style="color: rgb(170, 4, 249);">'k:'</span><span >);</span></span></div></div><div class="inlineWrapper"><div class = 'S7'><span style="white-space: pre"><span style="color: rgb(14, 0, 255);">end</span></span></div></div></div><h3 class = 'S2'><span>edteng()</span></h3><div class = 'S0'><span>The function edteng() receives a 4x4 full-turn matrix R as input and returns the 4x4 matrices </span><span texencoding="{\cal O}" style="vertical-align:-5px"><img src="data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAABsAAAAkCAYAAAB4+EEtAAACS0lEQVRYR+3WWaiNYRQG4OeglJKiKJSQKQpJUooIRRIZEkWZx2RIIi5MGYoMJZI7hRtREiWSSCFxoWTMhQsUpcy06tvaZ/v3+f99kqKzrnbfXmu933rXu9b31/mLVvcXsTSB/RG2/xkau2EWhqMLOuIFLmMTXlfS0ZjKemE/xuAjDuESHqAT5mACBuBNOWCtYOuxGS1xHkvwPKOh69AZyxsD1gyHMT8Fx++l+F5FOc1T9eHzy4pWdjzRE4F7sSplCKqi0qEJ+Co24nHWJYqAzcPRFHwB41PiGTiAbaln7bAd/dAfz2oVSIjhDlql3kSSd0l5dzEqCaOUtzUe4hYm1Qp2EtNSUFQUogjbg0/YkEHXQUxFh1rA+qRbhzhOYXpZcMzSWtzOAIvzXWiD90UFEjcsKa43HpUFxvDGZT5kgO1LIxHj8aMo2BN0xekyKkux9zEInzPAzqEvYsPUs2pqDGFEo8OG4XpF3AnsxL2K86A8hvwi5hYFm5IqepWUV48OjMRYxKYot9EJKJR4pijYsjRDMV8LMqgKRqLiaxX/hXqHoAe+FAXbkjbBigSagffb0WDcxGIcyQqo1rPVaZYm4mwBpBa4gW+p4q+1gAVIcJ7JfUaieHJmYyBCxZlWrbL2eIqt2JFTWezGoHtcRg8LST+cdqe107PKPMXijcEfgcmJxgbv1dDWj3VzBW+xKG2QmKPumImFaduvTD65rc17YmLbr0E8J23T8n2ZvjOOJapzQUoOeWCFExVxbAIrwlKuTxONuRQVcfh/afwJrmRuJVl019kAAAAASUVORK5CYII=" width="13.5" height="18" /></span><span>, </span><span texencoding="{\cal A}" style="vertical-align:-5px"><img src="data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAACAAAAAkCAYAAADo6zjiAAACeklEQVRYR+3WW6jNeRQH8M/hweXBpUTklksZhZkkDy4NE5Nc3uZpJCQp1GiieJBScnlyaUp4mEKuuSfXTIi8oESYFxGRO7lNMS2tU7vdPmfvczjOg72edv/9+6/fd33X+n7Xv0YzR00z368KoMrAd8FAa/yJaeiFcThbq76mZqAb9uEhJuBlgnj9LQAMxnbMxeRkYSmWF3pPUzHQFccwI6u+hvfoiedNDaANTmERzuEApmAlFhc7b1MwsBXXsQK/4CTeojcefQmAflnBcPTJZGewFpcz8VRMx3g+23w8H5Rn/ii1dyplYF5SuAE7cBOdc7IXYn1SfQKj8AC/I9j4gL6411gAUc3BYv0WJOuA4+iCOTiKFojB+wGbMLuurVsJA9HP/VhSz+pehlnogU/4Dbvyd4AIxkpGOQAd8RQjcb6OHD9iM27jr5z8KwgfOIJJ9X1zlAMQfQ4Xiwm+UyJRqwQ2MzU+AheTsTgethsqqDPKAYgXX2UV/5TIsgaPsSqHch1CLUNzBkIB9UYlAPbgLhYUZfo1ZTkWH1N20a4YyoiYiS1fA8BPOI1h+DcThqzC4YLikFxtvEC7ZCUG8l1jAASF83EJ2zJBeHpYa0jyGaId4Q0XCi4Ync/jUSycWDxlo1QLbqF/vjmgQEJxefQ8qgxgG4uyH8bENJ4Y2kJmKh7CmOo3ucGepI2Gj0fEQomq/sbAZCKs+Eb2/WrOwSHsTJcMpsKIKgYQB8NOg4Xw7v9SgmG33dPnowURbREMBVvBRvt8vju/ePbifrkelGpBJ6zGGLRM/ccFQXu4XHEMQRhPRFT8c7lLC/+vRIYNydfgs1UAVQaqDFQZ+B8hlIAlu5whRgAAAABJRU5ErkJggg==" width="16" height="18" /></span><span>, and </span><span style="font-family: STIXGeneral, STIXGeneral-webfont, serif; font-style: italic; font-weight: normal; color: rgb(0, 0, 0);">T</span><span> from Equation 3.104. The fourth output </span><span style=' font-family: monospace;'>para</span><span>, defined in the last line of the function, contains the eigentunes and beta functions as well as other parameters related to coupled beam lines. The algorithm is a straight implementation following reference [13] in the book.</span></div><div class="CodeBlock"><div class="inlineWrapper"><div class = 'S5'><span style="white-space: pre"><span style="color: rgb(14, 0, 255);">function </span><span >[O,A,T,para]=edteng(R)</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >TRMN=0.5*(R(1,1)-R(3,3)+R(2,2)-R(4,4));</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >DETM=R(3,1)*R(4,2)-R(3,2)*R(4,1);</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >TR=R(1,3)*R(3,1)+R(1,4)*R(4,1)+R(2,3)*R(3,2)+R(2,4)*R(4,2);</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >CCMU=sqrt(TRMN*TRMN+2*DETM+TR)*sign(TRMN);</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span style="color: rgb(14, 0, 255);">if </span><span >(abs(DETM) < 1E-10 & abs(TR)<1E-10) CCMU=TRMN; </span><span style="color: rgb(14, 0, 255);">end</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >QQ=TRMN/CCMU;</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span style="color: rgb(14, 0, 255);">if </span><span >(abs(QQ)>1) QQ=0; </span><span style="color: rgb(14, 0, 255);">end</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >phi=0.5*acos(QQ);</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >DENOM=CCMU*sin(2D0*phi);</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span style="color: rgb(14, 0, 255);">if </span><span >(abs(DENOM)>1E-10) </span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > D11=-(R(3,1)+R(2,4))/DENOM; D12=-(R(3,2)-R(1,4))/DENOM;</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > D21=-(R(4,1)-R(2,3))/DENOM; D22=-(R(4,2)+R(1,3))/DENOM;</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span style="color: rgb(14, 0, 255);">else</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span > D11=0; D12=0; D21=0; D22=0;</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span style="color: rgb(14, 0, 255);">end</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >A11=R(1,1)-(D22*R(3,1)-D12*R(4,1))*tan(phi);</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >A12=R(1,2)-(D22*R(3,2)-D12*R(4,2))*tan(phi);</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >A21=R(2,1)-(D11*R(4,1)-D21*R(3,1))*tan(phi);</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >A22=R(2,2)-(D11*R(4,2)-D21*R(3,2))*tan(phi);</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >B11=R(3,3)+(D11*R(1,3)+D12*R(2,3))*tan(phi);</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >B12=R(3,4)+(D11*R(1,4)+D12*R(2,4))*tan(phi);</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >B21=R(4,3)+(D21*R(1,3)+D22*R(2,3))*tan(phi);</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >B22=R(4,4)+(D21*R(1,4)+D22*R(2,4))*tan(phi);</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >[Q1,alpha1,beta1,gamma1]=R2beta([A11,A12;A21,A22]);</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >[Q2,alpha2,beta2,gamma2]=R2beta([B11,B12;B21,B22]);</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span style="color: rgb(14, 0, 255);">if </span><span >~isreal(Q1) disp(</span><span style="color: rgb(170, 4, 249);">'Mode 1 unstable'</span><span >); </span><span style="color: rgb(14, 0, 255);">end</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span style="color: rgb(14, 0, 255);">if </span><span >~isreal(Q2) disp(</span><span style="color: rgb(170, 4, 249);">'Mode 2 unstable'</span><span >); </span><span style="color: rgb(14, 0, 255);">end</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >A=zeros(4);</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >A(1,1)=1/sqrt(beta1); A(2,1)=alpha1/sqrt(beta1); A(2,2)=sqrt(beta1);</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >A(3,3)=1/sqrt(beta2); A(4,3)=alpha2/sqrt(beta2); A(4,4)=sqrt(beta2);</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >O=eye(4);</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >O(1,1)=cos(2*pi*Q1); O(1,2)=sin(2*pi*Q1); </span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >O(2,1)=-O(1,2); O(2,2)=O(1,1);</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >O(3,3)=cos(2*pi*Q2); O(3,4)=sin(2*pi*Q2);</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >O(4,3)=-O(3,4); O(4,4)=O(3,3);</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >T=eye(4)*cos(phi);</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >T(3:4,1:2)=[D11,D12;D21,D22]*sin(phi);</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >T(1:2,3:4)=[-D22,D12;D21,-D11]*sin(phi);</span></span></div></div><div class="inlineWrapper"><div class = 'S6'><span style="white-space: pre"><span >para=[Q1,alpha1,beta1,Q2,alpha2,beta2,phi,D11,D12,D21,D22];</span></span></div></div><div class="inlineWrapper"><div class = 'S7'><span style="white-space: pre"><span style="color: rgb(14, 0, 255);">end</span></span></div></div></div><div class = 'S9'></div>
<br>
<!--
##### SOURCE BEGIN #####
function [Racc,spos,nmat,nlines]=calcmat(beamline)
% CALCMAT Companion software for "Volker Ziemann, _Hands-on Accelerator physics using
% MATLAB, CRCPress, 2019_" (https://www.crcpress.com/9781138589940)
%% Beam optics support functions 4D (Section 3.5)
% Volker Ziemann, 211119
%
% In this live script we define the functions for the 4D beam optics calculations,
% such as |calcmat()| that a frequently used in other calculations. All described
% functions reside in the subdirectory |4D| that is contained in the archive |BeamOpticsSupportFile.zip|.
% Any scripts using these function need to include that subirectory with the command
% "|addpath ./4D|".
% The function |calcmat()| to calculate all transfer matrices
% The following function receives the |beamline| description as input and returns
%%
% * Racc(4,4,nmat): transfer matrices from the start to the each of each segment,
% such that R(:,:,end) is the transfer matrix from the start to the end of the
% beamline.
% * spos: position along the beamline after each segment, useful when plotting.
% * nmat: number of segments
% * nlines: number of lines in the beamline
ndim=size(DD(1),1);
nlines=size(beamline,1); % number of lines in beamline
nmat=sum(beamline(:,2))+1; % sum over repeat-count in column 2
Racc=zeros(ndim,ndim,nmat); % matrices from start to element-end
Racc(:,:,1)=eye(ndim); % initialize first with unit matrix
spos=zeros(nmat,1); % longitudinal position
ic=1; % element counter
for line=1:nlines % loop over input elements
for seg=1:beamline(line,2) % loop over repeat-count
ic=ic+1; % next element
Rcurr=eye(4); % matrix in next element
switch beamline(line,1)
case 1 % drift
Rcurr=DD(beamline(line,3));
case 2 % thin quadrupole
Rcurr=Q(beamline(line,4));
case 4 % sector dipole
phi=beamline(line,4)*pi/180; % convert to radians
rho=beamline(line,3)/phi;
Rcurr=SB(beamline(line,3),rho);
case 5 % thick quadrupole
Rcurr=QQ(beamline(line,3),beamline(line,4));
case 20 % coordinate roll
Rcurr=ROLL(beamline(line,4));
otherwise
disp('unsupported code')
end
Racc(:,:,ic)=Rcurr*Racc(:,:,ic-1); % concatenate
spos(ic)=spos(ic-1)+beamline(line,3); % position of element
end
end
end
% Transfer matrix for a drift space |DD(L)|
% The function |DD()| receives the length |L| of a drift space and resturns
% the 4x4 transfer matrix |out| for a drift space.
function out=DD(L)
out=eye(4);
out(1,2)=L;
out(3,4)=L;
end
% Transfer matrix for a thin-lens quadrupole |Q(F)|
% The function |Q()| receives the focal length |F| as input and returns the
% 4x4 transfer matrix |out| for a thin-lens quadrupole.
function out=Q(F)
out=eye(4);
if abs(F)<1e-8, return; end
out(2,1)=-1/F;
out(4,3)=1/F;
end
% Transfer matrix for a thick quadrupole |Q(F)|
% The function |QQ()| receives the length |L| and |k1| as input and returns
% the 4x4 transfer matrix |out| for a thick quadrupole.
function out=QQ(L,k)
ksq=sqrt(abs(k));
out=eye(4);
if abs(k) < 1e-6
out(1,2)=L;
out(3,4)=L;
else
A=[cos(ksq*L),sin(ksq*L)/ksq;-ksq*sin(ksq*L),cos(ksq*L)];
B=[cosh(ksq*L),sinh(ksq*L)/ksq;ksq*sinh(ksq*L),cosh(ksq*L)];
if k>0
out(1:2,1:2)=A;
out(3:4,3:4)=B;
else
out(1:2,1:2)=B;
out(3:4,3:4)=A;
end
end
end
% Transfer matrix for a sector dipole |SB(L,rho)|
% The function |SB()| receives the length |L| and bending radius |rho| of a
% horizontally deflecting sector dipole magnet and returns its 4x4 transfer matrix
% |out|.
function out=SB(L,rho)
phi=L/rho;
out=eye(4);
out(3,4)=L;
if abs(phi)<1e-8
out(1,2)=L;
else
out(1:2,1:2)=[cos(phi),rho*sin(phi); ...
-sin(phi)/rho,cos(phi)];
end
end
% Transfer matrix for coordinate rotation |ROLL(phi)|
% The function |ROLL()| receives the roll angle |phi| (in degree) around the
% s-direction as input and returns the corresponding 4x4 transfer matrix |out|.
function out=ROLL(phi) % phi in degree
c=cos(phi*pi/180); s=sin(phi*pi/180);
out=zeros(4);
out(1,1)=c; out(1,3)=s; out(2,2)=c; out(2,4)=s;
out(3,1)=-s; out(3,3)=c; out(4,2)=-s; out(4,4)=c;
end
% R2beta()
% The function |R2beta()| receives a transfer matrix |R| as input and returns
% the "tune" $Q=\mu/2\pi$ for the transfer matrix R, as well as the periodic Twiss
% parameters $\alpha$, $\beta$, and $\gamma$ following Equation 3.60.
function [Q,alpha,beta,gamma]=R2beta(R)
mu=acos(0.5*(R(1,1)+R(2,2)));
if (R(1,2)<0), mu=2*pi-mu; end
Q=mu/(2*pi);
beta=R(1,2)/sin(mu);
alpha=(0.5*(R(1,1)-R(2,2)))/sin(mu);
gamma=(1+alpha^2)/beta;
end
% plot_betas()
% The function |plot_betas()| receives the |beamline| description and the initial
% 4x4 beam matrix |sigma0| as input an produces a plot of the horizontal and the
% vertical beta function. This function assumes that the emittance of sigma0 is
% 1, or $\det\sigma_0=1$ in both planes, such that $\sigma_{11}=\beta_x$ and $\sigma_{33}=\beta_y$
% are the beta functions. It then uses Equation 3.43 to propagate $\sigma$.
function plot_betas(beamline,sigma0)
[Racc,spos,nmat,nlines]=calcmat(beamline);
betax=zeros(1,nmat); betay=betax;
for k=1:nmat
sigma=Racc(:,:,k)*sigma0*Racc(:,:,k)';
betax(k)=sigma(1,1); betay(k)=sigma(3,3);
end
plot(spos,betax,'k',spos,betay,'r-.');
xlabel(' s[m]'); ylabel('\beta_x,\beta_y [m]')
legend('\beta_x','\beta_y')
axis([0, max(spos), 0, 1.05*max([betax,betay])])
end
% plot_sigmas()
% The function |plot_sigmas()| receives the |beamline| description and the initial
% 4x4 beam matrix |sigma0| as input an produces a plot of the horizontal and the
% vertical beam sizes $\sigma_x$ and $\sigma_y$. It uses Equation 3.43 to propagate
% $\sigma$.
function plot_sigmas(beamline,sigma0)
[Racc,spos,nmat,nlines]=calcmat(beamline);
sigmax=zeros(nmat,1); sigmay=sigmax; % allocate space
for k=1:nmat
sigma=Racc(:,:,k)*sigma0*Racc(:,:,k)';
sigmax(k)=sqrt(sigma(1,1)); sigmay(k)=sqrt(sigma(3,3));
end
plot(spos,sigmax,'k',spos,sigmay,'k-.');
xlabel(' s[m]'); ylabel('\sigma_x,\sigma_y [m]')
legend('\sigma_x','\sigma_y')
axis([0, max(spos), 0, 1.05*max([sigmax,sigmay])])
end
% periodic_beammatrix()
% The function |periodic_beammatrix()| receives the 4x4 transfer matrix |Rend|
% and the emittances |epsx| and |epsy| as input and returns the 4x4 beam matrix
% $\sigma$ that obeys $\sigma=R_{end}\sigma R_{end}^t$. In other words, it is
% periodic.
function sigma=periodic_beammatrix(Rend,epsx,epsy)
[Qx,alphax,betax,gammax]=R2beta(Rend(1:2,1:2)); % eq. 3.60
[Qy,alphay,betay,gammay]=R2beta(Rend(3:4,3:4));
sigma=zeros(4,4);
sigma(1:2,1:2)=epsx*[betax,-alphax;-alphax,gammax]; % eq. 3.78
sigma(3:4,3:4)=epsy*[betay,-alphay;-alphay,gammay];
end
% tunes()
% The function tunes() receives the 4x4 transfer matrix |Rend| and returns the
% horizontal and vertical tunes, $Q_x$ and $Q_y$, respectively.
function Q=tunes(Rend);
[Qx,alphax,betax,gammax]=R2beta(Rend(1:2,1:2));
[Qy,alphay,betay,gammay]=R2beta(Rend(3:4,3:4));
Q=[Qx,Qy];
end
% drawmag()
% The function |drawmag()| receives the beamline description and the vertical
% position |vpos| and |height| of the magnets on the plot as input and produces
% a graphical rendition of the quadrupoles and dipoles on a plot.
function drawmag(beamline,vpos,height)
hold on
% legend('AutoUpdate','off') % avoids an extra entry for the magnet drawings
nlines=size(beamline,1);
nmat=sum(beamline(:,2))+1;
spos=zeros(nmat,1);
ic=1;
for line=1:nlines
for seg=1:beamline(line,2)
ic=ic+1;
switch beamline(line,1)
case 2
dv=0.15*height*sign(beamline(line,4));
rectangle('Position',[spos(ic-1),vpos+dv,0.1,height])
case 4
L=beamline(line,3);
rectangle('Position', ...
[spos(ic-1),vpos+0.25*height,L,0.5*height])
case 5
L=beamline(line,3);
dv=0.15*height*sign(beamline(line,4));
rectangle('Position',[spos(ic-1),vpos+dv,L,height])
end
spos(ic)=spos(ic-1)+beamline(line,3);
end
end
plot([spos(1),spos(end)],[vpos+0.5*height,vpos+0.5*height],'k:');
end
% edteng()
% The function edteng() receives a 4x4 full-turn matrix R as input and returns
% the 4x4 matrices ${\cal O}$, ${\cal A}$, and $T$ from Equation 3.104. The fourth
% output |para|, defined in the last line of the function, contains the eigentunes
% and beta functions as well as other parameters related to coupled beam lines.
% The algorithm is a straight implementation following reference [13] in the book.
function [O,A,T,para]=edteng(R)
TRMN=0.5*(R(1,1)-R(3,3)+R(2,2)-R(4,4));
DETM=R(3,1)*R(4,2)-R(3,2)*R(4,1);
TR=R(1,3)*R(3,1)+R(1,4)*R(4,1)+R(2,3)*R(3,2)+R(2,4)*R(4,2);
CCMU=sqrt(TRMN*TRMN+2*DETM+TR)*sign(TRMN);
if (abs(DETM) < 1E-10 & abs(TR)<1E-10) CCMU=TRMN; end
QQ=TRMN/CCMU;
if (abs(QQ)>1) QQ=0; end
phi=0.5*acos(QQ);
DENOM=CCMU*sin(2D0*phi);
if (abs(DENOM)>1E-10)
D11=-(R(3,1)+R(2,4))/DENOM; D12=-(R(3,2)-R(1,4))/DENOM;
D21=-(R(4,1)-R(2,3))/DENOM; D22=-(R(4,2)+R(1,3))/DENOM;
else
D11=0; D12=0; D21=0; D22=0;
end
A11=R(1,1)-(D22*R(3,1)-D12*R(4,1))*tan(phi);
A12=R(1,2)-(D22*R(3,2)-D12*R(4,2))*tan(phi);
A21=R(2,1)-(D11*R(4,1)-D21*R(3,1))*tan(phi);
A22=R(2,2)-(D11*R(4,2)-D21*R(3,2))*tan(phi);
B11=R(3,3)+(D11*R(1,3)+D12*R(2,3))*tan(phi);
B12=R(3,4)+(D11*R(1,4)+D12*R(2,4))*tan(phi);
B21=R(4,3)+(D21*R(1,3)+D22*R(2,3))*tan(phi);
B22=R(4,4)+(D21*R(1,4)+D22*R(2,4))*tan(phi);
[Q1,alpha1,beta1,gamma1]=R2beta([A11,A12;A21,A22]);
[Q2,alpha2,beta2,gamma2]=R2beta([B11,B12;B21,B22]);
if ~isreal(Q1) disp('Mode 1 unstable'); end
if ~isreal(Q2) disp('Mode 2 unstable'); end
A=zeros(4);
A(1,1)=1/sqrt(beta1); A(2,1)=alpha1/sqrt(beta1); A(2,2)=sqrt(beta1);
A(3,3)=1/sqrt(beta2); A(4,3)=alpha2/sqrt(beta2); A(4,4)=sqrt(beta2);
O=eye(4);
O(1,1)=cos(2*pi*Q1); O(1,2)=sin(2*pi*Q1);
O(2,1)=-O(1,2); O(2,2)=O(1,1);
O(3,3)=cos(2*pi*Q2); O(3,4)=sin(2*pi*Q2);
O(4,3)=-O(3,4); O(4,4)=O(3,3);
T=eye(4)*cos(phi);
T(3:4,1:2)=[D11,D12;D21,D22]*sin(phi);
T(1:2,3:4)=[-D22,D12;D21,-D11]*sin(phi);
para=[Q1,alpha1,beta1,Q2,alpha2,beta2,phi,D11,D12,D21,D22];
end
%%
%
##### SOURCE END #####
-->
</div></body></html>