HOW-TO Use Galphot in Astro-WISE

Introduction

Galphot is a surface photometry tool, which fits ellipses to isophotes in galaxy profiles. It was written by Marijn Franx and is available for download from his website. The version on his website only works by using IRAF; for Astro-WISE a number of changes were done to make it work outside of IRAF.

Astro-WISE implementation

The main classes in the Galphot object model for Astro-WISE (i.e. those classes that are stored in the database, and must be queried on to get results) are these:

• GalPhotModel: The main Galphot class. Conceptually this is the model of a single galaxy. The model image is stored on the dataservers and can be retrieved. GalPhotModel contains methods to create the residual image etc. The model also contains a list of GalPhotEllipses.
• GalPhotEllipse: This is the equivalent of one line in the output table produced by Galphot. A GalPhotModel contains a list of GalPhotEllipses.
• GalPhotParameters: This class stores the Galphot configuration parameters (see Table 1).
• GalPhotList: This class is used as tool to link any group of GalPhotModels together through an identfier and optionally a name.

Each GalPhotModel is identified by a unique number (GPID). This also connects the GalPhotEllipses to the GalPhotModel they belong to.

GalPhotModels are made based on SourceLists, so in order to run Galphot within Astro-WISE one first has to make a SourceList and determine which source one wants to model with Galphot. This source is then identified with the identifier combination SLID, SID, which uniquely defines one source in a SourceList.

First step: making a SourceList or querying existing SourceLists

First we must decide which source to run galphot on. This means either making a new SourceList (running Sextractor) or querying existing SourceLists. How to make SourceLists is described in HOW-TO SourceLists in the Astro-WISE System.

For now we assume you know which SourceList you have made, and now want to find the SIDs of suitable galaxies.

# Query for SourceList with identifier 57424
sl = (SourceList.SLID == 57424)[0]

# This will be the output dictionary
d = {'SID':[], 'MAG_ISO':[], 'A':[], 'B':[], 'E_WCS':[]}

# Do a query to find bright galaxy-like objects
r = sl.sources.sql_query(d, '"MAG_ISO"<20.0 AND "E_WCS"<0.5')

At this point “d” contains the list of identifiers, magnitudes and A, B of the sources in the SourceList that matched the query. The list “d[‘SID’]” can be given as the “sids” argument to the GalPhotTask (see the next section).

Running Galphot; using the GalPhotTask

Here is how to create a GalPhotModel for the sourcelist with ID 57424 and source with ID 71 within it:

task = GalPhotTask(instrument='WFI', slid=57424, sids=[71], commit=1)
task.execute()

or equivalently using the DPU:

dpu.run('GalPhot', i='WFI', slid=57424, sids=[71], C=1)

Configuring GalPhotModel

Table 1: Overview of GalPhotParameters

parameter	description	type	default
ngal	Number of galaxies to model, fixed to 1	int	1
iter1	Number of iterations without calculating residuals	int	3
iter2	Number of iterations with calculating residuals	int	3
hmax	Maximum number of harmonicals	int	6
nsammax	Maximum number of samples along the ellipse	int	200
debug	Debug parameter	int	0
npolres	Order in interpolation in intensity	int	4
rmin	Minimum radius	float	0.5
rmax	Maximum radius	float	75.0
radfac	Radial scaling factor	float	1.1
rshap	Maximum radius for modifying ellipticity and PA	float	10.0
rcen	Maximum radius for modifying centers	float	10.0
errshap	Maximum error in shape of ellipse	float	0.5
errcen	Maximum error in position of ellipse	float	0.1
dposmax	Maximum change in center per iteration	float	0.1
dellmax	Maximum change in ellipticity per iteration	float	0.1
dangmax	Maximum change in position angle per iteration	float	0.1
r1	Maximum radius for single annulus sampling	float	20.0
r2	Maximum radius for multiple annulus sampling	float	0.0
fracmin	Minimum fraction of good points along the ellipse	float	0.4
cliplow	Fraction of points to clip at low end	float	0.1
cliphigh	Fraction of points to clip at high end	float	0.1
linear	Linear scaling	int	0
extend	Extend radial interval after fit	int	1
outfill	Calculate residual outside of galaxy	int	1
xpos	X-position (in region) of galaxy to model	float	-1.0
ypos	Y-position (in region) of galaxy to model	float	-1.0
region_xmin	Minimum x-coordinate of cutout region in larger frame	int	-1 [1]
region_xmax	Maximum x-coordinate of cutout region in larger frame	int	-1 [1]
region_ymin	Minimum y-coordinate of cutout region in larger frame	int	-1 [1]
region_ymax	Maximum y-coordinate of cutout region in larger frame	int	-1 [1]
region_scale	Factor which determines the size of the region	float	12.0 [2]
mask_areas	List of GalPhotMask (circles, squares) objects	list	[] [3]
mask_type	“sourcelist”, “segmentation” or “”for masking	str	segmentation [3]
mask_scale	Semi-major axis scale factor while masking	float	3.0 [3]

