GNU Astronomy Utilities Coordinate and border operators

The operators here help you in defining or manipulating coordinates. For examples to define the “box” (a rectangular region) that surrounds an ellipse or to rotate a point around a reference point.


Rotate the given point (horizontal and vertical coordinates given in 5th and 4th popped operands) around a center/reference point (coordinates given in the 3rd and 2nd popped operands) by a given angle (first popped operand).

For example, if you want to trace the outer edge of a circle centered on (1.23,45.6) with a radius of 0.78, you can use this operator like below. The logic is that we assume a single point that is located on 0.78 units after the center on the horizontal axis (the point’s vertical axis position is the same as the center). We then rotate this point in each row by one degree to build the circle’s circumference.

$ cx=1.23
$ cy=45.6
$ rad=0.78
$ seq 0 360 \
      | awk '{print '$rad'+'$cx', '$cy', $1}' \
      | asttable -c'arith $1 $2 '$cx' '$cy' $3 rotate-coord' \

##  Within TOPCAT, after opening "Plane Plot", within "Axes" select
## "Aspect lock" so the steps in both axis is the same.
$ astscript-fits-view circle.fits

If you want the points to create a circle on the celestial sphere, you can use the eq-j2000-from-flat operator after this one (see Column arithmetic):

$ seq 0 360 \
      | awk '{print '$rad'+'$cx', '$cy', $1}' \
      | asttable -c'arith $1 $2 '$cx' '$cy' $3 rotate-coord \
                          '$cx' '$cy' TAN eq-j2000-from-flat' \

When you open TOPCAT, if you open the “Plane Plot”, you will see an ellipse. However, if you open “Sky Plot” (from the “Graphics” menu), and select the first and second columns respectively, you will see a circle.

The center coordinates and angle can be fixed for all the rows (as in the example above) or be different for every row. Recall that if you want these to change on every row, you should give the column name (or number followed by $) for these operands instead of the constant number above.


Return the width (along horizontal) and height (along vertical) of a box that encompasses an ellipse with the same center point. The top-popped operand is assumed to be the position angle (angle from the horizontal axis) in degrees. The second and third popped operands are the minor and major radii of the ellipse respectively. This operator outputs two operands on the general stack. The first one is the width and the second (which will be the top one when this operator finishes) is the height.

If the value to the second popped operand (minor axis) is larger than the third (major axis), a NaN value will be written for both the width and height of that element and a warning will be printed (the warning can be disabled with the --quiet option).

As an example, if your ellipse has a major axis radius of 10 units, a minor axis radius of 4 units and a position angle of 20 degrees, you can estimate the bounding box with this command:

$ echo "10 4 20" \
       | asttable -c'arith $1 $2 $3 box-around-ellipse'

Alternatively if your three values are in separate FITS arrays/images, you can use the command below to have the width and height in similarly sized fits arrays. In this example a.fits and b.fits are respectively the major and minor axis lengths and pa.fits is the position angle (in degrees). Also, in all three, we assume the first extension is used. After it is done, the height of the box will be put in h.fits and the width will be in w.fits. Just note that because this operator has two output datasets, you need to first write the height (top output operand) into a file and free it with the tofilefree- operator, then write the width in the file given to --output.

$ astarithmetic a.fits b.fits pa.fits box-around-ellipse \
                tofilefree-h.fits -ow.fits -g1

Finally, if you need to treat the width and height separately for further processing, you can call the set- operator two times afterwards like below. Recall that the set- operator will pop the top operand, and put it in memory with a certain name, bringing the next operand to the top of the stack.

For example, let’s assume catalog.fits has at least three columns MAJOR, MINOR and PA which specify the major axis, minor axis and position angle respectively. But you want the final width and height in 32-bit floating point numbers (not the default 64-bit, which may be too much precision in many scenarios). You can do this with the command below (note you can also break lines with \, within the single-quote environment)

$ asttable catalog.fits \
           -c'arith MAJOR MINOR PA box-around-ellipse \
                    set-height set-width \
                    width float32 height float32'

Convert a box center and width to the coordinates of the vertices of the box on a left-hand spherical coordinate system. In a left-handed spherical coordinate system, the longitude increases towards the left while north is up (as in the RA and Dec direction of the equatorial coordinate system used in astronomy). This operator therefore takes four input operands (the RA and Dec of the box’s center, as well as the width of the box in each direction).

After it is complete, this operator places 8 operands on the stack which contain the RA and Dec of the four vertices of the box in the following anti-clockwise order:

  1. Bottom-left vertice Longitude (RA)
  2. Bottom-left vertice Latitude (Dec)
  3. Bottom-right vertice Longitude (RA)
  4. Bottom-right vertice Latitude (Dec)
  5. Top-right vertice Longitude (RA)
  6. Top-right vertice Latitude (Dec)
  7. Top-left vertice Longitude (RA)
  8. Top-left vertice Latitude (Dec)

For example, with the command below, we will retrieve the vertice coordinates of a rectangle around a point with RA=20 and Dec=0 (on the equator). The rectangle will have a 1 degree edge along the RA direction and a 2 degree edge along the declination. In this example, we are using the -Afixed -B2 only for demonstration purposes here due to the round numbers! In general, it is best to write your outputs to a binary FITS table to preserve the full precision (see Printing floating point numbers).

$ echo "20 0 1 2" \
       | asttable -Afixed -B2 \
                  -c'arith $1 $2 $3 $4 box-vertices-on-sphere'
20.50  -1.00  19.50  -1.00  19.50  1.00   20.50  1.00

We see that the bottom-left vertice is at (RA,Dec) of \((20.50,-1.0)\) and the top-right vertice is at \((19.50,1.00)\). These could have easily been done by manually adding and subtracting! But you will see that the complexity arises at higher/lower declinations. For example, with the command below, let’s see how vertice coordinates of the same box, but after moving its center to (RA,Dec) of (20,85):

$ echo "20 85 1 2" \
       | asttable -Afixed -B2 \
                  -c'arith $1 $2 $3 $4 box-vertices-on-sphere'
24.78  84.00  15.22  84.00  12.83  86.00  27.17  86.00

Even though, we didn’t change the central RA (20) or the size of the box along the RA (1 degree), the RA of the bottom-left vertice is now at 24.78; almost 5 degrees away! This occurs because of the spherical coordinate system, we measure the longitude (e.g., RA) with the following way:

  1. Draw a meridian that passes your point. The meridian is half of a great-circle (which has a diameter that is equal to the sphere’s diameter) passes both poles.
  2. Find the intersection of that meridian with the equator.
  3. The distance of the intersection and the reference point (along the equator) defines the longitude angle.

As you get more distant from the equator (declination becomes non-zero), any change along the RA (towards the east; 1 degree in the example above) will on longer be on a great circle, but along a “small circle”. On a small circle that is defined by the fixed declination \(\delta\), the distance of two points is closer than the distances of their projection on the equator (as described in the definition of longitude above). It is smaller by a factor of \(\cos({\delta})\).

Therefore, an angular change (let’s call it \(\Delta_{lon}\)) along the small circle defined by the fixed declination of \(\delta\) corresponds to \(\Delta_{lon}/\cos(\delta)\) on the equator.