R. W. O'CONNELL
August 2013
IDL EXERCISES III:
BASIC IMAGE PROCESSING
The test data files (*.fits) for Part III are in the
IDL exercises data directory. Copy
those to your local directory.
0) From a UNIX shell, start X-windows if it is not already running.
From the shell in an X-window, "cd" into your "IDL" directory. Normally
that would be named "~/idl".
Choose a graphics display method for images
For making image displays during these exercises, you can use
either the "direct graphics" display techniques described
in Part II or the GUI-based ATV program written by Aaron Barth.
If you plan to use direct graphics, first check whether your
computer display supports "true-color" graphics (for background,
see the section titled "Setting a Visual Display Mode on Your
Computer" in the IDL Guide).
From the shell within an X-window, type: "xdpyinfo"
If the response includes "TrueColor" or "24 planes" then your display
is a 24-bit, true-color device. If not, it is an 8-bit, pseudo-color
device.
After you start IDL, set the visual class of your display to
pseudo-color or emulated pseudo-color as follows:
If you have an 8-bit display, type DEVICE,RETAIN=2,PSEUDO_COLOR=8
If you have a 24-bit display, type DEVICE,RETAIN=2,TRUE_COLOR=24,DECOMPOSED=0
If you plan to use ATV, check that it is in your IDL path.
If not download the software from this link and place it
in your IDL path. ATV performance is independent of display
bit depth. One limitation of ATV is that only one display
window can be opened at a time, although you can "blink"
between 2 or 3 pre-loaded images.
For plots, you should use the direct graphics techniques
described in the IDL Guide
(basic plotting commands are unaffected by the bit depth of
your display).
Start IDL from the X-window.
On most UNIX systems, just type "idl". This will start IDL
running in the command-line mode using the window from which
you called it.
Before starting the exercises, give the IDL command:
ON_ERROR,1
This will return control to the main program level if
an error occurs in a called routine.
1) CRYPTO-IMAGE
The image file "mystery-image.fits" contains an encrypted,
low-contrast, 2-D message. See if you can determine what it is by
using any of the tools you have tried so far: plotting, IMLIST,
contrast adjustment, interactive value inspection, different color
tables, histograms, flagging special values, and so forth. Once
you find the message, develop a method for making it stand out
(i.e. maximize the signal-to-noise) and make a hardcopy (from
TVLASER or the File/WriteImage menu of ATV).
2) IMAGE COMPARISONS
Read in the two images of M87 (courtesy of Brian McNamara).
These are roughly co-centered. One was taken in the B band, the
other in the I band.
Compare the general statistics for the two images: mean, min, max,
variance. One image contains some flagged bad pixels; which is
it?
Load a black-on-white color table. (See
the IDL Guide for a sample direct graphics script.
In ATV, invert the default color table.)
Using a 0.3 power law display to preserve the dynamic range,
display the images in two windows, with the [MIN,MAX] for each
adjusted to show as much of the range of scientifically
meaningful intensity values as possible. The central density of
the two displays should be roughly the same. What differences
between the two images are apparent?
Now use the AstUseLib BLINK command to blink between the two
windows (in ATV, use the BLINK menu). Are differences easier to
see?
In direct graphics, "roll" a color table using the ICROLL routine
from Part II in conjunction with one of the
following standard IDL color tables: 11, 25, 31,
or 33. Note how different color table settings affect the
appearance of the images and highlight different features.
3) UNSHARP MASKING
Reload the black-on-white color table.
M87 represents a modestly complicated scene with a diffuse
background, a structured extended source (the famous nonthermal
jet), and a system of unresolved globular star clusters, the
brightest of which can be seen in the outer parts of the image.
There are many image processing techniques available to extract
information on the structure of a scene like this by enhancing
lower contrast features. Here, we'll illustrate two simple ones
which suppress low-frequency structure (the diffuse background).
Use the M87 I band image for these:
i) Edge detection: this is fast, if crude. Results also
depend on the orientation of structures in the image.
Simply display IMAGE - SHIFT(IMAGE,K,K) , where K is a
small integer (2 is a good starting value). For best
appearance, adjust the display to show a small range of
values around zero (e.g. -3 to +3).
When properly adjusted, your display should readily show
both the jet (in the bright central parts of the galaxy)
and the faint globular clusters scattered throughout the
frame. You can try various color tables to adjust the
appearance, although B&W may actually be the most useful.
ii) Unsharp masking: this is a similar differencing technique,
except that the subtracted image is a smoothed version of
the original. This removes low-frequency structure and is
especially useful where there is a bright background, as
in the case of M87. The fastest approach with IDL is to
use the SMOOTH(IMAGE,K) function. However, the smoothed
image would include the effects of compact structures
(e.g. the jet here, or stars in the case of a star
cluster). A better, if somewhat slower, approach is to
use the MEDIAN(IMAGE,K) function, which is resistant to a
minority of bright pixels within the KxK smoothing box.
Make displays for several values of K. Good trial values
are in the range 9 to 21. Displays are best for a small
positive range of values, e.g. 0 to 2.
4) SURFACE PHOTOMETRY
Use the SKY routine to obtain an estimate of the sky background on
the M87 I-band frame. What is the result? Does it agree with
the approximate sky value you can infer by, for example, sampling
blank regions of the frame with CURVAL or other techniques?
Refer to the SKY header for an explanation of how the program
works.
Determine in which pixel the galaxy nucleus lies on the M87 I band
frame as follows. Adjust the display until you can easily see
the starlike nucleus. Measure an approximate location
interactively with the cursor. In direct graphics, you can use
the AstUseLib CURVAL routine or the built-in CURSOR routine.
Then, take the [xx,yy] location of the nucleus and use it as
input to the AstUseLib CNTRD centroiding routine. In ATV,
place the cursor near the center of the nucleus and press "p".
What are the centroid values for the nucleus?
Do aperture photometry on the M87 I-band frame as follows:
Read the information header for the AstUseLib routine DIST_CIRCLE,
and be sure you understand how it works. You may want to try
using it on a small test array and plotting slices across the
array.
Devise a method, similar to that given as a sample script in the
DIST_CIRCLE header, to obtain the integrated flux of the M87 image
within a circular aperture (centered on the galaxy nucleus).
Measure the total fluxes within circular apertures of 10, 20,
30, and 40 pixel radii.
Make up a simple table (printed to your screen and your journal file) of
the aperture radius, total flux in aperture, and net flux in aperture
after an appropriate correction for the sky background is made.
By what percentage did the sky correction reduce the final flux in
the 40 pixel aperture?
To illustrate the general shape of the galaxy brightness profile
and emphasize any asymmetries, make a "pixel plot" of the
galaxy's surface brightness as follows:
Extract a subarray of 71x71 pixels centered on the nuclear pixel
you identified above. Call this FLX.
Use DIST_CIRCLE to make a matching array containing the distance of
each pixel from the center. Call this RAD.
Now open a plotting window and make a plot of the radius of
each pixel in the subarray vs its flux by typing
PLOT,RAD,FLX,PSYM=5. (Setting PSYM = 5 plots points as open
triangles; alternatively, use PSYM = 1 for plus signs.) Most
data points should fall in a compact distribution showing
that the source is highly symmetrical. Have you reached the
sky background level at the edge of the subarray? Is there
evidence of an important asymmetry? [A linear-log plot using
PLOT_IO can also be useful, but use WHERE to exclude any
negative flux points.]
Test the sensitivity of this technique to asymmetries by
re-extracting the subarray but with the center offset by a
couple of pixels from the true nuclear cenrtroid.
In ATV, you can see a similar radial plot by using the "p"
function and pressing the "Show radial profile" button on
the pop-up menu.
Another way to visualize the source structure is to make a
surface plot with the SURFACE routine. Open a plotting window
and give the command SURFACE,FLX. You can explore the shape of
the surface by viewing it from other angles
(e.g. SURFACE,FLX,AX=15,AZ=70) or inverting it (SURFACE,-FLX).
In ATV, you can quickly make nice surface plots by pressing
the "s" key.
5) POINT SOURCE PHOTOMETRY
Read in the image of the globular cluster M13 from the file
"m13-ex8.fits".
Display the image. The center of the cluster is very crowded but
is not saturated in this exposure. Find a set of display
parameters which retains information in the center while still
showing the stars in the outer part of the cluster.
Extract a 100x100 subarray of the image near the cluster center.
Open a plotting window, and make a SURFACE plot of the
subarray. Explore several regions at different radii. How
far from the center do you have to be to avoid blended point-spread
functions?
The ATV "s" feature is more convenient for exploring the
overlap problem, since it automatically extracts a subarray
of the main image.
Extract a 400x400 subarray of the image in its upper left hand
quadrant, and use this for the remainder of exercise (5).
The AstUseLib FIND routine is part of the DAOPHOT photometry package
written by Peter Stetson. This version has been converted to IDL
from the original FORTRAN. FIND locates stellar objects within
a given image. It is efficient and employs only a few parameters,
the most important of which is HMIN, the minimum flux of a candidate
source above the background for a threshold detection.
Apply FIND to the subarray. Use an HMIN value of 120 and a
convolution FWHM of 4 pixels. Leave other parameters at their
default values. In the initial call to FIND, be sure to
specify output variables to store the [x,y] and flux values
produced by the routine.
Display the subarray on your terminal and use the AstUseLib
TVBOX (or ATVPLOT) routine to place small boxes around all
the stars chosen by FIND. (11 pixels is a good boxsize; use
the /DEVICE keyword in TVBOX.
APER (in AstUseLib) is the next routine used in a standard DAOPHOT
reduction. It obtains circular aperture photometry of stars
selected by FIND. It uses a small inner aperture to measure
each star and a slightly larger annulus to measure the local
sky.
Use APER to obtain brightnesses for the sources located in the
FIND run above. Reasonable values for this image are 5
pixels for the "first" star aperture radius (press the return
key when you are asked for "another aperture radius") and 6
and 9 pixels for the sky aperture. Set the number of "photons
per ADU" equal to 1. The output will be expressed in
magnitudes with an arbitrary zero point. Be sure to store
the magnitude values. [Note for non-astronomers: "magnitudes"
are related to fluxes measured in power per unit area by
m = -2.5*alog10(f) + const. ]
Make a histogram of the cluster luminosity function you obtained.
Note: the standard procedure with DAOPHOT to use a measuring
aperture sufficient to include a significant portion but not
all of the light from the point spread function. Extrapolation
to the "integrated" light of the source is then determined from
a PSF derived from well-exposed but isolated stellar images.
END OF IDL EXERCISES Part III
Part II 2D arrays and image displays.
Part I scalars, vectors, plotting.
IDL Tutorial Home Page
Translations of this page:
Latvian
Slovak
Copyright © 2000-2020 Robert W. O'Connell. All rights reserved.
Content last modified by RWO August 2013