Yoke Movement Plan

Code Preamble

In [1]:
cd ~/Workspace/gm2/studies/field-shimming/pole-geometry
/Users/matthias/Workspace/gm2/studies/field-shimming/pole-geometry
In [2]:
%matplotlib inline

import numpy as np
import matplotlib as mpl
import matplotlib.pyplot as plt
import ROOT as r
import seaborn as sns

# Set the figure size to span the notebook.
mpl.rcParams['figure.figsize'] = (18, 6)
mpl.style.use('fivethirtyeight')
mpl.rcParams['axes.grid'] = False

colors = sns.hls_palette(10, l=0.3, s=0.8)

def finish_plot():
    plt.xlabel(r'$\theta$ [deg]')
    plt.xlim(0, 360)
    plt.ylabel(r'z [$\mu$ m]')
    
    yoke_label = ord('A')
    
    for i in xrange(0, 37):

        x = (i * 10.0 - 15.0 - 0.1012) % 360.0

        if i % 3 == 0:
            plt.axvline(x, linestyle='-', color='k', alpha=0.2)

        else:
            plt.axvline(x, linestyle='--', color='k', alpha=0.2)
        
        if i is not 0:
            if i < 10:
                s = str(i) + ' '
            else:
                s = str(i)
                
            plt.figtext(0.052 + i * (0.868 / 36.0), 0.86, s, color='k', alpha=0.4)
            
        if (i % 3) == 1:
            plt.figtext(0.052 + i * (0.868 / 36.0), 0.82, chr(yoke_label), color='k', alpha=0.4)
            yoke_label += 1
            
    plt.show()


def tilt_plane(phi):
    return 832.1725450408629 * (np.cos((phi - 171.41451281) * np.pi / 180.0) + 1.0588562)
Welcome to ROOTaaS 6.06/00

Current Ring Geometry

We currently have two measurements that can be used to reconstruct the geometry of the ring (lower poles if you want to be specific). We can use a set of tilt measurements at the center of the poles in conjunction with the boundaries between poles as described in elog 432 and elog 430. We can also use the laser tracker device that follows our cart around the ring, and give a relative value that is calibrated against gravity (an example in elog 429). The two methods both agree within 100 $\mu m$ as shown below.

In [3]:
# Open the data for pole boundaries (elog 432) and azimuthal tilts (elog 430).
d_pb = np.genfromtxt('data-in/tilt_pole_boundaries_02172016.csv', skip_header=True)
d_pt = np.genfromtxt('data-in/PoleTilts_02_16_2016.csv', skip_header=True)

# Take the average of the two steps.
dz = (2 * d_pb[:, 3] - d_pb[:, 1] - d_pb[:, 5]) * 0.5

# Convert the value to microns.
dz *= 1.29 * 2.25 * 0.0254

# Extract the slope from dataset.
dzds = d_pt[:, 4]

# Allocate a vector to hold the height values.
z = np.empty(dz.shape[0] * 100)
phi = np.linspace(-5, 355, z.shape[0])

# And fill it using slopes and steps.
for i in xrange(z.shape[0]):
    z[i] = dz[:i/100].sum()
    z[i] += dzds[i/100] * 0.012 * (i % 100)
    z[i] += dzds[:i/100].sum() * 1.2

# Divide the different by the number of weight of pieces in between z_i and z_f.
dzds0 = (z[-1] - z[0] - dz[-1]) / (1.2 * 35.99)

# Add the offset to the slopes.
dzds = d_pt[:, 4] - dzds0

# And fill it again with the correction.
for i in xrange(z.shape[0]):
    z[i] = dz[:i/100].sum()
    z[i] += dzds[i/100] * 0.012 * (i % 100)
    z[i] += dzds[:i/100].sum() * 1.2


# Open the laser data from the last full scan.
f = r.TFile('data-in/full_scan_021.root')
t = f.Get('t')

z_laser_on = np.empty(t.GetEntries())
phi_laser_on = np.empty(t.GetEntries())

for i in xrange(t.GetEntries()):
    t.GetEntry(i)
    z_laser_on[i] = t.z_2 * 1e6 # Convert inches to microns
    phi_laser_on[i] = t.phi_2

# Remove points with bad data.    
indices = np.where(z_laser_on > 0.5)
z_laser_on = z_laser_on[indices]
phi_laser_on = phi_laser_on[indices]
    
# Subtract the height offset to put things in the same space as our reconstruction.
z_laser_on -= (z_laser_on[500] - z[int((phi_laser_on[500] + 5) * 10 + 0.5)])

f = r.TFile('data-in/full_scan_021.root')
t = f.Get('t')

z_laser_off = np.empty(t.GetEntries())
phi_laser_off = np.empty(t.GetEntries())

for i in xrange(t.GetEntries()):
    t.GetEntry(i)
    z_laser_off[i] = t.z_2 * 1e6 # Convert inches to microns
    phi_laser_off[i] = t.phi_2

# Subtract the height offset to put things in the same space as our reconstruction.
z_laser_off -= (z_laser_off[500] - z[int((phi_laser_off[500] + 5) * 10 + 0.5)])
    
