Premature broadcasting of coordinates in transforms creates slowdowns

@mhvk

Description

In summary, when transforming coordinate frames we broadcast before transforming, which ought not to be necessary.

In my PR to add a topocentric frame I noticed some odd behaviour. Suppose obstime is a (N, 1) array of times and we have M targets. Below I create an alt-az frame intended to broadcast against each other:

import numpy as np
from astropy.time import Time
import astropy.units as u
from astropy.coordinates import EarthLocation, CIRS, AltAz, SkyCoord
from astropy.coordinates.tests.utils import randomly_sample_sphere
from astropy.coordinates.erfa_astrom import ErfaAstromInterpolator, erfa_astrom

times = Time('2020-08-21T00:00') + np.linspace(-3*u.hour, 3*u.hour, 1000)
location = EarthLocation.of_site('lapalma')
aa_frame = AltAz(obstime=times[:, np.newaxis], location=location)

# 1000 random locations on the sky
ra, dec, _ = randomly_sample_sphere(200)
coos = SkyCoord(ra, dec)

Now, if I convert this to AltAz it is extremely slow. But it's much faster if I add CIRS as an intermediate step explicitly, even though CIRS is the intermediate frame between ICRS and AltAz...

coos.transform_to(aa_frame) # 18 seconds
coos.transform_to(CIRS()).transform_to(aa_frame) # 1 second

Note this isn't like this in the current master, but it illustrates the issue, which I believe is this. If you run coos.transform_to(aa_frame) the first step is to make a CIRS frame. The attributes of this frame are collected from the players in the transform, so it picks up the obstime from the AltAz frame and the positions from the ICRS frame, so the resulting CIRS frame is shape (N, M). This large frame is pushed through a CIRS->CIRS' which causes the slowdown.

In the second instance the CIRS->CIRS' transform is only for a coordinate of shape (M,) which is much faster.

Although this specific example doesn't affect astropy master, it will if my PR is merged, and @mhvk points out it probably also slows down a single ICRS position -> AltAz at many times transform.

Additional context

My question here is why is this CIRS->CIRS' transformation happening? It seems to me that if the source positions are picked up from the ICRS frame and the observer's location and the time of the observation are picked up from the AltAz frame, then, for each source position, the initial transform to CIRS will take into account the time of observation and the location of the observer so that the resulting CIRS frame for each source will be ready for transformation to AltAz.

It's intricately tied up with how astropy handles coordinate transforms - there's a big chain of transforms and the code doing the transformations "knows" the attributes of the various steps in the chain.

Ideally it would work as you suggest, but there are two issues:

The transformation only knows the attributes of the input coordinates and the destination frame. In this case it knows the times and the location - but it doesn't know how to use the location to create a topocentric CIRS frame (see Cannot directly transform output from GCRS to AltAz / need ITRS location? #10996). So the initial CIRS coord that gets made is geocentric.
The transformation code "knows" it will end up creating a (N, M) shaped coordinate, so it does that in the first step, making a (N, M) geocentric coordinate.

Then the code runs a check of the (N, M) obstime against the time in the altaz frame and the observatory location in the CIRS frame against the (M, 1) time of the altaz frame to see if it needs to run a conversion, sees that it does, and performs an extremely slow transform of the (N, M) geocentric CIRS frame.

There are clearly smarter ways of doing this, but actually working out how without breaking things is hard.

Assuming that it knows that its location attributes are the GCRS position and velocity of the observer, how is it that it doesn't know to call location.get_gcrs_posvel() to obtain these?

I don't think that you need to get overly general here. Given that CIRS is the designated frame for conversion to observed coordinates and, for observed coordinates, the location is specified to be an EarthLocation, CIRS should know how to obtain its location attributes from an EarthLocation.

Another possibility would be to add obsgeoloc and obsgeovel as attributes to the AltAz frame. These could be derived from the EarthLocation attribute when the frame is instantiated and passed directly to the CIRS frame instead of the EarthLocation itself.

@mkbrewer

@mkbrewer - yes, I have been thinking about adding obsgeoloc and obsgeovel as well - that seems somewhat cleaner. But I'm still thinking also about my suggestion of using a frame instance. Probably best to change the API only once!

Mostly just to keep track of where the problem originates: it is definitely inside the execution of the transformation graph, not in .transform_to(...):

from astropy.coordinates import frame_transform_graph
icrs_to_altaz = frame_transform_graph.get_transform(coos.frame.__class__, AltAz)
icrs_to_cirs = frame_transform_graph.get_transform(coos.frame.__class__, CIRS)
cirs_to_altaz = frame_transform_graph.get_transform(CIRS, AltAz)
assert icrs_to_altaz.transforms[0] is icrs_to_cirs.transforms[0]
assert icrs_to_altaz.transforms[1] is cirs_to_altaz.transforms[0]
%timeit icrs_to_altaz(coos.frame, aa_frame)
# slow
%timeit cirs_to_altaz(icrs_to_cirs(coos.frame, CIRS()), aa_frame)
# faster

Specifically,

astropy/astropy/coordinates/transformations.py

Lines 1404 to 1422 in ab1a1d3

    
           def __call__(self, fromcoord, toframe): 
        
               curr_coord = fromcoord 
        
               for t in self.transforms: 
        
                   # build an intermediate frame with attributes taken from either 
        
                   # `fromframe`, or if not there, `toframe`, or if not there, use 
        
                   # the defaults 
        
                   # TODO: caching this information when creating the transform may 
        
                   # speed things up a lot 
        
                   frattrs = {} 
        
                   for inter_frame_attr_nm in t.tosys.get_frame_attr_names(): 
        
                       if hasattr(toframe, inter_frame_attr_nm): 
        
                           attr = getattr(toframe, inter_frame_attr_nm) 
        
                           frattrs[inter_frame_attr_nm] = attr 
        
                       elif hasattr(fromcoord, inter_frame_attr_nm): 
        
                           attr = getattr(fromcoord, inter_frame_attr_nm) 
        
                           frattrs[inter_frame_attr_nm] = attr 
        
                   curr_toframe = t.tosys(**frattrs) 
        
                   curr_coord = t(curr_coord, curr_toframe)

