GNU Astronom y Utilities
Astronomical data manipulation and analysis programs and libraries
for v ersion 0.23, 13 July 2024
Mohammad Akhlaghi
This b o ok do cumen ts v ersion 0.23 of the GNU Astronom y Utilities (Gn uastro). Gn uastro
pro vides v arious programs and libraries for astronomical data manipulation and analysis.
Cop yrigh t c
2015-2024 F ree Soft w are F oundation, Inc.
P ermission is gran ted to cop y , distribute and/or mo dify this do cumen t under the
terms of the GNU F ree Do cumen tation License, V ersion 1.3 or an y later v ersion
published b y the F ree Soft w are F oundation; with no In v arian t Sections, no
F ron t-Co v er T exts, and no Bac k-Co v er T exts. A cop y of the license is included
in the section en titled “GNU F ree Do cumen tation License”.
Gn uastro (source co de, b o ok and w eb page) authors (sorted b y n um b er of commits):
Mohammad Akhlaghi ( [email protected] , 2275)
Raul Infan te-Sainz (infan [email protected] , 154)
Sepideh Esk andarlou (sepideh.esk [email protected] , 90)
P edram Ashofteh-Ardak ani (p edramardak [email protected] , 75)
F aezeh Bidjarc hian (fbidjarc [email protected] , 29)
Mos ` e Giordano (mose@gn u.org, 29)
Elham Saremi (saremi elham@y aho o.com, 21)
Vladimir Mark elo v ( [email protected] , 20)
Sac hin Kumar Singh (sac [email protected] , 18)
Zahra Sharbaf ( [email protected] , 13)
Boud Rouk ema (b [email protected] , 11)
Nat´ ali D. Anzanello ( [email protected] , 8)
Giacomo Lorenzetti ( [email protected] , 6)
Jash Shah ( [email protected] , 5)
Samane Ra ji (samanera [email protected] , 4)
Thorsten Alteholz ( [email protected] , 4)
Carlos Morales-So corro (cmorso [email protected] , 3)
Marjan Akbari ( [email protected] , 3)
Th ´ er ` ese Go defro y (go [email protected] , 3)
Zohreh Ghaffari ( [email protected] , 3)
F athma Mehno or (fathmamehno [email protected] , 2)
Joseph Putk o (josephputk [email protected] , 2)
S. Zahra Hosseini Shahisa v andi ( [email protected] , 2)
Alexey Dokuc haev ( [email protected] , 1)
Andreas Stieger ( [email protected] , 1)
Bharat Bhandari ( [email protected] , 1)
F ran¸ cois Oc hsen b ein (francois.o c hsen b [email protected] , 1)
Giulia Golini ( [email protected] , 1)
Kartik Ohri (k artik [email protected] , 1)
Lab eeb Asari (asari.r.lab [email protected] , 1)
Leindert Bo ogaard (leindertb o [email protected] , 1)
Lucas MacQuarrie ( [email protected] , 1)
Madha v Bansal (madha vbansal.cse18@itbh u.ac.in, 1)
Miguel de V al-Borro (miguel.dev [email protected] , 1)
Nafise Sedighi ( [email protected] , 1)
Rashid Y aaqib (r.y [email protected] , 1)
F or m yself, I am in terested in science and in philosoph y only b ecause I w an t to
learn something ab out the riddle of the w orld in whic h w e liv e, and the riddle
of man’s kno wledge of that w orld. And I b eliev e that only a reviv al of in terest
in these riddles can sa v e the sciences and philosoph y from narro w sp ecialization
and from an obscuran tist faith in the exp ert’s sp ecial skill, and in his p ersonal
kno wledge and authorit y; a faith that so w ell fits our ‘p ost-rationalist’ and ‘p ost-
critical’ age, proudly dedicated to the destruction of the tradition of rational
philosoph y , and of rational though t itself.
—Karl P opp er. The logic of scien tific disco v ery . 1959.
i
Short Con ten ts
1 In tro duction .................................. ........ 1
2 T utorials ...................... ...................... 21
3 Installation ................................. ........ 210
4 Common program b eha vior .. .......................... 246
5 Data con tainers ............................ .......... 294
6 Data manipulation ................................... 384
7 Data analysis .................................. ...... 509
8 Data mo deling .................................. ..... 631
9 High-lev el calculations ................................ 656
10 Installed scripts ............................... ....... 668
11 Mak efile extensions (for GNU Mak e) ..................... 720
12 Library .................................. ........... 729
13 Dev eloping ............................ .............. 933
A Other useful soft w are ................................. 964
B GNU F ree Do c. License ............................... 968
C GNU Gen. Pub. License v3 ............................ 976
Index: Macros, structures and functions ..................... 987
Index .................................. ............... 997
ii
T able of Con ten ts
1 In tro duction .................................... .. 1
1.1 Quic k start ..... ................................................ 1
1.2 Gn uastro programs list ......................................... 2
1.3 Gn uastro manifesto: Science and its to ols ....................... 6
1.4 Y our righ ts .................................................... 10
1.5 Logo of Gn uastro .............................................. 10
1.6 Naming con v en tion ............................................ 11
1.7 V ersion n um b ering ............................................ 11
1.7.1 GNU Astronom y Utilities 1.0 .............................. 12
1.8 New to GNU/Lin ux? .......................................... 12
1.8.1 Command-line in terface .................................... 13
1.9 Rep ort a bug .................................................. 15
1.10 Suggest new feature ........................................... 17
1.11 Announcemen ts ............................................ ... 18
1.12 Con v en tions .......... ......................................... 18
1.13 Ac kno wledgmen ts ............................................. 19
2 T utorials ............................ ............. 21
2.1 General program usage tutorial ................................ 22
2.1.1 Calling Gn uastro’s programs .............................. 23
2.1.2 Accessing do cumen tation ................................. 23
2.1.3 Setup and data do wnload ..... ............................ 25
2.1.4 Dataset insp ection and cropping .......................... 25
2.1.5 Angular co v erage on the sky .............................. 28
2.1.6 Cosmological co v erage and visualizing tables .............. 30
2.1.7 Building custom programs with the library ................ 33
2.1.8 Option managemen t and configuration files ............... 35
2.1.9 W arping to a new pixel grid ........................ ...... 37
2.1.10 NoiseChisel and Multi-Extension FITS files ............... 38
2.1.11 NoiseChisel optimization for detection .................... 41
2.1.12 NoiseChisel optimization for storage ...................... 46
2.1.13 Segmen tation and making a catalog ....................... 47
2.1.14 Measuring the dataset limits .............................. 49
2.1.15 W orking with catalogs (estimating colors) ................. 54
2.1.16 Column statistics (color-magnitude diagram) .............. 59
2.1.17 Ap erture photometry ..................................... 61
2.1.18 Matc hing catalogs ........................................ 62
2.1.19 Reddest clumps, cutouts and parallelization ............... 63
2.1.20 FITS images in a publication ............................. 65
2.1.21 Marking ob jects for publication ........................... 69
2.1.22 W riting scripts to automate the steps ............. ........ 73
2.1.23 Citing and ac kno wledging Gn uastro ....................... 80
2.2 Detecting large extended targets ............................... 80
iii
2.2.1 Do wnloading and v alidating input data .................... 81
2.2.2 NoiseChisel optimization ................................... 82
2.2.3 Sk ewness caused b y signal and its measuremen t ............ 88
2.2.4 Image surface brigh tness limit ............................. 92
2.2.5 Ac hiev ed surface brigh tness lev el ........................... 97
2.2.6 Extract clumps and ob jects (Segmen tation) ................ 99
2.3 Building the extended PSF ................................... 102
2.3.1 Preparing input for extended PSF ........................ 102
2.3.2 Saturated pixels and Segmen t’s clumps ................... 103
2.3.3 One ob ject for the whole detection ........................ 107
2.3.4 Building outer part of PSF ............................... 110
2.3.5 Inner part of the PSF .................................... 113
2.3.6 Uniting the differen t PSF comp onen ts .................... 114
2.3.7 Subtracting the PSF ...................................... 117
2.4 Sufi sim ulates a detection ................................. ... 123
2.5 Detecting lines and extracting sp ectra in 3D data ............. 134
2.5.1 Viewing sp ectra and redshifted lines ...................... 134
2.5.2 Sky lines in optical IFUs .................................. 138
2.5.3 Con tin uum subtraction ................................... 138
2.5.4 3D detection with NoiseChisel ............................ 140
2.5.5 3D measuremen ts and sp ectra ............................ 142
2.5.6 Extracting a single sp ectrum and plotting it .............. 146
2.5.7 Pseudo narro w-band images .............................. 148
2.6 Color images with full dynamic range ........................ . 150
2.6.1 Color c hannels in same pixel grid ......................... 151
2.6.2 Color image using linear transformation ................... 153
2.6.3 Color for brigh t regions and gra yscale for fain t ............ 155
2.6.4 Man ually setting color-blac k-gra y regions ................. 160
2.6.5 W eigh ts, con trast, mark ers and other customizations . . . . . . 162
2.7 Zero p oin t of an image ....................................... 165
2.7.1 Zero p oin t tutorial with reference image .................. 165
2.7.2 Zero p oin t tutorial with reference catalog ................. 174
2.8 P oin ting pattern design ...................................... 175
2.8.1 Preparing input and generating exp osure map ............ 176
2.8.2 Area of non-blank pixels on sky ........................... 179
2.8.3 Script with p oin ting sim ulation steps so far ............... 180
2.8.4 Larger steps sizes for b etter calibration ................... 182
2.8.5 P oin tings that accoun t for sky curv ature .................. 185
2.8.6 Accoun ting for non-exp osed pixels ........................ 188
2.9 Moir ´ e pattern in stac king and its correction .................. 190
2.10 Clipping outliers ............................................ . 194
2.10.1 Building inputs and analysis without clipping ............ 195
2.10.2 Sigma clipping .......................................... 199
2.10.3 MAD clipping ........................................... 205
2.10.4 Con tiguous outliers ...................................... 208
iv
3 Installation ..................................... 210
3.1 Dep endencies ............................................ ..... 210
3.1.1 Mandatory dep endencies ................................. 211
3.1.1.1 GNU Scien tific Library .............................. 211
3.1.1.2 CFITSIO ............................................ 211
3.1.1.3 W CSLIB ............................................ 212
3.1.2 Optional dep endencies .................................... 213
3.1.3 Bo otstrapping dep endencies .............................. 216
3.1.4 Dep endencies from pac k age managers ..................... 219
3.2 Do wnloading the source ....................................... 224
3.2.1 Release tarball ........................................... 224
3.2.2 V ersion con trolled source ................................. 225
3.2.2.1 Bo otstrapping ....................................... 226
3.2.2.2 Sync hronizing ....................................... 228
3.3 Build and install .............................................. 229
3.3.1 Configuring ............................................ .. 229
3.3.1.1 Gn uastro configure options .......................... 230
3.3.1.2 Installation directory ................................ 232
3.3.1.3 Executable names ................................... 237
3.3.1.4 Configure and build in RAM ......................... 238
3.3.2 Separate build and source directories ..................... 239
3.3.3 T ests ............................................ ......... 242
3.3.4 A4 prin t b o ok ............................................ 242
3.3.5 Kno wn issues ............................................. 243
4 Common program b eha vior .................. 246
4.1 Command-line ............................................ ... 246
4.1.1 Argumen ts and options ................................... 247
4.1.1.1 Argumen ts .......................................... 248
4.1.1.2 Options ............................................ . 248
4.1.2 Common options ......................................... 250
4.1.2.1 Input/Output options ............................... 251
4.1.2.2 Pro cessing options ................................... 254
4.1.2.3 Op erating mo de options ............................. 256
4.1.3 Shell T AB completion (highly customized) ................ 261
4.1.4 Standard input ........................................... 263
4.1.5 Shell tips ................................................. 264
4.1.5.1 Separate shell v ariables for m ultiple outputs .......... 264
4.1.5.2 T runcating start of long string FITS k eyw ord v alues . . 266
4.2 Configuration files ........................ .................... 267
4.2.1 Configuration file format ................................. 267
4.2.2 Configuration file precedence ............................. 268
4.2.3 Curren t directory and User wide .......................... 269
4.2.4 System wide .............................................. 269
4.3 Getting help ............................................ ..... 270
4.3.1 --usage .................................................. 270
v
4.3.2 --help ............................................ ....... 271
4.3.3 Man pages ............................................ ... 272
4.3.4 Info ............................................ .......... 272
4.3.5 help-gn uastro mailing list ................................. 273
4.4 Multi-threaded op erations .................................... 273
4.4.1 A note on threads ...................... .................. 274
4.4.2 Ho w to run sim ultaneous op erations ...................... 275
4.5 Numeric data t yp es .......................................... 276
4.6 Memory managemen t ........................................ 278
4.7 T ables ............................................ ........... 281
4.7.1 Recognized table formats ................................. 282
4.7.2 Gn uastro text table format ............................... 284
4.7.3 Selecting table columns ................................... 286
4.8 T essellation ............................................ ...... 288
4.9 Automatic output ............................................ 289
4.10 Output FITS files .......... .................................. 290
4.11 Numeric lo cale ............................................ ... 292
5 Data con tainers ................................ 294
5.1 Fits ............................................ ............... 294
5.1.1 In v oking Fits ............................................. 296
5.1.1.1 HDU information and manipulation .................. 298
5.1.1.2 Keyw ord insp ection and manipulation ............... 301
5.1.1.3 Pixel information images ............................. 312
5.2 Con v ertT yp e ............................................ ...... 313
5.2.1 Raster and V ector graphics ................ ............... 313
5.2.2 Recognized file formats ................................... 314
5.2.3 Color .............................................. ....... 317
5.2.3.1 Pixel colors .......................................... 317
5.2.3.2 Colormaps for single-c hannel pixels .................. 318
5.2.3.3 V ector graphics colors ............................... 319
5.2.4 Annotations for figure in pap er ............. .............. 319
5.2.4.1 F ull script of annotations on figure ................... 326
5.2.5 In v oking Con v ertT yp e .................................... 328
5.2.5.1 Con v ertT yp e input and output ...................... 328
5.2.5.2 Pixel visualization ................................... 330
5.2.5.3 Dra wing with v ector graphics ........................ 333
5.3 T able ............................................ ............. 339
5.3.1 Prin ting floating p oin t n um b ers ........................... 340
5.3.2 V ector columns ........................................... 342
5.3.3 Column arithmetic ....................................... 345
5.3.4 Op eration precedence in T able ............................ 353
5.3.5 In v oking T able ........................................... 358
5.4 Query ......................................................... 374
5.4.1 Av ailable databases ...................................... 375
5.4.2 In v oking Query ........................................... 378
xii
13.8.2 Implemen ting T AB completion in Gn uastro .............. 952
13.9 Dev elop er’s c hec klist ........................................ 953
13.10 Gn uastro pro ject w ebpage ................................... 954
13.11 Dev eloping mailing lists ..................................... 955
13.12 Con tributing to Gn uastro ................................... 956
13.12.1 Cop yrigh t assignmen t .................................. 956
13.12.2 Commit guidelines ..................................... 957
13.12.3 Pro duction w orkflo w ................................... 959
13.12.4 F orking tutorial ........................................ 960
App endix A Other useful soft w are ............. 964
A.1 SA O DS9 ......................................... ............ 964
A.2 TOPCA T ............................................ ......... 965
A.3 PGPLOT ..................................................... 966
App endix B GNU F ree Do c. License .......... 968
App endix C GNU Gen. Pub. License v3 ..... 976
Index: Macros, structures and functions ....... 987
Index .................................. ............ 997
1
1 In tro duction
GNU Astronom y Utilities (Gn uastro) is an official GNU pac k age consisting of separate
programs and libraries for the manipulation and analysis of astronomical data. All the pro-
grams share the same basic command-line user in terface for the comfort of b oth the users
and dev elop ers. Gn uastro is written to comply fully with the GNU co ding standards so it
in tegrates finely with the GNU/Lin ux op erating system. This also enables astronomers to
exp ect a fully familiar exp erience in the source co de, building, installing and command-line
user in teraction that they ha v e seen in all the other GNU soft w are that they use. The official
and alw a ys up to date v ersion of this b o ok (or man ual) is freely a v ailable under App endix B
[GNU F ree Do c. License], page 968, in v arious formats (PDF, HTML, plain text, info, and
as its T exinfo source) at http://www.gnu.org/software/gnuastro/manual/ .
F or users who are new to the GNU/Lin ux en vironmen t, unless otherwise sp ecified most of
the topics in Chapter 3 [Installation], page 210, and Chapter 4 [Common program b eha vior],
page 246, are common to all GNU soft w are, for example, installation, managing command-
line options or getting help (also see Section 1.8 [New to GNU/Lin ux?], page 12). So if y ou
are new to this emp o w ering en vironmen t, w e encourage you to go through these c hapters
carefully . They can b e a starting p oin t from whic h y ou can con tin ue to learn more from
eac h program’s o wn man ual and fully b enefit from and enjo y this w onderful en vironmen t.
Gn uastro also comes with a large set of libraries, so y ou can write y our o wn programs using
Gn uastro’s building blo c ks, see Section 12.1 [Review of library fundamen tals], page 729, for
an in tro duction.
In Gn uastro, no c hange to an y program or library will b e committed to its history , b efore
it has b een fully do cumen ted here first. As discussed in Section 1.3 [Gn uastro manifesto:
Science and its to ols], page 6, this is a founding principle of the Gn uastro.
1.1 Quic k start
The latest official release tarball is alw a ys a v ailable as gnuastro-latest.tar.lz ( http://
ftp.gnu.org/gnu/gnuastro/gnuastro-latest.tar.lz ). The Lzip ( http://www.
nongnu.org/lzip/lzip.html ) format is used for b etter compression (smaller output size,
th us faster do wnload), and robust archiv al features and standards. F or historical reasons
(those users that do not y et ha v e Lzip), the Gzip’d tarball 1 is a v ailable at the same URL
(just c hange the .lz suffix ab o v e to .gz ; ho w ev er, the Lzip’d file is recommended). See
Section 3.2.1 [Release tarball], page 224, for more details on the tarball release.
Let’s assume the do wnloaded tarball is in the TOPGNUASTRO directory . Y ou can follo w the
commands b elo w to do wnload and un-compress the Gn uastro source. Y ou need to ha v e the
lzip program for the decompression (see Section 3.1.4 [Dep endencies from pac k age man-
agers], page 219) If y our T ar implemen tation do es not recognize Lzip (the third command
fails), run the fourth command. Note that lines starting with ## do not need to b e t yp ed
(they are only a description of the follo wing command):
## Go into the download directory.
$ cd TOPGNUASTRO
1 The Gzip library and program are commonly a v ailable on most systems. Ho wev er, Gn uastro recommends
Lzip as describ ed ab o v e and the b eta-releases are also only distributed in tar.lz .
Chapter 1: In tro duction 2
## If you do not already have the tarball, you can download it:
$ wget http://ftp.gnu.org/gnu/gnuastro/gnuastro-latest.tar.lz
## If this fails, run the next command.
$ tar -xf gnuastro-latest.tar.lz
## Only when the previous command fails.
$ lzip -cd gnuastro-latest.tar.lz | tar -xf -
Gn uastro has three mandatory dep endencies and some optional dep endencies for extra
functionalit y , see Section 3.1 [Dep endencies], page 210, for the full list. In Section 3.1.4
[Dep endencies from pac k age managers], page 219, we ha v e prepared the command to easily
install Gn uastro’s dep endencies using the pac k age manager of some op erating systems.
When the mandatory dep endencies are ready , y ou can configure, compile, c hec k and install
Gn uastro on y our system with the follo wing commands. See Section 3.3.5 [Kno wn issues],
page 243, if y ou confron t an y complications.
$ cd gnuastro-X.X # Replace X.X with version number.
$ ./configure
$ make -j8 # Replace 8 with no. CPU threads.
$ make check -j8 # Replace 8 with no. CPU threads.
$ sudo make install
F or eac h program there is an ‘In v ok e ProgramName’ sub-section in this b o ok whic h ex-
plains ho w the programs should b e run on the command-line (for example, see Section 5.3.5
[In v oking T able], page 358).
In Chapter 2 [T utorials], page 21, w e ha v e prepared some complete tutorials with com-
mon Gn uastro usage scenarios in astronomical researc h. They ev en con tain links to do wn-
load the necessary data, and thoroughly describ e ev ery step of the pro cess (the science,
statistics and optimal usage of the command-line). W e therefore recommend to read (an
run the commands in) the tutorials b efore starting to use Gn uastro.
1.2 Gn uastro programs list
One of the most common w a ys to op erate Gn uastro is through its command-line programs.
F or some tutorials on sev eral real-w orld usage scenarios, see Chapter 2 [T utorials], page 21.
The list here is just pro vided as a general summary for those who are new to Gn uastro.
GNU Astronom y Utilities 0.23, con tains the follo wing programs. They are sorted in
alphab etical order and a short description is pro vided for eac h program. The description
starts with the executable names in thisfont follo w ed b y a p oin ter to the resp ectiv e section
in paren thesis. Throughout this b o ok, they are ordered based on their context, please see
the top-lev el con ten ts for con textual ordering (based on what they do).
Arithmetic
( astarithmetic , see Section 6.2 [Arithmetic], page 398) F or arithmetic op era-
tions on m ultiple (theoretically unlimited) n um b er of datasets (images). It has
a large and gro wing set of arithmetic, mathematical, and ev en statistical op era-
tors (for example, + , - , * , / , sqrt , log , min , average , median , see Section 6.2.4
[Arithmetic op erators], page 407).
Chapter 1: In tro duction 3
BuildProgram
( astbuildprog , see Section 12.2 [BuildProgram], page 737) Compile, link and
run custom C programs that dep end on the Gn uastro library (see Section 12.3
[Gn uastro library], page 741). This program will automatically link with the
libraries that Gn uastro dep ends on, so there is no need to explicitly mention
them ev ery time y ou are compiling a Gn uastro library dep enden t program.
Con v ertT yp e
( astconvertt , see Section 5.2 [Con v ertT yp e], page 313) Con v ert astronomical
data files (FITS or IMH) to and from sev eral other standard image and data
formats, for example, TXT, JPEG, EPS or PDF. Optionally , it is also p ossible
to add v ector graphics mark ers o v er the output image (for example, circles from
catalogs con taining RA or Dec).
Con v olv e ( astconvolve , see Section 6.3 [Conv olv e], page 470) Con v olv e (blur or smo oth)
data with a giv en k ernel in spatial and frequency domain on m ultiple threads.
Con v olv e can also do decon v olution to find the appropriate k ernel to PSF-matc h
t w o images.
CosmicCalculator
( astcosmiccal , see Section 9.1 [CosmicCalculator], page 656) Do cosmological
calculations, for example, the luminosity distance, distance mo dulus, como ving
v olume and man y more.
Crop ( astcrop , see Section 6.1 [Crop], page 384) Crop region(s) from one or man y
image(s) and stitc h sev eral images if necessary . Input co ordinates can b e in
pixel co ordinates or w orld co ordinates.
Fits ( astfits , see Section 5.1 [Fits], page 294) View and manipulate FITS file
extensions and header k eyw ords.
Mak eCatalog
( astmkcatalog , see Section 7.4 [Mak eCatalog], page 574) Mak e catalog of la-
b eled image (output of NoiseChisel). The catalogs are highly customizable and
adding new calculations/columns is v ery straigh tforw ard, see Akhlaghi 2019
( https://arxiv.org/abs/1611.06387 ).
Mak eProfiles
( astmkprof , see Section 8.1 [Mak eProfiles], page 631) Mak e mo c k 2D profiles
in an image. The cen tral regions of radial profiles are made with a configurable
2D Mon te Carlo in tegration. It can also build the profiles on an ov er-sampled
image.
Matc h ( astmatch , see Section 7.5 [Matc h], page 620) Given t w o input catalogs, find
the ro ws that matc h with eac h other within a giv en ap erture (ma y b e an ellipse).
NoiseChisel
( astnoisechisel , see Section 7.2 [NoiseChisel], page 544) Detect signal in
noise. It uses a tec hnique to detect v ery fain t and diffuse, irregularly shap ed
signal in noise (galaxies in the sky), using thresholds that are b elo w the Sky
v alue, see Akhlaghi and Ic hik aw a 2015 ( http://arxiv.org/abs/1505.01664 )
and Akhlaghi 2019 ( https://arxiv.org/abs/1909.11230 ).
Chapter 1: In tro duction 4
Query ( astquery , see Section 5.4 [Query], page 374) High-lev el in terface to query
pre-defined remote, or external databases, and directly do wnload the required
sub-tables on the command-line.
Segmen t ( astsegment , see Section 7.3 [Segmen t], page 563) Segmen t detected regions
based on the structure of signal and the input dataset’s noise prop erties.
Statistics ( aststatistics , see Section 7.1 [Statistics], page 509) Statistical calculations
on the input dataset (column in a table, image or data cub e). This includes
man op erations suc h as generating histogram, sigma clipping, and least squares
fitting.
T able ( asttable , Section 5.3 [T able], page 339) Con v ert FITS binary and ASCI I
tables in to other suc h tables, prin t them on the command-line, sa v e them in a
plain text file, do arithmetic on the columns or get the FITS table information.
F or a full list of op erations, see Section 5.3.4 [Op eration precedence in T able],
page 353.
W arp ( astwarp , see Section 6.4 [W arp], page 492) W arp image to new pixel grid. By
default it will align the pixel and W CS co ordinates, remo ving an y non-linear
W CS distortions. An y linear w arp (pro jectiv e transformation or Homograph y)
can also b e applied to the input images b y explicitly calling the resp ectiv e
op eration.
The programs listed ab o v e are designed to b e highly mo dular and generic. Hence, they
are naturally for lo w er-lev el op erations. In Gn uastro, higher-lev el op erations (com bining
m ultiple programs, or running a program in a sp ecial wa y), are done with installed Bash
scripts (all prefixed with astscript- ). They can b e run just lik e a program and b eha v e
v ery similarly (with minor differences, see Chapter 10 [Installed scripts], page 668).
astscript-ds9-region
(See Section 10.3 [SA O DS9 region files from table], page 680) Giv en a table
(either as a file or from standard input), create an SA O DS9 region file from
the requested p ositional columns (W CS or image co ordinates).
astscript-fits-view
(see Section 10.4 [Viewing FITS file con ten ts with DS9 or TOPCA T], page 683)
Giv en an y n um b er of FITS files, this script will either op en SA O DS9 (for images
or cub es) or TOPCA T (for tables) to view them in a graphic user in terface
(GUI).
astscript-pointing-simulate
(See Section 10.6 [P oin ting pattern sim ulation], page 692) Giv en a table of
p oin tings on the sky , create and a reference image that con tains y our camera’s
distortions and prop erties, generate a stack ed exp osure map. This is v ery useful
in testing the co v erage of dither patterns when designing y our observing strat-
egy and it is highly customizable. See Akhlaghi 2023 ( https://arxiv.org/
abs/2310.15006 ), or the dedicated tutorial in Section 2.8 [P oin ting pattern
design], page 175.
astscript-radial-profile
(See Section 10.2 [Generate radial profile], page 672) Calculate the 1D radial
profile or 2D p olar plot of an ob ject within an image. The ob ject can b e at any
Chapter 1: In tro duction 5
lo cation in the image, using v arious measures (median, sigma-clipp ed mean,
etc.), and the radial distance can also b e measured on an y general ellipse. See
Infan te-Sainz et al. 2024 ( https://arxiv.org/abs/2401.05303 ) and/or Es-
k andarlou and Akhlaghi 2024 ( https://arxiv.org/abs/2406.14619 ).
astscript-color-faint-gray
(see Section 10.7 [Color images with gra y fain t regions], page 697) Giv en three
images for the Red-Green-Blue (R GB) c hannels, this script will use the brigh t
pixels for color and will sho w the fain t/diffuse regions in gra yscale. This greatly
helps in visualizing the full dynamic range of astronomical data. See Infan te-
Sainz et al. 2024 ( https://arxiv.org/abs/2401.03814 ) or a dedicated tuto-
rial in Section 2.6 [Color images with full dynamic range], page 150.
astscript-sort-by-night
(See Section 10.1 [Sort FITS files b y nigh t], page 669) Giv en a list of FITS files,
and a HDU and k eyw ord name (for a date), this script separates the files in the
same nigh t (p ossibly o v er t w o calendar da ys).
astscript-zeropoint
(see Section 10.5 [Zero p oin t estimation], page 687) Estimate the zero p oin t
(to calibrate pixel v alues) of an input image using a reference image or a ref-
erence catalog. This is necessary to pro duce measuremen ts with ph ysical units
from new images. See Esk andarlou et al. 2023 ( https://arxiv.org/abs/
2312.04263 ), or a dedicated tutorial in Section 2.7 [Zero p oin t of an image],
page 165.
astscript-psf-*
The follo wing scripts are used to estimate the extended PSF estimation and
subtraction as describ ed in the tutorial Section 2.3 [Building the extended PSF],
page 102:
astscript-psf-select-stars
(see Section 10.8.2 [In v oking astscript-psf-select-stars], page 705)
Find all the stars within an image that are suitable for constructing
an extended PSF. If the image has W CS, this script can automat-
ically query Gaia to find the go o d stars.
astscript-psf-stamp
(see Section 10.8.3 [In v oking astscript-psf-stamp], page 708) build a
crop (stamp) of a certain width around a star at a certain co ordinate
in a larger image. This script will do sub-pixel re-p ositioning to
mak e sure the star is cen tered and can optionally mask all other
bac kground sources).
astscript-psf-scale-factor
(see Section 10.8.5 [In v oking astscript-psf-scale-factor], page 714)
Giv en a PSF mo del, and the cen tral co ordinates of a star in an
image, find the scale factor that has to b e multiplied b y the PSF
to scale it to that star.
Chapter 1: In tro duction 6
astscript-psf-unite
(see Section 10.8.4 [In v oking astscript-psf-unite], page 712) Unite
the v arious comp onen ts of a PSF in to one. Because of saturation
and non-linearit y , to get a go o d estimate of the extended PSF, it
is necessary to construct v arious parts from differen t magnitude
ranges.
astscript-psf-subtract
(see Section 10.8.6 [In v oking astscript-psf-subtract], page 717)
Giv en the mo del of a PSF and the cen tral co ordinates of a star in
the image, do sub-pixel re-p ositioning of the PSF, scale it to the
star and subtract it from the image.
1.3 Gn uastro manifesto: Science and its to ols
History of science indicates that there are alw a ys inevitably unseen faults, hidden assump-
tions, simplifications and appro ximations in all our theoretical mo dels, data acquisition and
analysis tec hniques. It is precisely these that will ultimately allo w future generations to
adv ance the existing exp erimen tal and theoretical kno wledge through their new solutions
and corrections.
In the past, scientists w ould gather d ata and pro cess them individually to ac hiev e an
analysis th us ha ving a m uc h more in tricate kno wledge of the data and analysis. The theo-
retical mo dels also required little (if an y) sim ulations to compare with the data. T o da y b oth
metho ds are b ecoming increasingly more dep enden t on pre-written soft w are. Scien tists are
disso ciating themselv es from the in tricacies of reducing ra w observ ational data in exp eri-
men tation or from bringing the theoretical mo dels to life in sim ulations. These ‘in tricacies’
are precisely those unseen faults, hidden assumptions, simplifications and appro ximations
that define scien tific progress.
Unfortunately , most p ersons who hav e recourse to a computer for statistical
analysis of data are not m uc h in terested either in computer programming or in
statistical metho d, b eing primarily concerned with their o wn prop er business.
Hence the common use of library programs and v arious statistical pac k ages. ...
It’s time that w as c hanged.
—F.J. Anscom b e. The American Statistician, V ol. 27, No. 1. 1973
Anscom b e’s quartet ( http://en.wikipedia.org/wiki/Anscombe%27s_quartet )
demonstrates ho w four data sets with widely differen t shap es (when plotted) giv e nearly
iden tical output from standard regression tec hniques. Anscom b e uses this (no w famous)
quartet, whic h w as in tro duced in the pap er quoted ab o v e, to argue that “ Go o d statistic al
analysis is not a pur ely r outine matter, and gener al ly c al ls for mor e than one p ass thr ough
the c omputer ”. Ec hoing Anscom b e’s concern after 44 y ears, some of the highly recognized
statisticians of our time (Leek, McShane, Gelman, Colquhoun, Nuijten and Go o dman),
wrote in Nature that:
W e need to appreciate that data analysis is not purely computational and al-
gorithmic – it is a h uman b eha vior....Researc hers who h un t hard enough will
turn up a result that fits statistical criteria – but their disco v ery will probably
b e a false p ositiv e.
—Fiv e w a ys to fix statistics, Nature, 551, No v 2017.
Chapter 1: In tro duction 7
Users of statistical (scien tific) metho ds (soft w are) are therefore not passiv e (ob jectiv e)
agen ts in their results. It is necessary to actually understand the metho d, not just use it
as a blac k b o x. The sub jectiv e exp erience gained b y frequen tly using a metho d/soft w are is
not sufficien t to claim an understanding of ho w the to ol/metho d w orks and ho w relev an t it
is to the data and analysis. This kind of sub jectiv e exp erience is prone to serious misunder-
standings ab out the data, what the soft w are/statistical-metho d really do es (esp ecially as it
gets more complicated), and th us the scien tific in terpretation of the result. This attitude
is further encouraged through non-free soft w are 2 , p o orly written (or non-existen t) scien tific
soft w are man uals, and non-repro ducible pap ers 3 . This approach to scien tific soft w are and
metho ds only helps in pro ducing dogmas and an “ obscur antist faith in the exp ert’s sp e cial
skil l, and in his p ersonal know le dge and authority ” 4 .
Program or b e programmed. Cho ose the former, and y ou gain access to the
con trol panel of civilization. Cho ose the latter, and it could b e the last real
c hoice y ou get to mak e.
—Douglas Rushk off. Program or b e programmed, O/R Bo oks (2010).
It is ob viously impractical for an y one h uman b eing to gain the in tricate kno wledge
explained ab o v e for ev ery step of an analysis. On the other hand, scien tific data can
b e large and n umerous, for example, images pro duced b y telescop es in astronom y . This
requires efficien t algorithms. T o mak e things w orse, natural scien tists ha v e generally not
b een trained in the adv anced soft w are tec hniques, paradigms and arc hitecture that are
taugh t in computer science or engineering courses and th us used in most soft w are. The
GNU Astronom y Utilities are an effort to tac kle this issue.
Gn uastro is not just a soft w are, this b o ok is as imp ortan t to the idea b ehind Gn uastro as
the source co de (soft w are). This b o ok has tried to learn from the success of the “Numerical
Recip es” b o ok in educating those who are not soft w are engineers and computer scien tists
but still hea vy users of computational algorithms, lik e astronomers. There are t w o ma jor
differences.
The first difference is that Gn uastro’s co de and the bac kground information are segre-
gated: the co de is mo v ed within the actual Gn uastro soft w are source co de and the under-
lying explanations are giv en here in this b o ok. In the source co de, ev ery non-trivial step is
hea vily commen ted and correlated with this b o ok, it follows the same logic of this bo ok, and
all the programs follo w a similar in ternal data, function and file structure, see Section 13.4
[Program source], page 940. Complemen ting the co de, this b o ok fo cuses on thoroughly
explaining the concepts b ehind those co des (history , mathematics, science, soft w are and
usage advice when necessary) along with detailed instructions on ho w to run the programs.
A t the exp ense of frustrating “professionals” or “exp erts”, this b o ok and the commen ts in
the co de also in ten tionally a v oid jargon and abbreviations. The source co de and this b o ok
2 https://www.gnu.org/philosophy/free-sw.html
3 Where the authors omit man y of the analysis/pro cessing “details” from the pap er b y arguing that they
w ould mak e the pap er to o long/unreadable. Ho w ever, soft w are engineers hav e been dealing with such
issues for a long time. There are th us soft ware managemen t solutions that allo w us to supplemen t pap ers
with all the details necessary to exactly repro duce the result. F or example, see Akhlaghi et al. 2021
( https://arxiv.org/abs/2006.03018 ).
4 Karl P opp er. The logic of scien tific disco very . 1959. Larger quote is giv en at the start of the PDF (for
prin t) v ersion of this b o ok.
Chapter 1: In tro duction 8
are th us in timately link ed, and when considered as a single entit y can be thought of as a
real (an actual soft w are accompan ying the algorithms) “Numerical Recip es” for astronom y .
The second ma jor, and arguably more imp ortan t, difference is that “Numerical Recip es”
do es not allo w y ou to distribute an y co de that y ou ha v e learned from it. In other w ords, it
do es not allo w y ou to release y our soft w are’s source co de if y ou ha v e used their co des, y ou can
only publicly release binaries (a blac k b o x) to the comm unit y . Therefore, while it emp o w ers
the privileged individual who has access to it, it exacerbates so cial ignorance. Exactly at
the opp osite end of the sp ectrum, Gnuastro’s source code is released under the GNU general
public license (GPL) and this b o ok is released under the GNU free do cumen tation license.
Y ou are therefore free to distribute an y soft w are y ou create using parts of Gn uastro’s source
co de or text, or figures from this b o ok, see Section 1.4 [Y our righ ts], page 10.
With these principles in mind, Gnuastro’s dev elopers aim to imp ose the minim um re-
quiremen ts on y ou (in computer science, engineering and ev en the mathematics b ehind the
to ols) to understand and mo dify an y step of Gn uastro if y ou feel the need to do so, see
Section 13.1 [Wh y C programming language?], page 933, and Section 13.2 [Program design
philosoph y], page 935.
Without prior familiarit y and exp erience with optics, it is hard to imagine ho w, Galileo
could ha v e come up with the idea of mo difying the Dutc h military telescop e optics to use in
astronom y . Astronomical ob jects could not b e seen with the Dutc h military design of the
telescop e. In other words, it is unlik ely that Galileo could ha v e ask ed a random optician
to mak e mo difications (not understo o d b y Galileo) to the Dutc h design, to do something
no astronomer of the time to ok seriously . In the paradigm of the da y , what could b e the
purp ose of enlarging geometric spheres (planets) or p oin ts (stars)? In that paradigm only
the p osition and mo v emen t of the hea v enly b o dies w as imp ortan t, and that had already
b een accurately studied (recen tly b y T yc ho Brahe).
In the b eginning of his “The Sidereal Messenger” (published in 1610) he cautions the
readers on this issue and b efor e describing his results/observ ations, Galileo instructs us on
ho w to build a suitable instrumen t. Without a detailed description of how he made his to ols
and done his observ ations, no reasonable p erson w ould b eliev e his results. Before he actually
sa w the mo ons of Jupiter, the moun tains on the Mo on or the crescent of V en us, Galileo w as
“ev asiv e” 5 to Kepler. Science is defined b y its to ols/metho ds, not its ra w results 6 .
The same is true to da y: science cannot progress with a blac k b o x, or p o orly released
co de. The source co de of a researc h is the new (abstractified) comm unication language
in science, understandable b y h umans and computers. Source co de (in an y programming
language) is a language/notation designed to express all the details that w ould b e to o
tedious/long/frustrating to rep ort in sp ok en languages lik e English, similar to mathematic
notation.
An article ab out computational science [almost all sciences to da y] ... is not
the sc holarship itself, it is merely advertising of the sc holarship. The Actual
5 Galileo G. (T ranslated by Maurice A. Finocchiaro). The essential Galile o .Hac k ett publishing compan y ,
first edition, 2008.
6 F or example, tak e the follo wing tw o results on the age of the univ erse: roughly 14 billion years (suggested
b y the curren t consensus of the standard mo del of cosmology) and less than 10,000 y ears (suggested from
some in terpretations of the Bible). Both these n umbers are r esults . What distinguishes these t w o results,
is the to ols/metho ds that w ere used to deriv e them. Therefore, as the term “Scientific method” also
signifies, a scien tific statemen t it defined b y its metho d , not its result.
Chapter 1: In tro duction 9
Sc holarship is the complete soft w are dev elopmen t en vironmen t and the complete
set of instructions whic h generated the figures.
—Buc kheit & Donoho, Lecture Notes in Statistics, V ol 103, 1996
T o da y , the quality of the source co de that go es in to a scien tific result (and the distri-
bution of that co de) is as critical to scien tific vitalit y and in tegrit y , as the qualit y of its
written language/English used in publishing/distributing its pap er. A scien tific pap er will
not ev en b e review ed b y an y resp ectable journal if its written in a p o or language/English.
A similar lev el of qualit y assessmen t is th us increasingly b ecoming necessary regarding the
co des/metho ds used to deriv e the results of a scien tific pap er. F or more on this, please see
Akhlaghi et al. 2021 ( https://arxiv.org/abs/2006.03018 )).
Bjarne Stroustrup (creator of the C ++ language) sa ys: “ Without understanding softwar e,
you ar e r e duc e d to b elieving in magic ”. Ken Thomson (the designer or the Unix op erating
system) sa ys “ I abhor a system designe d for the ‘user’ if that wor d is a c o de d p ejor ative
me aning ‘stupid and unsophistic ate d’ .” Certainly no scien tist (user of a scientific soft w are)
w ould w an t to b e considered a b eliev er in magic, or stupid and unsophisticated.
This can happ en when scien tists get to o distant from the ra w data and metho ds, and are
mainly discussing results. In other words, when they feel they ha ve tamed Nature in to their
o wn high-lev el (abstract) mo dels (creations), and are mainly concerned with scaling up, or
industrializing those results. Roughly fiv e y ears b efore sp ecial relativit y , and ab out t wo
decades b efore quan tum mec hanics fundamen tally c hanged Ph ysics, Lord Kelvin is quoted
as sa ying:
There is nothing new to b e disco vered in physics no w. All that remains is more
and more precise measuremen t.
—William Thomson (Lord Kelvin), 1900
A few y ears earlier Alb ert. A. Mic helson made the follo wing statemen t:
The more imp ortan t fundamental la ws and facts of ph ysical science ha v e all b een
disco v ered, and these are no w so firmly established that the p ossibilit y of their
ev er b eing supplan ted in consequence of new discov eries is exceedingly remote....
Our future disco v eries m ust b e lo ok ed for in the sixth place of decimals.
—Alb ert. A. Mic helson, dedication of Ry erson Ph ysics Lab, U. Chicago 1894
If scien tists are considered to b e more than mere puzzle solv ers 7 (simply adding to the
decimals of existing v alues or observing a feature in 10, 100, or 100000 more galaxies or
stars, as Kelvin and Mic helson clearly b eliev ed), they cannot just passiv ely sit bac k and
uncritically rep eat the previous (observ ational or theoretical) metho ds/to ols on new data.
T o da y there is a w ealth of ra w telescop e images ready (mostly for free) at the finger tips
of an y one who is in terested with a fast enough in ternet connection to do wnload them. The
only thing lac king is new w a ys to analyze this data and dig out the treasure that is lying
hidden in them to existing metho ds and tec hniques.
New data that w e insist on analyzing in terms of old ideas (that is, old mo dels
whic h are not questioned) cannot lead us out of the old ideas. Ho wev er man y
data w e record and analyze, we ma y just k eep rep eating the same old errors,
missing the same crucially imp ortan t things that the exp erimen t w as comp eten t
to find.
—Ja ynes, Probabilit y theory , the logic of science. Cam bridge U. Press (2003).
7 Thomas S. Kuhn. The Structur e of Scientific R evolutions , Universit y of Chicago Press, 1962.
Chapter 1: In tro duction 16
Be descriptiv e
Please pro vide as man y details as p ossible and b e v ery descriptiv e. Explain
what y ou exp ected and what the output w as: it migh t b e that y our exp ectation
w as wrong. Also please clearly state whic h sections of the Gn uastro b o ok (this
b o ok), or other references y ou ha v e studied to understand the problem. This
can b e useful in correcting the b o ok (adding links to lik ely places where users
will c hec k). But more imp ortan tly , it will b e encouraging for the dev elop ers,
since y ou are sho wing ho w serious y ou are ab out the problem and that y ou
ha v e actually put some though t in to it. “T o b e able to ask a question clearly
is t w o-thirds of the w a y to getting it answ ered.” – John Ruskin (1819-1900).
Individual and indep enden t bug rep orts
If y ou ha v e found m ultiple bugs, please send them as separate (and indep enden t)
bugs (as m uc h as p ossible). This will significan tly help us in managing and
resolving them so oner.
Repro ducible bug rep orts
If w e cannot exactly repro duce y our bug, then it is v ery hard to resolv e it. So
please send us a Minimal w orking example 17 along with the description. F or
example, in running a program, please send us the full command-line text and
the output with the -P option, see Section 4.1.2.3 [Op erating mo de options],
page 256. If it is caused only for a certain input, also send us that input file.
In case the input FITS is large, please use Crop to only crop the problematic
section and mak e it as small as p ossible so it can easily b e uploaded and do wn-
loaded and not w aste the arc hiv e’s storage, see Section 6.1 [Crop], page 384.
There are generally t w o w a ys to inform us of bugs:
• Send a mail to [email protected] . An y mail y ou send to this address will b e
distributed through the bug-gn uastro mailing list 18 . This is the simplest w ay to send
us bug rep orts. The dev elop ers will then register the bug in to the pro ject w eb page
(next c hoice) for y ou.
• Use the Gn uastro pro ject web page at https://savannah.gnu.org/projects/
gnuastro/ : There are tw o w a ys to get to the submission page as listed b elo w. Fill
in the form as describ ed b elo w and submit it (see Section 13.10 [Gn uastro pro ject
w ebpage], page 954, for more on the pro ject w eb page).
• Using the top horizon tal men u items, immediately under the top page title. Hov-
ering y our mouse on “Supp ort” will op en a drop-do wn list. Select “Submit new”.
Also if y ou hav e an accoun t in Sa v annah, y ou can c ho ose “Bugs” in the men u items
and then select “Submit new”.
• In the main b o dy of the page, under the “Comm unication to ols” section, clic k on
“Submit new item”.
Once the items ha v e b een registered in the mailing list or w eb page, the dev elop ers will
add it to either the “Bug T rac k er” or “T ask Manager” trac k ers of the Gn uastro pro ject
w eb page. These t w o trac k ers can only b e edited b y the Gn uastro pro ject dev elop ers, but
17 http://en.wikipedia.org/wiki/Minimal_Working_Example
18 https://lists.gnu.org/mailman/listinfo/bug-gnuastro
Chapter 1: In tro duction 17
they can b e bro wsed b y an y one, so y ou can follo w the progress on y our bug. Y ou are most
w elcome to join us in dev eloping Gn uastro and fixing the bug y ou ha v e found ma yb e a
go o d starting p oin t. Gn uastro is designed to b e easy for an y one to dev elop (see Section 1.3
[Gn uastro manifesto: Science and its to ols], page 6) and there is a full c hapter dev oted to
dev eloping it: Chapter 13 [Dev eloping], page 933.
Sa v annah’s Markup: When p osting to Sa v annah, it helps to hav e the co de displa y ed in
mono-space fon t and a differen t bac kground, y ou ma y also w an t to mak e a list of items or
mak e some w ords b old. F or features lik e these, y ou should use Sa v annah’s “Markup” guide
at https://savannah.gnu.org/markup-test.php . Y ou can access this page b y clic king
on the “F ull Markup” link that is just b eside the “Preview” button, near the b o x that y ou
write y our commen ts. As y ou see there, for example when you w an t to high-ligh t co de,
y ou should put it within a “ + v erbatim + ” and “-v erbatim-” en vironmen t lik e b elo w:
+verbatim+
astarithmetic image.fits image_arith.fits -h1 isblank nan where
-verbatim-
Unfortunately , Sa v annah do esn’t ha ve a w a y to edit submitted commen ts. Therefore b e
sure to press the “Preview” button and c hec k y our rep ort’s final format b efore the final
submission.
1.10 Suggest new feature
W e w ould alw a ys b e happ y to hear of suggested new features. F or every program, there are
already lists of features that w e are planning to add. Y ou can see the curren t list of plans
from the Gn uastro pro ject web page at https://savannah.gnu.org/projects/gnuastro/
and follo wing “T asks” → “Bro wse” on the horizontal men u at the top of the page immedi-
ately under the title, see Section 13.10 [Gn uastro pro ject w ebpage], page 954. If y ou w an t
to request a feature to an existing program, click on the “Displa y Criteria” ab o v e the list
and under “Category”, c ho ose that particular program. Under “Category” y ou can also see
the existing suggestions for new programs or other cases lik e installation, do cumen tation or
libraries. Also, b e sure to set the “Op en/Closed” v alue to “An y”.
If the feature y ou w an t to suggest is not already listed in the task manager, then follo w
the steps that are fully describ ed in Section 1.9 [Rep ort a bug], page 15. Please ha v e in mind
that the dev elop ers are all busy with their o wn astronomical researc h, and implementing
existing “task”s to add or resolv e bugs. Gn uastro is a v olun teer effort and none of the
dev elop ers are paid for their hard w ork. So, although w e will try our b est, please do
not exp ect for y our suggested feature to b e immediately included (for the next release of
Gn uastro).
The b est p erson to apply the exciting new feature y ou hav e in mind is y ou, since y ou
ha v e the motiv ation and need. In fact, Gn uastro is designed for making it as easy as p ossible
for y ou to hac k in to it (add new features, c hange existing ones and so on), see Section 1.3
[Gn uastro manifesto: Science and its to ols], page 6. Please ha v e a lo ok at the c hapter
dev oted to dev eloping (Chapter 13 [Dev eloping], page 933) and start applying y our desired
feature. Once you ha v e added it, you can use it for y our o wn w ork and if y ou feel you
w an t others to b enefit from y our w ork, y ou can request for it to b ecome part of Gn uastro.
Chapter 1: In tro duction 18
Y ou can then join the dev elop ers and start main taining y our o wn part of Gn uastro. If y ou
c ho ose to tak e this path of action please con tact us b eforehand (Section 1.9 [Rep ort a bug],
page 15) so w e can a v oid p ossible duplicate activities and get in terested p eople in con tact.
Gn uastro is a collection of lo w lev el programs: As describ ed in Section 13.2 [Program
design philosoph y], page 935, a founding principle of Gn uastro is that eac h library or
program should b e basic and lo w-lev el. High lev el jobs should b e done b y running the
separate programs or using separate functions in succession through a shell script or calling
the libraries b y higher level functions, see the examples in Chapter 2 [T utorials], page 21.
So when making the suggestions please consider ho w y our desired job can b est b e brok en
in to separate steps and mo dularized.
1.11 Announcemen ts
Gn uastro has a dedicated mailing list for making announcemen ts ( info-gnuastro ). An yone
can subscrib e to this mailing list. An ytime there is a new stable or test release, an email
will b e circulated there. The email con tains a summary of the o v erall c hanges along with
a detailed list (from the NEWS file). This mailing list is th us the b est w a y to sta y up to
date with new releases, easily learn ab out the up dated/new features, or dep endencies (see
Section 3.1 [Dep endencies], page 210).
T o subscrib e to this list, please visit https://lists.gnu.org/mailman/listinfo/
info-gnuastro . T raffic (n um b er of mails p er unit time) in this list is designed to b e lo w:
only a handful of mails p er y ear. Previous announcemen ts are a v ailable on its arc hiv e
( http://lists.gnu.org/archive/html/info-gnuastro/ ).
1.12 Con v en tions
In this b o ok w e ha v e the follo wing con v en tions:
• All commands that are to b e run on the shell (command-line) prompt as the user start
with a $ . In case they must be run as a sup eruser or system administrator, they will
start with a single # . If the command is in a separate line and next line is also in
the code type face , but do es not ha v e an y of the $ or # signs, then it is the output
of the command after it is run. As a user, y ou do not need to t yp e those lines. A line
that starts with ## is just a commen t for explaining the command to a h uman reader
and m ust not b e t yp ed.
• If the command b ecomes larger than the page width a \ is inserted in the co de. If y ou
are t yping the co de b y hand on the command-line, you do not need to use m ultiple
lines or add the extra space c haracters, so y ou can omit them. If y ou w an t to cop y and
paste these examples (highly discouraged!) then the \ should sta y .
The \ c haracter is a shell escap e c haracter whic h is used commonly to mak e characters
whic h ha v e sp ecial meaning for the shell, lose that sp ecial meaning (the shell will not
treat them esp ecially if there is a \ b ehind them). When \ is the last visible c haracter
in a line (the next c haracter is a new-line c haracter) the new-line c haracter loses its
meaning. Therefore, the shell sees it as a simple white-space c haracter not the end of
a command! This enables y ou to use m ultiple lines to write y our commands.
Chapter 1: In tro duction 19
This is not a con v en tion, but a bi-pro duct of the PDF building pro cess of the man ual: In
the PDF v ersion of this man ual, a single quote (or ap ostrophe) c haracter in the commands
or co des is sho wn like this: ' . Single quotes are sometimes necessary in com bination with
commands lik e awk or sed , or when using Column arithmetic in Gn uastro’s o wn T able
(see Section 5.3.3 [Column arithmetic], page 345). Therefore when t yping (recommended)
or cop y-pasting (not recommended) the commands that ha v e a ' , please correct it to the
single-quote (or ap ostrophe) c haracter, otherwise the command will fail.
1.13 Ac kno wledgmen ts
Gn uastro w ould not ha v e b een p ossible without sc holarships and gran ts from sev eral funding
institutions. W e th us ask that if y ou used Gn uastro in an y of y our pap ers/rep orts, please
add the prop er citation and ac kno wledge the funding agencies/pro jects. F or details of whic h
pap ers to cite (ma y b e differen t for differen t programs) and get the ac kno wledgmen t state-
men t to include in y our pap er, please run the relev an t programs with the common --cite
option lik e the example commands b elo w (for more on --cite , please see Section 4.1.2.3
[Op erating mo de options], page 256).
$ astnoisechisel --cite
$ astmkcatalog --cite
Here, w e will ac kno wledge all the institutions (and their gran ts) along with the p eople
who help ed mak e Gn uastro p ossible. The full list of Gn uastro authors is a v ailable at the
start of this b o ok and the AUTHORS file in the source co de (b oth are generated automatically
from the v ersion controlled history). The plain text file THANKS , whic h is also distributed
along with the source co de, con tains the list of p eople and institutions who pla y ed an indirect
role in Gn uastro (not committed an y co de in the Gn uastro v ersion con trolled history).
The Japanese Ministry of Education, Culture, Sp orts, Science, and T ec hnology (MEXT)
sc holarship for Mohammad Akhlaghi’s Masters and PhD degree in T ohoku Univ ersit y As-
tronomical Institute had an instrumen tal role in the long term learning and planning that
made the idea of Gn uastro p ossible. The very critical view points of Professor T ak ashi
Ic hik a w a (Mohammad’s adviser) w ere also instrumen tal in the initial ideas and creation
of Gn uastro. Afterw ards, the Europ ean Researc h Council (ER C) adv anced gran t 339659-
MUSICOS (Principal in v estigator: Roland Bacon) w as vital in the gro wth and expansion
of Gn uastro. W orking with Roland at the Cen tre de Rec herc he Astroph ysique de Ly on
(CRAL), enabled a thorough re-write of the core functionalit y of all libraries and programs,
turning Gn uastro in to the large collection of generic programs and libraries it is to da y . A t
the Instituto de Astrofisica de Canarias (IA C, and in particular in collab oration with Johan
Knap en and Ignacio T rujillo), Gn uastro matured and its user base significan tly grew. W ork
on impro ving Gn uastro is no w con tin uing primarily in the Cen tro de Estudios de F ´
isica del
Cosmos de Arag´ on (CEF CA), lo cated in T eruel, Spain.
In general, w e w ould lik e to gratefully thank the follo wing p eople for their useful and con-
structiv e commen ts and suggestions (in alphab etical order b y family name): V alen tina Abril-
melgarejo, Marjan Akbari, Carlos Allende Prieto, Hamed Altafi, Roland Bacon, Rob erto
Baena Gall ´ e, Zahra Bagheri, Karl Berry , F aezeh Bidjarc hian, Leindert Boogaard, Nicolas
Bouc h ´ e, Stefan Br ¨ uns, F ernando Buitrago Alonso, Adrian Bunk, Rosa Calvi, Mark Cal-
abretta, Juan Castillo Ram ´ ırez, Nushkia Cham ba, Sergio Ch ueca Urzay , T amara Civ era
Lorenzo, Benjamin Clement, Nima Dehdilani, Andr ´ es Del Pino Molina, An tonio Diaz
Chapter 1: In tro duction 20
Diaz, P aola Dimauro, Alexey Dokuc haev, Pierre-Alain Duc, Alessandro Edero clite, El-
ham Eftekhari, P aul Eggert, Sepideh Esk andarlou, S ´
ilvia F arras, Juan An tonio F ern´ andez
On tiv eros, Gaspar Galaz, Andr ´ es Garc ´ ıa-Serra Romero, Zohre Ghaffari, Th ´ er ` ese Go defro y ,
Giulia Golini, Craig Gordon, Martin Guerrero Roncel, Madusha Guna w ardhana, Bruno
Haible, Stephen Hamer, Siyang He, Zahra Hosseini, Leslie Hun t, T ak ashi Ic hik aw a, Ra ´ ul
Infan te Sainz, Brandon In v ergo, Oryna Iv ash tenk o, Aur ´ elien Jarno, Ooldo oz Kab o o d, Lee
Kelvin, Brandon Kelly , Mohammad-Reza Khellat, Johan Knap en, Geoffry Krouc hi, Martin
Kuemmel, T eet Kuutma, Clotilde Laigle, Floriane Leclercq, Alan Lefor, Ja vier Licandro,
Jerem y Lim, Giacomo Lorenzetti, Alejandro Lumbreras Calle, Sebasti´ an Luna V alero, Al-
b erto Madrigal, Guillaume Mahler, Juan Miro, Alireza Molaeinezhad, Ja vier Moldon, Juan
Molina T obar, F rancesco Montanari, Raphael Morales, Carlos Morales So corro, Sylv ain
Mottet, Dmitrii Oparin, F ran¸ cois Oc hsen b ein, Bertrand P ain, Rahna P a yy asseri Thandu-
parac k al, William P ence, Irene Pin tos Castro, Mam ta Pommier, Marcel P op escu, Bob
Proulx, Joseph Putk o, Samane Ra ji, Ignacio Ruiz Cejudo, T eymo or Saifollahi, Joanna
Sak o wsk a, Elham Saremi, Nafise Sedighi, Markus Sc haney , Y ah y a Sefidbakh t, Alejandro
Serrano Borlaff, Zahra Sharbaf, Da vid Sh up e, Leigh Smith, Jenn y Sorce, Manuel S´ anchez-
Bena v en te, Lee Spitler, Ric hard Stallman, Mic hael Stein, Ole Streic her, Alfred M. Szmidt,
Mic hel T allon, Juan C. T ello, Vincenzo T esta, ´
Eric Thi ´ ebaut, Ignacio T rujillo, P eter T eub en,
Mathias Urbano, Da vid V alls-Gabaud, Jes ´ us V arela, Jes ´ us V ega, Aaron W atkins, Ric hard
Wilbur, Phil Wy ett, Dennis Williamson, Mic hael H.F. Wilkinson, Christopher Willmer,
Greg W o oledge, Xiuqin W u, Sara Y ousefi T aemeh, Johannes Zabl. The GNU F renc h T rans-
lation T eam is also managing the F renc h v ersion of the top Gn uastro w eb page whic h w e
highly appreciate. Finally , w e should thank all the (sometimes anon ymous) p eople in v arious
online forums who patien tly answ ered all our small (but imp ortan t) tec hnical questions.
All w ork on Gn uastro has b een v olun tary , but the authors are most grateful to the
follo wing institutions (in c hronological order) for hosting/supp orting us in our researc h.
Where necessary , these institutions ha v e disclaimed an y o wnership of the parts of Gnuastro
that w ere dev elop ed there, th us insuring the freedom of Gn uastro for the future (see Sec-
tion 13.12.1 [Cop yrigh t assignmen t], page 956). W e highly appreciate their supp ort for free
soft w are, and th us free science, and therefore a free so ciet y .
T ohoku Univ ersit y Astronomical Institute, Sendai, Japan.
Univ ersit y of Salen to, Lecce, Italy .
Cen tre de Rec herc he Astroph ysique de Ly on (CRAL), Ly on, F rance.
Instituto de Astrofisica de Canarias (IA C), T enerife, Spain.
Cen tro de Estudios de F ´
isica del Cosmos de Arag´ on (CEF CA), T eruel, Spain.
Go ogle Summer of Co de 2020, 2021 and 2022
21
2 T utorials
T o help new users ha ve a smooth and easy start with Gn uastro, in this c hapter sev eral
thoroughly elab orated tutorials, or co okb o oks, are pro vided. These tutorials demonstrate
the capabilities of differen t Gn uastro programs and libraries, along with tips and guidelines
for the b est practices of using them in v arious realistic situations.
W e strongly recommend going through these tutorials to get a go o d feeling of ho w
the programs are related (built in a mo dular design to b e used together in a pip eline),
v ery similar to the core Unix-based programs that they w ere mo deled on. Therefore these
tutorials will help in optimally using Gn uastro’s programs (and generally , the Unix-lik e
command-line en vironmen t) effectiv ely for y our researc h.
The first three tutorials (Section 2.1 [General program usage tutorial], page 22, and Sec-
tion 2.2 [Detecting large extended targets], page 80, and Section 2.3 [Building the extended
PSF], page 102) use real input datasets from some of the deep Hubble Space T elescop e
(HST) images, the Sloan Digital Sky Surv ey (SDSS) and the Ja v alam bre Photometric Lo-
cal Univ erse Surv ey (J-PLUS) resp ectiv ely . Their aim is to demonstrate some real-w orld
problems that man y astronomers often face and how they can b e solv ed with Gn uastro’s
programs. The fourth tutorial (Section 2.4 [Sufi sim ulates a detection], page 123) fo cuses
on sim ulating astronomical images, whic h is another critical asp ect of an y analysis!
The ultimate aim of Section 2.1 [General program usage tutorial], page 22, is to detect
galaxies in a deep HST image, measure their p ositions, magnitude and select those with
the strongest colors. In the pro cess, it tak es man y detours to introduce you to the useful
capabilities of man y of the programs. So please b e patien t in reading it. If y ou do not ha v e
m uc h time and can only try one of the tutorials, w e recommend this one.
Section 2.2 [Detecting large extended targets], page 80, deals with a ma jor problem in
astronom y: effectiv ely detecting the faint outer wings of brigh t (and large) nearb y galaxies
to extremely lo w surface brigh tness levels (roughly one quart er of the lo cal noise lev el in the
example discussed). Besides the in teresting scien tific questions in these lo w-surface brigh t-
ness features, failure to prop erly detect them will bias the measuremen ts of the bac kground
ob jects and the surv ey’s noise estimates. This is an imp ortan t issue, esp ecially in wide
surv eys. Because brigh t/large galaxies and stars 1 , co v er a significan t fraction of the surv ey
area.
Section 2.3 [Building the extended PSF], page 102, tac kles an imp ortan t problem in
astronom y: ho w the extract the PSF of an image, to the largest p ossible exten t, without
assuming an y functional form. In Gn uastro w e ha v e m ultiple installed scripts for this job.
Their usage and logic b ehind b est tuning them for the particular step, is fully describ ed
in this tutorial, on a real dataset. The tutorial concludes with subtracting that extended
PSF from the science image; thus giving y ou a cleaner image (with no scattered ligh t of the
brigh ter stars) for y our higher-lev el analysis.
Section 2.4 [Sufi sim ulates a detection], page 123, has a fictional 2 setting! Sho wing ho w
Ab d al-rahman Sufi (903 – 986 A.D., the first recorded description of “nebulous” ob jects
1 Stars also ha v e similarly large and extended wings due to the p oin t spread function, see Section 8.1.1.2
[P oin t spread function], page 633.
2 The t w o historically motiv ated tutorials (Section 2.4 [Sufi sim ulates a detection], page 123, is not in tended
to b e a historical reference (the historical facts of this fictional tutorial used Wikip edia as a reference).)
This form of presen ting a tutorial w as influenced by the PGF/TikZ and Beamer man uals. They are b oth
Chapter 2: T utorials 22
in the hea v ens is attributed to him) could ha v e used some of Gn uastro’s programs for a
realistic sim ulation of his observ ations and see if his detection of nebulous ob jects w as trust-
able. Because all conditions are under con trol in a sim ulated/mo c k en vironmen t/dataset,
mo c k datasets can b e a v aluable to ol to insp ect the limitations of y our data analysis and
pro cessing. But they need to b e as realistic as p ossible, so this tutorial is dedicated to this
imp ortan t step of an analysis (sim ulations).
There are other tutorials also, on things that are commonly necessary in astronomi-
cal researc h: In Section 2.5 [Detecting lines and extracting sp ectra in 3D data], page 134,
w e use MUSE cub es (an IFU dataset) to sho w ho w y ou can subtract the con tin uum, de-
tect emission-line features, extract sp ectra and build pseudo narro w-band images. In Sec-
tion 2.6.1 [Color c hannels in same pixel grid], page 151, w e demonstrate ho w y ou can w arp
m ultiple images in to a single pixel grid (often necessary with m ulti-w a v elength data), and
build a single color image. In Section 2.9 [Moir ´ e pattern in stac king and its correction],
page 190, w e sho w how y ou can a v oid the un w an ted Moir ´ e pattern whic h happ ens when
w arping separate exp osures to build a stac k ed/co-add deep er image. In Section 2.7 [Zero
p oin t of an image], page 165, w e review the pro cess of estimating the zero p oin t of an
image using a reference image or catalog. Finally , in Section 2.8 [P ointing pattern design],
page 175, w e sho w the pro cess b y whic h y ou can sim ulate a dither pattern to find the b est
observing strategy for y our next exciting scien tific pro ject.
In these tutorials, w e ha v e in ten tionally a v oided to o man y cross references to mak e it
more easy to read. F or more information ab out a particular program, y ou can visit the
section with the same name as the program in this b o ok. Eac h program section in the
subsequen t c hapters starts b y explaining the general concepts b ehind what it do es, for ex-
ample, see Section 6.3 [Con v olv e], page 470. If y ou only w an t practical information on
running a program, for example, its options/configuration, input(s) and output(s), please
consult the subsection titled “In v oking ProgramName”, for example, see Section 7.2.2 [In-
v oking NoiseChisel], page 547. F or an explanation of the con ven tions w e use in the example
co des through the b o ok, please see Section 1.12 [Con v en tions], page 18.
2.1 General program usage tutorial
Measuring colors of astronomical ob jects in broad-band or narro w -band images is one of
the most basic and common steps in astronomical analysis. Here, w e will use Gn uastro’s
programs to get a ph ysical scale (area at certain redshifts) of the field w e are studying,
detect ob jects in a Hubble Space T elescop e (HST) image, measure their colors and iden tify
the ones with the strongest colors, do a visual insp ection of these ob jects and insp ect spatial
p osition in the image. After this tutorial, y ou can also try the Section 2.2 [Detecting large
extended targets], page 80, tutorial whic h go es in to a little more detail on detecting v ery
lo w surface brigh tness signal.
During the tutorial, w e will tak e man y detours to explain, and practically demonstrate,
the man y capabilities of Gn uastro’s programs. In the end y ou will see that the things y ou
learned during this tutorial are m uc h more generic than this particular problem and can b e
pac k ages in T
E X and L
A
T
E X, the first is a high-lev el vector graphic programming en vironmen t, while
with the second y ou can mak e presen tation slides. On a similar topic, there are also some nice words
of wisdom for Unix-lik e systems called Ro otless Ro ot ( http://catb.org/esr/writings/unix-koans ).
These also ha v e a similar style but they use a m ythical figure named Master F o o. If y ou already hav e
some exp erience in Unix-lik e systems, you will definitely find these Unix Koans en tertaining/educativ e.
Chapter 2: T utorials 23
used in solving a wide v ariet y of problems in v olving the analysis of data (images or tables).
So please do not rush, and go through the steps patien tly to optimally master Gn uastro.
In this tutorial, w e will use the HSTeXtreme Deep Field ( https://archive.stsci.edu/
prepds/xdf ) dataset. Lik e almost all astronomical surv eys, this dataset is free for do wnload
and usable b y the public. Y ou will need the follo wing to ols in this tutorial: Gn uastro, SA O
DS9 3 , GNU Wget 4 , and A WK (most common implementation is GNU A WK 5 ).
This tutorial w as first prepared for the “Exploring the Ultra-Lo w Surface Brigh tness
Univ erse” w orkshop (No v em b er 2017) at the ISSI in Bern, Switzerland. It w as further ex-
tended in the “4th Indo-F renc h Astronom y Sc ho ol” (July 2018) organized b y LIO, CRAL
CNRS UMR5574, UCBL, and IUCAA in Ly on, F rance. W e are v ery grateful to the orga-
nizers of these w orkshops and the attendees for the v ery fruitful discussions and suggestions
that made this tutorial p ossible.
W rite the example commands man ually: T ry to t yp e the example commands on y our
terminal man ually and use the history feature of y our command-line (by pressing the
“up” button to retriev e previous commands). Do not simply cop y and paste the commands
sho wn here. This will help simulate future situations when y ou are pro cessing y our o wn
datasets.
2.1.1 Calling Gn uastro’s programs
A handy feature of Gn uastro is that all program names start with ast . This will allo w y our
command-line pro cessor to easily list and auto-complete Gn uastro’s programs for y ou. T ry
t yping the follo wing command (press TAB k ey when y ou see <TAB> ) to see the list:
$ ast<TAB><TAB>
An y program that starts with ast (including all Gn uastro programs) will b e shown. By
c ho osing the subsequen t c haracters of y our desired program and pressing <TAB><TAB> again,
the list will narro w do wn and the program name will auto-complete once y our input c harac-
ters are unam biguous. In short, y ou often do not need to t yp e the full name of the program
y ou w an t to run.
2.1.2 Accessing do cumen tation
Gn uastro con tains a large n um b er of programs and it is natural to forget the details of eac h
program’s options or inputs and outputs. Therefore, b efore starting the analysis steps of
this tutorial, let’s review ho w y ou can access this b o ok to refresh y our memory an y time
y ou w an t, without ha ving to tak e y our hands off the k eyb oard.
When y ou install Gn uastro, this b o ok is also installed on y our system along with all
the programs and libraries, so y ou do not need an in ternet connection to access/read it.
Also, b y accessing this b o ok as describ ed b elow, y ou can b e sure that it corresp onds to your
installed v ersion of Gn uastro.
3 See Section A.1 [SA O DS9], page 964, a v ailable at http://ds9.si.edu/site/Home.html
4 https://www.gnu.org/software/wget
5 https://www.gnu.org/software/gawk
Chapter 2: T utorials 24
GNU Info 6 is the program in c harge of displa ying the man ual on the command-line
(for more, see Section 4.3.4 [Info], page 272). T o see this whole b o ok on y our command-
line, please run the follo wing command and press subsequen t k eys. Info has its o wn mini-
en vironmen t, therefore w e will sho w the k eys that must be pressed in the mini-environmen t
after a -> sign. Y ou can also ignore anything after th e # sign in the middle of the line, they
are only for y our information.
$ info gnuastro # Open the top of the manual.
-> <SPACE> # All the book chapters.
-> <SPACE> # Continue down: show sections.
-> <SPACE> ... # Keep pressing space to go down.
-> q # Quit Info, return to the command-line.
The thing that greatly simplifies na vigation in Info is the links (regions with an under-
line). Y ou can immediately go to the next link in the page with the <TAB> k ey and press
<ENTER> on it to go in to that part of the man ual. T ry the commands ab o v e again, but this
time also use <TAB> to go to the links and press <ENTER> on them to go to the resp ectiv e
section of the b o ok. Then follo w a few more links and go deep er in to the b o ok. T o return to
the previous page, press l (small L). If y ou are searching for a sp ecific phrase in the whole
b o ok (for example, an option name), press s and t yp e y our searc h phrase and end it with
an <ENTER> . Finally , you can return to the command line and quit Info b y pressing the q
k ey .
Y ou do not need to start from the top of the man ual ev ery time. F or example, to get to
Section 7.2.2 [In v oking NoiseChisel], page 547, run the follo wing command. In general, all
programs ha v e suc h an “In voking ProgramName” section in this b o ok. These sections are
sp ecifically for the description of inputs, outputs and configuration options of each program.
Y ou can access them directly for eac h program b y giving its executable name to Info.
$ info astnoisechisel
The other sections do not ha v e suc h shortcuts. T o directly access them from the
command-line, you need to tell Info to look in to Gn uastro’s man ual, then lo ok for the
sp ecific section (an unam biguous title is necessary). F or example, if y ou only w an t to re-
view/remem b er NoiseChisel’s Section 7.2.2.2 [Detection options], page 552), just run the
follo wing command. Note ho w case is irrelev an t for Info when calling a title in this manner.
$ info gnuastro "Detection options"
In general, Info is a p o w erful and conv enien t w a y to access this whole b o ok with detailed
information ab out the programs y ou are running. If y ou are not already familiar with it,
please run the follo wing command and just read along and do what it sa ys to learn it. Do
not stop un til y ou feel sufficien tly fluen t in it. Please in v est the half an hour’s time necessary
to start using Info comfortably . It will greatly impro v e y our pro ductivit y and y ou will start
reaping the rew ards of this in v estmen t v ery so on.
$ info info
As a go o d scien tist y ou need to feel comfortable to pla y with the features/options and
a v oid (b e critical to) using default v alues as m uc h as p ossible. On the other hand, our
h uman memory is limited, so it is imp ortan t to b e able to easily access an y part of this
b o ok fast and remem b er the option names, what they do and their acceptable v alues.
6 GNU Info is already a v ailable on almost all Unix-like operating systems.
Chapter 2: T utorials 25
If y ou just w ant the option names and a short d escription, calling the program with the
--help option migh t also b e a go o d solution lik e the first example b elo w. If y ou kno w a
few c haracters of the option name, y ou can feed the prin ted output to grep lik e the second
or third example commands.
$ astnoisechisel --help
$ astnoisechisel --help | grep quant
$ astnoisechisel --help | grep check
2.1.3 Setup and data do wnload
The first step in the analysis of the tutorial is to do wnload the necessary input datasets.
First, to k eep things clean, let’s create a gnuastro-tutorial directory and con tin ue all
future steps in it:
$ mkdir gnuastro-tutorial
$ cd gnuastro-tutorial
W e will b e using the near infra-red Wide Field Camera ( http://www.stsci.edu/hst/
wfc3 ) dataset. If y ou already ha v e them in another directory (for example, XDFDIR , with
the same FITS file names), y ou can set the download directory to b e a sym b olic link to
XDFDIR with a command lik e this:
$ ln -s XDFDIR download
Otherwise, when the following images are not already presen t on y our system, y ou can mak e
a download directory and do wnload them there.
$ mkdir download
$ cd download
$ xdfurl=http://archive.stsci.edu/pub/hlsp/xdf
$ wget $xdfurl/hlsp_xdf_hst_wfc3ir-60mas_hudf_f105w_v1_sci.fits
$ wget $xdfurl/hlsp_xdf_hst_wfc3ir-60mas_hudf_f125w_v1_sci.fits
$ wget $xdfurl/hlsp_xdf_hst_wfc3ir-60mas_hudf_f160w_v1_sci.fits
$ cd ..
In this tutorial, w e will just use these three filters. Later, y ou ma y need to do wnload more
filters. T o do that, y ou can use the shell’s for lo op to do wnload them all in series (one after
the other 7 ) with one command lik e the one b elo w for the WF C3 filters. Put this command
instead of the three wget commands ab o v e. Recall that all the extra spaces, bac kslashes
( \ ), and new lines can b e ignored if y ou are t yping on the lines on the terminal.
$ for f in f105w f125w f140w f160w; do \
wget $xdfurl/hlsp_xdf_hst_wfc3ir-60mas_hudf_"$f"_v1_sci.fits; \
done
2.1.4 Dataset insp ection and cropping
First, let’s visually insp ect the datasets w e do wnloaded in Section 2.1.3 [Setup and data
do wnload], page 25. Let’s tak e F160W image as an example. One of the most common
programs for viewing FITS images is SA O DS9, whic h is usually called through the ds9
7 Note that y ou only ha ve one port to the internet, so do wnloading in parallel will actually b e slo wer than
do wnloading in series.
Chapter 2: T utorials 32
done
F ortunately , the shell has a useful to ol/program to prin t a sequence of n um b ers that is nicely
called seq (short for “sequence”). Y ou can use it instead of typing all the differen t redshifts
in the lo op ab o v e. F or example, the lo op b elo w will calculate and prin t the tangen tial
co v erage of this field across a larger range of redshifts (0.1 to 5) and with finer incremen ts
of 0.1. F or more on the LC_NUMERIC command, see Section 4.11 [Numeric lo cale], page 292.
## If your system language uses ',' (not '.') as decimal separator.
$ export LC_NUMERIC=C
## The loop over the redshifts
$ for z in $(seq 0.1 0.1 5); do \
k=$(astcosmiccal -z$z --arcsectandist); \
echo $z $k $a | awk '{print $1, ($2*60)^2 * $3 / 1e6}'; \
done
Ha v e a lo ok at the t w o prin ted columns. The first is the redshift, and the second is
the area of this image at that redshift (in mega-parsecs squared). Redshift ( https://en.
wikipedia.org/wiki/Redshift ) ( z ) is often used as a pro xy for distance in galaxy ev olu-
tion and cosmology: a higher redshift corresp onds to larger line-of-sigh t como ving distance.
No w, ha v e a lo ok at the first few v alues. At z = 0 . 1 and z = 0 . 5, this image co v ers
0 . 05 M pc 2 and 0 . 57 M pc 2 resp ectiv ely . This increase of co v erage with redshift is exp ected
b ecause a fixed angle will co ver a larger tangen tial area at larger distances. Ho w ev er,
as y ou come do wn the list (to higher redshifts) y ou will notice that this relation do es
not hold! The largest co v erage is at z = 1 . 6: at higher redshifts, the area decreases,
and con tin ues decreasing!!! In flat FLR W cosmology (including ΛCDM), the only fac-
tor con tributing to this is the (1 + z )$ factor from the expansion of the univ erse, see
the Wikip edia page ( https://en.wikipedia.org/wiki/Angular_diameter_distance#
Angular_diameter_turnover_point ), with no curv ature effect.
In case y ou ha v e TOPCA T, y ou can visualize this as a plot (if y ou do not ha v e TOPCA T,
see Section A.2 [TOPCA T], page 965). T o do so, first y ou need to sa v e the output of the
lo op ab o v e in to a FITS table b y piping the output to Gn uastro’s T able program and giving
an output name:
$ for z in $(seq 0.1 0.1 5); do \
k=$(astcosmiccal -z$z --arcsectandist); \
echo $z $k $a | awk '{print $1, ($2*60)^2 * $3 / 1e6}'; \
done | asttable --output=z-vs-tandist.fits
Y ou can no w use Gn uastro’s astscript-fits-view to op en this table in TOPCA T with
the command b elo w. Do y ou remem b er this script from Section 2.1.4 [Dataset insp ection
and cropping], page 25? There, w e used it to view a FITS image with DS9! This script will
see if the first dataset in the image is a table or an image and will call TOPCA T or DS9
accordingly: making it a v ery con v enien t to ol to insp ect the con ten ts of all t yp es of FITS
data.
$ astscript-fits-view z-vs-tandist.fits
After TOPCA T op ens, y ou will see the name of the table z-vs-tandist.fits in the
left panel. On the top men u bar, select the “Graphics” men u, then select “Plain plot” to
Chapter 2: T utorials 33
visualize the t w o columns prin ted ab o v e as a plot and get a b etter impression of the turn
o v er p oin t of the image cosmological co v erage.
2.1.7 Building custom programs with the library
In Section 2.1.6 [Cosmological co v erage and visualizing tables], page 30, we repeated a
certain calculation/output of a program m ultiple times using the shell’s for lo op. This
simple w a y of rep eating a calculation is great when it is only necessary once. Ho w ev er, if
y ou commonly need this calculation and p ossibly for a larger n um b er of redshifts at higher
precision, the command ab o v e can b e slo w. Please try it out b y c hanging the sequence
command in the previous section to ‘ seq 0.1 0.01 10 ’. It will tak e ab out 11 seconds 13 !
This can b e impro v ed b y hundr e ds of times! This section will sho w y ou ho w.
Generally , rep eated calls to a generic program (lik e CosmicCalculator) are slo w, b ecause
a generic program can ha v e a lot of o verhead on eac h call. T o b e generic and easy to
op erate, CosmicCalculator has to parse the command-line and all configuration files (see
Section 2.1.8 [Option managemen t and configuration files], page 35) whic h con tain h uman-
readable c haracters and need a lot of pre-pro cessing to b e ready for pro cessing b y the
computer. Afterw ards, CosmicCalculator has to c hec k the sanit y of its inputs and c hec k
whic h of its man y options y ou ha v e ask ed for. All the this pre-pro cessing tak es as m uc h
time as the high-lev el calculation y ou are requesting, and it has to re-do all of these for
ev ery redshift in y our lo op.
T o greatly sp eed up the pro cessing, y ou can directly access the core work-horse of Cos-
micCalculator without all that o v erhead b y designing y our custom program for this job.
Using Gn uastro’s library , you can write y our o wn tin y program particularly designed for
this exact calculation (and nothing else!). T o do that, cop y and paste the follo wing C
program in a file called myprogram.c .
#include <math.h>
#include <stdio.h>
#include <stdlib.h>
#include <gnuastro/cosmology.h>
int
main(void)
{
double area=4.03817; /* Area of field (arcmin^2). */
double z, adist, tandist; /* Temporary variables. */
/* Constants from Plank 2018 (arXiv:1807.06209, Table 2) */
double H0=67.66, olambda=0.6889, omatter=0.3111, oradiation=0;
/* Do the same thing for all redshifts (z) between 0.1 and 5. */
for(z=0.1; z<10; z+=0.01)
{
13 T o measure how m uch time the loop of Section 2.1.6 [Cosmological cov erage and visualizing tables],
page 30, tak es on your system, y ou can use the time command. First put the whole lo op (and pipe) into
a plain-text file (to b e loaded as a shell script) called z-vs-tandist.sh . Then run this command: time
-p bash z-vs-tandist.sh . The relev ant time (in seconds) is sho wn after real .
Chapter 2: T utorials 34
/* Calculate the angular diameter distance. */
adist=gal_cosmology_angular_distance(z, H0, olambda,
omatter, oradiation);
/* Calculate the tangential distance of one arcsecond. */
tandist = adist * 1000 * M_PI / 3600 / 180;
/* Print the redshift and area. */
printf("%-5.2f %g\n", z, pow(tandist * 60,2) * area / 1e6);
}
/* Tell the system that everything finished successfully. */
return EXIT_SUCCESS;
}
Then run the follo wing command to compile y our program and run it.
$ astbuildprog myprogram.c
In the command ab o v e, y ou used Gn uastro’s BuildProgram program. Its job is to simplify
the compilation, linking and running of simple C programs that use Gn uastro’s library (lik e
this one). BuildProgram is designed to manage Gn uastro’s dep endencies, compile and link
y our custom program and then run it.
Did y ou notice ho w your custom program created the table almost instan taneously?
T ec hnically , it only to ok ab out 0.03 seconds! Recall that the for lo op of Section 2.1.6
[Cosmological co v erage and visualizing tables], page 30, to ok more than 11 seconds (or
∼ 367 times slo w er!).
Please run the ls command to see a listing of the files in the curren t directory . Y ou will
notice that a new file called myprogram has b een created. This is the compiled program
that w as created and run b y the command ab o v e (its in binary mac hine co de format, not
h uman-readable an y more). Y ou can run it again to get the same results b y executing it:
$ ./myprogram
The efficiency of y our custom myprogram compared to rep eated calls to CosmicCalculator
is b ecause in the latter, the requested pro cessing is comparable to the necessary o v erheads.
F or other programs that tak e large input datasets and do complicated pro cessing on them,
the o v erhead is usually negligible compared to the pro cessing. In suc h cases, the libraries
are only useful if y ou w an t a differen t/new pro cessing compared to the functionalities in
Gn uastro’s existing programs.
Gn uastro has a large library whic h is used extensively b y all the programs. In other
w ords, the library is like th e sk eleton of Gn uastro. F or the full list of a v ailable functions
classified b y con text, please see Section 12.3 [Gn uastro library], page 741. Gn uastro’s library
and BuildProgram are created to mak e it easy for you to use these p ow erful features as y ou
lik e. This giv es y ou a high lev el of creativit y , while also pro viding efficiency and robust-
ness. Sev eral other complete w orking examples (in v olving images and tables) of Gnuastro’s
libraries can b e see in Section 12.4 [Library demo programs], page 915.
But for this tutorial, let’s stop discussing the libraries here and get bac k to Gn uastro’s
already built programs (whic h do not need C programming). But b efore con tin uing, let’s
clean up the files w e do not need an y more:
Chapter 2: T utorials 35
$ rm myprogram* z-vs-tandist*
2.1.8 Option managemen t and configuration files
In the previous section (Section 2.1.6 [Cosmological co v erage and visualizing tables],
page 30), when y ou ran CosmicCalculator, y ou only sp ecified the redshfit with -z2 option.
Y ou did not sp ecify the cosmological parameters that are necessary for the calculations!
P arameters lik e the Hubble constan t ( H 0 ) and the matter densit y . In spite of this,
CosmicCalculator done its pro cessing and prin ted results.
None of Gn uastro’s programs keep a default v alue in ternally within their co de (they are
all set b y the user)! So where did the necessary cosmological parameters that are necessary
for its calculations come from? What w ere the v alues to those parameters? In short, they
come from a configuration file (see Section 4.2.2 [Configuration file precedence], page 268),
and the final used v alues can b e c hec k ed/edited on the command-line. In this section w e
will review this imp ortan t asp ect of all the programs in Gn uastro.
Configuration files are an imp ortan t part of all Gn uastro’s programs, esp ecially the ones
with a large n um b er of options, so it is imp ortan t to understand this part w ell. Once y ou
get comfortable with configuration files, you can mak e go o d use of them in all Gn uastro
programs (for example, NoiseChisel). F or example, to do optimal detection on v arious
datasets, y ou can ha v e configuration files for differen t noise prop erties. The configuration
of eac h program (b esides its v ersion) is vital for the repro ducibilit y of y our results, so it is
imp ortan t to manage them prop erly .
As w e sa w ab o v e, the full list of the options in all Gn uastro programs can b e seen with the
--help option. T ry calling it with CosmicCalculator as shown below. Note ho w options are
group ed b y con text to mak e it easier to find your desired option. Ho w ev er, in eac h group,
options are ordered alphab etically .
$ astcosmiccal --help
After running the command ab o ve, please scroll to the line that y ou ran this command
and read through the output (its the same format for all the programs). All options ha v e
a long format (starting with -- and a m ulti-c haracter name) and some ha v e a short format
(starting with - and a single c haracter), for more see Section 4.1.1.2 [Options], page 248.
The options that exp ect a v alue, ha v e an = sign after their long v ersion. The format of
their exp ected v alue is also sho wn as FLT , INT or STR for floating p oin t n um b ers, in teger
n um b ers, and strings (filenames for example) resp ectiv ely .
Y ou can see the v alues of all options that need one with the --printparams option (or
its short format: -P ). --printparams is common to all programs (see Section 4.1.2 [Com-
mon options], page 250). Y ou can see the default cosmological parameters, from the Plank
collab oration 2020 ( https://arxiv.org/abs/1807.06209 ), under the # Input: title:
$ astcosmiccal -P
# Input:
H0 67.66 # Current expansion rate (Hubble constant).
olambda 0.6889 # Current cosmological cst. dens. per crit. dens.
omatter 0.3111 # Current matter density per critical density.
oradiation 0 # Current radiation density per critical density.
Chapter 2: T utorials 36
Let’s sa y y ou w an t to do the calculation in the previous section using H 0 = 70 km/s/Mp c.
T o do this, just add --H0=70 after the command ab o v e (while k eeping the -P ). In the output,
y ou can see that the used Hubble constan t has also c hanged.
$ astcosmiccal -P --H0=70
Afterw ards, delete the -P and add a -z2 to see the calculations with the new cosmology (or
configuration).
$ astcosmiccal --H0=70 -z2
F rom the output of the --help option, note ho w the option for Hubble constan t has
b oth short ( -H ) and long ( --H0 ) formats. One final note is that the equal ( = ) sign is
not mandatory . In the short format, the v alue can stic k to the actual option (the short
option name is just one c haracter after-all, th us easily iden tifiable) and in the long format,
a white-space c haracter is also enough.
$ astcosmiccal -H70 -z2
$ astcosmiccal --H0 70 -z2 --arcsectandist
When an option do es not need a v alue, and has a short format (lik e --arcsectandist ),
y ou can easily app end it b efor e other short options. So the last command ab o ve can also
b e written as:
$ astcosmiccal --H0 70 -sz2
Let’s assume that in one pro ject, y ou w an t to only use rounded cosmological parameters
( H 0 of 70km/s/Mp c and matter densit y of 0.3). Y ou should therefore run CosmicCalculator
lik e this:
$ astcosmiccal --H0=70 --olambda=0.7 --omatter=0.3 -z2
But ha ving to t yp e these extra options ev ery time y ou run CosmicCalculator will b e
prone to errors (t yp os in particular), frustrating and slo w. Therefore in Gn uastro, y ou can
put all the options and their v alues in a “Configuration file” and tell the programs to read
the option v alues from there.
Let’s create a configuration file... With y our fa v orite text editor, mak e a file named
my-cosmology.conf (or my-cosmology.txt , the suffix do es not matter for Gn uastro, but a
more descriptiv e suffix lik e .conf is recommended for h umans reading y our co de and seeing
y our files: this includes y ou, lo oking in to y our o wn pro ject, in a couple of mon ths that y ou
ha v e forgot the details!). Then put the follo wing lines inside of the plain-text file. One
space b et w een the option v alue and name is enough, the v alues are just under eac h other to
help in readabilit y . Also note that y ou should only use long option names in configuration
files.
H0 70
olambda 0.7
omatter 0.3
Y ou can no w tell CosmicCalculator to read this file for option v alues immediately using the
--config option as sho wn b elo w. Do y ou see how the output of the follo wing command
corresp onds to the option v alues in my-cosmology.conf , and is therefore iden tical to the
previous command?
$ astcosmiccal --config=my-cosmology.conf -z2
Chapter 2: T utorials 37
But still, ha ving to t yp e --config=my-cosmology.conf ev ery time is anno ying, is not
it? If y ou need this cosmology ev ery time y ou are w orking in a sp ecific directory , y ou can
use Gn uastro’s default configuration file names and a v oid ha ving to t yp e it man ually .
The default configuration files (that are c hec k ed if they exist) m ust b e placed in the
hidden .gnuastro sub-directory (in the same directory y ou are running the program).
Their file name (within .gnuastro ) m ust also b e the same as the program’s executable
name. So in the case of CosmicCalculator, the default configuration file in a giv en directory
is .gnuastro/astcosmiccal.conf .
Let’s do this. W e will first make a directory for our custom cosmology , then build a
.gnuastro within it. Finally , w e will cop y the custom configuration file there:
$ mkdir my-cosmology
$ mkdir my-cosmology/.gnuastro
$ mv my-cosmology.conf my-cosmology/.gnuastro/astcosmiccal.conf
Once y ou run CosmicCalculator within my-cosmology (as sho wn b elo w), y ou will see
ho w y our custom cosmology has b een implemen ted without having to t yp e an ything extra
on the command-line.
$ cd my-cosmology
$ astcosmiccal -P # Your custom cosmology is printed.
$ cd ..
$ astcosmiccal -P # The default cosmology is printed.
T o further simplify the pro cess, you can use the --setdirconf option. If y ou are already
in y our desired w orking directory , calling this option with the others will automatically write
the final v alues (along with descriptions) in .gnuastro/astcosmiccal.conf . F or example,
try the commands b elo w:
$ mkdir my-cosmology2
$ cd my-cosmology2
$ astcosmiccal -P
$ astcosmiccal --H0 70 --olambda=0.7 --omatter=0.3 --setdirconf
$ astcosmiccal -P
$ cd ..
Gn uastro’s programs also ha v e default configuration files for a sp ecific user (when run in
an y directory). This allows y ou to set a sp ecial b eha vior every time a p rogram is run b y a
sp ecific user. Only the directory and filename differ from the ab o v e, the rest of the pro cess
is similar to b efore. Finally , there are also system-wide configuration files that can b e used
to define the option v alues for all users on a system. See Section 4.2.2 [Configuration file
precedence], page 268, for a more detailed discussion.
W e will stop the discussion on configuration files here, but y ou can alw a ys read ab out
them in Section 4.2 [Configuration files], page 267. Before con tin uing the tutorial, let’s
delete the t w o extra directories that w e do not need an y more:
$ rm -rf my-cosmology*
2.1.9 W arping to a new pixel grid
W e are no w ready to start pro cessing the deep HST images that w ere prepared in Sec-
tion 2.1.4 [Dataset insp ection and cropping], page 25. One of the most imp ortan t p oints
Chapter 2: T utorials 38
while using sev eral images for data pro cessing is that those images m ust ha v e the same pixel
grid. The pro cess of changing the pixel grid is named ‘w arp’. F ortunately , Gn uastro has
W arp program for w arping the pixel grid (see Section 6.4 [W arp], page 492).
W arping to a differen t/matched pixel grid is commonly needed b efore higher-lev el anal-
ysis esp ecially when y ou are using datasets from differen t instrumen ts. The XDF datasets
w e are using here are already aligned to the same pixel grid. But let’s ha v e a lo ok at some
of Gn uastro’s linear w arping features here. F or example, try rotating one of the images b y
20 degrees with the first command b elo w. With the second command, op en the output and
input to see ho w it is rotated.
$ astwarp flat-ir/xdf-f160w.fits --rotate=20
$ astscript-fits-view flat-ir/xdf-f160w.fits xdf-f160w_rotated.fits
W arp can generally b e used for man y kinds of pixel grid manipulation (w arping), not
just rotations. F or example, the outputs of the commands b elo w will ha v e larger pixels
resp ectiv ely (new resolution b eing one quarter the original resolution), get shifted b y 2.8
(b y sub-pixel), get a shear of 2, and b e tilted (pro jected). Run eac h of them and op en the
output file to see the effect, they will b ecome handy for y ou in the future.
$ astwarp flat-ir/xdf-f160w.fits --scale=0.25
$ astwarp flat-ir/xdf-f160w.fits --translate=2.8
$ astwarp flat-ir/xdf-f160w.fits --shear=0.2
$ astwarp flat-ir/xdf-f160w.fits --project=0.001,0.0005
$ astscript-fits-view flat-ir/xdf-f160w.fits *.fits
If y ou need to do m ultiple warps, y ou can com bine them in one call to W arp. F or example,
to first rotate the image, then scale it, run this command:
$ astwarp flat-ir/xdf-f160w.fits --rotate=20 --scale=0.25
If y ou ha v e m ultiple w arps, do them all in one command. Do not w arp them in separate
commands b ecause the correlated noise will b ecome to o strong. As y ou see in the matrix
that is prin ted when y ou run W arp, it merges all the warps in to a single w arping matrix
(see Section 6.4.2 [Merging m ultiple w arpings], page 495) and simply applies that (mixes the
pixel v alues) just once. How ev er, if y ou run W arp m ultiple times, the pixels will b e mixed
m ultiple times, creating a strong artificial blur/smo othing, or stronger correlated noise.
Recall that the merging of m ultiple warps is done through matrix m ultiplication, there-
fore order matters in the separate op erations. A t a lo w er lev el, through W arp’s --matrix
option, y ou can directly request y our desired final w arp and do not ha v e to break it up in to
differen t w arps lik e ab o v e (see Section 6.4.4 [In v oking W arp], page 498).
F ortunately these datasets are already aligned to the same pixel grid, so y ou do not
actually need the files that w ere just generated. Y ou can safely delete them all with the
follo wing command. Here, y ou see wh y w e put the pro cessed outputs that w e need later
in to a separate directory . In this w a y , the top directory can b e used for temp orary files for
testing that y ou can simply delete with a generic command lik e b elo w.
$ rm *.fits
2.1.10 NoiseChisel and Multi-Extension FITS files
In the previous sections, w e completed a review of the basics of Gn uastro’s programs.
W e are no w ready to do some more serious analysis on the do wnloaded images: extract
Chapter 2: T utorials 39
the pixels con taining signal from the image, find sub-structure of the extracted signal,
do measuremen ts o v er the extracted ob jects and analyze them (finding certain ob jects of
in terest in the image).
The first step is to separate the signal (galaxies or stars) from the bac kground noise in
the image. W e will b e using the results of Section 2.1.4 [Dataset insp ection and cropping],
page 25, so b e sure y ou already ha v e them. Gn uastro has NoiseChisel for this job. But
NoiseChisel’s output is a m ulti-extension FITS file, therefore to b etter understand ho w to
use NoiseChisel, let’s tak e a lo ok at m ulti-extension FITS files and ho w y ou can in teract
with them.
In the FITS format, eac h extension con tains a separate dataset (image in this case). Y ou
can get basic information ab out the extensions in a FITS file with Gn uastro’s Fits program
(see Section 5.1 [Fits], page 294). T o start with, let’s run NoiseChisel without an y options,
then use Gn uastro’s Fits program to insp ect the n um b er of extensions in this file.
$ astnoisechisel flat-ir/xdf-f160w.fits
$ astfits xdf-f160w_detected.fits
F rom the output list, w e see that NoiseChisel’s output contains 5 extensions. The zero-
th (coun ting from zero, with name NOISECHISEL-CONFIG ) is empt y: it has v alue of 0 in
the fourth column (whic h sho ws its size in pixels). Lik e NoiseChisel, in all of Gnuastro’s
programs, the first (or zero-th) extension of the output only con tains meta-data: data
ab out/describing the datasets within (all) the output’s extensions. This is recommended
b y the FITS standard, see Section 5.1 [Fits], page 294, for more. In the case of Gn uastro’s
programs, this generic zero-th/meta-data extension (for the whole file) con tains all the
configuration options of the program that created the file.
Metadata regarding ho w the analysis was done (or a dataset was created) is very impor-
tan t for higher-lev el analysis and repro ducibilit y . Therefore, Let’s first tak e a closer lo ok
at the NOISECHISEL-CONFIG extension. If y ou sp ecify a sp ecial header in the FITS file,
Gn uastro’s Fits program will prin t the header k eyw ords (metadata) of that extension. Y ou
can either sp ecify the HDU/extension coun ter (starting from 0), or name. Therefore, the
t w o commands b elo w are iden tical for this file. W e are usually tempted to use the first
(shorter format), but when putting y our commands in to a script, please use the second
format whic h is more h uman-friendly and understandable for readers of you r co de who ma y
not kno w what is in the 0-th extension (this includes y ourself in a few mon ths!):
$ astfits xdf-f160w_detected.fits -h0
$ astfits xdf-f160w_detected.fits -hNOISECHISEL-CONFIG
The first group of FITS header k eywords y ou see (con taining the SIMPLE and BITPIX
k eyw ords; b efore the first empt y line) are standard k eywords. They are required b y the FITS
standard and m ust b e present in an y FITS extension. The second group starts with the
input file name (v alue to the INPUT k eyw ord). The rest of the k eyw ords y ou see afterw ards
ha v e the same name as NoiseChisel’s options, and the v alue used b y NoiseChisel in this run
is sho wn after the = sign. Finally , the last group (starting with DATE ) con tains the date and
v ersion information of Gn uastro and its dep endencies that w ere used to generate this file.
Besides the option v alues, these are also critical for future repro ducibilit y of the result (y ou
ma y up date Gn uastro or its dep endencies, and they ma y b eha v e differen tly afterw ards).
The “v ersions and date” group of k eyw ords are presen t in all Gn uastro’s FITS extension
outputs, for more see Section 4.10 [Output FITS files], page 290.
Chapter 2: T utorials 40
Note that if a k eyw ord name is larger than 8 c haracters, it is preceded b y a HIERARCH
k eyw ord and that all k eyw ord names are in capital letters. These are all part of the FITS
standard and originate from its history . But in short, b oth can b e ignored! F or example,
with the commands b elo w, let’s see at first what the default v alues are, and then just c hec k
the v alue of --detgrowquant option (using the -P option describ ed in Section 2.1.8 [Option
managemen t and configuration files], page 35).
$ astnoisechisle -P
$ astnoisechisel -P | grep detgrowquant
T o confirm that NoiseChisel used this v alue when w e ran it ab o v e, let’s use grep to
extract the k eyword line with detgrowquant from the metadata extension. Ho wev er, as
y ou sa w ab o v e, k eyw ord names in the header is in all caps. So we need to ask grep to
ignore case with the -i option.
$ astfits xdf-f160w_detected.fits -h0 | grep -i detgrowquant
In the output of the ab o v e command, y ou see HIERARCH at the start of the line. According
to the FITS standard, HIERARCH is placed at the start of all keyw ords that ha ve a name
that is more than 8 c haracters long. Both the all-caps and the HIERARCH k eyw ord can b e
anno ying when y ou w an t to read/c hec k the v alue. Therefore, the b est solution is to use
the --keyvalue option of Gn uastro’s astfits program as sho wn b elo w. With it, y ou do
not ha v e to w orry ab out HIERARCH or the case of the name (FITS k eyw ord names are not
case-sensitiv e).
$ astfits xdf-f160w_detected.fits -h0 --keyvalue=detgrowquant -q
The metadata (that is stored in the output) can later b e used to exactly repro-
duce/understand y our result, ev en if y ou ha v e lost/forgot the command y ou used to create
the file. This feature is presen t in all of Gn uastro’s programs, not just NoiseChisel.
The rest of the HDUs in NoiseChisel ha v e data. So let’s op en them in a DS9 windo w
and then describ e eac h:
$ astscript-fits-view xdf-f160w_detected.fits
A “cub e” windo w op ens along with DS9’s main windo w. The buttons and horizon tal
scroll bar in this small new windo w can b e used to na vigate b et w een the extensions. In
this mo de, all DS9’s settings (for example, zo om or color-bar) will b e iden tical b et w een
the extensions. T ry zo oming in to one part and flipping through the extensions to see ho w
the galaxies w ere detected along with the Sky and Sky standard deviation v alues for that
region. Just ha v e in mind that NoiseChisel’s job is only detection (separating signal from
noise). W e will do segmen tation on this result later to find the individual galaxies/p eaks
o v er the detected pixels.
The second extension of NoiseChisel’s output (n um b ered 1, named INPUT-NO-SKY ) is
the Sky-subtracted input that y ou pro vided. The third ( DETECTIONS ) is NoiseChisel’s main
output whic h is a binary image with only t wo possible v alues for all pixels: 0 for noise and
1 for signal. Since it only has t w o v alues, to a v oid taking to o m uch space on y our computer,
its n umeric datat yp e an unsigned 8-bit in teger (or uint8 ) 14 . The fourth and fifth ( SKY and
SKY_STD ) extensions, ha v e the Sky and its standard deviation v alues for the input on a tile
grid and w ere calculated o v er the undetected regions (for more on the imp ortance of the
Sky v alue, see Section 7.1.4 [Sky v alue], page 520).
14 T o learn more ab out n umeric data types see Section 4.5 [Numeric data t yp es], page 276.
Chapter 2: T utorials 41
Eac h HDU/extension in a FITS file is an indep enden t dataset (image or table) whic h y ou
can delete from the FITS file, or cop y/cut to another file. F or example, with the command
b elo w, y ou can cop y NoiseChisel’s DETECTIONS HDU/extension to another file:
$ astfits xdf-f160w_detected.fits --copy=DETECTIONS -odetections.fits
There are similar options to con venien tly cut ( --cut , cop y , then remov e from the input)
or delete ( --remove ) HDUs from a FITS file also. See Section 5.1.1.1 [HDU information
and manipulation], page 298, for more.
2.1.11 NoiseChisel optimization for detection
In Section 2.1.10 [NoiseChisel and Multi-Extension FITS files], page 38, w e ran NoiseChisel
and review ed NoiseChisel’s output format. No w that y ou ha v e a b etter feeling for m ulti-
extension FITS files, let’s optimize NoiseChisel for this particular dataset.
One go o d w a y to see if y ou ha v e missed an y signal (small galaxies, or the wings of
brigh ter galaxies) is to mask all the detected pixels and insp ect the noise pixels. F or this, you
can use Gn uastro’s Arithmetic program (in particular its where op erator, see Section 6.2.4
[Arithmetic op erators], page 407). The command b elo w will pro duce mask-det.fits . In it,
all the pixels in the INPUT-NO-SKY extension that are flagged 1 in the DETECTIONS extension
(dominated b y signal, not noise) will b e set to NaN.
Since the v arious extensions are in the same file, for eac h dataset w e need the file and
extension name. T o mak e the command easier to read/write/understand, let’s use shell
v ariables: ‘ in ’ will b e used for the Sky-subtracted input image and ‘ det ’ will b e used for
the detection map. Recall that a shell v ariable’s v alue can b e retriev ed by addin g a $ b efore
its name, also note that the double quotations are necessary when w e hav e white-space
c haracters in a v ariable v alue (lik e this case).
$ in="xdf-f160w_detected.fits -hINPUT-NO-SKY"
$ det="xdf-f160w_detected.fits -hDETECTIONS"
$ astarithmetic $in $det nan where --output=mask-det.fits
T o in vert the result (only k eep the detected pixels), y ou can flip the detection map (from 0
to 1 and vice-v ersa) b y adding a ‘ not ’ after the second $det :
$ astarithmetic $in $det not nan where --output=mask-sky.fits
Lo ok again at the DETECTIONS extension, in particular the long w orm-like structure
around 15 pixel 1650 (X) and 1470 (Y). These t yp es of long wiggly structures sho w that w e
ha v e dug to o deep in to the noise, and are a signature of correlated noise. Correlated noise
is created when w e w arp (for example, rotate) individual exp osures (that are eac h sligh tly
offset compared to eac h other) in to the same pixel grid b efore adding them in to one deep er
image. During the w arping, nearb y pixels are mixed and the effect of this mixing on the
noise (whic h is in ev ery pixel) is called “correlated noise”. Correlated noise is a form of
con v olution and it sligh tly smo oths the image.
15 T o find a particular co ordian te easily in DS9, you can do this: Clic k on the “Edit” menu, and select
“Region”. Then click on any random part of the image to see a circle sho w up in that lo cation (this is
the “region”). Double-clic k on the region and a “Circle” window will open. If y ou hav e celestial co ordi-
nates, k eep the default “fk5” in the scroll-down men u after the “Cen ter”. But if y ou hav e pixel/image
co ordinates, click on the “fk5” and select “Image”. Now y ou can set the “Cen ter” co ordinates of the
region ( 1650 and 1470 in this case) b y man ually typing them in the t wo boxes in fron t of “Cen ter”. Fi-
nally , when ev erything is ready , clic k on the “Apply” button and y our region will go o ver y our requested
co ordinates. Y ou can zo om out (to see the whole image) and visually find it.
Chapter 2: T utorials 48
So let’s select the prop erties w e w an t to measure in this tutorial. First of all, we need
to kno w whic h measuremen t b elongs to whic h ob ject or clump, so w e will start with the
--ids (read as: IDs 16 ). W e also w an t to measure (in this order) the Righ t Ascension (with
--ra ), Declination ( --dec ), magnitude ( --magnitude ), and signal-to-noise ratio ( --sn ) of
the ob jects and clumps. F urthermore, as men tioned ab o v e, w e also w an t measuremen ts
on clumps, so w e also need to call --clumpscat . The follo wing command will mak e these
measuremen ts on Segmen t’s F160W output and write them in a catalog for eac h ob ject and
clump in a FITS table. F or more on the zero p oin t, see Section 7.4.2 [Brigh tness, Flux,
Magnitude and Surface brigh tness], page 576.
$ mkdir cat
$ astmkcatalog seg/xdf-f160w.fits --ids --ra --dec --magnitude --sn \
--zeropoint=25.94 --clumpscat --output=cat/xdf-f160w.fits
F rom the prin ted statemen ts on the command-line, y ou see that Mak eCatalog read all the
extensions in Segmen t’s output for the v arious measuremen ts it needed. Let’s lo ok at the
output of the command ab o v e:
$ astfits cat/xdf-f160w.fits
Y ou will see that the output of the Mak eCatalog has t w o extensions. The first extension
sho ws the measuremen ts o v er the OBJECTS , and the second extension shows the measure-
men ts o v er the clumps CLUMPS .
T o calculate colors, w e also need magnitude measuremen ts on the other filters. So let’s
rep eat the command ab o v e on them, just c hanging the file names and zero p oin t (whic h w e
got from the XDF surv ey w eb page):
$ astmkcatalog seg/xdf-f125w.fits --ids --ra --dec --magnitude --sn \
--zeropoint=26.23 --clumpscat --output=cat/xdf-f125w.fits
$ astmkcatalog seg/xdf-f105w.fits --ids --ra --dec --magnitude --sn \
--zeropoint=26.27 --clumpscat --output=cat/xdf-f105w.fits
Ho w ev er, the galaxy prop erties migh t differ b et w een the filters (whic h is the whole
purp ose b ehind observing in differen t filters!). Also, the noise prop erties and depth of the
datasets differ. Y ou can see the effect of these factors in the resulting clump catalogs, with
Gn uastro’s T able program. W e will go deep in to w orking with tables in the next section,
but in summary: the -i option will prin t information ab out the columns and n um b er of
ro ws. T o see the column v alues, just remo v e the -i option. In the output of eac h command
b elo w, lo ok at the Number of rows: , and note that they are differen t.
$ asttable cat/xdf-f105w.fits -hCLUMPS -i
$ asttable cat/xdf-f125w.fits -hCLUMPS -i
$ asttable cat/xdf-f160w.fits -hCLUMPS -i
Matc hing the catalogs is p ossible (for example, with Section 7.5 [Matc h], page 620).
Ho w ev er, the measuremen ts of each column are also done on differen t pixels: the clump
lab els can/will differ from one filter to another for one ob ject. Please op en them and fo cus
on one ob ject to see for y ourself. This can bias the result, if y ou matc h catalogs.
16 This option is plural b ecause w e need tw o ID columns for iden tifying “clumps” in the clumps cata-
log/table: the first column will b e the ID of the host “ob ject”, and the second one will b e the ID of the
clump within that ob ject. In the “ob jects” catalog/table, only a single column will b e returned for this
option.
Chapter 2: T utorials 49
An accurate color calculation can only b e done when magnitudes are measured from
the same pixels on all images and this can b e done easily with Mak eCatalog. In fact this
is one of the reasons that NoiseChisel or Segmen t do not generate a catalog like most
other detection/segmen tation soft w are. This giv es you the freedom of selecting the pixels
for measuremen t in an y w a y y ou lik e (from other filters, other soft w are, man ually , etc.).
F ortunately in these images, the P oin t spread function (PSF) is v ery similar, allo wing us to
use a single lab eled image output for all filters 17 .
The F160W image is deep er, th us pro viding b etter detection/segmen tation, and redder,
th us observing smaller/older stars and represen ting more of the mass in the galaxies. W e will
th us use the F160W filter as a reference and use its segmen t lab els to iden tify whic h pixels
to use for whic h ob jects/clumps. But we will do the measuremen ts on the sky-subtracted
F105W and F125W images (using Mak eCatalog’s --valuesfile option) as sho wn b elo w:
Notice that the only difference b et ween these calls and the call to generate the ra w F160W
catalog (excluding the zero p oin t and the output name) is the --valuesfile .
$ astmkcatalog seg/xdf-f160w.fits --ids --ra --dec --magnitude --sn \
--valuesfile=nc/xdf-f125w.fits --zeropoint=26.23 \
--clumpscat --output=cat/xdf-f125w-on-f160w-lab.fits
$ astmkcatalog seg/xdf-f160w.fits --ids --ra --dec --magnitude --sn \
--valuesfile=nc/xdf-f105w.fits --zeropoint=26.27 \
--clumpscat --output=cat/xdf-f105w-on-f160w-lab.fits
After running the commands ab o ve, lo ok in to what Mak eCatalog prin ted on the
command-line. Y ou can see that (as requested) the ob ject and clump pixel lab els in b oth
w ere tak en from the resp ectiv e extensions in seg/xdf-f160w.fits . Ho wev er, the pixel
v alues and pixel Sky standard deviation w ere resp ectiv ely tak en from nc/xdf-f105w.fits
and nc/xdf-f125w.fits . Since w e used the same lab eled image on all filters, the num b er
of ro ws in b oth catalogs are no w iden tical. Let’s ha v e a lo ok:
$ asttable cat/xdf-f105w-on-f160w-lab.fits -hCLUMPS -i
$ asttable cat/xdf-f125w-on-f160w-lab.fits -hCLUMPS -i
$ asttable cat/xdf-f160w.fits -hCLUMPS -i
Finally , MakeCatalog also does basic calculations on the full dataset (indep enden t of
eac h lab eled region but related to whole data), for example, pixel area or p er-pixel surface
brigh tness limit. They are stored as k eyw ords in the FITS headers (or lines starting with
# in plain text). This (and other w a ys to measure the limits of y our dataset) are discussed
in the next section: Section 2.1.14 [Measuring the dataset limits], page 49.
2.1.14 Measuring the dataset limits
In Section 2.1.13 [Segmen tation and making a catalog], page 47, w e created a catalog of the
differen t ob jects with the image. Before measuring colors, or doing an y other kind of analysis
on the catalogs (and detected ob jects), it is v ery imp ortant to understand the limitations
of the dataset. Without understanding the limitations of y our dataset, y ou cannot mak e
an y ph ysical in terpretation of y our results. The theory b ehind the calculations discussed
here is thoroughly in tro duced in Section 7.4.3 [Quan tifying measuremen t limits], page 581.
17 When the PSFs b et w een t wo images differ largely , you w ould ha ve to PSF-matc h the images b efore using
the same pixels for measuremen ts.
Chapter 2: T utorials 50
F or example, with the command b elo w, let’s sort all the detected clumps in the image
b y magnitude (with --sort=magnitude ) and and prin t the magnitude and signal-to-noise
ratio (S/N; with -cmagnitude,sn ):
$ asttable cat/xdf-f160w.fits -hclumps -cmagnitude,sn \
--sort=magnitude --noblank=magnitude
As y ou see, w e hav e clumps with a total magnitude of almost 32! This is extr emely
faint ! Are these things trustable? Let’s ha v e a lo ok at all of those with a magnitude
b et w een 31 and 32 with the command b elo w. W e are first using T able to only k eep the
relev an t columns ro ws, and using Gn uastro’s DS9 region file creation script ( astscript-
ds9-region ) to generate DS9 region files, and op en DS9:
$ asttable cat/xdf-f160w.fits -hclumps -cra,dec \
--range=magnitude,31:32 \
| astscript-ds9-region -c1,2 --radius=0.5 \
--command="ds9 -mecube seg/xdf-f160w.fits -zscale"
Zo om-out a little and y ou will see some green circles (DS9 region files) in some regions
of the image. There actually do es seem to b e a true p eak under the selected regions, but as
y ou see, they are v ery small, diffuse and noisy . Ho w reliable are the measured magnitudes?
Using the S/N column from the first command ab o v e, y ou can see that suc h ob jects only
ha v e a signal to noise of ab out 2.6 (whic h is indeed to o lo w for most analysis purp oses)
$ asttable cat/xdf-f160w.fits -hclumps -csn \
--range=magnitude,31:32 | aststatistics
This brings us to the first metho d of quan tifying y our dataset’s magnitude limit , whic h
is also sometimes called dete ction limit (see Section 7.4.3.6 [Magnitude limit of image],
page 587). T o estimate the 5 σ detection limit of y our dataset, y ou simply rep ort the
median magnitude of the ob jects that ha ve a signal to noise of (appro ximately) fiv e. This
is v ery easy to calculate with the command b elo w:
$ asttable cat/xdf-f160w.fits -hclumps --range=sn,4.8:5.2 -cmagnitude \
| aststatistics --median
29.9949
Let’s ha v e a lo ok at these ob jects, to get a feeling of what these clump lo oks lik e:
$ asttable cat/xdf-f160w.fits -hclumps --range=sn,4.8:5.2 \
-cra,dec,magnitude \
| astscript-ds9-region -c1,2 --namecol=3 \
--width=2 --radius=0.5 \
--command="ds9 -mecube seg/xdf-f160w.fits -zscale"
The n um b er y ou see on top of eac h region is the clump’s magnitude. Please go o v er the
ob jects and ha ve a close look at them! It is v ery imp ortan t to ha v e a feeling of what y our
dataset lo oks lik e, and ho w to in terpret the n um b ers to asso ciate an image with them.
Generally , they lo ok v ery small with differen t lev els of diffuse-ness! Those that are
sharp er mak e more visual sense (to b e 5 σ detections), but the more diffuse ones extend
o v er a larger area. F urthermore, the noise is measured on individual pixel measuremen ts.
Ho w ev er, during the reduction man y exp osures are co-added and stac k ed, mixing the pixels
lik e a small con v olution (creating “correlated noise”). Therefore you clearly see t w o main
issues with the detection limit as defined ab o v e: it dep ends on the morphology , and it do es
not tak e in to accoun t the correlated noise.
Chapter 2: T utorials 51
A more realistic w a y to estimate the significance of the detection is to tak e its fo otprint,
randomly place it in thousands of undetected regions of the image and use that distribution
as a reference. This is tec hnically kno wn as upp er-limit measuremen ts. F or a full discussion,
see Section 7.4.3.5 [Upp er limit magnitude of eac h detection], page 586).
Since it is for eac h separate ob ject, the upp er-limit measuremen ts should b e requested
as extra columns in Mak eCatalog’s output. F or example, with the command b elow, let’s
generate a new catalog of the F160W filter, but with t w o extra columns compared to the
one in cat/ : the upp er-limit magnitude and the upp er-limit m ultiple of sigma.
$ astmkcatalog seg/xdf-f160w.fits --ids --ra --dec --magnitude --sn \
--zeropoint=25.94 --clumpscat --upnsigma=3 \
--upperlimit-mag --upperlimit-sigma \
--output=xdf-f160w.fits
Let’s compare the upp er-limit magnitude with the measured magnitude of eac h clump:
$ asttable xdf-f160w.fits -hclumps -cmagnitude,upperlimit_mag
As y ou see, in almost all of the cases, the measured magnitude is sufficiently higher than
the upp er-limit magnitude. Let’s subtract the latter from the former to b etter see this
difference in a third column:
$ asttable xdf-f160w.fits -hclumps -cmagnitude,upperlimit_mag \
-c'arith upperlimit_mag magnitude -'
The ones with a p ositiv e third column (difference) show that the clump has sufficien tly
higher brigh tness than the noisy bac kground to b e usable. Let’s use T able’s Section 5.3.3
[Column arithmetic], page 345, to find only those that ha v e a negativ e difference:
$ asttable xdf-f160w.fits -hclumps -cra,dec --noblankend=3 \
-c'arith upperlimit_mag magnitude - set-d d d 0 gt nan where'
F rom more than 3500 clumps, this command only ga v e ∼ 150 ro ws (this n um b er ma y
sligh tly c hange on differen t runs due to the random nature of the upp er-limit sampling 18 )!
Let’s ha v e a lo ok at them:
$ asttable xdf-f160w.fits -hclumps -cra,dec --noblankend=3 \
-c'arith upperlimit_mag magnitude - set-d d d 0 gt nan where' \
| astscript-ds9-region -c1,2 --namecol=3 --width=2 \
--radius=0.5 \
--command="ds9 -mecube seg/xdf-f160w.fits -zscale"
Y ou see that they are all extremely fain t and diffuse/small p eaks. Therefore, if an ob ject’s
magnitude is fain ter than its upp er-limit magnitude, y ou should not use the magnitude: it
is not accurate! Y ou should use the upp er-limit magnitude instead (with an arro w in y our
plots to mark whic h ones are upp er-limits).
But the main p oin t (in relation to the magnitude limit) with the upp er-limit, is
the UPPERLIMIT_SIGMA column. y ou can think of this as a r e alistic S/N for extremely
fain t/diffuse/small ob jects). The ra w S/N column is simply calculated on a pixel-b y-pixel
basis, ho w ev er, the upp er-limit sigma is pro duced b y actually taking the lab el’s fo otprin t,
and randomly placing it thousands of time o v er undetected parts of the image and
measuring the brigh tness of the sky . The clump’s brigh tness is then divided b y the
18 Y ou can fix the random num b er generator seed, so you alw ays get the same sampling, see Section 6.2.3.4
[Generating random n um b ers], page 405.
Chapter 2: T utorials 52
standard deviation of the resulting distribution to giv e y ou exactly ho w significan t it is
(accoun ting for in ter-pixel issues lik e correlated noise, whic h are strong in this dataset).
Y ou can actually compare the t w o v alues with the command b elo w:
$ asttable xdf-f160w.fits -hclumps -csn,upperlimit_sigma
As y ou see, the second column (upp er-limit sigma) is almost alw a ys less than the S/N.
This clearly sho ws the effect of correlated noise! If you no w use this column as the reference
for deriving the magnitude limit, you will see that it will shift b y almost 0.5 magnitudes
brigh ter and is no w more reasonable:
$ asttable xdf-f160w.fits -hclumps --range=upperlimit_sigma,4.8:5.2 \
-cmagnitude | aststatistics --median
29.6257
W e see that the 5 σ detection limit is ∼ 29 . 6! This is extremely deep! F or example,
in the Legacy Surv ey 19 , the 5 σ detection limit for p oint sour c es is appro ximately 24.5 (5
magnitudes, or 100 times, shallo w er than this image).
As men tioned ab o v e, an imp ortant ca v eat in this simple calculation is that w e should
only b e lo oking at p oin t-like ob jects, not simply ev erything. This is b ecause the shap e
or radial slop e of the profile has an imp ortan t effect on this measuremen t: at the same
total magnitude, a sharp er ob ject will hav e a higher S/N. T o b e more precise, w e should
first p erform star-galaxy separation, then do this only for the ob jects that are classified as
stars. A crude, first-order, metho d is to use the --axis-ratio option so Mak eCatalog also
measures the axis ratio, then call T able with --range=upperlimit_sigma,,4.8:5.2 and
--range=axis_ratio,0.95:1 (in one command). Please do this for y ourself as an exercise
to see the difference with the result ab o v e.
Before con tin uing, let’s remo v e this temp orarily pro duced catalog:
$ rm xdf-f160w.fits
Another measure of the dataset’s limit is the completeness limit (Section 7.4.3.4 [Com-
pleteness limit of eac h detection], page 586). This is necessary when y ou are lo oking at
p opulations of ob jects o v er the image. Y ou w an t to kno w un til what magnitude you can be
sure that y ou ha v e detected an ob ject (if it w as presen t). As describ ed in Section 7.4.3.4
[Completeness limit of eac h detection], page 586, the b est w a y to do this is with mo c k
images. But a crude, first order result can b e obtained from the actual image: by simply
plotting the histogram of the magnitudes:
$ aststatistics cat/xdf-f160w.fits -hclumps -cmagnitude
...
Histogram:
| *
| ** ****
| ***********
| *************
| ****************
| *******************
| **********************
| **************************
19 https://www.legacysurvey.org/dr9/description
Chapter 2: T utorials 53
| *********************************
| *********************************************
|* * ** ** **********************************************************
|----------------------------------------------------------------------
This plot (the histogram of magnitudes; where fain ter magnitudes are to w ards the righ t)
is tec hnically called the dataset’s numb er c ount plot. Y ou see that the num b er of ob jects
increases with magnitude as the magnitudes get fain ter (to the righ t). Ho w ev er, b ey ond
a certain magnitude, y ou see it b ecomes flat, and so on afterw ards, the n um b ers suddenly
drop.
Once y ou ha v e y our catalog, you can easily find this point with the t w o commands b elo w.
First w e generate a histogram with few er bins (to ha v e more n um b ers in eac h bin). W e
then use A WK to find the magnitude bin where the n umber of p oints decrease compared
to the previous bin. But we only do this for bins that ha v e more than 50 items (to a v oid
scatter in the brigh t end). Finally , in Statistics, w e ha v e man ually set the magnitude range
and n um b er of bins so eac h bin is roughly 0.5 magnitudes thic k (with --greaterequal=20 ,
--lessthan=32 and --numbins=24 )
$ aststatistics cat/xdf-f160w.fits -hclumps -cmagnitude --histogram \
--greaterequal=20 --lessthan=32 --numbins=24 \
--output=f160w-hist.txt
$ asttable f160w-hist.txt \
| awk '$2>50 && $2<prev{print prevbin; exit} \
{prev=$2; prevbin=$1}'
28.932122667631
Therefore, to first order (and v ery crudely!) w e can sa y that if an ob ject is in our field
of view and has a magnitude of ∼ 29 or brigh ter, w e can b e highly confiden t that w e ha v e
detected it. But b efore con tin uing, let’s clean up b ehind ourselv es:
$ rm f160w-hist.txt
Another imp ortan t limiting parameter in a pro cessed dataset is the surface brigh tness
limit (Section 7.4.3.7 [Surface brigh tness limit of image], page 588). The surface brigh tness
limit of a dataset is an imp ortan t measure for extended structures (for example, when y ou
w an t to lo ok at the outskirts of galaxies). In the next tutorial, we ha v e thoroughly describ ed
the deriv ation of the surface brigh tness limit of a dataset. So we will just sho w the final
result here, and encourage y ou to follo w up with that tutorial after finishing this tutorial
(see Section 2.2.4 [Image surface brigh tness limit], page 92)
By default, MakeCatalog will estimate the surface brigh tness limit of a giv en dataset,
and put it in the k eyw ords of the output (all k eyw ords starting with SBL , whic h is short for
surface brigh tness limit):
$ astfits cat/xdf-f160w.fits -h1 | grep SBL
As y ou see, the only one with a unit of mag/arcsec^2 is SBLMAG . It con tains the surface
brigh tness limit of the input dataset o ver SBLAREA arcsec 2 with SBLNSIG m ultiples of σ .
In the curren t v ersion of Gn uastro, SBLAREA=100 and SBLNSIG=3 , so the surface brigh tness
limit of this image is 32.66 mag/arcsec 2 (3 σ , o ver 100 arcsec 2 ). Therefore, if this default area
and m ultiple of sigma are fine for y ou 20 (these are the most commonly used v alues), y ou can
20 Y ou can change these v alues with the --sfmagarea and --sfmagnsigma
Chapter 2: T utorials 54
simply read the image surface brigh tness limit from the catalogs pro duced b y Mak eCatalog
with this command:
$ astfits cat/*.fits -h1 --keyvalue=SBLMAG
2.1.15 W orking with catalogs (estimating colors)
In the previous step w e generated catalogs of ob jects and clumps o v er our dataset (see
Section 2.1.13 [Segmen tation and making a catalog], page 47). The catalogs are a v ailable
in the t w o extensions of the single FITS file 21 . Let’s see the extensions and their basic
prop erties with the Fits program:
$ astfits cat/xdf-f160w.fits # Extension information
Let’s insp ect the table in eac h extension with Gn uastro’s T able program (see Section 5.3
[T able], page 339). W e should ha v e used -hOBJECTS and -hCLUMPS instead of -h1 and -h2
resp ectiv ely . The n um b ers are just used here to con v ey that b oth names or n um b ers are
p ossible, in the next commands, w e will just use names.
$ asttable cat/xdf-f160w.fits -h1 --info # Objects catalog info.
$ asttable cat/xdf-f160w.fits -h1 # Objects catalog columns.
$ asttable cat/xdf-f160w.fits -h2 -i # Clumps catalog info.
$ asttable cat/xdf-f160w.fits -h2 # Clumps catalog columns.
As y ou see ab o ve, when giv en a sp ecific table (file name and extension), T able will prin t the
full con ten ts of all the columns. T o see the basic metadata ab out each column (for example,
name, units and commen ts), simply app end a --info (or -i ) to the command.
T o prin t the conten ts of sp ecial column(s), just giv e the column num b er(s) (coun ting
from 1 ) or the column name(s) (if they ha ve one) to the --column (or -c ) option. F or
example, if you just w an t the magnitude and signal-to-noise ratio of the clumps (in the
clumps catalog), y ou can get it with an y of the follo wing commands
$ asttable cat/xdf-f160w.fits -hCLUMPS --column=5,6
$ asttable cat/xdf-f160w.fits -hCLUMPS -c5,SN
$ asttable cat/xdf-f160w.fits -hCLUMPS -c5 -c6
$ asttable cat/xdf-f160w.fits -hCLUMPS -cMAGNITUDE -cSN
Similar to HDUs, when the columns ha ve names, alw a ys use the name: it is so common to
mis-write n um b ers or forget the order later! Using column names instead of n um b ers has
man y adv an tages:
1. Y ou do not ha v e to w orry ab out the order of columns in the table.
2. It acts as a do cumen tation in the script.
3. Column meta-data (including a name) are not just limited to FITS tables and can also
b e used in plain text tables, see Section 4.7.2 [Gn uastro text table format], page 284.
T able also has to ols to limit the displa y ed ro ws. F or example, with the first command
b elo w only ro ws with a magnitude in the range of 29 to 30 will b e sho wn. With the second
command, y ou can further limit the displa y ed ro ws to ro ws with an S/N larger than 10 (a
21 Mak eCatalog can also output plain text tables. Ho w ever, in the plain text format y ou can only ha ve
one table p er file. Therefore, if y ou also request measurements on clumps, t wo plain text tables will be
created (suffixed with _o.txt and _c.txt ).
Chapter 2: T utorials 55
range b et w een 10 to infinit y). Y ou can further sort the output ro ws, only sho w the top (or
b ottom) N ro ws, etc., see Section 5.3 [T able], page 339, for more.
$ asttable cat/xdf-f160w.fits -hCLUMPS --range=MAGNITUDE,28:29
$ asttable cat/xdf-f160w.fits -hCLUMPS \
--range=MAGNITUDE,28:29 --range=SN,10:inf
No w that y ou are comfortable in viewing table columns and ro ws, let’s lo ok in to merging
columns of m ultiple tables into one table (whic h is necessary for measuring the color of
the clumps). Since cat/xdf-f160w.fits and cat/xdf-f105w-on-f160w-lab.fits ha v e
exactly the same n um b er of ro ws and the ro ws corresp ond to the same clump, let’s merge
them to ha v e one table with magnitudes in b oth filters.
W e can merge columns with the --catcolumnfile option lik e b elow. Y ou giv e this
option a file name (whic h is assumed to b e a table that has the same n um b er of ro ws as
the main input), and all the table’s columns will b e concatenated/app ended to the main
table. No w, try it out with the commands b elo w. W e will first lo ok at the metadata of the
first table (only the CLUMPS extension). With the second command, w e will concatenate the
t w o tables and write them in, two-in-one.fits and finally , we will c hec k the new catalog’s
metadata.
$ asttable cat/xdf-f160w.fits -i -hCLUMPS
$ asttable cat/xdf-f160w.fits -hCLUMPS --output=two-in-one.fits \
--catcolumnfile=cat/xdf-f125w-on-f160w-lab.fits \
--catcolumnhdu=CLUMPS
$ asttable two-in-one.fits -i
By comparing the t wo metadata, w e see that b oth tables ha v e the same n um b er of ro ws.
But what migh t ha v e attracted y our atten tion more, is that two-in-one.fits has double
the n um b er of columns (as exp ected, after all, y ou merged b oth tables in to one file, and did
not ask for an y sp ecific column). In fact y ou can concatenate an y n um b er of other tables
in one command, for example:
$ asttable cat/xdf-f160w.fits -hCLUMPS --output=three-in-one.fits \
--catcolumnfile=cat/xdf-f125w-on-f160w-lab.fits \
--catcolumnfile=cat/xdf-f105w-on-f160w-lab.fits \
--catcolumnhdu=CLUMPS --catcolumnhdu=CLUMPS
$ asttable three-in-one.fits -i
As y ou see, to a v oid confusion in column names, T able has in ten tionally app ended a -1
to the column names of the first concatenated table if the column names are already presen t
in the original table. F or example, w e ha v e the original RA column, and another one called
RA-1 ). Similarly a -2 has b een added for the columns of the second concatenated table.
Ho w ev er, this example clearly sho ws a problem with this full concatenation: some
columns are iden tical (for example, HOST_OBJ_ID and HOST_OBJ_ID-1 ), or not needed
(for example, RA-1 and DEC-1 whic h are not necessary here). In suc h cases, y ou can use
--catcolumns to only concatenate certain columns, not the whole table. F or example, this
command:
$ asttable cat/xdf-f160w.fits -hCLUMPS --output=two-in-one-2.fits \
--catcolumnfile=cat/xdf-f125w-on-f160w-lab.fits \
--catcolumnhdu=CLUMPS --catcolumns=MAGNITUDE
$ asttable two-in-one-2.fits -i
Chapter 2: T utorials 56
Y ou see that w e ha v e no w only app ended the MAGNITUDE column of cat/xdf-f125w-
on-f160w-lab.fits . This is what w e needed to b e able to later subtract the magnitudes.
Let’s go ahead and add the F105W magnitudes also with the command b elo w. Note how
w e need to call --catcolumnhdu once for ev ery table that should b e app ended, but w e only
call --catcolumn once (assuming all the tables that should b e app ended ha v e this column).
$ asttable cat/xdf-f160w.fits -hCLUMPS --output=three-in-one-2.fits \
--catcolumnfile=cat/xdf-f125w-on-f160w-lab.fits \
--catcolumnfile=cat/xdf-f105w-on-f160w-lab.fits \
--catcolumnhdu=CLUMPS --catcolumnhdu=CLUMPS \
--catcolumns=MAGNITUDE
$ asttable three-in-one-2.fits -i
But w e are not finished y et! There is a v ery big problem: it is not immediately clear
whic h one of MAGNITUDE , MAGNITUDE-1 or MAGNITUDE-2 columns b elong to whic h filter!
Righ t no w, y ou kno w this b ecause y ou just ran this command. But in one hour, y ou’ll start
doubting y ourself and will b e forced to go through y our command history , trying to figure
out if y ou added F105W first, or F125W. Y ou should never torture y our future-self (or y our
colleagues) lik e this! So, let’s rename these confusing columns in the matc hed catalog.
F ortunately , with the --colmetadata option, y ou can correct the column metadata of
the final table (just b efore it is written). It tak es four v alues: 1) the original column name
or n um b er, 2) the new column name, 3) the column unit and 4) the column comments.
Since the commen ts are usually h uman-friendly sen tences and con tain space c haracters,
y ou should put them in double quotations lik e b elow. F or example, by adding three calls
of this option to the previous command, w e write the filter name in the magnitude column
name and description.
$ asttable cat/xdf-f160w.fits -hCLUMPS --output=three-in-one-3.fits \
--catcolumnfile=cat/xdf-f125w-on-f160w-lab.fits \
--catcolumnfile=cat/xdf-f105w-on-f160w-lab.fits \
--catcolumnhdu=CLUMPS --catcolumnhdu=CLUMPS \
--catcolumns=MAGNITUDE \
--colmetadata=MAGNITUDE,MAG-F160W,log,"Magnitude in F160W." \
--colmetadata=MAGNITUDE-1,MAG-F125W,log,"Magnitude in F125W." \
--colmetadata=MAGNITUDE-2,MAG-F105W,log,"Magnitude in F105W."
$ asttable three-in-one-3.fits -i
W e no w hav e all three magnitudes in one table and can start doing arithmetic on them
(to estimate colors, which are just a subtraction of magnitudes). T o use column arith-
metic, simply call the column selection option ( --column or -c ), put the v alue in single
quotations and start the v alue with arith (follo w ed b y a space) lik e the example b elow.
Column arithmetic uses the same “rev erse p olish notation” as the Arithmetic program (see
Section 6.2.1 [Rev erse p olish notation], page 398), with almost all the same op erators (see
Section 6.2.4 [Arithmetic op erators], page 407), and some column-sp ecific op erators (that
are not a v ailable for images). In column-arithmetic, y ou can iden tify columns b y n umber
(prefixed with a $ ) or name, for more see Section 5.3.3 [Column arithmetic], page 345.
So let’s estimate one color from three-in-one-3.fits using column arithmetic. All the
commands b elo w will pro duce the same output, try them eac h and fo cus on the differences.
Note that column arithmetic can b e mixed with other w a ys to c ho ose output columns (the
-c option).
Chapter 2: T utorials 57
$ asttable three-in-one-3.fits -ocolor-cat.fits \
-c1,2,3,4,'arith $5 $7 -'
$ asttable three-in-one-3.fits -ocolor-cat.fits \
-c1,2,RA,DEC,'arith MAG-F125W MAG-F160W -'
$ asttable three-in-one-3.fits -ocolor-cat.fits -c1,2 \
-cRA,DEC --column='arith MAG-F105W MAG-F160W -'
This example again highligh ts the imp ortant p oin t on using column names: if y ou do
not kno w the commands b efore, y ou ha v e no w a y of making sense of the first command:
what is in column 5 and 7? wh y not subtract columns 3 and 4 from each other? Do
y ou see ho w cryptic the first one is? Then lo ok at the last one: ev en if y ou ha v e no idea
ho w this table w as created, y ou immediately understand the desired op eration. When y ou
ha v e column names, please use them. If y our table do es not ha v e column names, giv e them
names with the --colmetadata (describ ed ab o ve) as y ou are creating them. But ho w ab out
the metadata for the column y ou just created with column arithmetic? Ha v e a lo ok at the
column metadata of the table pro duced ab o v e:
$ asttable color-cat.fits -i
The name of the column pro duced b y arithmetic column is ARITH_1 ! This is natural:
Arithmetic has no idea what the mo dified column is! Y ou could ha v e m ultiplied t w o columns,
or done m uc h more complex transformations with man y columns. Metadata c annot b e
set automatic al ly, your (the human) input is ne c essary. T o add metadata, y ou can use
--colmetadata lik e b efore:
$ asttable three-in-one-3.fits -ocolor-cat.fits -c1,2,RA,DEC \
--column='arith MAG-F105W MAG-F160W -' \
--colmetadata=ARITH_1,F105W-F160W,log,"Magnitude difference"
$ asttable color-cat.fits -i
Sometimes, b ecause of a particular w a y of storing data, y ou migh t need to tak e all input
columns. If there are man y columns (for example hundreds!), listing them (lik e ab o v e) will
b ecome anno ying, buggy and time-consuming. In suc h cases, y ou can giv e -c_all . Up on
execution, _all will b e replaced with a comma-separated list of all the input columns. This
allo ws y ou to add new columns easily , without ha ving to w orry ab out the n umber of input
columns that y ou w an t an yw a y . A lo w er-lev el but more customizable metho d is to use the
seq (sequence) command with the -s (separator) option set to ',' ). F or example, if you
ha v e 216 columns and only w an t to return columns 1 and 2 as w ell as all the columns
b et w een 12 to 58 (inclusiv e), y ou can use the command b elo w:
$ asttable table.fits -c1,2,$(seq -s',' 12 58)
W e are no w ready to make our final table. W e wan t it to ha v e the magnitudes in all
three filters, as well as the three possible colors. Recall that b y con v en tion in astronom y
colors are defined b y subtracting the bluer magnitude from the redder magnitude. In this
w a y a larger color v alue corresp onds to a redder ob ject. So from the three magnitudes,
w e can pro duce three colors (as sho wn b elo w). Also, b ecause this is the final table w e are
creating here and w ant to use it later, w e will store it in cat/ and w e will also giv e it a
clear name and use the --range option to only prin t columns with a signal-to-noise ratio
( SN column, from the F160W filter) ab o v e 5.
Chapter 2: T utorials 64
$ asttable cat/mags-with-color.fits --range=F105W-F160W,1.5,inf | wc -l
Let’s crop the F160W image around eac h of these ob jects, but w e first need a unique
iden tifier for them. W e will define this iden tifier using the ob ject and clump lab els (with an
underscore b et w een them) and feed the output of the command ab o ve t o A WK to generate
a catalog. Note that since w e are making a plain text table, w e will define the necessary
(for the string-t yp e first column) metadata man ually (see Section 4.7.2 [Gn uastro text table
format], page 284).
$ echo "# Column 1: ID [name, str10] Object ID" > cat/reddest.txt
$ asttable cat/mags-with-color.fits --range=F105W-F160W,1.5,inf \
| awk '{printf("%d_%-10d %f %f\n", $1, $2, $3, $4)}' \
>> cat/reddest.txt
Let’s see ho w these ob jects are p ositioned o ver the dataset. DS9 has the “Region”s
concept for this purp ose. And y ou build suc h regions easily from a table using Gn uastro’s
astscript-ds9-region installed script, using the command b elo w:
$ astscript-ds9-region cat/reddest.txt -c2,3 --mode=wcs \
--command="ds9 flat-ir/xdf-f160w.fits -zscale"
W e can no w feed cat/reddest.txt in to Gnuastro’s Crop program to get separate
p ostage stamps for eac h ob ject. T o k eep things clean, we will mak e a directory called
crop-red and ask Crop to sa v e the crops in this directory . W e will also add a -f160w.fits
suffix to the crops (to remind us whic h filter they came from). The width of the crops will
b e 15 arc-seconds (or 15/3600 degrees, whic h is the units of the W CS).
$ mkdir crop-red
$ astcrop flat-ir/xdf-f160w.fits --mode=wcs --namecol=ID \
--catalog=cat/reddest.txt --width=15/3600,15/3600 \
--suffix=-f160w.fits --output=crop-red
Lik e the Mak eProfiles command in Section 2.1.17 [Ap erture photometry], page 61, if
y ou lo ok at the order of the crops, y ou will notice that the crops are not made in order!
This is b ecause eac h crop is indep enden t of the rest, therefore crops are done in parallel,
and parallel op erations are async hronous. So the order can differ in each run, but the final
output is the same! In the command ab o v e, you can c hange f160w to f105w to mak e the
crops in b oth filters. Y ou can see all the cropp ed FITS files in the crop-red directory with
this command:
$ astscript-fits-view crop-red/*.fits
T o view the crops more easily (not ha ving to op en ds9 for eac h image), you can con v ert
the FITS crops in to the JPEG format with a shell lo op lik e b elo w.
$ cd crop-red
$ for f in *.fits; do \
astconvertt $f --fluxlow=-0.001 --fluxhigh=0.005 --invert -ojpg; \
done
$ cd ..
$ ls crop-red/
Y ou can no w use y our general graphic user in terface image view er to flip through the
images more easily , or imp ort them in to y our pap ers/rep orts.
The for lo op ab o v e to con v ert the images will do the job in series: each file is con-
v erted only after the previous one is complete. But lik e the crops, eac h JPEG image is
Chapter 2: T utorials 65
indep enden t, so let’s parallelize it. In other w ords, w e w an t to run more than one instance
of the command at an y momen t. T o do that, w e will use Mak e ( https://en.wikipedia.
org/wiki/Make_(software) ). Mak e is a v ery w onderful pip eline managemen t system, and
the most common and p o w erful implemen tation is GNU Mak e ( https://www.gnu.org/
software/make ), which has a complete man ual just lik e this one. W e cannot go in to the
details of Mak e here, for a hands-on video tutorial, se e this video in tro duction ( https://
peertube.stream/w/iJitjS3r232Z8UPMxKo6jq ). T o do the pro cess ab o v e in Make, please
cop y the con ten ts b elo w in to a plain-text file called Makefile . Just replace the __[TAB]__
part at the start of the line with a single ‘ TAB ’ button on y our k eyb oard.
jpgs=$(subst .fits,.jpg,$(wildcard *.fits))
all: $(jpgs)
$(jpgs): %.jpg: %.fits
__[TAB]__astconvertt $< --fluxlow=-0.001 --fluxhigh=0.005 \
__[TAB]__ --invert -o$
No w that the Makefile is ready , y ou can run Mak e on 12 threads using the commands
b elo w. F eel free to replace the 12 with an y n um b er of threads y ou ha v e on y our system
(y ou can find out b y running the nproc command on GNU/Lin ux op erating systems):
$ make -j12
Did y ou notice ho w m uc h faster this one w as? When p ossible, it is alw a ys v ery helpful to do
y our analysis in parallel. Y ou can build v ery complex w orkflo ws with Mak e, for example,
see Akhlaghi 2021 ( https://arxiv.org/abs/2006.03018 ) so it is w orth sp ending some
time to master.
2.1.20 FITS images in a publication
In the previous section (Section 2.1.19 [Reddest clumps, cutouts and parallelization],
page 63), w e visually insp ected the p ositions of the reddest ob jects using DS9. That is
v ery go o d for an in teractiv e insp ection of the ob jects: y ou can zo om-in and out, y ou can
do measuremen ts, etc. Once the exp erimen tation phase of y our pro ject is complete, y ou
w an t to sho w these ob jects o v er the whole image in a rep ort, pap er or slides.
One solution is to use DS9 itself ! F or example, run the astscript-fits-view command
of the previous section to op en DS9 with the regions o v er-plotted. Clic k on the “File” men u
and select “Sa v e Image”. In the side-men u that op ens, y ou ha v e m ultiple formats to select
from. Usually for publications, we w an t to sho w the regions and text (in the colorbar) in
v ector graphics, so it is b est to exp ort to EPS. Once y ou ha v e made the EPS, y ou can then
con v ert it to PDF with the epspdf command.
Another solution is to use Gn uastro’s Conv ertT yp e program. The main difference is that
DS9 is a Graphic User In terface (GUI) program, so it tak es relatively long (about a second)
to load, and it requires man y dep endencies. This will slo w-do wn automatic con v ersion
of man y files, and will mak e y our co de hard to mo v e to another op erating system. DS9
do es ha v e a command-line in terface that y ou can use to automate the creation of each file,
ho w ev er, it has a v ery p eculiar command-line in terface and formats (lik e the “region” files).
Ho w ev er, in Con v ertT yp e, there is no graphic interface, so it has v ery few dep endencies, it is
fast, and finally , it tak es normal tables (in plain-text or FITS) as input. So in this concluding
step of the analysis, let’s build a nice publication-ready plot, sho wing the p ositions of the
reddest ob jects in the image for our pap er.
Chapter 2: T utorials 66
In Section 2.1.19 [Reddest clumps, cutouts and parallelization], page 63, we already used
Con v ertT yp e to mak e JPEG p ostage stamps. Here, w e will use it to mak e a PDF image of
the whole deep region. T o start, let’s simply run Con v ertT yp e on the F160W image:
$ astconvertt flat-ir/xdf-f160w.fits -oxdf.pdf
Op en the output in a PDF view er. Y ou see that it is almost fully blac k! Let’s see wh y
this happ ens! First, with the t w o commands b elow, let’s calculate the maxim um v alue, and
the standard deviation of the sky in this image (using NoiseChisel’s output, which w e found
at the end of Section 2.1.11 [NoiseChisel optimization for detection], page 41). Note that
NoiseChisel writes the median sky standard deviation b efor e in terp olation in the MEDSTD
k eyw ord of the SKY_STD HDU. This is more robust than the median of the Sky standard
deviation image (whic h has gone through in terp olation).
$ max=$(aststatistics nc/xdf-f160w.fits -hINPUT-NO-SKY --maximum)
$ skystd=$(astfits nc/xdf-f160w.fits -hSKY_STD --keyvalue=MEDSTD -q)
$ echo $max $skystd
58.8292 0.000410282
$ echo $max $skystd | awk '{print $1/$2}'
143387
In the last command ab o ve, w e divided the maxim um b y the sky standard deviation. Y ou
see that the maxim um v alue is more than 140000 times larger than the noise lev el! On the
other hand common monitors or prin ters, usually hav e a maxim um dynamic range of 8-bits,
only allo wing for 2 8 = 256 la y ers. This is therefore the maxim um n um b er of “la y ers” y ou
can ha v e in a common displa y formats lik e JPEG, PDF or PNG! Dividing the result ab o v e
b y 256, w e get a la y er spacing of
$ echo $max $skystd | awk '{print $1/$2/256}'
560.106
In other w ords, the first la y er (whic h is blac k) will con tain all the pixel v alues b elo w
∼ 560! So all pixels with a signal-to-noise ratio lo w er than ∼ 560 will ha v e a blac k color
since they fall in the first la yer of an 8-bit PDF (or JPEG) image. This happ ens b ecause
b y default w e are assuming a linear mapping from floating p oin t to 8-bit in tegers.
T o fix this, w e should mo v e to a differen t mapping. A go o d, physically motiv ated,
mapping is Surface Brigh tness (whic h is in log-scale, see Section 7.4.2 [Brigh tness, Flux,
Magnitude and Surface brigh tness], page 576). F ortunately this is very easy to do with
Gn uastro’s Arithmetic program, as shown in the commands below (using the kno wn zero
p oin t 23 , and after calculating the pixel area in units of arcsec 2 ):
$ zeropoint=25.94
$ pixarcsec2=$(astfits nc/xdf-f160w.fits --pixelareaarcsec2)
$ astarithmetic nc/xdf-f160w.fits $zeropoint $pixarcsec2 counts-to-sb \
--output=xdf-f160w-sb.fits
With the t wo commands below, first, let’s lo ok at the dynamic range of the image no w
(dividing the maxim um b y the minim um), and then let’s op en the image and ha v e a lo ok
at it:
23 https://archive.stsci.edu/prepds/xdf/#science-images
Chapter 2: T utorials 67
$ aststatistics xdf-f160w-sb.fits --minimum --maximum
$ astscript-fits-view xdf-f160w-sb.fits
The go o d news is that the dynamic range has no w decreased to ab out 2! In other w ords,
w e can distribute the 256 la y ers of an 8-bit displa y o v er a m uc h smaller range of v alues, and
therefore b etter visualize the data. How ev er, there are t w o imp ortan t p oin ts to consider
from the output of the first command and a visual insp ection of the second.
• The largest pixel v alue (fain test surface brigh tness lev el) in the image is ∼ 43! This
is far to o lo w to b e realistic, and is just due to noise. As discussed in Section 2.1.14
[Measuring the dataset limits], page 49, the 3 σ surface brigh tness limit of this image,
o v er 100 arcsec 2 is roughly 32.66 mag/arcsec 2 .
• Y ou see man y NaN pixels in b et ween the galaxies! These are due to the fact that the
magnitude is defined on a logarithmic scale and the logarithm of a negativ e n um b er is
not defined.
In other w ords, w e should replace all NaN pixels, and pixels with a surface brigh tness
v alue fain ter than the image surface brigh tness limit to this limit. With the first command
b elo w, w e will first extract the surface brigh tness limit from the catalog headers that w e
calculated b efore, and then call Arithmetic to use this limit.
$ sblimit=$(astfits cat/xdf-f160w.fits --keyvalue=SBLMAG -q)
$ astarithmetic nc/xdf-f160w.fits $zeropoint $pixarcsec2 \
counts-to-sb set-sb \
sb sb $sblimit gt sb isblank or $sblimit where \
--output=xdf-f160w-sb.fits
Let’s con v ert this image in to a PDF with the command b elo w:
$ astconvertt xdf-f160w-sb.fits --output=xdf-f160w-sb.pdf
It is m uc h b etter no w and we can visualize man y features of the FITS file (from the
cen tral structures of the galaxies and stars, to a little in to the noise and their lo w surface
brigh tness features. Ho wev er, the image generally lo oks a little to o gra y! This is b ecause
of that brigh t star in the b ottom half of the image! Stars are v ery sharp! So let’s man ually
tell Con v ertT yp e to set an y pixel with a v alue less than (brigh ter than) 20 to blac k (and
not use the minim um). W e do this with the --fluxlow option:
$ astconvertt xdf-f160w-sb.fits --output=xdf-f160w-sb.pdf --fluxlow=20
W e are still missing some of the diffuse flux in this PDF. This is b ecause of those negativ e
pixels that w ere set to NaN. T o b etter sho w these structures, w e should w arp the image to
larger pixels. So let’s w arp it to a pixel grid where the new pixels are 4 × 4 larger than the
original pixels. But b e careful that w arping should b e done on the original image, not on
the surface brigh tness image. W e should re-calculate the surface brigh tness image after the
w arping is one. This is b ecause l og ( a + b ) 6 = l og ( a ) + l og ( b ). Recall that surface brigh tness
calculation in v olv es a logarithm, and w arping in v olv es addition of pixel v alues.
$ astwarp nc/xdf-f160w.fits --scale=1/4 --centeroncorner \
--output=xdf-f160w-warped.fits
$ pixarcsec2=$(astfits xdf-f160w-warped.fits --pixelareaarcsec2)
$ astarithmetic xdf-f160w-warped.fits $zeropoint $pixarcsec2 \
Chapter 2: T utorials 68
counts-to-sb set-sb \
sb sb $sblimit gt sb isblank or $sblimit where \
--output=xdf-f160w-sb.fits
$ astconvertt xdf-f160w-sb.fits --output=xdf-f160w-sb.pdf --fluxlow=20
Ab o v e, w e needed to re-calculate the pixel area of the w arpp ed image, but we did not
need to re-calculate the surface brigh tness limit! The reason is that the surface brigh tness
limit is indep enden t of the pixel area (in its deriv ation, the pixel area has b een accounted
for). As a side-effect of the warping, the n umber of pixels in the image also dramatically
decreased, therefore the v olume of the output PDF (in b ytes) is also smaller, making y our
pap er/rep ort easier to upload/do wnload or send by email. This visual resolution is still
more than enough for including on top of a column in y our pap er!
I do not ha v e the zero p oin t of m y image: The absolute v alue of the zero p oin t is irrelev an t
for the finally pro duced PDF. W e used it here b ecause it w as a v ailable and mak es the
n um b ers ph ysically understandable. If y ou do not hav e the zero p oin t, just set it to zero
(whic h is also the default zero p oin t used by Mak eCatalog when it estimates the surface
brigh tness limit). F or the v alue to --fluxlow ab o v e, y ou can simply subtract ∼ 10 from
the surface brigh tness limit.
T o summarize, and to k eep the image for the next section in a separate directory , here are
the necessary commands:
$ zeropoint=25.94
$ mkdir report-image
$ cd report-image
$ sblimit=$(astfits cat/xdf-f160w.fits --keyvalue=SBLMAG -q)
$ astwarp nc/xdf-f160w.fits --scale=1/4 --centeroncorner \
--output=warped.fits
$ pixarcsec2=$(astfits warped.fits --pixelareaarcsec2)
$ astarithmetic warped.fits $zeropoint $pixarcsec2 \
counts-to-sb set-sb \
sb sb $sblimit gt sb isblank or $sblimit where \
--output=sb.fits
$ astconvertt sb.fits --output=sb.pdf --fluxlow=20
Finally , let’s remo v e all the temp orary files w e built in the top-lev el tutorial directory:
$ rm *.fits *.pdf
Color images: In this tutorial w e just used one of the filters and sho w ed the surface
brigh tness image of that single filter as a gra yscale image. But the image can also b e
in color (using three filters) to b etter con v ey the ph ysical prop erties of the ob jects in
y our image. T o create an image that sho ws the full dynamic range of y our data, see this
dedicated tutorial Section 2.6 [Color images with full dynamic range], page 150.
Chapter 2: T utorials 69
2.1.21 Marking ob jects for publication
In Section 2.1.20 [FITS images in a publication], page 65, w e created a ready-to-prin t
visualization of the FITS image used in this tutorial. Ho w ev er, y ou rarely wan t to sho w
a nak ed image lik e that! Y ou usually w an t to highligh t some ob jects (that are the target
of y our science) o ver the image and sho w differen t marks for the v arious t yp es of ob jects
y ou are studying. In this tutorial, w e will do just that: select a sub-set of the full catalog
of clumps, and sho w them with differen t marks shap es and colors, while also adding some
text under eac h mark. T o add co ordinates on the edges of the figure in y our pap er, see
Section 5.2.4 [Annotations for figure in pap er], page 319.
T o start with, let’s put a red plus sign o v er the sub-sample of reddest clumps similar to
Section 2.1.19 [Reddest clumps, cutouts and parallelization], page 63. First, w e will need to
mak e the table of marks. W e will c ho ose those with a color stronger than 1.5 magnitudes
and a signal-to-noise ratio (in F160W) larger than 5. W e also only need the RA, Dec, color
and magnitude (in F160W) columns (recall that at the end of the previous section w e w ere
already in the report-image/ directory):
$ asttable cat/mags-with-color.fits --range=F105W-F160W,1.5:inf \
--range=sn-f160w,5:inf -cRA,DEC,MAG-F160w,F105W-F160W \
--output=reddest-cat.fits
Gn uastro’s Con v ertT yp e program also has features to add marks o v er the finally pro-
duced PDF. Belo w, w e will start with the same astconvertt command of the previous
section. The p ositions of the marks should b e giv en as a table to the --marks option.
Tw o other options are also mandatory: --markcoords iden tifies the columns that con tain
the co ordinates of eac h mark and --mode sp ecifies if the co ordinates are in image or W CS
co ordinates.
$ astconvertt sb.fits --output=reddest.pdf --fluxlow=20 \
--marks=reddest-cat.fits --mode=wcs \
--markcoords=RA,DEC
Op en the output reddest.pdf and see the result. Y ou will see relativ ely thic k red circles
placed o v er the giv en co ordinates. In y our PDF bro wser, zo om-in to one of the regions,
y ou will see that while the pixels of the bac kground image b ecome larger, the lines of these
regions do not degrade! This is the concept/p o w er of V ector Graphics: ideal for publication!
F or more on raster (pixelated) and v ector (infinite-resolution) graphics, see Section 5.2.1
[Raster and V ector graphics], page 313.
W e had planned to put a plus-sign on eac h ob ject. Ho w ev er, b ecause we did not explicitly
ask for a certain shap e, Con v ertT yp e put a circle. Each mark can ha v e its o wn separate
shap e. Shap es can b e giv en b y a name or a co de. The full list of av ailable shap es names
and co des is giv en in the description of --markshape option of Section 5.2.5.3 [Dra wing
with v ector graphics], page 333.
T o use a differen t shap e, w e need to add a new column to the base table, con taining
the iden tifier of the desired shap e for eac h mark. F or example, the co de for the plus sign
is 2 . With the commands b elo w, w e will add a new column with this fixed v alue. With
the first A WK command w e will mak e a single-column file, where all the ro ws hav e the
same v alue. W e pip e our base table in to A WK, so it has the same num b er of ro ws. With
the second command, w e concatenate (or app end) the new column with T able, and giv e
this new column the name SHAPE (to easily refer to it later and not ha v e to coun t). With
Chapter 2: T utorials 70
the third command, w e clean-up b ehind our selv es (deleting the extra params.txt file).
Finally , w e use the --markshape option to tell Conv ertT yp e whic h column to use for the
shap e iden tifier.
$ asttable reddest-cat.fits | awk '{print 2}' > params.txt
$ asttable reddest-cat.fits --catcolumnfile=params.txt \
--colmetadata=5,SHAPE,id,"Shape of mark" \
--output=reddest-marks.fits
$ rm params.txt
$ astconvertt sb.fits --output=reddest.pdf --fluxlow=20 \
--marks=reddest-marks.fits --mode=wcs \
--markcoords=RA,DEC --markshape=SHAPE
Op en the PDF and ha v e a lo ok! Y ou do see red signs o v er the co ordinates, but the thic k
plus-signs only b ecome visible after y ou zo om-in m ultiple times! T o mak e them larger, y ou
can giv e another column to sp ecify the size of eac h mark. Let’s set the full width of the
plus sign to extend 3 arcseconds. The commands are similar to ab o v e, try to follo w the
difference (in particular, ho w w e use --sizeinarcsec ).
$ asttable reddest-cat.fits | awk '{print 2, 3}' > params.txt
$ asttable reddest-cat.fits --catcolumnfile=params.txt \
--colmetadata=5,SHAPE,id,"Shape of mark" \
--colmetadata=6,SIZE,arcsec,"Size in arcseconds" \
--output=reddest-marks.fits
$ rm params.txt
$ astconvertt sb.fits --output=reddest.pdf --fluxlow=20 \
--marks=reddest-marks.fits --mode=wcs \
--markcoords=RA,DEC --markshape=SHAPE \
--marksize=SIZE --sizeinarcsec
The p o w er of this metho dology is that eac h mark can b e completely differen t! F or
example, let’s sho w the ob jects with a color less than 2 magnitudes with a circle, and those
with a stronger color with a plus (recall that the co de for a circle w as 1 and that of a plus
w as 2 ). Y ou only need to replace the first command ab ov e with the one b elo w. Afterw ards,
run the rest of the commands in the last co de-blo c k.
$ asttable reddest-cat.fits -cF105W-F160W \
| awk '{if($1<2) shape=1; else shape=2; print shape, 3}' \
> params.txt
Ha v e a lo ok at the resulting reddest.pdf . Y ou see that the circles are muc h larger than
the plus signs. This is b ecause the “size” of a cross is defined to b e its full width, but for a
circle, the v alue in the size column is the radius. The wa y eac h shap e in terprets the v alue
of the size column is fully describ ed under --markshape of Section 5.2.5.3 [Dra wing with
v ector graphics], page 333. T o mak e them more comparable, let’s set the circle sizes to b e
half of the cross sizes.
$ asttable reddest-cat.fits -cF105W-F160W \
Chapter 2: T utorials 71
| awk '{if($1<2) {shape=1; size=1.5} \
else {shape=2; size=3} \
print shape, size}' \
> params.txt
Let’s mak e things a little more complex (and sho w more information in the visualization)
b y using color. Gn uastro recognizes the full extended w eb colors ( https://en.wikipedia.
org/wiki/Web_colors#Extended_colors ), for their full list (con taining names and co des)
see Section 5.2.3.3 [V ector graphics colors], page 319. But lik e ev erything else, an ev en easier
w a y to view and select the color for y our figure is on the command-line! If y our terminal
supp orts 24-bit true-color, y ou can see all the colors b y running this command (supp orted
on mo dern GNU/Lin ux distributions):
$ astconvertt --listcolors
w e will giv e a “Sienna” color for the ob jects that are fain ter than 29th magnitude and a
“deeppink” color to the brigh ter ones (while k eeping the same shap es definition as b efore)
Since there are man y colors, using their co des can mak e the table hard to read b y a h uman!
So let’s use the color names instead of the color co des in the example b elo w (this is useful
in other columns require strings-only , lik e the fon t name).
The only in tricacy is in the making of params.txt . Recall that string columns need
column metadata (Section 4.7.2 [Gn uastro text table format], page 284). In this particular
case, since the string column is the last one, w e can safely use A WK’s print command.
But if y ou ha v e m ultiple string columns, to b e safe it is b etter to use A WK’s printf and
explicitly sp ecify the n um b er of c haracters in the string columns.
$ asttable reddest-cat.fits -cF105W-F160W,MAG-F160W \
| awk 'BEGIN{print "# Column 3: COLOR [name, str8]"}\
{if($1<2) {shape=1; size=1.5} \
else {shape=2; size=3} \
if($2>29) {color="sienna"} \
else {color="deeppink"} \
print shape, size, color}' \
> params.txt
$ asttable reddest-cat.fits --catcolumnfile=params.txt \
--colmetadata=5,SHAPE,id,"Shape of mark" \
--colmetadata=6,SIZE,arcsec,"Size in arcseconds" \
--output=reddest-marks.fits
$ rm params.txt
$ astconvertt sb.fits --output=reddest.pdf --fluxlow=20 \
--marks=reddest-marks.fits --mode=wcs \
--markcoords=RA,DEC --markshape=SHAPE \
--marksize=SIZE --sizeinarcsec --markcolor=COLOR
As one final example, let’s write the magnitude of each ob ject under it. Since the
magnitude is already in the marks.fits that w e pro duced ab o v e, it is v ery easy to add it
(just add --marktext option to Con v ertT yp e):
$ astconvertt sb.fits --output=reddest.pdf --fluxlow=20 \
Chapter 2: T utorials 72
--marks=reddest-marks.fits --mode=wcs \
--markcoords=RA,DEC --markshape=SHAPE \
--marksize=SIZE --sizeinarcsec \
--markcolor=COLOR --marktext=MAG-F160W
Op en the final PDF ( reddest.pdf ) and y ou will see the magnitudes written under each
mark in the same color. In the case of magnitudes (where the magnitude error is usually
m uc h larger than 0.01 magnitudes, four decimals is not to o meaningful. By default, for
prin ting floating p oin t columns, w e use the compiler’s default precision (whic h is ab out 4
digits for 32-bit floating p oin t n um b ers). But y ou can o v er-write this (to only sho w t w o
digits after the decimal p oin t) with the --marktextprecision=2 option.
Y ou can customize the written text b y sp ecifying a differen t line-width (for the text,
differen t from the main mark), or ev en sp ecifying a differen t fon t for each mark! Y ou can
see the full list of a v ailable fon ts for the text under a mark with the first command b elo w
and with the second, y ou can actually see them in a custom PDF (to only sho w the fon ts).
$ astconvertt --listfonts
$ astconvertt --showfonts
As y ou see, there are man y w a ys y ou can customize each mark! The ab o v e examples
w ere just the tip of the iceburg! But this section has already b ecome long so w e will stop
it here (see the b o x at the end of this section for y et another useful example). Lik e ab o v e,
eac h feature of a mark can b e con trolled with a column in the table of mark information.
Please see in Section 5.2.5.3 [Dra wing with v ector graphics], page 333, for the full list of
columns/features that y ou can use.
Dra wing ellipses: With the commands b elo w, you can measure the elliptical prop erties of
the ob jects and visualized them in a ready-to-publish PDF (w e will only sho w the ellipses
of the largest clumps):
$ astmkcatalog ../seg/xdf-f160w.fits --ra --dec --semi-major \
--axis-ratio --position-angle --clumpscat \
--output=ellipseinfo.fits
$ asttable ellipseinfo.fits -hCLUMPS | awk '{print 4}' > params.txt
$ asttable ellipseinfo.fits -hCLUMPS --catcolumnfile=params.txt \
--range=SEMI_MAJOR,10,inf -oellipse-marks.fits \
--colmetadata=6,SHAPE,id,"Shape of mark"
$ astconvertt sb.fits --output=ellipse.pdf --fluxlow=20 \
--marks=ellipse-marks.fits --mode=wcs \
--markcoords=RA,DEC --markshape=SHAPE \
--marksize=SEMI_MAJOR,AXIS_RATIO --sizeinpix \
--markrotate=POSITION_ANGLE
T o conclude this section, let us highlight an imp ortan t factor to consider in v ector graph-
ics. In Con v ertT yp e, things lik e line width or fon t size are defined in units of p oints . In
v ector graphics standards, 72 p oin ts corresp ond to one inc h. Therefore, one w a y y ou can
c hange these factors for all the ob jects is to assign a larger or smaller prin t size to the im-
age. The prin t size is just a meta-data en try , and will not affect the file’s v olume in b ytes!
Chapter 2: T utorials 73
Y ou can do this with the --widthincm option. T ry adding this option and giving it v ery
differen t v alues lik e 5 or 30 .
2.1.22 W riting scripts to automate the steps
In the previous sub-sections, w e w en t through a series of steps lik e do wnloading the necessary
datasets (in Section 2.1.3 [Setup and data do wnload], page 25), detecting the ob jects in the
image, and finally selecting a particular subset of them to insp ect visually (in Section 2.1.19
[Reddest clumps, cutouts and parallelization], page 63). T o b enefit most effectiv ely from
this subsection, please go through the previous sub-sections, and if y ou ha v e not actually
done them, w e recommended to do/run them b efore con tin uing here.
Eac h sub-section/step of the sub-sections ab o v e in v olv ed sev eral commands on the
command-line. Therefore, if y ou w an t to repro duce the previous results (for example, to
only c hange one part, and see its effect), you’ll ha v e to go through all the sections ab o ve
and read through them again. If y ou hav e ran the commands recen tly , y ou ma y also ha v e
them in the history of y our shell (command-line environ ment). Y ou can see man y of y our
previous commands on the shell (ev en if y ou ha v e closed the terminal) with the history
command, lik e this:
$ history
T ry it in y our teminal to see for y ourself. By default in GNU Bash, it shows the last
500 commands. Y ou can also sa v e this “history” of previous commands to a file using shell
redirection (to ha v e it after y our next 500 commands), with this command
$ history > my-previous-commands.txt
This is a go o d w ay to temporarily k eep trac k of ev ery single command y ou ran. But
in the middle of all the useful commands, y ou will ha v e man y extra commands, lik e tests
that y ou did b efore/after the go o d output of a step (that y ou decided to con tin ue w orking
on), or an unrelated job y ou had to do in the middle of this pro ject. Because of these
impurities, after a few days (that y ou ha v e forgot the con text: tests y ou did not end-up
using, or unrelated jobs) reading this full history will b e v ery frustrating.
Keeping the final commands that w ere used in eac h step of an analysis is a common
problem for an y one who is doing something serious with the computer. But simply k eeping
the most imp ortan t commands in a text file is not enough, the small steps in the middle
(lik e making a directory to k eep the outputs of one step) are also imp ortan t. In other w ords,
the only w ay y ou can b e sure that y ou are under control of your processing (and actually
understand ho w y ou pro duced y our final result) is to run the commands automatically .
F ortunately , t yping commands in teractiv ely with y our fingers is not the only w a y to
op erate the shell. The shell can also tak e its orders/commands from a plain-text file, whic h
is called a script . When giv en a script, the shell will read it line-b y-line as if you ha v e
actually t yp ed it man ually .
Let’s con tin ue with an example: try t yping the commands b elo w in y our shell. With these
commands w e are making a text file ( a.txt ) con taining a simple 3 × 3 matrix, con v erting
it to a FITS image and computing its basic statistics. After the first three commands op en
a.txt with a text editor to actually see the v alues w e wrote in it, and after the fourth, op en
the FITS file to see the matrix as an image. a.txt is created through the shell’s redirection
feature: ‘ > ’ o verwrites the existing con ten ts of a file, and ‘ >> ’ app ends the new con ten ts
after the old con ten ts.
Chapter 2: T utorials 80
if ! [ -f $f105w_flat ]; then
astcrop --mode=wcs -h0 --output=$f105w_flat \
--polygon=$deep_polygon $dldir/$f105w_in
fi
if ! [ -f $f160w_flat ]; then
astcrop --mode=wcs -h0 --output=$f160w_flat \
--polygon=$deep_polygon $dldir/$f160w_in
fi
2.1.23 Citing and ac kno wledging Gn uastro
In conclusion, w e hop e this extended tutorial has b een a go o d starting p oin t to help in y our
exciting researc h. If this b o ok or an y of the programs in Gn uastro ha v e b een useful for
y our researc h, please cite the resp ectiv e pap ers, and ackno wledge the funding agencies that
made all of this p ossible. Without citations, w e will not b e able to secure future funding
to con tin ue w orking on Gn uastro or impro ving it, so please tak e soft w are citation seriously
(for all the scien tific soft w are y ou use, not just Gn uastro).
T o help y ou in this, all Gnuastro programs ha v e a --cite option to facilitate the citation
and ac kno wledgmen t. Just note that it ma y b e necessary to cite additional pap ers for
differen t programs, so please try it out on all the programs that y ou used, for example:
$ astmkcatalog --cite
$ astnoisechisel --cite
2.2 Detecting large extended targets
The outer wings of large and extended ob jects can sink in to the noise v ery gradually and can
ha v e a large v ariet y of shap es (for example, due to tidal in teractions). Therefore separating
the outer b oundaries of the galaxies from the noise can b e particularly tric ky . Besides
causing an under-estimation in the total estimated brigh tness of the target, failure to detect
suc h fain t wings will also cause a bias in the noise measuremen ts, thereb y hamp ering the
accuracy of an y measuremen t on the dataset. Therefore even if they do not constitute a
significan t fraction of the target’s ligh t, or are not your primary target, these regions m ust
not b e ignored. In this tutorial, w e will w alk y ou through the strategy of detecting such
targets using Section 7.2 [NoiseChisel], page 544.
Do not start with this tutorial: If y ou ha v e not already completed Section 2.1 [General
program usage tutorial], page 22, w e strongly recommend going through that tutorial
b efore starting this one. Bas ic features lik e access to this b o ok on the command-line,
the configuration files of Gn uastro’s programs, b enefiting from the mo dular nature of the
programs, viewing m ulti-extension FITS files, or using NoiseChisel’s outputs are discussed
in more detail there.
W e will try to detect the fain t tidal wings of the b eautiful M51 group 25 in this tutorial.
W e will use a dataset/image from the public Sloan Digital Sky Surv ey ( http://www.sdss.
org/ ), or SDSS. Due to its more p eculiar lo w surface brigh tness structure/features, we will
fo cus on the dw arf companion galaxy of the group (or NGC 5195).
25 https://en.wikipedia.org/wiki/M51_Group
Chapter 2: T utorials 81
2.2.1 Do wnloading and v alidating input data
T o get the image, you can use the simple field searc h ( https://dr12.sdss.org/fields )
to ol of SDSS. As long as it is co v ered b y the SDSS, y ou can find an image con taining y our
desired target either b y pro viding a standard name (if it has one), or its co ordinates. T o
access the dataset w e will use here, write NGC5195 in the “Ob ject Name” field and press
“Submit” button.
T yp e the example commands: T ry to t yp e the example commands on y our terminal and
use the history feature of y our command-line (by pressing the “up” button to retriev e
previous commands). Do not simply cop y and paste the commands shown here. This will
help sim ulate future situations when y ou are pro cessing y our o wn datasets.
Y ou can see the list of a v ailable filters under the color image. F or this demonstration, w e
will use the r-band filter image. By clic king on the “r-band FITS” link, y ou can do wnload
the image. Alternativ ely , y ou can just run the following command to do wnload it with
GNU Wget 26 . T o keep things clean, let’s also put it in a directory called ngc5195 . With
the -O option, w e are asking Wget to sav e the do wnloaded file with a more manageable
name: r.fits.bz2 (this is an r-band image of NGC 5195, whic h w as the directory name).
$ mkdir ngc5195
$ cd ngc5195
$ topurl=https://dr12.sdss.org/sas/dr12/boss/photoObj/frames
$ wget $topurl/301/3716/6/frame-r-003716-6-0117.fits.bz2 -Or.fits.bz2
When y ou w an t to repro duce a previous result (a kno wn analysis, on a kno wn dataset,
to get a kno wn result: lik e the case here!) it is imp ortan t to verify that the file is correct:
that the input file has not c hanged (on the remote serv er, or in your o wn arc hiv e), or there
w as no do wnloading problem. Otherwise, if the data ha v e c hanged in y our serv er/arc hiv e,
and y ou use the same script, y ou will get a differen t result, causing a lot of confusion!
One go o d w a y to v erify the con ten ts of a file is to store its Che cksum in y our analysis
script and c heck it before an y other op eration. The Che cksum algorithms lo ok in to the
con ten ts of a file and calculate a fixed-length string from them. If an y c hange (ev en in a
bit or b yte) is made within the file, the resulting string will c hange, for more see Wikip edia
( https://en.wikipedia.org/wiki/Checksum ). There are man y common algorithms, but
a simple one is the SHA-1 algorithm ( https://en.wikipedia.org/wiki/SHA-1 ) (Secure
Hash Algorithm 1) that y ou can calculate easily with the command b elo w (the second line
is the output, and the c hec ksum is the first/long string: it is indep enden t of the file name)
$ sha1sum r.fits.bz2
5fb06a572c6107c72cbc5eb8a9329f536c7e7f65 r.fits.bz2
If the c hec ksum on y our computer is differen t from this, either the file has b een incorrectly
do wnloaded (most probable), or it has c hanged on SDSS serv ers (v ery unlik ely 27 ). T o get
a b etter feeling of c hec ksums op en y our fa v orite text editor and mak e a test file b y writing
26 T o make the command easier to view on screen or in a page, w e ha ve defined the top URL of the image
as the topurl shell v ariable. Y ou can just replace the v alue of this v ariable with $topurl in the wget
command.
27 If y our c hecksum is differen t, try uncompressing the file with the bunzip2 command after this, and op en
the resulting FITS file. If it op ens and y ou see the image of M51 and NGC5195, then there w as no
Chapter 2: T utorials 82
something in it. Sa v e it and calculate the text file’s SHA-1 c hec ksum with sha1sum . T ry
renaming that file, and y ou’ll see the c hec ksum has not c hanged (c hec ksums only lo ok in to
the con ten ts, not the name/lo cation of the file). Then op en the file with y our text editor
again, mak e a c hange and re-calculate its c hec ksum, y ou’ll see the c hec ksum string has
c hanged.
Its alw a ys go o d to k eep this short c hec ksum string with y our pro ject’s scripts and v alidate
y our input data b efore using them. Y ou can do this with a shell conditional lik e this:
filename=r.fits.bz2
expected=5fb06a572c6107c72cbc5eb8a9329f536c7e7f65
sum=$(sha1sum $filename | awk '{print $1}')
if [ $sum = $expected ]; then
echo "$filename: validated"
else
echo "$filename: wrong checksum!"
exit 1
fi
No w that w e kno w y ou ha v e the same data that w e wrote this tutorial with, let’s con tin ue.
The SDSS serv er k eeps the files in a Bzip2 compressed file format (that ha ve a .bz2 suffix).
So w e will first decompress it with the follo wing command to use it as a normal FITS file. By
con v en tion, compression programs delete the original file (compressed when uncompressing,
or uncompressed when compressing). T o k eep the original file, y ou can use the --keep or
-k option whic h is a v ailable in most compression programs for this job. Here, w e do not
need the compressed file an y more, so w e will just let bunzip delete it for us and k eep the
directory clean.
$ bunzip2 r.fits.bz2
2.2.2 NoiseChisel optimization
In Section 2.2.1 [Do wnloading and v alidating input data], page 81, w e do wnloaded the single
exp osure SDSS image. Let’s see how NoiseChisel operates on it with its default parameters:
$ astnoisechisel r.fits -h0
As describ ed in Section 2.1.10 [NoiseChisel and Multi-Extension FITS files],
page 38, NoiseChisel’s default output is a m ulti-extension FITS file. Op en the output
r_detected.fits file and ha v e a lo ok at the extensions, the 0-th extension is only
meta-data and con tains NoiseChisel’s configuration parameters. The rest are the
Sky-subtracted input, the detection map, Sky v alues and Sky standard deviation.
$ ds9 -mecube r_detected.fits -zscale -zoom to fit
Flipping through the extensions in a FITS view er, y ou will see that the first image (Sky-
subtracted image) lo oks reasonable: there are no ma jor artifacts due to bad Sky subtraction
compared to the input. The second extension also seems reasonable with a large detection
map that co v ers the whole of NGC5195, but also extends to w ards the b ottom of the image
where w e actually see fain t and diffuse signal in the input image.
No w try flipping b et w een the DETECTIONS and SKY extensions. In the SKY extension,
y ou’ll notice that there is still significan t signal b ey ond the detected pixels. Y ou can tell
do wnload problem, and the file has indeed c hanged on the SDSS serv ers! In this case, please con tact us
at [email protected] .
Chapter 2: T utorials 83
that this signal b elongs to the galaxy b ecause the far-righ t side of the image (a w a y from
M51) is dark (has lo wer v alues) and the brighter parts in the Sky image (with larger v alues)
are just under the detections and follo w a similar pattern.
The fact that signal from the galaxy remains in the SKY HDU sho ws that NoiseChisel
can b e optimized for a m uch better result. The SKY extension m ust not con tain an y ligh t
around the galaxy . Generally , any time y our target is m uc h larger than the tile size and
the signal is v ery diffuse and extended at lo w signal-to-noise v alues (like this case), this wil l
happ en. Therefore, when there are large ob jects in the dataset, the b est place to c hec k the
accuracy of y our detection is the estimated Sky image.
When dominated b y the background, noise has a symmetric distribution. Ho w ev er, sig-
nal is not symmetric (w e do not ha v e negativ e signal). Therefore when non-constan t 28 signal
is presen t in a noisy dataset, the distribution will b e p ositiv ely sk ew ed. F or a demonstra-
tion, see Figure 1 of Akhlaghi and Ic hik a w a 2015 ( https://arxiv.org/abs/1505.01664 ).
This sk ewness is a go o d measure of ho w m uc h fain t signal w e ha v e in the distribution. The
sk ewness can b e accurately measured b y the difference in the mean and median (assuming
no strong outliers): the more distan t they are, the more sk ew ed the dataset is. This imp or-
tan t concept will b e discussed more extensiv ely in the next section (Section 2.2.3 [Sk ewness
caused b y signal and its measuremen t], page 88).
Ho w ev er, sk ewness is only a pro xy for signal when the signal has structure (v aries p er
pixel). Therefore, when it is appro ximately constan t o v er a whole tile, or sub-set of the
image, the constan t signal’s effect is just to shift the symmetric cen ter of the noise distribu-
tion to the p ositiv e and there will not b e an y sk ewness (ma jor difference b et w een the mean
and median). This p ositive 29 shift that preserv es the symmetric distribution is the Sky
v alue. When there is a gradien t o v er the dataset, differen t tiles will ha v e differen t constan t
shifts/Sky-v alues, for example, see Figure 11 of Akhlaghi and Ic hik a w a 2015 ( https://
arxiv.org/abs/1505.01664 ).
T o mak e this v ery large diffuse/flat signal detectable, you will therefore need a larger tile
to con tain a larger c hange in the v alues within it (and impro v e n um b er statistics, for less
scatter when measuring the mean and median). So let’s play with the tessellation a little
to see ho w it affects the result. In Gn uastro, y ou can see the option v alues ( --tilesize in
this case) b y adding the -P option to y our last command. T ry running NoiseChisel with -P
to see its default tile size.
Y ou can clearly see that the default tile size is indeed m uc h smaller than this (h uge)
galaxy and its tidal features. As a result, NoiseChisel w as unable to iden tify the sk ewness
within the tiles under the outer parts of M51 and NGC 5159 and the threshold has b een
o v er-estimated on those tiles. T o see whic h tiles w ere used for estimating the quan tile
threshold (no sk ewness w as measured), y ou can use NoiseChisel’s --checkqthresh option:
$ astnoisechisel r.fits -h0 --checkqthresh
Did y ou see ho w NoiseChisel ab orted after finding and applying the quan tile thresholds?
When y ou call an y of NoiseChisel’s --check* options, b y default, it will ab ort as so on as
all the c hec k steps ha v e b een written in the c hec k file (a m ulti-extension FITS file). This
allo ws y ou to fo cus on the problem y ou w anted to c hec k as so on as p ossible (y ou can disable
this feature with the --continueaftercheck option).
28 b y constan t, w e mean that it has a single v alue in the region w e are measuring.
29 In pro cessed images, where the Sky v alue can b e o v er-estimated, this constan t shift can b e negativ e.
Chapter 2: T utorials 84
T o optimize the threshold-related settings for this image, let’s pla y with this quan tile
threshold c hec k image a little. Do not forget that “ Go o d statistic al analysis is not a pur ely
r outine matter, and gener al ly c al ls for mor e than one p ass thr ough the c omputer ” (Anscom b e
1973, see Section 1.3 [Gn uastro manifesto: Science and its to ols], page 6). A go o d scien tist
m ust ha v e a go o d understanding of her to ols to mak e a meaningful analysis. So do not
hesitate in pla ying with the default configuration and reviewing the man ual when y ou ha v e
a new dataset (from a new instrumen t) in fron t of y ou. Robust data analysis is an art,
therefore a go o d scien tist must first be a go o d artist. So let’s op en the c heck image as a
m ulti-extension cub e:
$ ds9 -mecube r_qthresh.fits -zscale -cmap sls -zoom to fit
The first extension (called CONVOLVED ) of r_qthresh.fits is the con v olv ed input im-
age where the threshold(s) is(are) defined (and later applied to). F or more on the effect
of con v olution and thresholding, see Sections 3.1.1 and 3.1.2 of Akhlaghi and Ic hik a w a
2015 ( https://arxiv.org/abs/1505.01664 ). The second extension ( QTHRESH_ERODE ) has
a blank/white v alue for all the pixels of an y tile that w as iden tified as ha ving significan t
signal. The other tiles ha v e the measured threshold o v er them. The next t w o extensions
( QTHRESH_NOERODE and QTHRESH_EXPAND ) are the other t wo quan tile thresholds that are
necessary in NoiseChisel’s later steps. Ev ery step in this file is rep eated on the three
thresholds.
Pla y a little with the color bar of the QTHRESH_ERODE extension, y ou clearly see ho w the
non-blank tiles around NGC 5195 ha v e a gradien t. As one line of attac k against discarding
to o m uc h signal b elo w the threshold, NoiseChisel rejects outlier tiles. Go forw ard b y three
extensions to VALUE1_NO_OUTLIER and y ou will see that man y of the tiles o v er the galaxy
ha v e b een remo v ed in this step. F or more on the outlier rejection algorithm, see the latter
half of Section 7.1.4.3 [Quan tifying signal in a tile], page 523.
Ev en though m uc h of the galaxy’s fo otprin t has b een rejected as outliers, there are still
tiles with signal remaining: pla y w ith the DS9 color-bar and y ou still see a gradien t near
the outer tidal feature of the galaxy . Before trying to correct this, let’s lo ok at the other
extensions of this c hec k image. W e will use a * as a wild-card that can b e 1, 2 or 3. In the
THRESH*_INTERP extensions, you see that all the blan k tiles hav e b een in terp olated using
their nearest neigh b ors (the relev an t option here is --interpnumngb ). In the follo wing
THRESH*_SMOOTH extensions, you can see the tile v alues after smo othing (configured with
--smoothwidth option). Finally , in QTHRESH-APPLIED , y ou see the thresholded image:
pixels with a v alue of 1 will b e ero ded later, but pixels with a v alue of 2 will pass the
erosion step un-touc hed.
Let’s get bac k to the problem of optimizing the result. Y ou ha v e t w o strategies for
detecting the outskirts of the merging galaxies: 1) Increase the tile size to get more accurate
measuremen ts of sk ewness. 2) Strengthen the outlier rejection parameters to discard more
of the tiles with signal (primarily b y increasing --outliernumngb ). F ortunately in this
image w e ha v e a sufficien tly large region on the righ t side of the image that the galaxy do es
not extend to. So w e can use the more robust first solution. In situations where this do es
not happ en (for example, if the field of view in this image was shifted to the left to ha v e
more of M51 and less sky) y ou are limited to a com bination of the t w o solutions or just to
the second solution.
Chapter 2: T utorials 85
Skipping con v olution for faster tests: The slo w est step of NoiseChisel is the con v olution of
the input dataset. Therefore when y our dataset is large (unlik e the one in this test), and
y ou are not c hanging the input dataset or k ernel in m ultiple runs (as in the tests of this
tutorial), it is faster to do the con v olution separately once (using Section 6.3 [Conv olv e],
page 470) and use NoiseChisel’s --convolved option to directly feed the con v olv ed image
and a v oid con v olution. F or more on --convolved , see Section 7.2.2.1 [NoiseChisel input],
page 549.
T o b etter iden tify the skewness caused b y the flat NGC 5195 and M51 tidal features on
the tiles under it, we ha v e to c ho ose a larger tile size. Let’s try a tile size of 100 b y 100
pixels and insp ect the c hec k image.
$ astnoisechisel r.fits -h0 --tilesize=100,100 --checkqthresh
$ ds9 -mecube r_qthresh.fits -zscale -cmap sls -zoom to fit
Y ou can clearly see the effect of this increased tile size: the tiles are m uc h larger and
when y ou lo ok in to VALUE1_NO_OUTLIER , y ou see that all the tiles are nicely group ed on the
righ t side of the image (the farthest from M51, where w e do not see a gradien t in QTHRESH_
ERODE ). Things lo ok go o d no w, so let’s remo v e --checkqthresh and let NoiseChisel pro ceed
with its detection.
$ astnoisechisel r.fits -h0 --tilesize=100,100
$ ds9 -mecube r_detected.fits -zscale -cmap sls -zoom to fit
The detected pixels of the DETECTIONS extension ha ve expanded a little, but not as
m uc h. Also, the gradien t in the SKY image is almost fully remo ved (and does not fall ov er
M51 an ymore). Ho w ev er, on the b ottom-right of the m51 detection, w e see man y holes
gradually increasing in size. This hin ts that there is still signal out there. Let’s chec k the
next series of detection steps b y adding the --checkdetection option this time:
$ astnoisechisel r.fits -h0 --tilesize=100,100 --checkdetection
$ ds9 -mecube r_detcheck.fits -zscale -cmap sls -zoom to fit
The output no w has 16 extensions, sho wing ev ery step that is tak en b y NoiseChisel.
The first and second ( INPUT and CONVOLVED ) are clear from their names. The third
( THRESHOLDED ) is the thresholded image after finding the quan tile threshold (last extension
of the output of --checkqthresh ). The fourth HDU ( ERODED ) is new: it is the name-stak e
of NoiseChisel, or ero ding pixels that are ab ov e the threshold. By erosion, w e mean that
all pixels with a v alue of 1 (ab o ve the threshold) that are touc hing a pixel with a v alue of
0 (b elo w the threshold) will b e flipp ed to zero (or “carv ed” out) 30 . Y ou can see its effect
directly b y going bac k and forth b et w een the THRESHOLDED and ERODED extensions.
In the fifth extension ( OPENED-AND-LABELED ) the image is “op ened”, whic h is a name
for ero ding once, then dilating (dilation is the in v erse of erosion). This is go o d to remo v e
thin connections that are only due to noise. Eac h separate connected group of pixels is also
giv en its unique lab el here. Do you see ho w just b ey ond the large M51 detection, there
are man y smaller detections that get smaller as y ou go more distan t? This hin ts at the
30 Pixels with a v alue of 2 are very high signal-to-noise pixels, they are not eroded, to preserv e sharp and
brigh t sources.
Chapter 2: T utorials 86
solution: the default num b er of erosions is to o m uch. Let’s see ho w man y erosions tak e
place b y default (b y adding -P | grep erode to the previous command)
$ astnoisechisel r.fits -h0 --tilesize=100,100 -P | grep erode
W e see that the v alue of erode is 2 . The default NoiseChisel parameters are primarily
targeted to pro cessed images (where there is correlated noise due to all the pro cessing
that has gone in to the warping and stac king of ra w images, see Section 2.1.11 [NoiseChisel
optimization for detection], page 41). In those scenarios 2 erosions are commonly necessary .
But here, w e ha v e a single-exp osure image where there is no correlated noise (the pixels are
not mixed). So let’s see ho w things c hange with only one erosion:
$ astnoisechisel r.fits -h0 --tilesize=100,100 --erode=1 \
--checkdetection
$ ds9 -mecube r_detcheck.fits -zscale -cmap sls -zoom to fit
Lo oking at the OPENED-AND-LABELED extension again, w e see that the main/large detec-
tion is no w m uc h larger than b efore. While the immediately-outer connected regions are
still presen t, they ha v e decreased dramatically , so w e can pass this step.
After the OPENED-AND-LABELED extension, NoiseChisel go es on to finding false detec-
tions using the undetected pixels. The pro cess is fully describ ed in Section 3.1.5. (Defining
and Remo ving F alse Detections) of Akhlaghi and Ic hik a w a 2015 ( https://arxiv.org/pdf/
1505.01664.pdf ). Please compare the extensions to what y ou read there and things will b e
v ery clear. In the last HDU ( DETECTION-FINAL ), w e hav e the final detected pixels that will
b e used to estimate the Sky and its Standard deviation. W e see that the main detection
has indeed b een detected v ery far out, so let’s see ho w the full NoiseChisel will estimate the
Sky and its standard deviation (b y remo ving --checkdetection ):
$ astnoisechisel r.fits -h0 --tilesize=100,100 --erode=1
$ ds9 -mecube r_detected.fits -zscale -cmap sls -zoom to fit
The DETECTIONS extension of r_detected.fits closely follo ws what the DETECTION-
FINAL of the c hec k image (lo oks go o d!). If y ou go ahead to the SKY extension, things still
lo ok go o d. But it can still b e impro v ed.
Lo ok at the DETECTIONS again, y ou will see the right-w ard edges of M51’s detected pixels
ha v e man y “holes” that are fully surrounded b y signal (v alue of 1 ) and the signal stretc hes
out in the noise v ery thinly (the size of the holes increases as w e go out). This suggests
that there is still undetected signal and that w e can still dig deep er in to the noise.
With the --detgrowquant option, NoiseChisel will “gro w” the detections in to the noise.
Its v alue is the ultimate limit of the gro wth in units of quantile (betw een 0 and 1). Therefore
--detgrowquant=1 means no gro wth and --detgrowquant=0.5 means an ultimate limit of
the Sky lev el (whic h is usually to o m uc h and will co v er the whole image!). See Figure 2 of
Akhlaghi 2019 ( https://arxiv.org/pdf/1909.11230.pdf ) for more on this option. T ry
running the previous command with v arious v alues (from 0.6 to higher v alues) to see this
option’s effect on this dataset. F or this particularly h uge galaxy (with signal that extends
v ery gradually in to the noise), w e will set it to 0.75 :
$ astnoisechisel r.fits -h0 --tilesize=100,100 --erode=1 \
--detgrowquant=0.75
$ ds9 -mecube r_detected.fits -zscale -cmap sls -zoom to fit
Bey ond this lev el (smaller --detgrowquant v alues), y ou see man y of the smaller bac k-
ground galaxies (to w ards the righ t side of the image) starting to create thin spider-leg-lik e
Chapter 2: T utorials 87
features, sho wing that w e are follo wing correlated noise for to o m uch. Please try it for
y ourself b y c hanging it to 0.6 for example.
When y ou lo ok at the DETECTIONS extension of the command sho wn ab o v e, y ou see the
wings of the galaxy b eing detected m uch farther out, But y ou also see man y holes whic h are
clearly just caused b y noise. After gro wing the ob jects, NoiseChisel also allo ws y ou to fill
suc h holes when they are smaller than a certain size through the --detgrowmaxholesize
option. In this case, a maxim um area/size of 10,000 pixels seems to b e go o d:
$ astnoisechisel r.fits -h0 --tilesize=100,100 --erode=1 \
--detgrowquant=0.75 --detgrowmaxholesize=10000
$ ds9 -mecube r_detected.fits -zscale -cmap sls -zoom to fit
When lo oking at the ra w input image (which is v ery “shallo w”: less than a min ute ex-
p osure!), y ou do not see an ything so far out of the galaxy . Y ou migh t just think to y ourself
that “this is all noise, I ha v e just dug to o deep and I’m follo wing systematics”! If y ou feel
lik e this, ha v e a lo ok at the deep images of this system in W atkins 2015 ( https://arxiv.
org/abs/1501.04599 ), or a 12 hour deep image of this system (with a 12-inc h telescop e):
https://i.redd.it/jfqgpqg0hfk11.jpg 31 . In these deep er images y ou clearly see ho w the
outer edges of the M51 group follo w this exact structure, b elo w in Section 2.2.5 [Ac hiev ed
surface brigh tness lev el], page 97, w e will measure the exact lev el.
As the gradien t in the SKY extension sho ws, and the deep images cited ab o v e confirm,
the galaxy’s signal extends ev en b eyond this. But this is already far deep er than what most
(if not all) other to ols can detect. Therefore, w e will stop configuring NoiseChisel at this
p oin t in the tutorial and let y ou pla y with the other options a little more, while reading
more ab out it in the pap ers: Akhlaghi and Ic hik a w a 2015 ( https://arxiv.org/abs/1505.
01664 ) and 2019 ( https://arxiv.org/abs/1909.11230 ) and Section 7.2 [NoiseChisel],
page 544. When y ou do find a b etter configuration feel free to con tact us for feedbac k. Do
not forget that go o d data analysis is an art, so lik e a sculptor, master y our c hisel for a go o d
result.
T o a v oid t yping all these options ev ery time y ou run NoiseChisel on this image, y ou can
use Gn uastro’s configuration files, see Section 4.2 [Configuration files], page 267. F or an
applied example of setting/using them, see Section 2.1.8 [Option managemen t and config-
uration files], page 35.
31 The image is tak en from this Reddit discussion: https://www.reddit.com/r/Astronomy/comments/
9d6x0q/12_hours_of_exposure_on_the_whirlpool_galaxy/
Chapter 2: T utorials 88
This NoiseChisel configuration is NOT GENERIC: Do not use the configuration deriv ed
ab o v e, on another instrumen t’s image blind ly . If y ou are unsure, just use the default
v alues. As y ou sa w ab o v e, the reason w e ch ose this particular configuration for NoiseChisel
to detect the wings of the M51 group w as strongly influenced b y the noise prop erties of
this particular image. Remem b er Section 2.1.11 [NoiseChisel optimization for detection],
page 41, where w e lo ok ed in to the v ery deep XDF image whic h had strong correlated
noise?
As long as y our other images ha ve similar noise prop erties (from the same data-
reduction step of the same instrumen t), you can use y our configuration on an y of them.
But for images from other instrumen ts, please follo w a similar logic to what w as presen ted
in these tutorials to find the optimal configuration.
Smart NoiseChisel: As you sa w during this section, there is a clear logic b ehind the optimal
parameter v alue for eac h dataset. Therefore, w e plan to add capabilities to (optionally)
automate some of the c hoices made here based on the actual dataset, please join us in
doing this if y ou are in terested. Ho w ev er, given the man y problems in existing “smart”
solutions, suc h automatic c hanging of the configuration ma y cause more problems than
they solv e. So ev en when they are implemen ted, w e w ould strongly recommend qualit y
c hec ks for a robust analysis.
2.2.3 Sk ewness caused b y signal and its measuremen t
In the previous section (Section 2.2.2 [NoiseChisel optimization], page 82) w e sho w ed ho w
to customize NoiseChisel for a single-exp osure SDSS image of the M51 group. During
the customization, w e also discussed the sk ewness caused b y signal. In the next section
(Section 2.2.4 [Image surface brigh tness limit], page 92), w e will use this to measure the
surface brigh tness limit of the image. Ho wev er, to b etter understand NoiseChisel and also,
the image surface brigh tness limit, understanding the sk ewness caused b y signal, and how
to measure it prop erly are v ery imp ortant. Therefore now that w e hav e separated signal
from noise, let’s pause for a momen t and lo ok in to sk ewness, ho w signal creates it, and find
the b est w a y to measure it.
Let’s start masking all the detected pixels found at the end of the previous section
(Section 2.2.2 [NoiseChisel optimization], page 82) and having a look at the noise distri-
bution with Gn uastro’s Arithmetic and Statistics programs as sho wn b elo w (while visually
insp ecting the mask ed image with DS9 in the middle).
$ astarithmetic r_detected.fits -hINPUT-NO-SKY set-in \
r_detected.fits -hDETECTIONS set-det \
in det nan where -odet-masked.fits
$ ds9 det-masked.fits
$ aststatistics det-masked.fits
Y ou will see that Gn uastro’s Statistics program prin ts an ASCI I histogram when no option
is giv en (it is sho wn b elo w). This is done to give y ou a fast and easy view of the distribution
of v alues in the dataset (pixels in an image, or ro ws in a table’s column).
Chapter 2: T utorials 89
-------
Input: det-masked.fits (hdu: 1)
-------
Number of elements: 903920
Minimum: -0.113543
Maximum: 0.130339
Median: -0.00216306
Mean: -0.0001893073877
Standard deviation: 0.02569057188
-------
Histogram:
| ** *
| * ** * *
| ** ** * *
| * ** ** ** *
| ** ** ** ** * **
| ** ** ** ** * ** *
| * ** ** ** ** * ** **
| ** ** ** ** **** ** ** *
| ** ** ** ** ** **** ** ** ** *
| ** ** ** ** ** ** ******* ** ** ** *
|*********** ** ** ** ******************* ** ** ** ** ***** ** ***** **
|----------------------------------------------------------------------
This histogram sho ws a roughly symmetric noise distribution, so let’s ha v e a lo ok at its
sk ewness. The most commonly used definition of sk ewness is known as the “P earson’s first
sk ewness co efficien t”. It measures the difference b et ween the mean and median, in units of
the standard deviation (STD):
Sk ewness ≡ (mean − median)
STD
The logic b ehind this definition is simple: as more signal is added to the same pixels
that originally only ha v e ra w noise (sk ewness is increased), the mean shifts to the p ositiv e
faster than the median, so the distance b et w een the mean and median should increase. Let’s
measure the sk ewness (as defined ab o v e) o v er the image without an y signal. Its v ery easy
with Gn uastro’s Statistics program (and piping the output to A WK):
$ aststatistics det-masked.fits --mean --median --std \
| awk '{print ($1-$2)/$3}'
0.0768279
W e see that the mean and median are only 0 . 08 σ (rounded) a w a y from eac h other (whic h is
v ery close)! All pixels with significan t signal are mask ed, so this is exp ected, and ev erything
is fine. No w, let’s c hec k the pixel distribution of the sky-subtracted input (where pixels with
significan t signal remain, and are not mask ed):
$ ds9 r_detected.fits
$ aststatistics r_detected.fits -hINPUT-NO-SKY
-------
Chapter 2: T utorials 96
--upmaskfile=r_detected.fits --upmaskhdu=DETECTIONS \
--upnsigma=3 --checkuplim=1 --upnum=1000 \
--ids --upperlimit-sb
The sbl.fits catalog no w con tains the upp er-limit surface brigh tness for a circle with
an area of 25 arcsec 2 . Y ou can c hec k the v alue with the command b elo w, but the great
thing is that no w y ou ha v e b oth of the surface brigh tness limiting magnitude in the headers
discussed ab o v e, and the upp er-limit surface brigh tness within the table. Y ou can also
add more profiles with differen t shap es and sizes if necessary . Of course, y ou can also use
--upperlimit-sb in y our actual science ob jects and clumps to get an ob ject-sp ecific or
clump-sp ecific v alue.
$ asttable sbl.fits -cUPPERLIMIT_SB
25.9119
Y ou will get a sligh tly different v alue from the command ab ov e. In fact, if y ou run the
Mak eCatalog command again and lo ok at the measured upp er-limit surface brigh tness, it
will b e sligh tly differen t with y our first trial! Please try exactly the same Mak eCatalog
command ab o v e a few times to see ho w it c hanges.
This is b ecause of the r andom factor in the upp er-limit measuremen ts: ev ery time y ou
run it, differen t random p oin ts will b e c hec k ed, resulting in a slightly differen t distribution.
Y ou can decrease the random scatter b y increasing the n um b er of random c hec ks (for
example, setting --upnum=100000 , compared to 1000 in the command ab o v e). But this
will b e slo w er and the results will not b e exactly repro ducible. The only w a y to ensure y ou
get an iden tical result later is to fix the random n um b er generator function and seed lik e the
command b elo w 34 . This is a v ery imp ortan t p oin t regarding an y statistical pro cess in v olving
random n um b ers, please see Section 6.2.3.4 [Generating random n um b ers], page 405.
export GSL_RNG_TYPE=ranlxs1
export GSL_RNG_SEED=1616493518
astmkcatalog lab.fits -h1 --zeropoint=$zeropoint -osbl.fits \
--sfmagarea=25 --sfmagnsigma=3 --forcereadstd \
--valuesfile=r_detected.fits --valueshdu=INPUT-NO-SKY \
--upmaskfile=r_detected.fits --upmaskhdu=DETECTIONS \
--upnsigma=3 --checkuplim=1 --upnum=1000 \
--ids --upperlimit-sb --envseed
But where do all the random ap ertures of the upp er-limit measuremen t fall on the image?
It is go o d to actually insp ect their lo cation to get a b etter understanding for the pro cess and
also detect p ossible bugs/biases. When MakeCatalog is run with the --checkuplim option,
it will prin t all the random lo cations and their measured brigh tnes s as a table in a file with
the suffix _upcheck.fits . With the first command b elo w y ou can use Gn uastro’s asttable
and astscript-ds9-region to con v ert the successful ap erture lo cations in to a DS9 region
file, and with the second can load the region file in to the detections and sky-subtracted
image to visually see where they are.
## Create a DS9 region file from the check table (activated
## with '--checkuplim')
asttable sbl_upcheck.fits --noblank=RANDOM_SUM \
34 Y ou can use any in teger for the seed. One recommendation is to run Mak eCatalog without --envseed
once and use the randomly generated seed that is prin ted on the terminal.
Chapter 2: T utorials 97
| astscript-ds9-region -c1,2 --mode=img \
--radius=$r_pixel
## Have a look at the regions in relation with NoiseChisel's
## detections.
ds9 r_detected.fits[INPUT-NO-SKY] -regions load ds9.reg
ds9 r_detected.fits[DETECTIONS] -regions load ds9.reg
In this example, w e w ere lo oking at a single-exp osure image that has no correlated
noise. Because of this, the surface brigh tness limit and the upp er-limit surface brigh tness
are v ery close. They will ha v e a bigger difference on deep datasets with stronger correlated
noise (that are the result of stac king man y individual exp osures). As an exercise, please
try measuring the upp er-limit surface brigh tness level and surface brigh tness limit for the
deep HST data that w e used in the previous tutorial (Section 2.1 [General program usage
tutorial], page 22).
2.2.5 Ac hiev ed surface brigh tness lev el
In Section 2.2.2 [NoiseChisel optimization], page 82, w e customized NoiseChisel for a single-
exp osure SDSS image of the M51 group and in Section 2.2.4 [Image surface brigh tness limit],
page 92, we measured the surface brigh tness limit and the upp er-limit surface brigh tness
lev el (whic h are b oth measures of the noise lev el). In this section, let’s do some mea-
suremen ts on the outer-most edges of the M51 group to see ho w they relate to the noise
measuremen ts found in the previous section.
F or this measuremen t, we will need to estimate the a v erage flux on the outer edges of the
detection. F ortunately all this can b e done with a few simple commands using Section 6.2
[Arithmetic], page 398, and Section 7.4 [Mak eCatalog], page 574. First, let’s separate eac h
detected region, or giv e a unique lab el/coun ter to all the connected pixels of NoiseChisel’s
detection map with the command b elo w. Recall that with the set- op erator, the p opp ed
op erand will b e giv en a name ( det in this case) for easy usage later.
$ astarithmetic r_detected.fits -hDETECTIONS set-det \
det 2 connected-components -olabeled.fits
Y ou can find the lab el of the main galaxy visually (b y op ening the image and ho v ering
y our mouse o v er the M51 group’s lab el). But to ha v e a little more fun, let’s do this auto-
matically (whic h is necessary in a general scenario). The M51 group detection is b y far the
largest detection in this image, this allo ws us to find its ID/lab el easily . W e will first run
Mak eCatalog to find the area of all the lab els, then w e will use T able to find the ID of the
largest ob ject and k eep it as a shell v ariable ( id ):
# Run MakeCatalog to find the area of each label.
$ astmkcatalog labeled.fits --ids --geo-area -h1 -ocat.fits
## Sort the table by the area column.
$ asttable cat.fits --sort=AREA_FULL
## The largest object, is the last one, so we will use '--tail'.
$ asttable cat.fits --sort=AREA_FULL --tail=1
## We only want the ID, so let's only ask for that column:
Chapter 2: T utorials 98
$ asttable cat.fits --sort=AREA_FULL --tail=1 --column=OBJ_ID
## Now, let's put this result in a variable (instead of printing)
$ id=$(asttable cat.fits --sort=AREA_FULL --tail=1 --column=OBJ_ID)
## Just to confirm everything is fine.
$ echo $id
W e can no w use the id v ariable to reject all other detections:
$ astarithmetic labeled.fits $id eq -oonly-m51.fits
Op en the image and ha ve a lo ok. T o separate the outer edges of the detections, w e will
need to “ero de” the M51 group detection. So in the same Arithmetic command as ab o v e,
w e will ero de three times (to ha v e more pixels and th us less scatter), using a maxim um
connectivit y of 2 (8-connected neigh b ors). W e will then sa v e the output in eroded.fits .
$ astarithmetic labeled.fits $id eq 2 erode 2 erode 2 erode \
-oeroded.fits
In labeled.fits , w e can now set all the 1-v alued pixels of eroded.fits to 0 using Arith-
metic’s where op erator added to the previous command. W e will need the pixels of the
M51 group in labeled.fits t wo times: once to do the erosion, another time to find the
outer pixel la y er. T o do this (and b e efficien t and more readable) w e will use the set-i
op erator (to giv e this image the name ‘ i ’). In this w a y w e can use it an y n um b er of times
afterw ards, while only reading it from disk and finding M51’s pixels once.
$ astarithmetic labeled.fits $id eq set-i i \
i 2 erode 2 erode 2 erode 0 where -oedge.fits
Op en the image and ha v e a lo ok. Y ou’ll see that the detected edge of the M51 group is
no w clearly visible. Y ou can use edge.fits to mark (set to blank) this b oundary on the
input image and get a visual feeling of ho w far it extends:
$ astarithmetic r.fits -h0 edge.fits nan where -oedge-masked.fits
T o quan tify ho w deep we ha v e detected the lo w-surface brigh tness regions (in units of
signal to-noise ratio), w e will use the command b elow. In short it just divides all the
non-zero pixels of edge.fits in the Sky subtracted input (first extension of NoiseChisel’s
output) b y the pixel standard deviation of the same pixel. This will giv e us a signal-to-noise
ratio image. The mean v alue of this image sho ws the lev el of surface brigh tness that w e
ha v e ac hiev ed. Y ou can also break the command b elo w in to m ultiple calls to Arithmetic and
create temp orary files to understand it b etter. Ho w ev er, if you ha v e a lo ok at Section 6.2.1
[Rev erse p olish notation], page 398, and Section 6.2.4 [Arithmetic op erators], page 407, y ou
should b e able to easily understand what y our computer do es when y ou run this command 35 .
$ astarithmetic edge.fits -h1 set-edge \
35 edge.fits (extension 1 ) is a binary (0 or 1 v alued) image. Applying the not op erator on it, just
flips all its pixels (from 0 to 1 and vice-v ersa). Using the where op erator, w e are then setting all
the newly 1-v alued pixels (pixels that are not on the edge) to NaN/blank in the sky-subtracted input
image ( r_detected.fits , extension INPUT-NO-SKY , whic h w e call skysub ). W e are then dividing all the
non-blank pixels (only those on the edge) b y the sky standard deviation ( r_detected.fits , extension
SKY_STD , which w e called skystd ). This gives the signal-to-noise ratio (S/N) for eac h of the pixels on
the b oundary . Finally , with the meanvalue op erator, w e are taking the mean v alue of all the non-blank
pixels and rep orting that as a single n um b er.
Chapter 2: T utorials 99
r_detected.fits -hSKY_STD set-skystd \
r_detected.fits -hINPUT-NO-SKY set-skysub \
skysub skystd / edge not nan where meanvalue --quiet
W e ha v e th us detected the wings of the M51 group do wn to roughly 1/3rd of the noise
lev el in this image whic h is a v ery go o d ac hiev emen t! But the p er-pixel S/N is a relative
measuremen t. Let’s also measure the depth of our detection in absolute surface brigh tness
units; or magnitudes p er square arc-seconds (see Section 7.4.2 [Brigh tness, Flux, Magnitude
and Surface brigh tness], page 576). W e will also ask for the S/N and magnitude of the full
edge w e ha v e defined. F ortunately doing this is v ery easy with Gn uastro’s Mak eCatalog:
$ astmkcatalog edge.fits -h1 --valuesfile=r_detected.fits \
--zeropoint=22.5 --ids --sb --sn --magnitude
$ asttable edge_cat.fits
1 25.6971 55.2406 15.8994
W e ha v e th us reac hed an outer surface brigh tness of 25 . 70 magnitudes/arcsec 2 (second
column in edge_cat.fits ) on this single exp osure SDSS image! This is v ery similar to the
surface brigh tness limit measured in Section 2.2.4 [Image surface brigh tness limit], page 92,
(whic h is a big ac hiev emen t!). But another p oint in the result abov e is v ery in teresting:
the total S/N of the edge is 55 . 24 with a total edge magnitude 36 of 15.90!!! This is v ery
large for suc h a fain t signal (recall that the mean S/N p er pixel w as 0.32) and sho ws a v ery
imp ortan t p oin t in the study of galaxies: While the p er-pixel signal in their outer edges
ma y b e v ery faint (and in visible to the ey e in noise), a lot of signal hides deeply buried in
the noise.
In in terpreting this v alue, y ou should just hav e in mind that NoiseChisel w orks based on
the con tiguit y of signal in the pixels. Therefore the larger the ob ject, the deep er NoiseChisel
can carv e it out of the noise (for the same outer surface brigh tness). In other w ords, this
rep orted depth, is the depth w e ha v e reac hed for this ob ject in this dataset, pro cessed with
this particular NoiseChisel configuration. If the M51 group in this image was larger/smaller
than this (the field of view w as smaller/larger), or if the image w as from a differen t instru-
men t, or if w e had used a differen t configuration, w e w ould go deep er/shallo w er.
2.2.6 Extract clumps and ob jects (Segmen tation)
In Section 2.2.2 [NoiseChisel optimization], page 82, we found a goo d detection map o v er the
image, so pixels harb oring signal hav e b een differen tiated from those that do not. F or noise-
related measuremen ts lik e the surface brigh tness limit, this is fine. Ho w ev er, after finding
the pixels with signal, y ou are most lik ely in terested in kno wing the sub-structure within
them. F or example, ho w many star forming regions (those brigh t dots along the spiral
arms) of M51 are within this image? What are the colors of eac h of these star forming
regions? In the outer most wings of M51, whic h pixels b elong to bac kground galaxies and
foreground stars? And man y more similar questions. T o address these questions, y ou can
use Section 7.3 [Segmen t], page 563, to iden tify all the “clumps” and “ob jects” o v er the
detection.
$ astsegment r_detected.fits --output=r_segmented.fits
$ ds9 -mecube r_segmented.fits -cmap sls -zoom to fit -scale limits 0 2
36 Y ou can run MakeCatalog on only-m51.fits instead of edge.fits to see the full magnitude of the M51
group in this image.
Chapter 2: T utorials 100
Op en the output r_segmented.fits as a m ulti-extension data cub e with the second
command ab o v e and flip through the first and second extensions, zo om-in to the spiral
arms of M51 and see the detected clumps (all pixels with a v alue larger than 1 in the second
extension). T o optimize the parameters and mak e sure y ou ha v e detected what y ou w anted,
w e recommend to visually insp ect the detected clumps on the input image.
F or visual insp ection, y ou can mak e a simple shell script lik e b elo w. It will first call
Mak eCatalog to estimate the p ositions of the clumps, then mak e an SA O DS9 region file
and op en ds9 with the image and region file. Recall that in a shell script, the n umeric
v ariables (lik e $1 , $2 , and $3 in the example b elo w) represen t the argumen ts giv en to the
script. But when used in the A WK argumen ts, they refer to column n um b ers.
T o create the shell script, using y our fa v orite text editor, put the con ten ts b elo w in to a
file called check-clumps.sh . Recall that ev erything after a # is just commen ts to help y ou
understand the command (so read them!). Also note that if y ou are copying from the PDF
v ersion of this b o ok, fix the single quotes in the A WK command.
#! /bin/bash
set -e # Stop execution when there is an error.
set -u # Stop execution when a variable is not initialized.
# Run MakeCatalog to write the coordinates into a FITS table.
# Default output is `$1_cat.fits'.
astmkcatalog $1.fits --clumpscat --ids --ra --dec
# Use Gnuastro's Table and astscript-ds9-region to build the DS9
# region file (a circle of radius 1 arcseconds on each point).
asttable $1"_cat.fits" -hCLUMPS -cRA,DEC \
| astscript-ds9-region -c1,2 --mode=wcs --radius=1 \
--output=$1.reg
# Show the image (with the requested color scale) and the region file.
ds9 -geometry 1800x3000 -mecube $1.fits -zoom to fit \
-scale limits $2 $3 -regions load all $1.reg
# Clean up (delete intermediate files).
rm $1"_cat.fits" $1.reg
Finally , you just ha v e to activ ate the script’s executable flag with the command b elo w. This
will enable y ou to directly/easily call the script as a command.
$ chmod +x check-clumps.sh
This script do es not exp ect the .fits suffix of the input’s filename as the first argumen t.
Because the script pro duces in termediate files (a catalog and DS9 region file, which are later
deleted). Ho wev er, w e do not w an t m ultiple instances of the script (on differen t files in the
same directory) to collide (read/write to the same in termediate files). Therefore, w e ha v e
used suffixes added to the input’s name to iden tify the in termediate files. Note ho w all the
$1 instances in the commands (not within the A WK command 37 ) are follo w ed b y a suffix.
If y ou w an t to k eep the in termediate files, put a # at the start of the last line.
37 In A WK, $1 refers to the first column, while in the shell script, it refers to the first argumen t.
Chapter 2: T utorials 101
The few, but high-v alued, bright pixels in the cen tral parts of the galaxies can hinder easy
visual insp ection of the fain ter parts of the image. With the second and third argumen ts
to this script, y ou can set the numerical v alues of the color map (first is minim um/blac k,
second is maxim um/white). Y ou can call this script with any 38 output of Segmen t (when
--rawoutput is not used) with a command lik e this:
$ ./check-clumps.sh r_segmented -0.1 2
Go ahead and run this command. Y ou will see the in termediate pro cessing b eing done
and finally it op ens SA O DS9 for y ou with the regions sup erimp osed on all the extensions
of Segmen t’s output. The script will only finish (and giv e y ou con trol of the command-line)
when y ou close DS9. If y ou need y our access to the command-line b efore closing DS9, add
a & after the end of the command ab o v e.
While DS9 is op en, slide the dynamic range (v alues for blac k and white, or mini-
m um/maxim um v alues in differen t color sc hemes) and zo om in to v arious regions of the
M51 group to see if y ou are satisfied with the detected clumps. Do not forget that through
the “Cub e” windo w that is op ened along with DS9, you can flip through the extensions and
see the actual clumps also. The questions you should be asking yourself are these: 1) Whic h
real clumps (as y ou visually fe el ) ha v e b een missed? In other words, is the c ompleteness
go o d? 2) Are there an y clumps whic h y ou fe el are false? In other words, is the purity go o d?
Note that completeness and purit y are not indep enden t of eac h other, they are an ti-
correlated: the higher your purit y , the lo w er y our completeness and vice-v ersa. Y ou can see
this b y pla ying with the purit y lev el using the --snquant option. Run Segmen t as sho wn
ab o v e again with -P and see its default v alue. Then increase/decrease it for higher/lo wer
purit y and c hec k the result as b efore. Y ou will see that if y ou w an t the b est purit y , y ou
ha v e to sacrifice completeness and vice v ersa.
One in teresting region to insp ect in this image is the man y brigh t p eaks around the
cen tral parts of M51. Zo om in to that region and insp ect ho w man y of them ha v e actually
b een detected as true clumps. Do y ou ha v e a go o d balance b et w een completeness and
purit y? Also lo ok out far in to the wings of the group and insp ect the completeness and
purit y there.
An easier w a y to insp ect completeness (and only completeness) is to mask all the pixels
detected as clumps and visually insp ecting the rest of the pixels. Y ou can do this using
Arithmetic in a command lik e b elo w. F or easy reading of the command, w e will define the
shell v ariable i for the image name and sa v e the output in masked.fits .
$ in="r_segmented.fits -hINPUT-NO-SKY"
$ clumps="r_segmented.fits -hCLUMPS"
$ astarithmetic $in $clumps 0 gt nan where -oclumps-masked.fits
Insp ecting clumps-masked.fits , y ou can see some v ery diffuse p eaks that ha v e b een
missed, esp ecially as y ou go farther a w a y from the group cen ter and in to the diffuse wings.
This is due to the fact that with this configuration, w e ha v e fo cused more on the sharp er
clumps. T o put the fo cus more on diffuse clumps, y ou can use a wider con v olution k ernel.
38 Some mo difications are necessary based on the input dataset: dep ending on the dynamic range, y ou
ha v e to adjust the second and third arguments. But more imp ortan tly , dep ending on the dataset’s w orld
co ordinate system, y ou hav e to change the region width , in the A WK command. Otherwise the circle
regions can b e to o small/large.
Chapter 2: T utorials 102
Using a larger k ernel can also help in detecting the existing clumps to fain ter lev els (th us
b etter separating them from the surrounding diffuse signal).
Y ou can mak e an y k ernel easily using the --kernel option in Section 8.1 [Mak eProfiles],
page 631. But note that a larger k ernel is also going to w ash-out man y of the sharp/small
clumps close to the cen ter of M51 and also some smaller p eaks on the wings. Please
con tin ue pla ying with Segmen t’s configuration to obtain a more complete result (while
k eeping reasonable purit y). W e will finish the discussion on finding true clumps at this
p oin t.
The prop erties of the clumps within M51, or the bac kground ob jects can then easily
b e measured using Section 7.4 [Mak eCatalog], page 574. T o measure the prop erties of the
bac kground ob jects (detected as clumps o v er the diffuse region), y ou should not mask the
diffuse region. When measuring clump prop erties with Section 7.4 [Mak eCatalog], page 574,
and using the --clumpscat , the am bien t flux (from the diffuse region) is calculated and
subtracted. If the diffuse region is mask ed, its effect on the clump brigh tness cannot b e
calculated and subtracted.
T o k eep this tutorial short, w e will stop here. See Section 2.1.13 [Segmen tation and
making a catalog], page 47, and Section 7.3 [Segmen t], page 563, for more on using Segmen t,
pro ducing catalogs with Mak eCatalog and using those catalogs.
2.3 Building the extended PSF
Deriving the extended PSF of an image is v ery imp ortan t in man y asp ects of the analysis of
the ob jects within it. Gn uastro has a set of installed scripts, designed to simplify the pro cess
follo wing the recip e of Infan te-Sainz et al. 2020 ( https://arxiv.org/abs/1911.01430 );
for more, see Section 10.8 [PSF construction and subtraction], page 703. An ov erview of
the pro cess is giv en in Section 10.8.1 [Ov erview of the PSF scripts], page 704.
2.3.1 Preparing input for extended PSF
W e will use an image of the M51 galaxy group in the r (SDSS) band of the Ja v alambre
Photometric Lo cal Univ erse Surv ey (J-PLUS) to extract its extended PSF. F or more infor-
mation on J-PLUS, and its unique features visit: http://www.j-plus.es .
First, let’s do wnload the image from the J-PLUS w eb page using wget . But to hav e a
generalize-able, and easy to read command, w e will define some base v ariables (in all-caps)
first. After the do wnload is complete, op en the image with SA O DS9 (or an y other FITS
view er y ou prefer!) to ha v e a feeling of the data (and of course, enjo y the b eaut y of M51 in
suc h a wide field of view):
$ urlend="jplus-dr2/get_fits?id=67510"
$ urlbase="http://archive.cefca.es/catalogues/vo/siap/"
$ mkdir jplus-dr2
$ wget $urlbase$urlend -O jplus-dr2/67510.fits.fz
$ astscript-fits-view jplus-dr2/67510.fits.fz
After enjo ying the large field of view, hav e a closer lo ok at the edges of the image. Please
zo om in to the corners. Y ou will see that on the edges, the pixel v alues are either zero or
with significan tly differen t v alues than the main b o dy of the image. This is due to the
Chapter 2: T utorials 103
dithering pattern that w as used to make this image and happ ens in all imaging surv eys 39 .
T o a v oid p oten tial issues or problems that these regions ma y cause, we will first crop out
the main b o dy of the image with the command b elo w. T o k eep the top-lev el directory clean,
let’s also put the crop in a directory called flat .
$ mkdir flat
$ astcrop jplus-dr2/67510.fits.fz --section=225:9275,150:9350 \
--mode=img -oflat/67510.fits
$ astscript-fits-view flat/67510.fits
Please zo om in to the edges again, y ou will see that they no w ha v e the same noise-lev el as
the rest of the image (the problematic parts are no w gone).
2.3.2 Saturated pixels and Segmen t’s clumps
A constan t-depth (flat) image w as created in the previous section (Section 2.3.1 [Preparing
input for extended PSF], page 102). As explained in Section 10.8.1 [Ov erview of the PSF
scripts], page 704, an imp ortan t step when building the PSF is to mask other sources in the
image. Therefore, b efore going onto selecting stars, let’s detect all significan t signal, and
iden tify the clumps of bac kground ob jects o v er the wings of the extended PSF.
There is a problem ho wev er: the saturated pixels of the brigh t stars are going to cause
problems in the segmen tation phase. T o see this problem, let’s make a 1000 × 1000 crop
around a brigh t star to sp eed up the test (and its solution). Afterwards w e will apply the
solution to the whole image.
$ astcrop flat/67510.fits --mode=wcs --widthinpix --width=1000 \
--center=203.3916736,46.7968652 --output=saturated.fits
$ astnoisechisel saturated.fits --output=sat-nc.fits
$ astsegment sat-nc.fits --output=sat-seg.fits
$ astscript-fits-view sat-seg.fits
Ha v e a lo ok at the CLUMPS extension. Y ou will see that instead of a single clump at the
cen ter of the brigh t star, w e ha v e man y clumps! This has happ ened b ecause of the saturated
pixels! When saturation o ccurs, the sharp p eak of the profile is lost (lik e cutting off the tip
of a moun tain to build a telescop e!) and all saturated pixels get a noisy v alue close to the
saturation lev el. T o see this saturation noise run the last command again and in SA O DS9,
set the “Scale” to “min max” and zo om in to the cen ter. Y ou will see the noisy saturation
pixels at the cen ter of the star in red.
This noise-at-the-p eak disrupts Segmen t’s assumption to expand clumps from a lo cal
maxima: eac h noisy p eak is b eing treated as a separate lo cal maxima and th us a separate
clump. F or more on ho w Segmen t defines clumps, see Section 3.2.1 and Figure 8 of Akhlaghi
and Ic hik a w a 2015 ( https://arxiv.org/abs/1505.01664 ). T o ha v e the cen ter iden tified
as a single clump, w e should mask these saturated pixels in a wa y that suites Segmen t’s
non-parametric metho dology .
First w e need to find the saturation lev el! The saturation lev el is usually fixed for an y
surv ey or input data that y ou receiv e from a certain database, so you will usually ha v e to do
this only once (the first time y ou get data from that database). Let’s mak e a smaller crop
of 50 × 50 pixels around the star with the first command b elo w. With the next command,
39 Recall the cropping in a previous tutorial for a similar reason (v arying “depth” across the image):
Section 2.1.4 [Dataset insp ection and cropping], page 25.
Chapter 2: T utorials 104
please lo ok at the crop with DS9 to visually understand the problem. Y ou will see the
saturated pixels as the noisy red pixels in the cen ter of the image. A non-saturated star will
ha v e a single pixel as the maxim um and will not ha v e suc h a large area co v ered b y a noisy
constan t v alue (find a few stars in the image and see for y ourself ). Visual and qualitativ e
insp ection of the pro cess is v ery imp ortan t for understanding the solution.
$ astcrop saturated.fits --mode=wcs --widthinpix --width=50 \
--center=203.3916736,46.7968652 --output=sat-center.fits
$ astscript-fits-view sat-center.fits --ds9scale=minmax
T o quan titativ ely iden tify the saturation lev el in this image, let’s ha v e a lo ok at the distri-
bution of pixels with a v alue larger than 100 (ab o v e the noise lev el):
$ aststatistics sat-center.fits --greaterequal=100
Histogram:
|*
|*
|*
|*
|* *
|** *
|*** **
|**** **
|****** ****
|********** * * * ******
|************************* ************ * *** ******* *** ************
|----------------------------------------------------------------------
The p eak y ou see in the righ t end (larger v alues) of the histogram sho ws the saturated
pixels (a constan t lev el, with some scatter due to the large P oisson noise). If there was no
saturation, the n um b er of pixels should ha v e decreased at increasing v alues; un til reac hing
the maxim um v alue of the profile in one pixel. But that is not the case here. Please try
this exp erimen t on a non-saturated (fain ter) star to see what w e mean.
If y ou still ha v e not exp erimen ted on a non-saturated star, please stop reading this
tutorial! Please op en flat/67510.fits in DS9, select a fain ter/smaller star and rep eat the
last three commands (with a differen t center). After y ou ha v e confirmed the p oin t ab o v e
(visually , and with the histogram), please con tin ue with the rest of this tutorial.
Finding the saturation lev el is easy with Statistics (by using the --lessthan option un til
the histogram b ecomes as exp ected: only decreasing). First, let’s try --lessthan=3000 :
$ aststatistics sat-center.fits --greaterequal=100 --lessthan=3000
-------
Histogram:
|*
|*
|*
|*
|*
|**
|*** *
Chapter 2: T utorials 105
|**** *
|******* **
|*********** * * * * * * * ****
|************************* * ***** ******* ***** ** ***** * ********
|----------------------------------------------------------------------
W e still see an increase in the histogram around 3000. Let’s try a threshold of 2500:
$ aststatistics sat-center.fits --greaterequal=100 --lessthan=2500
-------
Histogram:
|*
|*
|**
|**
|**
|**
|****
|*****
|*********
|************* * * * *
|********************************* ** ** ** *** ** * **** ** *****
|----------------------------------------------------------------------
The p eak at the large end of the histogram has gone! But let’s ha v e a closer lo ok at the
v alues (the resolution of an ASCI I histogram is limited!). T o do this, w e will ask Statistics
to sa v e the histogram in to a table with the --histogram option, then lo ok at the last 20
ro ws:
$ aststatistics sat-center.fits --greaterequal=100 --lessthan=2500 \
--histogram --output=sat-center-hist.fits
$ asttable sat-center-hist.fits --tail=20
2021.1849112701 1
2045.0495397186 0
2068.9141681671 1
2092.7787966156 1
2116.6434250641 0
2140.5080535126 0
2164.3726819611 0
2188.2373104095 0
2212.101938858 1
2235.9665673065 1
2259.831195755 2
2283.6958242035 0
2307.560452652 0
2331.4250811005 1
2355.289709549 1
2379.1543379974 1
2403.0189664459 2
2426.8835948944 1
Chapter 2: T utorials 112
First let’s ha v e a lo ok at all the mask ed p ostage stamps of the cropp ed stars. Once they
all op en, feel free to zo om-in, they are all matc hed and lo c ked. It is alw a ys go o d to c hec k
the differen t stamps to ensure the qualit y and p ossible t w o dimensional features that are
difficult to detect from the radial profiles (suc h as ghosts and in ternal reflections).
$ astscript-fits-view finding-normradii/cropped-masked*.fits
If ev erything lo oks go o d in the image, let’s op en all the radial profiles and visually c hec k
those with the command b elo w. Note that astscript-fits-view calls the topcat graphic
user in terface (GUI) program to visually insp ect (plot) tables. If y ou do not already hav e
it, see Section A.2 [TOPCA T], page 965.
$ astscript-fits-view finding-normradii/rprofile*.fits
After some study of this data, w e could sa y that a go o d normalization ring is those
pixels b et w een R=20 and R=30 pixels. Suc h a ring ensures ha ving a high n um b er of pixels
so the estimation of the flux normalization will b e robust. Also, at suc h distance from the
cen ter the signal to noise is high and there are not ob vious features that can affect the
normalization. Note that the profiles are differen t b ecause w e are considering a wide range
of magnitudes, so the fain ter stars are m uc h more noisy . Ho w ev er, in this tutorial w e will
k eep these stars in order to ha ve a higher n um b er of stars for the outer part. In a real case
scenario, w e should lo ok for stars with a m uc h more similar brigh tness (smaller range of
magnitudes) to not lose signal to noise as a consequence of the inclusion of fain ter stars.
$ rm -r finding-normradii
$ counter=1
$ mkdir outer/stamps
$ asttable outer/67510-6-10.fits \
| while read -r ra dec mag; do
astscript-psf-stamp label/67510-seg.fits \
--mode=wcs \
--nocentering \
--center=$ra,$dec \
--normradii=20,30 \
--widthinpix=1000,1000 \
--segment=label/67510-seg.fits \
--output=outer/stamps/67510-$counter.fits; \
counter=$((counter+1)); \
done
After the stamps are created, w e need to stac k them together with a simple Arithmetic
command (see Section 6.2.4.7 [Stac king op erators], page 421). The stac k is done using the
sigma-clipp ed mean op erator that will preserv e more of the signal, while rejecting outliers
(more than 3 σ with a tolerance of 0 . 2, for more on sigma-clipping see Section 2.10.2 [Sigma
clipping], page 199). Just recall that w e need to sp ecify the n um b er of inputs in to the
stac king op erators, so w e are reading the list of images and coun ting them as separate
v ariables b efore calling Arithmetic.
$ numimgs=$(echo outer/stamps/*.fits | wc -w)
$ astarithmetic outer/stamps/*.fits $numimgs 3 0.2 sigclip-mean \
-g1 --output=outer/stack.fits --wcsfile=none
Chapter 2: T utorials 113
Did y ou notice the --wcsfile=none option ab o v e? With it, the stac k ed image no longer
has an y W CS information. This is natural, b ecause the stac ked image does not corresp ond
to an y sp ecific region of the sky an y more.
Let’s compare this stac k ed PSF with the images of the individual stars that w ere used
to create it. Y ou can clearly see that the num b er of mask ed pixels is significan tly decreased
and the PSF is m uc h more “cleaner”.
$ astscript-fits-view outer/stack.fits outer/stamps/*.fits
Ho w ev er, the saturation in the cen ter still remains. Also, b ecause w e did not ha v e to o
man y images, some regions still are v ery noisy . If we had m ore brigh t stars in our selected
magnitude range, w e could ha v e filled those outer remaining patc hes. In a large surv ey lik e
J-PLUS (that w e are using here), y ou can simply lo ok in to other fields that w ere observ ed
so on b efore/after the image ID 67510 that w e used here (to ha v e a similar PSF) and get
more stars in those images to add to these. In fact, the J-PLUS DR2 image ID of the
field ab o v e w as in ten tionally preserv ed during the steps ab o ve to sho w ho w easy it is to use
images from other fields and blend them all in to the output PSF.
2.3.5 Inner part of the PSF
In Section 2.3.4 [Building outer part of PSF], page 110, w e w ere able to create a stac k of
the outer-most b eha vior of the PSF in a J-PLUS surv ey image. But the cen tral part that
w as affected b y saturation and non-linearit y is still remaining, and we still do not ha v e a
“complete” PSF! In this section, w e will use the same steps b efore to mak e stac ks of more
inner regions of the PSF to ultimately unite them all in to a single PSF in Section 2.3.6
[Uniting the differen t PSF comp onen ts], page 114.
F or the outer PSF, w e selected stars in the magnitude range of 6 to 10. So let’s hav e
a lo ok and see ho w many stars w e ha v e in the magnitude range of 12 to 13 with a more
relaxed condition on the minim um distance for neigh b ors, --mindistdeg=0.01 (36 arcsec,
since these stars are fain ter), and use the ds9 region script to visually insp ect them:
$ mkdir inner
$ astscript-psf-select-stars flat/67510.fits \
--magnituderange=12,13 --mindistdeg=0.01 \
--output=inner/67510-12-13.fits
$ astscript-ds9-region inner/67510-12-13.fits -cra,dec \
--namecol=phot_g_mean_mag \
--command="ds9 flat/67510.fits -zoom to fit -zscale"
W e ha v e 41 stars, but if y ou zo om in to their cen ters, y ou will see that they do not ha v e
an y ma jor bleeding-v ertical saturation any more. Only the v ery cen tral core of some of the
stars is saturated. W e can therefore use these stars to fill the strong bleeding fo otprin ts
that w ere presen t in the outer stac k of outer/stack.fits . Similar to b efore, let’s build
ready-to-stac k crops of these stars. T o get a b etter feeling of the normalization radii, follo w
the same steps of Section 2.3.4 [Building outer part of PSF], page 110, (setting --tmpdir
and --keeptmp ). In this case, since the stars are fain ter, we can set a smaller size for the
individual stamps, --widthinpix=500,500 , to sp eed up the calculations:
$ counter=1
$ mkdir inner/stamps
Chapter 2: T utorials 114
$ asttable inner/67510-12-13.fits \
| while read -r ra dec mag; do
astscript-psf-stamp label/67510-seg.fits \
--mode=wcs \
--normradii=5,10 \
--center=$ra,$dec \
--widthinpix=500,500 \
--segment=label/67510-seg.fits \
--output=inner/stamps/67510-$counter.fits; \
counter=$((counter+1)); \
done
$ numimgs=$(echo inner/stamps/*.fits | wc -w)
$ astarithmetic inner/stamps/*.fits $numimgs 3 0.2 sigclip-mean \
-g1 --output=inner/stack.fits --wcsfile=none
$ astscript-fits-view inner/stack.fits inner/stamps/*.fits
W e are no w ready to unite the t w o stac ks w e ha v e constructed: the outer and the inner
parts.
2.3.6 Uniting the differen t PSF comp onen ts
In Section 2.3.4 [Building outer part of PSF], page 110, w e built the outer part of the
extended PSF and the inner part w as built in Section 2.3.5 [Inner part of the PSF], page 113.
The outer part w as built with v ery brigh t stars, and the inner part using fain ter stars to
not ha v e saturation in the core of the PSF. The next step is to join these t w o parts in order
to ha v e a single PSF. First of all, let’s ha v e a lo ok at the t w o stac ks and also to their radial
profiles to ha v e a go o d feeling of the task. Note that y ou will need to ha v e TOPCA T to
run the last command and plot the radial profile (see Section A.2 [TOPCA T], page 965).
$ astscript-fits-view outer/stack.fits inner/stack.fits
$ astscript-radial-profile outer/stack.fits -o outer/profile.fits
$ astscript-radial-profile inner/stack.fits -o inner/profile.fits
$ astscript-fits-view outer/profile.fits inner/profile.fits
F rom the visual insp ection of the images and the radial profiles, it is clear that w e ha v e
saturation in the cen ter for the outer part. Note that the absolute flux v alues of the PSFs
are meaningless since they dep end on the normalization radii w e used to obtain them. The
uniting step consists in scaling up (or do wn) the inner part of the PSF to ha v e the same
flux at the junction radius, and then, use that flux-scaled inner part to fill the cen ter of the
outer PSF. T o get a feeling of the pro cess, first, let’s op en the t w o radial profiles and find
the factor man ually first:
1. Run this command to op en the t w o tables in Section A.2 [TOPCA T], page 965:
$ astscript-fits-view outer/profile.fits inner/profile.fits
2. On the left side of the screen, under “T able List”, y ou will see the t w o imp orted tables.
Clic k on the first one (profile of the outer part) so it is sho wn first.
3. Under the “Graphics” men u item, clic k on “Plane plot”. A new windo w will op en
with the plot of the first t w o columns: RADIUS on the horizontal axis and MEAN on the
v ertical. The rest of the steps are done in this windo w.
Chapter 2: T utorials 115
4. In the b ottom settings, within the left panel, clic k on the “Axes” item. This will allo w
customization of the plot axes.
5. In the b ottom-righ t panel, clic k on the b ox in fron t of “Y Log” to mak e the v ertical
axis logarithmic-scaled.
6. On the “La y ers” men u, select “Add P osition Con trol” to allo w adding the profile of
the inner region. After it, y ou will see that a new red-blue scatter plot icon op ened on
the b ottom-left men u (with a title of <no table> ).
7. On the b ottom-righ t panel, in the drop-do wn menu in fron t of Table: , select 2:
profile.fits . Afterwards, y ou will see the radial profile of the inner stac k as the
newly added blue plot. Our goal here is to find the factor that is necessary to m ultiply
with the inner profile so it matc hes the outer one.
8. On the b ottom-righ t panel, in fron t of Y: , you will see MEAN . Clic k in the white-space
after it, and t yp e this: *100 . This will displa y the MEAN column of the inner profile, after
m ultiplying it b y 100. Afterw ards, y ou will see that the inner profile (blue) matc hes
more cleanly with the outer (red); esp ecially in the smaller radii. A t larger radii, it
do es not drop lik e the red plot. This is b ecause of the extremely lo w signal-to-noise
ratio at those regions in the fain ter stars used to mak e this stac k.
9. T ak e y our mouse cursor o v er the profile, in particular ov er the bump around a radius
of 100 pixels. Scroll y our mouse do wn-w ard to zo om-in to the profile and up-w ard to
zo om-out. Y ou can clic k-and-hold an y part of the profile and if y ou mo v e y our cursor
(while still holding the mouse-button) to lo ok at differen t parts of the profile. This is
particular helpful when y ou ha v e zo omed-in to the profile.
10. Zo om-in to the bump around a radius of 100 pixels un til the horizon tal axis range
b ecomes around 50 to 130 pixels.
11. Y ou clearly see that the inner stac k (blue) is m uc h more noisy than the outer (red)
stac k. By “noisy”, w e mean that the scatter of the p oin ts is m uc h larger. If y ou further
zo om-out, y ou will see that the shallo w slop e at the larger radii of the inner (blue)
profile has also affected the heigh t of this bump in the inner profile. This is a very
imp ortant p oin t: this clearly sho ws that the inner profile is to o noisy at these radii.
12. Clic k-and-hold y our mouse to see the inner parts of the t w o profiles (in the range 0 to
80). Y ou will see that for radii less than 40 pixels, the inner profile (blue) p oin ts lo ose
their scatter (and th us ha v e a go o d signal-to-noise ratio).
13. Zo om-in to the plot and follo w the profiles un til smaller radii (for example, 10 pixels).
Y ou see that for eac h radius, the inner (blue) p oin ts are consisten tly ab o v e the outer
(red) p oin ts. This sho ws that the × 100 factor w e selected ab o v e w as to o m uc h.
14. In the b ottom-righ t panel, change the 100 to 80 and zo om-in to the same region. A t
eac h radius, the blue p oin ts are no w b elo w the red p oin ts, so the scale-factor 80 is not
enough. So let’s increase it and try 90 . After zo oming-in, you will notice that in the
inner-radii (less than 30 pixels), they are no w very similar. The ultimate aim of the
steps b elo w is to find this factor automatically .
15. But b efore con tin uing, let’s fo cus on another imp ortan t p oin t ab out the cen tral regions:
non-linearit y and saturation. While y ou are zo omed-in (from the step ab o v e), follo w
(clic k-and-drag) the profile to w ards smaller radii. Y ou will see that smaller than a
radius of 10, they start to div erge. But this time, the outer (red) profile is getting a
shallo w er slop e and div erges significan tly from ab out the radius of 8. W e had mask ed all
Chapter 2: T utorials 116
saturated pixels b efore, so this div ergence for radii smaller than 10 sho ws the effect of
the CCD’s non-linearit y (where the n um b er of electrons will not b e linearly correlated
with the n um b er of inciden t photons). This is presen t in all CCDs and pixels b eyond
this lev el should not b e used in measuremen ts (or prop erly corrected).
The items ab o ve w ere only listed so y ou get a go o d mental/visual understanding of the
logic b ehind the op eration of the next script (and to learn ho w to tune its parameters where
necessary): astscript-psf-scale-factor . This script is more general than this particular
problem, but can b e used for this sp ecial case also. Its job is to tak e a mo del of an ob ject
(PSF, or inner stac k in this case) and the p osition of an instance of that mo del (a star, or
the outer stac k in this case) in a larger image.
Instead of dealing with radial profiles (that enforce a certain shap e), this script will put
the cen ters of the inner and outer PSFs o ver eac h other and divide the outer b y the inner.
Let’s ha v e a lo ok with the command b elo w. Just note that w e are running it with --keeptmp
so the temp orary directory with all the in termediate files remain for further clarification:
$ astscript-psf-scale-factor outer/stack.fits \
--psf=inner/stack.fits --center=501,501 \
--mode=img --normradii=10,15 --keeptmp
$ astscript-fits-view stack_psfmodelscalefactor/cropped-*.fits \
stack_psfmodelscalefactor/for-factor-*.fits
With the second command, y ou see the four steps of the pro cess: the first t w o images
sho w the cropp ed outer and inner stac ks (to same width image). The third sho ws the radial
p osition of eac h pixel (which is used to only k eep the pixels within the desired radial range).
The fourth sho ws the p er-pixel division of the outer b y the inner within the requested radii.
The sigma-clipp ed median of these pixels is finally rep orted. Unlike the radial profile metho d
(whic h a v erages o v er a circular/elliptical ann ulus for eac h radius), this metho d imp oses no
a-priori shap e on the PSF. This mak es it v ery useful for complex PSFs (lik e the case here).
T o con tin ue, let’s remov e the temp orary directory and re-run the script but with --quiet
mo de so w e can put the output in a shell v ariable.
$ rm -r stack_psfmodelscalefactor
$ scale=$(astscript-psf-scale-factor outer/stack.fits \
--psf=inner/stack.fits --center=501,501 \
--mode=img --normradii=10,15 --quiet)
$ echo $scale
No w that w e kno w the scaling factor, w e are ready to unite the outer and the inner part
of the PSF. T o do that, we will use the script astscript-psf-unite with the command
b elo w (for more on this script, see Section 10.8.4 [In v oking astscript-psf-unite], page 712).
The basic parameters are the inner part of the PSF (giv en to --inner ), the inner part’s
scale factor ( --scale ), and the junction radius ( --radius ). The inner part is first scaled,
and all the pixels of the outer image within the giv en radius are replaced with the pixels of
the inner image. Since the flux factor w as computed for a ring of pixels b et w een 10 and 15
pixels, let’s set the junction radius to b e 12 pixels (roughly in b et w een 10 and 15):
$ astscript-psf-unite outer/stack.fits \
--inner=inner/stack.fits --radius=12 \
--scale=$scale --output=psf.fits
Chapter 2: T utorials 117
Let’s ha v e a lo ok at the outer stac k and the final PSF with the command b elo w. Since
w e w an t sev eral other DS9 settings to help y ou directly see the main p oin t, we are using
--ds9extra . After DS9 is op ened, y ou can see that the cen ter of the PSF has no w b een
nicely filled. Y ou can clic k on the “Edit” button and then the “Colorbar” and hold y our
cursor o v er the image and mo v e it. Y ou can see that b esides filling the inner regions nicely ,
there is also no ma jor discon tin uit y in the 2D image around the union radius of 12 pixels
around the cen ter.
$ astscript-fits-view outer/stack.fits psf.fits --ds9scale=minmax \
--ds9extra="-scale limits 0 22000 -match scale" \
--ds9extra="-lock scale yes -zoom 4 -scale log"
Nothing demonstrates the effect of a bad analysis than actually seeing a bad result! So
let’s c ho ose a bad normalization radial range (50 to 60 pixels) and unite the inner and outer
parts based on that. The last command will op en the t w o PSFs together in DS9, y ou should
b e able to immediately see the discon tin uit y in the union radius.
$ scale=$(astscript-psf-scale-factor outer/stack.fits \
--psf=inner/stack.fits --center=501,501 \
--mode=img --normradii=50,60 --quiet)
$ astscript-psf-unite outer/stack.fits \
--inner=inner/stack.fits --radius=55 \
--scale=$scale --output=psf-bad.fits
$ astscript-fits-view psf-bad.fits psf.fits --ds9scale=minmax \
--ds9extra="-scale limits 0 50 -match scale" \
--ds9extra="-lock scale yes -zoom 4 -scale log"
As y ou see, the selection of the normalization radii and the unite radius are v ery imp or-
tan t. The first time y ou are trying to build the PSF of a new dataset, it has to b e explored
with a visual insp ection of the images and radial profiles. Once y ou ha v e found a go o d
normalization radius for a certain part of the PSF in a surv ey , y ou can generally use it
comfortably without c hange. But for a new surv ey , or a differen t part of the PSF, b e sure
to rep eat the visual c hecks abov e to c ho ose the b est radii. As a summary , a go o d junction
radius is one that:
• Is large enough to not let saturation and non-linearit y (from the outer profile) in to the
inner region.
• Is small enough to ha v e a sufficien tly high signal to noise ratio (from the inner profile)
to a v oid adding noise in the union radius.
No w that the complete PSF has b een obtained, let’s remo v e that bad-lo oking PSF, and
stic k with the nice and clean PSF for the next step in Section 2.3.7 [Subtracting the PSF],
page 117.
$ rm -rf psf-bad.fits
2.3.7 Subtracting the PSF
Previously (in Section 2.3.6 [Uniting the differen t PSF comp onents], page 114) w e con-
structed a full PSF, from the cen tral pixel to a radius of 500 pixels. No w, let’s use the PSF
to subtract the scattered ligh t from eac h individual star in the image.
Chapter 2: T utorials 118
By construction, the pixel v alues of the PSF came from the normalization of the indi-
vidual stamps (that w ere created for stars of differen t magnitudes). As a consequence, it is
necessary to compute a scale factor to fit that PSF image to eac h star. This is done with
the same astscript-psf-scale-factor command that w e used previously in Section 2.3.6
[Uniting the differen t PSF comp onen ts], page 114. The difference is that no w w e are not
aiming to join t w o differen t PSF parts but lo oking for the necessary scale factor to matc h
the star with the PSF. Afterw ards, w e will use astscript-psf-subtract for placing the
PSF image at the desired co ordinates within the same pixel grid as the image. Finally , once
the stars ha ve been mo deled b y the PSF, w e will subtract it.
First, let’s start with a single star. Later, when the basic idea has b een explained, we
will generalize the metho d for an y n um b er of stars. With the follo wing command w e obtain
the co ordinates (RA and DEC) and magnitude of the brigh test star in the image (whic h is
on the top edge of the image):
$ mkdir single-star
$ center=$(asttable flat/67510-bright.fits --sort phot_g_mean_mag \
--column=ra,dec --head 1 \
| awk '{printf "%s,%s", $1, $2}')
$ echo $center
With the cen ter p osition of that star, let’s obtain the flux factor using the same normal-
ization ring w e used for the creation of the outer part of the PSF:
$ scale=$(astscript-psf-scale-factor label/67510-seg.fits \
--mode=wcs --quiet \
--psf=psf.fits \
--center=$center \
--normradii=10,15 \
--segment=label/67510-seg.fits)
No w w e ha v e all the information necessary to mo del the star using the PSF: the p osition
on the sky and the flux factor. Let’s use this data with the script astscript-psf-subtract
for mo deling this star and ha v e a lo ok with DS9.
$ astscript-psf-subtract label/67510-seg.fits \
--mode=wcs \
--psf=psf.fits \
--scale=$scale \
--center=$center \
--output=single-star/subtracted.fits
$ astscript-fits-view label/67510-seg.fits single-star/subtracted.fits \
--ds9center=$center --ds9mode=wcs --ds9extra="-zoom 4"
Y ou will notice that there is something wrong with this “subtraction”! The b o x of the
extended PSF is clearly visible! The sky noise under the b o x is clearly larger than the rest
of the noise in the image. Before reading on, please try to think ab out the cause of this
y ourself.
T o understand the cause, let’s lo ok at the scale factor, the n um b er of stamps used to
build the outer part (and its square ro ot):
$ echo $scale
Chapter 2: T utorials 119
$ ls outer/stamps/*.fits | wc -l
$ ls outer/stamps/*.fits | wc -l | awk '{print sqrt($1)}'
Y ou see that the scale is almost 19! As a result, the PSF has b een multiplied b y 19 b efore
b eing subtracted. How ev er, the outer part of the PSF w as created with only a handful of
star stamps. When y ou stac k N images, the stac k’s signal-to-noise ratio (S/N) impro v es b y
√ N . W e had 8 images for the outer part, so the S/N has only impro v ed b y a factor of just
under 3! When w e m ultiply the final stac k ed PSF with 19, w e are also scaling up the noise
b y that same factor (most imp ortan tly: in the outer most regions where there is almost no
signal). So the stac ked image’s noise-lev el is 19 / 3=6 . 3 times larger than the noise of the
input image. This terrible noise-lev el is what y ou clearly see as the fo otprin t of the PSF.
T o confirm this, let’s use the commands b elo w to subtract the fain test of the brigh t-stars
catalog (note the use of --tail when finding the cen tral p osition). Y ou will notice that the
scale factor ( ∼ 1 . 3) is no w smaller than 3. So when w e multiply the PSF with th is factor,
the PSF’s noise lev el is low er than our input image and w e should not see any fo otprin t lik e
b efore. Note also that w e are using a larger zo om factor, b ecause this star is smaller in the
image.
$ center=$(asttable flat/67510-bright.fits --sort phot_g_mean_mag \
--column=ra,dec --tail 1 \
| awk '{printf "%s,%s", $1, $2}')
$ scale=$(astscript-psf-scale-factor label/67510-seg.fits \
--mode=wcs --quiet \
--psf=psf.fits \
--center=$center \
--normradii=10,15 \
--segment=label/67510-seg.fits)
$ echo $scale
$ astscript-psf-subtract label/67510-seg.fits \
--mode=wcs \
--psf=psf.fits \
--scale=$scale \
--center=$center \
--output=single-star/subtracted.fits
$ astscript-fits-view label/67510-seg.fits single-star/subtracted.fits \
--ds9center=$center --ds9mode=wcs --ds9extra="-zoom 10"
In a large surv ey like J-PLUS, it is easy to use more and more brigh t stars from differen t
p oin tings (ideally with similar FWHM and similar telescop e prop erties 41 ) to impro v e the
S/N of the PSF. As explained b efore, w e designed the output files of this tutorial with the
67510 (whic h is this image’s p oin ting lab el in J-PLUS) where necessary so y ou see ho w easy
it is to add more p oin tings to use in the creation of the PSF.
41 F or example, in J-PLUS, the baffle of the secondary mirror was adjusted in 2017 because it pro duced
extra spik es in the PSF. So all images after that date ha v e a PSF with 4 spikes (like this one), while
those b efore it ha ve man y more spikes.
Chapter 2: T utorials 120
Let’s consider no w more than one single star. W e should ha v e t w o things in mind:
• The brigh test (subtract-able, see the p oint below) star should b e the first star to b e
subtracted. This is b ecause of its extended wings whic h ma y affect the scale factor of
nearb y stars. So w e should sort the catalog b y magnitude and come do wn from the
brigh test.
• W e should only subtract stars where the scale factor is less than the S/N of the PSF
(in relation to the data).
Since it can get a little complex, it is easier to implemen t this step as a script (that is
hea vily commen ted for y ou to easily understand ev e ry step; esp ecially if you put it in a
go o d text editor with color-co ding!). Y ou will notice that script also creates a .log file,
whic h sho ws whic h star w as subtracted and whic h one w as not (this is imp ortan t, and will
b e used b elo w!).
#!/bin/bash
# Abort the script on first error.
set -e
# ID of image to subtract stars from.
imageid=67510
# Get S/N level of the final PSF in relation to the actual data:
snlevel=$(ls outer/stamps/*.fits | wc -l | awk '{print sqrt($1)}')
# Put a copy of the image we want to subtract the PSF from in the
# final file (this will be over-written after each subtraction).
subtracted=subtracted/$imageid.fits
cp label/$imageid-seg.fits $subtracted
# Name of log-file to keep status of the subtraction of each star.
logname=subtracted/$imageid.log
echo "# Column 1: RA [deg, f64] Right ascension of star." > $logname
echo "# Column 2: Dec [deg, f64] Declination of star." >> $logname
echo "# Column 3: Stat [deg, f64] Status (1: subtracted)" >> $logname
# Go over each item in the bright star catalog:
asttable flat/67510-bright.fits -cra,dec --sort phot_g_mean_mag \
| while read -r ra dec; do
# Put a comma between the RA/Dec to pass to options.
center=$(echo $ra $dec | awk '{printf "%s,%s", $1, $2}')
# Calculate the scale value
scale=$(astscript-psf-scale-factor label/67510-seg.fits \
--mode=wcs --quiet\
--psf=psf.fits \
--center=$center \
Chapter 2: T utorials 121
--normradii=10,15 \
--segment=label/67510-seg.fits)
# Subtract this star if the scale factor is less than the S/N
# level calculated above.
check=$(echo $snlevel $scale \
| awk '{if($1>$2) c="good"; else c="bad"; print c}')
if [ $check = good ]; then
# A temporary file to subtract this star.
subtmp=subtracted/$imageid-tmp.fits
# Subtract this star from the image where all previous stars
# were subtracted.
astscript-psf-subtract $subtracted \
--mode=wcs \
--psf=psf.fits \
--scale=$scale \
--center=$center \
--output=$subtmp
# Rename the temporary subtracted file to the final one:
mv $subtmp $subtracted
# Keep the status for this star.
status=1
else
# Let the user know this star did not work, and keep the status
# for this star.
echo "$center: $scale is larger than $snlevel"
status=0
fi
# Keep the status in a log file.
echo "$ra $dec $status" >> $logname
done
Cop y the con tents abov e in to a file called subtract-psf-from-cat.sh and run the
follo wing commands. Just note that in the script ab o v e, w e assumed the output is written
in the subtracted/ , directory , so w e will first mak e that.
$ mkdir subtracted
$ chmod +x subtract-psf-from-cat.sh
$ ./subtract-psf-from-cat.sh
$ astscript-fits-view label/67510-seg.fits subtracted/67510.fits
Can y ou visually find the stars that ha ve been subtracted? Its a little hard, is not it?
This sho ws that y ou done a go o d job this time (the sky-noise is not significan tly affected)!
Chapter 2: T utorials 128
>> cat.txt
$ echo "0 0.000 0.000 2 5 4.7 0.0 1.0 30.0 5.0" >> cat.txt
$ echo "1 250.0 250.0 1 40 1.0 -25 0.4 3.44 5.0" >> cat.txt
$ echo "2 50.00 50.00 4 0 0.0 0.0 0.0 6.00 0.0" >> cat.txt
$ echo "3 450.0 50.00 4 0 0.0 0.0 0.0 6.50 0.0" >> cat.txt
$ echo "4 50.00 450.0 4 0 0.0 0.0 0.0 7.00 0.0" >> cat.txt
$ echo "5 450.0 450.0 4 0 0.0 0.0 0.0 7.50 0.0" >> cat.txt
T o mak e sure if the catalog’s con ten t is correct (and there w as no t yp o for example!), Sufi
runs ‘ cat cat.txt ’, and confirms that it is correct.
No w that the catalog is created, Sufi is ready to call Mak eProfiles to build the image
con taining these ob jects. He lo oks in to his records and finds that the zero p oin t magnitude
for that nigh t, and that particular detector, w as 18 magnitudes. The studen t w as a little
confused on the concept of zero p oin t, so Sufi p oin ted him to Section 7.4.2 [Brigh tness,
Flux, Magnitude and Surface brigh tness], page 576, whic h the studen t can study in detail
later. Sufi therefore runs Mak eProfiles with the command b elo w:
$ astmkprof --prepforconv --mergedsize=499,499 --zeropoint=18.0 cat.txt
MakeProfiles 0.23 started on Sat Oct 6 16:26:56 953
- 6 profiles read from cat.txt
- Random number generator (RNG) type: ranlxs1
- Basic RNG seed: 1652884540
- Using 12 threads.
---- row 3 complete, 5 left to go
---- row 4 complete, 4 left to go
---- row 6 complete, 3 left to go
---- row 5 complete, 2 left to go
---- ./0_cat_profiles.fits created.
---- row 1 complete, 1 left to go
---- row 2 complete, 0 left to go
- ./cat_profiles.fits created. 0.092573 seconds
-- Output: ./cat_profiles.fits
MakeProfiles finished in 0.293644 seconds
Sufi encourages the studen t to read through the prin ted output. As the statemen ts sa y ,
t w o FITS files should ha v e b een created in the running directory . So Sufi ran the command
b elo w to confirm:
$ ls
0_cat_profiles.fits cat_profiles.fits cat.txt
The file 0_cat_profiles.fits is the PSF Sufi had ask ed for, and cat_profiles.fits is
the image con taining the main ob jects in the catalog. Sufi op ened the main image with the
command b elo w (using SA O DS9):
$ astscript-fits-view cat_profiles.fits --ds9scale=95
The studen t could clearly see the main elliptical structure in the cen ter. Ho w ev er, the
size of cat_profiles.fits w as surprising for the studen t, instead of 499 b y 499 (as w e had
requested), it w as 2615 b y 2615 pixels (from the command b elo w):
$ astfits cat_profiles.fits
Fits (GNU Astronomy Utilities) 0.23
Chapter 2: T utorials 129
Run on Sat Oct 6 16:26:58 953
-----
HDU (extension) information: 'cat_profiles.fits'.
Column 1: Index (counting from 0, usable with '--hdu').
Column 2: Name ('EXTNAME' in FITS standard, usable with '--hdu').
Column 3: Image data type or 'table' format (ASCII or binary).
Column 4: Size of data in HDU.
-----
0 MKPROF-CONFIG no-data 0
1 Mock profiles float32 2615x2615
So Sufi explained wh y ov ersampling is imp ortan t in mo deling, esp ecially for parts of the
image where the flux c hange is significan t o v er a pixel. Recall that when y ou o v ersample
the mo del (for example, b y 5 times), for ev ery desired pixel, y ou get 25 pixels (5 × 5). Sufi
then explained that after con v olving (next step b elo w) w e will do wn-sample the image to
get our originally desired size/resolution.
After seeing the image, the studen t complained that only the large elliptical mo del for
the Andromeda nebula can b e seen in the cen ter. He could not see the four stars that w e had
also requested in the catalog. So Sufi had to explain that the stars are there in the image,
but the reason that they are not visible when lo oking at the whole image at once, is that
they only co v er a single pixel! T o prov e it, he cen tered the image around the co ordinates
2308 and 2308, where one of the stars is lo cated in the o v er-sampled image [y ou can do this
in ds9 b y selecting “P an” in the “Edit” men u, then clic king around that p osition]. Sufi
then zo omed in to that region and so on, the star’s non-zero pixel could b e clearly seen.
Sufi explained that the stars will tak e the shap e of the PSF (cov er an area of more than
one pixel) after con v olution. If w e did not ha v e an atmosphere and w e did not need an
ap erture, then stars w ould only co v er a single pixel with normal CCD resolutions. So Sufi
con v olv ed the image with this command:
$ astconvolve --kernel=0_cat_profiles.fits cat_profiles.fits \
--output=cat_convolved.fits
Convolve started on Sat Oct 6 16:35:32 953
- Using 8 CPU threads.
- Input: cat_profiles.fits (hdu: 1)
- Kernel: 0_cat_profiles.fits (hdu: 1)
- Input and Kernel images padded. 0.075541 seconds
- Images converted to frequency domain. 6.728407 seconds
- Multiplied in the frequency domain. 0.040659 seconds
- Converted back to the spatial domain. 3.465344 seconds
- Padded parts removed. 0.016767 seconds
- Output: cat_convolved.fits
Convolve finished in: 10.422161 seconds
When con v olution finished, Sufi op ened cat_convolved.fits and the four stars could b e
easily seen no w:
$ astscript-fits-view cat_convolved.fits --ds9scale=95
It w as in teresting for the studen t that all the flux in that single pixel is no w distributed
o v er so man y pixels (the sum of all the pixels in eac h conv olv ed star is actually equal to
Chapter 2: T utorials 130
the v alue of the single pixel b efore con v olution). Sufi explained ho w a PSF with a larger
FWHM w ould mak e the p oin ts ev en wider than this (distributing their flux in a larger area).
With the con v olv ed image ready , they w ere prepared to resample it to the original pixel
scale Sufi had planned [from the $ astmkprof -P command ab o ve, recall that Mak eProfiles
had o v er-sampled the image b y 5 times]. Sufi explained the basic concepts of w arping the
image to his studen t and ran W arp with the follo wing command:
$ astwarp --scale=1/5 --centeroncorner cat_convolved.fits
Warp started on Sat Oct 6 16:51:59 953
Using 8 CPU threads.
Input: cat_convolved.fits (hdu: 1)
matrix:
0.2000 0.0000 0.4000
0.0000 0.2000 0.4000
0.0000 0.0000 1.0000
$ astfits cat_convolved_scaled.fits --quiet
0 WARP-CONFIG no-data 0
1 Warped float32 523x523
cat_convolved_scaled.fits no w has the correct pixel scale. Ho wev er, the image is still
larger than what w e had w an ted, it is 523 × 523 pixels (not our desired 499 × 499). The
studen t is sligh tly confused, so Sufi also resamples the PSF with the same scale b y running
$ astwarp --scale=1/5 --centeroncorner 0_cat_profiles.fits
$ astfits 0_cat_profiles_scaled.fits --quiet
0 WARP-CONFIG no-data 0
1 Warped float32 25x25
Sufi notes that 25 = 12 + 12 + 1 and that 523 = 499 + 12 + 12. He go es on to explain that
frequency space con volution will dim the edges and that is wh y he added the --prepforconv
option to Mak eProfiles ab ov e. No w that con v olution is done, Sufi can remo v e those extra
pixels using Crop with the command b elo w. Crop’s --section option accepts co ordinates
inclusiv ely and coun ting from 1 (according to the FITS standard), so the crop region’s first
pixel has to b e 13, not 12.
$ astcrop cat_convolved_scaled.fits --section=13:*-12,13:*-12 \
--mode=img --zeroisnotblank
Crop started on Sat Oct 6 17:03:24 953
- Read metadata of 1 image. 0.001304 seconds
---- ...nvolved_scaled_cropped.fits created: 1 input.
Crop finished in: 0.027204 seconds
T o fully con vince the studen t, Sufi chec ks the size of the output of the crop command ab o v e:
$ astfits cat_convolved_scaled_cropped.fits --quiet
0 n/a no-data 0
1 n/a float32 499x499
Finally , cat_convolved_scaled_cropped.fits is 499 × 499 pixels and the mo c k
Andromeda galaxy is cen tered on the cen tral pixel. This is the same dimensions as Sufi
had desired in the b eginning. All this trouble w as certainly w orth it b ecause no w there is
no dimming on the edges of the image and the profile cen ters are more accurately sampled.
Chapter 2: T utorials 131
The final step to sim ulate a real observ ation w ould b e to add noise to the image. Sufi
set the zero p oin t magnitude to the same v alue that he set when making the mo c k profiles
and lo oking again at his observ ation log, he had measured the bac kground flux near the
nebula had a p er-pixel magnitude of 7 that nigh t. F or more on ho w the bac kground v alue
determines the noise, see Section 6.2.3 [Noise basics], page 402. So using these v alues
he ran Arithmetic’s mknoise-sigma-from-mean op erator, and with the second command,
he visually insp ected the image. The mknoise-sigma-from-mean op erator tak es the noise
standard deviation in linear units, not magnitudes (whic h are logarithmic). Therefore within
the same Arithmetic command, he has con v erted the sky bac kground magnitude to coun ts
using Arithmetic’s counts-to-mag op erator.
$ astarithmetic cat_convolved_scaled_cropped.fits \
7 18 mag-to-counts mknoise-sigma-from-mean \
--output=out.fits
$ astscript-fits-view out.fits
The out.fits file no w con tains the noised image of the mo c k catalog Sufi had ask ed for.
The studen t had not observ ed the nebula in the sky , so when he sa w the mo c k image in
SA O DS9 (with the second command ab o ve), he understo od why Sufi w as dubious: it was
v ery diffuse!
Seeing ho w the --output option allo ws the user to sp ecify the name of the output file,
the studen t w as confused and w an ted to kno w wh y Sufi had not used it more regularly
b efore? Sufi explained that for in termediate steps, y ou can rely on the automatic output
of the programs (see Section 4.9 [Automatic output], page 289). Doing so will giv e all the
in termediate files a similar basic name structure, so in the end y ou can simply remo v e them
all with the Shell’s capabilities, and it will b e familiar for other users. So Sufi decided to
sho w this to the studen t b y making a shell script from the commands he had used b efore.
The command-line shell has the capabilit y to read all the separate input commands from
a file. This is useful when y ou w an t to do the same thing m ultiple times, with only the
names of the files or minor parameters c hanging b et w een the differen t instances. Using the
shell’s history (b y pressing the up k eyb oard k ey) Sufi review ed all the commands and then
he retriev ed the last 5 commands with the $ history 5 command. He selected all those
lines he had input and put them in a text file named mymock.sh . Then he defined the
edge and base shell v ariables for easier customization later, and b efore ev ery command, he
added some commen ts (lines starting with # ) for future readabilit y . Finally , Sufi p oin ted
the studen t to Gn uastro’s Section 2.1 [General program usage tutorial], page 22, whic h has
a full section on Section 2.1.22 [W riting scripts to automate the steps], page 73.
#!/bin/bash
edge=12
base=cat
# Stop running next commands if one fails.
set -e
# Remove any (possibly) existing output (from previous runs)
# before starting.
Chapter 2: T utorials 132
rm -f out.fits
# Run MakeProfiles to create an oversampled FITS image.
astmkprof --prepforconv --mergedsize=499,499 --zeropoint=18.0 \
"$base".txt
# Convolve the created image with the kernel.
astconvolve "$base"_profiles.fits \
--kernel=0_"$base"_profiles.fits \
--output="$base"_convolved.fits
# Scale the image back to the intended resolution.
astwarp --scale=1/5 --centeroncorner "$base"_convolved.fits
# Crop the edges out (dimmed during convolution). '--section'
# accepts inclusive coordinates, so the start of the section
# must be one pixel larger than its end.
st_edge=$(( edge + 1 ))
astcrop "$base"_convolved_scaled.fits --zeroisnotblank \
--mode=img --section=$st_edge:*-$edge,$st_edge:*-$edge
# Add noise to the image.
$ astarithmetic "$base"_convolved_scaled_cropped.fits \
7 18 mag-to-counts mknoise-sigma-from-mean \
--output=out.fits
# Remove all the temporary files.
rm 0*.fits "$base"*.fits
He used this c hance to remind the studen t of the imp ortance of commen ts in co de or
shell scripts! Just like metadata in a dataset, when writing the co de, y ou ha v e a go o d
men tal picture of what y ou are doing, so writing commen ts migh t seem sup erfluous and
excessiv e. Ho w ev er, in one mon th when y ou w an t to re-use the script, you ha v e lost that
men tal picture and remem b ering it can b e time-consuming and frustrating. The imp ortance
of commen ts is further amplified when y ou w an t to share the script with a friend/colleague.
So it is go o d to accompan y an y step of a script, or co de, with useful commen ts while y ou
are writing it (create a go o d men tal picture of wh y y ou are doing something: do not just
describ e the command, but its purp ose).
Sufi then explained to the eager studen t that y ou define a v ariable b y giving it a name,
follo w ed b y an = sign and the v alue y ou w an t. Then y ou can reference that v ariable from
an ywhere in the script b y calling its name with a $ prefix. So in the script whenev er y ou
see $base , the v alue w e defined for it ab o v e is used. If y ou use adv anced editors lik e GNU
Emacs or ev en simpler ones lik e Gedit (part of the GNOME graphical user in terface) the
v ariables will b ecome a differen t color whic h can really help in understanding the script.
W e ha v e put all the $base v ariables in double quotation marks ( " ) so the v ariable name
and the follo wing text do not get mixed, the shell is going to ignore the " after replacing
the v ariable v alue. T o mak e the script executable, Sufi ran the follo wing command:
Chapter 2: T utorials 133
$ chmod +x mymock.sh
Then finally , Sufi ran the script, simply b y calling its file name:
$ ./mymock.sh
After the script finished, the only file remaining is the out.fits file that Sufi had w an ted
in the b eginning. Sufi then explained to the studen t how he could run this script an ywhere
that he has a catalog if the script is in the same directory . The only thing the studen t had
to mo dify in the script w as the name of the catalog (the v alue of the base v ariable in the
start of the script) and the v alue to the edge v ariable if he c hanged the PSF size. The
studen t w as also happ y to hear that he will not need to mak e it executable again when he
mak es c hanges later, it will remain executable unless he explicitly c hanges the executable
flag with chmod .
The studen t w as really excited, since no w, through simple shell scripting, he could really
sp eed up his w ork and run any command in an y fashion he lik es allo wing him to b e m uch
more creativ e in his w orks. Un til no w he w as using the graphical user in terface whic h do es
not ha v e suc h a facilit y and doing rep etitiv e things on it w as really frustrating and some
times he w ould mak e mistak es. So he left to go and try scripting on his o wn computer.
He later reminded Sufi that the second tutorial in the Gn uastro b o ok as more complex
commands in data analysis, and a more adv anced in tro duction to scripting (see Section 2.1
[General program usage tutorial], page 22).
Sufi could no w get bac k to his o wn w ork and see if the sim ulated nebula whic h re-
sem bled the one in the Andromeda constellation could b e detected or not. Although
it w as extremely fain t 43 . Therefore, Sufi ran Gn uastro’s detection soft w are (Section 7.2
[NoiseChisel], page 544) to see if this ob ject is detectable or not. NoiseChisel’s output
( out_detected.fits ) is a m ulti-extension FITS file, so he used Gn uastro’s astscript-
fits-view program in the second command to see the output:
$ astnoisechisel out.fits
$ astscript-fits-view out_detected.fits
In the “Cub e” windo w (that w as op ened with DS9), if Sufi click ed on the “Next” button
to see the pixels that w ere detected to con tain significan t signal. F ortunately the nebula’s
shap e w as detectable and he could finally confirm that the nebula he k ept in his noteb o ok
w as actually observ able. He wrote this result in the draft man uscript that w ould later
b ecome “Bo ok of fixed stars” 44 .
He still had to c hec k the other nebula he sa w from Y emen and sev eral other such ob jects,
but they could w ait un til tomorro w (thanks to the shell script, he only has to define a new
catalog). It w as nearly sunset and they had to b egin preparing for the nigh t’s measuremen ts
on the ecliptic.
43 The brigh tness of a diffuse ob ject is added o v er all its pixels to giv e its final magnitude, see Section 7.4.2
[Brigh tness, Flux, Magnitude and Surface brigh tness], page 576. So although the magnitude 3.44 (of the
mo c k nebula) is orders of magnitude brigh ter than 6 (of the stars), the central galaxy is m uc h fainter.
Put another w a y , the brigh tness is distributed o v er a large area in the case of a nebula.
44 https://en.wikipedia.org/wiki/Book_of_Fixed_Stars
Chapter 2: T utorials 134
2.5 Detecting lines and extracting sp ectra in 3D data
3D data cub es are an increasingly common format of data pro ducts in observ ational as-
tronom y . As opp osed to 2D images (where eac h 2D “picture elemen t” or “pixel” co vers an
infinitesimal area on the surface of the sky), 3D data cub es con tain “v olume elemen ts” or
“v o xels” that are also connected in a third dimension.
The most common case of 3D data in observ ational astroph ysics is when the first t wo
dimensions are spatial (RA and Dec on the sky), and the third dimension is w a v elength. This
t yp e of data is generically (also outside of astronom y) kno wn as Hyp ersp ectral imaging 45 .
F or example high-lev el data pro ducts of In tegral Field Units (IFUs) like MUSE 46 in the
optical, A CIS 47 in the X-ra y , or in the radio where most data are 3D cub es.
In this tutorial, w e’ll use a small crop of a reduced deep MUSE cub e cen tered on the
Ab ell 370 ( https://en.wikipedia.org/wiki/Abell_370 ) galaxy cluster from the Pilot-
WINGS surv ey; see Lagattuta et al. 2022 ( https://arxiv.org/abs/2202.04663 ). Ab ell
370 has a spiral galaxy in its bac kground that is stretc hed due to the cluster’s gra vitational
p oten tial to create a b eautiful arc h. If y ou ha v en’t seen it y et, ha v e a lo ok at some of its
images in the Wikip edia link ab o v e b efore con tin uing.
The Pilot-WINGS surv ey data are a v ailable in its w ebpage 48 . The cub e of the c or e
region is 10.2GBs. This can b e prohibitiv ely large to download (and later pro cess) on
man y net w orks and smaller computers. Therefore, in this demonstration w e w on’t b e using
the full cub e. W e ha v e prepared a small crop 49 of the full cub e that you can do wnload
with the first command b elo w. The randomly selected crop is cen tered on (RA,Dec) of
(39.96769,-1.58930), with a width of ab out 27 arcseconds.
$ mkdir tutorial-3d
$ cd tutorial-3d
$ wget http://akhlaghi.org/data/a370-crop.fits # Downloads 287 MB
In the sections b elo w, w e will first review ho w y ou can visually insp ect a 3D data cub e
in DS9 and in teractiv ely see the sp ectra of an y region. W e will then subtract the con tin uum
emission, detect the emission-lines within this cub e and extract their sp ectra. W e will finish
with creating pseudo narro w-band images optimized for some of the emission lines.
2.5.1 Viewing sp ectra and redshifted lines
In Section 2.5 [Detecting lines and extracting sp ectra in 3D data], page 134, w e do wnloaded
a small crop from the Pilot-WINGS surv ey of Ab ell 370 cluster; observ ed with MUSE. In
45 https://en.wikipedia.org/wiki/Hyperspectral_imaging
46 https://en.wikipedia.org/wiki/Multi-unit_spectroscopic_explorer
47 https://en.wikipedia.org/wiki/Advanced_CCD_Imaging_Spectrometer
48 https://astro.dur.ac.uk/~hbpn39/pilot-wings.html
49 Y ou can download the full cube and create the crop your self with the commands below. Due to the
decompression of the + 10GB file that is necessary for the compressed do wnloaded file (note that its suffix
is .fits.gz ), the Crop command will tak e a little long.
$ wget https://astro.dur.ac.uk/~hbpn39/pilotWINGS/A370_PilotWINGS_CORE.fits.gz
$ astcrop A370_PilotWINGS_CORE.fits.gz -hDATA --mode=img \
--section=200:300,100:200 -oa370-crop.fits --metaname=DATA
$ astcrop A370_PilotWINGS_CORE.fits.gz -hSTAT --mode=img --append \
--section=200:300,100:200 -oa370-crop.fits --metaname=STAT
Chapter 2: T utorials 135
this section, w e will review ho w y ou can visualize/insp ect a data cub e using that example.
With the first command b elo w, we’ll open DS9 suc h that eac h 2D slice of the cub e (at a
fixed w a v elength) is seen as a single image. If y ou mov e the slider in the “Cub e” windo w
(that also op ens), y ou can view the same field at differen t w a v elengths. W e are ending the
first command with a ‘ & ’ so y ou can contin ue viewing DS9 while using the command-line
(press one extra ENTER to see the prompt). With the second command, y ou can see that
the spacing b et w een eac h slice is 1 . 25 × 10 − 10 meters (or 1.25 Angstroms).
$ astscript-fits-view a370-crop.fits -h1 --ds9scale="limits -5 20" &
$ astfits a370-crop.fits --pixelscale
Basic info. for --pixelscale (remove info with '--quiet' or '-q')
Input: a370-crop.fits (hdu 1) has 3 dimensions.
Pixel scale in each FITS dimension:
1: 5.55556e-05 (deg/pixel) = 0.2 (arcsec/pixel)
2: 5.55556e-05 (deg/pixel) = 0.2 (arcsec/pixel)
3: 1.25e-10 (m/slice)
Pixel area (on each 2D slice) :
3.08642e-09 (deg^2) = 0.04 (arcsec^2)
Voxel volume:
3.85802e-19 (deg^2*m) = 5e-12 (arcsec^2*m) = 0.05 (arcsec^2*A)
In the DS9 “Cub e” windo w, y ou will see t w o num b ers on the t w o sides of the scroller.
The left n um b er is the w a v elength in meters (W CS co ordinate in 3rd dimension) and the
righ t n um b er is the slice n umber (slice num b er or arra y co ordinates in 3rd dimension).
Y ou can man ually edit an y of these num b ers and press ENTER to go to that slice in an y
co ordinate system. If y ou w an t to go one-b y-one, simply press the “Next” button. The first
few slides are v ery noisy , but in the rest the noise lev el decreases and the galaxies are more
ob vious.
As y ou slide b et w een the differen t w a v elengths, y ou see that the noise-lev el is not constan t
and in some slices, the sky noise is v ery strong (for example, go to slice 3201 and press the
“Next” button). W e will discuss these issues b elo w (in Section 2.5.2 [Sky lines in optical
IFUs], page 138). T o view the sp ectra of a region in DS9 tak e the follo wing steps:
1. Clic k somewhere on the image (to mak e sure DS9 receiv es y our k eyb oard inputs),
then press Ctrl+R to activ ate regions and clic k on the brigh test galaxy of this cub e
(cen ter-righ t, at RA, Dec of 39.9659175 and -1.5893075).
2. A thin green circle will sho w up; this is called a “region” in DS9.
3. Double-clic k on the region, and y ou will see a “Circle” windo w.
4. Within the “Circle” windo w, clic k on the “Analysis” men u and select “Plot 3D”.
5. A second “Circle” windo w will op en that sho ws the sp ectra within y our selected region.
This is just the sum of v alues on eac h slice within the region.
6. Don’t close the second “circle” windo w (that sho ws the sp ectrum). Clic k and hold the
region in DS9, and mo ve it to other ob jects within the cub e. Y ou will see that the
sp ectrum c hanges as y ou mo v e the region, and y ou can see that differen t ob jects ha v e
v ery differen t sp ectra. Y ou can even see the spectra of only one part of a galaxy , not
the whole galaxy .
7. T ak e the region bac k to the first (brigh test) galaxy that w e originally started with.
Chapter 2: T utorials 136
8. Slide o v er differen t w a v elengths in the “Cub e” windo w, you will see the ligh t-blue line
mo ving through the sp ectrum as y ou slide to differen t w a v elengths. This line sho ws
the w a v elength of the displa y ed image in the main windo w o v er the sp ectra.
9. The strongest emission line in this galaxy app ears to b e around 8500 Angstroms or
8 . 5 × 10 − 7 meters. F rom the p osition of the Balmer break ( https://en.wikipedia.
org/wiki/Balmer_jump ) (blue-w ard of 5000 Angstroms for this galaxy), the strong
seems to b e H-alpha.
10. T o confirm that this is H-alpha, y ou can select the “Edit” men u in the sp ectrum windo w
and select “Zo om”.
11. Double-click and hold (for next step also) somewhere before the strongest line
and sligh tly ab o v e the con tin uum (for example at 8E-07 in the horizon tal and
60 × 10 − 20 erg/Angstrom/cm 2 /s on the v ertical). As y ou mo v e y our cursor (while
holding), y ou will see a rectangular b o x getting created.
12. Mo v e the b ottom-left corner of the b o x to somewhere after the strongest line and b elo w
the con tin uum. F or example at 9E-07 and 20 × 10 − 20 erg/Angstrom/cm 2 /s.
13. Once y ou remo v e y our finger from the mouse/touc hpad, it will zo om-in to that part of
the sp ectrum.
14. T o zo om out to the full sp ectrum, just press the right mouse button o v er the sp ectra
(or tap with t w o fingers on a touc hpad).
15. Select that zo om-b o x again to see the brigh test line m uc h more clearly . Y ou can also
see the t wo lines of the Nitrogen I I doublet that sandwic h H-alpha. Beside its relativ e
p osition to the Balmer break, this is further evidence that the strongest line is H-alpha.
16. Let’s ha v e a lo ok at the galaxy in its b est glory: righ t o v er the H-alpha line: Mo v e
the w a v elength slider accurately (b y pressing the “Previous” or “Next” buttons) suc h
that the blue line falls in the middle of the H-alpha line. W e see that the w a v elength
at this slice is 8.56593e-07 meters or 8565.93 Angstroms. Please compare the image
of the galaxy at this w a v elength with the w a v elengths b efore (b y pressing “Next” or
“Previous”). Y ou will also see that it is m uc h more extended and brigh ter than other
w a v elengths! H-alpha sho ws the un-obscured star formation of the galaxy!
Automaticly going to next slice: When y ou w an t to get a general feeling of the cub e,
pressing the “Next” button man y times is anno ying and slo w. T o automatically shift
b et w een the slices, y ou can press the “Pla y” button in the DS9 “Cub e” windo w. Y ou can
adjust the time it sta ys on eac h slice b y clic king on the “In terv al” men u and selecting
lo w er v alues.
Kno wing that this is H-alpha at 8565.93 Angstroms, y ou can get the redshift of the galaxy
with the first command b elo w and the lo cation of all other exp ected lines in Gn uastro’s
sp ectral line database with the second command. Because there are many lines in the
second command (more than 200!), with the third command, w e’ll only limit it to the
Balmer series (that start with H- ) using grep . The output of the second command prin ts
the metadata on the top (that is not sho wn any more in the third command due to the
grep call). T o b e complete, the first column is the observ ed w a v elength of the giv en line in
the giv en redshift and the second column is the name of the line.
# Redshift where H-alpha falls on 8565.93.
Chapter 2: T utorials 137
$ astcosmiccal --obsline=H-alpha,8565.93 --usedredshift
0.305221
# Wavelength of all lines in Gnuastro's database at this redshift
$ astcosmiccal --obsline=H-alpha,8565.93 --listlinesatz
# Only the Balmer series (Lines starting with 'H-'; given to Grep).
$ astcosmiccal --obsline=H-alpha,8565.93 --listlinesatz | grep H-
4812.13 H-19
4818.29 H-18
4825.61 H-17
4834.36 H-16
4844.95 H-15
4857.96 H-14
4874.18 H-13
4894.79 H-12
4921.52 H-11
4957.1 H-10
5006.03 H-9
5076.09 H-8
5181.83 H-epsilon
5353.68 H-delta
5665.27 H-gamma
6345.11 H-beta
8565.93 H-alpha
4758.84 H-limit
Zo om-out to the full sp ectrum and mo v e the displa y ed slice to the lo cation of the
first emission line that is blue-w ard (at shorter w av elengths) of H-alpha (at around 6300
Angstroms) and follo w the previous steps to confirm that y ou are on its cen ter. Y ou will
see that it falls exactly on 6 . 34468 × 10 − 7 m or 6344.68 Angstroms. No w, ha v e a lo ok at
the Balmer lines ab o v e. Y ou ha v e found the H-b eta line!
The rest of the Balmer series ( https://en.wikipedia.org/wiki/Balmer_series ) that
y ou see in the list ab o v e (lik e H-gamma, H-delta and H-epsilon) are visible only as absorp-
tion lines. Please c hec k their lo cation b y mo ving the blue line on the w a v elengths ab o v e
and confirm the sp ectral absorption lines with the ones ab o v e. The Balmer break is caused
b y the fact that these stronger Balmer absorption lines b ecome to o close to eac h other.
Lo oking bac k at the full sp ectrum, y ou can also confirm that the only other relatively
strong emission line in this galaxy , that is on the blue side of the sp ectrum is the w eak est
OI I line that is appro ximately lo cated at 4864 Angstroms in the observed spectra of this
galaxy . The n um b ers after the v arious OI I emission lines sho w their rest-frame w a v elengths
(“OI I” can corresp ond to man y electron transitions, so w e should b e clear ab out whic h one
w e are talking ab out).
$ astcosmiccal --obsline=H-alpha,8565.93 --listlinesatz | grep O-II-
4863.3 O-II-3726
4866.93 O-II-3728
5634.82 O-II-4317
Chapter 3: Installation 240
White space b et w een option and v alue: developer-build do es not accept an = sign
b et w een the options and their v alues. It also needs at least one c haracter b et w een the
option and its v alue. Therefore -n 4 or --numthreads 4 are acceptable, while -n4 , -n=4 ,
or --numthreads=4 are not. Finally m ultiple short option names cannot b e merged: for
example, y ou can sa y -c -n 4 , but unlik e Gn uastro’s programs, -cn4 is not acceptable.
Reusable for other pac k ages: This script can b e used in any soft w are whic h is configured
and built using the GNU Build System. Just copy it in the top source directory of that
soft w are and run it from there.
Example usage: See Section 13.12.4 [F orking tutorial], page 960, for an example usage of
this script in some scenarios.
-b STR
--top-build-dir STR
The top build directory to mak e a directory for the build. If this option is
not called, the top build directory is /dev/shm (only a v ailable in GNU/Lin ux
op erating systems, see Section 3.3.1.4 [Configure and build in RAM], page 238).
-V
--version
Prin t the v ersion string of Gn uastro that will b e used in the build. This string
will b e app ended to the directory name con taining the built files.
-a
--autoreconf
Run autoreconf -f b efore building the pac k age. In Gn uastro, this is necessary
when a new commit has b een made to the pro ject history . In Gn uastro’s build
system, the Git description will b e used as the v ersion, see Section 1.7 [V ersion
n um b ering], page 11, and Section 3.2.2.2 [Sync hronizing], page 228.
-c
--clean Delete the con ten ts of the build directory (clean it) b efore starting the config-
uration and building of this run.
This is useful when y ou ha v e recen tly pulled c hanges from the main Git rep osi-
tory , or committed a c hange y ourself and ran autoreconf -f , see Section 3.2.2.2
[Sync hronizing], page 228. After running GNU Auto conf, the v ersion will b e
up dated and y ou need to do a clean build.
-d
--debug Build with debugging flags (for example, to use in GNU Debugger, also kno wn
as GDB, or V algrind), disable optimization and also the building of shared
libraries. Similar to running the configure script of b elo w
$ ./configure --enable-debug
Chapter 3: Installation 241
Besides all the debugging adv an tages of building with this option, it will also
b e significan tly sp eed up the build (at the cost of slo w er built programs). So
when y ou are testing something small or w orking on the build system itself, it
will b e m uc h faster to test y our w ork with this option.
-v
--valgrind
Build all make check tests within V algrind. F or more, see the description
of --enable-check-with-valgrind in Section 3.3.1.1 [Gn uastro configure op-
tions], page 230.
-j INT
--jobs INT
The maxim um n um b er of threads/jobs for Mak e to build at an y momen t. As
the name suggests (Mak e has an iden tical option), the n um b er giv en to this
option is directly passed on to an y call of Mak e with its -j option.
-C
--check After finishing the build, also run make check . By default, make check is not
run b ecause the dev elop er usually has their own c hec ks to w ork on (for example,
defined in tests/during-dev.sh ).
-i
--install
After finishing the build, also run make install .
-D
--dist Run make dist-lzip pdf to build a distribution tarball (in .tar.lz format)
and a PDF man ual. This can b e useful for archiving, or sending to colleagues
who do not use Git for an easy build and man ual.
-u STR
--upload STR
Activ ate the --dist ( -D ) option, then use secure cop y ( scp , part of the SSH
to ols) to cop y the tarball and PDF to the src and pdf sub-directories of the
sp ecified serv er and its directory (v alue to this option). F or example, --upload
my-server:dir , will cop y the tarball in the dir/src , and the PDF man ual
in dir/pdf of my-server serv er. It will then mak e a sym b olic link in the top
serv er directory to the tarball that is called gnuastro-latest.tar.lz .
-p STR
--publish=STR
Clean, b o otstrap, build, c hec k and upload the c hec k ed tarball and PDF of the
b o ok to the URL giv en as STR . This option is just a wrapp er for --autoreconf
--clean --debug --check --upload STR . --debug is added b ecause it will
greatly sp eed up the build. --debug will ha v e no effect on the pro duced tarball
(p eople who later do wnload will b e building with the default optimized, and
non-debug mo de). This option is go o d when y ou ha v e made a commit and are
ready to publish it on y our server (if nothing crashes). Recall that if an y of the
previous steps fail the script ab orts.
Chapter 3: Installation 242
-I
--install-archive
Short for --autoreconf --clean --check --install --dist . This is useful
when y ou actually w an t to install the commit y ou just made (if the build and
c hec ks succeed). It will also produce a distribution tarball and PDF manual for
easy access to the installed tarball on y our system at a later time.
Ideally , Gn uastro’s Git v ersion history mak es it easy for a prepared system to
rev ert bac k to a differen t p oin t in history . But Gn uastro also needs to b o otstrap
files and also y our collab orators migh t (usually do!) find it to o m uc h of a burden
to do the b o otstrapping themselv es. So it is conv enien t to ha v e a tarball and
PDF man ual of the v ersion y ou ha v e installed (and are using in y our researc h)
handily a v ailable.
-h
--help
-P
--printparams
Prin t a description of this script along with all the options and their curren t
v alues.
3.3.3 T ests
After successfully building (compiling) the programs with the $ make command y ou can
c hec k the installation b efore installing. T o run the tests, run
$ make check
F or ev ery program some tests are designed to c hec k some p ossible op erations. Running
the command ab o v e will run those tests and giv e y ou a final rep ort. If ev erything is OK
and y ou ha v e built all the programs, all the tests should pass. In case an y of the tests fail,
please ha v e a lo ok at Section 3.3.5 [Kno wn issues], page 243, and if that still do es not fix
y our problem, lo ok that the ./tests/test-suite.log file to see if the source of the error is
something particular to y our system or more general. If you feel it is general, please con tact
us b ecause it migh t b e a bug. Note that the tests of some programs dep end on the outputs
of other program’s tests, so if you ha v e not installed them they migh t b e skipp ed or fail.
Prior to releasing ev ery distribution all these tests are c hec k ed. If you ha v e a reasonably
mo dern terminal, the outputs of the successful tests will b e colored green and the failed
ones will b e colored red.
These scripts can also act as a go o d set of examples for y ou to see ho w the programs are
run. All the tests are in the tests/ directory . The tests for eac h program are shell scripts
(ending with .sh ) in a sub-directory of this directory with the same name as the program.
See Section 13.7 [T est scripts], page 947, for more detailed information ab out these scripts
in case y ou w an t to insp ect them.
3.3.4 A4 prin t b o ok
The default prin t v ersion of this b o ok is pro vided in the letter pap er size. If y ou w ould lik e
to ha v e the prin t v ersion of this b o ok on pap er and y ou are living in a coun try whic h uses
A4, then y ou can rebuild the b o ok. The great thing ab out the GNU build system is that
Chapter 3: Installation 243
the b o ok source co de whic h is in T exinfo is also distributed with the program source co de,
enabling y ou to do suc h customization (hac king).
In order to c hange the pap er size, y ou will need to hav e GNU T exinfo installed. Op en
doc/gnuastro.texi with an y text editor. This is the source file that created this b o ok. In
the first few lines y ou will see this line:
@c@afourpaper
In T exinfo, a line is commen ted with @c . Therefore, un-commen t this line b y deleting the
first t w o c haracters suc h that it c hanges to:
@afourpaper
Sa v e the file and close it. Y ou can no w run the follo wing command
$ make pdf
and the new PDF b o ok will b e a v ailable in SRCdir/doc/gnuastro.pdf . By c hanging the
pdf in $ make pdf to ps or dvi y ou can ha v e the b o ok in those formats. Note that y ou can
do this for an y b o ok that is in T exinfo format, they migh t not ha v e @afourpaper line, so
y ou can add it close to the top of the T exinfo source file.
3.3.5 Kno wn issues
Dep ending on y our op erating system and the v ersion of the compiler you are using, y ou
migh t confron t some kno wn problems during the configuration ( $ ./configure ), compila-
tion ( $ make ) and tests ( $ make check ). Here, their solutions are discussed.
• $ ./configure : Configur e c omplains ab out not finding a libr ary even though you have
instal le d it. The p ossible solution is based on ho w y ou installed the pac k age:
• F rom y our distribution’s pac k age manager. Most probably this is b ecause y our
distribution has separated the header files of a library from the library parts.
Please also install the ‘dev elopment’ pac k ages for those libraries to o. Just add a
-dev or -devel to the end of the pac k age name and re-run the pac k age manager.
This will not happ en if y ou install the libraries from source. When installed from
source, the headers are also installed.
• F rom source. Then y our linker is not looking where y ou installed the library . If
y ou follo w ed the instructions in this c hapter, all the libraries will b e installed in
/usr/local/lib . So y ou ha v e to tell y our link er to lo ok in this directory . T o do
so, configure Gn uastro lik e this:
$ ./configure LDFLAGS="-L/usr/local/lib"
If y ou w an t to use the libraries for y our other programming pro jects, then exp ort
this en vironmen t v ariable in a start-up script similar to the case for LD_LIBRARY_
PATH explained b elo w, also see Section 3.3.1.2 [Installation directory], page 232.
• $ make : Complains ab out an unknown function on a non-GNU b ase d op er ating system.
In this case, please run $ ./configure with the --enable-gnulibcheck option to see if
the problem is from the GNU P ortabilit y Library (Gn ulib) not supp orting y our system
or if there is a problem in Gn uastro, see Section 3.3.1.1 [Gn uastro configure options],
page 230. If the problem is not in Gn ulib and after all its tests y ou get the same
complain t from make , then please con tact us at [email protected] . The cause is
probably that a function that w e ha v e used is not supp orted b y y our op erating system
Chapter 3: Installation 244
and w e did not included it along with the source tarball. If the function is a v ailable in
Gn ulib, it can b e fixed immediately .
• $ make : Cannot find the he aders (.h files) of instal le d libr aries. Y our C prepro cessor
(CPP) is not lo oking in the righ t place. T o fix this, configure Gn uastro with an addi-
tional CPPFLAGS lik e b elo w (assuming the library is installed in /usr/local/include :
$ ./configure CPPFLAGS="-I/usr/local/include"
If y ou w an t to use the libraries for y our other programming pro jects, then exp ort
this en vironmen t v ariable in a start-up script similar to the case for LD_LIBRARY_PATH
explained b elo w, also see Section 3.3.1.2 [Installation directory], page 232.
• $ make check : Only the first c ouple of tests p ass, al l the r est fail or get skipp e d. It is
highly lik ely that when searc hing for shared libraries, y our system do es not lo ok in to the
/usr/local/lib directory (or wherev er y ou installed Gn uastro or its dep endencies). T o
mak e sure it is added to the list of directories, add the follo wing line to y our ~/.bashrc
file and restart y our terminal. Do not forget to c hange /usr/local/lib if the libraries
are installed in other (non-standard) directories.
export LD_LIBRARY_PATH="$LD_LIBRARY_PATH:/usr/local/lib"
Y ou can also add more directories b y using a colon ‘ : ’ to separate them. See Sec-
tion 3.3.1.2 [Installation directory], page 232, and Section 12.1.2 [Linking], page 733,
to learn more on the PATH v ariables and dynamic linking resp ectiv ely .
• $ make check : The tests r elying on external pr o gr ams (for example, fitstopdf.sh
fail.) This is probably due to the fact that the v ersion num b er of the external programs
is to o old for the tests w e hav e preformed. Please update the program to a more recen t
v ersion. F or example, to create a PDF image, y ou will need GPL Ghostscript, but older
v ersions do not w ork, w e ha v e successfully tested it on v ersion 9.15. Older v ersions
migh t cause a failure in the test result.
• $ make pdf : The PDF b o ok c annot b e made. T o mak e a PDF b o ok, y ou need to ha v e
the GNU T exinfo program (lik e an y program, the more recen t the b etter). A w orking
T
E X program is also necessary , whic h you can get from T ex Liv e 25 .
• After make check : do not cop y the programs’ executables to another (for example,
the installation) directory man ually (using cp , or mv for example). In the default
configuration 26 , the program binaries need to link with Gn uastro’s shared library whic h
is also built and installed with the programs. Therefore, to run successfully b efore and
after installation, linking mo difications need to b e made b y GNU Libto ol at installation
time. make install do es this in ternally , but a simple cop y migh t giv e linking errors
when y ou run it. If y ou need to cop y the executables, y ou can do so after installation.
• $ make (when b o otstrapping): After y ou ha v e b o otstrapp ed Gn uastro from the v ersion-
con trolled source, y ou ma y confron t the follo wing (or a similar) error when con v erting
images (for more on b o otstrapping, see Section 3.2.2.1 [Bo otstrapping], page 226):
convert: attempt to perform an operation not allowed by the
security policy `gs' error/delegate.c/ExternalDelegateCommand/378.
25 https://www.tug.org/texlive/
26 If y ou configure Gn uastro with the --disable-shared option, then the libraries will b e statically link ed
to the programs and this problem will not exist, see Section 12.1.2 [Linking], page 733.
Chapter 3: Installation 245
This error is a kno wn issue 27 with ImageMagick security policies in some op erating
systems. In short, imagemagick uses Ghostscript for PDF, EPS, PS and XPS parsing.
Ho w ev er, b ecause some securit y vulnerabilities ha v e b een found in Ghostscript 28 , b y
default, ImageMagic k ma y b e compiled without Ghostscript library . In suc h cases, if
allo w ed, ImageMagic k will fall bac k to the external gs command instead of the library .
But this ma y b e disabled with the follo wing (or a similar) lines in /etc/ImageMagick-
7/policy.xml (an ything related to PDF, PS, or Ghostscript).
<policy domain="delegate" rights="none" pattern="gs" />
<policy domain="module" rights="none" pattern="{PS,PDF,XPS}" />
T o fix this problem, simply commen t suc h lines (b y placing a <!-- b efore eac h state-
men t/line and --> at the end of that statemen t/line).
If y our problem w as not listed ab o v e, please file a bug rep ort (Section 1.9 [Rep ort a bug],
page 15).
27 https://wiki.archlinux.org/title/ImageMagick
28 https://security.archlinux.org/package/ghostscript
[Document text truncated for crawler view.]