# Remove points with bad data.    
indices = np.where(z_laser_off - tilt_plane(phi_laser_off) > -1000)
z_laser_off = z_laser_off[indices]
phi_laser_off = phi_laser_off[indices]
                
plt.clf()
plt.scatter(phi_laser_off, z_laser_off, color=colors[2])
plt.scatter(phi_laser_on, z_laser_on, color=colors[1])
plt.scatter(phi, z, color=colors[0])
finish_plot()

plt.clf()
plt.scatter(phi_laser_off, z_laser_off - tilt_plane(phi_laser_off), color=colors[2])
plt.scatter(phi_laser_on, z_laser_on - tilt_plane(phi_laser_on), color=colors[1])
plt.scatter(phi, z - tilt_plane(phi), color=colors[0])
finish_plot()

Yoke Movements

Taking confidence from the tilt sensor reconstruction, we will actually use the laser tracker data taken with the field on to advise movements on each individual yoke. Each yoke has 4 jacks supported it, so we can adjust the height, pitch and roll individually. The plan from a meeting between field team members, engineers and management was to dividide and conquer into yokes D-I and A-C plus J-L. The plan accomadates a finite supply of dynamic support jacks, and has the advantage that one set only needs to be raised while the other set only needs to be lowered. Intially we treat each yoke movement as a height adjustment and a pitch adjustment (azimuthal tilt).

In [4]:
from scipy.stats import linregress

yoke_z = np.zeros(12)
yoke_m = np.zeros(12)

for i in xrange(12):

    if i == 0:
        indices1 = np.where((phi_laser_on + 180.0) % 360.0 < 195.0)
        indices2 = np.where((phi_laser_on[indices1] + 180.0) % 360.0 >= 155.0)
    
    else:
        indices1 = np.where(phi_laser_on >= i * 30.0 - 15.0)
        indices2 = np.where(phi_laser_on[indices1] < i * 30.0 + 15.0)
    
    temp_z = z_laser_on[indices1][indices2]
    temp_phi = phi_laser_on[indices1][indices2]

    yoke_m[i], x, x, x, x = linregress(temp_phi, temp_z)
    yoke_z[i] = temp_z.mean()
    
    # Can convert to microradians like this.
    # yoke_m[i] *= 360.0 / (2 * np.pi * 6.99) 

plt.clf()
plt.scatter(range(0, 360, 30), yoke_z, s=150.0, color=colors[0])

z0 = 734.2
yoke_label = ord('A')
plt.axhline(z0, linestyle='--', color=colors[1])
print "Yoke\tUpstream\tDownstream\t"

for i in xrange(12):
    x = np.linspace(i * 30.0 - 15.0, i * 30.0 + 15.0, 100)
    y = np.ones(100) * yoke_z[i]
    y += (x - i * 30) * yoke_m[i]
    print "{}\t{: 8.2f}\t{: 8.2f}\t".format(chr(yoke_label), z0 - y[0], z0 - y[-1])
    yoke_label += 1
    plt.scatter(x % 360, y, color=colors[0])

finish_plot()
Yoke	Upstream	Downstream	
A	  808.70	  809.62	
B	  777.51	  454.39	
C	  349.42	  -78.84	
D	 -200.00	 -529.98	
E	 -488.23	 -826.46	
F	 -836.35	 -817.72	
G	 -865.18	 -750.55	
H	 -732.88	 -519.03	
I	 -511.76	 -311.91	
J	 -276.89	   67.23	
K	   58.35	  461.27	
L	  396.46	  843.82	
In [5]:
# And make a plot with laser data to compare.
for i in xrange(12):
    x = np.linspace(i * 30.0 - 15.0, i * 30.0 + 15.0, 100)
    y = np.ones(100) * yoke_z[i]
    y += (x - i * 30) * yoke_m[i]
    plt.scatter(x % 360, y, color=colors[0])

plt.scatter(phi_laser_on, z_laser_on, color=colors[1])
finish_plot()

# And make a plot with laser data to compare.
z_laser_on_corr = np.zeros(z_laser_on.shape)

for i in xrange(z_laser_on.shape[0]):

    for j in xrange(12):

        phi_i = j * 30.0 - 15.0
        phi_f = j * 30.0 + 15.0

        if (phi_i <= phi_laser_on[i] < phi_f):
    
            z_laser_on_corr[i] = (phi_laser_on[i] - j * 30) * yoke_m[j] + yoke_z[j]

plt.scatter(phi_laser_on, z_laser_on - z_laser_on_corr, color=colors[1])
plt.axhline(0, color='k', alpha=0.2)
finish_plot()
In [9]:
# At second glance, we need to make a different model that moves the pole interfaces.

yoke_if_iz = np.zeros(12)
yoke_if_oz = np.zeros(12)

from scipy.optimize import curve_fit

def pol(x, p0, p1, p2):
    return p0 + p1 * x + p2 * x**2


