Global Error Minimization In Image Mosaicking Using Graph Laplacians And Its Applications In Microscopy

Abstract

An image mosaicking method includes performing pairwise registration of a plurality of tiles (101), determining absolute homographies for each of the plurality of tiles according to the pairwise registration (102B), and performing a blending of the plurality of tiles to obtain a stitched image according to the absolute homographies (103)

Description

BACKGROUND

1. Technical Field

The present disclosure relates to image editing, and more particularly to an image mosaicking method.

2. Discussion of Related Art

Several applications, such as underwater sea-floor mapping, microscopy, etc., require the mosaicking of images or tiles acquired by a camera as the camera traverses an unknown trajectory in three dimensional (3D) space.

A homography relates the image coordinates of a point in each tile to those of a reference tile provided the image surface is planar. One exemplary approach in such applications is to first perform pairwise alignment of the tiles that have imaged common features (points) in order to recover a homography relating the tile pair and to then find the global set of homographies relating each individual tile to a reference tile such that the homographies relating all tiles pairs are kept as consistent as possible.

Using these global homographies, one can generate a mosaic of the entire scene.

BRIEF SUMMARY

According to an embodiment of the present disclosure, an image mosaicking method includes performing pairwise registration of a plurality of tiles, determining absolute homographies for each of the plurality of tiles according to the pairwise registration, and performing a blending of the plurality of tiles to obtain a stitched image according to the absolute homographies.

According to an embodiment of the present disclosure, an image mosaicking method includes performing pairwise registration of a plurality of tiles, determining absolute locations for each of the plurality of tiles according to the pairwise registration, and performing a blending of the plurality of tiles to obtain a stitched image according to the absolute locations.

According to an embodiment of the present disclosure, an image mosaicking method includes performing pairwise registration of a plurality of tiles, receiving prior homography information for each of the plurality of tiles, determining absolute homographies for each of the plurality of tiles according to the pairwise registration and the prior homography information of each of the plurality of tiles, and performing a blending of the plurality of tiles to obtain a stitched image according to the absolute homographies.

BRIEF DESCRIPTION OF THE DRAWINGS

Preferred embodiments of the present disclosure will be described below in more detail, with reference to the accompanying drawings:

FIGS. 1A-B are flow diagrams of an image mosaicking method according to an embodiment of the present disclosure;

FIG. 2 shows exemplary X-Y coordinates for a component of a specimen comprising 79×34 tiles according to an embodiment of the present disclosure;

FIG. 3 is an exemplary stitched image from the 79×34 tiles according to an embodiment of the present disclosure; and

FIG. 4 is a system for executing an image mosaicking method according to an embodiment of the present disclosure.

DETAILED DESCRIPTION OF PREFERRED EMBODIMENTS

According to an embodiment of the present disclosure, a global homography may be determined from the pairwise homographies via global error minimization under an assumption that a 3D object being imaged is planar.

According to an embodiment of the present disclosure, an image mosaicking method (see FIGS. 1A-B) includes performing pairwise registration for all neighboring tiles (101). The pairwise registration may include maximizing a similarity measure such as a normalized cross-correlation (NCC) function. By way of example, a normalized cross-correlation may be used as a similarity measure, which normalizes the image and feature vectors to unit length, yielding a cosine-like correlation coefficient. This operation can be performed in a transform domain.

The method includes obtaining or determining absolute locations for each tile (102A—FIG. 1A) or absolute homographies for each tile (102B—FIG. 1B) (e.g., in the case of a general transformation). Referring to FIG. 1A, the absolute locations may be obtained by solving a linear system such that the global measure of registration in (6) is minimized. For example, the linear system can be solved using conjugate gradients with an incomplete Cholesky preconditioner.

The method further includes performing a blending to obtain a stitched image with non-overlapping tiles (103). The blending may be performed using techniques including an alpha-blending approach to assign different tiles various transparency settings. A weight α can be found by raising the normalized distance of a pixel from the edge of its tile to an exponent.

The following describes exemplary methods corresponding to blocks 102A and 102B for determining absolute locations and absolute homographies for each tile.

