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.