.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "auto_examples/tutorials/coordinates_example.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_tutorials_coordinates_example.py: Coordinates usage in ctapipe ============================ .. GENERATED FROM PYTHON SOURCE LINES 6-28 .. code-block:: Python :lineno-start: 7 import copy import astropy.units as u import matplotlib.pyplot as plt from astropy.coordinates import AltAz, EarthLocation, SkyCoord from astropy.time import Time from ctapipe.coordinates import CameraFrame, NominalFrame, TelescopeFrame from ctapipe.instrument import SubarrayDescription from ctapipe.io import EventSource from ctapipe.utils import get_dataset_path from ctapipe.visualization import CameraDisplay # %matplotlib inline # make plots and fonts larger plt.rcParams["figure.figsize"] = (12, 8) plt.rcParams["font.size"] = 16 .. GENERATED FROM PYTHON SOURCE LINES 29-32 Open test dataset ----------------- .. GENERATED FROM PYTHON SOURCE LINES 32-42 .. code-block:: Python :lineno-start: 33 filename = get_dataset_path("gamma_prod5.simtel.zst") source = EventSource(filename) events = [copy.deepcopy(event) for event in source] event = events[4] layout = set(source.subarray.tel_ids) .. GENERATED FROM PYTHON SOURCE LINES 43-46 Choose event with LST ~~~~~~~~~~~~~~~~~~~~~ .. GENERATED FROM PYTHON SOURCE LINES 49-52 This ensures that the telescope is not “parked” (as it would be in an event where it is not triggered) but is actually pointing to a source. .. GENERATED FROM PYTHON SOURCE LINES 52-57 .. code-block:: Python :lineno-start: 53 print(f"Telescope with data: {event.r1.tel.keys()}") tel_id = 3 .. rst-class:: sphx-glr-script-out .. code-block:: none Telescope with data: dict_keys([3, 4, 8, 9, 14, 15, 23, 24, 25, 29, 32, 36, 42, 46, 52, 56, 103, 104, 109, 110, 118, 119, 120, 126, 127, 129, 130, 137, 139, 143, 161]) .. GENERATED FROM PYTHON SOURCE LINES 58-71 AltAz ----- See `Astropy Docs on AltAz `__. Pointing direction of telescopes or the origin of a simulated shower are described in the ``AltAz`` frame. This is a local, angular coordinate frame, with angles ``altitude`` and ``azimuth``. Altitude is the measured from the Horizon (0°) to the Zenith (90°). For the azimuth, there are different conventions. In Astropy und thus ctapipe, Azimuth is oriented East of North (i.e., N=0°, E=90°). .. GENERATED FROM PYTHON SOURCE LINES 71-87 .. code-block:: Python :lineno-start: 73 obstime = Time("2013-11-01T03:00") location = EarthLocation.of_site("Roque de los Muchachos") altaz = AltAz(location=location, obstime=obstime) array_pointing = SkyCoord( alt=event.monitoring.pointing.array_azimuth, az=event.monitoring.pointing.array_altitude, frame=altaz, ) print(array_pointing) .. rst-class:: sphx-glr-script-out .. code-block:: none .. GENERATED FROM PYTHON SOURCE LINES 88-108 CameraFrame ----------- Camera coordinate frame. The camera frame is a 2d cartesian frame, describing position of objects in the focal plane of the telescope. The frame is defined as in H.E.S.S., starting at the horizon, the telescope is pointed to magnetic north in azimuth and then up to zenith. Now, x points north and y points west, so in this orientation, the camera coordinates line up with the CORSIKA ground coordinate system. MAGIC and FACT use a different camera coordinate system: Standing at the dish, looking at the camera, x points right, y points up. To transform MAGIC/FACT to ctapipe, do x’ = -y, y’ = -x. **Typical usage**: Position of pixels in the focal plane. .. GENERATED FROM PYTHON SOURCE LINES 108-114 .. code-block:: Python :lineno-start: 109 geometry = source.subarray.tel[tel_id].camera.geometry pix_x = geometry.pix_x pix_y = geometry.pix_y focal_length = source.subarray.tel[tel_id].optics.equivalent_focal_length .. GENERATED FROM PYTHON SOURCE LINES 115-131 .. code-block:: Python :lineno-start: 115 telescope_pointing = SkyCoord( alt=event.monitoring.tel[tel_id].pointing.altitude, az=event.monitoring.tel[tel_id].pointing.azimuth, frame=altaz, ) camera_frame = CameraFrame( focal_length=focal_length, rotation=0 * u.deg, telescope_pointing=telescope_pointing, ) cam_coords = SkyCoord(x=pix_x, y=pix_y, frame=camera_frame) print(cam_coords) .. rst-class:: sphx-glr-script-out .. code-block:: none , obstime=None, location=None): (x, y) in m [( 0. , 0. ), (-0.0377967 , -0.03273243), (-0.04724547, 0.01636666), ..., ( 0.67088033, -0.96561883), (-0.45356484, -1.08017546), (-0.50081024, -1.06380875)]> .. GENERATED FROM PYTHON SOURCE LINES 132-139 .. code-block:: Python :lineno-start: 132 plt.scatter(cam_coords.x, cam_coords.y) plt.title(f"Camera type: {geometry.name}") plt.xlabel(f"x / {cam_coords.x.unit}") plt.ylabel(f"y / {cam_coords.y.unit}") plt.axis("square") .. image-sg:: /auto_examples/tutorials/images/sphx_glr_coordinates_example_001.png :alt: Camera type: LSTCam :srcset: /auto_examples/tutorials/images/sphx_glr_coordinates_example_001.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none (np.float64(-1.2888595582149416), np.float64(1.2888595582149416), np.float64(-1.1881999735124904), np.float64(1.3895191429173928)) .. GENERATED FROM PYTHON SOURCE LINES 140-145 The implementation of the coordinate system with astropy makes it easier to use time of the observation and location of the observing site, to understand, for example which stars are visible during a certain night and how they might be visible in the camera. .. GENERATED FROM PYTHON SOURCE LINES 145-192 .. code-block:: Python :lineno-start: 147 location = EarthLocation.of_site("Roque de los Muchachos") obstime = Time("2018-11-01T04:00") crab = SkyCoord.from_name("crab nebula") altaz = AltAz(location=location, obstime=obstime) pointing = crab.transform_to(altaz) camera_frame = CameraFrame( telescope_pointing=pointing, focal_length=focal_length, obstime=obstime, location=location, ) subarray = SubarrayDescription.read("dataset://gamma_prod5.simtel.zst") cam = subarray.tel[1].camera.geometry fig, ax = plt.subplots() display = CameraDisplay(cam, ax=ax) ax.set_title( f"La Palma, {obstime}, az={pointing.az.deg:.1f}°, zenith={pointing.zen.deg:.1f}°, camera={geometry.name}" ) for i, name in enumerate(["crab nebula", "o tau", "zet tau"]): star = SkyCoord.from_name(name) star_cam = star.transform_to(camera_frame) x = star_cam.x.to_value(u.m) y = star_cam.y.to_value(u.m) ax.plot(x, y, marker="*", color=f"C{i}") ax.annotate( name, xy=(x, y), xytext=(5, 5), textcoords="offset points", color=f"C{i}", ) plt.show() .. image-sg:: /auto_examples/tutorials/images/sphx_glr_coordinates_example_002.png :alt: La Palma, 2018-11-01T04:00:00.000, az=169.1°, zenith=6.8°, camera=LSTCam :srcset: /auto_examples/tutorials/images/sphx_glr_coordinates_example_002.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 193-209 TelescopeFrame -------------- Telescope coordinate frame. A ``Frame`` using a ``UnitSphericalRepresentation``. This is basically the same as a ``HorizonCoordinate``, but the origin is at the telescope’s pointing direction. This is what astropy calls a ``SkyOffsetFrame``. The axis of the telescope frame, ``fov_lon`` and ``fov_lat``, are aligned with the horizontal system’s azimuth and altitude respectively. Pointing corrections should applied to the transformation between this frame and the camera frame. .. GENERATED FROM PYTHON SOURCE LINES 209-217 .. code-block:: Python :lineno-start: 210 telescope_frame = TelescopeFrame( telescope_pointing=pointing, obstime=pointing.obstime, location=pointing.location, ) telescope_coords = cam_coords.transform_to(telescope_frame) .. GENERATED FROM PYTHON SOURCE LINES 218-243 .. code-block:: Python :lineno-start: 218 wrap_angle = telescope_pointing.az + 180 * u.deg plt.axis("equal") plt.scatter( telescope_coords.fov_lon.deg, telescope_coords.fov_lat.deg, alpha=0.2, color="gray" ) for i, name in enumerate(["crab nebula", "o tau", "zet tau"]): star = SkyCoord.from_name(name) star_tel = star.transform_to(telescope_frame) plt.plot(star_tel.fov_lon.deg, star_tel.fov_lat.deg, "*", ms=10) plt.annotate( name, xy=(star_tel.fov_lon.deg, star_tel.fov_lat.deg), xytext=(5, 5), textcoords="offset points", color=f"C{i}", ) plt.xlabel("fov_lon / {}".format(telescope_coords.altaz.az.unit)) plt.ylabel("fov_lat / {}".format(telescope_coords.altaz.alt.unit)) .. image-sg:: /auto_examples/tutorials/images/sphx_glr_coordinates_example_003.png :alt: coordinates example :srcset: /auto_examples/tutorials/images/sphx_glr_coordinates_example_003.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none Text(101.97222222222221, 0.5, 'fov_lat / deg') .. GENERATED FROM PYTHON SOURCE LINES 244-256 NominalFrame ------------ Nominal coordinate frame. A Frame using a ``UnitSphericalRepresentation``. This is basically the same as a ``HorizonCoordinate``, but the origin is at an arbitrary position in the sky. This is what astropy calls a ``SkyOffsetFrame`` If the telescopes are in divergent pointing, this ``Frame`` can be used to transform to a common system. - 2D reconstruction (``HillasIntersector``) is performed in this frame - 3D reconstruction (``HillasReconstructor``) doesn’t need this frame .. GENERATED FROM PYTHON SOURCE LINES 259-261 Let’s play a bit with 3 LSTs with divergent pointing .. GENERATED FROM PYTHON SOURCE LINES 261-292 .. code-block:: Python :lineno-start: 262 location = EarthLocation.of_site("Roque de los Muchachos") obstime = Time("2018-11-01T02:00") altaz = AltAz(location=location, obstime=obstime) crab = SkyCoord.from_name("crab nebula") # let's observe crab array_pointing = crab.transform_to(altaz) # let the telescopes point to different positions alt_offsets = u.Quantity([1, -1, -1], u.deg) az_offsets = u.Quantity([0, -2, +2], u.deg) tel_pointings = SkyCoord( alt=array_pointing.alt + alt_offsets, az=array_pointing.az + az_offsets, frame=altaz, ) camera_frames = CameraFrame( telescope_pointing=tel_pointings, # multiple pointings, so we get multiple frames focal_length=focal_length, obstime=obstime, location=location, ) nom_frame = NominalFrame(origin=array_pointing, obstime=obstime, location=location) .. GENERATED FROM PYTHON SOURCE LINES 293-334 .. code-block:: Python :lineno-start: 293 fig, ax = plt.subplots(figsize=(15, 10)) ax.set_aspect(1) for i in range(3): cam_coord = SkyCoord(x=pix_x, y=pix_y, frame=camera_frames[i]) nom_coord = cam_coord.transform_to(nom_frame) ax.scatter( x=nom_coord.fov_lon.deg, y=nom_coord.fov_lat.deg, label=f"Telescope {i + 1}", s=30, alpha=0.15, ) for i, name in enumerate(["Crab", "o tau", "zet tau"]): s = SkyCoord.from_name(name) s_nom = s.transform_to(nom_frame) ax.plot( s_nom.fov_lon.deg, s_nom.fov_lat.deg, "*", ms=10, ) ax.annotate( name, xy=(s_nom.fov_lon.deg, s_nom.fov_lat.deg), xytext=(5, 5), textcoords="offset points", color=f"C{i}", ) ax.set_xlabel("fov_lon / deg") ax.set_ylabel("fov_lat / deg") ax.legend() plt.show() .. image-sg:: /auto_examples/tutorials/images/sphx_glr_coordinates_example_004.png :alt: coordinates example :srcset: /auto_examples/tutorials/images/sphx_glr_coordinates_example_004.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 335-346 GroundFrame ----------- Ground coordinate frame. The ground coordinate frame is a simple cartesian frame describing the 3 dimensional position of objects compared to the array ground level in relation to the nomial centre of the array. Typically this frame will be used for describing the position on telescopes and equipment **Typical usage**: positions of telescopes on the ground (x, y, z) .. GENERATED FROM PYTHON SOURCE LINES 346-350 .. code-block:: Python :lineno-start: 347 source.subarray.peek() .. image-sg:: /auto_examples/tutorials/images/sphx_glr_coordinates_example_005.png :alt: MonteCarloArray :srcset: /auto_examples/tutorials/images/sphx_glr_coordinates_example_005.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none .. GENERATED FROM PYTHON SOURCE LINES 351-354 In case a layout is selected, the following line will produce a different output from the picture above. .. GENERATED FROM PYTHON SOURCE LINES 354-358 .. code-block:: Python :lineno-start: 355 source.subarray.select_subarray(layout, name="Prod3b layout").peek() .. image-sg:: /auto_examples/tutorials/images/sphx_glr_coordinates_example_006.png :alt: Prod3b layout :srcset: /auto_examples/tutorials/images/sphx_glr_coordinates_example_006.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none .. GENERATED FROM PYTHON SOURCE LINES 359-363 .. figure:: ground_frame.png :alt: Ground Frame Ground Frame .. GENERATED FROM PYTHON SOURCE LINES 366-369 In this image all the telescope from the ``gamma_test.simtel.gz`` file are plotted as spheres in the GroundFrame. .. GENERATED FROM PYTHON SOURCE LINES 372-375 TiltedGroundFrame ----------------- .. GENERATED FROM PYTHON SOURCE LINES 378-388 Tilted ground coordinate frame. The tilted ground coordinate frame is a cartesian system describing the 2 dimensional projected positions of objects in a tilted plane described by pointing_direction. The plane is rotated along the z_axis by the azimuth of the ``pointing_direction`` and then it is inclined with an angle equal to the zenith angle of the ``pointing_direction``. This frame is used for the reconstruction of the shower core position. .. GENERATED FROM PYTHON SOURCE LINES 391-395 .. figure:: tilted_ground_frame.png :alt: Tilted Ground Frame Tilted Ground Frame .. GENERATED FROM PYTHON SOURCE LINES 398-407 This image picture both the telescopes in the GroundFrame (red) and in the TiltedGroundFrame (green) are displayed: in this case since the azimuth of the ``pointing_direction`` is 0 degrees, then the plane is just tilted according to the zenith angle. For playing with these and with more 3D models of the telescopes themselves, have a look at the `CREED_VTK `__ library. .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 16.840 seconds) .. _sphx_glr_download_auto_examples_tutorials_coordinates_example.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: coordinates_example.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: coordinates_example.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: coordinates_example.zip ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_