According to an embodiment of the present disclosure, it may be assumed that N images or tiles (herein after tiles) have been acquired by a camera. A global homography for tile i may be expressed as

$M_i=K\left(R_i-\frac{1}{d_i}p_i,n^T\right)K^{-1},$

where K denotes the intrinsic camera matrix, R_i denotes the rotation of the camera axis of tile i with respect to a reference tile, p_i denotes the translation of the tile i's center with respect to a reference tile's center, n and d_i are the normal vector of the plane and the distance to the plane, respectively.

A connectivity graph may be formed between all acquired tiles by planning an edge between tiles that share a border. Let Nbd(i) denote the set of tiles that share any overlap with tile i. Let {tilde over (M)}_ij denote the pairwise homography between tiles i and j found via pairwise alignment. It should be understood that a homography is an invertible transformation between projective planes that maps straight lines to straight lines. To determine global homographies, the following objective function may be minimized:

$\begin{array}{cc}E\left(M_1,\dots ,M_N\right)=\frac{1}{2}\sum _{i=1}^N\sum _{j\in \mathrm{Nbd}\left(i\right)}{\mathrm{dist}\left(M_i,M_{\mathrm{ij}},M_j\right)}^2& \left(1\right)\end{array}$

The squared distance between any two homography matrices A and B may be taken as Σ_mΣ_n(A_mn−B_mn)²=Trace((A−B)(A−B)^T), which is a squared Frobenius norm of A−B. Alternative matrix distance measures such as those based upon the Log-Euclidean transform or the Baker-Campbell-Hausdorff formula may also be used. Note that it may be desirable to use two symmetric distance measures subject to dist(M_i,{tilde over (M)}_ij,M_j)²=dist(M_j,{tilde over (M)}_ji,M_i)², a condition that may be satisfied by the squared Frobenius norm as shown below:

$\begin{array}{cc}\begin{array}{c}{\left[\mathrm{Fro}\left(M_j-M_{\mathrm{ji}}M_i\right)\right]}^2={\left[\mathrm{Fro}\left({\tilde{M}}_{\mathrm{ji}}\left(M_i-{\tilde{M}}_{\mathrm{ji}}M_i\right)\right)\right]}^2\\ =\mathrm{Trace}\left({\tilde{M}}_{\mathrm{ji}}\left(M_i-{\tilde{M}}_{\mathrm{ij}}M_j\right){\left(M_i-{\tilde{M}}_{\mathrm{ij}}M_j\right)}^TM_{\mathrm{ji}}^T\right)\\ =\mathrm{Trace}\left(\left(M_i-{\tilde{M}}_{\mathrm{ij}}M_j\right){\left(M_i-{\tilde{M}}_{\mathrm{ij}}M_j\right)}^T\right)\\ ={\left[\mathrm{Fro}\left(M_i-{\tilde{M}}_{\mathrm{ij}}M_j\right)\right]}^2.\end{array}& \left(2\right)\end{array}$

The first equality in (2) assumes inverse consistency in the pairwise registration, i.e., {tilde over (M)}_ij={tilde over (M)}_ji⁻¹. The second-to-last equality follows from the fact that a similarity transformation does not affect the eigenvalues or the trace of a matrix.

The minimum of the objective function in (1) can be found by setting its gradient to zero. The gradient of the objective function in (1) with respect to homography matrix M_i may be determine as:

$\begin{array}{cc}{\nabla }_{M_i}E\left(M_1,\dots ,M_N\right)=\sum _{j\in \mathrm{Nbd}\left(i\right)}{\nabla }_{M_i}{\left[\mathrm{Fro}\left(M_i-{\tilde{M}}_{\mathrm{ij}}M_j\right)\right]}^2& \left(3\right)\end{array}$

Note that

∇_Mi[Fro(M_i−{tilde over (M)}_ijM_j)]=∇_MiTrace((M_i−{tilde over (M)}_ijM_j)(M_i−{tilde over (M)}_ijM_j)^T)=2(M_i−{tilde over (M)}_ijM_j).

Assuming that a reference homography M_N is fixed, the remaining N−1 homographies may be solved for by solving the following set of N−1 simultaneous matrix equations:

