.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "auto_examples/core/InstrumentDescription.py" .. LINE NUMBERS ARE GIVEN BELOW. .. only:: html .. note:: :class: sphx-glr-download-link-note :ref:`Go to the end ` to download the full example code. .. rst-class:: sphx-glr-example-title .. _sphx_glr_auto_examples_core_InstrumentDescription.py: Working with Instrumental Descriptions ====================================== the instrumental description is loaded by the event source, and consists of a hierarchy of classes in the ctapipe.instrument module, the base of which is the ``SubarrayDescription`` .. GENERATED FROM PYTHON SOURCE LINES 10-23 .. code-block:: Python :lineno-start: 11 from astropy.coordinates import SkyCoord from ctapipe.io import EventSource from ctapipe.utils.datasets import get_dataset_path from ctapipe.visualization import CameraDisplay filename = get_dataset_path("gamma_prod5.simtel.zst") with EventSource(filename, max_events=1) as source: subarray = source.subarray .. GENERATED FROM PYTHON SOURCE LINES 24-27 the SubarrayDescription: ------------------------ .. GENERATED FROM PYTHON SOURCE LINES 27-30 .. code-block:: Python :lineno-start: 28 subarray.info() .. rst-class:: sphx-glr-script-out .. code-block:: none Subarray : MonteCarloArray Num Tels : 180 Footprint: 4.92 km2 Height : 2147.00 m Lon/Lat : 0.0 deg, 0.0 deg Type Count Tel IDs ----------------- ----- --------------- SST_ASTRI_CHEC 120 30-99,131-180 LST_LST_LSTCam 4 1-4 MST_MST_NectarCam 28 100-124,128-130 MST_MST_FlashCam 28 5-29,125-127 .. GENERATED FROM PYTHON SOURCE LINES 31-34 .. code-block:: Python :lineno-start: 31 subarray.to_table() .. raw:: html
Table length=180
tel_idnametypepos_xpos_ypos_zcamera_nameoptics_namecamera_indexoptics_indextel_description
mmm
uint16str5str3float32float32float32str9str5int64int64str17
1LSTLST-20.643-64.81734.3LSTCamLST21LST_LST_LSTCam
2LSTLST79.993996-0.76829.4LSTCamLST21LST_LST_LSTCam
3LSTLST-19.39665.231.0LSTCamLST21LST_LST_LSTCam
4LSTLST-120.0331.15133.1LSTCamLST21LST_LST_LSTCam
5MSTMST-0.017-0.00124.35FlashCamMST12MST_MST_FlashCam
6MSTMST-1.468-151.22131.0FlashCamMST12MST_MST_FlashCam
7MSTMST-3.1379998-325.24539.0FlashCamMST12MST_MST_FlashCam
8MSTMST1.4339999151.2225.0FlashCamMST12MST_MST_FlashCam
9MSTMST3.1039999325.24323.5FlashCamMST12MST_MST_FlashCam
.................................
172ASTRISST260.0-920.065.0CHECASTRI00SST_ASTRI_CHEC
173ASTRISST-500.0815.015.0CHECASTRI00SST_ASTRI_CHEC
174ASTRISST-500.0-815.075.0CHECASTRI00SST_ASTRI_CHEC
175ASTRISST500.0815.045.0CHECASTRI00SST_ASTRI_CHEC
176ASTRISST500.0-815.053.0CHECASTRI00SST_ASTRI_CHEC
177ASTRISST-810.0655.012.0CHECASTRI00SST_ASTRI_CHEC
178ASTRISST-810.0-655.068.0CHECASTRI00SST_ASTRI_CHEC
179ASTRISST810.0655.020.0CHECASTRI00SST_ASTRI_CHEC
180ASTRISST810.0-655.041.0CHECASTRI00SST_ASTRI_CHEC


.. GENERATED FROM PYTHON SOURCE LINES 35-39 You can also get a table of just the ``OpticsDescriptions`` (``CameraGeometry`` is more complex and can’t be stored on a single table row, so each one can be converted to a table separately) .. GENERATED FROM PYTHON SOURCE LINES 39-43 .. code-block:: Python :lineno-start: 40 subarray.to_table(kind="optics") .. raw:: html
Table length=3
optics_namesize_typereflector_shapemirror_arean_mirrorsn_mirror_tilesequivalent_focal_lengtheffective_focal_length
m2mm
str5str3str20float64int64int64float64float64
ASTRISSTSCHWARZSCHILD_COUDER14.126235008239746222.15000009536743162.1519100666046143
LSTLSTPARABOLIC386.7332458496094119828.029.30565071105957
MSTMSTHYBRID106.241355895996118616.016.445049285888672


