diff options
Diffstat (limited to 'noao/nproto/ace/doc')
| -rw-r--r-- | noao/nproto/ace/doc/detect.hlp | 470 | ||||
| -rw-r--r-- | noao/nproto/ace/doc/installation.hlp | 208 | ||||
| -rw-r--r-- | noao/nproto/ace/doc/objmasks.hlp | 710 |
3 files changed, 1388 insertions, 0 deletions
diff --git a/noao/nproto/ace/doc/detect.hlp b/noao/nproto/ace/doc/detect.hlp new file mode 100644 index 00000000..ac18c675 --- /dev/null +++ b/noao/nproto/ace/doc/detect.hlp @@ -0,0 +1,470 @@ +.help detect Sep00 ace +.ih +NAME +detect -- detect and catalog objects in images +.ih +SYNOPSIS +.ih +USAGE +detect images objmasks catalogs +.ih +PARAMETERS +.ls images +List of images containing objects to be detected. The images should generally +have read and write permission to allow addition of header information. +However, the task will still run without write access with the consequence +that the header will not be updated. +.le +.ls masks = "!BPM" +List of bad pixel masks for the images. This may consist of no bad pixel +mask specified as the empty string "", a single bad pixel mask to apply to +all images or a list of bad pixel masks which must match the images list. +Mask names beginning with "!" are image header keywords which point to the +bad pixel mask. +.le +.ls skys = "SKYFIT" +List of sky images, constant values, sky fit names, or keyword +indirection. If only one value is specified then it applies to all input +images otherwise the list must match the images list. Values beginning +with "!" specify image header keywords containing the image name, constant +value, or sky fit name to be used. The value is first checked to see +if an image with that name exists, then if sky fit keywords are in the +header, and finally if it is a number. Sky fit keywords are formed from +the sky fit name with two digit sequence numbers and are interpreted as +surface fit coefficients. + +If none of these are found then the value is treated as the sky fit name +to be used to save sky fitting performed by this task. +.le +.ls sigmas = "SKYSIG" +List of sky sigma images, constant values, sigma fit names, or keyword +indirection. If only one value is specified then it applies to all input +images otherwise the list must match the images list. Values beginning +with "!" specify image header keywords containing the image name, constant +value, or sigma fit name to be used. The value is first checked to see +if an image with that name exists, then if sigma fit keywords are in the +header, and finally if it is a number. Sigma fit keywords are formed from +the sigma fit name with two digit sequence numbers and are interpreted as +surface fit coefficients. + +If none of these are found then the value is treated as the sigma fit name +to be used to save sky fitting performed by this task. +.le + +The following parameters specify the output. +.ls objmasks +List of output object masks. If no list is given then no object masks +will be created. Otherwise there must be one object mask name for each +input image. The object mask name will be recorded in the input image +header and in any output catalog. +.le +.ls catalogs +List of output catalogs. If no list is given then no catalogs will be +created. Otherwise there must be one catalog name for each input image. +The catalog name will be recorded in the input image header and in any +object mask. The catalog is created as a "table" (see \fBtables\fR +for information about the tables and general tools to interact with the +tables). If the name has an explicit ".fits" extension then a FITS binary +table is created otherwise an IRAF table (".tab" extension) is created. +.le +.ls logfiles = "STDOUT" +List of output log files. If no list is given then no output log information +will be produced. If only one file is specified it applies to all input +images otherwise the list of files must match the images list. Note that +the special name "STDOUT" corresponds to terminal output. +.le + +The following parameters define the initial sky fit determination. This is +only done if no sky image or sky constant value and sigma image or sigma +constant value are specified. +# Sky +.ls newsky = no +Determine new sky fit if one already exists? When the specified sky +corresponds to an existing sky fit (the sky fit coefficients are in the +image header) then this parameter is used to override that fit with a new +fit. Otherwise the fit is used and the initial sky fitting is skipped. +The sky fitting is also skipped if the specified sky is an image or +constant. +.le +.ls nskylines = 100 +Number of sky sample lines to use. This number of lines spread evenly +through the image are used to determine the initial sky fit. +.le +.ls skyblk1d = 10 +Sky block size for 1D sky estimation. +.le +.ls skyhclip = 2. +High sky clipping during 1D sky estimation +.le +.ls skylclip = 3. +Low sky clippling during 1D sky estimation +.le +.ls skyxorder = 4 +Sky fitting x order +.le +.ls skyyorder = 4 +Sky fitting y order +.le +.ls skyxterms = "half" (none|half|full) +Sky fitting y order +.le + +# Iterated Sky +.ls skyupdate = no +Update sky after detection iterations? +.le +.ls niterate = 1 +Maximum number of sky iterations +.le +.ls skyblk2d = 50 +Sky block size during detection +.le +.ls maxskyres = 0.2 +Maximum sky residual for iteration +.le + +# Detection +.ls convolve = "block 3 3" +Convolution kernel +.le +.ls hsigma = 3. +Sigma threshold above sky +.le +.ls lsigma = 10. +Sigma threshold below sky +.le +.ls hdetect = yes +Detect objects above sky? +.le +.ls ldetect = yes +Detect objects below sky? +.le +.ls minpix = 10 +Minimum number of pixels in detected objects +.le +.ls sigavg = 4. +Sigma of mean flux cutoff +.le +.ls sigmax = 4. +Sigma of maximum pixel +.le +.ls bpval = 1 +Output bad pixel value +.le + +# Splitting" +.ls split = yes +Split objects? +.le +.ls splitmax = INDEF +Maximum sigma above sky for splitting +.le +.ls splitstep = 0.4 +Splitting steps in convolved sigma +.le +.ls splitthresh = 5. +Splitting threshold in sigma +.le +.ls sminpix = 10 +Minimum number of pixels in split objects +.le +.ls ssigavg = 10. +Sigma of mean flux cutoff +.le +.ls ssigmax = 5. +Sigma of maximum pixel +.le + + +# Growing" +.ls ngrow = 2 +Number of grow rings +.le +.ls agrow = 2. +Area grow factor +.le +.ih +DESCRIPTION + +SKY DETERMINATION + +A critical part of detecting objects in astronomical images is determining +the background sky and sky sigma at each point in the image. In the +following discussion sky means both the mean sky level and the sky sigma. +\fBDetect\fR provides for either the user to specify the sky or for the +task to use a sky fitting algorithm. The user may specify a sky either as +another image or as a constant value. Note that the image name or +value may be specified either explicitly or with a keyword associated +with the image. + +If the sky is not specified by an image or constant value then a surface +fit to the sky is used. The surface fit is recorded in the image header as +a sequence of keywords with a specified name (the keyword prefix which may +be up to six characters) and two digit sequence number. The values of the +keywords contain the coefficients of the fit. The the surface fit +coefficients are defined in the SURFACE FIT section. + +Note that it is possible to specify the mean sky and the sky sigma in +different ways. When one is given as an image or constant and the other +as a fit. The one given as an image or constant will be kept fixed and +the fit determination and updating will be done only on the other. + +The sky surface fit is computed in two stages. There is an initial +determination using a subsample of image lines. Then there is an +optional update of the sky sample during the object detection step. +The detection step with sky updating may be iterated a specified number +of times until the maximum difference in the mean sky is less than some +amount. + +INITIAL SKY DETERMINATION + +If an existing surface fit is specified then the parameter \fInewsky\fR +selects whether a new surface fit is to be computed. If the value is "no" +then the initial sky determination is skipped though the detection update +may still be selected. + +The initial sky fit uses a combination of block averaging to reduce the +number of points in the fitting, one dimensional line fitting with sigma +clipping rejection to eliminate objects, and finally fitting a two +dimensional surface to the set of block averages over all the sample lines +which cover the image. + +The parameter \fInskylines\fR defines the number of sample lines across +the image to be used. The lines are evenly spaced starting with the +first line and ending with the last line. The number of lines affects +how fast the sky estimation is done. + +The pixels from the input line are initially all given unit weight. Bad +pixels identified by the input bad pixel mask are excluded by setting their +weights to zero. A weighted block average, with the weight of each block +being the sum of the weights, is computed. The size of the blocks is given +by the \fIskyblk1d\fR parameter. This is done to speed the fitting by +reducing the number of points. Note that when all pixels in a block have +zero weight due to the bad pixel mask or subsequent rejection the weight of +the composite block average point is zero. + +If only one of sky mean and sky sigma quantities is being determined with +the other quantity given by an input image, constant, or previous fit +then those values are simple block averaged with the same block size +to produce sample points for the mean sky or sky sigma. Note that the +sky sigma of the sample points also requires division by the square root +of the block size to give the sky sigma per block average point. The +line fitting described next is then skipped for this quantity. + +The weighted one dimensional line fitting to the block averages uses +Chebyshev polynomials of order given by the \fIskyxorder\fR. Note that +this order is the number of polynomial terms, which is one higher than the +maximum power of the polynomial so that a value of 3 corresponds to a +quadratic polynomial. + +When the mean sky is being determined, the line fitting is performed and +the fitted values at the block centers are evaluated. + +When the sky sigma is being determined, the absolute value of the residuals +relative to the mean sky divided by 0.7979 are computed. A gaussian noise +distribution will have a mean value of this quantity equal to the sigma of +the distribution. In other words, the mean of the absolute deviations of a +gaussian distribution is 0.7979 times sigma. By fitting a function to +these residual values a position variable estimate of the sky sigma is +obtained without needing to compute standard deviations over some set of +points. The fitted values at the block centers are evaluated to give the +sky sigmas for the block averaged data. + +With the set of block averaged data points and estimated mean skys and sky +sigmas points that deviate by more than the number of sigma given by the +\fIskyhclip\fR and \fIskylclip\fR parameters are rejected by setting their +weights to zero. The line fitting is then repeated until no points are +rejected with a maximum of 10 iterations. + +When the iteration completes the block average points for that image line +are accumulated for a two dimensional surface fit. Note that the weights +are used to exclude rejected averages and to weight blocks that had fewer +points due to bad pixels. The surface fit is a two dimensional Chebyshev +polynomial of orders given by the \fIskyxorder\fR and \fIskyyorder\fR. The +orders have the same meaning as in the one dimensional polynomial, namely +the number of terms in powers of x and y. There are also cross terms which +are a mixture of powers of both x and y. The \fIskyxterms\fR select +whether to use any cross terms, only cross terms whose total power does not +exceed the maximum of the pure x and y terms, or all combinations of +powers. + +After all the sample lines are completed the final surface fits are +computed. The coefficients of the fits are written to the image header +under the specified sky fit names and the fits are passed on to the +detection phase. Note that if the input image is read only then the +fit will not be written to the header but the task continues. + +UPDATED TO SKY DURING DETECTION + + +DETECTION + +The detection of objects in an image is conceptually quite simple. Each +pixel is compared against the expected sky at that point and if it is +more that a specified number of sky sigma above the sky it is a candidate +object pixels. Candidate object pixels are grouped into objects on the basis +of being connected along the eight neighboring directions. The candidate +object is then accepted if it satisfies the criteria of a minimum +number of pixels, a sufficiently significant maximum pixel, and a sufficiently +significant flux above sky. + +To detect faint objects where individual pixels are not significantly above +the sky but all pixels taken together are significant a detection filter is +applied. This consists of applying a convolution function to the image and +performing the detection described in the previous paragraph on the +convolved pixels with the sky sigma suitable adjusted for the convolution. +The convolution acts as an optimizing filter for objects with shapes +corresponding to the convolution weights. The remaining discussion +is in terms of the convolved pixel values. The case of no convolution +can be thought of as a convolution with a delta function though the +implementation is not done as a convolution for efficiency. + +Two other options to the detection are to also find pixels that are +significantly below sky (using an independent threshold to that used for +detecting pixels above sky) and form them into "dark" objects and to +take the remaining pixels that are not significantly above or below the +sky and use them to define a sky sample for output or for updating the +initial sky. + +We now go into more detail. The background sky and sky sigma against which +the detection is performed is initially set as described earlier. If desired +the sky pixels may be accumulated to update the sky. After updating the +sky the detection step may be repeated using the new sky. This is +discussed futher when we reach the end of the detection step description. + +The convolution is specified by the \fIconvolve\fR parameter. The values for +this parameter and the definition of the convolution are given in the +CONVOLUTION DETECTION FILTER section. The input pixel data is convolved +and the sky sigma is appropriately adjusted. + +When the central pixel in the convolution is flagged as a bad pixel by the +bad pixel mask (any non-zero value is a bad pixels) then the convolved +value is considered to be a bad pixel. If an output object masks is +specified the pixel will be marked with the value specified by the +\fIbpval\fR parameter. The value may be set to not show the bad pixel in +the object mask, to set all input bad pixels to some value, or to pass the +input bad pixel value to the object mask. Note that bad pixel masks in the +object mask must be between 1 and 10 to avoid confusion with the values +used to identify objects. If other pixels in the convolution are flagged +as bad pixels they are excluded from the convolution and the +convolved sky sigma is adjusted but the convolution value is still used +as a valid image pixel for detection. + +The sigma threshold for pixels to be detected as part of an object above +sky is given by the \fIhsigma\fR. This number is multiplied by the sky +sigma to get the deviation from sky. As noted earlier the sky sigma is +for the convolved pixels and the + + +CONVOLUTION DETECTION FILTER + + +The convolution detection filter is specified with the \fIconvolve\fR +parameter. There is only one convolution that can be specified and it applies +to all input images in a list. If a null string ("") is specified +then no convolution is performed. The task has been optimizations for +this case to avoid treating this as a 1x1 convolution and to avoid extra +memory allocations required when a convolution is done. + +The convolved value at pixel (i,j), denoted I(i,j), within an image of size +CxL is defined by + +.nf + I_convolved(i,j) = sum_kl{I_unconvolved(m,n)*W(k,l)} / sum_kl{W(k,l)} +.fi + +where I(m,n) is the unconvolved value at pixel (m,n), W(k,l) are the NX x +NY (both must be odd) convolution weights, sum_kl is the double sum over k +and l, and + +.nf + m' = i + k - (NX+1)/2 for k = 1 to NX + n' = j + l - (NY+1)/2 for l = 1 to NY + + m = m' (1<=m'<=C) m = 1-m' (m'<1) m = 2C-m' (m'>C) + n = n' (1<=n'<=L) n = 1-n' (n'<1) n = 2L-n' (m'>L) +.fi + +The last two lines represent boundary reflection at the edges of the image. + +The sky sigma of a convolved pixel is approximated by + +.nf + sigma_convolved(i,j) = sigma_unconvolved(i,j) / sum_kl{W(k,l)} +.fi + +In the presence of bad pixels identified by a bad pixel mask the convolution +weight applied to a bad pixel is set to zero. The sum of the weights +used to normalize the convolution is then modified from the situation with +no bad pixels. This will correct the convolved pixel value for the missing +data and the estimated sky sigma is appropriately larger. + +A convolution can be computational slow, especially for larger sizes. +The implementation of the convolution has been optimized to recognize +bilinear symmetries or lines which are scaled versions of other lines. +So if possible such symmetries should be used. The "block", "bilinear", +and "gauss" special convolutions described below have such symmetries. + +There is also an overhead in checking for bad pixels. The convolution +has an optimization to avoid such checks in the case where no bad pixel +mask is specified. + +The \fIconvolve\fR parameter is a string which can take one of the +following forms. + +.ls "" +There is no convolution or, equivalently, NX=1, NY=1. +.le +.ls @[filename] +The weights are given in the specified file. The format consists of lines +of whitespace separated values. The number of values on each line must be +the same and defines NX and the number of lines defines NY. +.le +.ls block [NX] [NY] +The weights are all the same and the convolution size is given by the +two numbers following the word "block". +.le +.ls bilinear [NX] [NY] +The weights are the bilinear matrix product of triangular one dimensional +matrices of sizes given by the two numbers following the word "bilinear". +The weights are described by the matrix product relation + +.nf + [1 ... (NX+1)/2 ... 1] * Transpose{[1 ... (NY+2)/2 ... 1]} +.fi + +For example for NX=5, and NY=3 the weights would be + +.nf + 1 2 3 2 1 + 2 4 6 4 2 + 1 2 3 2 1 +.fi +.le +.ls gauss [NX] [NY] [SX] [SY] +The weights are bidimensional gaussian values on a grid of size NX by NY +with sigma values SX and SY (real numbers) in units of pixel spacing. +.le +.ls [W(1,1)] ... [W(NX,1)], ..., [W(1,NY)] ... [W(NX,NY)] +The weights are specified as a string of real values. The values are +whitespace separated within each line and the lines are delimited by +comma. For example + +.nf + 1 2 1 + 1 2 1, 2 3 2, 1 2 1 ==> 2 3 2 + 1 2 1 +.le + +When a logfile is defined the weights are included in the log output. + + +OBJECT MASKS + +.ih +EXAMPLES +.ih +REVISIONS +.ih +SEE ALSO +.endhelp diff --git a/noao/nproto/ace/doc/installation.hlp b/noao/nproto/ace/doc/installation.hlp new file mode 100644 index 00000000..c399ad4c --- /dev/null +++ b/noao/nproto/ace/doc/installation.hlp @@ -0,0 +1,208 @@ +.help installation Jan01 ace + +.ce +\fBACE: Astronomical Cataloging Environment\fR +.ce +Release Notes and Installation Instructions + +.sh +SUMMARY +The ACE external package is used to catalog objects in images and manipulate +the catalogs. + +.sh +RELEASE INFORMATION +The following summary only highlights the major changes. There will also +be minor changes and bug fixes. + +.ls V0.2: January 27, 2001 +Alpha test version. +.le +.sh +INSTALLATION INSTRUCTIONS +Installation of this external package consists of obtaining the files, +creating a directory containing the package, compiling the executables or +installing precompiled executables, and defining the environment to load +and run the package. The package may be +installed for a site or as a personal installation. If you need help with +these installation instructions contact iraf@noao.edu or call the IRAF +HOTLINE at 520-318-8160. +.ls [arch] +In the following steps you will need to know the IRAF architecture +identifier for your IRAF installation. This identifier is similar to the +host operating system type. The identifiers are things like "ssun" for +Solaris, "alpha" for Dec Alpha, and "linux" or "redhat" for most Linux +systems. The IRAF architecture identifier is defined when you run IRAF. +Start the CL and then type + +.nf + cl> show arch + .ssun +.fi + +This is the value you need to know without the leading '.'; i.e. the +IRAF architecture is "ssun" in the above example. +.le +.ls [1-site] +If you are installing the package for site use, login as IRAF +and edit the IRAF file defining the packages. + +.nf + % cd $hlib +.fi + +Define the environment variable ace to be the pathname to +the ace package root directory. The '$' +character must be escaped in the VMS pathname and UNIX pathnames must be +terminated with a '/'. Edit extern.pkg to include the following. + +.nf + reset ace = /local/ace/ + task ace.pkg = ace$ace.cl +.fi + +Near the end of the hlib$extern.pkg file, update the definition of +helpdb so it includes the ace help database, copying the syntax +already used in the string. Add this line before the line +containing a closing quote: + +.nf + ,ace$lib/helpdb.mip\ +.fi +.le +.ls [1-personal] +If you are installing the package for personal use define a host +environment variable with the pathname of the directory where the package +will be located (needed in order to build the package from the source +code). Note that pathnames must end with '/'. For example: + +.nf + % setenv ace /local/ace/ +.fi + +In your login.cl or loginuser.cl file make the following definitions +somewhere before the "keep" statement. + +.nf + reset ace = /local/ace/ + task ace.pkg = ace$ace.cl + printf ("reset helpdb=%s,ace$lib/helpdb.mip\nkeep\n", + envget("helpdb")) | cl + flpr +.fi + +If you will be compiling the package, as opposed to installing a binary +distribution, then you need to define various environment variables. +The following is for Unix/csh which is the main supported environment. + +.nf + # Example + % setenv iraf /iraf/iraf/ # Path to IRAF root (example) + % source $iraf/unix/hlib/irafuser.csh # Define rest of environment + % setenv IRAFARCH ssun # IRAF architecture +.fi + +where you need to supply the appropriate path to the IRAF installation root +in the first step and the IRAF architecture identifier for your machine +in the last step. +.le +.ls [2] +Login into IRAF. Create a directory to contain the package files and the +instrument database files. These directory should be outside the standard +IRAF directory tree. + +.nf + cl> mkdir ace$ + cl> cd ace +.fi +.le +.ls [3] +The package is distributed as a tar archive for the +sources and, as an optional convenience, a tar archive of the executables +for select host computers. Note that IRAF includes a tar reader. The tar +file(s) are most commonly obtained via anonymous ftp. Below is an example +from a Unix machine where the compressed files have the ".Z" extension. +Files with ".gz" or ".tgz" can be handled similarly. + +.nf + cl> ftp iraf.noao.edu (140.252.1.1) + login: anonymous + password: [your email address] + ftp> cd iraf/extern + ftp> get ace.readme + ftp> binary + ftp> get ace.tar.Z + ftp> get ace-bin.<arch>.Z (optional) + ftp> quit + cl> !uncompress ace.tar + cl> !uncompress ace-bin.<arch> (optional) +.fi + +The readme file contains these instructions. The <arch> in the +optional executable distribution is replaced by the IRAF architecture +identification for your computer. + +Upon request the tar file(s) may be otained on tape for a service +charge. In this case you would mount the tape use rtar to extract +the tar files. +.le +.ls [4] +Extract the source files from the tar archive using 'rtar". + +.nf + cl> softools + so> rtar -xrf ace.tar + so> bye +.fi + +On some systems, an error message will appear ("Copy 'bin.generic' +to './bin fails") which can be ignored. +Sites should leave the symbolic link 'bin' in the package root +directory pointing to 'bin.generic' but can delete any of the +bin.<arch> directories that won't be used. If there is no binary +directory for the system you are installing it will be created +when the package is compiled later or when the binaries are installed. + +If the binary executables have been obtained these are now extracted +into the appropriate bin.<arch> directory. + +.nf + # Example of sparc installation. + cl> cd ace + cl> rtar -xrf ace-bin.sparc # Creates bin.sparc directory +.fi + +The various tar file can be deleted once they have been +successfully installed. +.ls [5] +For a source installation you now have to build the package +executable(s). The "tables" package must be installed first if not +already available. First you configure the package for the particular +architecture. + +.nf + cl> cd ace + cl> mkpkg <arch> # Substitute sparc, ssun, alpha, etc. +.fi + +This will change the bin link from bin.generic to bin.<arch>. The binary +directory will be created if not present. If an error occurs in setting +the architecture then you may need to add an entry to the file "mkpkg". +Just follow the examples in the file. + +To create the executables and move them to the binary directory + +.nf + cl> mkpkg -p ace # build executables + cl> mkpkg generic # optionally restore generic setting +.fi + +Check for errors. If the executables are not moved to the binary directory +then step [1] to define the path for the package was not done correctly. +The last step restores the package to a generic configuration. This is not +necessary if you will only have one architecture for the package. +.le + +This should complete the installation. You can now load the package +and begin testing and use. +.endhelp diff --git a/noao/nproto/ace/doc/objmasks.hlp b/noao/nproto/ace/doc/objmasks.hlp new file mode 100644 index 00000000..1c35c4c9 --- /dev/null +++ b/noao/nproto/ace/doc/objmasks.hlp @@ -0,0 +1,710 @@ +.help objmasks Jan02 nproto +.ih +NAME +objmasks -- detect objects in images and create masks and sky maps +.ih +SYNOPSIS +.ih +USAGE +objmasks images objmasks skys +.ih +PARAMETERS +.ls images +List of images or multiextension files for which object masks are desired. +.le +.ls objmasks +List of object masks to be created. This list must match the input list. +Multiextension input files will produce multiextension mask files. If the +input image is writable, the name of the created mask will recorded in the +image header. Note that it is possible to specify a null image to +not produce an output mask. This might be done if the background sky +or sky sigma maps are desired or to just see the log information. +.le + +.ls omtype = "numbers" (boolean|numbers|colors|all) +The type of encoding for the object mask values. In all cases non-object pixels +(that is background) have mask values of zero. The choices for the mask +values are "boolean", "numbers", "colors", and "all". These are described +in the \fIOutput Data\fR section. +.le +.ls skys = "", sigmas = "" +Optional lists of input or output sky and sigma maps. Maps are either +constant values or images which are interpolated to the size of the input +images. If a list is given it must match the input \fIimages\fR list. +If constant values or existing maps are specified then those are used +without change. If a new filename is given then an output file is created +with the values computed by the task. Multiextension input images create +or apply the same extension names to the specified sky or sigma files. +Constant input values apply to all extensions. The sigma values are +per single input image pixel. +.le +.ls masks = "!BPM" +List of bad pixel masks for the input images. Non-zero masks values are +ignored in the object detection and are passed on to the output object +masks based on the \fIomtype\fR parameter. An empty list applies no bad +pixel mask, a single mask applies to all input images, and a matching +list matches the masks with the input image. A mask is specified by a +filename or by reference to a filename given by the value of a header +keyword in the input image. A header keyword reference is made with the +syntax "!<keyword>" where <keyword> is the desired keyword with case +ignored. For multiextension files the input masks may be either a +multiextension file with matching extension names or a directory of +pixel list files with the extension names as filenames. +.le +.ls extnames = "" +Extensions to select from multiextension files. A null string matches all +extension names. Otherwise the parameter is a comma separated list of +patterns that match the entire extension name. Thus, an explicit list of +extension names may be specified or the pattern matching characters '?' for +any character or '[]' for a set of characters may be used. The set may +include ranges in ascii order by using hyphens; i.e. 1-3 matches the +characters 1, 2, and 3. +.le +.ls logfiles = "STDOUT" +List of output log files. If no list is given then no output log information +will be produced. If only one file is specified it applies to all input +images otherwise the list of files must match the images list. Note that +the special name "STDOUT" corresponds to terminal output. +.le + +.ls blkstep = 1 +The mean and sigma of the background or sky pixels are determined in a +first pass through the image. If \fIblkstep\fR is one all lines are used. +To skip lines in order to speed up this computation, the parameter may be +set to a larger value to define the increment between lines. However, the +task will enforce a preset minimum number to insure a sufficient sample. +.le +.ls blksize = -10 +The background mean sky and sky sigma are determined in a set of square +blocks from which the values are linearly interpolated to each point in the +input image. The size of the blocks may be specified as a number of blocks +spanning the smaller image dimension by using a negative integer value. +Or the size may be specified as the number of pixels across a block. +The task will enforce a preset minimum number of pixels per block which may +require using bigger blocks than specified. The background determination +algorithm is described further in the "Background Determination" section. +.le + +.ls convolve = "block 3 3" +Convolution filter to be applied prior to threshold detection. The +convolution filter is defined by a set of weights in a 2D array. These +may be specified in files or with certain forms given by special strings. +The options are described in the "Convolution Filter" section. +.le +.ls hsigma = 3., lsigma = 10. +Object pixels are identified by sigma thresholds about the mean background +based on the estimated background sigma at each point in the image. +The sigma factors are specified in terms of the "per pixel" sigma before +convolution. The \fIhsigma\fR value is the "high" or above background +limit and the \fIlsigma\fR value is the "low" or below background limit. +Typically detections are one-sided, such as detecting objects above +the background, and so the thresholds need not be equal. +.le +.ls hdetect = yes, ldetect = no +Identify objects as pixels which are above the background (\fIhdetect\fR) +and below the background (\fIldetect\fR)? If objects are detected but the +corresponding parameter is no then the output mask will not include those +objects. +.le +.ls neighbors = "8" (8|4) +The threshold selected pixels are associated with other neighboring pixels to +form an object. The criterion for a neighbor being part of the +same object is defined by this parameter. The choices are "8" for +pixels touching in any of the 8 directions or "4" to identify neighbors +as only horizontal or vertically adjacent. +.le +.ls minpix = 6 +The minimum number of neighboring pixels which define an acceptable object. +.le +.ls ngrow = 2, agrow = 2. +After an object is identified as a set of threshold detected pixels, +additional neighboring pixels may be added to the object. This allows +expanding the object into the faint wings of the light distribution. The +additional pixels are those which touch the boundary pixels. Pixels are +added in multiple passes, each time extending the previous boundary. The +parameter \fIngrow\fR (an integer value) defines the maximum number of +boundary extensions. The parameter \fIagrow\fR (a real value) specifies +the maximum increase in area (number of pixels) from the original +detection. +.le +.ih +DESCRIPTION +\fBOBJMASKS\fR is a task for creating masks covering objects in images. +An optional secondary product of this task is to produce background +and sigma maps. Objects are identified by threshold sigma detection. +These object masks may be used by other applications to exclude the object +data or focus on the objects. The detection consists of determining a +smooth, spatially variable mean background and background sigma (if no +input maps are provided), convolving the data by an optional filter to +optimize detection of faint sources, collecting pixels satisfying the +detection thresholds, assigning neighboring pixels to a common object, +applying a minimum number of pixels test to the objects, and growing +objects to extend into the wings of the object light distribution. +The last step is writing out the identified object pixels as a mask. + +1. Input Data + +The input data consists of one or more 2D images. The images are assumed +to contain a moderately smooth background and multiple sources or +objects. This task is most useful for images with large numbers of small +sources rather than one large object such as a nearby galaxy. The input +images, specified by the \fIimages\fR parameter, may be individual images +(which includes images selected from multiextension files as explicit +image extensions) or multiextension files specified by a root filename. In +the latter case the image extension names selected by the \fIextnames\fR +parameter are used. + +Background means and sigmas (specified per image pixels) may be specified +by "maps". These may be constant numerical values or images. The map +images will be linearly interpolated to the size of the input images. +For multi-extension input data, constant map values apply to all extensions +and maps are also multiextension files with map images having the same +extension names. + +Bad pixel masks may be associated with the input images to +exclude pixels from the background and object determinations. These +bad pixels are also included in the output object masks. The bad pixel +masks are specified by the \fImasks\fR parameter. This parameter may +identify a mask by a filename or a keyword. A single mask may be +specified to apply to all images or a matching list of masks may be +given. + +The masks are in one of the supported mask formats. As of IRAF V2.12 this +includes pixel list (.pl) files and FITS "type=mask" extensions. When the +input files are multiextension files, the selected extension names are +appended to the specified mask filename to select masks with the same +extension name. If a mask file of the form "name[ext]" is not found +the task will treat the filename as a directory of pixel list files and +select the pixel list file with the extension name; i.e. "name/ext.pl". + +2. Output Data + +The output of this task are object masks, sky maps, sigma maps, and log +information. The output object masks default to mask type extensions. If an +extension name is not specified explicitly the default extension name +"pl" is created. To select a pixel list output format an explicit ".pl" +extension must be used. + +When the input data are multiextension files, the output masks, mean sky +maps, and sky sigma maps will be multiextension files with the specified +rootnames and the same extension name as the input. + +The output mask values identify non-object pixels with zero. The non-zero +values are encoded as selected by the \fIomtype\fR parameter. The choices +are: + +.ls "boolean" +All object and bad pixels have a mask value of one; i.e. the output masks +consists only of the values 0 and 1. +.le +.ls "numbers" +Input bad pixels values between 1 and 10 preserve their value and all +other input mask values are mapped to 10. The object mask pixels have +object numbers starting with 11. The object numbers are assigned by +the task (roughly in order from the first line to the last line) and +all pixels from a single object have the same unique object number. +.le +.ls "colors" +Input bad pixels are mapped to output values of one. The object numbers +are modulo 8 plus 2; i.e. values between 2 and 9. The purpose of this +numbering is to allow mapping to the nine standard display colors for an +interesting overlay with the \fBdisplay\fR task and "ocolors='+203'". +.le +.ls "all" +This is the same as "numbers" except that bits 24 to 27 in the mask values +are used for various purposes. In particular bit 24 is set for the boundary +pixels. This numbering will be used in the future by special tasks. +.le + +Output mean sky and sky sigma maps consist of the mean and sigma values +in blocks as described in the "Background Determination" section. +Therefore, the size of the map images are smaller than the input data images. +These maps need to be interpolated to the size of the input image +to obtain the values used for particular pixels in the data images. +This interpolation expansion is done automatically by some tasks such +as \fBmscred.rmfringe\fR. + +The log output provides information about the files, the phase of the +processing, some of the parameters, and the convolution filter weights. +The output begins with the task identifier ACE. This is because this +prototype task is a first release piece of a major package called ACE +(Astronomical Cataloging Environment), which is under development. + +3. Background Determination + +Detection of sources in an image begins with determining the background. +By this we mean estimating the probability distribution of the background +pixel values at every pixel in the image. In practice we only estimate +the central value and width and assume a normal distribution for evaluating +the significance of deviations from the central value. Since we normally +won't have a sample of values at each pixel the distribution is +determined from a sample of nearby pixels. + +In this discussion the central value of a distribution is denoted by <I>. +It is estimated by the mean or mode of the sample. The width of the +distribution about <I> is denoted by <S> and is estimated by the absolute +mean residual converted to the standard deviation of a normal distribution +with the same absolute mean residual. The normal deviation of a value I +from the distribution is defined as R = (I - <I>) / <S>. + +The background may be specified by input maps for one or both of the +background quantities. The maps may be constant values which apply +to all pixels or a grid of values given in an image which are linearly +interpolated to the full size of the input data. For those quantities +which are not input the following algorithm is used for computing +a map. The maps may be output and used as a product of this task. + +The background and/or sigma are estimated in two initial passes through the +data. The first pass algorithm fits linear functions to a subsample of +lines using sigma clipping iteration to eliminate objects. The subsample +is used to speed up the algorithm and is reasonable since only linear +functions are used. Each sample line is block averaged in blocks of 10 +pixels and a linear function is fit by least squares to obtain an estimate +for <I> along the line. The fitting weights are the number of good pixels +in each block average after elimination of bad pixels specified by the +user in a bad pixel mask. The absolute values of the residuals are also +fit to produce a constant function for <S>. + +To exclude objects from affecting these estimates the fitting is iterated +using sigma clipping rejection on the normal deviations R. In the +first iteration the fitting function for <S> is a constant and in +subsequent steps a linear fit is used. When the sigma clipping iteration +rejects no more data, the remaining block averages, absolute residuals, and +weights are used to fit a 2D plane for both <I> and <S>. The <S> surface +is a constant in order to avoid potential negative sigma values. + +This first pass algorithm is fast and produces good estimates for the +planar approximation to the background. The second pass divides the image +into large, equal sized blocks, as specified by the \fIblksize\fR +parameter, and estimates <I> and <S> in each block. The size of the blocks +needs to be large enough to give good estimates of the statistics though +small enough to handle the scale of variations in the sky. Each block is +divided into four subblocks for independent estimates which are then +combined into a final value for the block. As with the first pass, the +second pass can be speeded up by using a subsample of lines (parameter +\fBblkstep\fR) provided some minimum number of lines per subblock is +maintained. + +The background estimates in each subblock are made using histograms of the +normal deviations R computed relative to the first pass estimates of <I> +and <S>. When pixels are added into the histogram the <I> and <S> used to +compute R are accumulated into means of these quantities in order +to convert estimates from the normalized deviation histogram back into data +values. The histograms are truncated at +/-2.5 and have bin widths +determined by requiring a specified average bin population based on the +number of pixels in the block. Typically the bin population is of order +500. The histogram truncation is essentially an object-background +discrimination. + +When all the pixels in a subblock have been accumulated, new estimates of +<I> and <S> are computed. If the number of pixels in the histogram is +less than two-thirds of the subblock pixels the estimates are set to be +indefinite. This flags the subblock as too contaminated by objects to be +used. All subblock neighbors, which may cross the full block boundaries, +are also rejected to minimize contamination by the wings of big galaxies +and very bright stars. + +If the histogram has enough pixels, the bin populations are squared to +emphasize the peak of the distribution and reduce the effects of the +truncated edges of the histogram. Because of noise and the fine binning of +the histogram, a simple mode cannot be used and squaring the bin numbers +helps to approach the mode with a centroid. Squaring the bin values and +then computing the centroid can also be thought of as a weighted centroid. + +Generally a mode is considered the best estimate to use for the central +value <I> of the sky distribution. But it is unclear how to best estimate +the mode without an infinite number of pixels. One could do something like +fit a parabola to the histogram peak. But instead we use the empirical +relation for a skewed distribution between the mean, mode, and median; +<I>=mean-3*(mean-median). The mean is the weighted centroid and the median +is obtained numerically from the histogram using linear interpolation to +get a subbin value. + +The <S> values are obtained from the absolute mean residual of the +unweighted histogram about the previously derived central value <I> of the +histogram. The conversion to a standard deviation is made by computing the +ratio between the standard deviation and mean absolute deviation of a +Gaussian distribution. The standard value over the entire distribution +cannot be used because the histogram is truncated. However, it is easy to +numerically compute the ratio with the same truncation. + +Once <I> and <S> are obtained in bin numbers it is converted to data +values by using the mean and sigma of the input pixel values used +to create the histogram. + +The averages of the subblock <I> and <S> values which are not indeterminate +in each block are computed. If any of the full blocks are indeterminate +when all the subblocks have been eliminated as contaminated, values are +obtained for them by interpolation from nearby blocks. The block values +are then linearly interpolated to get background values for every +pixel in the input image. + +Note that the background pixels used in the block algorithm before +detection are derived by simple sigma clipping of the histogram values +around the planar background. If an output map for either the mean +values or the sigmas is specified then during the object detection stage +the background and sigmas are updated using the detected sky pixels about +the initial block sampled background. This is a more sensitive selection +of sky pixels since convolution filtering can exclude pixels from faint +objects and the wings of all objects. The new set of sky pixels are +accumulated and used in the same way as described earlier. + +4. Convolution Filters + +In order to improve the detection of faint sources dominated by the +background noise, the input data may be convolved to produce filtered +values in which the noise has been suppressed. The threshold detection +is then performed on the filtered data values. + +The convolution detection filter is specified with the \fIconvolve\fR +parameter. There is only one convolution that can be specified and it +applies to all input images in a list. If a null string ("") is specified +then no convolution is performed. The task has been optimizations for this +case to avoid treating this as a 1x1 convolution and to avoid extra memory +allocations required when a convolution is done. + +The convolved value at pixel (i,j), denoted I'(i,j), is defined by + +.nf + I'(i,j) = sum_kl{I(m,n)*W(k,l)} / sum_kl{W(k,l)} +.fi + +where I(m,n) is the unconvolved value at pixel (m,n), W(k,l) are the NX x +NY (both must be odd) convolution weights, sum_kl is the double sum over k +and l, and + +.nf + m' = i + k - (NX+1)/2 for k = 1 to NX + n' = j + l - (NY+1)/2 for l = 1 to NY + + m = m' (1<=m'<=C) m = 1-m' (m'<1) m = 2C-m' (m'>C) + n = n' (1<=n'<=L) n = 1-n' (n'<1) n = 2L-n' (m'>L) +.fi + +The size of the image is C x L. The last two lines represent boundary +reflection at the edges of the image. + +The sky sigma of a convolved pixel is approximated by + +.nf + sigma'(i,j) = sigma(i,j) / sum_kl{W(k,l)} +.fi + +In the presence of bad pixels specified in the bad pixel mask the +convolution weight applied to a bad pixel is set to zero. If the central +pixel is bad then the convolved value is also considered to be bad. The +sum of the weights used to normalize the convolution is then modified from +the situation with no bad pixels. This will correct the convolved pixel +value for the missing data and the estimated sky sigma is appropriately +larger. Since there is an overhead in checking for bad pixels the +convolution has an optimization to avoid such checks in the case where no +bad pixel mask is specified. + +A convolution can be computational slow, especially for larger convolution +kernel sizes. The implementation of the convolution has been optimized to +recognize bilinear symmetries or lines which are scaled versions of other +lines. So if possible users should chose convolutions with such symmetries +to be most efficient. The "block", "bilinear", and "gauss" special +convolutions described below all have such symmetries. + +The \fIconvolve\fR parameter is a string with one of the following forms. + +.ls "" +There is no convolution or, equivalently, NX=1, NY=1. +.le +.ls @[filename] +The weights are given in the specified file. The format consists of lines +of whitespace separated values. The number of values on each line must be +the same and defines NX and the number of lines defines NY. +.le +.ls block [NX] [NY] +The weights are all the same and the convolution size is given by the +two numbers following the word "block". This is a moving block average +filter. +.le +.ls bilinear [NX] [NY] +The weights are the bilinear matrix product of triangular one dimensional +matrices of sizes given by the two numbers following the word "bilinear". +The weights are described by the matrix product relation + +.nf + [1 ... (NX+1)/2 ... 1] * Transpose{[1 ... (NY+2)/2 ... 1]} +.fi + +For example for NX=5, and NY=3 the weights would be + +.nf + 1 2 3 2 1 + 2 4 6 4 2 + 1 2 3 2 1 +.fi +.le +.ls gauss [NX] [NY] [SX] [SY] +The weights are bidimensional gaussian values on a grid of size NX by NY +with sigma values SX and SY (real numbers) in units of pixel spacing. +.le +.ls [W(1,1)] ... [W(NX,1)], ..., [W(1,NY)] ... [W(NX,NY)] +The weights are specified as a string of real values. The values are +whitespace separated within each line and the lines are delimited by +comma. For example + +.nf + 1 2 1 + 1 2 1, 2 3 2, 1 2 1 ==> 2 3 2 + 1 2 1 +.fi +.le + +When a logfile is defined the convolution weights are included in the +output. + +5. Object Detection + +The detection of objects in an image is conceptually quite simple once the +background is known. If an input pixel, before any convolution, is +identified in the bad pixel mask the output object mask pixel is also +identified as bad. Otherwise the input data is convolved as described +previously. + +Each convolved pixel is compared against the expected background at that +point and, if it is more that a specified number of convolution adjusted +background sigma above (\fIhsigma\fR) or below (\fIlsigma\fR) the +background, it is identified as a candidate object pixel. Candidate object +pixels, with the same sense of deviation, are grouped into objects on +the basis of being connected along the four or eight neighboring directions +as specified by the \fIneighbor\fR parameter. The candidate object is then +accepted if it satisfies the minimum number of pixels (\fIminpix\fR) in +an object and the \fIhdetect\fR or \fIldetect\fR parameter selects that +type of object. The accepted objects are assigned sequential numbers +beginning with 11. The object numbers are used, as described in the +section on the output data, to set the output object mask values. + +If an output mean sky or sigma map is requested, the output is that +updated by the sky pixels identified during the detection. + +6. Object Growing + +Astronomical objects do not have sharp edges but have light distributions +that merge into the background. This is due not only to the nature of +extended sources but to the atmospheric and instrument point spread function +effects on unresolved sources. In order to include pixels which extend +away from the threshold detection and contain some amount of light +apart from the background, the task provides options to extend or grow +the object boundaries. This is done by making multiple passes where +pixels which have not been identified as object pixels but which neighbor +object pixels are assigned to the object which they neighbor in any of +the eight directions. Each pass can be thought of as adding a ring +of new pixels following the boundary of the object from the previous +pass. + +When a non-object pixel neighbors two or more object pixels it is +assigned to the object with the greater "flux". The flux is the sum +of the pixel value deviations from the background. + +The parameter \fIngrow\fR selects the maximum number of growing iterations. +The parameter \fIagrow\fR selects the maximum fractional increase in +the number of original detected object pixels. The number of pixels +is called the "area" of the object. The growing of an object stops +when either maximum is exceedd at the end of a growing iteration. +.ih +EXAMPLES +1. The following is a test example with default parameters that can be run +by anyone. An artificial galaxy field image is generated with the task +\fBmkexample\fR (the \fBartdata\fR package is assumed to already be loaded) +and a mask is created with \fBobjmasks\fR. The image is displayed with +the object mask overlayed in colors. + +.nf + np> mkexample galfield galfield + Creating example galfield in image galfield ... + np> objmasks omtype=color + List of images or MEF files: galfield + List of output object masks: gfmask + ACE: + Image: galfield - Example artificial galaxy field + Set sky and sigma: + Determine sky and sigma by surface fits: + start line = 1, end line = 512, step = 51.1 + xorder = 2, yorder = 2, xterms = half + hclip = 2., lclip = 3. + Determine sky and sigma by block statistics: + Number of blocks: 5 5 + Number of pixels per block: 100 100 + Number of subblocks: 10 10 + Number of pixels per subblock: 50 50 + Detect objects: + Convolution: + 1. 1. 1. + 1. 1. 1. + 1. 1. 1. + 422 objects detected + Grow objects: ngrow = 2, agrow = 2. + Write object mask: gfmask[pl,type=mask] + np> display galfield 1 + z1=371.5644 z2=455.8792 + np> display galfield 2 overlay=gfmask[pl] ocolors="+203" + z1=371.5644 z2=455.8792 +.fi + +2. In the first example there was no input mask. The next example +creates a new object mask using the first object mask as an input +"bad pixel mask". While this is not the usual usage of the bad pixel +mask it does illustrate an interesting option. Note that the mask +values in the input mask are mapped to an output value of 1 in the +"colors" output. In this example the output is forced to be a pl +file by using the explicit extension. + +.nf + np> objmasks omtype=colors mask=gfmask[pl] + List of images or MEF files (galfield): + List of output object masks (gfmask): gfmask1.pl + ACE: + Image: galfield - Example artificial galaxy field + Bad pixel mask: gfmask.pl + Set sky and sigma: + Determine sky and sigma by surface fits: + start line = 1, end line = 512, step = 51.1 + xorder = 2, yorder = 2, xterms = half + hclip = 2., lclip = 3. + Determine sky and sigma by block statistics: + Number of blocks: 5 5 + Number of pixels per block: 100 100 + Number of subblocks: 10 10 + Number of pixels per subblock: 50 50 + Detect objects: + Convolution: + 1. 1. 1. + 1. 1. 1. + 1. 1. 1. + 44 objects detected + Grow objects: ngrow = 2, agrow = 2. + Write object mask: gfmask1.pl + np> display galfield 2 overlay=gfmask1 ocolors="+203" + z1=371.5644 z2=455.8792 +.fi + +3. The next example illustrates use with a multiextension file. The +example is two realizations of the galfield artificial data. + +.nf + np> mkexamples galfield mef.fits[im1] + Creating example galfield in image mef[im1] ... + np> mkexamples galfield mef[im2,append] oseed=2 + Creating example galfield in image mef[im2,append] ... + np> objmasks + List of images or MEF files (galfield): mef + List of output object masks (gfmask1.pl): mefmask + ACE: + Image: mef[im1] - Example artificial galaxy field + Set sky and sigma: + Determine sky and sigma by surface fits: + start line = 1, end line = 512, step = 51.1 + xorder = 2, yorder = 2, xterms = half + hclip = 2., lclip = 3. + Determine sky and sigma by block statistics: + Number of blocks: 5 5 + Number of pixels per block: 100 100 + Number of subblocks: 10 10 + Number of pixels per subblock: 50 50 + Detect objects: + Convolution: + 1. 1. 1. + 1. 1. 1. + 1. 1. 1. + 422 objects detected + Grow objects: ngrow = 2, agrow = 2. + Write object mask: mefmask[im1,append,type=mask] + ACE: + Image: mef[im2] - Example artificial galaxy field + Set sky and sigma: + Determine sky and sigma by surface fits: + start line = 1, end line = 512, step = 51.1 + xorder = 2, yorder = 2, xterms = half + hclip = 2., lclip = 3. + Determine sky and sigma by block statistics: + Number of blocks: 5 5 + Number of pixels per block: 100 100 + Number of subblocks: 10 10 + Number of pixels per subblock: 50 50 + Detect objects: + Convolution: + 1. 1. 1. + 1. 1. 1. + 1. 1. 1. + 410 objects detected + Grow objects: ngrow = 2, agrow = 2. + Write object mask: mefmask[im2,append,type=mask] + np> display mef[im1] 1 over=mefmask[im1] + z1=371.5644 z2=455.8792 + np> display mef[im2] 2 over=mefmask[im2] + z1=371.5666 z2=455.7844 +.fi + +4. This example shows outputing the sky information. + +.nf + np> objmasks galfield gfmask2 sky=gfsky2 + ACE: + Image: galfield - Example artificial galaxy field + Set sky and sigma: + Determine sky and sigma by surface fits: + start line = 1, end line = 512, step = 51.1 + xorder = 2, yorder = 2, xterms = half + hclip = 2., lclip = 3. + Determine sky and sigma by block statistics: + Number of blocks: 5 5 + Number of pixels per block: 100 100 + Number of subblocks: 10 10 + Number of pixels per subblock: 50 50 + Write sky map: gfsky2 + Detect objects: + Convolution: + 1. 1. 1. + 1. 1. 1. + 1. 1. 1. + 422 objects detected + Update sky map: gfsky2 + Grow objects: ngrow = 2, agrow = 2. + Write object mask: gfmask2[pl,append,type=mask] + np> imstat gfsky2 + # IMAGE NPIX MEAN STDDEV MIN MAX + gfsky2 25 401.1 0.4397 400.3 401.9 +.fi + +5. This examples shows specifying the sky information as constant values. +In this case we already know that the artificial image has a +constant background of 400 and a sigma of 10. + +.nf + np> objmasks galfield gfmask3 sky=400 sigma=10 + ACE: + Image: galfield - Example artificial galaxy field + Set sky and sigma: + Use constant input sky: 400. + Use constant input sigma: 10. + Detect objects: + Convolution: + 1. 1. 1. + 1. 1. 1. + 1. 1. 1. + 432 objects detected + Grow objects: ngrow = 2, agrow = 2. + Write object mask: gfmask3[pl,append,type=mask] +.fi + +.ih +REVISIONS +.le +.ih +SEE ALSO +.endhelp +266c266 +< fit to produce a function for <S>. +--- +> fit to produce a constant function for <S>. +273c273,274 +< weights are used to fit a 2D plane for both <I> and <S>. +--- +> weights are used to fit a 2D plane for both <I> and <S>. The <S> surface +> is a constant in order to avoid potential negative sigma values. + |