$\begin{array}{cc}\sum _{j\in \mathrm{Nbd}\left(i\right)}\left(M_i-{\tilde{M}}_{\mathrm{ij}}M_j\right)=0& \left(4\right)\end{array}$

Each such matrix equation is equivalent to a set of 9 scalar equations for the individual matrix elements. Other than the reference homography M_N, the (a,b)^th elements M_i^(a,b) of the set of the remaining 3×3 matrices {M_i} may be concatenated into a vector x and then this vector may be determined by solving the simultaneous equations. This is equivalent to solving the linear system Ax=y, where A is a 9(N−1)×9(N−1) matrix. The exact form of A and y is given below:

$\begin{array}{cc}\begin{array}{cc}{A}_{\mathrm{iab};\mathrm{iab}}=\sum _{j\in \mathrm{Nbd}\left(i\right)}^{}1& \forall i\in 1,\dots ,N-1\forall a,b\in \left\{1,2,3\right\}\\ A_{\mathrm{iab};\mathrm{iab}}0,& \forall a,b,c,d\in \left\{1,2,3\right\}\mathrm{and}\left(c\ne a\mathrm{or}d\ne b\right)\\ A_{\mathrm{iab};\mathrm{iab}}=-{\tilde{M}}_{\mathrm{ij}}^{ac},& \mathrm{if}j\in \mathrm{Nbd}\left(i\right)\mathrm{and}j\ne i\forall a,b,c\in \left\{1,2,3\right\}\\ A_{\mathrm{iab};\mathrm{iab}}=0,& \mathrm{if}j\in \mathrm{Nbd}\left(i\right)\mathrm{and}j\ne i\forall a,b,c,d\in \left\{1,2,3\right\}\mathrm{and}d\ne b\\ A_{\mathrm{iab};\mathrm{iab}}=0,& \mathrm{if}j\notin \mathrm{Nbd}\left(i\right)\mathrm{and}j\ne i\forall a,b,c,d\in \left\{1,2,3\right\}\\ y_{\mathrm{iab}}=\sum _{c=1}^3{\tilde{M}}_{\mathrm{iN}}^{ac}M_N^{\mathrm{cb}}& \mathrm{if}N\in \mathrm{Nbd}\left(i\right)\\ y_{\mathrm{iab}}=0& \mathrm{if}N\notin \mathrm{Nbd}\left(i\right)\end{array}& \left(5\right)\end{array}$

Note that the above derivation assumes a general 3×3 homography matrix without any constraints on its form. However, general homography matrices have at most 8 degrees of freedom (up to scale), with 3 degrees of freedom each for the rotations and translations, and additional degrees of freedom arising from the variations in the intrinsic camera matrix during the scanning process, due to a changing focus, for example. This 8-degree homography problem may be solved by setting M_i^(3,3)=1 for all homographies and dropping the corresponding rows and columns of A to get a 8(N−1)×8(N−1) matrix. The vector y may also be updated with the appropriate terms transferred from the original left-hand side.

In case the intrinsic matrix is fixed, but both translations and rotations are allowed, the homography matrix has the form M_i=R_i−t_in^T and has only six degrees of freedom (the scalar factor d_i has been absorbed into t_i). In this case, R_i and t_i may be recovered from the M_i found by the 8(N−1)×8(N−1) linear system solver.

For other special cases of homographies such as pure rotations, which may occur when the camera is rotated from a single vantage point to cover a large wall painting, these constraints may be imposed on the solution found via the 8(N−1)×8(N−1) linear system solver using the singular value decomposition (SVD) based rotation extraction approach.

For the case of 2D translations, only two independent equations need be obtained per matrix equation and additional constraints on the resulting linear system solution may not be needed, as demonstrated by the first-principles. An exemplary case of multiple connected components in the graph connectivity is also described herein.

Exemplary image mosaicking methods can be used in histological microscopic imaging. Since the motion reported by a slide-scanning microscopy X-Y stage in histophology imaging may be unreliable and an image-based mosaicking method may be useful in determining an accurate composite image. Assuming that the microscope moves in a plane parallel to the planar specimen (assuming this plane to be substantially normal to the Z-axis) and K=I, the homography M_i may be reduced to

