Module Contents

This documentation is automatically generated by scanning all the source code. Parts may be incomplete.

pymad8 - python tools for working with MAD8 output and input.

Authors:

  • Stewart Boogert

  • Laurie Nevay

  • Andrey Abramov

  • William Shields

  • Jochem Snuverink

  • Stuart Walker

  • Marin Deniaud

Dependencies: (package - minimum version required)

  • numpy - 1.7.1

  • matplotlib - 1.3.0

  • pylab - 1.3.0 (dependancy of matplotlib)

  • pandas - 1.4.3

  • fortranformat - 1.2.0

Modules: (script name - usage)

  • Input - Tidy Mad8 input

  • Output - Load Mad8 files into dataframes

  • Plot - Draw machine lattice and quick optics plots

  • Sim - Perform simulations on a machine, like particle tracking

  • Track - Old particle tracking code

  • Visualisation - Old survey plotting code

Copyright Royal Holloway, University of London 2019.

pymad8.Input module

pymad8.Input.decodeFileLine(input)
Decode line for each element type
input : string of a mad8 line
pymad8.Input.removeComments(input)
Remove comment lines
input : list of file lines
pymad8.Input.removeContinuationSymbols(input)
Remove continuation symbols from input
input : list of file lines
pymad8.Input.tidy(input)
Tidy input, remove EOL, remove empty lines
input : list of file lines

pymad8.Output module

class pymad8.Output.Output(filename, filetype='twiss')

Bases: object

Class to load different Mad8 output files in a Pandas DataFrame
>>> twiss = pymad8.Output(‘/twiss.tape’)
>>> rmat = pymad8.Output(‘/rmat.tape’,’rmat’)
>>> chrom = pymad8.Output(‘/chrom.tape’,’chrom’)
>>> envel = pymad8.Output(‘/envel.tape’,’envel’)
>>> survey = pymad8.Output(‘/survey.tape’,’survey’)
>>> struct = pymad8.Output(‘/struct.tape’,’struct’)
>>> track = pymad8.Output(‘/track.tape’,’track’)
>>> optics = pymad8.Output(‘/optics.tape’,’optics’)

By default the filetype is twiss
Clear()

Empties all data structures in this instance

calcBeamSize(EmitX, EmitY, Esprd, BunchLen=0)

Calculate the beam sizes and beam divergences in both planes for all elements Then the four columns are added at the end of the DataFrame Works only if a twiss file was loaded previously

getAperture(index, defaultAperSize=0.1)
Get aperture of an element using corresponding index
If none provided or is 0 : set to a default aperture of 0.1m
getColumnsByKeys(keylist)

Return Columns that ar in keylist

getElementByIndex(indexlist, keylist)

Return value of elements given their indices and for chosen columns

getElementByNames(namelist, keylist)

Return value of elements given their names and for chosen columns

getElementByNearestS(s, keylist)

Return value of elements at a given s position and for chosen columns

getIndexByNames(namelist)

Return Index or Index list of elements that are in namelist

getIndexByNearestS(s)

Return Index of the closest element in the beamline

getIndexByTypes(typelist)

Return Index or Index list of elements that are in typelist

getIndexByValues(**args)
Return Index or Index list of elements that have certain values for a chosen column
Same arguments as getRowsByValues
getNameByNearestS(s)

Return Name of the closest element in the beamline

getNamesByIndex(indexlist)

Return Name or Name list of elements that are in indexlist

getNamesByTypes(typelist)

Return Name or Name list of elements that are in typelist

getNamesByValues(**args)
Return Name or Name list of elements that have certain values for a chosen column
Same arguments as getRowsByValues
getRowByNearestS(s)

Return Rows of the closest element in the beamline

getRowsByFunction(f)

Return a sub-datafarme using a boolean function

getRowsByIndex(indexlist)

Return Rows of elements that are in indexlist

getRowsByNames(namelist)

Return Rows of elements that are in namelist

getRowsByTypes(typelist)

Return Rows of elements that are in typelist

getRowsByValues(key=None, minValue=-inf, maxValue=inf, equalValues=None)

Return Rows of elements that have certain values for a chosen column

getTypeByNearestS(s)

Return Type of the closest element in the beamline

getTypesByIndex(indexlist)

Return Type or Type list of elements that are in indexlist

getTypesByNames(namelist)

