Rowland Grating Spectrometer

In the previous sections, we either used pre-defined instrument configurations or placed all components “by hand”, giving their absolute coordinates in mm. But how do we know the best places for gratings, detectors etc.? In some cases, analytical work can tell us where e.g. diffraction gratings and detectors need to go to achieve the best possible resolving power. The module marxs.design.rowland implements those formulas for a spectrometer in a Rowland Torus design (see e.g. Paerels 1999 and Paerels 2010 for a derivation).

In the following example, we define a Rowland torus based spectrometer. In this case, we select the radius of the Rowland Circle to be very small (1 m) and also make the grating constant of the gratings very small. This way we get a much larger diffraction angle than usual, and we end up with a model that does real ray-tracing (albeit with some unrealistic input values) and delivers very clear images to demonstrate the operating principles for teaching. We begin by importing the relevant packages, and then we define the Rowland torus. In this case, both radii of the torus are set to 500 mm, so the Rowland circle in the x, y plane touches the focal point at (0,0) and has a diameter of 1 m, so grating on the optical axis would be 1 m away from the focal point:

>>> import numpy as np
>>> import astropy.units as u
>>> from marxs import optics
>>> from marxs import simulator
>>> from marxs.design import rowland
>>> my_torus = rowland.RowlandTorus(500, 500)

Our instrument approximates an imaging optic by using a circular aperture that illuminates a marxs.optics.PerfectLens with some additional scatter (marxs.optics.RadialMirrorScatter). The reflection on the mirror and the scatter always happen at the same position, so those two functions are summarized in a single element (a marxs.optics.FlatStack) as explained in Complex designs:

>>> aperture = optics.CircleAperture(position=[1200, 0, 0], zoom=[1, 100, 100])
>>> mirror = optics.FlatStack(position=[1100, 0, 0], zoom=[20, 100, 100],
...                            elements=[optics.PerfectLens, optics.RadialMirrorScatter],
...                            keywords = [{'focallength': 1100},
...                                        {'inplanescatter': 3 * u.arcmin, 'perpplanescatter': 0.3 * u.arcmin}])

Next, behind the mirror, we place a GratingArrayStructure, which is simply a collection of gratings. The class of the grating is defined in the parameter elem_class and any keyword arguments needed to initialize every individual grating are passed in elem_args. d is the grating period, the size of the gratings set by the zoom and the order_selector that gives the probabilities to decide into which diffraction order a particular photon is diffracted. Note that position and orientation are missing from this list, because they will be calculated by the GratingArrayStructure.

The remaining parameters set how the GratingArrayStructure places the gratings. Looking along the optical axis, it will select y and z positions for the gratings such that they lay in a ring, where the inner and outer radius is given by the radius parameter. The x-position of each grating is adjusted to make sure it touches the Rowland torus. To do so, we need to pass GratingArrayStructure an instance of the Rowland Torus, the distance d_element between gratings (measured center to center - note that this is typically larger than the size of each element to account for the frame the gratings are mounted on), and the range along the optical axis (the x-axis of the torus) in which the gratings are expected to fall (the equation of the torus is solved numerically, and the solver needs bounded regions with a unique solution to work):

>>> gas = rowland.GratingArrayStructure(rowland=my_torus, d_element=[25, 25],
...                                     radius=[20, 100], elem_class=optics.FlatGrating,
...                                     elem_args={'d': 1e-5, 'zoom':[1, 10, 10],
...                                                'order_selector': optics.OrderSelector([-1, 0, 1])})

Last, we place detectors on the inner part of the Rowland circle. This is very similar to the GratingArrayStructure, except that the CCD detectors are not placed in a 2D pattern, but just in a 1D line along the Rowland Circle. All other parameters work the same way as above:

>>> det = rowland.RectangularGrid(rowland=my_torus, y_range=[-50, 50],
...                               d_element=[10.5, 10.5], elem_class=optics.FlatDetector,
...                               elem_args={'zoom': [1, 5, 5], 'pixsize': 0.01})

Last, all the elements above as combined into a single element:

>>> demo = simulator.Sequence(elements=[aperture, mirror, gas, det])

Reference / API

marxs.design.rowland Module

Tools for setting up instruments in the Rowland Torus geometry.

This includes an object that represents the Rowland torus itself (imaginatively called RowlandTorus), some helper functions used to set the right parameters for the torus and a few classes that are derived from marxs.simulator.ParallelCalculated, each placing elements such as gratings or detectors on the Rowland torus. There are many ways to do that, e.g. the limits can be defined in x,y,z coordinates or in theta, phi coordinates on the torus, gratings can be ordered in concentric circles or pack densely in a rectangular area etc. At this point, the different classes do not exploit all these possibilities, they merely give a few of many possible ways to set up an instrument.

These classes may be generalized in the future.

Functions

find_radius_of_photon_shell(photons, ...[, ...])

Find the radius the photons coming from a single mirror shell have.

design_tilted_torus(f, alpha, beta)

Design a torus with specifications similar to Heilmann et al. 2010.

double_rowland_from_channel_distance(d_BF, R, f)

This prescription is borne out of the engineering needs where it turns out that the spacing between channels d_BF is actually an important parameter so it makes sense to use that here to parameterize the torus, too.

add_offset_double_rowland_channels(geometry)

Generate Rowland tori for multiple spectral channels in a double-tori design.

Classes

RowlandTorus(R, r, **kwargs)

Torus with y axis as symmetry axis.

ElementsOnTorus(**kwargs)

A collection of elements on a Rowland torus.

GratingArrayStructure(**kwargs)

Collection of diffraction gratings on the Rowland Torus

RectangularGrid(**kwargs)

A collection of diffraction gratings on the Rowland torus.

CircularMeshGrid(**kwargs)

A collection of diffraction gratings on the Rowland torus filling a circle.

Class Inheritance Diagram

Inheritance diagram of marxs.design.rowland.RowlandTorus, marxs.design.rowland.ElementsOnTorus, marxs.design.rowland.GratingArrayStructure, marxs.design.rowland.RectangularGrid, marxs.design.rowland.CircularMeshGrid