PURPOSE: GEOMA is a VICAR applications program which makes geometric changes in pictures. It can be used for correcting geometric distortion, or performing other types of geometric transformations. The geometric changes (transformations) are defined as a mapping from a set of locations (called tiepoints or control points) in the input image to a corresponding set of locations in the output image. Between the tiepoints, the mapping is defined by bilinear interpolation. GEOMA is similar in function to VICAR program LGEOM. EXECUTION: The following TAE command line formats show the most common usages: GEOMA INP=(a,a2) OUT=b SIZE=(1,1,nl,ns) optional parameters GEOMA INP=(a,a2) OUT=b NL=nl NS=ns optional parameters GEOMA (a,a2) b (1,1,nl,ns) optional parameters GEOMA (a,a2) b optional parameters GEOMA INP=a OUT=b SIZE=(1,1,nl,ns) optional parameters GEOMA INP=a OUT=b NL=nl NS=ns optional parameters GEOMA a b (1,1,nl,ns) optional parameters GEOMA a b optional parameters GEOMA INP=a OUT=b SIZE=(1,1,nl,ns) PARMS=c optional parameters GEOMA INP=a OUT=b NL=nl NS=ns PARMS=c optional parameters GEOMA a b (1,1,nl,ns) PARMS=c optional parameters GEOMA a b PARMS=c optional parameters Here 'a' represents the input image file name, 'a2' represents an IBIS input file containing tiepoint information, 'b' represents the output image file name, and 'c' represents the parameter file name. Note: GEOMA no longer supports a non-IBIS second input file. Old non-IBIS 'GEOMA' files (which were often used to specify the geometric transformation needed to correct for camera distortion) can be converted on the VAX or Alpha using the program OLDGEOMA2IBIS. USER INFORMATION: In both GEOMA and LGEOM the picture to be geometrically transformed is subdivided into a number of four-sided areas (quadrilaterals). The user must specify to the program the exact mapping from input picture to transformed picture for each of the four corner points of each quadrilateral. Within the quadrilaterals, the programs use a bilinear interpolation to map the interior points from input picture to transformed picture. The LGEOM program restricts the specification of the mapping of the quadrilateral corners (control points) such that, in the transformed (output) picture, they form a rectangular grid, (i.e. the quadrilaterals are rectangles in the output image.) The GEOMA program does not impose this restriction, and allows more arbitrary specification of the control point mapping. In fact, two of the four control points forming the corners of any given quadrilateral may coincide, causing the quadrilateral to degenerate into a triangle (a quadrilateral having one side of zero length). The potential user is WARNED that GEOMA can produce some undesirable side effects due to the basic difference from the LGEOM program. This is explained further in the OPERATION section. As with LGEOM, GEOMA will accept its control parameters as the second input file or as a parameter file instead of or in addition to the conventional VICAR parameters. The size of the output image is determined by the number of lines and number of samples in the SIZE field if the SIZE field is entered by the user. If the SIZE field is not entered, the output file is the same size as the input file. The data type of the input image may either be byte or halfword data. The data type is obtained from the VICAR label of the input image or from the FORMAT parameter. The output image has the same data format (byte or halfword) as the input image. The tiepoint grid must be specified as a set of parameters in the TAE command line or in the parameter file or in a second input file. The tiepoint grid is specified using the parameters NAH, NAV, and TIEPOINT. NAH=nah NAV=nav TIEPOINT=(nl1,ns1,ol1,os1 nl2,ns2,ol2,os2... ...nlk,nsk,olk,osk) where nah is the number of grid cells horizontally (across the top in the sample direction), nav is the number of grid cells vertically in the output image space, the point pairs (nli,nsi,oli,osi) are line-sample coordinates in the output (new) and input (old) spaces respectively. The number of pairs k must equal (nah+1)*(nav+1). The RESTART parameter has been removed since it was deemed not maintainable. (See FR 16750 of February 1986.) Applications of GEOMA: GEOMA and LGEOM are both VICAR programs for performing general geometric transformations. Both of them require a grid of tiepoints for defining the transformation. In LGEOM, the grid must be rectangular (i.e. the tiepoints are arranged in straight rows and columns.) In GEOMA the grid does not need to be rectangular. This allows GEOMA to be more flexible or general purpose. GEOMA and LGEOM are almost always used following a VICAR program that generates the control information (tiepoints) in either a parameter file or an input file. This keeps the user from having to enter the tiepoints manually. This also makes most of the details regarding the tiepoint grid transparent to the user. GEOMA is often used following VICAR program RESLOC and sometimes following VICAR program PICREG. GEOMA's greater flexibility is not often needed. GEOMA is generally slower than LGEOM, and it is less commonly used than LGEOM. The most common use of GEOMA is in conjunction with RESLOC. In this application GEOMA is used to remove geometric distortion in Voyager pictures based on the reseau locations. RESLOC writes the control information in a file which is used as the second input file for GEOMA. Sometimes VICAR program FARENC is used between RESLOC and GEOMA to add a translation to the geometric transformation when dealing with a limb of a planet or satelite. GEOMA is used in this application because the reseau do not lie in a rectangular grid. In general, program LGEOM is recommended over GEOMA because of speed unless the greater flexibility of GEOMA is needed. For instance, it is more common to use LGEOM than GEOMA following PICREG. Where GEOMA is used, it is common to use a grid of degenerate quadrilaterals (triangles) as shown below in the diagrams in the OPERATION section. Users can also use VICAR program MGEOM as an alternative to GEOMA when the output tiepoints lie in a rectangular grid. The advantages of MGEOM are listed in the help for MGEOM. EXAMPLES 1. RESLOC (IMG,TF) (R,G) GEOMA (IMG,G) IMGB (1,1,800,800) In this example GEOMA is used to remove geometric distortion in a Voyager image named IMG based on the reseau locations. RESLOC writes the control information in file G which is used as the second input file for GEOMA. GEOMA produces an 800 by 800 output image named IMGB. 2. PICREG (J,J2) (T,G) ... numerous steps to interactively select tiepoints NHOR=10 NVER=10 TPFORM=3 'GEOMA 'EXIT GEOMA INP=J OUT=JB PARMS=G (1,1,800,800) In this example GEOMA is used to register image J to image J2 using tiepoints selected interactively in program PICREG. The tiepoints selected by the user are converted into a rectangular grid with 10+1 rows and columns using a first-order surface fit (TPFORM=3). Each rectangle is divided into two triangles to make the final grid of degenerate quadrilaterals (triangles) for GEOMA. The control information is passed to GEOMA via parameter file G. GEOMA produces the registered image JB. 3. GEOMA INP=T3 OUT=T4 SIZE=(1,1,20,20) NAH=1 NAV=1 TIEPOINT=(1,1,1,1 1,20,1,15, 20,1,10,1 20,1,10,1) In this example GEOMA is used to perform an affine (first order) transformation defined by three tiepoints. (1,1) is mapped to (1,1). (1,15) is mapped to (1,20). (10,1) is mapped to (20,1). The rest of the mapping is defined by linear interpolation. 4. The last example is the test procedure for GEOMA. This is a complete example that could be run by the user and that demonstrates uses of the possible parameters. !THIS IS A TEST FILE FOR GEOMA !THIS TAKES A 10 X 10 MATRIX AND ENLARGES IT TO A 20 X 20 MATRIX gen tgen1 10 10 linc=5 sinc=10 geoma inp=tgen1 out=tgen2 size=(1,1,20,20) nah=1 nav=1 tiepoint=(1,1,1,1, 1,20,1,10,20,1,10,1,20,20,10,10) list tgen2 !THIS TAKES A 10 X 10 HALFWORD MATRIX AND ENLARGES IT TO A 20 X 20 HALFWORD !MATRIX gen tgen3 10 20 'half sinc=-1 linc=2 geoma inp=tgen3 out=tgen4 size=(1,1,20,20) 'half nah=1 nav=1 tiepoint=(1,1,1,1,1,20,1,20,20,1,10,1,20,20,10,20) list tgen4 !REPEAT WITH PARAMETER NOINTERP gen tgen1 10 10 ival=10 sinc=2 linc=5 geoma inp=tgen1 out=tgen6 size=(1,1,100,100) nah=1 nav=1 tiepoint=(1,1,1,1, 1,100,1,10,100,1,10,1,100,100,10,10) 'nointerp list tgen6 !THIS ROTATES A 10 X 10 MATRIX 90 DEGREES !ALSO TEST PARAMETER NOTIFY gen tgen1 10 10 geoma inp=tgen1 out=tgen7 size=(1,1,10,10) nah=1 nav=1 tiepoint=(1,1,10,10, 1,10,10,1,10,1,1,10,10,10,1,1) notify=20 list tgen7 !GEOMETRICALLY CORRECTS AN IO IMAGE USING SPECIFIED PARAMETER FILE gen tgen1 10 10 ival=10 sinc=2 linc=5 DCL DEFINE VGR MIPLDISK:[MIPL.VGR] geoma inp=(VGR:F1636832.FIC,VGR:F1636832.IGPR) out=tgen5 label-list tgen5 OPERATION: GEOMA computes the DN value of each pixel in the output picture as follows. First, decide within which quadrilateral area in the output picture the pixel lies. Let the coordinates of the four control points defining that quadrilateral in the output picture be called (x(j),y(j),j=1,2,3,4), and in the input picture let them be called (x'(j), y'(j),j=1,2,3,4). (x may be considered the line coordinate and y the sample, although they could as well be reversed for this discussion.) The values of all 16 of these numbers are part of the control parameters for GEOMA. Define a transformation by x' = ax + by + cxy + d (1) y' = ex + fy + gxy + h where the values of the coefficients a,b,..,h are determined by requiring that x'(j) = ax(j) + by(j) + cx(j)y(j) + d (2) y'(j) = ex(j) + fy(j) + gx(j)y(j) + h for j = 1,2,3,4 This is the condition that the defined transformation exactly maps the control point coordinates as specified by the parameters. The transformation is used to transform the line-sample coordinates of the output picture pixel being processed into its corresponding coordinates in the input picture. In general, although the output coordinates are always integer, the input coordinates are not, so the input coordinates do not correspond to any one pixel but instead fall in between the input pixels. If intensity interpolation is NOT to be performed, as specified by the NOINTERP parameter, then the DN value assigned to the output pixel is that of the pixel in the input picture nearest the transformed coordinates. The use of this mode considerably reduces the CPU time required by GEOMA. The quality of the result is in general degraded but still may be quite acceptable, depending on the nature of the transformation. Assuming interpolation is to be performed to compute the DN value of the output picture, the four pixels nearest to the transformed coordinates in the input picture are selected. The interpolation formula is DN = px' + qy' + rx'y' + s (3) where x' and y' are the transformed coordinates of the output picture pixel, DN is the data number to be assigned to the output pixel, and p,q,r and s are constants chosen such that DN(k) = px'(k) + qy'(k) + rx'(k)y'(k) + s (4) for k = 1,2,3,4 where x'(k), y'(k) and DN(k), k = 1,2,3,4 are the coordinates and DN values of the four selected input pixels nearest to the transformed coordinates. In order for GEOMA to work properly, the tiepoint coordinates in the output picture must satisfy certain constraints. They must be organizable into a rectangular matrix of points. That is, each tiepoint must belong to one row and one column of points, each row must have the same number of points as every other row, and each column must have the same number as every other column. It should be possible to number the rows sequentially from top to bottom, and the columns form left to right. Quadrilaterals used by the program for performing the geometric transformation are formed by connecting each tiepoint to the ad- jacent tiepoints in the same column and in the same row. The quadrilaterals should all have the property that they are not concave (no interior angle of any quadrilateral shoud be greater than 180 degrees). In addition, each tie- point should be below any tiepoint in its column and in any preceding row, and above any tiepoint in its column in any subsequent row. Similarly, each tie- point should be to the right of any tiepoint in its row and in a preceding column, and to the left of any tiepoint in its row and in a subsequent column. The number of areas vertically is the number of rows of tiepoints less 1, and the number of areas horizontally is the number of columns of tiepoints less 1. The following caution is issued to potential users of GEOMA. Although the transformation in equation (1) is very well behaved within the quadrilateral in which it applies, there is no assurance that it is continuous across the boundary between adjacent quadrilaterals unless the boundary is precisely vertical or horizontal. The degree of discontinuity depends on the details of the tiepoint specifications. In some cases the discontinuity may be so small that the transformation picture has no visible defect, but in others it may be quite visible. The user has two choices to minimize the discontinuity. He may carefully choose the transformation such that discontinuities are minimized, or he may specify the tiepoints in such a way that the quadrilaterals degenerate into triangles. A quadrilateral degenerates into a triangle if two adjacent tiepoints defining it have the same coordinates. By suitable choice of tiepoints, the user can specify the coordinates in such a way that all quadrilaterals degenerate into triangles. This, in fact, is commonly done. (See programs RESLOC and PICREG.) The figures below show sample arrays of tiepoints satisfying this condition in output-image views. (The rows and columns are straight below because it was easier to draw, but they are not required to be straight.) 1,2 3,4 5 +---------+--------- | \ | \ | | \ | \ | | \ | \ | | \7|8 \ | 6 +---------+---------+ 9,10 | / | / | | / | / | | / | / | | / | / | +---------+--------- 11,12 13,14 15 Sample Triangular Tiepoint Array: 4 areas horizontally, 2 areas vertically. 1,2 3,4 5,6 7,8 +---------+---------+--------- | \ | \ | \ | \ | \ | \ | \ | \ | \ | \ | \ | \ | \ |10 \ |12 \ |14 \ 9 +---------+---------+---------+---------+ 16 | / |11 / |13 / |15 / | / | / | / | / | / | / | / | / | / | / | / | / +---------+---------+--------- 17,18 19,20 21,22 23,24 Sample Triangular Tiepoint Array: 7 areas horizontally, 2 areas vertically. When an area is a triangle, the cross terms in equation (1) are eliminated and the transformation becomes x' = ax + by + d (5) y' = ex + fy + h where the values of the coefficients are determined by requiring that x'(j) = ax(j) + by(j) + d (6) y'(j) = ex(j) + fy(j) + h for j = 1,2,3,4 The values of j designate the three distinct tiepoints defining the triangle. Because equation (5) is linear, the transformation is guaranteed to be continuous at the triangle boundaries and discontinuities in the output picture cannot occur. There is a question about how to transform picture samples which fall outside all the quadrilaterals defined by parameters. One answer is to avoid this situation by specifying the control point locations so that every sample in the output picture falls within some quadrilateral. This requires that some control points fall exactly on or outside the border of the output picture. The program is not affected in any way by using control points outside the picture. If some picture samples do fall outside the defined quadrilaterals, the program will process each such sample by assigning it to a nearby quadrila- teral. This is done by "extending" the boundaries of all "edge" quadrilaterals to the picture borders. The details of the precise extrapolation algorithm will not be given. The user should understand that the algorithm can lead to discontinuities in the processed picture in the regions outside the defined quadrilaterals, even when the control areas are triangles. RESTRICTIONS: 1. The input and output images must be byte or halfword data. 2. There are essentially no restrictions in GEOMA on picture size or number of tiepoints. (The input image may not have >32767 lines.) 3. The maximum number of tiepoints that can be passed via a parameter file is currently 2047. 4. The maximum number of tiepoints that can be specified as part of the TAE command line is 125. 5. The maximum number of tiepoints that can be passed via an IBIS tiepoint file is set to 10,000. (This could be increased, but you would probably want to increase the size of BIGBUF proportionally.) TIMING: Under certain circumstances, such as 90 degree rotations, GEOMA will be much slower than LGEOM. (An "LGEOMA" could be written, using the internal algorithm of LGEOM, which would be much faster in those circumstances. However, such a program does not currently exist.) For small rotations GEOMA will be faster than LGEOM. As with most VICAR programs, the elapsed time for execution will be shortened if the input and output files are on separate disk drives. Under VMS the execution time and the number of page faults can be affected by the working set size that was set by the system manager. If a large number of page faults are observed, the working set size should be increased. Alternatively, it may be possible to reduce the NBBUF value in the program. PRECISION: This portable version of GEOMA will generally produce identical output on all of the MIPS-supported platforms. This output will frequently not be identical to the output of the unported version of GEOMA, but the differences will be trivial and should be improvements. This section is required by the MSTP SRD for all programs that do not meet the precision requirements given in the SRD. The two requirements not met by GEOMA are: 1) Integer data output by a ported program must be the same as that produced by the unported program. 2) Integer data output by a ported program must be the same on all of the MIPS-supported platforms. All of the programs that do geometric transformations are going to have problems meeting one or both of these requirements. This includes MAP3, LGEOM, MGEOM, GEOMA, and several others. GEOMA was modified to increase the use of Double Precision floating point arithmetic so that the first requirement is met "virtually always" and the second requirement is met "well enough". The definition of the above qouted terms is given towards the end of this section. In the port, GEOMA was modified to retain the full precision returned from subroutine DGELG. This was done to eliminate or minimize the differences produced on different platforms. The results are somewhat more accurate than those of the unported GEOMA, but since the output DN values are rounded to the nearest integer, the output images are essentially the same as for the unported GEOMA. In computing the output pixel value (DN), GEOMA usually interpolates between four input pixel values, but if the input line coordinate is close enough to an integer value, it interpolates between two input pixel values. For byte images, the two pixel interpolation and the four pixel interpolation yield the same result virtually always once the result is rounded to the nearest integer. For halfword images, where the pixel values are often a hundred times greater, the two pixel interpolation and the four pixel interpolation occassionally yield slightly different results. The image f1636832.fic used in the GEOMA test pdf shows what type of differences are likely to occur (with the ported GEOMA) between machines, and what type of differences are likely to occur between the ported GEOMA and the unported GEOMA. The output file from the ported GEOMA is identical on all of the MIPS supported machines. The output file from the ported GEOMA is very slightly different from the output file from the ported GEOMA. 98% of the pixels are identical. Almost all of the rest differ in value by 1. The minimum difference is -10. The maximum difference is 17. Program DIFPIC shows the AVE VAL OF DIFFS=-0.125000E-03, or approximately .0001. The cause of the differences is that the single precision computations (unported GEOMA) round differently than the double precision computations (ported GEOMA), resulting in differences in when two pixel interpolation is used and when four pixel interpolation is used. The above halfword image has a number of locations where the difference between adjacent pixels is greater than 10,000. This combines to produce a slight but insignificant difference. This difference is estimated at less than 1% of the difference between using 'NOINTERP and using regular interpolation, which is the definition I will use for "well enough". I believe that such differences will not be humanly observable in an image display. To complete this discussion, I want to define the other quoted term mentioned above. "Virtually always" is always minus the probability of a difference in output between different MIPS supported machines. Differences are possible when 1) one machine does a two pixel interpolation while another machine does a four pixel interpolation. 2) one machine produces a double-precision output pixel value very slightly less than halfway between two integers and another machine produces an output pixel value equal to or slightly greater than halfway between two integers. 3) one machine determines that an output pixel lies outside the geometrically transformed boundaries of the input image while another machine decides the output pixel lies within the boundaries. Possibilty 2 is approximately the probability of a double precision result in the range between 0.0 and 1.0 being close enough to 0.5 to round differently on two machines. Double-precision rounding differences start at about 1.E-14 on the MIPS-supported machines and might grow to a maximum of 1.E-8 going through the computations in GEOMA. So the above probability is less than the probability of the result being in the range 0.49999999 and 0.50000000. Assuming a fairly uniform distribution, this is probability is 1.E-8. Possibility 1 can be shown to have a similar probability. Possibility 3 has an even lower probability because the border pixels are such a small fraction of the total pixels. However, possibily 3 could cause one machine to generate a pixel value of 0 while another machine could generate a relatively large pixel value from the edge of the input image. The other possibilities only generate a difference of 1 (possibility 2) or something small relative to the difference between adjacent pixels (possibility 1). In an image display, none of these differences should be humanly observable. Based on this discussion, "virtually always" is defined as no more than 1 different pixel in a set of 1.E8 pixels. This amounts to no more than one difference in a hundred 1000 by 1000 pixel images. WRITTEN BY: J. SEIDMAN 22 DECEMBER 1970 CONVERTED TO VAX BY: HELEN B. MORTENSON 15 SEPTEMBER 1983 CONVERTED TO VICAR2 BY: HELEN B. MORTENSON 8 APRIL 1985 Made portable for UNIX ... Jim Turner (CRI) 5 September 1994 COGNIZANT PROGRAMMER: Steve Pohorsky REVISION HISTORY 7-95 SP Added support for an IBIS tiepoint file as a second input. 6-96 OAM Deleted default value of parameter OUT in geoma.pdf to make program compatible with SAGE (DFR).