[1]	(1, 2, 3, 4) This is not a Galphot parameter, it is used in Astro-WISE to define the region extracted from SourceList.frame

[2]	This is not a Galphot parameter, instead it is used to define the region_* parameters

[3]	(1, 2, 3) This is not a Galphot parameter, instead it is used to write the ``deletion’’ file

Galphot can be configured with the Pars class as conventional in Astro-WISE :

p = Pars(GalPhotModel)
p.show()

this results in a list of configurable parameters (see Table 1), now set one of them:

p.GalPhotModel.process_params.r1 = 10.0

Feed the configuration to the task using the “get()” method of Pars:

task = GalPhotTask(instrument='WFI', slid=57424, sids=[71], pars=p.get(),
                   commit=1)
task.execute()

The process parameter dictionary can also be specified by hand directly:

d = {'GalPhotModel.process_params.r1': 10.0,
     'GalPhotModel.process_params.r2': 20.0}

task = GalPhotTask(instrument='WFI', slid=57424, sids=[71], pars=d, commit=1)
task.execute()

Masking other sources in the field

The following pixels/sources may be masked:

• Pixels with a weight of 0 (always).

• Rectangles and/or circles specified in a mask file.

  The GalPhotTask can be given a list of mask file filenames. This only works when working locally; no files are uploaded to the DPU. See below.

• Sources other than the primary source that are in the cutout may be automatically masked. This is done on the basis of either the SourceList or the SExtractor “segmentation” check image.

  This is controlled with the process parameters “mask_type” and “mask_scale”. See below.

Automatic masking

The automatic masking process is controlled by the following process parameters:

• mask_type: Values: ‘segmentation’, ‘sourcelist’ or “. If the value is segmentation the segmentation check image created along with each SourceList will be used to identify which pixels belong to other sources than the modelled source and they will be masked. If the value is sourcelist, the sourcelist is scanned for nearby and bright sources, and a circle of size mask_scale*A is used to mask such sources. If the value is ” (empty string) no masking of nearby sources will be performed.
• mask_scale: See above.

Manual masking

A mask file (“.del” file, in galphot convention) can be created by the user. This text file contains lines of either 3 or 4 numbers, describing circles and rectangles respectively. Example:

25 28 10      # Describes a circle at position (25,28) with radius 10 pixels.
20 30 35 45   # Describes the rectangle [20:30, 35:45]

Note that the pixel coordinates are in the coordinate frame of the original image, rather than the cutout region. If such a file is created it can be specified when running the task:

task = GalPhotTask(instrument='WFI', slid=155211, sids=[17044],
                   mask_files=['masks.txt'], commit=0)

If a mask file is specified, one must be given for each SID specified in the “sids” list.

Using an existing model as initial values

It is possible to use an existing GalPhotModel as input, i.e. as an “intable”, in the Galphot configuration. This is done by specifying the existing GalPhotModel’s GPID in the “gpid_in” option of the task:

task = GalPhotTask(instrument='WFI', slid=57424, sids=[71], gpid_in=6721,
                   commit=1)
task.execute()

or equivalently in the dpu call:

