1# -*- coding: utf-8 -*- 2# Licensed under a 3-clause BSD style license - see LICENSE.rst 3 4from astropy import units as u 5from astropy.utils.decorators import format_doc 6from astropy.time import Time 7from astropy.coordinates import representation as r 8from astropy.coordinates.baseframe import (BaseCoordinateFrame, 9 RepresentationMapping, 10 frame_transform_graph, base_doc) 11from astropy.coordinates.transformations import AffineTransform 12from astropy.coordinates.attributes import DifferentialAttribute 13 14from .baseradec import BaseRADecFrame, doc_components as doc_components_radec 15from .icrs import ICRS 16from .galactic import Galactic 17from .fk4 import FK4 18 19# For speed 20J2000 = Time('J2000') 21 22v_bary_Schoenrich2010 = r.CartesianDifferential([11.1, 12.24, 7.25]*u.km/u.s) 23 24__all__ = ['LSR', 'GalacticLSR', 'LSRK', 'LSRD'] 25 26 27doc_footer_lsr = """ 28 Other parameters 29 ---------------- 30 v_bary : `~astropy.coordinates.representation.CartesianDifferential` 31 The velocity of the solar system barycenter with respect to the LSR, in 32 Galactic cartesian velocity components. 33""" 34 35 36@format_doc(base_doc, components=doc_components_radec, footer=doc_footer_lsr) 37class LSR(BaseRADecFrame): 38 r"""A coordinate or frame in the Local Standard of Rest (LSR). 39 40 This coordinate frame is axis-aligned and co-spatial with `ICRS`, but has 41 a velocity offset relative to the solar system barycenter to remove the 42 peculiar motion of the sun relative to the LSR. Roughly, the LSR is the mean 43 velocity of the stars in the solar neighborhood, but the precise definition 44 of which depends on the study. As defined in Schönrich et al. (2010): 45 "The LSR is the rest frame at the location of the Sun of a star that would 46 be on a circular orbit in the gravitational potential one would obtain by 47 azimuthally averaging away non-axisymmetric features in the actual Galactic 48 potential." No such orbit truly exists, but it is still a commonly used 49 velocity frame. 50 51 We use default values from Schönrich et al. (2010) for the barycentric 52 velocity relative to the LSR, which is defined in Galactic (right-handed) 53 cartesian velocity components 54 :math:`(U, V, W) = (11.1, 12.24, 7.25)~{{\rm km}}~{{\rm s}}^{{-1}}`. These 55 values are customizable via the ``v_bary`` argument which specifies the 56 velocity of the solar system barycenter with respect to the LSR. 57 58 The frame attributes are listed under **Other Parameters**. 59 """ 60 61 # frame attributes: 62 v_bary = DifferentialAttribute(default=v_bary_Schoenrich2010, 63 allowed_classes=[r.CartesianDifferential]) 64 65 66@frame_transform_graph.transform(AffineTransform, ICRS, LSR) 67def icrs_to_lsr(icrs_coord, lsr_frame): 68 v_bary_gal = Galactic(lsr_frame.v_bary.to_cartesian()) 69 v_bary_icrs = v_bary_gal.transform_to(icrs_coord) 70 v_offset = v_bary_icrs.data.represent_as(r.CartesianDifferential) 71 offset = r.CartesianRepresentation([0, 0, 0]*u.au, differentials=v_offset) 72 return None, offset 73 74 75@frame_transform_graph.transform(AffineTransform, LSR, ICRS) 76def lsr_to_icrs(lsr_coord, icrs_frame): 77 v_bary_gal = Galactic(lsr_coord.v_bary.to_cartesian()) 78 v_bary_icrs = v_bary_gal.transform_to(icrs_frame) 79 v_offset = v_bary_icrs.data.represent_as(r.CartesianDifferential) 80 offset = r.CartesianRepresentation([0, 0, 0]*u.au, differentials=-v_offset) 81 return None, offset 82 83 84# ------------------------------------------------------------------------------ 85 86 87doc_components_gal = """ 88 l : `~astropy.coordinates.Angle`, optional, keyword-only 89 The Galactic longitude for this object (``b`` must also be given and 90 ``representation`` must be None). 91 b : `~astropy.coordinates.Angle`, optional, keyword-only 92 The Galactic latitude for this object (``l`` must also be given and 93 ``representation`` must be None). 94 distance : `~astropy.units.Quantity` ['length'], optional, keyword-only 95 The Distance for this object along the line-of-sight. 96 (``representation`` must be None). 97 98 pm_l_cosb : `~astropy.units.Quantity` ['angular speed'], optional, keyword-only 99 The proper motion in Galactic longitude (including the ``cos(b)`` term) 100 for this object (``pm_b`` must also be given). 101 pm_b : `~astropy.units.Quantity` ['angular speed'], optional, keyword-only 102 The proper motion in Galactic latitude for this object (``pm_l_cosb`` 103 must also be given). 104 radial_velocity : `~astropy.units.Quantity` ['speed'], optional, keyword-only 105 The radial velocity of this object. 106""" 107 108 109@format_doc(base_doc, components=doc_components_gal, footer=doc_footer_lsr) 110class GalacticLSR(BaseCoordinateFrame): 111 r"""A coordinate or frame in the Local Standard of Rest (LSR), axis-aligned 112 to the `Galactic` frame. 113 114 This coordinate frame is axis-aligned and co-spatial with `ICRS`, but has 115 a velocity offset relative to the solar system barycenter to remove the 116 peculiar motion of the sun relative to the LSR. Roughly, the LSR is the mean 117 velocity of the stars in the solar neighborhood, but the precise definition 118 of which depends on the study. As defined in Schönrich et al. (2010): 119 "The LSR is the rest frame at the location of the Sun of a star that would 120 be on a circular orbit in the gravitational potential one would obtain by 121 azimuthally averaging away non-axisymmetric features in the actual Galactic 122 potential." No such orbit truly exists, but it is still a commonly used 123 velocity frame. 124 125 We use default values from Schönrich et al. (2010) for the barycentric 126 velocity relative to the LSR, which is defined in Galactic (right-handed) 127 cartesian velocity components 128 :math:`(U, V, W) = (11.1, 12.24, 7.25)~{{\rm km}}~{{\rm s}}^{{-1}}`. These 129 values are customizable via the ``v_bary`` argument which specifies the 130 velocity of the solar system barycenter with respect to the LSR. 131 132 The frame attributes are listed under **Other Parameters**. 133 """ 134 135 frame_specific_representation_info = { 136 r.SphericalRepresentation: [ 137 RepresentationMapping('lon', 'l'), 138 RepresentationMapping('lat', 'b') 139 ] 140 } 141 142 default_representation = r.SphericalRepresentation 143 default_differential = r.SphericalCosLatDifferential 144 145 # frame attributes: 146 v_bary = DifferentialAttribute(default=v_bary_Schoenrich2010) 147 148 149@frame_transform_graph.transform(AffineTransform, Galactic, GalacticLSR) 150def galactic_to_galacticlsr(galactic_coord, lsr_frame): 151 v_bary_gal = Galactic(lsr_frame.v_bary.to_cartesian()) 152 v_offset = v_bary_gal.data.represent_as(r.CartesianDifferential) 153 offset = r.CartesianRepresentation([0, 0, 0]*u.au, differentials=v_offset) 154 return None, offset 155 156 157@frame_transform_graph.transform(AffineTransform, GalacticLSR, Galactic) 158def galacticlsr_to_galactic(lsr_coord, galactic_frame): 159 v_bary_gal = Galactic(lsr_coord.v_bary.to_cartesian()) 160 v_offset = v_bary_gal.data.represent_as(r.CartesianDifferential) 161 offset = r.CartesianRepresentation([0, 0, 0]*u.au, differentials=-v_offset) 162 return None, offset 163 164 165# ------------------------------------------------------------------------------ 166 167# The LSRK velocity frame, defined as having a velocity of 20 km/s towards 168# RA=270 Dec=30 (B1900) relative to the solar system Barycenter. This is defined 169# in: 170# 171# Gordon 1975, Methods of Experimental Physics: Volume 12: 172# Astrophysics, Part C: Radio Observations - Section 6.1.5. 173 174 175class LSRK(BaseRADecFrame): 176 r""" 177 A coordinate or frame in the Kinematic Local Standard of Rest (LSR). 178 179 This frame is defined as having a velocity of 20 km/s towards RA=270 Dec=30 180 (B1900) relative to the solar system Barycenter. This is defined in: 181 182 Gordon 1975, Methods of Experimental Physics: Volume 12: 183 Astrophysics, Part C: Radio Observations - Section 6.1.5. 184 185 This coordinate frame is axis-aligned and co-spatial with `ICRS`, but has 186 a velocity offset relative to the solar system barycenter to remove the 187 peculiar motion of the sun relative to the LSRK. 188 """ 189 190 191# NOTE: To avoid a performance penalty at import time, we hard-code the ICRS 192# offsets here. The code to generate the offsets is provided for reproducibility. 193# GORDON1975_V_BARY = 20*u.km/u.s 194# GORDON1975_DIRECTION = FK4(ra=270*u.deg, dec=30*u.deg, equinox='B1900') 195# V_OFFSET_LSRK = ((GORDON1975_V_BARY * GORDON1975_DIRECTION.transform_to(ICRS()).data) 196# .represent_as(r.CartesianDifferential)) 197 198V_OFFSET_LSRK = r.CartesianDifferential([0.28999706839034606, 199 -17.317264789717928, 200 10.00141199546947]*u.km/u.s) 201 202ICRS_LSRK_OFFSET = r.CartesianRepresentation([0, 0, 0]*u.au, differentials=V_OFFSET_LSRK) 203LSRK_ICRS_OFFSET = r.CartesianRepresentation([0, 0, 0]*u.au, differentials=-V_OFFSET_LSRK) 204 205 206@frame_transform_graph.transform(AffineTransform, ICRS, LSRK) 207def icrs_to_lsrk(icrs_coord, lsr_frame): 208 return None, ICRS_LSRK_OFFSET 209 210 211@frame_transform_graph.transform(AffineTransform, LSRK, ICRS) 212def lsrk_to_icrs(lsr_coord, icrs_frame): 213 return None, LSRK_ICRS_OFFSET 214 215 216# ------------------------------------------------------------------------------ 217 218# The LSRD velocity frame, defined as a velocity of U=9 km/s, V=12 km/s, 219# and W=7 km/s in Galactic coordinates or 16.552945 km/s 220# towards l=53.13 b=25.02. This is defined in: 221# 222# Delhaye 1965, Solar Motion and Velocity Distribution of 223# Common Stars. 224 225 226class LSRD(BaseRADecFrame): 227 r""" 228 A coordinate or frame in the Dynamical Local Standard of Rest (LSRD) 229 230 This frame is defined as a velocity of U=9 km/s, V=12 km/s, 231 and W=7 km/s in Galactic coordinates or 16.552945 km/s 232 towards l=53.13 b=25.02. This is defined in: 233 234 Delhaye 1965, Solar Motion and Velocity Distribution of 235 Common Stars. 236 237 This coordinate frame is axis-aligned and co-spatial with `ICRS`, but has 238 a velocity offset relative to the solar system barycenter to remove the 239 peculiar motion of the sun relative to the LSRD. 240 """ 241 242 243# NOTE: To avoid a performance penalty at import time, we hard-code the ICRS 244# offsets here. The code to generate the offsets is provided for reproducibility. 245# V_BARY_DELHAYE1965 = r.CartesianDifferential([9, 12, 7] * u.km/u.s) 246# V_OFFSET_LSRD = (Galactic(V_BARY_DELHAYE1965.to_cartesian()).transform_to(ICRS()).data 247# .represent_as(r.CartesianDifferential)) 248 249V_OFFSET_LSRD = r.CartesianDifferential([-0.6382306360182073, 250 -14.585424483191094, 251 7.8011572411006815]*u.km/u.s) 252 253ICRS_LSRD_OFFSET = r.CartesianRepresentation([0, 0, 0]*u.au, differentials=V_OFFSET_LSRD) 254LSRD_ICRS_OFFSET = r.CartesianRepresentation([0, 0, 0]*u.au, differentials=-V_OFFSET_LSRD) 255 256 257@frame_transform_graph.transform(AffineTransform, ICRS, LSRD) 258def icrs_to_lsrd(icrs_coord, lsr_frame): 259 return None, ICRS_LSRD_OFFSET 260 261 262@frame_transform_graph.transform(AffineTransform, LSRD, ICRS) 263def lsrd_to_icrs(lsr_coord, icrs_frame): 264 return None, LSRD_ICRS_OFFSET 265 266 267# ------------------------------------------------------------------------------ 268 269# Create loopback transformations 270frame_transform_graph._add_merged_transform(LSR, ICRS, LSR) 271frame_transform_graph._add_merged_transform(GalacticLSR, Galactic, GalacticLSR) 272