p.s. Oddly, the coordinates are not exactly identical either (off by ~4e-13)

Since this doesn't cause a problem in current master, it is the new presence of the observer location that is actually causing the trouble. Possibly, we should keep the geocentric CIRS so that transformations can be written differently?

EDIT: am still a bit confused: the position does not get passed on (I thought), so why for 100 positions and 100 times, is the difference like this:

In [35]: %timeit coos.frame.transform_to(CIRS(obstime=times)).transform_to(aa_frame)
120 ms ± 2.58 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)

In [36]: %timeit coos.frame.transform_to(CIRS()).transform_to(aa_frame)
88.6 ms ± 1.14 ms per loop (mean ± std. dev. of 7 runs, 10 loops each)

In [37]: %timeit coos.frame.transform_to(aa_frame)
1.44 s ± 17.5 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)

@mhvk

Hi @mhvk

I found it illuminating to print the shape of cirs_coo as the first line of cirs_to_altaz.

Maybe there’s a smarter way of writing the transforms in my topocentric CIRS branch, but you can’t handle aberration properly unless you make a topocentric CIRS coordinate.

You can keep a geocentric only CIRS frame in astropy if you write the ICRS->AltAz transform and remove the CIRS->AltAz one (because it’s innacurate). That avoids this issue.

It does rather beg the question of why have the CIRS frame at all. You couldn’t use a geocentric CIRS frame to point a telescope accurately for example.

For example, if I add print statements to the conversion routines in play I get this:

In [48]: _ = coos.transform_to(CIRS()).transform_to(aa_frame)                                                                                    
icrs_to_cirs:  (200,)
cirs_to_altaz:  (200,)
cirs_to_icrs:  (200,)
icrs_to_cirs:  (200,)

In [49]: _ = coos.transform_to(aa_frame)                                                                                                         
icrs_to_cirs:  (200,)
cirs_to_altaz:  (100, 200)
cirs_to_icrs:  (100, 200)
icrs_to_cirs:  (100, 200)

In terms of where the time is being spent:

    24                                           def cirs_to_altaz(cirs_coo, altaz_frame):
    25         1         53.0     53.0      0.0      print('cirs_to_altaz: ', cirs_coo.shape)
    26         1     202211.0 202211.0      4.5      obsgeoloc, obsgeovel = altaz_frame.location.get_gcrs_posvel(altaz_frame.obstime)
    27         1       1461.0   1461.0      0.0      if (np.any(cirs_coo.obstime != altaz_frame.obstime) or
    28         1      57662.0  57662.0      1.3              np.any(cirs_coo.obsgeoloc != obsgeoloc)):
    29         1         89.0     89.0      0.0          cirs_coo = cirs_coo.transform_to(CIRS(obstime=altaz_frame.obstime,
    30         1          2.0      2.0      0.0                                                obsgeoloc=obsgeoloc,
    31         1    4205034.0 4205034.0     93.8                                                obsgeovel=obsgeovel))

That transform calls cirs_to_icrs and then icrs_to_cirs, the bulk of the time is spent in the first call:

    67                                           def cirs_to_icrs(cirs_coo, icrs_frame):
    68                                               # set up the astrometry context for ICRS<->cirs and then convert to
    69                                               # astrometric coordinate direction
    70         1         44.0     44.0      0.0      print('cirs_to_icrs: ', cirs_coo.shape)
    71         1    4264289.0 4264289.0     99.3      astrom = erfa_astrom.get().apci(cirs_coo)

Navigation Menu

Repository navigation

Premature broadcasting of coordinates in transforms creates slowdowns #10997

Description

Additional context

Metadata

Assignees

Labels

Type

Projects

Milestone

Relationships

Development

	def __call__(self, fromcoord, toframe):
	curr_coord = fromcoord
	for t in self.transforms:
	# build an intermediate frame with attributes taken from either
	# `fromframe`, or if not there, `toframe`, or if not there, use
	# the defaults
	# TODO: caching this information when creating the transform may
	# speed things up a lot
	frattrs = {}
	for inter_frame_attr_nm in t.tosys.get_frame_attr_names():
	if hasattr(toframe, inter_frame_attr_nm):
	attr = getattr(toframe, inter_frame_attr_nm)
	frattrs[inter_frame_attr_nm] = attr
	elif hasattr(fromcoord, inter_frame_attr_nm):
	attr = getattr(fromcoord, inter_frame_attr_nm)
	frattrs[inter_frame_attr_nm] = attr

	curr_toframe = t.tosys(**frattrs)
	curr_coord = t(curr_coord, curr_toframe)

Premature broadcasting of coordinates in transforms creates slowdowns #10997

Description

Description

Additional context

Activity

mkbrewer commented on Nov 6, 2020

StuartLittlefair commented on Nov 6, 2020

mkbrewer commented on Nov 6, 2020

mkbrewer commented on Nov 6, 2020

mkbrewer commented on Nov 6, 2020

mhvk commented on Nov 6, 2020

mhvk commented on Nov 6, 2020

mhvk commented on Nov 6, 2020

StuartLittlefair commented on Nov 6, 2020

StuartLittlefair commented on Nov 6, 2020

StuartLittlefair commented on Nov 6, 2020

Metadata

Metadata

Assignees

Labels

Type

Projects

Milestone

Relationships

Development

Issue actions