$M_i=I-\frac{1}{d}t_ie_z^T.$

Then

$M_i-{\tilde{M}}_{\mathrm{ij}}M_j=\frac{-1}{d}\left(p_i-p_j-{\tilde{p}}_{\mathrm{ij}}\right)e_z^T$

and the objective function to be minimized can be re-phrased as

$\begin{array}{cc}E\left(p_1,\dots ,p_N\right)=\frac{1}{2}\sum _i\sum _{j\in \mathrm{Nbd}\left(i\right)}{p_i-p_j-{\tilde{p}}_{\mathrm{ij}}}^2& \left(6\right)\end{array}$

An exemplary linear system is described herein, which when solved minimizes eq. (6), in lieu of eq. (5).

First-principles derivation; an exemplary first-principles derivation may be based on a Graph Laplacian. For example, let r,c denote the row and column number for any acquired image tile and I_r,c denote the acquired tile. Let {tilde over (p)}_r1,c1;r2,c2 denote the translation obtained by applying pairwise registration between neighbors I_r1,c1 and I_r2,c2. Note that {tilde over (p)}_r1,c1;r2,c2=−{tilde over (p)}_r2,c2;r1,c1. Let p_r,c denote the global location of the time at r,c. Note that the indices for the N acquired image (N=RC) can be obtained by lexicographic ordering so that both the row and column number is encoded in a single index i . Based upon the 8-connected neighborhood information, a directed graph may be formed with N vertices and with the edge weights determined by the registration-based translation recovered between the neighboring images.

The location of the last tile p_R,C (equivalently p_N) is fixed. Revised locations ({circumflex over (p)}_1,1, . . . , {circumflex over (p)}_R,C,C-1) of the remaining RC−1 tiles may be found be minimizing the objective:

$\begin{array}{cc}E\left(p_1,\dots ,p_N\right)=\sum _{r=1}^{R-1}\sum _{c=1}^{C-1}\frac{1}{2}\left({\left[\left(p_{r,c}-p_{r+1,c}\right)-{\tilde{p}}_{r,c;r+1,c}\right]}^2+{\left[\left(p_{r,c}-p_{r,c+1}\right)-{\tilde{p}}_{r,c;r,c+1}\right]}^2+{\left[\left(p_{r,c}-p_{r+1,c+1}\right)-{\tilde{p}}_{r,c;r+1,c+1}\right]}^2+{\left[\left(p_{r,c+1}-p_{r+1,c}\right)-{\tilde{p}}_{r,c+1;r+1,c}\right]}^2\right)+\sum _{r=1}^{R-1}\frac{1}{2}\left[\left(p_{r,C}-p_{r+1,C}\right)-{\tilde{p}}_{r,C;r+1,C}\right]+\sum _{c=1}^{C-1}{\frac{1}{2}\left[\left(p_{R,c}-p_{R,c+1}\right)-{\tilde{p}}_{R,c;R,c+1}\right]}^2=\sum _{i=1}^N\sum _{j\in \mathrm{Nbd}\left(i\right),i<j}{\frac{1}{2}\left[\left(p_i=p_j\right)-{\tilde{p}}_{i;j}\right]}^2& \left(7\right)\end{array}$

The partial derivative of the above objective with respect to p_i, i=1, . . . , N−1 may be given by:

$\begin{array}{cc}\frac{\partial E}{\partial p_i}=p_i\left(\sum _{j\in \mathrm{Nbd}\left(i\right)}1\right)-\sum _{j\in \mathrm{Nbd}\left(i\right)}p_j-\sum _{j\in \mathrm{Nbd}\left(i\right)}{\tilde{p}}_{i;j}& \left(8\right)\end{array}$

Therefore, if x=[p₁^x, . . . , p_N-1^x] (note that the dimension of x is (N−1) and it is obtained by concatenating all the unknown X-coordinates), the revised X-coordinates may be obtained by solving the linear system Ax=y, where:

$\begin{array}{cc}\begin{array}{cc}{A}_{\mathrm{ii}}=\sum _{j\in \mathrm{Nbd}\left(i\right)}1& \\ A_{\mathrm{ij}}=-1,& \mathrm{if}j\in \mathrm{Nbd}\left(i\right)\mathrm{and}j\ne i\\ A_{\mathrm{ij}}=0,& \mathrm{if}j\notin \mathrm{Nbd}\left(i\right)\mathrm{and}j\ne i\\ y_i=\sum _{j\in \mathrm{Nbd}\left(i\right)}{\tilde{p}}_{i;j}^X& \mathrm{if}N\notin \mathrm{Nbd}\left(i\right)\\ y_i=p_N^X+\sum _{j\in \mathrm{Nbd}\left(i\right)}{\tilde{p}}_{i;j}^X& \mathrm{if}N\in \mathrm{Nbd}\left(i\right)\end{array}& \left(9\right)\end{array}$

Revised Y-coordinates may also be obtained by solving a similar linear equation.

Note that with full 8NN connectivity (as described above), the above matrix A of dimension (N−1)×(N−1) is invertible because the N^th tile is included in the neighborhood connectivity. The matrix A is obtained by dropping the last row and column of an N×N connectivity graph Laplacian L=D−W , where D is a diagonal matrix such that D_ii=Σ_j∈Ndb(i) 1 and W is a binary adjacency matrix such that W_ij=1 if j∈Nbd(1) and j≠i.

Note that smaller weights may be assigned to the diagonal edges within the 8NN connectivity. In the example above the smaller weights are not used because the diagonal input translations appear to be as reliable as the vertical or horizontal translations.

Also, the above linear system can be modified to accommodate holes within the specimen by dropping the corresponding images, appropriately changing the 8NN connectivity and modifying the matrix A and the column vector b. If the holes end up disconnecting the specimen, then reference pixel coordinates need to be supplied for each connected component while solving one linear system per component. Note that partial holes that only remove a part of the 8NN connectivity may be accommodated.

Referring to FIG. 1A, all connected components may be found in N steps using depth-first-search or breadth-first-search using the binary adjacency matrix W (block 102Ai) and determining the matrices A^c corresponding to each component (block 102Aii) until all tiles have been searched (block 102Aiii). According to an exemplary embodiment of the present disclosure, where each matrix A^c is positive-definite, a Cholesky decomposition or a preconditioning conjugate gradient method may be used to solve a corresponding linear system.

Existing approaches do not consider a graph Laplacian and only work with an E×N edge-incidence matrix E with non-zero row entries E_ki=1 and E_kj=−1 if i<j and j∈E Nbd(i). Note that L=E^TE, where E is the edge-incidence matrix and A=F^TF, where F is obtained by dropping the last column of E.

Further, if a column vector t of all the E pairwise translations (X-components) and a column vector a of all the tile positions (X-components) are formed, then t−Ea. The vector x may be obtained by dropping the last entry of a and adding the X-component of the N^th tile position to the appropriate elements of t resulting in the linear system {tilde over (t)}=Fx.

According to an exemplary embodiment of the present disclosure, the linear system {tilde over (t)}=Fx may be obtained with no spectral SVD decompositions.

FIG. 2 shows global positions, e.g., 201, determined for the largest component of a specimen comprising 79×34 tiles that contained both partial and full holes. On account of the holes, only 7942 pairwise translation out of a total possible 11207 8NN-connectivity edges were fed to a system embodiment a method according to the present disclosure and of these the system eliminated more than 50 components, e.g., relatively small components. Note the missing tiles (holes), e.g., 202, in the center of the specimen in FIG. 2. Also note the slight tilt in the X-coordinates possible on account of the backlash that occurs when the X-Y stage returns to its left-most X-coordinate after scanning an entire row. The complete stitched image 301 is shown in FIG. 3.

For speed-up, the matrix inverse of A may be pre-determined for different specimen sizes (e.g., different R and C) before acquiring any images. A run-time of about 1 minute has been observed on 80×40 images (tiles) for the global registration and a run-time of a few hours to obtain all pairwise translations.

Variations and Extensions

As described above, only one reference pixel position is used per connected component and this reference position may be the corresponding coordinates reported by the X-Y image. In the case of multiple connected components, extrapolation may be used in order to obtain reference coordinates that are a corrected version of the coordinates reported by the X-Y image.

An exemplary alternative to the connected components analysis (e.g., shown in FIG. 1A, block 102A) including using a prior term, e.g., location information reported by the slide-scanning microscopy X-Y stage, added to the objective function in (7) and enforcing an agreement between the unknown global coordinates and the coordinates reported by the X-Y image stage (block 102Bi). The prior term may have a weight of 0.5, though the value may be adjusted according to an application. The revised objective function for a quadratic prior term may be given by:

$\begin{array}{cc}E\left(p_1,\dots ,p_N\right)=\frac{1}{2}\sum _{i=1}^N\sum _{j\in \mathrm{Nbd}\left(i\right)}{p_i-p_j-{\tilde{p}}_{\mathrm{ij}}}^2+\frac{\lambda }{2}\sum _{i=1}^N{p_i-{\tilde{p}}_i}^2& \left(10\right)\end{array}$

where {tilde over (p)}_i denotes the prior coordinates, e.g., those reported by the X-Y image stage (block 102Bii).

It should be understood that the connected components analysis shown in FIG. 1A and the objective function using prior information shown in FIG. 1B are alternatives, such that the connected components analysis may be implemented in the homography based example of FIG. 1B and the objective function using prior information may be implemented in the location based example of FIG. 1A.

The above microscopy mosaicking method determines 2D global positions for a planar specimen from 2D translations obtained via pairwise registration. However, it is possible to extend it to determine 3D global positions from 3D translations as in the confocal microscopy volumetric stitching application.

The determination of global coordinates may be made by L1 or L2 metrics. The L1 method is relatively outlier-resistant. For the L1 method, instead of a linear system solver, an iterative revision (e.g., optimization) may be used.

Moreover, the application scope of an exemplary general mosaicking method can be extended beyond microscopy to other applications such as geographical information systems or sea bed exploration wherein the overlapping parts of the acquired neighboring images (tiles) differ by a translation or by more general transformations. Pairwise transformations other than translations would use the objective function in (1). A solution may be obtained by zeroing the gradient of (1). Note that the same graph connectivity plays a role in the solutions for both the generic and the microscopy scenarios, as seen in equations (5) and (9).

It is to be understood that embodiments of the present disclosure may be implemented in various forms of hardware, software, firmware, special purpose processors, or a combination thereof. In one embodiment, a software application program is tangibly embodied on a non-transitory computer-readable storage medium, such as a program storage device or program product, with an executable program stored thereon. The application program may be uploaded to, and executed by, a machine comprising any suitable architecture.

Referring now to FIG. 4, according to an embodiment of the present disclosure, a computer system (block 401) for an image mosaicking method includes, inter alia, a central processing unit (CPU) (block 402), a memory (block 403) and an input/output (I/O) interface (block 404). The computer system (block 401) is generally coupled through the I/O interface (block 404) to a display (block 405) and various input devices (block 406) such as a mouse, keyboard and a slide-scanning microscopy X-Y stage. The support circuits can include circuits such as cache, power supplies, clock circuits, and a communications bus. The memory (block 403) can include random access memory (RAM), read only memory (ROM), disk drive, tape drive, etc., or a combination thereof. The present invention can be implemented as a routine (block 407) that is stored in memory (block 403) and executed by the CPU (block 402) to process the signal from the signal source (block 408). For example, the signal may include image information from a camera, which may be stored to memory (block 403) As such, the computer system (block 401) is a general purpose computer system that becomes a specific purpose computer system when executing the routine (block 407) of the present disclosure.

The computer platform (block 401) also includes an operating system and micro instruction code. The various processes and functions described herein may either be part of the micro instruction code or part of the application program (or a combination thereof) which is executed via the operating system. In addition, various other peripheral devices may be connected to the computer platform such as an additional data storage device and a printing device.

