Morphing in 2-D and 3-D
Where image processing meets computer graphics
Valerie Hall
Valerie is a PhD student in computer science at Curtin University of Technology in Western Australia. Her main area of research is computer graphics, and her PhD topic is speech-driven facial animation. She can be reached at val@ marsh.cs.curtin.edu.au.
Special effects in the movies have always captured our imagination, but morphing seems to have attracted special attention. The average moviegoer is familiar with shape-changing sequences in films such as Willow, Terminator II, The Abyss, and Lawnmower Man; MTV fans regularly see this technique in music videos such as Michael Jackson's "Black or White;" and couch potatoes view it daily in everything from commercials pitching Chrysler vans and Schick razor blades to prime-time TV shows like "Deep Space 9." Even PC graphics packages--Autodesk's 3D Studio comes to mind--now offer morphing as part of their feature set.
It's often difficult to get precise information about the methods various commercial animators use because of the proprietary nature of the software in high-end animation shops. Still, morphing source code is publicly available for workstations; Mark Hall's morphine program for X Windows workstations and Tim Heidmann's demo program for the Silicon Graphics environment are good examples of this. For PC users, shareware morphing apps such as Richard Goedeken's Rmorf (see the accompanying textbox entitled "Rmorf:A Shareware Morphing Program for MS-DOS") are available on CompuServe and similar sources.
In this article, I'll provide a technical overview of this fascinating technique. But first, let's start with definitions and a bit of history.
What is Morphing?
The term morphing comes from the Greek word morphé, which means form or shape. The study of shapes is known as "morphology," and morphing has come to mean shape-changing via digital techniques. Here, I broadly use the definition to cover both two- and three-dimensional techniques.
The term originated in the late '80s at George Lucas's Industrial Light and Magic (ILM). Douglas Smythe and others at ILM developed Morf, a program for interpolating sequences between two images. Morf was originally written to handle the transformation scenes in the movie Willow and has since been used on other projects. (See the accompanying textbox, "How Do They Do It?".) Although the term "morphing" may be recently coined, many of the basic algorithms have been around for a decade in the fields of computer graphics and digital image processing.
While it's difficult to single out individuals, certain milestones are evident as you look back in time. In 1990, George Wolberg wrote Digital Image Warping (IEEE Computer Society Press, 1990), the definitive compendium of 2-D image-warping techniques. (For more details on Wolberg's approach, see the accompanying textbox entitled, "The Canonical Implementation in C," which implements in C the algorithm used by ILM.) I've already mentioned Doug Smythe's work at ILM in the late '80s. In the early '80s, Tom Brigham and Paul Heckbert developed some interesting image-transformation sequences at NYIT's Computer Graphics Lab. In 1980, Ed Catmull and Alvy Ray Smith published a seminal paper on computing efficient geometric transformations. In the late '70s Julian Gomez developed a "tweening" program at Ohio State (and later at NYIT). Prior to that, in 1976, Jim Blinn and Martin Newell published a key paper on texture mapping.
The earliest work in geometric transformation of digital images came from the field of remote sensing, which gained wide attention in the mid-'60s, when NASA undertook projects to observe Earth from space. Photographic instruments on Landsat and Skylab produced multiple overlapping views of the same terrestrial region. To align these images with each other, it's necessary to compensate for distortions such as lens aberrations and differences in viewing angles. Similar image-correction techniques have been used in medical imaging and digital radiology.
But even before digital processors and computer graphics, Hollywood moviemakers accomplished rudimentary morphing by cross-dissolving images. One memorable metamorphosis sequence appeared in the 1941 horror film The Wolf Man. Cross-dissolving is still used as a special effect in conjunction with the two kinds of digital morphing.
The Two Faces of Morphing
The morphing process varies, depending on whether the morph is 2-D or 3-D. 2-D morphing gives the visual effect of a 3-D change of shape by warping a 2-D image from an initial shape to a final shape. Using digital image-warping algorithms, an initial image is stretched and deformed to conform to the shape of the target. At the same time, the textures for each image are gradually blended from the initial texture to the final one. For greater control, source and target images are broken up into small regions that map onto each other.
In 3-D morphing, a 3-D geometric model of the object is transformed from one shape into another. At each stage in the metamorphosis, the 3-D model is rendered and texture-mapped to produce a 2-D screen representation. The rendering and texture-mapping techniques are standard in 3-D computer graphics. For 3-D morphing, the animator must set up a correspondence between different components of the initial and final shapes. Difficulties arise when trying to morph 3-D objects that are structurally different. For example, morphing a torus (doughnut shape) to a rectangle poses problems with the doughnut's hole and the rectangle's edges. But before getting into this, let's first consider 2-D morphing in more detail.
2-D Morphing and Texture Mapping
The techniques for 2-D morphing come from digital image warping and texture mapping. Texture mapping is widely used in computer graphics. I distinguish between computer graphics (the generation of synthetic images) and image processing (the manipulation of captured-image data). Although texture mapping is used in the visualization of synthetic images, its counterpart in the domain of image processing is digital image warping. Both texture mapping and digital image warping rely on geometric transformations to redefine the spatial relationship between sets of points in an image.
For texture mapping, the basic technique is a two-step process of first mapping a 2-D texture plane onto a 3-D surface and then projecting the surface onto the 2-D screen display. An oft-used analogy to explain texture mapping is that of clothing: You start with an object (body) and then wrap an image (article of clothing) around it. Texture mapping serves to create the appearance of complexity by applying elaborate image detail to relatively simple surfaces. In computer graphics, textures can be used to perturb surface normals and thus allow simulation of bumps and wrinkles without the effort of modeling intricate 3-D geometries.
Digital image warping leaves out the intermediate step of mapping to 3-D object space, and instead maps directly from one 2-D space (the input image) to another (the output image). Both digital image warping and texture mapping need an efficient way to accomplish a spatial transformation--a mapping that establishes a spatial correspondence between the two 2-D coordinate spaces, of input image and output image. Although digital image warping predates texture mapping, it is in the field of computer graphics that efficient two-pass algorithms have been developed to compute arbitrary geometric transformations, beginning with Catmull and Smith's 1980 paper.
Geometric Mappings
There are two ways to calculate the geometric correspondence between points in each of these two spaces: from input to output image (known as "forward mapping"), or from output back to input (known as "inverse mapping" or "screen-order traversal"). Screen order is the more common. Your program traverses the output image scanline by scanline, pixel by pixel; for each pixel in the output grid, the corresponding value in the input image is derived. This method is most useful when your program is required to write to the screen sequentially, when the mapping is readily invertible, and when the texture allows random access.
Various kinds of geometric transformations are possible between two coordinate systems, including affine, projective, bilinear, and polynomial transformations. The equations for calculating basic 2-D transforms use 3x3 matrix multiplications based on homogenous coordinates. A full discussion of this subject is beyond the scope of this article, but you can find a good explanation in George Wolberg's book, as well as in standard computer-graphics texts such as Foley and van Dam's Computer Graphics (Addison-Wesley, 1990).
Regardless of the particular geometric transformation, a problem arises when mapping from one grid to another: "Holes" and overlaps can occur. This problem occurs whether you are mapping from output image to input, or vice versa. It results from the fact that you are working with discrete, rather than continuous, images. Consequently, you'll end up with pixels in one space that don't have an exact correspondent in the alternate space. Deriving the appropriate value in the alternate space requires interpolation and sampling. Interpolation is a process for arriving at a continuous surface that passes through the discrete points in your image data. This continuous surface can then be sampled at arbitrary positions, not just the ones on the coarse coordinate grid.
There are a host of interpolation functions possible: cubic convolution, bilinear, cubic spline, and sinc-function convolution. The easiest approach is "nearest neighbor"--taking the pixel closest to the one wanted. Regardless of approach, the basic issue with image interpolation and resampling is still how to arrive at the most representative value for a given pixel. Arriving at the "wrong" value results in aliasing in which unexpected pixel values can make the image look chunky or motley. You can bring various techniques into play to reduce these aliasing artifacts. For more on these techniques, see Wolberg's book or a signal-processing text.
Many of these techniques are computationally intensive. Catmull and Smith's contribution in their 1980 paper is a method for decomposing a spatial transform into a sequence of computationally cheaper mapping operations. Specifically, they show how a 2-D resampling problem can be replaced with two orthogonal 1-D resampling stages. That is, your program first transforms the x-coordinates, then the y-coordinates. Although the basic idea is simple, there are many subtleties beyond the scope of this discussion.
Subdividing an Image
To specify a morph, the animator defines a correspondence between the initial and final images. This is commonly done using points, triangles, or meshes. (Figure 1 shows the meshing process the Rmorf program uses.) Once the relationship has been defined, the textures in the images are blended from the initial image to the final image to produce the morph sequence.
The method of subdividing the images varies between different morphing programs. Triangulations can be specified by hand or generated by the morphing program. For example, in the morphing programs written by Mark Hall, Michael Kaufman, and Jay Windley, the user begins by specifying the triangles to be morphed in each image. Hall's program warps the triangular mesh of the input image to match the mesh of the output image. As the animation goes from the initial to the final image, the colors change from that of the first image to that of the second. Likewise, the user of Kaufman's program begins by defining common points in two images. For a human face, these points include the eyes, eyebrows, nose, and mouth. The program builds triangles from these points to create corresponding areas on the images. Windley's program does a linear blend of the contents of the triangles using barycentric geometry. One of the most automatic triangulation algorithms is the Delaunay triangulation, also known as the Voronoi tessellation. (See "Spatial Data and the Voronoi Tessellation," by Hrvoje Lukatela and John Russell, DDJ, December 1992.)
Morphing isn't always done with triangular patches. Douglas Smythe's program fits a Catmull-Rom spline through the x-coordinates of the control points (in the first pass) and the y-coordinates (second pass) to realize a piecewise continuous mapping function. Likewise, the Montage program for the Amiga (developed by Thomas Krehbiel and Kermit Woodall at Nova Design) uses cubic splines in a deformable mesh to define the input and output images. Cubic splines are better for defining curves in the image. Triangular patches have to be small and well chosen to give a similarly smooth result, but on the positive side, they are simpler to work with.
2-D morphing is a very manual process. It relies heavily on the experience and know-how of the animator. The choice of subject matter can make all the difference. Care should also be taken to use images of the subject captured from similar angles. Similarly, the position of parts of the subject should be closely matched.
3-D Morphing
In 3-D morphing, the transformation is from one 3-D object model to another. The objects can be similar in type (people's faces, for example) or they can be as geometrically different as a cube and a sphere. When the topology is similar, the morphing will be a point-to-point mapping of some type of mesh. When the topology differs, finding corresponding points is more difficult.
The methods used for 3-D morphing depend on the type of objects the animator is working with. Convex shapes are the easiest, because they have no inward curves on their surface. The simplest 3-D morphs are when the initial and final shapes have similar topologies. For example, regions of one facial model can be mapped onto the same region on another. In their paper on human prototyping, Magnenat-Thalmann et al. use morphing to create new synthetic actors. They reorganize the input faces by creating a new structure of facets. These facets need to line up with reference points on each face. Corresponding regions are defined on each face, and then the insides of the regions are divided up using an automated algorithm. Once the correspondence between the faces is defined, a morph can be done using linear interpolation.
In a 1989 article on shape distortion, Wes Bethel and Sam Uselton outlined their 3-D morphing system, which can perform a warp between two 3-D objects. Although the initial version required objects to be polyhedrons, the principles allow the system to work with B-spline surfaces as well. The first step is to construct a B* tree for each object. The nodes of the tree represent faces, and its branches represent adjacency relationships between faces. The user then selects one face and one vertex from that face, for each object, to set up the correspondence between the objects.
Once the trees and the correspondences have been defined for the objects, a union of their topologies has to be created. A new B* tree is created to hold the information about the union. The goal is to preserve the adjacency relationships among the faces. If one of the objects has more faces than the other, then an extra face is added to the union tree. Missing faces and vertices are added as needed. If one of the objects has a hole in it, then the union will also have a hole. The result is a master object that has all the topological characteristics of each of the input objects.
The most difficult part of the processing is assigning new coordinate values to the master object. These coordinates need to be defined for the initial and final keyframes. Convex and star-shaped objects are fairly straightforward, but objects with holes are more difficult. If there's a hole in the initial frame but not in the final, then the hole shouldn't be visible in the initial frame. The hole is really there, but it has been geometrically collapsed to hide it. The resultant system can semi-automatically match the components of two objects. Using a simple interpolation scheme, the user is able to do keyframe animation of 3-D scenes of shape-changing objects.
In their paper on topological merging, Kent, Parent, and Carlson present a technique for computing smooth transformations between two polyhedral models. By taking topological and geometric information into account, their system gives results that maintain the connectivity of the polyhedrons at intermediate steps of the transformation, a desirable quality in this type of system. Their system displays far less distortion than that obtained by previous shape-transformation methods.
The algorithm, when given two 3-D models, generates two new models that have the same shapes as the original ones. The new models share the same topology, which is a merger of the original objects' topologies. The merger allows the transformations from one model to the other to be easily computed. Kent's system is currently restricted to "star-shaped" (topologically nonconcave) models without holes. This restriction is only temporary, as the concepts involved are applicable for arbitrary polyhedral models.
There are several issues to address when developing a shape-transforming system. Kent evaluates previous transformation systems thusly:
- Is face connectivity maintained for all intermediate shapes?
- Are the intermediate shapes overly distorted?
- What restrictions exist on the original models?
Kent's system establishes a correspondence between the objects in three stages. First, each model is projected onto the surface of a unit sphere. Then the union of the two models' projections is taken using a modified version of the Weiler polygon-clipping algorithm. The merged topology is then projected back onto the original models. The results of Kent's approach are highly encouraging when working with convex and star-shaped objects. Future efforts will look at developing projections for general concave polyhedrons and polyhedrons with holes. They are also looking at giving the animator more control over the transformations. By using both topological and geometric information in their system, Kent maintains the integrity of the objects (in the Eulerian sense) while the output is intuitive in appearance. There's much less distortion of the shapes in this method.
Dave Bock, a visualization programmer at the North Carolina Supercomputing Center, takes a different approach to morphing by transforming objects using their mathematical descriptions. Bock generates volume data files that contain varying percentages of the two objects' data sets. The initial and final frames have 100 percent of one object and 0 percent of the other. He interpolates both volume data sets to define their respective objects at a common data threshold to use with an isosurface generator. Then a series of data files is created to hold the combined percentages of the initial data files at different stages of the morph. The number of data files depends on the number of in-between frames the user specifies. For each frame, volume data files are created containing increasing percentages of the final image and decreasing percentages of the initial object. Then the isosurface generator is used to create the geometric object for each frame, based on the volume data in the files. Bock is enhancing his system to accept geometrically defined objects as well as mathematically defined ones.
In all these 3-D morphing systems, the emphasis has been on getting the 3-D models defined--rendering and texture mapping is done later. This is one of the main differences between 2-D and 3-D morphing. Another difference is that 3-D morphing systems tend to be more automated than their 2-D counterparts.
How to Improve Results
Morphing software won't necessarily produce compelling results on its own. The animator must work hard to make the morph smooth and believable. One way of improving results is to define a larger number of points. The use of a smaller mesh gives a higher resolution to the morph so it can stand close scrutiny. This is particularly important if the sequence is to be shown on a large screen.
Key areas on the images need to be separated to stop them from dissolving during morphing. Animators need to pay particular attention to features such as the eyes and mouth on a face. Therefore, a good strategy is to choose images that are fairly closely matched on key features. If the features don't match up, you can warp the images before use so that reference points line up. Even with something as familiar as a human face, observers will rarely notice a well-executed warp. Humans remember the shapes of certain features rather than an exact image of the whole face.
A good idea for a more-deceptive morphing sequence is to stagger the morphs as the image changes; that is, to morph different parts of the image at different times. A good example of this is Michael Jackson's "Black or White" film clip. In part of the dance sequence, as the faces changed to the different dancers, each frame would show some of the features from the initial image while other parts were from the final image. This mix distracts the viewer's eye so that the viewer can't predict where the next feature change will appear. In some parts of "Black or White," there were up to seven planes of morphing going on independently.
Finally, as in all animation, there's no law against doing touch-ups after the computer animation. Cleaning up rough edges to give an improved finish is very common. The need for touching up can be minimized by careful planning before animation begins. Unquestionably, patience on the part of the animator is a necessity.
References
Bethel, E.W. and S.P. Uselton. "Shape Distortion in Computer-Assisted Keyframe Animation." State-of-the-Art in Computer Animation (Proc. of Computer Animation '89), Magnenat-Thalmann & Thalmann, eds. New York, NY: Springer-Verlag, 1989.
Blinn, J.F. and M.E. Newell. "Texture and reflection in computer generated images." Communications of the ACM (October, 1976).
Heckbert, P.S. "Digital Image Warping: A Review." IEEE Computer Graphics and Applications (January, 1991).
Heckbert, P.S. "Survey of texture mapping." IEEE Computer Graphics and Applications (November, 1986).
Kent, J., R. Parent, and W. Carlson. "Establishing Correspondences by Topological Merging: A New Approach to 3D Shape Transformation." Proceedings of Graphics Interface '91.
Magnenat-Thalmann N., H.T. Minh, M. de Angelis, and D. Thalmann. "Human Prototyping." New Trends in Computer Graphics: Proceedings of CG International '88.
Sorenson, P. "Morphing Magic." Computer Graphics World (January, 1992).
Wolberg, G. Digital Image Warping. Los Alamitos, CA: IEEE Computer Society Press, 1990.
How Do They Do It?
There are numerous examples of morphing in recent TV commercials and movies. Here are a few. Much of this information comes from Peter Sorensen's article "Morphing Magic" (Computer Graphics World, January 1992). The rest is from various Internet discussions and e-mail conversations with morphing practitioners.
Willow (Industrial Light and Magic, 1987). In a number of scenes, the character Willow tries to transform a sorceress back to her true form, using a wand with which he is not very familiar. In one scene, the sorceress starts as a goat, later an ostrich, a peacock, a turtle, a tiger, and finally her true human form. For this scene, ILM used deformable puppets of each animal so that they could stretch them into the correct shape for each change. The actual changeover between the puppets was done using 2-D morphing.
Indiana Jones and the Last Crusade (Industrial Light and Magic, 1988). The villain believes he has found the grail and drinks from the cup. He has chosen the wrong cup and proceeds to age rapidly until he dies, shriveling up like a mummy. ILM used 2-D morphing to blend a series of three increasingly grotesque masks depicting the villain's face as he died.
The Abyss (Industrial Light and Magic, 1989). In the movie, the pseudopod--made out of water and shaped like a worm--comes into the human habitat and explores. ILM animated the body of the pseudopod using traditional computer-animation techniques. The face was morphed in 3-D using data from complete 3-D digitizations of the required facial expressions. Doug Smythe's Morf program was used to do the interpolations between facial expressions.
Terminator II (Industrial Light and Magic, 1991). Nearly every scene with the T1000 character included morphing, as he reconstituted himself in response to various attacks. For scenes with relatively little movement and fairly close initial and final images, 2-D morphing was used. The scene in which Sarah changes to the T-1000 is an example of the 2-D morphing. To distract the viewer and make the scene look more realistic, different parts of the actors were morphed at different times. In other scenes, 3-D morphing was used--for example, when the T-1000 comes up out of the floor, and when he slips into the helicopter.
Star Trek 6: The Undiscovered Country (Industrial Light and Magic, 1991). The shape-changer, played by Iman, was morphed in 2-D to change into many different forms, including: adult female, young girl, monster (her true form), and Captain Kirk. The actors didn't move much in their morphing scenes, making life easier for the animators. The difference in sizes between the actors increased the difficulty of the morphing, requiring extreme warps between the images.
Plymouth Voyager Commercial (Pacific Data Images, 1991). In this advertisement, PDI used 2-D morphing to transform the 1990 Plymouth Voyager van into a 1991 model. Different parts of the car were morphed at different times. To reduce distortion of the background during the morph, parts of the vehicle were photographed separately for easy manipulation. For other elements, painting software was used to separate the elements and generate mattes.
Exxon Commercial (Pacific Data Images, 1991). In this advertisement, a moving car turns into a tiger. (Editor's Note: Refer to the images on the cover of this issue.) The tiger was filmed on a stage while the car was filmed on a mountain. Both shots were filmed using motion-control systems to allow exact retakes. The 2-D morph was done with flair; the car ripples as it morphs, the ripples becoming the stripes of the tiger.
Michael Jackson's "Black or White" Music Video (Pacific Data Images, 1991). A series of dancers are morphed into each other. Jackson turns into and out of a black panther while walking. Both were 2-D morphs. The dancers were edited in many orders to create the smoothest result possible. Then the end of each dancer's part was morphed onto the beginning of the next. The morphs were staggered so that different parts of each dancer changed at different times.
--V.H.
Rmorf
A Shareware Morphing Program for MS-DOS
While Valerie Hall's article and George Wolberg's source code provide a solid overview of morphing techniques, it's often necessary to work hands-on with an actual program to get a visceral understanding of otherwise abstract concepts.
To this end, we're providing Rmorf, a shareware program for MS-DOS written by Richard Goedeken (see "Availability," page 5). Figure 1 shows the Rmorf user interface, which brings up two images, side by side, and allows you to specify corresponding mesh points in each image. The program then produces a series of in-between images, as many as the frame count you specify. You'll need a way to display the resulting image sequence, however. One possibility is to use Dave Mason's DTA or Autodesk's shareware TGAFLI. EXE program, both available on CompuServe, to chain together a sequence of Targa files into a Flic file that can then be displayed with a player such as the READFLIC program Jim Kent included with his article, "The Flic File Format" (DDJ, March 1993). In addition to Rmorf and the sample image files, we are providing TGAFLIC.EXEand READFLIC.EXE electronically (see "Availability," page 5).
Although industrial-strength morphing programs rely on interpolated splines and sophisticated antialiasing, Richard opts for a straightforward implementation and raw speed. Here's his description of the genesis of Rmorf:
Recently, I was browsing through CompuServe Magazine and noticed a small item about a shareware morphing program, so I logged on and downloaded it. The sparse documentation described the basics of image morphing, of which I was previously unaware. In running the program, it took me several attempts to produce a moderately decent morph.
After spending several hours with the software, I got to know its problems: It crashed frequently and was painfully slow.
I've done some work in computer graphics, so, over the following week, I considered the possibility of writing a better morphing program. I sketched out, in my mind and on paper, how two images could be morphed together. It wasn't as difficult as I had first imagined.
All I had to do was make a transfer mesh for each frame, which would be somewhere in between the first two meshes (depending upon the frame number), warp each image to the transfer mesh, and mix the colors in a given ratio (again depending upon the frame number). After approximately nine days of working in my spare time, I completed a rough copy of Rmorf. Two days later the finishing touches were done, and I released the program. The end result of these 11 days of intense work are about 3700 lines of source code--2400 in assembly language and the remaining 1300 in C.
Because the morphing portion is in assembler--and because it uses integer rather than floating-point math--Rmorf can calculate a frame in 7.2 seconds on my 33-MHz 386, instead of the 30--50 minutes required by other programs.
Rmorf currently supports only Targa files, but upcoming versions will include .FLI and .GIF support. The unregistered version of Rmorf supports 320x200 resolution; the registered version handles 1024x768.
Richard, a junior in high school, has written a number of other programs (a product-inventory system, a personal accounting package, and an SVGA game) which are available from his software company. He can be reached on CompuServe at 70304,1065.
--editors
Figure 1: Running the RMORF program. Note the mesh lines on the images which define a correspondence between one image and another.
The Canonical Implementation in C
The two-pass mesh warping algorithm implemented by Douglas Smythe at Industrial Light and Magic in 1989 has become the canonical implementation of morphing. In his classic book on digital image warping (Digital Image Warping, IEEE Computer Society Press, 1990), George Wolberg provides a lucid and thorough explanation of how the algorithm works, as well presenting his own version implemented in C. For this issue of DDJ, however, Wolberg has shared with us a new implementation that's more self-contained and optimized than the one in his book; see Listing One (page 92). For a complete description, refer to Wolberg's book and the documentation that accompanies the electronic version of this article; see "Availability," page 5. Here's a quick summary of Wolberg's description.
The algorithm relies on the fact that an arbitrary, one-pass, spatial transformation can be decomposed into a computationally cheaper two-pass operation. This approach stems from the seminal 1980 paper by Ed Catmull and Alvy Ray Smith, which is generally applicable to affine and perspective transformations on planar and nonplanar surfaces. Here, a 2-D resampling problem is replaced by two orthogonal 1-D resampling stages.
The user must specify two sets of control points. Both input and output images are thereby partitioned into a mesh of patches. Each patch delimits an image region over which a continuous mapping function applies. Mapping between both images now becomes a matter of transforming each patch onto its counterpart in the second image--known as mesh warping. Since the mapping function is defined only at these discrete points, it's necessary to determine the mapping function over all points to perform the warp. The patches can be fitted with a bivariate function to realize a piecewise continuous mapping function. Wolberg adds:
The benefit of using a mesh derives from the simplicity in interpolating the new positions of intermediate points (between the mesh points). A bilinear or bicubic function can be used. We use a Catmull-Rom cubic spline to implement bicubic interpolation here.
The algorithm requires that all four edges of each mesh be frozen. This means that the first and last rows and columns all remain intact throughout the warp. Wolberg continues:
The input includes a source image I1 and two meshes, M1 and M2. Mesh M1 is used to select landmark positions in I1, and M2 identifies their corresponding positions in the output image. In this manner, arbitrary points in I1 can be "pulled" to new positions. Although the use of a parametric mesh might seem to place unnecessary constraints on the positions of these points, a large class of useful transformations is possible. It is important, though, that the mesh not self-intersect in order to avoid the image from folding upon itself.
The algorithm's first pass is responsible for resampling each row independently. It maps all (u,v) points in the source image I1 to their (x,v) coordinates in the intermediate image, thereby positioning each input point into its proper output column. The intermediate image is that whose x-coordinates are the same as those in I2 and whose y-coordinates are taken from I1. The second pass then resamples each column in the intermediate image, mapping every (x,v) to its final (x,y) position. Each point now lies in its proper row as well as its column.
The implementation presented here is for gray-scale images only. It is straightforward to extend the program to handle three-channel images (like RGB color images), by handling each channel separately.
Wolberg's implementation is in the classic UNIX style of a rudimentary command-line user interface. Wolberg explains:
The code is missing a program to help the user create and edit meshes interactively. A good mesh editor is a critical component to any mesh warping program. Such code falls outside of the scope of this presentation. A more full-featured implementation would allow the user to control the cross-dissolve schedule at each mesh point, as well as its position. This permits the intensities in different regions of the image to interpolate at different rates.
--editors
Copyright © 1993, Dr. Dobb's Journal
[LISTING ONE]<a name="01b8_0010">
_MORPHING IN 2D AND 3D_
by Valerie Hall
Source code written by George Wolberg and accompanies the
sidebar entitled "The Canonical Implementation in C"
/* ========================================================================
* meshwarp.h -- Header file for meshwarp.c. (C) 1993 by George Wolberg.
* ======================================================================*/
#include <stdio.h>
#define BW 0
#define MESH 1
#define MAX(A,B) ((A) > (B) ? (A) : (B))
#define MIN(A,B) ((A) < (B) ? (A) : (B))
typedef unsigned char uchar;
typedef struct { /* image data structure */
int width; /* image width (# cols) */
int height; /* image height (# rows) */
void *ch[2]; /* pointers to channels */
} imageS, *imageP;
extern void meshWarp(); /* extern decls for funcs in meshwarp.c */
extern void resample();
extern imageP readImage(); /* extern decls for funcs in util.c */
extern int saveImage();
extern imageP allocImage();
extern void freeImage();
/* =======================================================================
* morph.c - Generate a metamorphosis sequence. (C) 1993 by George Wolberg.
* =======================================================================*/
#include "meshwarp.h"
/*------- main: Collect user parameters and pass them to morph() ------- */
main(argc, argv)
int argc; char **argv;
{ int nframes;
char name[10];
imageP I1, I2;
imageP M1, M2;
if(argc != 7) /* make sure user invokes this program properly */
{ fprintf(stderr,
"Usage: morph src.bw dst.bw src.XY dst.XY frames name\n");
exit(1);
}
/*----------- read input image and meshes --------------*/
I1 = readImage(argv[1], BW); /* source image */
I2 = readImage(argv[2], BW); /* target image */
M1 = readImage(argv[3], MESH); /* source mesh */
M2 = readImage(argv[4], MESH); /* target mesh */
nframes = atoi(argv[5]); /* # frames */
strcpy(name, argv[6]); /* out basename */
/*----------- call morph -------------------------------*/
morph(I1, I2, M1, M2, nframes, name);
}
/* ------------------------------------------------------------------------
* morph: Generate a morph sequence of frames between images I1 and I2.
* Correspondence points among I1 and I2 are given in meshes M1 and M2.
* nframes frames are generated (including I1 and I2). The output is stored
* in files "basename_xxx.bw" where xxx are sequential 3-digit frame numbers.
*------------------------------------------------------------------------*/
void morph(I1, I2, M1, M2, nframes, basename)
imageP I1, I2, M1, M2; int nframes; char * basename;
{ int i, j, totalI, totalM;
double w1, w2;
char name[20];
uchar *p1, *p2, *p3;
float *x1, *y1, *x2, *y2, *x3, *y3;
imageP I3, Iw1, Iw2, M3;
/* allocate space for tmp images and mesh */
M3 = allocImage(M1->width, M1->height, MESH);
I3 = allocImage(I1->width, I1->height, BW);
Iw1 = allocImage(I1->width, I1->height, BW);
Iw2 = allocImage(I1->width, I1->height, BW);
/* eval total number of points in mesh (totalM) and image (totalI) */
totalM = M1->width * M1->height;
totalI = I1->width * I1->height;
/* copy 1st frame to basename_000.bw*/
sprintf(name, "%s_000.bw", basename);
saveImage(I1, name, BW);
printf("Finished Frame 0\n");
for(i=1; i<nframes-1; i++)
{ /* M3 <- linearly interpolate between M1 and M2 */
w2 = (double) i / (nframes-1);
w1 = 1. - w2;
/* linearly interpolate M3 grid */
x1 = (float *) M1->ch[0]; y1 = (float *) M1->ch[1];
x2 = (float *) M2->ch[0]; y2 = (float *) M2->ch[1];
x3 = (float *) M3->ch[0]; y3 = (float *) M3->ch[1];
for(j=0; j<totalM; j++)
{
x3[j] = x1[j]*w1 + x2[j]*w2;
y3[j] = y1[j]*w1 + y2[j]*w2;
}
/* warp I1 and I2 according to grid M3 */
meshWarp(I1, M1, M3, Iw1);
meshWarp(I2, M2, M3, Iw2);
/* cross-dissolve warped images Iw1 and Iw2 */
p1 = (uchar *) Iw1->ch[0];
p2 = (uchar *) Iw2->ch[0];
p3 = (uchar *) I3->ch[0];
for(j=0; j<totalI; j++)
p3[j] = p1[j]*w1 + p2[j]*w2;
/* save frame into file */
sprintf(name, "%s_%03d.bw", basename, i);
saveImage(I3, name, BW);
printf("Finished Frame %d\n", i);
}
/* copy last frame to basename_xxx.bw */
sprintf(name, "%s_%03d.bw", basename, i);
saveImage(I2, name, BW);
printf("Finished Frame %d\n", i);
}
/* ===========================================================================
* catmullRom.c - Catmull-Rom interpolating spline. (C) 1993 by George Wolberg
* =======================================================================*/
#include "meshwarp.h"
/* ------------------------------------------------------------------------
* catmullRom: Compute a Catmull-Rom spline passing thru the len1 points in
* arrays x1, y1, where y1 = f(x1). len2 positions on the spline are to be
* Their positions are given in x2. The spline values are stored in y2.
*------------------------------------------------------------------------*/
void catmullRom(x1, y1, len1, x2, y2, len2)
float *x1, *y1, *x2, *y2;
int len1, len2;
{ int i, j, dir, j1, j2;
double x, dx1, dx2;
double dx, dy, yd1, yd2, p1, p2, p3;
double a0y, a1y, a2y, a3y;
/* find direction of monotonic x1; skip ends */
if(x1[0] < x1[1]) /* increasing */
{ if(x2[0]<x1[0] || x2[len2-1]>x1[len1-1]) dir=0;
else dir = 1;
}
else /* decreasing */
{ if(x2[0]>x1[0] || x2[len2-1]<x1[len1-1]) dir=0;
else dir = -1;
}
if(dir == 0) /* error */
{ printf("catmullRom: Output x-coord out of range of input\n");
return;
}
/* p1 is first endpoint of interval
* p2 is resampling position
* p3 is second endpoint of interval
* j is input index for current interval
*/
if(dir==1) p3 = x2[0] - 1; /* force coefficient initialization */
else p3 = x2[0] + 1;
for(i=0; i<len2; i++)
{ p2 = x2[i]; /* check if in new interval */
if( (dir== 1 && p2>p3 ) ||
(dir== -1 && p2<p3 ))
{ if(dir) /* find the interval which contains p2 */
{ for(j=0; j<len1 && p2>x1[j]; j++) ;
if(p2 < x1[j]) j--;
}
else
{ for(j=0; j<len1 && p2<x1[j]; j++) ;
if(p2 > x1[j]) j--;
}
p1 = x1[j]; /* update 1st endpt */
p3 = x1[j+1]; /* update 2nd endpt */
/* clamp indices for endpoint interpolation */
j1 = MAX(j-1, 0);
j2 = MIN(j+2, len1-1);
/* compute spline coefficients */
dx = 1.0 / (p3 - p1);
dx1 = 1.0 / (p3 - x1[j1]);
dx2 = 1.0 / (x1[j2] - p1);
dy = (y1[j+1] - y1[ j ]) * dx;
yd1 = (y1[j+1] - y1[ j1]) * dx1;
yd2 = (y1[j2 ] - y1[ j ]) * dx2;
a0y = y1[j];
a1y = yd1;
a2y = dx * ( 3*dy - 2*yd1 - yd2);
a3y = dx*dx*(-2*dy + yd1 + yd2);
}
/* use Horner's rule to calculate cubic polynomial */
x = p2 - p1;
y2[i] = ((a3y*x + a2y)*x + a1y)*x + a0y;
}
}
/* ========================================================================
* meshwarp.c -- Mesh warping program. Copyright (C) 1993 by George Wolberg.
* =======================================================================*/
#include "meshwarp.h"
/* -------------------------------------------------------------------------
* meshWarp: Warp I1 with correspondence points given in meshes M1 and M2.
* Result goes in I2. Based on Douglas Smythe's algorithm at ILM.
*------------------------------------------------------------------------*/
void meshWarp(I1, M1, M2, I2)
imageP I1, I2, M1, M2;
{ int I_w, I_h, M_w, M_h;
int x, y, u, v, n;
float *x1, *y1, *x2, *y2;
float *xrow, *yrow, *xcol, *ycol, *coll, *indx, *map;
uchar *src, *dst;
imageP Mx, My, I3;
I_w = I1->width; I_h = I1->height;
M_w = M1->width; M_h = M1->height;
/* alloc enough memory for a scanline along the longest dimension */
n = MAX(I_w, I_h);
indx = (float *) malloc(n * sizeof(float)); /* should check if err */
xrow = (float *) malloc(n * sizeof(float));
yrow = (float *) malloc(n * sizeof(float));
map = (float *) malloc(n * sizeof(float));
/* create table of x-intercepts for source mesh's vert splines */
Mx = allocImage(M_w, I_h, MESH);
for(y=0; y < I_h; y++) indx[y] = y;
for(u=0; u < M_w; u++) /* visit each vert spline */
{ /* store col as row for spline function */
xcol = (float *) M1->ch[0] + u;
ycol = (float *) M1->ch[1] + u;
coll = (float *) Mx->ch[0] + u;
/* scan-convert vert splines */
for(v=0; v < M_h; v++, xcol+=M_w) xrow[v] = *xcol;
for(v=0; v < M_h; v++, ycol+=M_w) yrow[v] = *ycol;
catmullRom(yrow, xrow, M_h, indx, map, I_h);
/* store resampled row back into column */
for(y=0; y < I_h; y++, coll+=M_w) *coll = map[y];
}
/* create table of x-intercepts for dst mesh's vert splines */
for(u=0; u < M_w; u++) /* visit each vert spline */
{ /* store column as row for spline fct */
xcol = (float *) M2->ch[0] + u;
ycol = (float *) M2->ch[1] + u;
coll = (float *) Mx->ch[1] + u;
/* scan-convert vert splines */
for(v=0; v < M_h; v++, xcol+=M_w) xrow[v] = *xcol;
for(v=0; v < M_h; v++, ycol+=M_w) yrow[v] = *ycol;
catmullRom(yrow, xrow, M_h, indx, map, I_h);
/* store resampled row back into column */
for(y=0; y < I_h; y++, coll+=M_w) *coll = map[y];
}
/*------------ first pass: warp x using tables in Mx --------*/
I3 = allocImage(I_w, I_h, BW);
x1 = (float *) Mx->ch[0];
x2 = (float *) Mx->ch[1];
src = (uchar *) I1->ch[0];
dst = (uchar *) I3->ch[0];
for(x=0; x < I_w; x++) indx[x] = x;
for(y=0; y < I_h; y++)
{ /* fit spline to x-intercepts; resample over all cols */
catmullRom(x1, x2, M_w, indx, map, I_w);
/* resample source row based on map */
resample(src, I_w, 1, map, dst);
/* advance pointers to next row */
src += I_w;
dst += I_w;
x1 += M_w;
x2 += M_w;
}
freeImage(Mx);
/* create table of y-intercepts for intermediate mesh's hor splines */
My = allocImage(I_w, M_h, MESH);
x1 = (float *) M2->ch[0];
y1 = (float *) M1->ch[1];
y2 = (float *) My->ch[0];
for(x=0; x < I_w; x++) indx[x] = x;
for(v=0; v < M_h; v++) /* visit each horz spline */
{ /* scan-convert horz splines */
catmullRom(x1, y1, M_w, indx, y2, I_w);
x1 += M_w; /* advance pointers to next row */
y1 += M_w;
y2 += I_w;
}
/* create table of y-intercepts for dst mesh's horz splines */
x1 = (float *) M2->ch[0];
y1 = (float *) M2->ch[1];
y2 = (float *) My->ch[1];
for(v=0; v < M_h; v++) /* visit each horz spline */
{ /* scan-convert horz splines */
catmullRom(x1, y1, M_w, indx, y2, I_w);
x1 += M_w; /* advance pointers to next row */
y1 += M_w;
y2 += I_w;
}
/*----------------------- second pass: warp y ------------------*/
src = (uchar *) I3->ch[0];
dst = (uchar *) I2->ch[0];
for(y=0; y < I_h; y++) indx[y] = y;
for(x=0; x < I_w; x++)
{ /* store column as row for spline fct */
xcol = (float *) My->ch[0] + x;
ycol = (float *) My->ch[1] + x;
for(v=0; v < M_h; v++, xcol+=I_w) xrow[v] = *xcol;
for(v=0; v < M_h; v++, ycol+=I_w) yrow[v] = *ycol;
/* fit spline to y-intercepts; resample over all rows */
catmullRom(xrow, yrow, M_h, indx, map, I_h);
/* resample source column based on map */
resample(src, I_h, I_w, map, dst);
/* advance pointers to next column */
src++; dst++;
}
freeImage(My); freeImage(I3); free((char *) indx);
free((char *) xrow); free((char *) yrow); free((char *) map);
}
/* ------------------------------------------------------------------------
* resample: Resample the len elements of src (with stride offst) into dst
* according to the spatial mapping given in xmap. Perform linear interpola-
* tion for magnification and box filtering (unweighted averaging) for
* minification. Based on Fant's algorithm (IEEE Comp. Graphics & Appl. 1/86)
*------------------------------------------------------------------------*/
void resample(src, len, offst, xmap, dst)
uchar *src, *dst; float *xmap; int len, offst;
{ int u, x, v0, v1;
double val, sizfac, inseg, outseg, acc, inpos[1024];
/* precompute input index for each output pixel */
for(u=x=0; x<len; x++)
{ while(xmap[u+1]<x) u++;
inpos[x] = u + (double) (x-xmap[u]) / (xmap[u+1]-xmap[u]);
}
inseg = 1.0;
outseg = inpos[1];
sizfac = outseg;
acc = 0.;
v0 = *src; src += offst;
v1 = *src; src += offst;
for(u=1; u<len; )
{ val = inseg*v0 + (1-inseg)*v1;
if(inseg < outseg)
{ acc += (val * inseg);
outseg -= inseg;
inseg = 1.0;
v0 = v1;
v1 = *src;
src += offst;
}
else
{ acc += (val * outseg);
acc /= sizfac;
*dst = (int) MIN(acc, 0xff);
dst += offst;
acc = 0.;
inseg -= outseg;
outseg = inpos[u+1] - inpos[u];
sizfac = outseg;
u++;
}
}
}
