It is the year 953 A.D. and Abd al-rahman Sufi (903 – 986 A.D.)62 is in Shiraz as a guest astronomer. He had come there to use the advanced 123 centimeter astrolabe for his studies on the ecliptic. However, something was bothering him for a long time. While mapping the constellations, there were several non-stellar objects that he had detected in the sky, one of them was in the Andromeda constellation. During a trip he had to Yemen, Sufi had seen another such object in the southern skies looking over the Indian ocean. He was not sure if such cloud-like non-stellar objects (which he was the first to call ‘Sahābi’ in Arabic or ‘nebulous’) were real astronomical objects or if they were only the result of some bias in his observations. Could such diffuse objects actually be detected at all with his detection technique?
He still had a few hours left until nightfall (when he would continue his studies on the ecliptic) so he decided to find an answer to this question. He had thoroughly studied Claudius Ptolemy’s (90 – 168 A.D) Almagest and had made lots of corrections to it, in particular in measuring the brightness. Using his same experience, he was able to measure a magnitude for the objects and wanted to simulate his observation to see if a simulated object with the same brightness and size could be detected in simulated noise with the same detection technique. The general outline of the steps he wants to take are:
Fortunately Sufi had heard of GNU Astronomy Utilities from a colleague in Isfahan (where he worked) and had installed it on his computer a year before. It had tools to do all the steps above. He had used MakeProfiles before, but was not sure which columns he had chosen in his user or system-wide configuration files for which parameters, see Configuration files. So to start his simulation, Sufi runs MakeProfiles with the -P option to make sure what columns in a catalog MakeProfiles currently recognizes, and confirm the output image parameters. In particular, Sufi is interested in the recognized columns (shown below).
$ astmkprof -P [[[ ... Truncated lines ... ]]] # Output: type float32 # Type of output: e.g., int16, float32, etc. mergedsize 1000,1000 # Number of pixels along first FITS axis. oversample 5 # Scale of oversampling (>0 and odd). [[[ ... Truncated lines ... ]]] # Columns, by info (see `--searchin'), or number (starting from 1): ccol 2 # Coord. columns (one call for each dim.). ccol 3 # Coord. columns (one call for each dim.). fcol 4 # sersic (1), moffat (2), gaussian (3), point # (4), flat (5), circumference (6), distance # (7), custom-prof (8), azimuth (9), # custom-img (10). rcol 5 # Effective radius or FWHM in pixels. ncol 6 # Sersic index or Moffat beta. pcol 7 # Position angle. qcol 8 # Axis ratio. mcol 9 # Magnitude. tcol 10 # Truncation in units of radius or pixels. [[[ ... Truncated lines ... ]]]
In Gnuastro, column counting starts from 1, so the columns are ordered such that the first column (number 1) can be an ID he specifies for each object (and MakeProfiles ignores), each subsequent column is used for another property of the profile. It is also possible to use column names for the values of these options and change these defaults, but Sufi preferred to stick to the defaults. Fortunately MakeProfiles has the capability to also make the PSF which is to be used on the mock image and using the --prepforconv option, he can also make the mock image to be larger by the correct amount and all the sources to be shifted by the correct amount.
For his initial check he decides to simulate the nebula in the Andromeda constellation.
The night he was observing, the PSF had roughly a FWHM of about 5 pixels, so as the first row (profile) in the table below, he defines the PSF parameters.
Sufi sets the radius column (
rcol above, fifth column) to
5.000, he also chooses a Moffat function for its functional form.
Remembering how diffuse the nebula in the Andromeda constellation was, he decides to simulate it with a mock Sérsic index 1.0 profile.
He wants the output to be 499 pixels by 499 pixels, so he can put the center of the mock profile in the central pixel of the image which is the 250th pixel along both dimensions (note that an even number does not have a “central” pixel).
Looking at his drawings of it, he decides a reasonable effective radius for it would be 40 pixels on this image pixel scale (second row, 5th column below). He also sets the axis ratio (0.4) and position angle (-25 degrees) to approximately correct values too, and finally he sets the total magnitude of the profile to 3.44 which he had measured. Sufi also decides to truncate both the mock profile and PSF at 5 times the respective radius parameters. In the end he decides to put four stars on the four corners of the image at very low magnitudes as a visual scale. While he was preparing the catalog, one of his students approached him and was also following the steps.
As described above, the catalog of profiles to build will be a table (multiple columns of numbers) like below:
0 0.000 0.000 2 5 4.7 0.0 1.0 30.0 5.0 1 250.0 250.0 1 40 1.0 -25 0.4 3.44 5.0 2 50.00 50.00 4 0 0.0 0.0 0.0 6.00 0.0 3 450.0 50.00 4 0 0.0 0.0 0.0 6.50 0.0 4 50.00 450.0 4 0 0.0 0.0 0.0 7.00 0.0 5 450.0 450.0 4 0 0.0 0.0 0.0 7.50 0.0
This contains all the “data” to build the profile, and you can easily pass it to Gnuastro’s MakeProfiles: since Sufi already knows the columns and expected values very good, he has placed the information in the proper columns. However, when the student sees this, he just sees a mumble-jumble of numbers! Generally, Sufi explains to the student that even if you know the number positions and values very nicely today, in a couple of months you will forget! It will then be very hard for you to interpret the numbers properly. So you should never use naked data (or data without any extra information).
Data (or information) that describes other data is called “metadata”!
One common example is column names (the name of a column is itself a data element, but data that describes the lower-level data within that column: how to interpret the numbers within it).
Sufi explains to his student that Gnuastro has a convention for adding metadata within a plain-text file; and guides him to Gnuastro text table format.
Because we do not want metadata to be confused with the actual data, in a plain-text file, we start lines containing metadata with a ‘
For example, see the same data above, but this time with metadata for every column:
# Column 1: ID [counter, u8] Identifier # Column 2: X [pix, f32] Horizontal position # Column 3: Y [pix, f32] Vertical position # Column 4: PROFILE [name, u8] Radial profile function # Column 5: R [pix, f32] Effective radius # Column 6: N [n/a, f32] Sersic index # Column 7: PA [deg, f32] Position angle # Column 8: Q [n/a, f32] Axis ratio # Column 9: MAG [log, f32] Magnitude # Column 10: TRUNC [n/a, f32] Truncation (multiple of R) 0 0.000 0.000 2 5 4.7 0.0 1.0 30.0 5.0 1 250.0 250.0 1 40 1.0 -25 0.4 3.44 5.0 2 50.00 50.00 4 0 0.0 0.0 0.0 6.00 0.0 3 450.0 50.00 4 0 0.0 0.0 0.0 6.50 0.0 4 50.00 450.0 4 0 0.0 0.0 0.0 7.00 0.0 5 450.0 450.0 4 0 0.0 0.0 0.0 7.50 0.0
The numbers now make much more sense for the student! Before continuing, Sufi reminded the student that even though metadata may not be strictly/technically necessary (for the computer programs), metadata are critical for human readers! Therefore, a good scientist should never forget to keep metadata with any data that they create, use or archive.
To start simulating the nebula, Sufi creates a directory named simulationtest in his home directory.
Note that the
pwd command will print the “parent working directory” of the current directory (its a good way to confirm/check your current location in the full file system: it always starts from the root, or ‘
$ mkdir ~/simulationtest $ cd ~/simulationtest $ pwd /home/rahman/simulationtest
It is possible to use a plain-text editor to manually put the catalog contents above into a plain-text file. But to easily automate catalog production (in later trials), Sufi decides to fill the input catalog with the redirection features of the command-line (or shell). Sufi’s student was not familiar with this feature of the shell! So Sufi decided to do a fast demo; giving the following explanations while running the commands:
Shell redirection allows you to “re-direct” the “standard output” of a program (which is usually printed by the program on the command-line during its execution; like the output of
pwd above) into a file.
For example, let’s simply “echo” (or print to standard output) the line “This is a test.”:
$ echo "This is a test." This is a test.
As you see, our statement was simply “echo”-ed to the standard output!
To redirect this sentence into a file (instead of simply printing it on the standard output), we can simply use the
> character, followed by the name of the file we want it to be dumped in.
$ echo "This is a test." > test.txt
This time, the
echo command did not print anything in the terminal.
Instead, the shell (command-line environment) took the output, and “re-directed” it into a file called test.txt.
Let’s confirm this with the
ls command (
ls is short for “list” and will list all the files in the current directory):
$ ls test.txt
Now that you confirm the existence of test.txt, you can see its contents with the
cat command (short for “concatenation”; because it can also merge multiple files together):
$ cat test.txt This is a test.
Now that we have written our first line in test.txt, let’s try adding a second line (do not forget that our final catalog of objects to simulate will have multiple lines):
$ echo "This is my second line." > test.txt $ cat test.txt This is my second line.
As you see, the first line that you put in the file is no longer present!
This happens because ‘
>’ always starts dumping content to a file from the start of the file.
In effect, this means that any possibly pre-existing content is over-written by the new content!
To append new lines (or dumping new content at the end of existing content), you can use ‘
for example, with the commands below, first we will write the first sentence (using ‘
>’), then use ‘
>>’ to add the second and third sentences.
Finally, we will print the contents of test.txt to confirm that all three lines are preserved.
$ echo "My first sentence." > test.txt $ echo "My second sentence." >> test.txt $ echo "My third sentence." >> test.txt $ cat test.txt My first sentence. My second sentence. My third sentence.
The student thanked Sufi for this explanation and now feels more comfortable with redirection. Therefore Sufi continues with the main project. But before that, he deletes the temporary test file:
$ rm test.txt
To put the catalog of profile data and their metadata (that was described above) into a file, Sufi uses the commands below. While Sufi was writing these commands, the student complained that “I could have done in this in a text editor”. Sufi reminded the student that it is indeed possible; but it requires manual intervention. The advantage of a solution like below is that it can be automated (for example, adding more rows; for more profiles in the final image).
$ echo "# Column 1: ID [counter, u8] Identifier" > cat.txt $ echo "# Column 2: X [pix, f32] Horizontal position" >> cat.txt $ echo "# Column 3: Y [pix, f32] Vertical position" >> cat.txt $ echo "# Column 4: PROF [name, u8] Radial profile function" \ >> cat.txt $ echo "# Column 5: R [pix, f32] Effective radius" >> cat.txt $ echo "# Column 6: N [n/a, f32] Sersic index" >> cat.txt $ echo "# Column 7: PA [deg, f32] Position angle" >> cat.txt $ echo "# Column 8: Q [n/a, f32] Axis ratio" >> cat.txt $ echo "# Column 9: MAG [log, f32] Magnitude" >> cat.txt $ echo "# Column 10: TRUNC [n/a, f32] Truncation (multiple of R)" \ >> 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
To make sure if the catalog’s content is correct (and there was no typo for example!), Sufi runs ‘
cat cat.txt’, and confirms that it is correct.
Now that the catalog is created, Sufi is ready to call MakeProfiles to build the image containing these objects. He looks into his records and finds that the zero point magnitude for that night, and that particular detector, was 18 magnitudes. The student was a little confused on the concept of zero point, so Sufi pointed him to Brightness, Flux, Magnitude and Surface brightness, which the student can study in detail later. Sufi therefore runs MakeProfiles with the command below:
$ astmkprof --prepforconv --mergedsize=499,499 --zeropoint=18.0 cat.txt MakeProfiles 0.19 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 student to read through the printed output. As the statements say, two FITS files should have been created in the running directory. So Sufi ran the command below to confirm:
$ ls 0_cat_profiles.fits cat_profiles.fits cat.txt
The file 0_cat_profiles.fits is the PSF Sufi had asked for, and cat_profiles.fits is the image containing the main objects in the catalog. Sufi opened the main image with the command below (using SAO DS9):
$ astscript-fits-view cat_profiles.fits --ds9scale=95
The student could clearly see the main elliptical structure in the center. However, the size of cat_profiles.fits was surprising for the student, instead of 499 by 499 (as we had requested), it was 2615 by 2615 pixels (from the command below):
$ astfits cat_profiles.fits Fits (GNU Astronomy Utilities) 0.19 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 why oversampling is important in modeling, especially for parts of the image where the flux change is significant over a pixel. Recall that when you oversample the model (for example, by 5 times), for every desired pixel, you get 25 pixels (\(5\times5\)). Sufi then explained that after convolving (next step below) we will down-sample the image to get our originally desired size/resolution.
After seeing the image, the student complained that only the large elliptical model for the Andromeda nebula can be seen in the center.
He could not see the four stars that we 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 looking at the whole image at once, is that they only cover a single pixel!
To prove it, he centered the image around the coordinates 2308 and 2308, where one of the stars is located in the over-sampled image [you can do this in
ds9 by selecting “Pan” in the “Edit” menu, then clicking around that position].
Sufi then zoomed in to that region and soon, the star’s non-zero pixel could be clearly seen.
Sufi explained that the stars will take the shape of the PSF (cover an area of more than one pixel) after convolution. If we did not have an atmosphere and we did not need an aperture, then stars would only cover a single pixel with normal CCD resolutions. So Sufi convolved 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 convolution finished, Sufi opened cat_convolved.fits and the four stars could be easily seen now:
$ astscript-fits-view cat_convolved.fits --ds9scale=95
It was interesting for the student that all the flux in that single pixel is now distributed over so many pixels (the sum of all the pixels in each convolved star is actually equal to the value of the single pixel before convolution).
Sufi explained how a PSF with a larger FWHM would make the points even wider than this (distributing their flux in a larger area).
With the convolved image ready, they were prepared to re-sample it to the original pixel scale Sufi had planned [from the
$ astmkprof -P command above, recall that MakeProfiles had over-sampled the image by 5 times].
Sufi explained the basic concepts of warping the image to his student and ran Warp with the following 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 now has the correct pixel scale. However, the image is still larger than what we had wanted, it is \(523\times523\) pixels (not our desired \(499\times499\)). The student is slightly confused, so Sufi also re-samples the PSF with the same scale by 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 goes on to explain that frequency space convolution will dim the edges and that is why he added the --prepforconv option to MakeProfiles above. Now that convolution is done, Sufi can remove those extra pixels using Crop with the command below. Crop’s --section option accepts coordinates inclusively and counting from 1 (according to the FITS standard), so the crop region’s first pixel has to be 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
To fully convince the student, Sufi checks the size of the output of the crop command above:
$ 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\times499\) pixels and the mock Andromeda galaxy is centered on the central pixel. This is the same dimensions as Sufi had desired in the beginning. All this trouble was certainly worth it because now there is no dimming on the edges of the image and the profile centers are more accurately sampled.
The final step to simulate a real observation would be to add noise to the image. Sufi set the zero point magnitude to the same value that he set when making the mock profiles and looking again at his observation log, he had measured the background flux near the nebula had a per-pixel magnitude of 7 that night. For more on how the background value determines the noise, see Noise basics. So using these values he ran MakeNoise, and with the second command, he visually inspected the image.
$ astmknoise --zeropoint=18 --background=7 --output=out.fits \ cat_convolved_scaled_cropped.fits MakeNoise 0.19 started on Sat Oct 6 17:05:06 953 - Generator type: ranlxs1 - Generator seed: 1428318100 MakeNoise finished in: 0.033491 (seconds) $ astscript-fits-view out.fits
The out.fits file now contains the noised image of the mock catalog Sufi had asked for. The student had not observed the nebula in the sky, so when he saw the mock image in SAO DS9 (with the second command above), he understood why Sufi was dubious: it was very diffuse!
Seeing how the --output option allows the user to specify the name of the output file, the student was confused and wanted to know why Sufi had not used it more regularly before? Sufi explained that for intermediate steps, you can rely on the automatic output of the programs (see Automatic output). Doing so will give all the intermediate files a similar basic name structure, so in the end you can simply remove them all with the Shell’s capabilities, and it will be familiar for other users. So Sufi decided to show this to the student by making a shell script from the commands he had used before.
The command-line shell has the capability to read all the separate input commands from a file.
This is useful when you want to do the same thing multiple times, with only the names of the files or minor parameters changing between the different instances.
Using the shell’s history (by pressing the up keyboard key) Sufi reviewed all the commands and then he retrieved 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
base shell variables for easier customization later, and before every command, he added some comments (lines starting with #) for future readability.
Finally, Sufi pointed the student to Gnuastro’s General program usage tutorial, which has a full section on Writing scripts to automate the steps.
#!/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. 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. astmknoise --zeropoint=18 --background=7 --output=out.fits \ "$base"_convolved_scaled_cropped.fits # Remove all the temporary files. rm 0*.fits "$base"*.fits
He used this chance to remind the student of the importance of comments in code or shell scripts! Just like metadata in a dataset, when writing the code, you have a good mental picture of what you are doing, so writing comments might seem superfluous and excessive. However, in one month when you want to re-use the script, you have lost that mental picture and remembering it can be time-consuming and frustrating. The importance of comments is further amplified when you want to share the script with a friend/colleague. So it is good to accompany any step of a script, or code, with useful comments while you are writing it (create a good mental picture of why you are doing something: do not just describe the command, but its purpose).
Sufi then explained to the eager student that you define a variable by giving it a name, followed by an
= sign and the value you want.
Then you can reference that variable from anywhere in the script by calling its name with a
So in the script whenever you see
$base, the value we defined for it above is used.
If you use advanced editors like GNU Emacs or even simpler ones like Gedit (part of the GNOME graphical user interface) the variables will become a different color which can really help in understanding the script.
We have put all the
$base variables in double quotation marks (
") so the variable name and the following text do not get mixed, the shell is going to ignore the
" after replacing the variable value.
To make the script executable, Sufi ran the following command:
$ chmod +x mymock.sh
Then finally, Sufi ran the script, simply by calling its file name:
After the script finished, the only file remaining is the out.fits file that Sufi had wanted in the beginning.
Sufi then explained to the student how he could run this script anywhere that he has a catalog if the script is in the same directory.
The only thing the student had to modify in the script was the name of the catalog (the value of the
base variable in the start of the script) and the value to the
edge variable if he changed the PSF size.
The student was also happy to hear that he will not need to make it executable again when he makes changes later, it will remain executable unless he explicitly changes the executable flag with
The student was really excited, since now, through simple shell scripting, he could really speed up his work and run any command in any fashion he likes allowing him to be much more creative in his works. Until now he was using the graphical user interface which does not have such a facility and doing repetitive things on it was really frustrating and some times he would make mistakes. So he left to go and try scripting on his own computer. He later reminded Sufi that the second tutorial in the Gnuastro book as more complex commands in data analysis, and a more advanced introduction to scripting (see General program usage tutorial).
Sufi could now get back to his own work and see if the simulated nebula which resembled the one in the Andromeda constellation could be detected or not.
Although it was extremely faint63.
Therefore, Sufi ran Gnuastro’s detection software (NoiseChisel) to see if this object is detectable or not.
NoiseChisel’s output (out_detected.fits) is a multi-extension FITS file, so he used Gnuastro’s
astscript-fits-view program in the second command to see the output:
$ astnoisechisel out.fits $ astscript-fits-view out_detected.fits
In the “Cube” window (that was opened with DS9), if Sufi clicked on the “Next” button to see the pixels that were detected to contain significant signal. Fortunately the nebula’s shape was detectable and he could finally confirm that the nebula he kept in his notebook was actually observable. He wrote this result in the draft manuscript that would later become “Book of fixed stars”64.
He still had to check the other nebula he saw from Yemen and several other such objects, but they could wait until tomorrow (thanks to the shell script, he only has to define a new catalog). It was nearly sunset and they had to begin preparing for the night’s measurements on the ecliptic.
In Latin Sufi is known as Azophi. He was an Iranian astronomer. His manuscript “Book of fixed stars” contains the first recorded observations of the Andromeda galaxy, the Large Magellanic Cloud and seven other non-stellar or ‘nebulous’ objects.
The brightness of a diffuse object is added over all its pixels to give its final magnitude, see Brightness, Flux, Magnitude and Surface brightness. So although the magnitude 3.44 (of the mock nebula) is orders of magnitude brighter than 6 (of the stars), the central galaxy is much fainter. Put another way, the brightness is distributed over a large area in the case of a nebula.