First, let’s create a directory named tutorial-zeropoint to keep things clean and work in that. Then, with the commands below, you can download an image from J-PLUS and SDSS. To speed up the analysis, the image is cropped to have a smaller region around its center.
$ mkdir tutorial-zeropoint
$ cd tutorial-zeropoint
$ jplusdr2=http://archive.cefca.es/catalogues/vo/siap/jplus-dr2/reduced
$ wget $jplusdr2/get_fits?id=771463 -O jplus.fits.fz
$ astcrop jplus.fits.fz --center=107.7263,40.1754 \
--width=0.6 --output=jplus-crop.fits
Although we cropped the J-PLUS image, it is still very large in comparison with the SDSS image (the J-PLUS field of view is almost \(1.5\times1.5\) deg\(^2\), while the field of view of SDSS in each filter is almost \(0.3\times0.5\) deg\(^2\)). Therefore, let’s download two SDSS images (and then decompress them) in the region of the cropped J-PLUS image to have a more accurate result compared to a single SDSS footprint: generally, your zero point estimation will have less scatter with more overlap between your reference image(s) and your input image.
$ sdssbase=https://dr12.sdss.org/sas/dr12/boss/photoObj/frames
$ wget $sdssbase/301/6509/5/frame-r-006509-5-0115.fits.bz2 \
-O sdss1.fits.bz2
$ wget $sdssbase/301/6573/5/frame-r-006573-5-0174.fits.bz2 \
-O sdss2.fits.bz2
$ bunzip2 sdss1.fits.bz2
$ bunzip2 sdss2.fits.bz2
To have a feeling of the data, let’s open the three images with astscript-fits-view using the command below.
Wait a few seconds to see the three images “blinking” one after another.
The largest one is the J-PLUS crop and the two smaller ones that partially cover it in different regions are from SDSS.
$ astscript-fits-view sdss1.fits sdss2.fits jplus-crop.fits \
--ds9extra="-lock frame wcs -single -zoom to fit -blink yes"
The test above showed that the three images are already astrometrically calibrated (the coverage of the pixel positions on the sky is correct in both). To confirm, you can zoom-in to a certain object and confirm it on a pixel level. It is always good to do the visual check above when you are confronted with new images (and may not be confident about the accuracy of the astrometry). Do not forget that the goal here is to find the calibration of pixel values; and that we assume pixel positions are already calibrated (the image already has a good astrometry).
The SDSS images are Sky subtracted, while this single-exposure J-PLUS image still contains the counts related to the Sky emission within them.
In the J-PLUS survey, the sky-level in each pixel is kept in a separate BACKGROUND_MODEL HDU of jplus.fits.fz; this allows you to use a different sky if you like.
The SDSS image FITS files also have multiple extensions.
To understand our inputs, let’s have a fast look at the basic info of each:
$ astfits sdss1.fits
Fits (GNU Astronomy Utilities) 0.23
Run on Fri Apr 14 11:24:03 2023
-----
HDU (extension) information: 'sdss1.fits'.
Column 1: Index (counting from 0, usable with '--hdu').
Column 2: Name ('EXTNAME' in FITS standard, usable with '--hdu').
('n/a': no name in HDU metadata)
Column 3: Image data type or 'table' format (ASCII or binary).
Column 4: Size of data in HDU.
Column 5: Units of data in HDU (only images).
('n/a': no unit in HDU metadata, or HDU is a table)
-----
0 n/a float32 2048x1489 nanomaggy
1 n/a float32 2048 n/a
2 n/a table_binary 1x3 n/a
3 n/a table_binary 1x31 n/a
$ astfits jplus.fits.fz
Fits (GNU Astronomy Utilities) 0.23
Run on Fri Apr 14 11:21:30 2023
-----
HDU (extension) information: 'jplus.fits.fz'.
Column 1: Index (counting from 0, usable with '--hdu').
Column 2: Name ('EXTNAME' in FITS standard, usable with '--hdu').
('n/a': no name in HDU metadata)
Column 3: Image data type or 'table' format (ASCII or binary).
Column 4: Size of data in HDU.
Column 5: Units of data in HDU (only images).
('n/a': no unit in HDU metadata, or HDU is a table)
-----
0 n/a no-data 0 n/a
1 IMAGE float32 9216x9232 adu
2 MASKED_PIXELS int16 9216x9232 n/a
3 BACKGROUND_MODEL float32 9216x9232 n/a
4 MASK_MODEL uint8 9216x9232 n/a
Therefore, in order to be able to compare the SDSS and J-PLUS images, we should first subtract the sky from the J-PLUS image.
To do that, we can either subtract the BACKGROUND_MODEL HDU from the IMAGE HDU using Arithmetic, or we can use NoiseChisel to find a good sky ourselves.
As scientists we like to tweak and be creative, so let’s estimate it ourselves with the command below.
Generally, you may not have a pre-estimated Sky estimation like above, so you should be prepared to subtract the sky yourself.
$ astnoisechisel jplus-crop.fits --output=jplus-nc.fits $ astscript-fits-view jplus-nc.fits
Notice that there is a relatively bright star in the center-bottom of the image.
In the “Cube” window, click on the “Next” button to see the DETECTIONS HDU.
The large footprint of the bright star is obvious.
Press the “Next” button one more time to get to the SKY HDU.
You see that in the center-bottom, the footprint of the large star is clearly visible in the measured Sky level.
This is not good!
With Sky values above 54 ADU in the center of the star (the white pixels).
This over-subtracted Sky level in part of the image will affect your magnitude measurements and thus the zero point!
In General program usage tutorial, we have a section on NoiseChisel optimization for detection, there is also a full tutorial on this in Detecting large extended targets.
Therefore, we will not go into the details of NoiseChisel optimization here.
Given the large images of J-PLUS, we will increase the tile-size to \(100\times100\) pixels and the number of neighbors to identify outlying tiles to 50 (these are usually the first parameters you should start editing when you are confronted with a new image).
After the second command, check the SKY extension to confirm that there is no footprint of any bright object there.
You will still see a gradient, but note the minimum and maximum values of the Sky level: their difference is more than 26 times smaller than the noise standard deviation (so statistically speaking, it is pretty flat!)
$ astnoisechisel jplus-crop.fits --output=jplus-nc.fits \
--tilesize=100,100 --outliernumngb=50
$ astscript-fits-view jplus-nc.fits
## Check that the gradient in the sky is statistically negligible.
$ aststatistics jplus-nc.fits -hSKY --minimum --maximum \
| awk '{print $2-$1}'
0.32809
$ aststatistics jplus-nc.fits -hSKY_STD --median
8.377977e+00
We are now ready to find the zero point!
First, let’s run the astscript-zeropoint with --help to see the option names (recall that you can see more details of each option in Invoking astscript-zeropoint).
For the first time, let’s use the script in the most simple state possible.
We will keep only the essential options: the names of the input and reference images (and their HDUs), the name of the output, and also two apertures with radii of 3 arcsec to start with:
$ astscript-zeropoint --help
$ astscript-zeropoint jplus-nc.fits --hdu=INPUT-NO-SKY \
--refimgs=sdss1.fits,sdss2.fits \
--output=jplus-zeropoint.fits \
--refimgszp=22.5,22.5 \
--refimgshdu=0,0 \
--aperarcsec=3
The output is a FITS table (because generally, you will give more apertures and choose the best one based on a higher-level analysis).
Let’s check the output’s internal structure with Gnuastro’s astfits program.
$ astfits jplus-zeropoint.fits ----- 0 n/a no-data 0 n/a 1 ZEROPOINTS table_binary 1x3 n/a 2 APER-3 table_binary 321x2 n/a
You can see that there are two HDUs in this file.
The HDU names give a hint, so let’s have a look at each extension with Gnuastro’s asttable program:
$ asttable jplus-zeropoint.fits --hdu=1 -i -------- jplus-zeropoint.fits (hdu: 1) ------- ----- ---- ------- No.Name Units Type Comment ------- ----- ---- ------- 1 APERTURE arcsec float32 n/a 2 ZEROPOINT mag float32 n/a 3 ZPSTD mag float32 n/a -------- Number of rows: 1 --------
As you can see, in the first extension, for each of the apertures you requested (APERTURE), there is a zero point (ZEROPOINT) and the standard deviation of the measurements on the apertures (ZPSTD).
In this case, we only requested one aperture, so it only has one row.
Now, let’s have a look at the next extension:
$ asttable jplus-zeropoint.fits --hdu=2 -i -------- jplus-zeropoint.fits (hdu: 2) ------- ----- ---- ------- No.Name Units Type Comment ------- ----- ---- ------- 1 MAG-REF f32 float32 Magnitude of reference. 2 MAG-DIFF f32 float32 Magnitude diff with input. -------- Number of rows: 321 --------
It contains a table of measurements for the aperture with the least scatter.
In this case, we only gave one aperture, so it is the same.
If you give multiple apertures, only the one with least scatter will be present by default.
In the MAG-REF column you see the magnitudes within each aperture on the reference (SDSS) image(s).
The MAG-DIFF column contains the difference of the input (J-PLUS) and reference (SDSS) magnitudes for each aperture (see Zero point estimation).
The two catalogs, created by the aperture photometry from the SDSS images, are merged into one so that there are more stars to compare.
Therefore, no matter how many reference images you provide, there will only be a single table here.
If the two SDSS images overlapped, each object in the overlap region would have two rows (one row for the measurement from one SDSS image, and another from the measurement from the other).
Now that we have obtained the zero point of the J-PLUS image, let’s go a little deeper into lower-level details of how this script operates. This will help you better understand what happened and how to interpret and improve the outputs when you are confronted with a new image and strange outputs.
To keep intermediate results the astscript-zeropoint script keeps temporary files in a temporary directory and later deletes it (and all the intermediate products).
If you like to check the temporary files of the intermediate steps, you can use --keeptmp option to not remove them.
Let’s take a closer look into the contents of each HDU.
First, we’ll use Gnuastro’s asttable to see the measured zero point for this aperture.
We are using -Y to have human-friendly (non-scientific!) numbers (which are sufficient here) and -O to also show the metadata of each column at the start.
$ asttable jplus-zeropoint.fits -Y -O # Column 1: APERTURE [arcsec,f32,] Aperture used. # Column 2: ZEROPOINT [mag ,f32,] Zero point (sig-clip median). # Column 3: ZPSTD [mag ,f32,] Zero point Standard deviation. 3.000 26.435 0.057
Now, let’s have a look at the first 10 rows of the second (APER-3) extension.
From the previous check we did above, we see that it contains 321 rows!
$ asttable jplus-zeropoint.fits -Y -O --hdu=APER-3 --head=10 # Column 1: MAG-REF [f32,f32,] Magnitude of reference. # Column 2: MAG-DIFF [f32,f32,] Magnitude diff with input. 16.461 30.035 16.243 28.209 15.427 26.427 20.064 26.459 17.334 26.425 20.518 26.504 17.100 26.400 16.919 26.428 17.654 26.373 15.392 26.429
But the table above is hard to interpret, so let’s plot it.
To do this, we’ll use the same astscript-fits-view command above that we used for images.
It detects if the file has a image or table HDU and will call DS9 or TOPCAT respectively.
You can also use any other plotter you like (TOPCAT is not part of Gnuastro), this script just calls it.
$ astscript-fits-view jplus-zeropoint.fits --hdu=APER-3
After TOPCAT opens, you can select the “Graphics” menu and then “Plain plot”.
This will show a plot with the SDSS (reference image) magnitude on the horizontal axis and the difference of magnitudes between the the input and reference (the zero point) on the vertical axis.
In an ideal world, the zero point should be independent of the magnitude of the different stars that were used. Therefore, this plot should be a horizontal line (with some scatter as we go to fainter stars). But as you can see in the plot, in the real world, this expected behavior is seen only for stars with magnitudes about 16 to 19 in the reference SDSS images. The stars that are brighter than 16 are saturated in one (or both) surveys71. Therefore, they do not have the correct magnitude or mag-diff. You can check some of these stars visually by using the blinking command above and zooming into some of the brighter stars in the SDSS images.
On the other hand, it is natural that we cannot measure accurate magnitudes for the fainter stars because the noise level (or “depth”) of each image is limited.
As a result, the horizontal line becomes wider (scattered) as we go to the right (fainter magnitudes on the horizontal axis).
So, let’s limit the range of used magnitudes from the SDSS catalog to calculate a more accurate zero point for the J-PLUS image.
For this reason, we have the --magnituderange option in astscript-zeropoint.
|
Necessity of sky subtraction:
To obtain this horizontal line, it is very important that both your images have been sky subtracted.
Please, repeat the last Remember: if you see a sloped line instead of a horizontal line, the input or reference image(s) are not sky subtracted. |
Another key parameter of this script is the aperture size (--aperarcsec) for the aperture photometry of images.
On one hand, if the selected aperture is too small, you will be at the mercy of the differing PSFs between your input and reference image(s): part of the light of the star will be lost in the image with the worse PSF.
On the other hand, with large aperture size, the light of neighboring objects (stars/galaxies) can affect the photometry.
We should select an aperture radius of the same order than the one used in the reference image, typically 2 to 3 times the PSF FWHM of the images.
For now, let’s assume the values 2, 3, 4, 5, and 6 arcsec for the aperture sizes parameter.
The script will compare the result for several aperture sizes and choose the one with least standard deviation value, ZPSTD column of the ZEROPOINTS HDU.
Let’s re-run the script with the following changes:
$ astscript-zeropoint jplus-nc.fits --hdu=INPUT-NO-SKY \
--refimgs=sdss1.fits,sdss2.fits \
--output=jplus-zeropoint.fits \
--refimgszp=22.5,22.5 \
--aperarcsec=2,3,4,5,6 \
--magnituderange=16,18 \
--refimgshdu=0,0 \
--keepzpap
Now, check number of HDU extensions by astfits.
$ astfits jplus-zeropoint.fits ----- 0 n/a no-data 0 n/a 1 ZEROPOINTS table_binary 5x3 n/a 2 APER-2 table_binary 319x2 n/a 3 APER-3 table_binary 321x2 n/a 4 APER-4 table_binary 323x2 n/a 5 APER-5 table_binary 323x2 n/a 6 APER-6 table_binary 325x2 n/a
You can see that the output file now has a separate HDU for each aperture (thanks to --keepzpap.)
The ZEROPOINTS hdu contains the final zero point values for each aperture and their error.
The best zero point value belongs to the aperture that has the least scatter (has the lowest standard deviation).
The rest of extensions contain the zero point value computed within each aperture (as discussed above).
Let’s check the different tables by plotting all magnitude tables at the same time with TOPCAT.
$ astscript-fits-view jplus-zeropoint.fits
After TOPCAT has opened take the following steps:
APER-6, with red points).
The Bottom-left panel has the logo of a red-blue scatter plot that has written 6:jplus-zeropoint.fits in front of it (showing that this is the 6th HDU of this file).
In the bottom-right panel, you see the names of the columns that are being displayed.
<no table>.
In the bottom-right panel, in front of “Table:”, select any other extension.
This will plot the same two columns of that extension as blue points.
Zoom-in to the region of the horizontal line to see/compare the different scatters.
Change the HDU given to “Table:” and see the distribution of zero points for the different apertures.
The manual/visual operation above is critical if this is your first time with a new dataset (it shows all kinds of systematic biases (like the Sky issue above)!
But once you know your data has no systematic biases, choosing between the different apertures is not easy visually!
Let’s have a look at the table the ZEROPOINTS HDU (we don’t need to explicitly call this HDU since it is the first one):
$ asttable jplus-zeropoint.fits -O -Y # Column 1: APERTURE [arcsec,f32,] Aperture used. # Column 2: ZEROPOINT [mag ,f32,] Zero point (sig-clip median). # Column 3: ZPSTD [mag ,f32,] Zero point Standard deviation. 2.000 26.405 0.028 3.000 26.436 0.030 4.000 26.448 0.035 5.000 26.458 0.042 6.000 26.466 0.056
The most accurate zero point is the one where ZPSTD is the smallest.
In this case, minimum of ZPSTD is with radii of 2 and 3 arcseconds.
Run the astscript-fits-view command above again to open TOPCAT.
Let’s focus on the magnitude plots in these two apertures and determine a more accurate range of magnitude.
The more reliable option is the range between 16.4 (where we have no saturated stars) and 18.5 mag (fainter than this, the scatter becomes too strong).
Finally, let’s set some more apertures between 2 and 3 arcseconds radius:
$ astscript-zeropoint jplus-nc.fits --hdu=INPUT-NO-SKY \
--refimgs=sdss1.fits,sdss2.fits \
--output=jplus-zeropoint.fits \
--magnituderange=16.4,18.5 \
--refimgszp=22.5,22.5 \
--aperarcsec=2,2.5,3,3.5,4 \
--refimgshdu=0,0 \
--keepzpap
$ asttable jplus-zeropoint.fits -Y
2.000 26.405 0.037
2.500 26.425 0.033
3.000 26.436 0.034
3.500 26.442 0.039
4.000 26.449 0.044
The aperture with the least scatter is therefore the 2.5 arcsec radius aperture, giving a zero point of 26.425 magnitudes for this image.
However, you can see that the scatter for the 3 arcsec aperture is also acceptable.
Actually, the ZPSTD for of the 2.5 and 3 arcsec apertures only have a difference of \(3\%\) (\(= (0.034−0.0333)/0.033\times100\)).
So simply choosing the minimum is just a first-order approximation (which is accurate within \(26.436−26.425=0.011\) magnitudes)
Note that in aperture photometry, the PSF plays an important role (because the aperture is fixed but the two images can have very different PSFs). The aperture with the least scatter should also account for the differing PSFs. Overall, please, always check the different and intermediate steps to make sure the parameters are the good so the estimation of the zero point is correct.
If you are happy with the minimum, you don’t have to search for the minimum aperture or its corresponding zero point yourself.
This script has written it in ZPVALUE keyword of the table.
With the first command, we also see the name of the file also, (you can use this on many files for example).
With the second command, we are only printing the number by adding the -q (or --quiet) option (this is useful in a script where you want to write the value in a shell variable to use later).
$ astfits jplus-zeropoint.fits --keyvalue=ZPVALUE jplus-zeropoint.fits 2.642512e+01 $ astfits jplus-zeropoint.fits --keyvalue=ZPVALUE -q 2.642512e+01
Generally, this script will write the following FITS keywords (all starting with ZP) for your future reference in its output:
$ astfits jplus-zeropoint.fits -h1 | grep ^ZP ZPAPER = 2.5 / Best aperture. ZPVALUE = 26.42512 / Best zero point. ZPSTD = 0.03276644 / Best std. dev. of zeropoint. ZPMAGMIN= 16.4 / Min mag for obtaining zeropoint. ZPMAGMAX= 18.5 / Max mag for obtaining zeropoint.
Using the --keyvalue option of the Fits program, you can easily get multiple of the values in one run (where necessary):
$ astfits jplus-zeropoint.fits --hdu=1 --quiet \
--keyvalue=ZPAPER,ZPVALUE,ZPSTD
2.500000e+00 2.642512e+01 3.276644e-02
To learn more about saturated pixels and recognition of the saturated level of the image, please see Saturated pixels and Segment’s clumps
JavaScript license information
GNU Astronomy Utilities 0.23 manual, July 2024.