Alignment Check

Look at misalignment of quadrupole magnets in a long beam line, and its effect on final trajectory at the end of the line.

Nqds = 99
L    = 10
s    = c(0,1:Nqds)*L
F    = 13.8379
q = (-1)^((0:Nqds)+1)/F

Let’s check the beam line:

plot(s,q,typ="h",ylim=c(-1,1)*0.25)
abline(h=0)

# assign random misalignments to the quadrupoles...
d = rnorm(Nqds+1,0,0.5)

plot(s,d,typ="h",ylim=c(-1,1)*2)
abline(h=0)

Write a function to track a particle through the system, given a set of misalignments:

trk = function(V,d){
   i = 0
   while(i < Nqds){
      i = i+1
      xp = xp-x*q[i]  + d[i]*q[i]
      x  = x+xp*L
      V  = rbind(V,c(x,xp))
   }
   V
}

Test it out with one vector:

x  = 0
xp = 0
V = c(x,xp)  # values at s=0
V = trk(V,d)
xt  = V[,1]  # position through the line
plot(s, xt,typ="l",ylim=c(-1,1)*20,col="red",lwd=4)
abline(h=0,col="blue")

x  = 0       # again, with more lines this time
xp = 0
V = c(x,xp)  # values at s=0
V = trk(V,d)
xt  = V[,1]  # position through the line
plot(s, xt,typ="l",ylim=c(-1,1)*20,col="red",lwd=4)
abline(h=0,col="blue")

# draw a nearby trajectory...
x  =  rnorm(1,0, 2)
xp = (rnorm(1,0, 2) - 0.685*x)/100
V = c(x,xp)  # values at s=0
V = trk(V,d)
xt  = V[,1]  # position through the line
points(s, xt,typ="l",lty=2)

# draw another nearby trajectory
x  =  rnorm(1,0, 2)
xp = (rnorm(1,0, 2) - 0.685*x)/100
V = c(x,xp)  # values at s=0
V = trk(V,d)
xt  = V[,1]  # position through the line
points(s, xt,typ="l",lty=2)

# draw another nearby trajectory
x  =  rnorm(1,0, 2)
xp = (rnorm(1,0, 2) - 0.685*x)/100
V = c(x,xp)  # values at s=0
V = trk(V,d)
xt  = V[,1]  # position through the line
points(s, xt,typ="l",lty=2)

x  = 0       # once again, with 100 lines ...
xp = 0
V = c(x,xp)  # values at s=0
V = trk(V,d)
xt  = V[,1]  # position through the line
plot(s, xt,typ="l",ylim=c(-1,1)*20,col="red",lwd=4)
abline(h=0,col="blue")

# draw a series of nearby trajectories
Ntr = 100
i = 0
while(i<Ntr){
   i = i+1
   x  =  rnorm(1,0, 2)
   xp = (rnorm(1,0, 2) - 0.685*x)/100
   V = c(x,xp)  # values at s=0
   V = trk(V,d)
   xt  = V[,1]  # position through the line
   points(s, xt,typ="l",lty=2)
}