Return Type or Type list of elements that are in namelist

getTypesByValues(**args)
Return Type or Type list of elements that have certain values for a chosen column
Same arguments as getRowsByValues
plotXY(Xkey, Ykey, color=None, label=None)

Quick plot of one colums of our dataframe w.r.t. another

sMax()

Return maximal S value

sMin()

Return minimal S value

subline(start, end)

Select only a portion of the inital lattice

pymad8.Plot module

pymad8.Plot.AddMachineLatticeToFigure(figure, mad8opt, tightLayout=True)
Add a diagram above the current graph in the figure that represents the accelerator based on a mad8 twiss file.

Note you can use matplotlib’s gcf() ‘get current figure’ as an argument.

>>> pymad8.Plot.AddMachineLatticeToFigure(gcf(), ‘/mad8_twiss_tape’)
class pymad8.Plot.Optics(mad8file, beamParams=None)

Bases: object

Class to load pymad8 DataFrames and make optics plots
>>> plot_data = pymad8.Plot.Optics(‘/mad8_twiss_tape’)
>>> plot_data.Betas()

Plot avaliable are plotBetas, plotAlphas, plotMus, plotDisp and plotSigmas
Alpha()

Plot the Alpha functions for both planes and both Mad

Beta()

Plot the Beta functions for both planes and both Mad

Disp()

Plot the Dispertion functions for both planes and both Mad

Mu()

Plot the Mu functions for both planes and both Mad

Sigma()

Plot the beam size and beam divergence functions for both planes and both Mad

pymad8.Sim module

pymad8.Sim.CheckUnits(coord)
class pymad8.Sim.Track_Collection(E_0)

Bases: object

Class to create and store different initial particles for future tracking
Requires only the reference enegy in GeV
>>> track_collection = pymad8.Sim.Track_Collection(14)
AddTrack(x, xp, y, yp, z, DE)
Add a new track to the collection with given position, angle and energy difference from the reference one
GenerateNtracks(nb, paramdict)
Add multiple track for which all parameters follow a gaussian disstribution
WriteBdsimTrack(outputfile)
Write the track collection in a file that can be read by BDSIM
WriteMad8Track(outputfile)
Write the track collection in a file that can be read by Mad8
class pymad8.Sim.Tracking(twiss, rmat)

Bases: object

Class to run tracking for a given track collection using Mad8 Rmat files
Requires the Mad8 twiss and rmat files
>>> tracking = pymad8.Sim.Tracking(‘twiss_tape’, ‘rmat_tape’)

Include a set of plotting function for tracjectory, phase space, …
One can also load tracking data generated by BDSIM in order to make comparison plots
AddAllElementsAsSamplers()
Register all lattice elements as samplers for the tracking
AddSamplers(value, select='index')
Add one or multiple samplers which will be used when runing the tracking
By default taking the indices of the elements we want as a samplers but one can change the select parameter to find samplers with names or types
GenerateSamplers(nb)
Add multiple samplers spread evenly along the lattice
LoadBdsimTrack(inputfilename)
Load Bdsim root file and generate orbit files for each particle in the track collection
Then store the data in a structure similar to the pymad8 generated one
LoadMad8Track(inputfilename)
Load the Mad8 track files and store it in a structure similar to the pymad8 generated data
/!/ NOT WORKING /!/
PlotCorrelation(S, coord, ref_S, ref_coord, linFit=False, partlimit=200, bdsimCompare=False)
Correlation plot for a given coordinate at the closest sampler from given S
By default the reference sampler is LUXE IP
PlotHist(S, coord, bins=50, calcSigma=False, bdsimCompare=False)
Plot histogram for a given coordinate at the closest sampler from given S
PlotPhaseSpace(S, coord, linFit=False, partlimit=200, bdsimCompare=False)
Phase space plot for a given coordinate at the closest sampler from given S
PlotRelatTrajSTD(coord)
STD Trajectory difference plot between Mad8 and BDSIM for a given coordinate
PlotRelatTrajectory(particle, coord)
Trajectory difference plot between Mad8 and BDSIM for a given coordinate and for a given particle number
PlotTrajectory(particle, coord, plotLegend=True, bdsimCompare=False)
Trajectory plot for a given coordinate and for a given particle number
RunPymad8Tracking(track_collection, turns=1)
Run tracking for all particles using the rmat from Mad8
>>> track.RunPymad8Tracking(track_collection)
getPlaneCoord(coord)
getTheoryAndFit(coord, ref_S, ref_coord, initial_fit)
class pymad8.Sim.Tracking_data(dataframe)

Bases: object

Class that stores tracking data from Mad8 or from BDSIM
Comes with a bunch of useful function to extrat data
CalcCorrelChi2(S, coord, ref_S, ref_coord, partlimit=200)

Calculate the chi2, slope and intercept of the correlation fit between two positions and two coordinates

CalcResolution(ref_coord, ref_S, BPM_list, noise=1e-05)
SVD(M, ref_Vect)

Return the correlation coefficients from a given matrix M using a Singular Value Decomposition method

SortCoeff(ref_coord, ref_S, BPM_list, noise=1e-05)
buildBPMmatrix(ref_S, ref_coord, BPM_list=None, BPM_list_type='pos', s_range=[-inf, inf], noise=None, mean_sub=False)

Build the BPM matrix M with respect to a given list of BPMs

getChi2(coord, ref_S, ref_coord)

Return the chi2 vector along the lattice

getColumnsByKeys(keylist)

Return Columns that ar in keylist

getFirstSamplerByNearestS(s)

Return the first sampler name that correspond to the closest s

getRowsByNearestS(s)

Return Rows of the closest element in the beamline

getRowsByParticles(particle)

Return Rows of a given partcile

getRowsByS(s)

Return Rows of elements with the correspondig s position

getRowsBySamplerAndNearestS(sampler, s)

Return Rows of elements with the corresponding sampler name and closest s position

getRowsBySamplerAndS(sampler, s)

Return Rows of elements with the corresponding sampler name and s position

getRowsBySamplers(sampler)

Return Rows of elements with the correspondig sampler name

getRowsByValues(key=None, minValue=-inf, maxValue=inf, equalValues=None)

Return Rows of elements that have certain values for a chosen column

getSByNearestS(s)

Return the s position of the closest element in the beamline

getSBySamplers(sampler)

Retrun the list of s positions that correspond to a given sampler name

getSamplersByNearestS(s)

Return the list of sampler name that correspond to the closest s

getSamplersByS(s)

Return the list of sampler name that correspond to a given s

getVectsByNearestS(keys, s)

Return vector of a given column and for the closest s

getVectsByParticle(keys, particle)

Return vector of a given column and for a given particle

getVectsByS(keys, s)

Return vector of a given column and for a given s

getVectsBySampler(keys, sampler)

Return vector of a given column and for a given sampler name

getVectsBySamplerAndNearestS(keys, sampler, s)

Return vector of a given column for a given sampler name and for the closest s

getVectsBySamplerAndS(keys, sampler, s)

Return vector of a given column for a given sampler name and for a given s

sMax()

Return maximal S value

sMin()

Return minimal S value

pymad8.Sim.setSamplersAndTrack(twissfile, rmatfile, nb_sampl)
pymad8.Sim.setTrackCollection(nb_part, energy, paramdict)
pymad8.Sim.testTrack(twissfile, rmatfile)

pymad8.Track module

pymad8.Track.MakeObserveFile(elementlist, filename)
pymad8.Track.MakeTableArchiveFile(elementlist, filename)
pymad8.Track.MakeTableMapFile(observeElements, observeIndex, filename)
pymad8.Track.MakeTrackCallingFile(fileNameStub)
pymad8.Track.MakeTrackFiles(savelineFileName, line, outputFileNameStub)

pymad8.Visualisation module

pymad8.Visualisation.MakeCombinedSurveyPlot(name, QUAD=True, RBEN=True, SBEN=True, MONI=True, MARK=True)

Takes a list of Survey filenames, plots them all on the same 2D plot. For branching machines or segmented models. Elements selectable via booleans, default to true

class pymad8.Visualisation.OneDim(survey, debug)

Bases: object

plot(colour=True)
class pymad8.Visualisation.TwoDim(survey, debug=False, annotate=False, fancy=False)

Bases: object

plot(event=None)
plotUpdate(event)
pymad8.Visualisation.testOneDim()
pymad8.Visualisation.testTwoDim()
pymad8.Visualisation.transformedPoly(xy, xyc, theta)
pymad8.Visualisation.transformedRect(xyc, dx, dy, theta)