Finds ellipse for specified fraction
dmellipse infile outfile fraction [shape] [x_centroid] [y_centroid] [angle] [ellipticity] [xcol] [ycol] [zcol] [fix_centroid] [fix_angle] [fix_ellipticity] [tolerance] [minstep] [maxwalk] [step] [normalize] [clobber] [verbose]
`dmellipse' finds an ellipse that encloses the user specified fraction of the flux in an image or table to within some specified tolerance. The default behavior is to find an ellipse whose center is the centroid of the data, and whose radii align with the orientation (moments) of the data.
The tool works by computing the centroid and moments of the entire dataset. The radii of the ellipse and angle of the ellipse are computed from the 1st and 2nd moments. The resulting ellipse is then used to filter the data. If the sum of the pixel values in the ellipse is too small, the radii are increased by 1 pixel. If the sum is too larger, the radii is decreased by 1 pixel. A new set of centroid, radii, and angle are computed from the data in just the filtered region. Again if the sum of the pixels exceeds the fraction by more than the tolerance, the radii are shrunk by 1. If the enclosed fraction is too small, the radii are increased by 1 pixel. The process then repeats. The step size is reduced by half each time the previous iteration and the current iteration differ in whether the fraction is too big or too small. This provides a means for the algorithm to converge onto the solution. The iterations stop when the fraction of the flux in the image is with 'tolerance' of the specified fraction.
Due to finite pixel sizes, there may not be an ellipse that can meet the fraction criteria. In that case an warning is printed to the screen and the closest value of the ellipse parameter are saved.
In some rare situations changing the radii does not result in any change of the image moments (centroid, fraction or radii). For example a single bright pixel containing 50% of the image flux with a very large region of pixels with 0 flux surrounding it. dmellipse may get "stuck" in this case trying to increase or decrease the radii by fractions of a pixel without effecting any change. There is a safeguard built into the tool to terminate any such searches after "maxwalk" iterations. If this occurs users may need to set the 'step' parameter to a larger initial value (or maybe rebin the image).
dmellipse can also be used to find the ellipse the encloses a weighted fraction of rows table. If the infile is a table the xcol and ycol parameters are used to specify the names of the columns to use for the coordinates. The optional zcol can be used to specify the column name to used to weight the value for the rows. If no zcol is specified, all the rows are given the same unit weight. The rest of the algorithm as discussed above is then used to find the ellipse that encloses the requested fraction.
dmellipse img.fits 50perc.reg 0.5
Finds the ellipse that encloses 50% of the energy/flux in the image.
dmellipse img.fits 50perc.reg 0.5,0.9
Finds the ellipses that encloses 50% and 90% of the energy/flux in the image.
dmellipse img.fits 50perc.reg lgrid(0.5:0.91:0.1)
Finds the ellipse that encloses 50%, 60%, 70%, 80% and 90% of the energy/flux in the image.
Note that to get the 90% radii we used an upper limit of 0.91 in the grid specification. This is because the stack grid computes values x<value for floating-point-value not x<=value.
dmellipse img.fits 50perc.reg 0.5 fix_centroid=yes
Finds the ellipses that encloses 50% and 90% of the energy/flux in the image. The centroid is fixed to be the centroid of the entire image; rather than being left to vary.
dmellipse img.fits 50perc.reg 0.5 fix_centroid=yes x_centroid=4096 y_centroid=4096
Similar to above. Finds the ellipses that encloses 50% and 90% of the energy/flux in the image. The centroid is fixed to be the user supplied value (4096,4096).
dmellipse img.fits 50perc.reg 0.5 shape=rotbox
Finds the rotated-box (rectangle) that encloses 50% of the energy/flux in the image.
dmellipse img.fits[sky=circle(4096,4096,100)] 50perc.reg 0.5 normalize=no
Find ellipse that encloses 50% of the flux; setting normalize=no lets the tool know that the data have already been normalized (so 'no' the tool should not normalize the data).
dmellipse img.fits eef.fits fraction=0.9 shape=ellipse fix_ellipticity=yes ellipticity=1
By freezing the ellipticity to a value=1, the ellipses are forced to be circular. Since the angle of a circle is undefined, it could also be frozen to be 0 by setting fix_angle=yes angle=0. The output file will still contain "shape=ellipse" however both radii will be the same.
dmellipse marx_output_evts.fits outfile=90percent_ellipse.reg xcol=x ycol=y fraction=0.9 norm=yes
This example shows how a table can be used. The input table must have the column specified by the xcol and ycol parameters. The values in those columns are used to compute an ellipse that encloses 90% of the data.
dmellipse "sherpa_mcmc.dat[#row=1000:]" out=1sigma.reg frac=0.68 xcol=g1_xpos ycol=g1_ypos
Another common use may be to find an error ellipse for fitted parameters computed from the MCMC draws output. In this example the 1st 1000 rows are skipped (taken to be the 'burn in' interval), and then the position error ellipse is computed from the draws values for the fitted Gaussian postion values: g1_xpos and g1_ypos.
Detailed Parameter Descriptions
Parameter=infile (file required filetype=input)
Input image or table
Input image. It must have all positive values (so if background subtracted, make sure to threshold).
If the input is a table, then the xcol and ycol must also be specified. A per-row weighting column, the same as an image pixel value, can also be supplied using the zcol.
Parameter=outfile (file required filetype=output)
Output region file
ASC-FITS region file with additional columns indicating which fraction and status information. It can be used directly in ds9 to visualize the ellipses or used as any other CIAO region filter.
The "fraction" column gives the fraction of the flux found in the region. The requested input fraction is stored in the "input_fraction" column. The tolerance on the fraction is stored in the FRACTOL keyword.
If the fraction is not with tolerance of the specified value, then the "state" and "converge" columns will indicate
- 'OK' : state=0, fraction is within tolerance of requested value
- 'CANNOT_CONVERGE' : state=1, due to finite pixel sizes cannot make ellipse that meets tolerance but value should be "close".
- 'ELLIPSE_TOO_SMALL' : state=2, more than specified fraction is contained in a single pixel; cannot make ellipse any smaller.
- 'ELLIPSE_TOO_BIG' : state=3, specified fraction is not contained in the entire image (usually when normalize=no). Cannot make an ellipse big enough.
- 'WALK_TOO_MANY_STEPS' : state=4, too many iterations where centroid and radii did not change.
Parameter=fraction (string required stacks=yes)
Enclosed counts/flux fraction.
The value(s) to find the ellipse. This can be specified using any of the stack syntax including lgrid() or just a comma separated list (or just a single value).
Parameter=shape (string not required default=ellipse)
shape of region
can be either 'ellipse' or 'rotbox'
Parameter=x_centroid (real not required default=INDEF)
The default x-centroid can be input via this parameter or if the value is set to INDEF the value will be computed over the entire image.
Parameter=y_centroid (real not required default=INDEF)
The default y-centroid can be input via this parameter or if the value is set to INDEF the value will be computed over the entire image.
Parameter=angle (real not required default=INDEF)
The default angle can be input via this parameter or if the value is set to INDEF the value will be computed over the entire image.
Parameter=ellipticity (real not required default=INDEF)
The default ellipticity can be input via this parameter or if the value is set to INDEF the value will be computed over the entire image.
Parameter=xcol (string not required default=)
X-coordinate column name in the input table
The column name in the input table to use for the X coordinate of the dataset.
This parameter is ignored when the infile is an image.
Parameter=ycol (string not required default=)
Y-coordinate column name in the input table
The column name in the input table to use for the Y coordinate of the dataset.
This parameter is ignored when the infile is an image.
Parameter=zcol (string not required default=)
Column name of the optional row weights in the input table.
If zcol is left blank, then each row in the input table is weighted equally. If not, then each row in the input table is weighted by the values in the zcol column. The easiest interpretation of this is that xcol and ycol are the coordinates, zcol is the pixel value in the image.
This parameter is ignored when the infile is an image.
Parameter=fix_centroid (boolean not required default=no)
allow centroid to vary
during the iterations the centroid can be allowed to vary or it can be held fixed at the default value (see x_ and y_centroid for default value)
Parameter=fix_angle (boolean not required default=no)
allow angle to vary
during the iterations the angle can be allowed to vary or it can be held fixed at the default value (see angle parameter for the default value)
Parameter=fix_ellipticity (boolean not required default=no)
allow ellipticity to vary
during the iterations the ellipticity can be allowed to vary or it can be held fixed at the default value (see ellipticity parameter for the default value).
Parameter=tolerance (real not required default=0.001 min=0 max=1)
How close to desired fraction is close enough?
The iterations will continue until the fraction is within fraction-tolerance to fraction+tolerance (so total-width is 2*tolerance).
Parameter=minstep (real not required default=0.001 min=0 max=1)
minimum step size
As the iterations proceed and we change from increasing to decreasing the radii we make the step size increasing smaller; this is the minimum step size to use. After reaching this the iterations will stop.
Parameter=maxwalk (integer not required default=10 min=1)
maximum number of equally applicable iterations.
When we have a candidate radii we use the ellipticity, angle, and centroid to do one last check. If the ellipse defined by those is still good (fraction is within tolerance) the we stop. If not, we 'walk' to the next iteration. This could lead to a situation where we oscillate between good and bad without changing step-size. This parameter limits this condition.
Parameter=step (real not required default=1 min=0)
Initial step size
The intial step size. If the program seems to be taking a long time it may be because the algorithm is creaping slowly towards the solution. Setting verbose greater than 2 will show the progress of the tool. Increasing the step size will allow the program to get close to the optimal solution quicky and then more quickly refine the solution until the tolerance (or other exit conditions) is (are) met.
Parameter=normalize (boolean not required default=yes)
Normalize input image or not?
Is the fraction really a fraction (% of total flux) or is it the integral?
Parameter=clobber (boolean default=no)
Remove output if it exists?
Used to specify whether or not to clobber existing file that has the same name as the specified output file
Parameter=verbose (integer default=0 min=0 max=5)
The tool chatter level
Verbose can be from 0 to 5, generating different amounts of debugging output.
Changes in CIAO 4.14
- Copies the generalized world coordinate system information into the output file.
There are no known bugs for this tool.
- dmappend, dmcontour, dmfilth, dmimg2jpg, dmimgadapt, dmimgblob, dmimgcalc, dmimgdist, dmimgfilt, dmimghist, dmimghull, dmimglasso, dmimgpick, dmimgpm, dmimgproject, dmimgreproject, dmimgthresh, dmmaskbin, dmmaskfill, dmnautilus, dmradar, dmregrid, dmregrid2, dmstat, evalpos, imgmoment, mean_energy_map, pileup_map