Skip to content

Add utu.rotation, differential rotation for named arrays - #8

Merged
roytsmart merged 5 commits into
mainfrom
feature/rotation
Sep 7, 2026
Merged

roytsmart merged 5 commits into
mainfrom
feature/rotation

Conversation

@roytsmart

Copy link
Copy Markdown
Contributor

Summary

Adds utu.rotation.rotate(), which differentially rotates helioprojective coordinates to a new time. It is sunpy.coordinates.propagate_with_solar_surface() expressed in named arrays, and it accepts a time per point, so a raster can give one per raster step.

rotated = utu.rotation.rotate(
    position=position,      # na.Cartesian2dVectorArray of helioprojective coordinates
    time=time,              # scalar, or one per raster step
    time_out=epoch,
)

Why it is not a thin wrapper

Calling sunpy directly for this is both incorrect and slow, and most of the module exists to avoid that.

  • A multidimensional obstime silently returns NaN. A 4x3 grid with a matching time array comes back 33% finite; flattened, it is correct throughout. get_earth raises outright on a multidimensional time.
  • A scalar observer against an array obstime gives NaN, so the observer is evaluated at the same times as the frame.
  • An array obstime is far slower than a scalar one, and not because of how many distinct times it holds. _rotation_matrix_reprs_to_reprs takes one vectorised branch for a scalar angle and a Python list comprehension building a 3x3 matrix per point for an array. Measured: 8000 points at one time is 20 ms, and at 8000 times 2.6 s, including when every one of those times is identical.

The interpolation

Since the transform is only ever handed a scalar time, num anchor times span the range of time and each point is interpolated between the two which surround it. Rotation is smooth and nearly linear over the hour a raster takes.

Against an exact per-step calculation for a 400-step raster:

Anchors Max error Cost
exact, per step reference 8.2 s
2 0.010 arcsec 0.22 s
3 (default) 0.004 arcsec 0.25 s
5 0.002 arcsec 0.32 s

An IRIS pixel is 0.33 arcsec, so the default is about one percent of a pixel. A full 216,540-point raster rotates in 0.5 s rather than the roughly three minutes sunpy would spend.

Off-limb points

Material beyond the limb has no surface to be carried along, so sunpy returns NaN and off_disk="nan" passes that on. off_disk="static" returns those points where they were observed, which keeps spicules and prominences in a mosaic instead of discarding them. The two kinds of point then mean different things and nothing marks which is which, so the default stays NaN.

The seam is smaller than it sounds: foreshortening takes the apparent shift from 57 arcsec at disk centre to 6 arcsec at the radius where sunpy stops answering. SphericalScreen is deliberately not offered, since it moves on-disk results by tens of arcseconds.

Tests

Twelve tests covering shape and axis broadcasting, the identity rotation, westward motion with the equator outrunning mid latitudes, a time per raster step giving monotonic shifts (the case a naive wrapper returns NaN for), both off_disk modes, num=1 reducing to the midpoint time, and an accuracy test asserting two anchors stay within a tenth of an IRIS pixel of the exact result with three and five doing no worse.

Adds sunpy~=8.0 as a dependency.

🤖 Generated with Claude Code

`utu.rotation.rotate` carries helioprojective coordinates along the solar
surface to a new time, which is
`sunpy.coordinates.propagate_with_solar_surface` said in named arrays. It
takes a time per point, so a raster may give one per step, and it returns
the vector it was given the shape of.

Calling `sunpy` the obvious way is both wrong and slow here, and most of
this module is about not doing that. An array `obstime` silently returns
NaN, and `get_earth` raises outright on a multidimensional one, so
everything is flattened on the way in. An array `obstime` is also far
slower than a scalar one, and not because of how many distinct times it
holds: `_rotation_matrix_reprs_to_reprs` takes one vectorised branch for a
scalar angle and a Python loop building a 3x3 matrix per point for an
array. Eight thousand points at one time takes 20 ms; at eight thousand
times, 2.6 s, and just as long when every one of those times is equal.

So the transform is only ever handed a scalar time. Rotation is smooth and
nearly linear over the hour a raster takes, so `num` of them span the range
and each point is interpolated between the two which surround it. Against
an exact per-step calculation, two anchors are worth 0.010 arcsec, three
0.004, and five 0.002, where an IRIS pixel is 0.33. A 400-step raster
takes 0.5 s rather than the three minutes `sunpy` would spend on it.

