Tracking¶
A. Petrenko (Novosibirsk, 2019)
In [1]:
import numpy as np
import pandas as pd
import holoviews as hv
import os
from IPython.display import Latex
from IPython.display import Image
hv.extension("bokeh")
In [2]:
%opts Curve Spread [width=700 height=300 show_grid=True \
default_tools=['box_zoom','pan','wheel_zoom','reset']]
Elegant lattice file:
Plot sdds-beam¶
In [3]:
w = "../p-linac/FODO_and_e-p_sep/results/w3.sdds"
In [4]:
#!plaindata2sdds
In [5]:
png = "results/image.png"
!sddsplot -lay=2,2 -device=lpng -output=$png \
-col=x,y -graph=dot,type=2 $w -endPanel \
-col=yp,y -graph=dot,type=2 $w -endPanel \
-col=x,xp -graph=dot,type=2 $w -endPanel \
-col=yp,xp -graph=dot,type=2 $w
Image(png)
Out[5]:
In [6]:
!sddsplot -device=lpng -output=$png -lay=2,2 "-title=" \
-col=t,p -factor=yMult=0.511,xMult=-3e11 \
"-xlabel=z (mm)" "-ylabel=p (MeV/c)" -graph=dot,type=2 $w -endPanel \
-col=xp,p -factor=yMult=0.511,xMult=1e3 \
"-xlabel=xp (mrad)" "-ylabel=p (MeV/c)" -graph=dot,type=2 $w -endPanel \
-col=t,x -factor=xMult=-3e11,yMult=1e3 \
"-xlabel=z (mm)" "-ylabel=x (mm)" -graph=dot,type=2 $w -endPanel \
-col=xp,x -factor=yMult=1e3,xMult=1e3 \
"-xlabel=xp (mrad)" "-ylabel=x (mm)" -graph=dot,type=2 $w -endPanel
Image(png)
Out[6]:
Run sddsanalyzebeam to check the beam parameters:
In [7]:
beam_analyzed = 'results/beam_analyzed.sdds'
sdds = "../p-linac/FODO_and_e-p_sep/results/w3.sdds"
!sddsanalyzebeam $sdds $beam_analyzed
In [8]:
!sddsprintout -par=PassLength $sdds
In [9]:
items = np.array([
# name in file symbol factor units
['enx', '\\epsilon_{nx}', 1.0e6, 'mm \\cdot mrad'],
['betax', '\\beta_x', 1.0, 'm' ],
['alphax', '\\alpha_x', 1.0, '' ],
['Sx', '\\sigma_x', 1.0e3, 'mm' ],
['Sxp', "\\sigma_{x'}", 1.0e3, 'mrad' ],
['eny', '\\epsilon_{ny}', 1.0e6, 'mm \\cdot mrad'],
['betay', '\\beta_y', 1.0, 'm' ],
['alphay', '\\alpha_y', 1.0, '' ],
# ['Sy', '\\sigma_y', 1.0e3, 'mm' ],
# ['Syp', "\\sigma_{y'}", 1.0e3, 'mrad' ],
['St', '\\sigma_t', 1e12, 'ps' ],
['St', '\\sigma_z', 3e11, 'mm' ],
['Sdelta', '\\sigma_{\\Delta p/p}', 100.0, '\\%' ],
['pAverage', 'p_{average}', 0.511, 'MeV/c' ],
['Charge', 'Q', 1.0e9, 'nC' ],
])
names = items[:,0]
cols_str = "-col=" + ",".join(names)
cols_str
Out[9]:
In [10]:
out = !sdds2stream $beam_analyzed $cols_str -pipe=out
values = out[0].split()
values = np.array(values)
values = values.astype(np.float)
values = dict(zip(names, values))
In [11]:
latex_str = r"\begin{aligned}"
for name,symbol,factor,units in items:
value = values[name]*factor.astype(float)
latex_str = latex_str + r"%s &= %.3f~\rm{%s} \\" % (symbol, value, units)
latex_str = latex_str + r"\end{aligned}"
In [12]:
Latex(latex_str)
Out[12]:
In [13]:
with open("machine.lte", "r") as f:
print(f.read())
In [14]:
with open("track.ele", "r") as f:
print(f.read())
Run Elegant
In [15]:
#!elegant track.ele
In [16]:
import subprocess
import sys
def run_cmd(cmd):
sub_process = subprocess.Popen(cmd, shell=True,
stdout=subprocess.PIPE, stderr=subprocess.PIPE)
while sub_process.poll() is None:
out = sub_process.stdout.read(10)
sys.stdout.write(out)
sys.stdout.flush()
In [17]:
run_cmd("elegant track.ele")
In [18]:
w = "results/w0.sdds"
!sddsplot -lay=2,2 -device=lpng -output=$png \
-col=x,y -graph=dot,type=2 $w -endPanel \
-col=yp,y -graph=dot,type=2 $w -endPanel \
-col=x,xp -graph=dot,type=2 $w -endPanel \
-col=yp,xp -graph=dot,type=2 $w
Image(png)
Out[18]:
In [19]:
!sddsplot -device=lpng -output=$png -lay=2,2 "-title=" \
-col=dt,p -factor=yMult=0.511,xMult=-3e11 \
"-xlabel=z (mm)" "-ylabel=p (MeV/c)" -graph=dot,type=2 $w -endPanel \
-col=xp,p -factor=yMult=0.511,xMult=1e3 \
"-xlabel=xp (mrad)" "-ylabel=p (MeV/c)" -graph=dot,type=2 $w -endPanel \
-col=dt,x -factor=xMult=-3e11,yMult=1e3 \
"-xlabel=z (mm)" "-ylabel=x (mm)" -graph=dot,type=2 $w -endPanel \
-col=xp,x -factor=yMult=1e3,xMult=1e3 \
"-xlabel=xp (mrad)" "-ylabel=x (mm)" -graph=dot,type=2 $w -endPanel
Image(png)
Out[19]:
In [20]:
w = "results/w3.sdds"
!sddsplot -lay=2,2 -device=lpng -output=$png \
-col=x,y -graph=dot,type=2 $w -endPanel \
-col=yp,y -graph=dot,type=2 $w -endPanel \
-col=x,xp -graph=dot,type=2 $w -endPanel \
-col=yp,xp -graph=dot,type=2 $w
Image(png)
Out[20]:
In [21]:
!sddsplot -device=lpng -output=$png -lay=2,2 "-title=" \
-col=dt,p -factor=yMult=0.511,xMult=-3e11 \
"-xlabel=z (mm)" "-ylabel=p (MeV/c)" -graph=dot,type=2 $w -endPanel \
-col=xp,p -factor=yMult=0.511,xMult=1e3 \
"-xlabel=xp (mrad)" "-ylabel=p (MeV/c)" -graph=dot,type=2 $w -endPanel \
-col=dt,x -factor=xMult=-3e11,yMult=1e3 \
"-xlabel=z (mm)" "-ylabel=x (mm)" -graph=dot,type=2 $w -endPanel \
-col=xp,x -factor=yMult=1e3,xMult=1e3 \
"-xlabel=xp (mrad)" "-ylabel=x (mm)" -graph=dot,type=2 $w -endPanel
Image(png)
Out[21]:
Now let's plot the results with holoviews¶
In [22]:
os.listdir("results")
Out[22]:
Reading beamline.mag
In [23]:
out = !sdds2stream results/beamline.mag -col=ElementName,s,Profile -pipe=out
In [24]:
out[:3]
Out[24]:
In [25]:
print("\n".join(out[:3]))
In [26]:
from io import StringIO
In [27]:
DATA = StringIO("\n".join(out))
In [28]:
df = pd.read_csv(DATA, names=['ElementName', 's', 'Profile'], delim_whitespace=True)
In [29]:
df.head(3)
Out[29]:
In [30]:
df.tail(3)
Out[30]:
In [31]:
dim_s = hv.Dimension('s', unit='m', label="s") # range=(-80,0))
In [32]:
#hv.help(hv.Curve)
In [33]:
from bokeh.models import HoverTool
#hover = HoverTool(tooltips=[("Name", "@ElementName")])
hover = HoverTool(tooltips="@ElementName")
%opts Curve.mag [height=70 show_frame=False show_title=False \
xaxis=None yaxis=None tools=['xbox_zoom, xpan', hover]] (color='black', alpha=0.3)
mag = hv.Curve(df, kdims=dim_s, vdims=['Profile', 'ElementName'], group='mag')
mag
Out[33]:
Reading track.sig file
In [34]:
#!sddsquery results/track.sig
In [35]:
out = !sdds2stream results/track.sig -col=ElementName,s,minimum1,maximum1,minimum3,maximum3,Sx,Sy -pipe=out
DATA = StringIO("\n".join(out))
df = pd.read_csv(DATA, names=['ElementName','s','xmin','xmax','ymin','ymax','sigma_x','sigma_y'],
delim_whitespace=True)
In [36]:
df.tail(3)
Out[36]:
In [37]:
df['xmin'] = df['xmin']*1e3; df['xmax'] = df['xmax']*1e3 # mm
df['ymin'] = df['ymin']*1e3; df['ymax'] = df['ymax']*1e3 # mm
df['sigma_x'] = df['sigma_x']*1e3; df['sigma_y'] = df['sigma_y']*1e3 # mm
In [38]:
#df.tail(3)
In [39]:
dim_x = hv.Dimension('x', unit='mm', label="x", range=(0,None))
dim_y = hv.Dimension('y', unit='mm', label="y", range=(0,None))
%opts Curve.fx (color='red', alpha=0.7, line_width=3)
%opts Curve.fy (color='blue', alpha=0.7, line_width=3)
In [40]:
Sx = hv.Curve((df.s, df.sigma_x), label='σx', kdims=dim_s, vdims=dim_x, group='fx')
Sy = hv.Curve((df.s, df.sigma_y), label='σy', kdims=dim_s, vdims=dim_y, group='fy')
In [41]:
(Sx*Sy + mag).cols(1)
Out[41]:
Reading track.cen file
In [42]:
#!sddsquery results/track.cen
In [43]:
out = !sdds2stream results/track.cen -col=ElementName,s,Particles,pCentral,Cx,Cy,Charge -pipe=out
DATA = StringIO("\n".join(out))
df_cen = pd.read_csv(DATA, names=['ElementName','s','Particles','p','Cx','Cy','Charge'],
delim_whitespace=True)
df_cen['p'] = 0.511*df_cen['p'] # MeV/c
df_cen['Cx'] = 1e3*df_cen['Cx'] # mm
df_cen['Cy'] = 1e3*df_cen['Cy'] # mm
In [44]:
df_cen.tail(3)
Out[44]:
In [45]:
dim_Q = hv.Dimension('Q', range=(0,None), unit='nC') # range=(-80,0))
Charge = hv.Curve((df_cen.s, df_cen.Charge*1e9), kdims=dim_s, vdims=dim_Q)
In [46]:
(Charge + mag).cols(1)
Out[46]:
Reading twiss.sdds file
In [47]:
out = !sdds2stream results/twiss.twi \
-col=ElementName,s,betax,betay,alphax,alphay,etax,etay,pCentral0,xAperture,yAperture -pipe=out
DATA = StringIO("\n".join(out))
df_twi = pd.read_csv(DATA, names=['ElementName','s','betax','betay','alphax','alphay',
'etax','etay','p','xAperture','yAperture'],
delim_whitespace=True)
In [48]:
df_twi['xAperture'] = 1e3*df_twi['xAperture'] # mm
df_twi['yAperture'] = 1e3*df_twi['yAperture'] # mm
In [49]:
df_twi.tail(3)
Out[49]:
In [50]:
dim_x = hv.Dimension('x', unit='mm', label="x", range=(-30,+30))
dim_y = hv.Dimension('y', unit='mm', label="y", range=(-30,+30))
dim_p = hv.Dimension('p', unit='MeV/c', label="p", range=(0,None))
In [51]:
%opts Curve.Aper (color='gray', alpha=0.5, line_width=3)
dim_sigma = hv.Dimension('sigma', unit='mm', label="2σ")#, range=(0,None))
dim_aper = hv.Dimension('aper', unit='mm', label="Aperture")#, range=(0,None))
Ax = hv.Curve((df_twi.s, df_twi.xAperture), label='Aperture', kdims=dim_s, vdims=dim_sigma, group='Aper')
Ax_pos = Ax
Ax_neg = hv.Curve((df_twi.s, -df_twi.xAperture), label='Aperture', kdims=dim_s, vdims=dim_sigma, group='Aper')
In [52]:
#Ax_pos*Ax_neg
In [53]:
#s_p = hv.Curve((df_cen.s, df_cen.p), kdims=dim_s, vdims=dim_p)
Cx = hv.Spread((df_cen.s, df_cen.Cx, 2*df.sigma_x), label='x', \
kdims=dim_s, vdims=[dim_x,dim_aper])
Cy = hv.Spread((df_cen.s, df_cen.Cy, 2*df.sigma_y), label='y', \
kdims=dim_s, vdims=[dim_y,dim_aper])
In [54]:
(Ax_pos*Ax_neg*Cy*Cx + mag).cols(1)
Out[54]:
In [55]:
s_Cx = hv.Curve((df_cen.s, df_cen.Cx), kdims=dim_s, vdims=dim_x, group='fx', label='x')
s_Cy = hv.Curve((df_cen.s, df_cen.Cy), kdims=dim_s, vdims=dim_x, group='fy', label='y')
In [56]:
(s_Cx*s_Cy*Ax_pos*Ax_neg + mag).cols(1)
Out[56]:
Resulting beam paramaters:
In [57]:
beam_analyzed = 'results/beam_analyzed.sdds'
w = "results/w3.sdds"
!sddsanalyzebeam $w $beam_analyzed
In [58]:
items = np.array([
# name in file symbol factor units
['enx', '\\epsilon_{nx}', 1.0e6, 'mm \\cdot mrad'],
['betax', '\\beta_x', 1.0, 'm' ],
['alphax', '\\alpha_x', 1.0, '' ],
['Sx', '\\sigma_x', 1.0e3, 'mm' ],
['Sxp', "\\sigma_{x'}", 1.0e3, 'mrad' ],
['eny', '\\epsilon_{ny}', 1.0e6, 'mm \\cdot mrad'],
['betay', '\\beta_y', 1.0, 'm' ],
['alphay', '\\alpha_y', 1.0, '' ],
# ['Sy', '\\sigma_y', 1.0e3, 'mm' ],
# ['Syp', "\\sigma_{y'}", 1.0e3, 'mrad' ],
['St', '\\sigma_t', 1e12, 'ps' ],
['St', '\\sigma_z', 3e11, 'mm' ],
['Sdelta', '\\sigma_{\\Delta p/p}', 100.0, '\\%' ],
['pAverage', 'p_{average}', 0.511, 'MeV/c' ],
['Charge', 'Q', 1.0e9, 'nC' ],
])
names = items[:,0]
cols_str = "-col=" + ",".join(names)
cols_str
Out[58]:
In [59]:
out = !sdds2stream $beam_analyzed $cols_str -pipe=out
values = out[0].split()
values = np.array(values)
values = values.astype(np.float)
values = dict(zip(names, values))
In [60]:
latex_str = r"\begin{aligned}"
for name,symbol,factor,units in items:
value = values[name]*factor.astype(float)
latex_str = latex_str + r"%s &= %.3f~\rm{%s} \\" % (symbol, value, units)
latex_str = latex_str + r"\end{aligned}"
In [61]:
from IPython.display import Latex
In [62]:
Latex(latex_str)
Out[62]:
In [ ]: