Field Maps

The pybdsim.Field module has several functions and classes to help create, convert, and visualise BDSIM-format field maps easier. Various useful workflows are described below, but every function and class is documented in: pybdsim.Field module.

Loading

Existing BDSIM-format field maps can be loaded with:

fm = pybdsim.Field.Load("nameOfFieldMap.dat")

This will automatically detect the number of dimensions in the field map and load it into a numpy array inside a pybdsim class. The class is one of:

  • pybdsim.Field.Field1D

  • pybdsim.Field.Field2D

  • pybdsim.Field.Field3D

  • pybdsim.Field.Field4D

These all have member variables such as:

  • data - the numpy array of the field map data

  • columns - a list of column names

  • header - a dictionary of all information in the header

Writing

All of the pybdsim classes described in Loading (e.g. Field3D) have a method called Write. This will write out the instance of that class (i.e. that particular field map data) to a BDSIM-format field map file.

Therefore, these classes represent a a very good route for creation and conversion of field maps. The recommended strategy is to get the information into one of those classes then it can be written out to file, loaded, plotted with ease.

Creation

A field map can be created from scratch in Python from knowledge of a field, such as an equation. A numpy array of the right shape should be created, then an instance of pybdsim.Field.Field[N]D created (where N is one of (1,2,3,4).

Creating the array is entirely up to the user as this depends on how they get the information. However, a small example script is shown below here (in Python language):

import numpy as np
import pybdsim

fx, fy, fz = 0.0, 1.0, 0.0
xmax = 30  # cm is the default units for bdsim field maps
ymax = 30
data = []
# loop over and build up 3d lists of lists of lists
for xi in [-xmax, xmax]:
   v = []
   for yi in [-ymax, ymax]:
       # here, fx,fy,fz could be from a function
       v.append([xi, yi, fx, fy, fz])
       data.append(v)

# convert to numpy array
data = np.array(data)

# construct a BDSIM format field object and write it out
f = pybdsim.Field.Field2D(data)

# optionally add some strings that will appear as comments at the
# start of the file
f.AddComment("I = 200A")
f.AddComment("Generated by fictional software v1.2.3")
f.Write("box-field-map.dat')

Note

BDSIM’s default convention for a field map is to loop over the lowest dimension first, i.e. x. Therefore, the loop above is with x on the outside. See below for purposely handling the other order of looping.

This minimal example creates a 2D ‘box’ of 4 points in space each with the same field value of [0,1,0] (unit Y direction with magnitude 1). The box is \(\pm 30 cm\) in size. Units aren’t used explicitly - just numbers - but the units of BDSIM field map files are cm and s.

Note

The file written out by the class at the end can be loaded again with pybdsim.Field.Load to get a (new) object with the same contents as the original Field2D instance.

The array should be structured to contain the coordinates and field components at each point. So for a 2D field map, this would be: \((x, y, F_x, F_y, F_z)\). The array would then have 2 dimensions representing its structure. So if the field map had 10 points in x and 20 in y, it would ultimately be a numpy array with the shape (10,20,5).

Generally, the numpy.shape(array) would look like:

  • 1D: (N_x, 4)

  • 2D: (N_x, N_y, 5)

  • 3D: (N_x, N_y, N_z, 6)

  • 4D: (N_x, N_y, N_z, N_t, 7)

Warning

A key check is looking at the field map, the higher dimension coordinates (e.g. Y, not X) should change first. So for a given X value we should see the Y values cycle through a range, then the X should increment then the Y values cycle again. If this is not the case, then the loop order of dimensions is backwards. You can use “loopOrder” in the header or rewrite the field map correctly.

Alternative Dimensions

In the case of alternative dimension (e.g. a 2D field map with x and z dimensions but no y), the construction is the same but we can label the dimensions differently. The dimensions must be in order (e.g. x, y, z, then t for whichever ones are used).

Example:

fm = pybdsim.Field.Field2D(arrayData, firstColumn='X', secondColumn='Z')

Alternative Loop Order

It is possible for BDSIM to read a file where the right-most coordinate column varies first. However, for each value, the coordinate columns must still be in x,y,z,t order left to right. Below is an example similar to above but writing out the file the other way (note the write function). This will also write the line loopOrder> tzyx in the header so BDSIM can load the field map equivalently.

import numpy as np
import pybdsim

fx, fy, fz = 0.0, 1.0, 0.0
xmax = 30  # cm is the default units for bdsim field maps
ymax = 30
data = []
# loop over and build up 3d lists of lists of lists
for yi in [-ymax, ymax]:
   v = []
   for xi in [-xmax, xmax]:
       # here, fx,fy,fz could be from a function
       # note, xi and yi coordinates must still be in that order for the value in the array
       v.append([xi, yi, fx, fy, fz])
       data.append(v)

# convert to numpy array
data = np.array(data)

# construct a BDSIM format field object and write it out
# this will be written out in the BDSIM conventional looping order
f = pybdsim.Field.Field2D(data, flip=True)
f.Write("box-field-map.dat')

Below is a script included with BDSIM (bdsim/examples/features/maps_bdsim/Generate2DLoopOrder.py) that shows 4 ways to write a field map with the same information. Ultimately, they convey the exact same field map to BDSIM although the file contents differ (2 sets of possible contents).

import numpy as _np
import pybdsim

B = 2.0

# LOOP METHOD 1
data = []
# loop over and build up 3d lists of lists of lists
for x in [-1,0,1]:
    v = []
    for z in [3,4]:
        v.append([x, z, B*x, B*x*z, B*z])
    data.append(v)

# convert to numpy array
data = _np.array(data)

# we looped in x first as per bdsim, so we need only tell it that
# the 2nd column is Z and not Y
f = pybdsim.Field.Field2D(data, secondColumn='Z')
f.Write('2dexample_loopOrder_for_xz.dat')
# but we can purposively write it out the other loop way for testing purposes
# note the header keys are still the same apart from loopOrder> tzyx
f.Write('2dexample_loopOrder_for_xz_tzyx.dat', writeLoopOrderReversed=True)


# LOOP METHOD 2
data2 = []
# loop over other way - outer dimension first
# this isn't the bdsim way, but we may get a field map from some other source that's
# structured like this - so even if you're not creating it in a loop, it may have this
# structure already.
for z in [3,4]:
    v = []
    for x in [-1,0,1]:
        v.append([x, z, B*x, B*x*z, B*z]) # values must still be in xyzt order
    data2.append(v)

# convert to numpy array
data2 = _np.array(data2)

# array structure is z is outer dimension, then x - we need it the other way
# around, so we use flip=True when constructing the field instance
g = pybdsim.Field.Field2D(data2, flip=True, secondColumn='Z')
# this will write out a file identical to the very first one
g.Write('2dexample_loopOrder_for_zx.dat')
# this will write out a file identical to the second one
g.Write('2dexample_loopOrder_for_zx_tzyx.dat', writeLoopOrderReversed=True)

Visualisation and Plotting

To visualise a field map, it is possible to do so in BDSIM / Geant4. See the BDSIM manual for this information. This draws a selection of arrows in the 3D model and gives a rough indication that the field map is as intended.

An alternative way is to load the data in pybdsim in Python and plot it, either fully or in slices (for 3D or 4D maps).

Any library desired can be used in Python and the classes described above in Loading provide an excellent way to get a numpy array, that is ubiquitous in Python programming and libraries.

pybdsim provides a variety of small plotting functions mostly for 1D and 2D field maps using matplotlib. These functions are inside the pybdsim.Field module and all start with Plot. A list is:

  • pybdsim.Field.Plot1DFxFyFz

  • pybdsim.Field.Plot2D

  • pybdsim.Field.Plot2DXY

  • pybdsim.Field.Plot2DXYMagnitude

  • pybdsim.Field.Plot2DXYConnectionOrder

  • pybdsim.Field.Plot2DXYStream

  • pybdsim.Field.Plot2DXYComponent

  • pybdsim.Field.Plot2DXYFxFyFz

  • pybdsim.Field.Plot2DXYBx

  • pybdsim.Field.Plot2DXYBy

  • pybdsim.Field.Plot2DXYBz

  • pybdsim.Field.Plot3DXY

  • pybdsim.Field.Plot3DXZ

Warning

Plots that use arrows or stream plots do not depend on the order of the points so they cannot be relied upon to tell if the field map being prepared is in the correct order. Use Plot2DXYMagnitude or Plot2DXYConnectionOrder to verify the order of the points.

A (guaranteed) complete list can be found in pybdsim.Field module.

Each can be inspected (in IPython, which is recommended) with a question mark to see its description:

>>> import pybdsim
>>> pybdsim.Field.Plot2DXY?
Signature: pybdsim.Field.Plot2DXY(filename, scale=None)
Docstring:
Plot a bdsim field map file using the X,Y plane.

:param filename: name of field map file or object
:type filename: str, pybdsim.Field._Field.Field2D instance
:param scale: numerical scaling for quiver plot arrow lengths.
:type scale: float
>>>

Conversion

your data -> numpy array -> pybdsim.Field.FieldND(data) class

To convert a field map, you should first write a loader from your own format to the field map into a numpy array with a structure described in Creation. Then, this array can be wrapped in an instance of one of the pybdsim Field classes. This class can then be used to write out the field map in BDSIM’s format. This would look something like:

def LoadMyFormatFieldMap(filename):
    # ... some implementation...
    # assume variable 'data' of type numpy.array
    return data

def Convert(inputfilename, outputfilename):
    d = LoadMyFormatFieldMap(inputfilename)
    # assume here it's a 2D field map... need to know which class to use
    bd = pybdsim.Field2D(d)
    bd.Write(outputfilename)

The source of a field map data should represent an equally spaced grid of points and provided in order, such that it can be converted easily to BDSIM’s format with the various classes.

Note

See Importance of Order. Make sure to validate the order with plots before using in a simulation in BDSIM. You can also visualise the fields in BDSIM to check. It is recommended to do this with a single component before using in a bigger model.

Sorting Points

If the data points for the field map correspond to a rectilinear grid but are not provided in order (sometimes can happen from finite-element programs), you should ideally try to get a field map in order. Failing that, you can try to sort the data into an ordered array. An example implementation is given in pybdsim.Field.SortUnorderedFieldMap2D. Although, this is provided there is no guarantee the implementation will work depending on the numerical precision of the coordinates. It is still recommended to go back to the origin field program and get a correct grid of points.

Importance of Order

In a BDSIM field map file, the coordinates at each point are written but BDSIM itself does not use these. BDSIM reads the header information and loops over the data assuming the number of points specified in the header. Therefore if the data is provided in the wrong order the field map will appear scrambled in BDSIM. This can happen with hand-preparation and editing of files.

It is recommended to use the pybdsim.Field classes as these are guaranteed to write the data out correctly.

Plots that use arrows or stream plots do not depend on the order of the points so they cannot be relied upon to tell if the field map being prepared is in the correct order. Use Plot2DXYMagnitude or Plot2DXYConnectionOrder to verify the order of the points.