GNU Astronomy Utilities


Next: , Previous: , Up: General program usage tutorial   [Contents][Index]


2.2.5 Angular coverage on the sky

This is the deepest image we currently have of the sky. The first thing that comes to mind may be this: “How large is this field on the sky?”. You can get a fast and crude answer with Gnuastro’s Fits program using this command:

astfits flat-ir/xdf-f160w.fits --skycoverage

It will print the sky coverage in two formats (all numbers are in units of degrees for this image): 1) the image’s central RA and Dec and full width around that center, 2) the range of RA and Dec covered by this image. You can use these values in various online query systems. You can also use this option to automatically calculate the area covered by this image. With the --quiet option, the printed output of --skycoverage will not contain human-readable text, making it easier for further processing:

astfits flat-ir/xdf-f160w.fits --skycoverage --quiet

The second row is the coverage range along RA and Dec (compare with the outputs before using --quiet). We can thus simply subtract the second from the first column and multiply it with the difference of the fourth and third columns to calculate the image area. We’ll also multiply each by 60 to have the area in arc-minutes squared.

astfits flat-ir/xdf-f160w.fits --skycoverage --quiet \
        | awk 'NR==2{print ($2-$1)*60*($4-$3)*60}'

The returned value is \(9.06711\) arcmin\(^2\). However, this method ignores the fact many of the image pixels are blank! In other words, the image does cover this area, but there is no data in more than half of the pixels. So let’s calculate the area coverage over-which we actually have data.

The FITS world coordinate system (WCS) meta data standard contains the key to answering this question. Run the following command to see all the FITS keywords (metadata) for one of the images (almost identical with the other images because they were are scaled to the same region of Sky):

astfits flat-ir/xdf-f160w.fits -h1

Look into the keywords grouped under the ‘World Coordinate System (WCS)’ title. These keywords define how the image relates to the outside world. In particular, the CDELT* keywords (or CDELT1 and CDELT2 in this 2D image) contain the “Coordinate DELTa” (or change in coordinate units) with a change in one pixel. But what is the units of each “world” coordinate? The CUNIT* keywords (for “Coordinate UNIT”) have the answer. In this case, both CUNIT1 and CUNIT1 have a value of deg, so both “world” coordinates are in units of degrees. We can thus conclude that the value of CDELT* is in units of degrees-per-pixel31.

With the commands below, we’ll use CDELT (along with the image size) to find the answer of our initial question: “how much of the sky does this image cover?”. The lines starting with ## are just comments for you to read and understand each command. Don’t type them on the terminal. The commands are intentionally repetitive in some places to better understand each step and also to demonstrate the beauty of command-line features like history, variables, pipes and loops (which you will commonly use as you master the command-line).

Use shell history: Don’t forget to make effective use of your shell’s history: you don’t have to re-type previous command to add something to them. This is especially convenient when you just want to make a small change to your previous command. Press the “up” key on your keyboard (possibly multiple times) to see your previous command(s) and modify them accordingly.

## If your system language uses ',' (not '.') as decimal separator.
$ export LANG=C

## See the general statistics of non-blank pixel values.
$ aststatistics flat-ir/xdf-f160w.fits

## We only want the number of non-blank pixels.
$ aststatistics flat-ir/xdf-f160w.fits --number

## Keep the result of the command above in the shell variable `n'.
$ n=$(aststatistics flat-ir/xdf-f160w.fits --number)

## See what is stored the shell variable `n'.
$ echo $n

## Show all the FITS keywords of this image.
$ astfits flat-ir/xdf-f160w.fits -h1

## The resolution (in degrees/pixel) is in the `CDELT' keywords.
## Only show lines that contain these characters, by feeding
## the output of the previous command to the `grep' program.
$ astfits flat-ir/xdf-f160w.fits -h1 | grep CDELT

## Since the resolution of both dimensions is (approximately) equal,
## we'll only use one of them (CDELT1).
$ astfits flat-ir/xdf-f160w.fits -h1 | grep CDELT1

## To extract the value (third token in the line above), we'll
## feed the output to AWK. Note that the first two tokens are
## `CDELT1' and `='.
$ astfits flat-ir/xdf-f160w.fits -h1 | grep CDELT1 | awk '{print $3}'

## Save it as the shell variable `r'.
$ r=$(astfits flat-ir/xdf-f160w.fits -h1 | grep CDELT1   \
              | awk '{print $3}')

## Print the values of `n' and `r'.
$ echo $n $r

## Use the number of pixels (first number passed to AWK) and
## length of each pixel's edge (second number passed to AWK)
## to estimate the area of the field in arc-minutes squared.
$ echo $n $r | awk '{print $1 * ($2*60)^2}'

The output of the last command (area of this field) is 4.03817 (or approximately 4.04) arc-minutes squared. Just for comparison, this is roughly 175 times smaller than the average moon’s angular area (with a diameter of 30arc-minutes or half a degree).

Some FITS writers don’t use the CDELT convention, making it hard to use the steps above. In such cases, you can extract the pixel scale with the --pixelscale option of Gnuastro’s Fits program like the command below. Like the --skycoverage option above, you can also use the --quiet option to allow easy usage of the values in scripts.

$ astfits flat-ir/xdf-f160w.fits --pixelscale

AWK for table/value processing: As you saw above AWK is a powerful and simple tool for text processing. You will see it often in shell scripts. GNU AWK (the most common implementation) comes with a free and wonderful book in the same format as this book which will allow you to master it nicely. Just like this manual, you can also access GNU AWK’s manual on the command-line whenever necessary without taking your hands off the keyboard. Just run info awk.

Your locale doesn’t use ‘.’ as decimal separator: the input/output of some core operating system tools like awk or seq depend on the system locale. For example in Spanish and some other languages the decimal separator (symbol used to separate the integer and fractional part of a number), is a comma. Therefore in systems that have Spanish as their default Locale, seq will print half of unity as ‘0,5’ (instead of ‘0.5’). This can cause problems for parsing the printed numbers in other programs. You can check your current locale with the locale command. You can test your default decimal separator with this command:

seq 0.5 1

To avoid these kinds of locale-specific problems (for example another program not being able to read ‘0,5’ as half of unity), you can change the locale by setting the LANG environment variable (or the lower-level/generic LC_ALL). You can do it only for a single command (the first one below), or all commands within the running session (the second command below):

## Change the locale to the standard, only for this 'seq' command.
$ LANG=C seq 0.5 1

## Change the locale to the standard, for all commands after it.
$ export LANG=C

If you want to change it generally for all future sessions, you can put the second command in your shell’s startup file. For more on startup files, please see Installation directory.


Footnotes

(31)

With the FITS CDELT convention, rotation (PC or CD keywords) and scales (CDELT) are separated. In the FITS standard the CDELT keywords are optional. When CDELT keywords aren’t present, the PC matrix is assumed to contain both the coordinate rotation and scales. Note that not all FITS writers use the CDELT convention. So you might not find the CDELT keywords in the WCS meta data of some FITS files. However, all Gnuastro programs (which use the default FITS keyword writing format of WCSLIB) write their output WCS with the CDELT convention, even if the input doesn’t have it. If your dataset doesn’t use the CDELT convention, you can feed it to any (simple) Gnuastro program (for example Arithmetic) and the output will have the CDELT keyword. See Section 8 of the FITS standard for more


Next: Cosmological coverage, Previous: Dataset inspection and cropping, Up: General program usage tutorial   [Contents][Index]