Points beyond the limb have no surface to be carried along. `off_disk`
says whether they come back as NaN, which is what `sunpy` does and what
this does by default, or where they were seen, which keeps the spicules in
a mosaic. The seam that leaves is smaller than it sounds, since
foreshortening takes the apparent shift from 57 arcsec at disk centre to 6
at the radius where `sunpy` stops answering.

Co-Authored-By: Claude Fable 5.1 <noreply@anthropic.com>
@codecov

codecov Bot commented Sep 7, 2026

Copy link
Copy Markdown

Codecov Report

✅ All modified and coverable lines are covered by tests.
✅ Project coverage is 100.00%. Comparing base (43ca22b) to head (3cd1662).

Additional details and impacted files
@@            Coverage Diff             @@
##              main        #8    +/-   ##
==========================================
  Coverage   100.00%   100.00%            
==========================================
  Files           10        13     +3     
  Lines          314       501   +187     
==========================================
+ Hits           314       501   +187     

☔ View full report in Codecov by Harness.
📢 Have feedback on the report? Share it here.

🚀 New features to boost your workflow:
  • ❄️ Test Analytics: Detect flaky tests, report on failures, and find test suite problems.

@read-the-docs-community

read-the-docs-community Bot commented Sep 7, 2026

Copy link
Copy Markdown

…n is.

`num` was three whatever it was given, which is only ever right for one
length of observation. What the interpolation costs in accuracy is set by
how far apart those times are and grows as the square of it, so three of
them are worth 0.004 arcsec over an hour, 0.22 over a day, and 163 across
a disk transit. Nothing says so at the time; the answer just quietly stops
being right.

`num` is now taken from the span of `time` unless it is given, and taken to
be an hour apart, which holds the error between 0.001 and 0.010 arcsec from
one hour to thirteen days. That is inside a pixel of anything likely to ask,
DKIST included. A raster is unaffected, wanting the same three times and the
same half second it wanted before; a day wants twenty-five and a transit
three hundred, which is what that answer costs.

The span is measured by subtracting two Julian dates, which are numbers near
two and a half million and so only good to about forty microseconds between
them. An hour comes out of that as a shade over an hour, and would take a
time it does not need, so the ratio is rounded before the ceiling.
The examples returned a coordinate and a list of numbers, which say that
something moved but not what the motion looks like. They are now the two
pictures the numbers were standing in for, drawn on the disk the points
sit on.

The first takes a grid across the disk a day forward and draws a line from
each point to where it lands. Both of the things which make the motion
awkward are in it: the lines are shortest near the limb, where the movement
is mostly toward the observer and hardly changes where a point appears, and
shorter toward the poles, which turn more slowly than the equator.

The second is the reason the rotation is called differential. A meridian is
straight to begin with and is not straight a day later, and the picture is
that bend, once a day for five days.
`utu.rotation` refers to `sunpy` throughout, and to
`propagate_with_solar_surface` and `SphericalScreen` by name, none of which
were links because `sunpy` was not among the inventories the documentation
is built against. Sphinx does not fail a build over a reference it cannot
resolve, so nothing said so.

Building with `-n` now leaves the module without a warning, and the
examples pick up a link to `angular_radius` from the code as well.
Two things were left over once `sunpy` was added to the inventories.

`matplotlib` and `adjustText` were not among them either, so the axes
`utu.spectrum.stem` draws on, the two methods it draws with, and the solver
it moves the labels apart with were all named without being linked.

The `Returns` of that function was a sentence with no type above it, which
`napoleon` reads as the type rather than the description, and then splits on
its commas: `The label`, `each line that was labelled` and `in the order
they were drawn.` were each looked up as a class. The sentence is now the
description of a return type, which is the one the annotation already gives.

Building with `-n` leaves no reference unresolved anywhere in the
documentation.
@roytsmart
roytsmart merged commit 3121255 into main Sep 7, 2026
13 checks passed
@roytsmart
roytsmart deleted the feature/rotation branch September 7, 2026 19:41
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

1 participant