Getting Started Documentation Glish Learn More Programming Contact Us
Version 1.9 Build 1488
News FAQ
Search Home


next up previous contents index
Next: regionmanager - Constructor Up: images - Module Previous: imagetest - Function


regionmanager - Tool



Package general
Module images


Postscript file available

Create and manipulate regions of interest

include "regionmanager.g"

Constructors
regionmanager Construct a regionmanager
Functions
absreltype Convert region type value to a string
box Create a pixel box region
complement Create the complement of a world region
concatenation Concatenate world regions along a new axis
copyregions Copy regions from one Table to another
def Default pixel coordinate value
deletefromtable Delete regions from a Table
difference Create the difference of two world regions
displayedplane Find a region describing the plane displayed by the viewer
done Destroy this regionmanager
extension Extend a world region to extra dimensions
extractsimpleregions Extract all simple regions into a record
fromglobaltotable Save regions into a Table
fromrecordtotable Save regions stored in a record into a Table
fromtabletoglobal Restore a region from a Table to Glish
fromtabletonakedrecord Restore all regions from a Table to a naked record
fromtabletorecord Restore regions from a Table to a record
getselectcallback Get the select callback function
gui Start the GUI interface
intersection Create the intersection of some world regions
ispixelregion Is this region a pixel region ?
isworldregion Is this region a world region ?
namesintable Find the names of the regions stored in a Table
pixeltoworldregion Convert a simple pixel region to a world region
pseudotoworldregion Convert a Viewer pseudoregion to a world region
quarter Create a quarter pixel region
setcoordinates Set new default Coordinate System
setselectcallback Set the callback function when GUI region selection is made
type Return the type of this tool
union Create a union of world regions
wbox Create a world box region
wmask Create a world mask region
wpolygon Create a world polygon region with quantities
wrange Create a world region for a pixel range on one axis


Overview of Regionmanager functionality


Regions and the Regionmanager