It is to be further understood that, because some of the constituent system components and method steps depicted in the accompanying figures may be implemented in software, the actual connections between the system components (or the process steps) may differ depending upon the manner in which the system is programmed. Given the teachings of the present disclosure provided herein, one of ordinary skill in the related art will be able to contemplate these and similar implementations or configurations of the present disclosure.

Having described embodiments for an image mosaicking method, it is noted that modifications and variations can be made by persons skilled in the art in light of the above teachings. It is therefore to be understood that changes may be made in embodiments of the present disclosure that are within the scope and spirit thereof.

Claims

1. A computer program product embodying instructions executed by a processor to perform an image mosaicking method, the method comprising: performing pairwise registration of a plurality of tiles;determining absolute homographies for each of the plurality of tiles according to the pairwise registration; andperforming a blending of the plurality of tiles to obtain a stitched image according to the absolute homographies.

2. The computer program product of claim 1, wherein performing the pairwise registration comprises maximizing a similarity measure.

3. The computer program product of claim 1, wherein determining the absolute homographies for each tile comprises solving a linear system such that a global measure of registration is minimized.

4. The computer program product of claim 1, wherein determining the absolute homographies for each tile comprises determining connected components given the pairwise registration of the plurality of tiles.

5. The computer program product of claim 1, wherein determining the absolute homographies for each tile is performed using a weight enforcing a relationship between unknown global coordinates and prior homography information for each of the plurality of tiles.

6. The computer program product of claim 1, wherein determining the absolute homographies for each tile comprises: receiving prior homography information of each of the plurality of tiles; andrevising the pairwise registration according to an objective function for an edge incidence matrix of the plurality of tiles and the prior homography information.

7. The computer program product of claim 1, wherein performing the blending comprises an alpha-blending approach wherein a weight alpha is determined by raising a normalized distance of a pixel from an edge of its tile to an exponent.

8. A computer program product embodying instructions executed by a processor to perform an image mosaicking method, the method comprising: performing pairwise registration of a plurality of tiles;determining absolute locations for each of the plurality of tiles according to the pairwise registration; andperforming a blending of the plurality of tiles to obtain a stitched image according to the absolute locations.

9. The computer program product of claim 8, wherein performing the pairwise registration comprises maximizing a similarity measure.

10. The computer program product of claim 8, wherein determining the absolute locations for each tile comprises solving a linear system such that a global measure of registration is minimized.

11. The computer program product of claim 8, wherein determining the absolute locations for each tile is performed using a weight enforcing a relationship between unknown global coordinates and prior location information for each of the plurality of tiles.

12. The computer program product of claim 8, wherein determining the absolute locations for each tile comprises: receiving prior location information of each of the plurality of tiles; andrevising the pairwise registration according to an objective function for an edge incidence matrix of the plurality of tiles and the prior location information.

13. The computer program product of claim 8, wherein determining the absolute locations for each tile comprises an iterative revision of the pairwise registration.

14. The computer program product of claim 8, wherein performing the blending comprises an alpha-blending, wherein a weight alpha is determined by raising a normalized distance of a pixel from an edge of its tile to an exponent.

15. A non-transitory computer-readable storage medium with an executable program stored thereon, wherein the program instructs a microprocessor to perform an image mosaicking method, the method comprising: performing pairwise registration of a plurality of tiles;receiving prior homography information for each of the plurality of tiles;determining absolute homographies for each of the plurality of tiles according to the pairwise registration and the prior homography information of each of the plurality of tiles; andperforming a blending of the plurality of tiles to obtain a stitched image according to the absolute homographies.

16. The computer-readable storage medium of claim 15, wherein performing the pairwise registration comprises maximizing a similarity measure.

17. The computer-readable storage medium of claim 15, wherein the location of each of the tiles is associated with a weight.

18. The computer-readable storage medium of claim 15, wherein determining the absolute homographies for each of the plurality of tiles comprises revising the pairwise registration according to an objective function for an edge incidence matrix of the plurality of tiles.

19. The computer-readable storage medium of claim 15, wherein performing the blending comprises an alpha-blending approach wherein a weight alpha is determined by raising a normalized distance of a pixel from an edge of its tile to an exponent.