#!/bin/bash

# Make main beam file
elegant run1.ele

# Get the number of particles in the main beam so we can 
# compute the starting particleID for the probe particles. 
# This allows us to identify them after tracking
np=`sdds2stream run1.bun -rows=bare`

# Make probe particles
sddssequence -pipe=out -define=delta -sequence=begin=-0.01,end=0.01,number=11 \
    | sddsprocess -pipe=in probeParticles.sdds \
    -define=parameter,pCentral,10.7382941763 \
    "-define=col,x,0,units=m" \
    "-define=col,y,0,units=m" \
    "-define=col,xp,0" \
    "-define=col,yp,0" \
    "-define=col,t,0,units=s" \
    "-define=col,p,pCentral 1 delta + *" \
    "-define=col,particleID,$np i_row + 1 +,type=long" 

# Combine
sddscombine run1.bun probeParticles.sdds beam.sdds -merge -overwrite

# Track the combined beam + probe particles
elegant run2.ele 

# Get the final data for the probe particles only
sddsxref run2.out probeParticles.sdds -pipe=out -equate=particleID -take=* -nowarning -edit=col,*,ei/0/ \
    | sddsprocess -pipe -define=col,delta,"p pCentral / 1 -" \
    | sddsmpfit -pipe=in run2.fit \
    -terms=2 -indep=delta0 -depen=x,xp,y,yp,t,delta 

# plot
sddsplot -groupby=nameindex -sep=nameindex \
         -column=delta0,'(x,xp,y,yp,t,delta)' run2.fit -graph=sym,scale=2,fill \
         -column=delta0,'(xFit,xpFit,yFit,ypFit,tFit,deltaFit)' run2.fit 

# Display slopes (Ri6 elements)
sddscollapse run2.fit -pipe=out \
    | sddsprintout -pipe -column=*Slope