.. GENERATED FROM PYTHON SOURCE LINES 44-46 Make a sub-array with only SC-type telescopes: .. GENERATED FROM PYTHON SOURCE LINES 46-52 .. code-block:: Python :lineno-start: 47 sc_tels = [tel_id for tel_id, tel in subarray.tel.items() if tel.optics.n_mirrors == 2] newsub = subarray.select_subarray(sc_tels, name="SCTels") newsub.info() .. rst-class:: sphx-glr-script-out .. code-block:: none Subarray : SCTels Num Tels : 120 Footprint: 4.92 km2 Height : 2147.00 m Lon/Lat : 0.0 deg, 0.0 deg Type Count Tel IDs -------------- ----- ------------- SST_ASTRI_CHEC 120 30-99,131-180 .. GENERATED FROM PYTHON SOURCE LINES 53-55 can also do this by using ``Table.group_by`` .. GENERATED FROM PYTHON SOURCE LINES 58-61 Explore some of the details of the telescopes --------------------------------------------- .. GENERATED FROM PYTHON SOURCE LINES 61-65 .. code-block:: Python :lineno-start: 62 tel = subarray.tel[1] tel .. rst-class:: sphx-glr-script-out .. code-block:: none TelescopeDescription(type='LST', optics_name='LST', camera_name='LSTCam') .. GENERATED FROM PYTHON SOURCE LINES 66-68 .. code-block:: Python :lineno-start: 66 tel.optics.mirror_area .. rst-class:: sphx-glr-script-out .. code-block:: none .. GENERATED FROM PYTHON SOURCE LINES 69-71 .. code-block:: Python :lineno-start: 69 tel.optics.n_mirror_tiles .. rst-class:: sphx-glr-script-out .. code-block:: none 198 .. GENERATED FROM PYTHON SOURCE LINES 72-74 .. code-block:: Python :lineno-start: 72 tel.optics.equivalent_focal_length .. rst-class:: sphx-glr-script-out .. code-block:: none .. GENERATED FROM PYTHON SOURCE LINES 75-77 .. code-block:: Python :lineno-start: 75 tel.camera .. rst-class:: sphx-glr-script-out .. code-block:: none CameraDescription(name=LSTCam, geometry=LSTCam, readout=LSTCam) .. GENERATED FROM PYTHON SOURCE LINES 78-80 .. code-block:: Python :lineno-start: 78 tel.camera.geometry.pix_x .. rst-class:: sphx-glr-script-out .. code-block:: none .. GENERATED FROM PYTHON SOURCE LINES 81-82 %matplotlib inline .. GENERATED FROM PYTHON SOURCE LINES 82-85 .. code-block:: Python :lineno-start: 83 CameraDisplay(tel.camera.geometry) .. image-sg:: /auto_examples/core/images/sphx_glr_InstrumentDescription_001.png :alt: LSTCam :srcset: /auto_examples/core/images/sphx_glr_InstrumentDescription_001.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none .. GENERATED FROM PYTHON SOURCE LINES 86-89 .. code-block:: Python :lineno-start: 86 CameraDisplay(subarray.tel[98].camera.geometry) .. image-sg:: /auto_examples/core/images/sphx_glr_InstrumentDescription_002.png :alt: CHEC :srcset: /auto_examples/core/images/sphx_glr_InstrumentDescription_002.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none .. GENERATED FROM PYTHON SOURCE LINES 90-100 Plot the subarray ----------------- We’ll make a subarray by telescope type and plot each separately, so they appear in different colors. We also calculate the radius using the mirror area (and exaggerate it a bit). This is just for debugging and info, for any “real” use, a ``visualization.ArrayDisplay`` should be used .. GENERATED FROM PYTHON SOURCE LINES 100-103 .. code-block:: Python :lineno-start: 101 subarray.peek() .. image-sg:: /auto_examples/core/images/sphx_glr_InstrumentDescription_003.png :alt: MonteCarloArray :srcset: /auto_examples/core/images/sphx_glr_InstrumentDescription_003.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none .. GENERATED FROM PYTHON SOURCE LINES 104-107 .. code-block:: Python :lineno-start: 104 subarray.footprint .. rst-class:: sphx-glr-script-out .. code-block:: none .. GENERATED FROM PYTHON SOURCE LINES 108-111 Get info about the subarray in general -------------------------------------- .. GENERATED FROM PYTHON SOURCE LINES 111-114 .. code-block:: Python :lineno-start: 112 subarray.telescope_types .. rst-class:: sphx-glr-script-out .. code-block:: none (TelescopeDescription(type='SST', optics_name='ASTRI', camera_name='CHEC'), TelescopeDescription(type='LST', optics_name='LST', camera_name='LSTCam'), TelescopeDescription(type='MST', optics_name='MST', camera_name='NectarCam'), TelescopeDescription(type='MST', optics_name='MST', camera_name='FlashCam')) .. GENERATED FROM PYTHON SOURCE LINES 115-117 .. code-block:: Python :lineno-start: 115 subarray.camera_types .. rst-class:: sphx-glr-script-out .. code-block:: none (CameraDescription(name=CHEC, geometry=CHEC, readout=CHEC), CameraDescription(name=FlashCam, geometry=FlashCam, readout=FlashCam), CameraDescription(name=LSTCam, geometry=LSTCam, readout=LSTCam), CameraDescription(name=NectarCam, geometry=NectarCam, readout=NectarCam)) .. GENERATED FROM PYTHON SOURCE LINES 118-121 .. code-block:: Python :lineno-start: 118 subarray.optics_types .. rst-class:: sphx-glr-script-out .. code-block:: none (OpticsDescription(name=ASTRI, size_type=SST, reflector_shape=SCHWARZSCHILD_COUDER, equivalent_focal_length=2.15 m, effective_focal_length=2.15 m, n_mirrors=2, mirror_area=14.13 m2), OpticsDescription(name=LST, size_type=LST, reflector_shape=PARABOLIC, equivalent_focal_length=28.00 m, effective_focal_length=29.31 m, n_mirrors=1, mirror_area=386.73 m2), OpticsDescription(name=MST, size_type=MST, reflector_shape=HYBRID, equivalent_focal_length=16.00 m, effective_focal_length=16.45 m, n_mirrors=1, mirror_area=106.24 m2)) .. GENERATED FROM PYTHON SOURCE LINES 122-127 .. code-block:: Python :lineno-start: 122 center = SkyCoord("10.0 m", "2.0 m", "0.0 m", frame="groundframe") coords = subarray.tel_coords # a flat list of coordinates by tel_index coords.separation(center) .. rst-class:: sphx-glr-script-out .. code-block:: none /home/docs/checkouts/readthedocs.org/user_builds/ctapipe/envs/2868/lib/python3.10/site-packages/astropy/coordinates/baseframe.py:1784: NonRotationTransformationWarning: transforming other coordinates from to . Angular separation can depend on the direction of the transformation. warnings.warn(NonRotationTransformationWarning(self, other_frame)) .. GENERATED FROM PYTHON SOURCE LINES 128-139 Telescope IDs vs Indices ------------------------ Note that ``subarray.tel`` is a dict mapped by ``tel_id`` (the identifying number of a telescope). It is possible to have telescope IDs that do not start at 0, are not contiguouous (e.g. if a subarray is selected). Some functions and properties like ``tel_coords`` are numpy arrays (not dicts) so they are not mapped to the telescope ID, but rather the *index* within this SubarrayDescription. To convert between the two concepts you can do: .. GENERATED FROM PYTHON SOURCE LINES 139-143 .. code-block:: Python :lineno-start: 140 subarray.tel_ids_to_indices([1, 5, 23]) .. rst-class:: sphx-glr-script-out .. code-block:: none array([ 0, 4, 22]) .. GENERATED FROM PYTHON SOURCE LINES 144-146 or you can get the indexing array directly in numpy or dict form: .. GENERATED FROM PYTHON SOURCE LINES 146-149 .. code-block:: Python :lineno-start: 147 subarray.tel_index_array .. rst-class:: sphx-glr-script-out .. code-block:: none array([-9223372036854775808, 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, 40, 41, 42, 43, 44, 45, 46, 47, 48, 49, 50, 51, 52, 53, 54, 55, 56, 57, 58, 59, 60, 61, 62, 63, 64, 65, 66, 67, 68, 69, 70, 71, 72, 73, 74, 75, 76, 77, 78, 79, 80, 81, 82, 83, 84, 85, 86, 87, 88, 89, 90, 91, 92, 93, 94, 95, 96, 97, 98, 99, 100, 101, 102, 103, 104, 105, 106, 107, 108, 109, 110, 111, 112, 113, 114, 115, 116, 117, 118, 119, 120, 121, 122, 123, 124, 125, 126, 127, 128, 129, 130, 131, 132, 133, 134, 135, 136, 137, 138, 139, 140, 141, 142, 143, 144, 145, 146, 147, 148, 149, 150, 151, 152, 153, 154, 155, 156, 157, 158, 159, 160, 161, 162, 163, 164, 165, 166, 167, 168, 169, 170, 171, 172, 173, 174, 175, 176, 177, 178, 179]) .. GENERATED FROM PYTHON SOURCE LINES 150-152 .. code-block:: Python :lineno-start: 150 subarray.tel_index_array[[1, 5, 23]] .. rst-class:: sphx-glr-script-out .. code-block:: none array([ 0, 4, 22]) .. GENERATED FROM PYTHON SOURCE LINES 153-160 .. code-block:: Python :lineno-start: 153 subarray.tel_indices[ 1 ] # this is a dict of tel_id -> tel_index, so we can only do one at once ids = subarray.get_tel_ids_for_type(subarray.telescope_types[0]) ids .. rst-class:: sphx-glr-script-out .. code-block:: none (30, 31, 32, 33, 34, 35, 36, 37, 38, 39, 40, 41, 42, 43, 44, 45, 46, 47, 48, 49, 50, 51, 52, 53, 54, 55, 56, 57, 58, 59, 60, 61, 62, 63, 64, 65, 66, 67, 68, 69, 70, 71, 72, 73, 74, 75, 76, 77, 78, 79, 80, 81, 82, 83, 84, 85, 86, 87, 88, 89, 90, 91, 92, 93, 94, 95, 96, 97, 98, 99, 131, 132, 133, 134, 135, 136, 137, 138, 139, 140, 141, 142, 143, 144, 145, 146, 147, 148, 149, 150, 151, 152, 153, 154, 155, 156, 157, 158, 159, 160, 161, 162, 163, 164, 165, 166, 167, 168, 169, 170, 171, 172, 173, 174, 175, 176, 177, 178, 179, 180) .. GENERATED FROM PYTHON SOURCE LINES 161-164 .. code-block:: Python :lineno-start: 161 idx = subarray.tel_ids_to_indices(ids) idx .. rst-class:: sphx-glr-script-out .. code-block:: none array([ 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, 40, 41, 42, 43, 44, 45, 46, 47, 48, 49, 50, 51, 52, 53, 54, 55, 56, 57, 58, 59, 60, 61, 62, 63, 64, 65, 66, 67, 68, 69, 70, 71, 72, 73, 74, 75, 76, 77, 78, 79, 80, 81, 82, 83, 84, 85, 86, 87, 88, 89, 90, 91, 92, 93, 94, 95, 96, 97, 98, 130, 131, 132, 133, 134, 135, 136, 137, 138, 139, 140, 141, 142, 143, 144, 145, 146, 147, 148, 149, 150, 151, 152, 153, 154, 155, 156, 157, 158, 159, 160, 161, 162, 163, 164, 165, 166, 167, 168, 169, 170, 171, 172, 173, 174, 175, 176, 177, 178, 179]) .. GENERATED FROM PYTHON SOURCE LINES 165-168 .. code-block:: Python :lineno-start: 165 subarray.tel_coords[idx] .. rst-class:: sphx-glr-script-out .. code-block:: none .. GENERATED FROM PYTHON SOURCE LINES 169-173 so, with that method you can quickly get many telescope positions at once (the alternative is to use the dict ``positions`` which maps ``tel_id`` to a position on the ground .. GENERATED FROM PYTHON SOURCE LINES 173-175 .. code-block:: Python :lineno-start: 174 subarray.positions[1] .. rst-class:: sphx-glr-script-out .. code-block:: none .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 1.830 seconds) .. _sphx_glr_download_auto_examples_core_InstrumentDescription.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: InstrumentDescription.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: InstrumentDescription.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: InstrumentDescription.zip ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_