.. _sect-optics: ******************************************* Optical elements that make up an instrument ******************************************* A simulation is defined by creating objects for all optical elements such as mirrors and gratings and sorting these optical elements in the order in which a photon encounters them. Marxs works under the assumption that the order of elements is fixed. A photon does not have to interact with each element (e.g. it can miss a grating, but still be detected by another focal plane instrument), but under no circumstance will a photon ever go back to an element it already passed. There are ways to describe layouts where photons split between different paths, such as in the design of XMM-Newton where some photons pass the gratings and get reflected into the RGS detectors, while others continue on their path without ever encountering a grating until they hit the imaging detectors. Such complex setups are described in :ref:`complexgeometry`. .. _args for optical elements: Define an optical element for the simulation ============================================ All optical elements have to be derived from :obj:`marxs.optics.base.OpticalElement` and they all share a common interface. Optical elements are generated with the following keywords: - Optical elements are placed in 3-d space with either the ``pos4d`` keyword, or the ``zoom``, ``orientation`` and ``position`` keywords as explained in :ref:`coordsys`: >>> from marxs.optics import FlatDetector >>> det1 = FlatDetector(position=[5., 1., 0.], zoom=[1, 50.2, 50.2], pixsize=0.1) - Optical elements can have a name: >>> from marxs.optics import FlatDetector >>> det1 = FlatDetector(name='back-illuminated CCD No 1') - Individual optical elements might add more keywords (e.g. the pixel size for a detector or the filename of a parameter file); those are listed in the documentation for individual optical elements. Here is an example for a detector: >>> from marxs.optics import FlatDetector >>> det1 = FlatDetector(pixsize=0.04) .. _sect-optics-output: Output columns for the photon list ================================== Each MARXS run starts with a source and a pointing. The source generates the physical properties of each photon (energy, polarization) and places it in the universe (RA, dec, time). The pointing transforms position and polarization angle into vectors in spacecraft coordiantes (these vectors are called "dir" and "polarization") - see :ref:`sect-results-output`. Next, we need to assign an initial position in spacecraft coordinates to each photon. This is done by an aperture. One common approach in ray-tracing is to subdivide the entrace aperture into a regular grid and trace one ray per grid cell; alternatively, a large number of photons can be randomly distributed on the aperture. MARXS takes this second approach because (i) the result does not only depend on the position in the entrance aperture, but on a number of other parmeters (energy, polarization, time for time variable sources or pointings) and (ii) the photon list produced by the simulation closely resembles the photons lists from real observations. The column "pos" in the table holds the photon position as assigned on the aperture in homogeneous coordinates. After passing through the aperture, the photon list will have at least the following columns: pos Position where the photon last interacted with an optical element in homogeneous coordinates. dir Direction of the phootn ray in homogeneous coordiantes. polarization Polarization vector in homogeneous coordinates (can be complex for cicularly polarized light). energy in keV probability Probability that the photon is not absorbed anywhere until this point Optical elements may add any number of diagnostic columns in their code. The values of these columns can change with every interaction. All other columns are just passed along. They carry information that helps interpret the results of a ray-trace, e.g. the (RA, Dec) sky coordinate where the photons was originally emitted. Most optical elements in MARXS support at least two ways to add columns: Which optical element did a photon pass? Instruments often have several identical CCDs or several gratings. If an element has a property ``id_col`` set to a string and ``id_num`` set to a number, then all photons passing it will have ``id_num`` set in the column ``id_col``. Photons that do not pass through any element of this collection have ``photons[id_col]==-1``. If a photon passes through multiple elements with the same ``id_col`` the ``id_num`` recorded is the last element passed. Where did a photon hit an element? Flat optical elements (those inheriting from `marxs.optics.base.FlatOpticalElement`) have an attribute called ``loc_coos_name``. Setting this to a list of two strings, will cause the element to add two columns with those names that hold the local coordiantes (normalized to the range -1 to 1 for both coordinates). If the column names already exisit, values will be overwritten. The default are the generic names "x" and "y". The following example demonstrates these two ways of adding columns with diangostic information where exactly the photon passed which element. First, we generate a bunch of photons all parallel to the x-axis. >>> import numpy as np >>> from marxs.utils import generate_test_photons >>> photons = generate_test_photons(5) >>> photons['pos'].data[:, 1] = np.arange(5) Then, we add two CCD detectors:: >>> det1 = FlatDetector(pixsize=0.1) >>> det1.loc_coos_name ['det_x', 'det_y'] >>> det1.id_col = 'CCD_ID' >>> det1.id_num = 1 >>> det2 = FlatDetector(pixsize=0.2, position=[0, 2.5, 0.]) >>> det2.id_col = 'CCD_ID' >>> det2.id_num = 2 Last, we process the photons through both detectors:: >>> photons = det1(photons) >>> photons = det2(photons) >>> photons # doctest: +IGNORE_OUTPUT