The above plot is for one set of misalignments. The red curve is for a beam coming in to the beamline with initial (\(x\),\(x'\)) = (0,0). We see that trajectories in the vicinity of the “ideal” beam will follow along, with potentially larger excursions.

Prior to construction, we do not know which misalignments we will have; so, we look at a distribution of misalignments. Suppose we believe we can align with an accuracy of 0.25 mm (typical of beam line alignment at our working example – the Fermilab Muon Campus). We make a variety of quad misalignment distributions, each with the same rms accuracy, and look at the distribution of final displacements at the end:

# Introduce various sets of random errors;
# look at final trajectory

drms = 0.25
d = rnorm(Nqds+1,0,drms)

V = c(0,0)  # values at s=0
V = trk(V,d)
xt  = V[,1]  # position through the line
plot(s, xt,typ="l",ylim=c(-1,1)*20,
   xlab="path length [m]",ylab="x [mm]")
abline(h=0,col="blue")

# change the quad offsets and track again
Ntrials = 1000
xtL = c(1:Ntrials)*0
i = 1
while(i<Ntrials){
   i = i+1
   d = rnorm(Nqds+1,0,drms)
   V = c(0,0)  # values at s=0
   V = trk(V,d)
   xt  = V[,1]  # position through the line
   points(s, xt,typ="l",lty=2)
   xtL[i] = xt[length(xt)]
}

This shows that – without any intervening corrections – we could expect the final transverse position of the initial ideal beam to grow, up to some maximum amount in the range of \(\pm\) 17 mm or so by the end of the line.

Since an offset quadrupole affects the trajectory by mis-steering the beam: \[ \Delta x' = \frac{d}{F} \] and a particle with initial \(x_0 = 0\) and initial angle \(\Delta x'_0\) will perform a betatron oscillation downstream, \[ \Delta x (s) = \Delta x'_0 \sqrt{\beta(s)\beta(s_0)}\sin[\psi(s)-\psi(s_0)] \] then \[ \Delta x (s) = \sum_n \frac{d_n}{F_n} \sqrt{\beta(s)\beta{s_n}}\sin[\psi(s)-\psi(s_n)] ~~~~~~~~ s_n < s \] So, through a random walk process, we would expect \[ \Delta x(s) \approx \frac{d_{rms}}{|F|}\sqrt{\beta(s)\langle\beta\rangle}\cdot\sqrt{N/2} \]

hist(xtL,breaks=25,xlab="position at end of line [mm]")
text(10,200,paste("d_rms =",drms," mm"))

HistXend = hist(abs(xtL),breaks=100)

HistXend

$breaks
 [1]  0.0  0.2  0.4  0.6  0.8  1.0  1.2  1.4  1.6  1.8  2.0  2.2  2.4  2.6  2.8  3.0  3.2  3.4
[19]  3.6  3.8  4.0  4.2  4.4  4.6  4.8  5.0  5.2  5.4  5.6  5.8  6.0  6.2  6.4  6.6  6.8  7.0
[37]  7.2  7.4  7.6  7.8  8.0  8.2  8.4  8.6  8.8  9.0  9.2  9.4  9.6  9.8 10.0 10.2 10.4 10.6
[55] 10.8 11.0 11.2 11.4 11.6 11.8 12.0 12.2 12.4 12.6 12.8 13.0 13.2 13.4 13.6 13.8 14.0 14.2
[73] 14.4 14.6 14.8 15.0 15.2 15.4 15.6 15.8 16.0 16.2 16.4 16.6 16.8 17.0 17.2 17.4 17.6

$counts
 [1] 104 107  88 113 104  93 103 128  85 109  94  85  99  86  85  95  82  76  74  78  72  62
[23]  52  51  76  54  47  40  51  43  54  34  42  32  32  28  34  25  27  27  25  21  18  14
[45]  11  17  11  12  11  13   7  12   3   4   6   6   3   5   4   1   3   5   5   0   1   3
[67]   0   0   1   2   0   0   0   1   0   1   1   0   0   0   0   0   1   0   0   0   0   1

$density
 [1] 0.173333333 0.178333333 0.146666667 0.188333333 0.173333333 0.155000000 0.171666667
 [8] 0.213333333 0.141666667 0.181666667 0.156666667 0.141666667 0.165000000 0.143333333
[15] 0.141666667 0.158333333 0.136666667 0.126666667 0.123333333 0.130000000 0.120000000
[22] 0.103333333 0.086666667 0.085000000 0.126666667 0.090000000 0.078333333 0.066666667
[29] 0.085000000 0.071666667 0.090000000 0.056666667 0.070000000 0.053333333 0.053333333
[36] 0.046666667 0.056666667 0.041666667 0.045000000 0.045000000 0.041666667 0.035000000
[43] 0.030000000 0.023333333 0.018333333 0.028333333 0.018333333 0.020000000 0.018333333
[50] 0.021666667 0.011666667 0.020000000 0.005000000 0.006666667 0.010000000 0.010000000
[57] 0.005000000 0.008333333 0.006666667 0.001666667 0.005000000 0.008333333 0.008333333
[64] 0.000000000 0.001666667 0.005000000 0.000000000 0.000000000 0.001666667 0.003333333
[71] 0.000000000 0.000000000 0.000000000 0.001666667 0.000000000 0.001666667 0.001666667
[78] 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.001666667 0.000000000
[85] 0.000000000 0.000000000 0.000000000 0.001666667

$mids
 [1]  0.1  0.3  0.5  0.7  0.9  1.1  1.3  1.5  1.7  1.9  2.1  2.3  2.5  2.7  2.9  3.1  3.3  3.5
[19]  3.7  3.9  4.1  4.3  4.5  4.7  4.9  5.1  5.3  5.5  5.7  5.9  6.1  6.3  6.5  6.7  6.9  7.1
[37]  7.3  7.5  7.7  7.9  8.1  8.3  8.5  8.7  8.9  9.1  9.3  9.5  9.7  9.9 10.1 10.3 10.5 10.7
[55] 10.9 11.1 11.3 11.5 11.7 11.9 12.1 12.3 12.5 12.7 12.9 13.1 13.3 13.5 13.7 13.9 14.1 14.3
[73] 14.5 14.7 14.9 15.1 15.3 15.5 15.7 15.9 16.1 16.3 16.5 16.7 16.9 17.1 17.3 17.5

$xname
[1] "abs(xtL)"

$equidist
[1] TRUE

attr(,"class")
[1] "histogram"

Prob = cumsum(HistXend$counts)/sum(HistXend$counts)
plot(HistXend$mids,Prob,typ="l",
   xlab="x [mm]",
   ylab="Probability (|max. displacement| < x)")
abline(h=0.9, lty=2)

abline(v=HistXend$mids[min(which(Prob>0.9))], col="red",lty=3)
text(12,0.6,paste("d_rms =",drms," mm"))

Thus for an rms alignment accuracy of 0.25 mm, there is a 90% chance that the final displacement will be below 7.5 mm.

As a check:

x90 = 1.645      #  +\- 1.645 sigma contains 90% of cases
pnorm(x90)-pnorm(-x90)  

[1] 0.9000302

betaMax = 2*F*sqrt((1+L/2/F)/(1-L/2/F))
betaMin = 2*F*sqrt((1-L/2/F)/(1+L/2/F))
sd(xtL)

[1] 4.469646

drms/F*sqrt((betaMax+betaMin)/2*betaMax)*sqrt(Nqds/2)

[1] 4.401829

sd(xtL)*x90

[1] 7.352568