dpu.run('GalPhot', i='WFI', slid=57424, sids=[71], gpid_in=6721, C=1)

Using GalPhotList

The GalPhotList class is intended as a simple way to group GalPhotModels (and GalPhotEllipses). The way to use this class is to create it first, and only then run Galphot. I.e. first create and commit a GalPhotList, and then specify its GalPhotList identifier GPLID when running the GalPhot task on the DPU or locally:

# Create the GalPhotList object
l = GalPhotList()
l.name = 'test-run-1'
l.make()
l.commit()

# [schmidt] 16:27:12 - Set GalPhotList identifier GPLID to 100231

# Refer to the GalPhotList object by specifying its GPLID, as reported
# after committing the GalPhotList (see above).
dpu.run('GalPhot', i='WFI', slid=75637, sids=range(10,20), gplid=100231, C=1)

Now it is easy to query on the group of GalPhotModels you called “test-run-1” and inspect their residual images:

query = GalPhotModel.GPLID == 100211
for model in query: model.get_residual()

This will create all residual images, which can then be inspected with a FITS viewer.

Querying for results

When you have made a GalPhotModel, the ellipse parameters are stored in the database as a list of GalPhotEllipse objects, and the residual image is stored on a dataserver. How does one get to this information?

# Query for GalPhotModel based on its GPID

q = GalPhotModel.GPID == 20
model = q[0]

# Query for GalPhotModel based on SLID/SID

q = (GalPhotModel.SLID == 57424) & (GalPhotModel.SID == 71)
model = q.max('GPID')

# Find all GalPhotEllipses (radii) for the source defined by SLID=57424 and
# SID=71, for which the fitted intensity is larger than 50000 (ADU)

q = (GalPhotEllipse.SLID == 57424) & (GalPhotEllipse.SID == 71) &\
    (GalPhotEllipse.i > 50000)
GPIDs = [ellipse.GPID for ellipse in q]
gpids = []
for gpid in GPIDs:
    if not gpid in gpids:
        gpids.append(gpid)

models = []
for gpid in gpids:
    m = (GalPhotModel.GPID == gpid)[0]
    models.append(m)

# Now we can have a look at the residual images:

for m in models:
    m.retrieve()
    m.display()

# And plot ellipse parameters against eachother:

for m in models:
    x = [ellipse.i for ellipse in m.ellipses]
    y = [ellipse.r for ellipse in m.ellipses]
    pylab.scatter(x,y)

Table 2: Overview of GalPhotEllipse parameters

parameter	description	unit
r, dr	The root of the product of the major and minor axis of the ellipse	pixel
i, di	Intensity at the ellipse	counts/pixel\(^2\)
s, ds	The slope of the intensity, the derivative of i wrt r	None
x, dx	Central x position in the region	pixel
x_orig	Central x position in the original image	pixel
y, dy	Central y position in the region	pixel
y_orig	Central y position in the original image	pixel
eps, deps	The ellipticity of the ellipse, 1-b/a	None
pos, dpos	The position angle of the ellipse, measured with respect to the	None
	axis y=0, counter clockwise
c1, c2, …, c6	The cos(n*theta) term of the residuals along the ellipse	None
dc1, dc2, …, dc6	Errors in c1, c2, …, c6	None
s1, s2, …, s6	The sin(n*theta) term of the residuals along the ellipse	None
ds1, ds2, …, ds6	Errors in ds1, ds2, …, ds6	None
f1, f2, f3, f4	Flag?	None

Description of useful public methods of GalPhotModel

• get_model()

  Downloads the model image and returns it as a BaseFrame object.

• get_residual()

  Creates the residual image and returns it as a BaseFrame object.

• get_science()

  Extracts and downloads the region in the science image for which the model was derived, and returns it as a BaseFrame object.

• get_weight()

  Extracts and downloads the region in the weight image for which the model was derived and returns it as a BaseFrame object.

• show_model_parameters()

  Display a list of all ellipse parameters.

• get_model_parameters()

  Returns a list of dictionaries of all ellipse parameters. I.e. each item of the list is a dictionary which contains the description of one ellipse. See Table 2 for a description of the parameters.