for i in xrange(12):

    if i == 11:
        indices1 = np.where((phi_laser_on + 180.0) % 360.0 < 175.0)
        indices2 = np.where((phi_laser_on[indices1] + 180.0) % 360.0 >= 155.0)
    
    else:
        indices1 = np.where(phi_laser_on >= i * 30.0 + 5.0)
        indices2 = np.where(phi_laser_on[indices1] < i * 30.0 + 25.0)
    
    temp_z = z_laser_on[indices1][indices2]
    temp_phi = phi_laser_on[indices1][indices2]

    p0 = [temp_z.mean(), 0.0, 0.0]
    
    fit, corr = curve_fit(pol, temp_phi - i * 30.0 - 15.0, temp_z, p0=p0)

    plt.scatter(temp_phi, temp_z, color=colors[1])
    plt.scatter(temp_phi, pol(temp_phi - i * 30.0 - 15.0, fit[0], fit[1], fit[2]), color=colors[2])

    yoke_if_iz[i] = fit[0]
    
    plt.scatter((i * 30 + 15) % 360, yoke_if_iz[i], s=150.0, color=colors[0])

finish_plot()

# Print the values.
z0 = 734.2
yoke_label = ord('A')

print "Yoke-IF\tInner-dZ"

for i in xrange(12):
    if i == 11:
        print "{}/{}\t{: 8.2f}".format(chr(yoke_label), chr(yoke_label - 11), z0 - yoke_if_iz[i])
    else:
        print "{}/{}\t{: 8.2f}".format(chr(yoke_label), chr(yoke_label + 1), z0 - yoke_if_iz[i])
    yoke_label += 1

# The residuals after adjusting with this model.
yoke_label = ord('A')
z_laser_on_corr = np.zeros(z_laser_on.shape)

for i in xrange(z_laser_on.shape[0]):

    for j in xrange(12):

        phi_i = j * 30.0 - 15.0
        phi_f = j * 30.0 + 15.0
        
        if j == 2:
            m = (z0 + 45.0 - yoke_if_iz[j-1]) / 30.0
            y0 = yoke_if_iz[j-1]
        elif j == 3:
            m = (yoke_if_iz[j] - (z0 + 45.0)) / 30.0
            y0 = z0 + 45.0
        else:
            m = (yoke_if_iz[j] - yoke_if_iz[j-1]) / 30.0
            y0 = yoke_if_iz[j-1]

        if (phi_i <= phi_laser_on[i] < phi_f):

            z_laser_on_corr[i] = (phi_laser_on[i] - j * 30 + 15) * m + y0

plt.scatter(phi_laser_on, z_laser_on - z_laser_on_corr, color=colors[1])
plt.axhline(0, color='k', alpha=0.2)
finish_plot()
Yoke-IF	Inner-dZ
A/B	  733.69
B/C	  378.18
C/D	 -105.88
D/E	 -480.08
E/F	 -815.55
F/G	 -838.74
G/H	 -751.30
H/I	 -542.78
I/J	 -307.00
J/K	   50.78
K/L	  428.53
L/A	  819.99
In [7]:
# Calculate the outer ring
if_phi = np.linspace(15.0, 345.0, 12)

plt.plot(if_phi % 360, yoke_if_iz - tilt_plane(if_phi % 360))
finish_plot()

yoke_if_oz = yoke_if_iz + (1.042 / 7.005) * (tilt_plane(if_phi) - 832.1725450408629 * 1.0588562)

print "Yoke-IF\tInner-dZ\tOuter-dZ"

for i in xrange(12):
    if i == 11:
        print "{}/{}\t{: 8.2f}\t{: 8.2f}".format(chr(yoke_label), chr(yoke_label - 11), z0 - yoke_if_iz[i], z0 - yoke_if_oz[i])
    else:
        print "{}/{}\t{: 8.2f}\t{: 8.2f}".format(chr(yoke_label), chr(yoke_label + 1), z0 - yoke_if_iz[i], z0 - yoke_if_oz[i])
    yoke_label += 1
Yoke-IF	Inner-dZ	Outer-dZ
A/B	  733.69	  847.14
B/C	  378.18	  451.67
C/D	 -105.88	  -92.05
D/E	 -480.08	 -529.61
E/F	 -815.55	 -915.17
F/G	 -838.74	 -961.75
G/H	 -751.30	 -864.75
H/I	 -542.78	 -616.26
I/J	 -307.00	 -320.83
J/K	   50.78	  100.31
K/L	  428.53	  528.15
L/A	  819.99	  943.00

Bottom Line

The current plan is to adjust ring elevation and remove the tilt plane first. This involves a differential movement on the C/D yoke interface and fixed height adjustments at all other interfaces. Therefore, the second model should be used for all yoke movements except interfaces C/D. That yoke can be taken from the first model's values for yoke C downstream and yoke D upstream. The final plan for removing the tilt plane and eleveation changes is then:

Yoke Inner Change Outer Change
A/B 733.69 847.14
B/C 378.18 451.67
C/D -45.00 -31.17
D/E -480.08 -529.61
E/F -815.55 -915.17
F/G -838.74 -961.75
G/H -751.30 -864.75
H/I -542.78 -616.26
I/J -307.00 -320.83
J/K 50.78 100.31
K/L 428.53 528.15
L/A 819.99 943.00
In [ ]: