1# -*- coding: utf-8 -*- 2# Licensed under a 3-clause BSD style license - see LICENSE.rst 3 4import copy 5from collections.abc import MappingView 6from types import MappingProxyType 7 8import numpy as np 9 10from astropy import units as u 11from astropy.utils.state import ScienceState 12from astropy.utils.decorators import format_doc, classproperty, deprecated 13from astropy.coordinates.angles import Angle 14from astropy.coordinates.matrix_utilities import rotation_matrix, matrix_product, matrix_transpose 15from astropy.coordinates import representation as r 16from astropy.coordinates.baseframe import (BaseCoordinateFrame, 17 frame_transform_graph, 18 base_doc) 19from astropy.coordinates.attributes import (CoordinateAttribute, 20 QuantityAttribute, 21 DifferentialAttribute) 22from astropy.coordinates.transformations import AffineTransform 23from astropy.coordinates.errors import ConvertError 24 25from .icrs import ICRS 26 27__all__ = ['Galactocentric'] 28 29 30# Measured by minimizing the difference between a plane of coordinates along 31# l=0, b=[-90,90] and the Galactocentric x-z plane 32# This is not used directly, but accessed via `get_roll0`. We define it here to 33# prevent having to create new Angle objects every time `get_roll0` is called. 34_ROLL0 = Angle(58.5986320306*u.degree) 35 36 37class _StateProxy(MappingView): 38 """ 39 `~collections.abc.MappingView` with a read-only ``getitem`` through 40 `~types.MappingProxyType`. 41 42 """ 43 44 def __init__(self, mapping): 45 super().__init__(mapping) 46 self._mappingproxy = MappingProxyType(self._mapping) # read-only 47 48 def __getitem__(self, key): 49 """Read-only ``getitem``.""" 50 return self._mappingproxy[key] 51 52 def __deepcopy__(self, memo): 53 return copy.deepcopy(self._mapping, memo=memo) 54 55 56class galactocentric_frame_defaults(ScienceState): 57 """This class controls the global setting of default values for the frame 58 attributes in the `~astropy.coordinates.Galactocentric` frame, which may be 59 updated in future versions of ``astropy``. Note that when using 60 `~astropy.coordinates.Galactocentric`, changing values here will not affect 61 any attributes that are set explicitly by passing values in to the 62 `~astropy.coordinates.Galactocentric` initializer. Modifying these defaults 63 will only affect the frame attribute values when using the frame as, e.g., 64 ``Galactocentric`` or ``Galactocentric()`` with no explicit arguments. 65 66 This class controls the parameter settings by specifying a string name, 67 with the following pre-specified options: 68 69 - 'pre-v4.0': The current default value, which sets the default frame 70 attribute values to their original (pre-astropy-v4.0) values. 71 - 'v4.0': The attribute values as updated in Astropy version 4.0. 72 - 'latest': An alias of the most recent parameter set (currently: 'v4.0') 73 74 Alternatively, user-defined parameter settings may be registered, with 75 :meth:`~astropy.coordinates.galactocentric_frame_defaults.register`, 76 and used identically as pre-specified parameter sets. At minimum, 77 registrations must have unique names and a dictionary of parameters 78 with keys "galcen_coord", "galcen_distance", "galcen_v_sun", "z_sun", 79 "roll". See examples below. 80 81 This class also tracks the references for all parameter values in the 82 attribute ``references``, as well as any further information the registry. 83 The pre-specified options can be extended to include similar 84 state information as user-defined parameter settings -- for example, to add 85 parameter uncertainties. 86 87 The preferred method for getting a parameter set and metadata, by name, is 88 :meth:`~galactocentric_frame_defaults.get_from_registry` since 89 it ensures the immutability of the registry. 90 91 See :ref:`astropy:astropy-coordinates-galactocentric-defaults` for more 92 information. 93 94 Examples 95 -------- 96 The default `~astropy.coordinates.Galactocentric` frame parameters can be 97 modified globally:: 98 99 >>> from astropy.coordinates import galactocentric_frame_defaults 100 >>> _ = galactocentric_frame_defaults.set('v4.0') # doctest: +SKIP 101 >>> Galactocentric() # doctest: +SKIP 102 <Galactocentric Frame (galcen_coord=<ICRS Coordinate: (ra, dec) in deg 103 (266.4051, -28.936175)>, galcen_distance=8.122 kpc, galcen_v_sun=(12.9, 245.6, 7.78) km / s, z_sun=20.8 pc, roll=0.0 deg)> 104 >>> _ = galactocentric_frame_defaults.set('pre-v4.0') # doctest: +SKIP 105 >>> Galactocentric() # doctest: +SKIP 106 <Galactocentric Frame (galcen_coord=<ICRS Coordinate: (ra, dec) in deg 107 (266.4051, -28.936175)>, galcen_distance=8.3 kpc, galcen_v_sun=(11.1, 232.24, 7.25) km / s, z_sun=27.0 pc, roll=0.0 deg)> 108 109 The default parameters can also be updated by using this class as a context 110 manager:: 111 112 >>> with galactocentric_frame_defaults.set('pre-v4.0'): 113 ... print(Galactocentric()) # doctest: +FLOAT_CMP 114 <Galactocentric Frame (galcen_coord=<ICRS Coordinate: (ra, dec) in deg 115 (266.4051, -28.936175)>, galcen_distance=8.3 kpc, galcen_v_sun=(11.1, 232.24, 7.25) km / s, z_sun=27.0 pc, roll=0.0 deg)> 116 117 Again, changing the default parameter values will not affect frame 118 attributes that are explicitly specified:: 119 120 >>> import astropy.units as u 121 >>> with galactocentric_frame_defaults.set('pre-v4.0'): 122 ... print(Galactocentric(galcen_distance=8.0*u.kpc)) # doctest: +FLOAT_CMP 123 <Galactocentric Frame (galcen_coord=<ICRS Coordinate: (ra, dec) in deg 124 (266.4051, -28.936175)>, galcen_distance=8.0 kpc, galcen_v_sun=(11.1, 232.24, 7.25) km / s, z_sun=27.0 pc, roll=0.0 deg)> 125 126 Additional parameter sets may be registered, for instance to use the 127 Dehnen & Binney (1998) measurements of the solar motion. We can also 128 add metadata, such as the 1-sigma errors. In this example we will modify 129 the required key "parameters", change the recommended key "references" to 130 match "parameters", and add the extra key "error" (any key can be added):: 131 132 >>> state = galactocentric_frame_defaults.get_from_registry("v4.0") 133 >>> state["parameters"]["galcen_v_sun"] = (10.00, 225.25, 7.17) * (u.km / u.s) 134 >>> state["references"]["galcen_v_sun"] = "https://ui.adsabs.harvard.edu/full/1998MNRAS.298..387D" 135 >>> state["error"] = {"galcen_v_sun": (0.36, 0.62, 0.38) * (u.km / u.s)} 136 >>> galactocentric_frame_defaults.register(name="DB1998", **state) 137 138 Just as in the previous examples, the new parameter set can be retrieved with:: 139 140 >>> state = galactocentric_frame_defaults.get_from_registry("DB1998") 141 >>> print(state["error"]["galcen_v_sun"]) # doctest: +FLOAT_CMP 142 [0.36 0.62 0.38] km / s 143 144 """ 145 146 _latest_value = 'v4.0' 147 _value = None 148 _references = None 149 _state = dict() # all other data 150 151 # Note: _StateProxy() produces read-only view of enclosed mapping. 152 _registry = { 153 "v4.0": { 154 "parameters": _StateProxy( 155 { 156 "galcen_coord": ICRS( 157 ra=266.4051 * u.degree, dec=-28.936175 * u.degree 158 ), 159 "galcen_distance": 8.122 * u.kpc, 160 "galcen_v_sun": r.CartesianDifferential( 161 [12.9, 245.6, 7.78] * (u.km / u.s) 162 ), 163 "z_sun": 20.8 * u.pc, 164 "roll": 0 * u.deg, 165 } 166 ), 167 "references": _StateProxy( 168 { 169 "galcen_coord": "https://ui.adsabs.harvard.edu/abs/2004ApJ...616..872R", 170 "galcen_distance": "https://ui.adsabs.harvard.edu/abs/2018A%26A...615L..15G", 171 "galcen_v_sun": [ 172 "https://ui.adsabs.harvard.edu/abs/2018RNAAS...2..210D", 173 "https://ui.adsabs.harvard.edu/abs/2018A%26A...615L..15G", 174 "https://ui.adsabs.harvard.edu/abs/2004ApJ...616..872R", 175 ], 176 "z_sun": "https://ui.adsabs.harvard.edu/abs/2019MNRAS.482.1417B", 177 "roll": None, 178 } 179 ), 180 }, 181 "pre-v4.0": { 182 "parameters": _StateProxy( 183 { 184 "galcen_coord": ICRS( 185 ra=266.4051 * u.degree, dec=-28.936175 * u.degree 186 ), 187 "galcen_distance": 8.3 * u.kpc, 188 "galcen_v_sun": r.CartesianDifferential( 189 [11.1, 220 + 12.24, 7.25] * (u.km / u.s) 190 ), 191 "z_sun": 27.0 * u.pc, 192 "roll": 0 * u.deg, 193 } 194 ), 195 "references": _StateProxy( 196 { 197 "galcen_coord": "https://ui.adsabs.harvard.edu/abs/2004ApJ...616..872R", 198 "galcen_distance": "https://ui.adsabs.harvard.edu/#abs/2009ApJ...692.1075G", 199 "galcen_v_sun": [ 200 "https://ui.adsabs.harvard.edu/#abs/2010MNRAS.403.1829S", 201 "https://ui.adsabs.harvard.edu/#abs/2015ApJS..216...29B", 202 ], 203 "z_sun": "https://ui.adsabs.harvard.edu/#abs/2001ApJ...553..184C", 204 "roll": None, 205 } 206 ), 207 }, 208 } 209 210 @classproperty # read-only 211 def parameters(cls): 212 return cls._value 213 214 @classproperty # read-only 215 def references(cls): 216 return cls._references 217 218 @classmethod 219 def get_from_registry(cls, name: str): 220 """ 221 Return Galactocentric solar parameters and metadata given string names 222 for the parameter sets. This method ensures the returned state is a 223 mutable copy, so any changes made do not affect the registry state. 224 225 Returns 226 ------- 227 state : dict 228 Copy of the registry for the string name. 229 Should contain, at minimum: 230 231 - "parameters": dict 232 Galactocentric solar parameters 233 - "references" : Dict[str, Union[str, Sequence[str]]] 234 References for "parameters". 235 Fields are str or sequence of str. 236 237 Raises 238 ------ 239 KeyError 240 If invalid string input to registry 241 to retrieve solar parameters for Galactocentric frame. 242 243 """ 244 # Resolve the meaning of 'latest': latest parameter set is from v4.0 245 # - update this as newer parameter choices are added 246 if name == 'latest': 247 name = cls._latest_value 248 249 # Get the state from the registry. 250 # Copy to ensure registry is immutable to modifications of "_value". 251 # Raises KeyError if `name` is invalid string input to registry 252 # to retrieve solar parameters for Galactocentric frame. 253 state = copy.deepcopy(cls._registry[name]) # ensure mutable 254 255 return state 256 257 @deprecated("v4.2", alternative="`get_from_registry`") 258 @classmethod 259 def get_solar_params_from_string(cls, arg): 260 """ 261 Return Galactocentric solar parameters given string names 262 for the parameter sets. 263 264 Returns 265 ------- 266 parameters : dict 267 Copy of Galactocentric solar parameters from registry 268 269 Raises 270 ------ 271 KeyError 272 If invalid string input to registry 273 to retrieve solar parameters for Galactocentric frame. 274 275 """ 276 return cls.get_from_registry(arg)["parameters"] 277 278 @classmethod 279 def validate(cls, value): 280 if value is None: 281 value = cls._latest_value 282 283 if isinstance(value, str): 284 state = cls.get_from_registry(value) 285 cls._references = state["references"] 286 cls._state = state 287 parameters = state["parameters"] 288 289 elif isinstance(value, dict): 290 parameters = value 291 292 elif isinstance(value, Galactocentric): 293 # turn the frame instance into a dict of frame attributes 294 parameters = dict() 295 for k in value.frame_attributes: 296 parameters[k] = getattr(value, k) 297 cls._references = value.frame_attribute_references.copy() 298 cls._state = dict(parameters=parameters, 299 references=cls._references) 300 301 else: 302 raise ValueError("Invalid input to retrieve solar parameters for " 303 "Galactocentric frame: input must be a string, " 304 "dict, or Galactocentric instance") 305 306 return parameters 307 308 @classmethod 309 def register(cls, name: str, parameters: dict, references=None, 310 **meta: dict): 311 """Register a set of parameters. 312 313 Parameters 314 ---------- 315 name : str 316 The registration name for the parameter and metadata set. 317 parameters : dict 318 The solar parameters for Galactocentric frame. 319 references : dict or None, optional 320 References for contents of `parameters`. 321 None becomes empty dict. 322 **meta : dict, optional 323 Any other properties to register. 324 325 """ 326 # check on contents of `parameters` 327 must_have = {"galcen_coord", "galcen_distance", "galcen_v_sun", 328 "z_sun", "roll"} 329 missing = must_have.difference(parameters) 330 if missing: 331 raise ValueError(f"Missing parameters: {missing}") 332 333 references = references or {} # None -> {} 334 335 state = dict(parameters=parameters, references=references) 336 state.update(meta) # meta never has keys "parameters" or "references" 337 338 cls._registry[name] = state 339 340 341doc_components = """ 342 x : `~astropy.units.Quantity`, optional 343 Cartesian, Galactocentric :math:`x` position component. 344 y : `~astropy.units.Quantity`, optional 345 Cartesian, Galactocentric :math:`y` position component. 346 z : `~astropy.units.Quantity`, optional 347 Cartesian, Galactocentric :math:`z` position component. 348 349 v_x : `~astropy.units.Quantity`, optional 350 Cartesian, Galactocentric :math:`v_x` velocity component. 351 v_y : `~astropy.units.Quantity`, optional 352 Cartesian, Galactocentric :math:`v_y` velocity component. 353 v_z : `~astropy.units.Quantity`, optional 354 Cartesian, Galactocentric :math:`v_z` velocity component. 355""" 356 357doc_footer = """ 358 Other parameters 359 ---------------- 360 galcen_coord : `ICRS`, optional, keyword-only 361 The ICRS coordinates of the Galactic center. 362 galcen_distance : `~astropy.units.Quantity`, optional, keyword-only 363 The distance from the sun to the Galactic center. 364 galcen_v_sun : `~astropy.coordinates.representation.CartesianDifferential`, `~astropy.units.Quantity` ['speed'], optional, keyword-only 365 The velocity of the sun *in the Galactocentric frame* as Cartesian 366 velocity components. 367 z_sun : `~astropy.units.Quantity` ['length'], optional, keyword-only 368 The distance from the sun to the Galactic midplane. 369 roll : `~astropy.coordinates.Angle`, optional, keyword-only 370 The angle to rotate about the final x-axis, relative to the 371 orientation for Galactic. For example, if this roll angle is 0, 372 the final x-z plane will align with the Galactic coordinates x-z 373 plane. Unless you really know what this means, you probably should 374 not change this! 375 376 Examples 377 -------- 378 379 To transform to the Galactocentric frame with the default 380 frame attributes, pass the uninstantiated class name to the 381 ``transform_to()`` method of a `~astropy.coordinates.SkyCoord` object:: 382 383 >>> import astropy.units as u 384 >>> import astropy.coordinates as coord 385 >>> c = coord.SkyCoord(ra=[158.3122, 24.5] * u.degree, 386 ... dec=[-17.3, 81.52] * u.degree, 387 ... distance=[11.5, 24.12] * u.kpc, 388 ... frame='icrs') 389 >>> c.transform_to(coord.Galactocentric) # doctest: +FLOAT_CMP 390 <SkyCoord (Galactocentric: galcen_coord=<ICRS Coordinate: (ra, dec) in deg 391 (266.4051, -28.936175)>, galcen_distance=8.122 kpc, galcen_v_sun=(12.9, 245.6, 7.78) km / s, z_sun=20.8 pc, roll=0.0 deg): (x, y, z) in kpc 392 [( -9.43489286, -9.40062188, 6.51345359), 393 (-21.11044918, 18.76334013, 7.83175149)]> 394 395 396 To specify a custom set of parameters, you have to include extra keyword 397 arguments when initializing the Galactocentric frame object:: 398 399 >>> c.transform_to(coord.Galactocentric(galcen_distance=8.1*u.kpc)) # doctest: +FLOAT_CMP 400 <SkyCoord (Galactocentric: galcen_coord=<ICRS Coordinate: (ra, dec) in deg 401 (266.4051, -28.936175)>, galcen_distance=8.1 kpc, galcen_v_sun=(12.9, 245.6, 7.78) km / s, z_sun=20.8 pc, roll=0.0 deg): (x, y, z) in kpc 402 [( -9.41284763, -9.40062188, 6.51346272), 403 (-21.08839478, 18.76334013, 7.83184184)]> 404 405 Similarly, transforming from the Galactocentric frame to another coordinate frame:: 406 407 >>> c = coord.SkyCoord(x=[-8.3, 4.5] * u.kpc, 408 ... y=[0., 81.52] * u.kpc, 409 ... z=[0.027, 24.12] * u.kpc, 410 ... frame=coord.Galactocentric) 411 >>> c.transform_to(coord.ICRS) # doctest: +FLOAT_CMP 412 <SkyCoord (ICRS): (ra, dec, distance) in (deg, deg, kpc) 413 [( 88.22423301, 29.88672864, 0.17813456), 414 (289.72864549, 49.9865043 , 85.93949064)]> 415 416 Or, with custom specification of the Galactic center:: 417 418 >>> c = coord.SkyCoord(x=[-8.0, 4.5] * u.kpc, 419 ... y=[0., 81.52] * u.kpc, 420 ... z=[21.0, 24120.0] * u.pc, 421 ... frame=coord.Galactocentric, 422 ... z_sun=21 * u.pc, galcen_distance=8. * u.kpc) 423 >>> c.transform_to(coord.ICRS) # doctest: +FLOAT_CMP 424 <SkyCoord (ICRS): (ra, dec, distance) in (deg, deg, kpc) 425 [( 86.2585249 , 28.85773187, 2.75625475e-05), 426 (289.77285255, 50.06290457, 8.59216010e+01)]> 427 428""" 429 430 431@format_doc(base_doc, components=doc_components, footer=doc_footer) 432class Galactocentric(BaseCoordinateFrame): 433 r""" 434 A coordinate or frame in the Galactocentric system. 435 436 This frame allows specifying the Sun-Galactic center distance, the height of 437 the Sun above the Galactic midplane, and the solar motion relative to the 438 Galactic center. However, as there is no modern standard definition of a 439 Galactocentric reference frame, it is important to pay attention to the 440 default values used in this class if precision is important in your code. 441 The default values of the parameters of this frame are taken from the 442 original definition of the frame in 2014. As such, the defaults are somewhat 443 out of date relative to recent measurements made possible by, e.g., Gaia. 444 The defaults can, however, be changed at runtime by setting the parameter 445 set name in `~astropy.coordinates.galactocentric_frame_defaults`. 446 447 The current default parameter set is ``"pre-v4.0"``, indicating that the 448 parameters were adopted before ``astropy`` version 4.0. A regularly-updated 449 parameter set can instead be used by setting 450 ``galactocentric_frame_defaults.set ('latest')``, and other parameter set 451 names may be added in future versions. To find out the scientific papers 452 that the current default parameters are derived from, use 453 ``galcen.frame_attribute_references`` (where ``galcen`` is an instance of 454 this frame), which will update even if the default parameter set is changed. 455 456 The position of the Sun is assumed to be on the x axis of the final, 457 right-handed system. That is, the x axis points from the position of 458 the Sun projected to the Galactic midplane to the Galactic center -- 459 roughly towards :math:`(l,b) = (0^\circ,0^\circ)`. For the default 460 transformation (:math:`{\rm roll}=0^\circ`), the y axis points roughly 461 towards Galactic longitude :math:`l=90^\circ`, and the z axis points 462 roughly towards the North Galactic Pole (:math:`b=90^\circ`). 463 464 For a more detailed look at the math behind this transformation, see 465 the document :ref:`astropy:coordinates-galactocentric`. 466 467 The frame attributes are listed under **Other Parameters**. 468 """ 469 470 default_representation = r.CartesianRepresentation 471 default_differential = r.CartesianDifferential 472 473 # frame attributes 474 galcen_coord = CoordinateAttribute(frame=ICRS) 475 galcen_distance = QuantityAttribute(unit=u.kpc) 476 477 galcen_v_sun = DifferentialAttribute( 478 allowed_classes=[r.CartesianDifferential]) 479 480 z_sun = QuantityAttribute(unit=u.pc) 481 roll = QuantityAttribute(unit=u.deg) 482 483 def __init__(self, *args, **kwargs): 484 # Set default frame attribute values based on the ScienceState instance 485 # for the solar parameters defined above 486 default_params = galactocentric_frame_defaults.get() 487 self.frame_attribute_references = \ 488 galactocentric_frame_defaults.references.copy() 489 490 for k in default_params: 491 if k in kwargs: 492 # If a frame attribute is set by the user, remove its reference 493 self.frame_attribute_references.pop(k, None) 494 495 # Keep the frame attribute if it is set by the user, otherwise use 496 # the default value 497 kwargs[k] = kwargs.get(k, default_params[k]) 498 499 super().__init__(*args, **kwargs) 500 501 @classmethod 502 def get_roll0(cls): 503 """ 504 The additional roll angle (about the final x axis) necessary to align 505 the final z axis to match the Galactic yz-plane. Setting the ``roll`` 506 frame attribute to -this method's return value removes this rotation, 507 allowing the use of the `Galactocentric` frame in more general contexts. 508 """ 509 # note that the actual value is defined at the module level. We make at 510 # a property here because this module isn't actually part of the public 511 # API, so it's better for it to be accessible from Galactocentric 512 return _ROLL0 513 514# ICRS to/from Galactocentric -----------------------> 515 516 517def get_matrix_vectors(galactocentric_frame, inverse=False): 518 """ 519 Use the ``inverse`` argument to get the inverse transformation, matrix and 520 offsets to go from Galactocentric to ICRS. 521 """ 522 # shorthand 523 gcf = galactocentric_frame 524 525 # rotation matrix to align x(ICRS) with the vector to the Galactic center 526 mat1 = rotation_matrix(-gcf.galcen_coord.dec, 'y') 527 mat2 = rotation_matrix(gcf.galcen_coord.ra, 'z') 528 # extra roll away from the Galactic x-z plane 529 mat0 = rotation_matrix(gcf.get_roll0() - gcf.roll, 'x') 530 531 # construct transformation matrix and use it 532 R = matrix_product(mat0, mat1, mat2) 533 534 # Now need to translate by Sun-Galactic center distance around x' and 535 # rotate about y' to account for tilt due to Sun's height above the plane 536 translation = r.CartesianRepresentation(gcf.galcen_distance * [1., 0., 0.]) 537 z_d = gcf.z_sun / gcf.galcen_distance 538 H = rotation_matrix(-np.arcsin(z_d), 'y') 539 540 # compute total matrices 541 A = matrix_product(H, R) 542 543 # Now we re-align the translation vector to account for the Sun's height 544 # above the midplane 545 offset = -translation.transform(H) 546 547 if inverse: 548 # the inverse of a rotation matrix is a transpose, which is much faster 549 # and more stable to compute 550 A = matrix_transpose(A) 551 offset = (-offset).transform(A) 552 offset_v = r.CartesianDifferential.from_cartesian( 553 (-gcf.galcen_v_sun).to_cartesian().transform(A)) 554 offset = offset.with_differentials(offset_v) 555 556 else: 557 offset = offset.with_differentials(gcf.galcen_v_sun) 558 559 return A, offset 560 561 562def _check_coord_repr_diff_types(c): 563 if isinstance(c.data, r.UnitSphericalRepresentation): 564 raise ConvertError("Transforming to/from a Galactocentric frame " 565 "requires a 3D coordinate, e.g. (angle, angle, " 566 "distance) or (x, y, z).") 567 568 if ('s' in c.data.differentials and 569 isinstance(c.data.differentials['s'], 570 (r.UnitSphericalDifferential, 571 r.UnitSphericalCosLatDifferential, 572 r.RadialDifferential))): 573 raise ConvertError("Transforming to/from a Galactocentric frame " 574 "requires a 3D velocity, e.g., proper motion " 575 "components and radial velocity.") 576 577 578@frame_transform_graph.transform(AffineTransform, ICRS, Galactocentric) 579def icrs_to_galactocentric(icrs_coord, galactocentric_frame): 580 _check_coord_repr_diff_types(icrs_coord) 581 return get_matrix_vectors(galactocentric_frame) 582 583 584@frame_transform_graph.transform(AffineTransform, Galactocentric, ICRS) 585def galactocentric_to_icrs(galactocentric_coord, icrs_frame): 586 _check_coord_repr_diff_types(galactocentric_coord) 587 return get_matrix_vectors(galactocentric_coord, inverse=True) 588 589 590# Create loopback transformation 591frame_transform_graph._add_merged_transform(Galactocentric, ICRS, Galactocentric) 592