When working with an image, one is generally interested in some part of that image on which astrophysical analysis is performed. This region-of-interest (or more generically and simply, the `region') might be the whole image, or some subset of it.

Regions come in a few varieties. There are simple regular shapes (box, ellipsoid), simple irregular shapes (polygon), as well as compound regions made by combining other regions (simple or compound). For example unions or intersections of regions. In addition, the simple regions can be pixel-coordinate based or world-coordinate based. However, a compound regions must always comprise world-coordinate based regions.

It is the task of the Regionmanager (tool) to manage your regions; it creates, lists and manipulates them. Apart from a Regionmanager, the only other way to create a region-of-interest is with the view function of the Image tool. This allows you to interactively make, with the cursor and an image display, a box or polygonal region. The region so made may be collected by the Regionmanager (in fact the complete process can be initiated from the Regionmanager).

Like many other AIPS++ tools, a default Regionmanager is created for you when you load its definition script (in this case, regionmanager.g). This pre-made tool is called drm (the defaultregionmanager). If it doesn't exist for some reason (you might have boldly deleted it maybe), then you can make another one with the regionmanager constructor. Although a Regionmanager tool can contain a little bit of state, in general, you should not need to make another one. The drm Regionmanager should be all you will ever need.

The Regionmanager has a custom GUI and a command line interface. It is probably fair to say that for the simplest regions, such as a pixel box, there is little to choose between making a region with the GUI or the command line interface. However, for the world regions, the GUI is significantly better; it hides the details of handling the coordinate information. Of course, the GUI always invokes the Regionmanager command line for you when it actually creates a region, but it removes the need for you to learn the Regionmanager command-line syntax.



Example

- include 'regionmanager.g'   
- drm.type()
regionmanager
- drm.gui()

This loads the Regionmanager code and creates the default Regionmanager tool called drm - we ask it to tell us its type. We then start the GUI interface from which you can create all sorts of regions.


Simple Pixel and World regions

Pixel regions are specified purely in terms of pixel coordinates. Pixel coordinates are considered to run from 1 at the bottom left corner (blc) of an image to the image shape at the top-right corner (trc).

For example, a pixel region might be box that runs from a bottom blc of [20,20] to a trc of [64,90] in a 2D image. Such a pixel region isn't very portable though. If you are interested in a region about an astronomical source, pixel coordinates are only useful for a particular image. In some other image, the source may well have different pixel coordinates. However, pixel regions are easy to make and use, and they are valuable.

So far, we have been talking about absolute pixel coordinates. However, pixel regions can also be specified in relative coordinates. This is controlled through the argument absrel which may take on one of the values 'abs' (absolute coordinates) 'relref' (relative to the reference pixel of the image) or 'relcen' (relative to the center of the image). You can use the summary function of the image module to see the reference pixel. The sense of the offset is rel = abs - ref.

You may also define pixel regions in terms of fractional coordinates that are in the range [0,1] (blc to trc) - look for argument frac which can be T or F.



Example

From the command line, let's make a simple pixel box region with

- include 'image.g'
- im := image('myimage')
- const r1 := drm.box([10,20], [30,40])
- im.statistics(region=r1);
Number points = 441           
Flux density  = -1.728895e-03  Jy
Sum           = -9.522969e-02        Mean     = -2.159403e-04 
Variance      = 1.518106e-05         Sigma    = 3.896288e-03  
Rms           = 3.897854e-03  

Minimum value is -1.208747e-02  at [13, 32]
Maximum value is 1.287862e-02   at [21, 37]

You have now created a region tool called r1 which describes a 2D box from [10,20] to [30,40] in absolute pixel coordinates. The region is made `const' as well so you can't overwrite it by mistake. We then apply it to an image and evaluate some statistics in the specified region.

World regions are more portable, because the region is specified in terms of the world coordinates (e.g. RA and DEC, or Frequency etc). You can create a region with one image, and then apply it to another. For example, you might make an RA/DEC world box (i.e. blc and trc specified as RA and DEC) from an optical image. You can then apply it to your radio image of the same source. The software will look for the axes of the region (in this case RA and DEC), and then attempt to match them with the image to which the region is being applied. For each matching axis, the region will be applied.

When you make a world region, you must supply a coordinate system from the appropriate image to the Regionmanager. The coordinate system comes from the coordsys function of the Image tool. This associates a coordinate with an image pixel axis.



Example

Now let's do an example with a world box. It is better to make a world box with the GUI interface, but the command line is still manageable. First, let's summarise the header of our image (the one from the previous example).

- im.summary()
Image name       : myimage
Image mask       : Absent
Image units      : JY/BEAM

Direction system : J2000

Name             Proj Shape Tile   Coord value at pixel    Coord incr Units
--------------------------------------------------------------------------- 
Right Ascension   SIN   155  155  17:42:29.303    90.00 -1.000000e+00 arcsec
Declination       SIN   178  178 -28.59.18.600    90.00  1.000000e+00 arcsec
T

And now we will make a world box with the blc at the reference pixel, and the trc offset somewhat. We use the quanta module for creating associations of values and units; in particular, the quantity function.

- blc := "17:42:29.303 -28.59.18.600";
- trc := "17:42:28.303 -28.59.10.600";
- cs := im.coordsys();
- r2 := drm.wbox(blc=blc,trc=trc,pixelaxes=[1,2],csys=cs)
- im.boundingbox(r2);
[blc=[90 90], trc=[103 98], inc=[1 1], bbShape=[14 9], regionShape=[14 9], imageShape=[155 178] ]

As well as passing in the blc and trc quantities, we have also told the Regionmanager which axes these values pertain to (they are numbered in the order in which they list in the summary function), and we also supplied it with the coordinate system of the image. Note that with the GUI interface, both the axes and coordinates handling is hidden. If we wished, we could have left out the pixelaxes argument, and it would have defaulted to [1,2]. We then find the bounding box of the region when it's applied to the image. You can see from the summary listing that it is correct.

We generally think about our images as being in some order; RA/DEC/Frequency (for the X/Y/Z axes) or whatever. In the application of world regions, this order is unimportant. If you created a world region from an RA/DEC/Frequency image, and applied it to a Frequency/DEC/RA image, that image order change would be accounted for.

Note that for pixel regions however, because the region is not tagged with a coordinate system, the order of the image is relevant. So if you made a pixel box from an RA/DEC image and applied it to a DEC/RA image, the RA range would be applied to the DEC axis and vice versa. The pixel regions are present for simple usage only.


Simple pixel regions as world regions

It is possible to create a world region that is specified in pixel coordinates, and maintains the coordinate axis information. This can be useful because compound regions must be world-coordinate based (see below). We can do this because we have defined in the Regionmanager, ``new'' world units, called ``pix'' (which are then known to the quanta module) and these are recognized in the world regions. Similarly, we have also defined a unit called ``frac''. This allows you to specify world coordinates in units of the fraction of the image shape. For example, the blc of a 2D image has a ``frac'' coordinate of [0,0] and a trc of [1,1]. But be careful, such regions are still not portable to another image, because they are still effectively pixel-coordinate based (although masquerading as world regions).



Example
- include 'image.g'
- im := imagemaketestimage()
- cs := im.coordsys()
- local statsout
-
- r1 := drm.wbox(blc="10pix 10pix", trc="20pix 20pix", pixelaxes=[1,2], csys=cs)     
- im.stats(statsout, region=r1, list=F, verbose=F)
- print 'Region 1: number points = ', statsout.npts
Region 1: number points =  121
-
- r2 := drm.wbox(blc="30pix 30pix", trc="40pix 40pix", pixelaxes=[1,2], csys=cs)  
- im.stats(statsout,region=r2, list=F, verbose=F)
- print 'Region 2: number points = ', statsout.npts
Region 2: number points =  121

- r3 := drm.union(r1, r2)
- im.stats(statsout,region=r3, list=F, verbose=F)
- print 'Region 3: number points = ', statsout.npts
Region 3: number points =  242
In this example we make a union from two pixel boxes masquerading as world regions.


Compound Regions

It is often desirable to make a region which combines others. These are called compound regions. For example, give me the union of 3 regions and evaluate the statistics in that union. Or intersect this region with that one and extract the image data from that region. Compound regions are fully recursive; you can nest them as deeply as you like.

You must be aware that compound regions can only be made from world-coordinated based regions. You can convert a pixel region to a world region with the pixeltoworldregion function (at some point this function will move from the Regionmanager tool to the Image tool.



Example
- include 'image.g'
- im := image('myimage')
- cs := im.coordsys();
-
- blc := "17:42:29.303 -28.59.18.600";
- trc := "17:42:28.303 -28.59.10.600";
- r1 := drm.wbox(blc=blc,trc=trc,pixelaxes=[1,2],csys=cs);
- im.boundingbox(r1);
[blc=[90 90] , trc=[103 98] , inc=[1 1], bbShape=[14 9] , regionShape=[14 9] , imageShape=[155 178] ] 
-
- blc := "17:42:27.303 -28.59.18.600"
- trc := "17:42:26.303 -28.59.10.600"
- r2 := drm.wbox(blc=blc,trc=trc,pixelaxes=[1,2],csys=cs);
- im.boundingbox(r2);
[blc=[116 90] , trc=[129 98] , inc=[1 1], bbShape=[14 9] , regionShape=[14 9] , imageShape=[155 178] ] 
-
- r3 := drm.union(r1, r2)
- bb := im.boundingbox(r3)
- bb
[blc=[90 90] , trc=[129 98] , regionShape=[40 9] , imageShape=[155 178] ] 
- 
- im.statistics(statsout=s1, region=r1, list=F)
- im.statistics(statsout=s2, region=r2, list=F)
- im.statistics(statsout=s3, region=r3, list=F)
-
- np := prod((bb.trc-bb.blc) + 1)
- print 'No. pts in r1, r2, r3, bb(r3) = ', s1.npts, s2.npts, s3.npts, np
No. pts in r1, r2, r3, bb(r3) =  126, 126, 252, 360
-
-
- r4 := drm.intersection(r1, r2)
- im.boundingbox(r4)
SEVERE: Caught an exception! Event type=run exception=LCIntersection::LCIntersection - regions do not overlap
<fail>: Caught an exception! Event type=run exception=LCIntersection::LCIntersection - regions do not overlap
        File:   servers.g, Line 711
        Stack:  .()
                .()

We make two world boxes and find the bounding box of each. Then we make the union of the two and find the bounding box which we store in a record called bb. We then find some statistics over the union (supressing the logger output and storing the results in records). You can see that the number of points found in the region reflects the union, not the bounding box of the union.

Finally, we make an intersection of our two regions. Because regions r1 and r2 don't actually intersect, a fail is generated.


Automatic Extension

One general philosophy behind the regions software is that if you apply a region to an image, but don't explicitly specify all axes of the image in that region, then for the unspecified axes, the full image is taken. For example, imagine that from an RA/DEC optical image you made an RA/DEC region. You then applied it to an RA/DEC/Frequency radio spectral-line cube. For the frequency axis, which was not specified in the region, all frequency pixels would be selected. Another philosophy is that if, on application, a region extends beyond an image, the overhang is silently discarded.



Example
- include 'image.g'
- im.boundingbox(drm.box())
[blc=[1 1] , trc=[155 178] , inc=[1 1] , bbShape=[155 178] , regionShape=[155 178] , imageShape=[155 178] ] 
-
- r1 := drm.box(trc=[20])
- im.boundingbox(r1)
[blc=[1 1] , trc=[20 178] , inc=[1 1] , bbShape=[20 178] , regionShape=[20 178] , imageShape=[155 178] ] 
-
- r2 := drm.box(trc=[9000,20])
- im.boundingbox(r2)
[blc=[1 1] , trc=[155 20] , inc=[1 1] , bbShape=[155 20] , regionShape=[155 20] , imageShape=[155 178] ]
In the first example, none of the axes are specified in the region, so it defaults to the full image on all axes. In the second example only the first axis of the trc is specified. In the third example, the trc overhang is silently discarded.


Masks

Some comment about the combination of masks and regions-of-interest is useful here. Consider a simple polygonal region. This region-of-interest is defined by a bounding box, the polygonal vertices, and a mask called a region mask. The region mask specifies whether a pixel within the bounding box is inside or outside the polygon. For a simple box region-of-interest, there is obviously no need for a region mask.

Now imagine that you wish to recover the mask of an image from a polygonal region-of-interest. Now, necessarily, the mask is returned to you in regular Boolean array. Thus, the array shape reflects the bounding-box of the polygonal region. If the actual pixel mask that you apply is all good, then the retrieved mask would be good inside of the polygonal region and bad outside of it. If the actual pixel mask had some bad values in it as well, the retrieved mask would be bad outside of the polygonal region. Inside the polygonal region it would be bad if the pixel mask was bad.

More simply put, the mask that you recover is just a logical ``and'' of the pixel mask and the region mask; if the pixel mask is T and the region mask is T then the retrieved mask is T (good), else it is F (bad).


Vector Inputs Many of the functions of a Regionmanager tool take vectors (numeric or strings) as their arguments. We take advantage of the weak-typing of Glish so that whenever you see a vector argument, you can assume that

  • For numeric vectors you can enter as an actual numeric vector (e.g. [1.2,2.5,3]), a vector of strings (e.g. "1.2 2.5 3") or even a string with white space and/or comma delimiters (e.g. '1.2 2.5 3').
  • For string vectors you can enter as an actual vector of strings (e.g. "1.2 2.5 3") or a string with white space or comma delimiters (e.g. '1.2 2.5 3').


Default Values

When specifying blc and trc vectors for some regions (usually boxes) it can often be that you wish to default one axis but not others.

For example, you may like to specify a value for axis 2 but leave axis 1 at its defaults. This is trivial with the GUI interface (you just don't fill in the ones you are not interested in), but is a little harder with the command line interface.

There is a function called def which provides a magic value which can be used with pixel boxes for this purpose. For example,

- include 'image.g'
- im := imagefromshape('x', [20,30])
- r := drm.box(trc=[drm.def(),10])
- im.boundingbox(region=r)
[blc=[1 1] , trc=[20 10] , regionShape=[20 10] , imageShape=[20 30] ]

You can see that trc for axis 1 has defaulted to the shape of the image.

With world boxes, we have defined a magic unit called `default' instead. Thus

- include 'image.g'
- im := imagefromshape('x', [20,30])
- cs := im.coordsys()
- drm.setcoordinates(cs)
- r := drm.wbox(trc="0default 10pix")
- im.boundingbox(region=r)
[blc=[1 1] , trc=[20 10] , regionShape=[20 10] , imageShape=[20 30] ]


Error Checking

Although some error checking is done when the region is created (in Glish), much of it does not occur until the region is applied to an image. In particular, when creating a compound region, it is not checked that the compound region contains only world-coordinate based regions. It is only when you apply the region to an image that these checks are made. Any errors that occur will cause exceptions to be generated. Thus, when writing scripts, you should always check for a fail status for when you create or use a region.



Example
- include 'image.g'
- const stats := function(imageobject, blc, trc, csys)
{
  r1 := drm.wbox(blc=blc, trc=trc, csys=csys)
  if (is_fail(r1)) {
     fail
  } else {
     ok := imageobject.statistics(r1)
     if (is_fail(ok)) fail
  }
  return T
}

This fairly silly function detects if the region creation failed (function given invalid arguments) and then returns the fail condition. Statistics are then evaluated from the region and the fail status checked in case the region was invalid (or for some other reason).


How a region is stored in AIPS++

In AIPS++, a region is itself a tool. This means it can be transmitted about with Glish, saved to Tables and restored from Tables. A region is actually made with a generic container tool called an Itemcontainer. The regions that you make have no intelligence regarding their greater purpose in life (although of course an Itemcontainer does have some functions so they can be manipulated). All of the knowledge about regions lies with the Regionmanager. It should not be necessary for you to know anything about the insides of regions. You just need to know how to make them and how to use them.



Example
- include 'regionmanager.g'
- r1 := drm.box([10,20], [30,40])
- r1
[gui=<function>, done=<function>, get=<function>, has_item=<function>, 
makeconst=<function>, length=<function>, names=<function>, 
set=<function>, add=<function>, fromrecord=<function>, 
torecord=<function>, type=<function>] 
-
- r1.torecord()
[name=LCSlicer, isRegion=3, blc=[10 20] , trc=[30 40] , inc=, 
fracblc=[F F] , fractrc=[F F] , fracinc=, 
arblc=[1 1] , artrc=[1 1] , oneRel=T, comment=]

First we make a simple pixel-box region. We then try and print it out. You can see that this isn't very enlightening. What you see are the Itemcontainer's functions. To really see inside the region, you can use the torecord function of Itemcontainer. This converts the region to a record which you can view if you really must (the record is of no functional use when dealing with regions).

It's perhaps also useful to know that once you have made a region you can't break it. You can delete it, maybe overwrite it, but you can't mistakenly mess up its contents (even with the Itemcontainer functions with which it was created). It's user proof - there's a challenge for you !




next up previous contents index
Next: regionmanager - Constructor Up: images - Module Previous: imagetest - Function   Contents   Index
Please send questions or comments about AIPS++ to aips2-request@nrao.edu.
Copyright © 1995-2000 Associated Universities Inc., Washington, D.C.

Return to AIPS++ Home Page
2006-08-01