Sensors 2009, 9, 9629-9665; doi10 3390/s91209629

O PEN ACCESS

sensors

ISSN 1424-8220 www m dpicom /pumal/sensors

Article

Bundle Block Adjustm ent with 3D Natural Cubic Splines

W on H ee Lee and K iyun Yu *

Departm entof Civil and Environm ental Engineering, Seoul NatonalUniveisity, 599 G w anak-ro, Gwanak-gu, Seoul 151-742, Korea; E-M ail: wlee2432@ snuackr

* Authorto whom correspondenceshould be addressed; E-Mail: kiyun@ snuackr; Tel:: +82-2-880-1355; Fax: +82-2-889-0032.

Received: 12 October 2009 /Accepted: 15 November 2009 / Published: 2 D ecember 2009

Abstract: Point-based methods undertaken by experienced human operators arevery effective tor traditional photogram m etric activities, but they are not appropriate in the autonomous environment of digital photogrammetry. To develop more reliable and accurate techniques, higher level objects w ith linear features accom m odating elem ents other than points are alternatively adopted for aerial triangulaton. Even though recent advanced algorithm s provide accurate and reliable linear feature extraction, the use of such features that can consistof com plex curve form s is m ore difficult than extracting a dia:rete set of points. Control points that are the initial input data, andbreak points that are end points of ^gmented curves, are readily obtained. Employment of high level features increases the feasibility of using geom etric inform aton and provides access to appropriate analytical solutions for advanced com puter technology.

K eyw ords: bundle block adjustm ent 3D natural cubic splines; arc-length param eterization; linear features; line photogram m etry

1. introduction

O ne of the m ajor tasks in digital photogram m etry is to determ ine the orientation param eters of aerial im ages quickly and accurately, which involves the two prim arysteps of -interior and exterior orientation. W hile the original aerial photography provides the interior orientation param eters, the problem rem ains to determ ine the exterior orientationw ihrespect to the object coordinate ^stern . Exterior orientation establishes the position of the cam era projection center in the ground coordinate ^stem and the three rotation angles of the cam era axis represent the transform aton betw een the im age

and the object coordinate system . Exterior orientation param eters (EO Ps) fora stereo m odel consisting of two aerialim ages can be obtained using relative and absolute orientation. This is a fundam entaltask in m any applications such as surface reconstruction, orthophoto generation, im age registration, and object. recognition.. The EO Ps of m ultiple overlapping aerial im ages can be com puted using a bundle block adjustm entt The position and oren^ta^tion of each exposure station are obtained by bundle block adjustm ents using collinearity equations that are linearized as having an unknow n position and orientation with the object space coordinate ^stem .

The program for bundle block adjustm ent in m ost sfttopy workstations em ploys point features as the control inform atton. Phottogram m etric trangulation losing (digital phottogram m etric workstations is m ore auttom ated tthanaerial triangulationusihg analog instrum ents because the stereo m odel can be directyset using analytical triangulaton outputs. Bundle bUockadjustm ent reduces the cost of field anveying in difficult areas and verifies the accuracy of field observations during the adjustm ent process. Even though each stereo m odel requires at least two horizontal and three vertical control points, this m ethod can reduce the num ber of control points w ith accurate orientation param eters. EO Ps of all the photographs in the target area are determ ined by the adjustm ent, which im proves the accuracy and reliabilityof photogram m etric tasks. B ecause object reconstructionis processedbyan intersection employing m ore than two im ages, bundle block adjustm ent provides the redundancyfor the intersection geom etry and contributes to the elim ination of the gross errorin the recovery of EO Ps.

A stereo model consisting of two images with 12 EOPs is a common orientation unit. The m echanism of object reconstruction from a stereo m odel is com parable w ith that of an anim al or hum an visual ^stem . The principle aspects of the hum an vision system , icluding its neurophysiology, anatom y, and visual perception, are well described in Schenk [1]. A point-based procedure relationship betw een point prim itives is widely developed in traditional phottogram m etry, such that one m easured point on an im age is identified in another image. Even for linear features, data for a stereo m odel in a sfttopy workstation is collected as points so that further application and analysis relies on points as prim ary input data units. The coefficients of interior, relative, and absolute orientation are com puted from the point relationship. ]terbrbrienLtatioh compensates for lens distortion, film ±irinkage, scanner error, andatmosphere refraction. Relative orientation makes the stereoscopic view possible, and the relationship between a model coordinate ^stem and an object space coordinate system is reconstructed byabsolute orientation. G roundcontrol points (GCPs) are w idelyem ployedto com pute orientation param eters. Although the use of m any GCPs is a tim e-consum ing procedure and inhibits the robust and accurate auttom ationthat research into digital phottogram m etryaim s to achieve, the deploym ent of a computer, storage capacity, phottogrammetric software, and a digital camera can reduce the com putational and tim e com plexity.

Employing high level features increases the feasibility of gaining geometric information and provides a aiitable analytical situation for advanced computer technology. W ith advancing developm en.t in the extraction, segm entation, classification, and recognition of features, the input data for feature-based phottogram m etry has been expanded at the expense of a redundancy in the application of aerial triangulaton. Because the identification, formulation, and application of reasonable linear features is a crucial procedure for autonom ous phottogram m etry, higher order geom etric feature-based m odeling plays an im portantrole in m odern digital phottogram m etry. The digitalim age form at is suited

to tthis purpose, especially in feature extraction and m easurem ent, and it is uœful for précisé and rigorous modeling of ieatures frm images.

2. Line Photogram m etry

21. O verview of Line Photogram m etry

Line photogram m etry refers to applications such as single photo resection, relative orientation, triangulation, im age m atthing, im age registration, and audace reœnstuiction, which are im plem ented using linear ffeatures and tthe correspondence between linear ffeatures rather tthan points. Interest conjugate points such as edge points, corner points, and points on parking lanes operate w ell ibr determ ining EO Ps w itthrespect tothe object space cœrdinatefram e intraditional photogram m etry. The most wel-known edge and interest point detectors arethe Canny [2], Fôrstner [3], Hairis, otherwise wel-known as the Plessy detector [4], M oraavec [5], Prew itt [6], Sobel [7], and SU SAN [8] detectors. The Canny, Prew itt andSobel operators are edge detectors and the Fôrstner, H arris, and SU SAN operators are corner detectors. Other wel-known corner detection algorithms are the Laplacian of G aussian, the différence of G aussians, and the determ inant of H essian. Interest point operators that detect w el-defined points, edges, and corners play an im portant role in autom ated triangulationandstereo m atching. For exam ple, the H arris operator is definedas a m easurem ent of corner strength as:

H (x,y) = detM )-a(traae(M ))2 (1)

where M isthe local structure matrix and a isa paraametera) that 0 <a< 0 25 . A defaultvalue is0.04. The gradients of the x and y directions are :

g = — , g = — (2)

with I as an image. The local structure matrix M is:

where A = gX,B = g^ and C = gxgy. A corner is detected when:

H (x,y)> Hthr (4)

where H thr is the thresholdparam eter onco^ier strength. The H arris operator searches points where variations in two orthogonal directions are large using the local autocorrelatbn function and provides good repeatability under varying rotation, scale, and illum ination. The Fôrstner corner detector is alto based on the covarance m atrix for the gradientat a targetpoint.

M arr [9] proposes the zero-crowing edge point detector'utilizing second order rather than firstorder directional derivatives. The m axim um of first order derivatives indicates the location of an edge whereas it is the zero of second orderderivatives that indicates an edge. Physical boundaries of objects

are easily detected because the gray levels cHange abruptHy in boundaries. B ecause no Single operator exists in edge det^ctio^, ^veial crieria are required ibr each specific application. M atching point features present large percEn^tag^ of m atch errors because point features are am biguous and an analytical Solution for point m atching is not yet developed. B ecause of the geom etric inform aton and sym bolic m eaning of linear features, m atching them is m ore relable tthan m atching point features in the autonom ous environm ent of digital phottogram m etry. As the use of linear features does not require a point-to-point Gorrespondence, the m atching of linear features is m ore flexible than tthat forpoints.

A num ber of reseaüchers have published studies on autom atc feature extraction and iis applicatbn for varixis phottogramm etric tasks: Former [10], Hannah [11], Schenk et al. [12], Schenk and Toth [13], Ebner and Ohlhof [14], Ackernan and Tsingas [15], Haala and Vo^eGman [16], Drewniok and Rohr [17,18], Zahran [19], and Schenk [20]. However, point-based phottogram m etry based on m anual m easurem ent and the identifiGatbnof interest points is not com patible w iththe autonom ous environm ent of digital phottogram m etry, but is a labor-intensive interpretation w ith the lim itations of occlusion, am bigu^ityr, and sem antic inform ation when com pared with appropria^te robbst auttom ation.

B ecause point features do not provide explicit geom etric inform ation, geom etrical know ledge is achieved by perceptual organization [21-24]. Perceptual brgah:кätihderives features andstructures from im agery w itthout a prior know ledge of the geom etric, ^pectral or radiom etric properties of features and is a required step inobjsct reoognition. Perceptual organization is an interm edite level process for various vision taiäks aich as target-to-background discrim ination, obj^recog^ittio^, target cueing, motion-based grouping, auface recbhst:ruciioh, im age interpretation, and chahge det^ctio^. B ecause objecto camot be distinguished by one gray level pixel an im age m ust be investigated entiely to obta:! perc^tual inform ation. The m ost recent r^^arch related to perceptual organization cbncems the 2D im age im plem entation atsgnal prim itive, and structural levels.

In general grouping or ^gm entation has the sam e m eaning as perceptual organization in com puter vision. This segm entation is typicaly aiddre^ed by two approaches, a m odel-based m etthod (top-dow n approach^) and a data-based m ethod Obottom -up approach^), and m any reseaüchers have em ployed edges and regions in segm entation. In the edge-based approach^, edges are likened to general form s of linear features w ithout diarontinuities, and in region-based approach^, iterative region grow ing techniques using seed points are preferred for auface fitting. M odel-based m etthods require dom ain know ledge for each spec^iSc äpplication in a m anner sm ilar to the hum an visual ^stem , whereas data-based m etthods em ploy data properties for data recognition in a global faähLion. In data-based m etho^s, the sam e invariant properties :hdiffprentpbsitiohs andbr:ier^ta^tio^s are cbmbinedi^tb the ^m e r^io^s or the ^m e features. One approach albhe, how ever, canhbt guaran^t^ cb^si£^^t quality ^ com bined approach^ are im plem ented to m inim ize ^rror ^gm ^ta^tio^.

Sym bblic r^r^eh^ta^tioh u^sing d^ii^ct points is ^ifEi::ullt becau^se points co^tain no explc^it

in^form ationäOout ph^sical i^ality^. W hile traditional phottogram m ^tric tech^iques bbtainthe cam era param ^t^ from the correspondence betw ^^ 2D and 3D points, a m ore general and relable proce^ is req^ir^ for ädvähced computer techhblogy such as the adoption of lin^jr f^tur^. Line phottogram m etry is aiperior in higher level taek^s ^ch as object recog^iiioh and au^ttom ation as com par^ with point-based phottogram m ^try^, b^t selectoii of the corr^ cän^idäte linear features is a com plicated proce^. The developm e^t of the allgbrithm from poi^t- to lne-based phottogram m ^try uses the advantages of both approaches. The ^lection of ^itable features is easier than the f^traG:iioh

of straight linear features andthe candidate feature can be used inhigher applications. A reason for developing curve features is that they will be prior to, and a fundamental aspect of, the next highest features arch as airfaces, areas, and 3D volum es that consist of free-form linear features. Line-based photogram m etry is m ost suitable in the developm ent of robust and accurate techniques forauttom ation. If linear features are em ployed as control features, they provide advantages over points in the autom ation of aerial trangulation. Photogram m etry based on the manual measurement and ientification of conjugate points is less reliable thanline-basedphotogram m etrybecause it has the lim itations of occlusion (visibility), am biguity (repetitive patterns), and sem antic inform ation w hen considering the need for reliable and effective autom aton. The m anual ientifration of corresponding entities w ithin two im ages is crucial in the autom ation of point based photogram m etric tasks. No know ledge of the point-to-point correspondence is required in line-based photogram m etry. In addition, point features do not caryinform ationabout the scene w hereas linear features cbntaintne sem antic information related to real object features. Additional information concerning linear features can increase the redundancy of the point ^stem .

2 2. Literature Review

A review of related works begins with those using m ethods of pose estim ation in im agery based on linear features that appear in m ost m an-m ade objects such as buildings, roads, and parking lots. Over the years, a num ber of researchers in photogram m etry and com puter vision have used line instead of point features; for example, M asry [25], Heikkila [26], Kubik [27], Petsaand Patias [28], Gulch [29], W iman and Axelsson [30], Chen and Shibasaki [31], Habib [32], Heuvel [33], Tommaselli [34], Vosselm an and Veldhuis [35], Forstner [36], Sm ith and Park [37], Schenk [38], Tangelder et al. [39], and Parian and Gruen [40]. Mulawa and M ikhail [41] originally proved the feasibility of linear features for close-range photogram m etric applications such as space interaction and resection, and relative and absolute orientation. This was the first stEp in employing linear feature-based methods in close-range photogram m etric applications. M ulaw a [42] later developed linear feature-based m ethods for different sensors.

W hereas straight linear features and conic sections can be represented as unique m athem atcal expressions, free-form lines in nature cannot be described by algebraic equations. Hence, M ikhail and W eerawong [43] used splines and polylines to represent ftee-form lines as analytical expressions. Tom m aselli and Tozzi [44] proposed that the accuracy of the straight line param eterbe a subpixelwith the representation of four degrees of freedom in an infinite line. M any researchers in photogram m etry have described straight lines as infinite lines using m inim al representation to reduce unknow n param eters. The m ain consideration in straight line expression is in the singularities. Habib et al. [45] m ade a bundle block adjustm ent using a 3D point set lying on control linear features instead of traditional control points. EO Ps w ere reconstructed hierarchicallyem ploying autom atc single photo resection (SPR).

Habib et al. [46] aimmarized linear features extracted from a mobile mapping ^stem , a GIS database, and maps for various photogrammetric applications such as SPR, trangulation, digital cam era calibration, im age m atching, 3D reconstruction, im age to im age registration, andsurface to aufface registration. In theirwork, m atched linear feature prim itives were utilized in space intersection

for the reœnsttuctonof object space features, andlinear features inthe object space w ere used as control features in triangulation and digital cam era calibration.

M ikhail [47] and Habib et al. [48] accomplished the geometrical modeling and the perspective transform atton of linear features within a triangulation process. Linear features were used to recover relative orientation parameters. Habib et al. proposed a fteee-form line in object space by a sequence of 3D points along the object space line.

Schenk [49] extended the concept of aerial triangulation from point features to linear features. The line equation of six dependentparam eters replaced the point-based colinearity equation :

X - X A + t- a Y - Ya + t- b Z = ZA + t- c

where a real variable ist,the startpoint (XA,YA,ZA), and the direction vector (abc). The traditional point-based œlinearity equation was extended to line features:

(X A + t- a - XJr + Ya + t- b - YC)Î2 + (Za + t- c - Zc)r

Xp = -f-

yp = -f-

(Xa + t- a - X c )r3i + Ya + t- b - Yjr^ + (Za + t- c - Zc (Xa + t- a - X c )r2i + Ya + t- b - Yç)Ç2 + (Za + t- c ~ Zc )^r3 (Xa + t- a - X c )r3i + Ya + t- b - Yjr^ + (Za + t- c - Zc )r^

with x yp as photo coordinates, f the focal length, Xc,Yc,Zc the camera perspective center, and rj the elem entsof the 3D orthogonal rotation matrix. The extended colinearity equation with six param eters was derived as the line expression of four param eters tyfixy0) because a 3D straight line has only four independent param eters. Tw o constraints are required to solve a com m on form of the 3D straight equations using six param eters determ ined by two vectors :

cos#cos^ • x - sii^ - y0 + sin #cos^ • z

cos^sii^- x + cos^- y0 + sin0 sinz 7)

- sin 0 ■ x + cos# • z

where z is a real variable. The advantage of the 3D straight line using four independent param eters is that it reduces the com putation and tim e com plexity in processes such as bundle block adjustm entt The collinearity equation as the straight line function of fourparam eters was developed :

(X - Xc )rii + Y - Yc)ç2 + (Z - Zc)Î3

xp =" f yp = -f

(X - XJri + Y - YJr + (Z - Zc )r^3 (X - Xc )rTi + Y ~ Yc )rr2 + (Z ~ Zc)3T3 (X - Xc)ri + Y - YJr + (Z - Zc)3T3

where X, Y, and Z were defined in equation (7).

Zalm anson [50] updated EOPs using the correspondence betw een the param ettric control free-form line in object space and the projected 2D freee-form line in im age space. The hierarchical approach, the m odified iterativelyclose point CcP)m ethod, was devetpedto estim ate curve param eters. The ray lies on the freee-form line w hose param ettric equation represented by one param eter follow s. Besl and M cKay [5i] employed the ccp algorithm to solve a matching problem of point sets, free-form curves,

auffBœs, and terrain models in 2D and 3D space. In thBirworrk, an ICP algorithm was executed without prior knowledge of the correspondence between points. The ICP method affected the Zalm anson' s dis^rtatdon on the developn ent of the recovery of EO Ps using 3D free-form lines in photogram m etry. Euclidean 3D transform aton was then em ployed in в search for the closest entity in the geom etric data set RBbbani et al. 52] utilized the IC P m ethod in the registration of Lidar point clouds to divide then into four categories (spheres, planes, cylinders, and tori) by direct and indirectm etthods.

S l) =

" X l) " Xk " P1

Yl) - Y0k + P2

_ z l)_ zk P3 _

where X0lY0lZ0la>l(plK are the EOPs and p1,p2,p3 are the direction vector-.

The paranetric curve T(t = IX (t) Y (t) Z (t)]T was obtainedbym inin izing the Euclidian distance between two paran etric curves:

Ф (l) = ||Г t) - S l))2 - (X (t) - X0 - Al)2 +Y(t)-Y0 - p2l)2 + (Z (t) - Z0 - p3l)2 (10)

Ф (tl had в mininum value at5Ф /<5l= 5Ф /<5t= 0 with two indepBndentvBriablBs l and t as in 11). 0Ф /öl- -2p1 (X (t) - X0 - pl - 2p2 Y (t) - Y0 -p2l)~ 2p3 (Z (t) - Z0 - p3l) = 0

0Ф /at= 2X (t) (X (t)- X0 -p1l)+ 2Y (t) Y (t)-Y0-p2l)+ 2Z (t) (Z (t)- Z0~p3D = 0

Akav etal. [53] em ployedplanar free-form curves for aerial trangulationw iththe ICP method. B ecause the effBat of the Z paran eter as con pared w itth that of X or Y w as large in в norm al plane equation aX + bY + cZ = 1, в different plane representation was developed to avoid nun erial problem s:

11 sii в cos^ n. - sii в sii ç n3 aos^

n. (X - X0) + n, Y -Y0) + n3(Z - Z0)= 0 (12)

^X + nY + n3Z - D

with в the angle fron the XY plane, (p the angle around the Z axis, n the unit vector of plane norm al and D the distance between the plane and the origin. Five relative orBntationpBrBm eters andtthrBB planar paran eters w егв obtainedbyusing the hon ographym apping ^stem , w hich marched for the con.jugate point in an in age corresponding to в point in the other in age.

Lin [54] proposed the m ethod of autonom ous recovery of exterior orBntation paran eters by an extension of the traditional point-based modified iterated Hough transform (M IHT) to the 3D free-form linear íBaturB-basedM IHT. Straight polylines w ere generalizó for m atthing prim ittves in the pose estdn atom because the m athem atdcal representation of straight lines is m uch clearer than the algebraic expre^ion of conic sections and splines.

Gruen and Akca [55] matched 3D curves whose forms were defined by a cubic spline using m atching least squares. Subpixels were localized by the m atching, and the quality of the localization was decided by the geom etry of im age patches. Tw o free-form lines were defined as:

f (u) = [x(u) y(u) z(u)]T = a0 + a1u + a2u2 + aju3 g(u) = [x'li) y'(u) z'(u)f -b0 + bu + b2u2 + tyu3

where ue 0>1],a0,a1 ,a2,a3,b0,ki,b2,bj are variables and the f(u),g(u)e^3 . The Taylorexpansion was em ployed to adopt the Gauss-M arkov m odel:

f (a)-e(u) - g (u)

f(u)-e(u) - g0 (u) + 5g u) du (14)

— < i . . 0 , N ög0 (u) du ^ dg0 (u) du dg0 (u) du ^

f (u)- e(u) - g0 (u) + —--dx +—--dy+ —--dz

öu dx öu öy öu dz

3. Bundle Block Adjustment with 3D Natural Cubic Splines

31. 3D Natural Cubic Splines

The choice of the right feature m odel is im portant in the developm ent of a feature-based approach because an am biguous feature rEpresentation leads to unstable adjustm entt A spline is a piecew ise polynomialfunction in the n of vector graphics. Splines are widely used for data fitting in computer aoience because of the resultant smplicityin curve reconstruction. Comply figures are well approximated through curve fitting and a spline lends strength to the accuracy evaluation, data interpolation, and curve sm oothing. One of the im porant properties of a goline is that it can easily be m orphed. A goline represents a 2D or 3D continuous line within a sequence of pixels and ^gm entation. The relationship betw een pixels and lines is applied to a bundle block adjustm ent or a functional

representation. A spline of degree 0 is the sim plest goline, a linear goline has degree 1, a quadratic

■ ■ ■ ■ ■ ■ 2

golrne has degree 2, and a com m on natural cubic spline has degree 3 w ith continuity C . The

geom etrical m eaning of continuity C is that the first and second derivatives are proportonal at joint pouts and the param etrc im portance of continuity C is that the first and second derivatives are equal atconnected points.

The num ber of break points that are the determ ination of a ^t of piecew cubic functions varies depending uponthe spline param eters. A natural cubic degree guarantees a second-order continuity, which m eans that the first andsecondorder derivatives of two consecutive natural cubic splines are continuous at the break point. The intervals for a natural cubic spline do not need to be the sam e as the distance of every two consecutive data points. The best intervals are chosen by a least squares m etthod. In general, the total num berof break points is less than that of original inputpoints. The algorithm of a natural cubic spline is as below .

Generate the break point (control point) set for the spline of the original input data. Calculate the maximum distance betw een the approximated spline and the original input data while the maximum distance > the threshold of the maximum distance. Add the break point to the break point setat the location of the maximum distance. Compare the maximum distance with the threshold.

A larger thresholdm akes for m orebreak points w itha m ore aœurate spline to the original input datai. N piecew ise cubic polynom ial fonctions betw een two adjacent break points aie defined from the N + 1 break points. There isa sepatate cubic polynom ial for each segmentwih itsown rre^fficipntp!

Figure 1. Natal cubi: splines.

X o (t) - aoo + aoit+ ao2tf + ao3t3 ,te [01]

XI (t) - aio + a11t+ a12t2 + al3t3 ,te [oi]

X i (t) = aio + ai1t+ a12t2 + a13t3 ,te [o1] Y, (t) - bo + bi1t+ bi2t2 + bi3t3 ,te [o1] Zi (t) - cio + ci1t+ ci2tf + ci3t? ,te [o1]

The strngth of tth^is appioach is that segm ented lines represent a fre-form line w ith anaLytical param eters. The num ber of break points is ieduœd and the input error should be absorbed by a m athem atcal m odel, espedally in the expression of points on a straight line. A natuualcubic spline is a data-independent curve fitting. The disadvantage is that the w hole curve shape depends on all of the passing poir^ts>, and changing any one of tthese alters the entire curve.

The cor]tppsponden<cp betw een the 3D curve in the object space aoordinate system and its projected 2D curve in the im age coordinate system is im plem ented using an accom m odating natuual cubic spline curve featur because of itsboundary conditions thatretain zero second derivatives attthe end points. A natural cubic spline is corn posed of a sequence of cubi: polynom ial segm ents as in Figure 1 w ith xo ix1,^,xn as the n + 1 controlpoints and Xo X1, — ,Xn-1 as the ground œordinates of n segm ents.

3 2. Extended CoHinearity E quation M odel for Splines

Collinearty equations are the com m only used condition equations to determ ine relative orientation. The space intersection calculates a point lxation in object. space using the projection ray intersection from two or m ore im ages, and the space resection determ ines the coordinates of a point on an im age and the EOPs with respect to the object space coordinate system . The space intersection and the space

resection are the fundam ental operations in photogram m etry for further prooes^s such as triangulaton. The basic concept of the collinearity equation is that all points on the im age, the perspective center, and the corresponding point in the object space are all on a straight line. The relationship between the im age andobject coordinate ^stem s is expressedbytthree positionparam eters andtthree orientation parameters. Collinearity equations play an important role in photogrammetry because each control point in object space produces two collinearity equations for every photograph in which the point appears. If m points appear in n im ages, then 2m n collinearity equations can be em ployed in the bundle blockadjustm entt The extendedcoliearityequations relating a natural cubic spline inobject space with ground coordinates Xi(ttYi(tj,Zi(tj) with image space having photo coordinates (xpi,ypi) are ^en as (16). A natural cubic spline allow s the utilization of the collinearity m odel forexpre^ing orientation param eters and curve param eters as below :

(Xi (t) - Xc )ri + Yi (t) - Yc )r2 + (Zi (t) - Zc )ri3

Xpi =" f"

ypi = - f

(Xi (t) - Xc )rri + Yi (t) - Yc )rr2 + (Zi (t) - Zc )rr3 (X i (t) - Xc )r2i + Yi (t) - Yc )rr2 + (Zi (t) - Zc )rr3 (Xi (t) - Xc )ri + Yi (t) - Yc )r32 + (Zi (t) - Zc

with (xpi,ypi) as the photo coordinate, f the focallength, XC ,YC ,ZC the cam era perspective center, and rij the elementsof the 3D orthogonal rotation matrix RT by the angularelements (a>,(p,K) of EOPs.

Figure 2. The projection of a pointon a spline.

The extended collinearity equations can be written as follow s:

x = -fw' yp = -f77

= RT ko,(p,K

Xi(t)- Xc Yi(t)- Yc

Zi(t)" Zc

To recover the 3D natural cubic spline param eters and the exterior orientation param eters in a bundle block adjustm ent, a nonlinear m athem atical m odel of the extended collinearity equation is differentiated The models of exterior orientation recovery are classified into linear and nonlinear

m ethods. W hereas linear m etthods decrease the com putation load, the accuracy and reliability of linear algorithms are not excellent-. Nonlinear methods are more accurate and predictable. However, nonlinear methods requite initial estimates and they increase the computational complexity. The relationship betw een a point in im age goace and a corresponding point in object space is establined by the extended collinearity equation. Prior know ledge of the correspondences betw een individual points in the 3D object space and their projected features in the 2D im age space is not required in extended collinearity equations with 3D natural golines. One point on a cubic spline has 19 parameters Xc,Yc,Zc,®,^,A:aoaiA2A3bobi,b2,b3,cociAA,tt. The differentials of (17) are derived by 18):

f fu f

dxp =--du H—j dw, dyp =--dvH—j dw

with differentials of du, dv, and dw (19).

du dv dw

X±(t)- Xc

Yi(t)- Yc

Z±(t)- Zc

X±t)- Xc

Yi(t)" Yc

Z±t)- Zc

X±t)- Xc

Yi(t)-Yc

Z±t)- Zc

6K- RT

"0" "0" "1" " t t2 t3 "0

RT 1 dYc - RT 0 dZc + RT 0 da0 + RT 0 da, + RT 0 da2 + RT 0 da3 + RT 1

0 1 0 0 0 0 0

"0" "0 " " 0" "0" "0" 0" 0"

RT t dbL - RT t2 db2 + RT t3 db, + RT 0 dc0 + RT 0 dq + RT 0 dq> + RT 0

0 0 0 1 t t2 t3

a + 2a2t+ 3a3t2 b + 2b2t+ 3b3t2 q + 2qt+ 3qt2

Substituting du, dv, and dw in (18) by the expressions found in (19) leads to:

dxp =M 1dXc + M 2dYc + M 3dZc + M 4d« + M 5d<p + M 6d^ + M 7da0 + M 8da + M 9da2 + M ^da

+ M 11db0 + M 12db1 + M 13db, + M 14db3 + M 15dco + M 16dq + M 17dc + M 18dc3 + M 19dt dyp = N1dXc + N2dYc + N3dZc + N4d« + N5d^ + N6d^ + N7da0 + N8da1 + N9da2 + N10da3 + N11db0 + N12db1 + N13db, + N14db3 + N15dc0 + N16dq + N17dc, + N18dc3 + N19dt

M 1v. .M 19jN1,—N19 denote the partial derivatives of the extended collinearity equation forcurves. The linearize extended collinearity equations by Taylor expansion, ignoring the 2nd and higher order term s, can be written as follow s:

x + f —- =M idXC + M 2dYC + M 3dZC + M 4d® + M 5d<p + M 6d^ + M 7da0 + M 8da + M 9da w

+ M i0da3 + M iidb0 + M i2dbi + M i3db, + M i4db3 + M i5dc0 + M i6dq + M i7dc^ + M i8dq + M i9dt+ ex

yp + f—0 = NidXC + N2dYC + N3dZC + N4d® + N5d^ + N6d^ + N7da0 + N8da + N9da2 w

+ Ni0da3 + Niidb0 + Ni2dbi + Ni3db, + N^db, + Ni5dc0 + Ni6dq + Ni7dcL, + N^dq + Ni9dt+ ey

withu0vr°w° being the approximate param eteisby (x0 ,yc0 ,zc0 V y ,a00 a ,a0° ,a30 ,b00 b ,b20 b ,c00 ,q0 ,c0 ,c0,t0), and exey. the stochastic errors of xp,yp, the observed photo coordinates with zero expectation, respectively. Orientation parameters including the 3D natural cubic spline parametersare expected. to recover correctly because the extended collinearity equations with these golines increase redundancy.

3 3. Arc-Length Param eterization of 3D Natural Cubic Splines

The assum ption m ade in bundle block adjustm ent by the G auss-M arkov m odel is that all the estim attd param eters are uincorrelated. Hence, the design m atrix of the adjustm ent m ust be full rank, nonsingular, and norm al. However, because the goline param enters are not independent of their location param eters, additional observations are required to obtain param eter estim ations. In a point-based approach, the point locatonrelationship betw een im age andobject space is established for ihe pose estim ation to include the fundam ental cam era position and orientation, the rem ote sensing, and the com putervision. The coordinates of a point are necessary for the space interaction and resection. To rem ove any rank deficiency caused by datum defects in point-based photogram m etry, som. e constraints are adopted to estimate the unknown param eters. The most common constraints are coplanarity, ^mm etry, perpendicularity, and parallelism . The m inimum number of constraints is equal to the rank deficiency of the ^stem . Inner constraints are often used in a photogram m etric network:, which can be applied to both the object features and the cam era orientation param eters. A ngle or distance condition equations provide inform ation on the relativity betw een observations in object space and points in im age space. Absolute inform ation can be obtained from the fixed control points.

In this research, an arc-length param eterization is applied as an additional condition equationtto solve the rank deficientproblem in extended collineariy equations using 3D natural cubic splines. The concept of different-able param eterization is that the arc length of a curve can be divided into m inute pieces and these can be sum m ed such that each piece w illbe approxim ateHylinear. The sum of ihe squares of derivatives is the sam e with a velocity vectorbecause a param etric curve can be considered as a point trajectory. A velocity vector describes the path of a curve and m ovem ent characteristics. If a particle on a curve moves at a constant rate, the curve is param eterized by the arc length. W hile the extended collineariy equation provides the only information, curves have additional geometric constraints aich as arc length, tangent of location, and curvalue. These support space resection under the a^Lim ption of properly accounting for additional independent observations in both the im age and object space. Arc length in object space is determ ined by a geom etric integration using a construction from the differentiable param denization of a spline. A rc length in im age space is calculated by a

geometric integration of a construction from the differentiable parameterizai-Qn of the photo coordinates derived from a spline in the objsct space:

Arc (t) =

t)=fV (xpt))2 +(ypt)): f f

u'(t)w (t)- u(t)w '(t)

v'(t)w (t)- v(t)w'(t)

where f is the foc^l^gth and:

u(t) v(t)

w (t)_

u'(t) v'(t) w'(t)

= RT ko,(p,K

= RT ko,(p,K

a0 + aLt+ a2t + ajt3 - X C b0 + tit+ b2t + bjt3 - YC co + qt+ c2t2 + c3t3 - Zc

aL + 2a2t+ 3a3t2 b + 2b2t+ 3b3t2 q + 2qt+ 3qt

B ecause the problem of the arc-length param ^te^riza^tion of golines has no anaLytical solution, several num erical approx;im atons of reparam ^^riiza^tion tBch^;iques for splines or other curve r^r^sen^ta^tio^s have been developed. W hie most curves are not param^te^riz^ for arc length, the arc length of a B-pline can be reparam ^te^rized by adjusting its knots. W ang et al. [56] approx;im ated the param e^te^riz^ arc length of spline curves by generating a new curve that accurately approxin attd the original spline curve to reduce the com putaeorqom p]ex:;tyof the arc-lengthparam . They

show ed that the approx;im ation of the arc-length param ^t^riiza^tion works wBI in a variety of real-tim e applications including driving sim nations.

Guвnter and Parent [57] employed a hierarchical approach algorihm to develop a linear search arc-length subd;iv;ision table for param eterizeri curves to reduce the arc-length qomputationtime. A table of the correspondence betw een param eter t and the arc length can be established to accelerate the arc-length computation. After dividing the parameter range into intervals, the arc length of each interval -s com puted for m apping param eters to the arc length. A reference table for various intervals of arc length can be developed. Another method of arc-length approx^ination -s to use explicit functions such as the Bez;er curve, which has advantages in fast function ^/aljatio^s. Adaptive G aussian integrations em ploy a recursive m ethod that starts from a few sam p]es and adds on m ore as necessary. Adaptive G aussian integration also uses a table that m aps curves or spline param eter values according to the arc-length values.

Nasri et al. [58] proposed an arc-]ergthapprox;im ation m ethodof circles and p-bcbw ise circular splines generated by control polygons or points using a recursive subd;iv;ision algorihm . W hie B-plies have various tangents over the curve depending upon the arc-length param eterizatiori, circular splines have constant tangents w hose vectors are useful in arc-length com p^ta^tio^.

Simpson's rule is the num ericalapproxim attonof definite integrals. The geom etric integration of the arc length in the im age space can be calculated by this rule as follow s:

Arc (t) = jf t)dt t f t)dt=

f (t) + 4 f

+ f (tt:

f (t) ^(xp (t))2 + (yp (t))2

df= - (ft))2 f 2

, w' , 1 , w' , 1

2x (t)— du -2xp (t) — du^2yp (t) — dv- 2yp (t) — dv'

w2 w w2 w

I , ,%u'w2 - UW -uw')2w , ,%v'w2 - Vw - vw')2w I n + i 2xp (t)-4--2yp (t)-4-^dw

+ ^2xp (t)^ + 2yp, (t) [ w

tt0 -10

Arc (t) -

f (t°)+ 4 f

+ f (t0

=A1dXC + A2dYC + AjdZC + A4d« + Agd^ + A6dK + A7da0 + A8da1 + A9da2 + A10da3 + A11db0 + A^db + A13db, + A14db! + A15dq3 + AL6dcL + A^dc, + AL8dc3 + A^dt + A^dt; + ea

with t0,^, f (t0) being the approxinate parameterusing the following:

(x C ,YC0 ,z0 ,co° ,(p° ,k° ,a0 í^0 'a30 £>0 ,]q0 ,b0 'b0 /c0 ,q0 ,c0 ,cc0 ,tf) and w ith ea being the stochastic error of the arc length between two locations with zero expectation. A1,_A20 denote the partial derivatives of the arc-length param eterizatbn of a 3D natural cubic spline.

3.4. M odel Integration

The objective of bundle blockadjustm ent is twofold, nam elyto calculatetthe exterior orientation param etersof a block of images and also the coordinates of the ground features in objectspace. In the determ ination of orientation param eters, additional interior conditions such as lens distortion, atm ospheric refraction, and principal point offset can be obtained by self-calibration. In general, orentationparam eters are determ ined by bundle block adjustm ent losing a large num ber of control points. This establishm ent of control points, how ever, m eans expensive fieldwork, so an econom ical andaacurate adjustm ent m ethodis required. L inear features have several advantages to com plem ent points in that they are useful for higher level tasks and they are easily extracted in m an-m ade

environm ents. The line photogram m etric bundle adjustm ent in this research aim s at the estim ation of exterior orientation param eters and 3D natural cubic spline param eters using the correspondence between splines in object space and spline observations of multiple images in image space. Nonlinear functions of orientation param eters, spline param eters, and spline location param eters are represented by extended collineariy and arc-length param eterization equations. Five observation equations are produced by each two points, and these are four extended collineariy equations 21) and one arc-length param eterization equation (24). An integrated m odel provides not only for the recovery of the im age orientation param eters but also enables surface reconstruction using 3D curves. Of course, as the equation system of the integrated m odel has seven datum defects, control inform ation about the coordinate system is required to obtain param eters. This is a step tow ard higher level vision tasks such as object recognition and surface reconstruction. in the case of straight lines and conic sections, tangents are additional observations in the integrated m odel. C onic sections, like points, provide good m athem atical constraints because such sections provide nonsingular second degree equations. Such equations provide inform ation for reconstruction and transform ation and conic sections are divided by the eccentricity e. Because such sections can adopt more constraints than points and straight line features, they are useful for close range photogram m etric applications. in addition, conic sections have strengthincorrespondence establishment between 3D sections inobject space and theircounterpart features in 2D projected image space.

Ji etal. [59] em ployed conic sections for the recovery of EO Ps, and Heikkila used them forcam era calibration. A H ough transform ation reduces the tim e com plexity of conic section extraction using five param eterspaces fora SPR, cam era calibration, and trangulation.

Param eters are linearized in the previous sections and the Gauss-M arkov m odel is em ployed for the unknow n param eter esttim ation. The equation system of the integrated m odel is described as!:

ASP Ati

7EOP ^SP Й

m k1 m k

m km 1 1 i1

Nk1 Nkk

M 1 M 1

1 1 i7 1 1 i8

1 1 i7

1 1 is

1 1 i18

M 1 M 1

1 1 i19 1 1 i20

M m M m

1 1 i19 1 1 i2

Nm Nm 14 ilo J-^nO

1 1 i1Sf n

Ai1 Ak1

д km д km

Ai1 Ai2

Akm Ai20

= [dxk dYck dzk d^k drk f

C =[dai0 dai1 dai2 dai3 dbi0 dbi1 dbi2 dbi3 dci0 dci1 dci2 dci3 ]T

£ = [ct dt2 - dtn ]T

xP1 + ypi + Arct)pi - Arc t)

with Arc t)0 =

f0 t?)+ 4 f2

tL±t 2

+ f010:

, m as the num ber of im ages, n the num ber of

points on a spline segm ent, k the kth im age, and i the ith spline segm ent. B ecause the equation system of the integrated m odel has seven datum defects, the control inform atton for the coordinate system is required to obtainseventransform ationparam eters. Iha general phottogram m etric network, the rank deficiency referred to as datum defects is seven. Estm attes of the unknow n param eters are obtained by the least squares solution, which m inim izes the sum of squared deviations. A nonlinear least squares system is required in a conventional nonlinear phottogram m etric solution to obtain orientation param eters. M any observations in phottogram m etry are random variable that are considered as different values inthe case of rt^p^at^d observations such as the im age coordinates of points. Each m easured observation represents a random varable estim ate. If im age point coordinates are m easured using a digital phottogrammetric workstation, the values are measured slightly differently. The integrated and linearized Gauss-M arkov m odel and the least squares estim attd param eter vector with its dispersion matrix are:

7"= Aim + e

Ak a1 A1 ^op -^SP "t

=|£eop ^p 'tha = (Aim PAim ) am Pytl

D (|M )=CT02(AiaPAim ^

with e~ N (2,o"0;P_1) being the error vector with zero mean and cofactor m atrix P, a variance component^2, which can be known or not, is the least squares estimated param eter vector, and D )

is the disperson m atrix.

If one or m ore of the three estim attd param eter sets ^Eop are considered as stochastic

constraints, the reduction of the normal equation matrix can be applied. Control information is im plem ented as stochastic constraints in a bundle blockadjustm entt The distrfcutinandqualityof control features depend on the num ber and the density of control features, the num ber of tie features, and the degree of overlap of the tie features. If adding stochastic constraints rem oves the rank deficiency of the G auss-M arkovm odel, bundle adjustm ent can be im plem ented em ploying onlythe extended collinearity equations fortthe 3D natural cubic splines. Fixed exteriororentation param eters, control splines, or control spline location param eters can be stochastic constraints.

3 5. Evaluation of Bundle Block Adjustment

B undle block adjustm ent m ust be follow ed by an evaluation postadjustm ent analysis to check the suitability of project specifications and requiem ents. leratively rew eighed least squares and least m edian of squares are the appropriate im plem entation of a statistical evaluation that rem oves poor observations. The important elem ent affecting bundle block adjustm ent is the geom etry of aerial im ages. Generally, the previous flight plan is adopted to obtain suitable results. A simulation bundle block adjustm ent is im plem ented before em ploying a fightplan within the new projectdesign because such a sim ulation can reduce the effect of errorm easurem ents.

A qualitative evaluation that allow s the operator to recognize the adjustm ent characteristics is often used after bundle block adjustm ent. The sizes of the residuals in im ages are drawn for the evaluation. The im age residuals can be points or long lines and if all im age residuals have the sam e orientation, then the im age has a system atic error such as atm ospheric refraction or an orientation param eter error:. In addition, a lack of flatness in the focal plane m ay cause system atic errors in the im age space, which affects the accuracy of a bundle block adjustm ent. D istortions are different from one location to another in the entire im age space:. The topographic m easurem ent of the focalplane can correct the lack of focal plane flatness. Im age coordinate errors are correlated in the case of system atic im age errors. A poorm easurem ent can result in an indicated opposite residual direction or an exaggerated residual.

The three m ain elem ents in the statistical evaluation of bundle block adjustm ents are precision, accuracy, and reliability. Precision is calculated employing param eter variances and covariances, because a sm all variance indicates that the estim ated values have a sm all range and a large variance m eans that the estim ated values are not calculated properly. The range of the param eter variance is from zero, in the case of error free param eters, to infinity, in the case of completely unknown param eters. A dispersbnm atrixm aycontaindiagonal elem ents that are param eter variances. These and any off-diagonal elements are covariances between two param eters. Accuracy can be verified using check points that are not contained in bundle block adjustm ent like control points. Reliability can be confirmed from other redundant observations. The extended collinearity equations are a m athem atical m odel forbundle block adjustm ent: The m athem atical m odel consists of both functional and stochastic m odels. The functional one represents the geom etrcal properties and the stochastic one describes the statistical properties. Repeated m easurem ents at the sam e location in the im age space are represented with respect to the functional m odel and the redundant observations of im age locations in the im age space are expressed with respect to the stochastic m odel. W hile the Gauss-M arkov m odel uses indirect observations, condition equations such as coordinate transform ations and the coplanarty condition can be em ployed in the ad justm ent.

The G auss-M arkov m odel and the condition equation can be com bined into the G auss-Helm ert m odel. maddition, functional constraints such as points having the sam e height or straight railroad segm ents can be added into the block ad justm ent.

The difference betw een condition and constraint equations is that condition equations consist of observations and param eters, and constraint equations consist of only param eters. W ith the advance of technology, the photogram m etrical input data has increased so adequate form ulation of adjustm ent is required. All the variables are involvedinthe m athem atical equations and the weight matrixof the variables changes from zero to infinity depending upon the variances. Variables w ith near to zero

weight are considered as unknown parameters and variables with close to infinite weight are considered as constants. M ost actual observatibns exist between the two boundary cases. A ssessm ent by postadjustmene analysis is imports in phottogrammetry to evaluate the resuls. One of the assessm entm etods is to com pare the estim ated variance with the two -ailed confidence interval based on the norm al distrjbution. The two-tailed confidence interval is com puted by a reference variance a02

with %2 distribution as!:

where r is degrees of frEedom and a is a confidence coefficiEnt (or a confidence level). If has a

value outside of the interval, we can assume that the m atthem atcal m odel of adjustmentis incorrect through the wrong form ulation orlinear:izatibn, blunders!, or system atic errors.

3.6. Pose Estimation with an ICP Algorithm

In the previous spline segm ent case, the correspondence betw e:n spline segm ents in the im age and the object space was assum ed. In the present consideration, it is not know n which im age points belong to which spline segm eat The ICP algorithm can be utilized for the recovery of EO Ps because the initial estim ated param eters of the relative pose can be obtained from the orientation data for general phottogram m etric tasks. The original K P algorithm stEps are as follow s. The closEstpoint operators search the associate point using the nearest neighboring algorithm and then the transform ation param eters are estimated using a mean square cbsefUnction. The point is transformed by the estimated param eters andtthis step is iterativelyestablished tow ards convergence into a local m inim um of the m ean square distance. The transform ation, which includes translation and rotation betw.en two clouds of points, is estim ated iteratively towards convergence into a global m inimum . In other words, the iterative calculation of the m ean square errors is term inated when a local m inim um falls below a predefined threshold. A sm all global m inim um or a fluctuated curve requires m ore m em ory-intensive and tim e-consum ing com putation. In every iteration step, a local m inim um is calculated with varying transform ation param eters, but convergence into a global m inim um w ith the correct transform ation param eters is not always the result.

By the definition of a natural cubic spline, each param etric equation of a spline segm ent (Si t)) can be expressed as:

with Xi t),Yi t),Zi t) as the object space coordinates and a^b^q as the coefficients of the ith spline segm ent.

The ray from the perspective center XC ,YC ,ZC) to the im age point (xp yp ,-f) is:

" X i (t)~| aio + a±Lt+ a^e + a^t3 S, (t) - Y, (t) = b0 + bi1t+ bi2t2 + bi3t3 ,te [21] _ Z± (t)J cio + ci1t+ c12t2 + c^t3

" X l) " Xk "

Yl) - Yck + d2

_ z l)_ Zk

where:

d. d2 d

- RT (fOk ,«pk ,Kk

yp - f

with XCk,YCk,Z,k,iyk,^k,^k EOPs atthe kth itEration.

A point on the ray searches the closest to a natural cubic spline by m inim izing the follow ing target function forevery spline segm ent. Transform ation param eters relate to an im age point and its closest spline segm entcan be established using the leastsquares m ethod:

O l,t) = ||S l) - Si (t))2 - stationary l, t

The global minimum of ® l,t can be calculated by V 0QL,t) = 2 or50 /5l= 50 /5t= 2 . Substituting (28) and (29) into (31) and taking the derivatives with respectta l and t leads to:

1 ^ = (xc + ^l- aio - a^t- a12t2 - a^t3+ (yc + d2l- bio - bi1t- b12t2 - b^t3^

+ (zC + d31" Gio - cut- c^ - ci3t3 )d3 - 0

1 IT = (Xc + dxl- ai0 - ai1t- a^e - a^t3)(- afi - 2a12t- 3ai3t2) (32)

+ (Yc + ^l- bio - but- b12t - b^t3 X- b - 2b12t— 3bi3t2)

+ (Zc + ^l- cio - ci1t- c^t - c^t3X" Oi! - 2c12t- 3^) = 0

C onvergence into a global m inim um does not exist because equation (32) is not a linear system in l and t The relationship betw een an im age space point and its corresponding spline segm ent cannot be established with the minim ization method.

4. Experim entsand Results

This section dem onstrates the feasibility and the perform ance of the proposed m odel for the acquisitonof spline param eters, spline locationparam eters, andim age orientatbnparam eters based on control and tie splines in the object space within the simulated and real data sets. In general photogram m etric tasks, the correspondence betw een im age edge features m ust be established either autom atcaly or m anually, but in this study correspondence betw een im age edge features is not required. In a series of sax experim ents with the synthetic data set, the first test recovers spline param eters and spline locationparam eters in an error free EO Ps case. The secondtest recovers the partial spline param eters related to the spline shape. The third procedure estim atts the spline location param eters w ith error free EO Ps. The fourth step calculates EO Ps and spline location param eters,

ioIOw ed bytthe fifth step that estm ates EO Ps w ith full controlled splines in which the param etric curves used as control features are assum ed to be errorftee. in the Istexperm ent, EO Ps and tie spline param eters are obtained using the control spline.

Object space knowledge concerning splines, their relationships, and the orientation informatton of im ages can be considered as control inform atton. Spline param eters in a partial control spline or orientation param eters can be considered as stochastic constraints in the integrated adjustm ent m odel. The starting point of a spline is considered tobe a know nparam eter intthe partial control spline in which a0b>0, and c of the X, Y, and Z coordinates of a spline are known. The number of unknowns is displayed inTable 1 and Figure 3, where n is the num ber of points intthe object space, t show s the num ber of spline location param eters, and m represents the num ber of overlapped im ages in the target area.

Table 1. Num berof unknow ns.

EO P Spline N um ber of unknow ns

Known EO P Tie spline Partial control spline Full control spline 12 (n-1) + t 9n-1) + t t

Unknown Partial control spline 6m + 9 n-1) + t

EO P Full control spline 6m + t

Figure 3. D ifferent examples. (a) KnownEOPs withtie splines, b>) Known EOPs with partial control splines, (c) Known EO Ps with full control splines, d) Unknown EOPs with partial control splines, and (e) UnknownEO Ps w ithfull control splines. Red: Unknown param eters,Green: Partially fixed param eters,Blue: Fixed param eters).

(c) (d) (e)

Four points on a spline segm ent in one im age are the only independent observations so additional points on the sam e segm ent do not provide nonredundant inform aton to reduce the overall deficiency of the EO P and spline param eter recovery. To verify the inform aton content of an im age spline, we dem onstrate that any five points on a spline segm ent generate a dependent set of extended colinearity equations. Any combination of four points yielding eight colinearity equations are independent observations, but five points bearing 10 colinearity equations produce a dependent set of observations

related to the correspondence between a natural cubic spline in the im age and the object space. M ore than four point observations on an im age spline segm ent increase the redundancy related to the accuracy but do not decrease the overall rankdeficiency of the proposed adjustm ent system . Inthe sam e fashion, the case using a polynom ial of degree 2 can be im plem ented. Three points on a quadratic polynom ial curve in one im age are the only independent sets, so additional points on the sam e curve segm ent are a dependent observation. M ore than the independent point observations on a polynom ial increase the redundancy related to the accuracy, but they do notprovide nonredundant inform ation.

The amount of inform ation carried by a natural cubic spline can be calculated with the redundancy budget. Every spline segm ent has 12 param eters and every point m easured on a spline segm ent contributes one additional param eter. Let n be the num berof points m easured on one spline segm ent in the im age space and m be the num ber of im ages that contain a tie spline. 2nm collinearity equations and m (n - 1), the arc-length parametpriai-inns, are equations and 12 (the number of one spline segm ent param eters) + nm (the num ber of spline location param eters) are unknow ns. The redundancy is 2nm - m - 12 for one spline segment, so thatiftwo images (m = 2) are used for bundle block adjustm ent, the redundancy is 4n - 14. Fourpoints are required to determ ine spline and spline location param eters, inw hich case one spline segm ent andone degree of freedom to the overall redundancy budget is solved by each point m easurem ent w ith the extended œllinearity equation. Arc-length param eterizatbnalso contributes one degree of freedom to the overall redundancybudget. The fifth point does not provide additional inform ation to reduce the overall deficiency but only strengthens the spline param eters. This m eans it increases the overall precision of the estim ated param eters.

This fact show s the advantage of adopting splines inw hich the num ber of degrees of freedom is four because in straight tie lines only two points per line are independent. Independent inform ation, the num ber of degrees of freedom of a straight line, is two from two points or a point w ithits tangent direction. A redundancy isr= 2m - 4 with a line expression of fourparam eters because there are 2 nm œliearityequiations and the unknowns are 4 + nm [49]. Onlytwo points (n= 2) are available to determine four line param eters with two images (m = 2) so a^tl^stthr^ images mustoontaun a tie line. The information contentof t tie lines on m images ist 2m - 4). One straight line adds two degrees of freedom to the redundancy budget and at least three lines are required in the space resection. An additional point on a straight line does notprovide additional inform ation to reduce the rank deficiency of the recovery of EO Ps but only contributes im age line cœfficients. If spline location param eters or spline param eters enter the integrated adjustm ent m odel through stochastic constraints, em ploying extended œllneariy equations i enough to solve the system without the arc-length param eterization.

The redundancy budget of a tie point isr=2m-3sotie points provide one m ore independent equation than the tie lines. However, using tie points requires a ssm iautom atic m atching procedure to identify the tie points on all the im ages, and using linear features provides a m ore reliable and accurate basis for object recognition, pose determ ination, and other higher phottogram m etric activities than using point features.

41. Synthetic Data D escrp>tion

To evaluate the new bundle block adjustm ent m odel using natural cubic splines, an analysis of the sensitivity and robustness of the model is required. The model suitability can be verified by using the

estim ated param eters with a dispersion m atrix that includes standard deviations and correlations. The accuracy of bundle block adjustm ent is determ ined by the geom etry of a com plete block of im ages and the quality of the position and attitude inform aton of a cam era. A novel approach is a sim ulation of the bundle block adjustm ent. This is required prior to an actual experim ent w ith real data in order to evaluate the perform ance of the proposed algorithm s. Such a sim ulation can control the m easunem ent errors to m inim ize random noise affecting the overall geom etry of a block. Individual observations are generated based on the general situation of bundle block ad justm ent in order to estim ate the properties of the proposed algorithm s. A sim ulation allow s adjustm ent forgeom etric problem s orconditions with various experiments. A gpline isderived via three ground controlpoints (3232, 4261, 18), (3335, 4343, 52), and (3373, 4387, 34). Several factors that affect the estim atts of exterororientation param eters, spline param eters, and spline location param eters are discerned using the proposed bundle block adjustm ent m odel togetherwith both the sm ulated im age and the real im age blocks.

2500 3000 3500 4000 4500

Table 2. EO Ps of six bundle block images for simulation.

Param eter XC hi YC hi zc mi co degi (p deg] k Icieg]

In a je 1 3000 00 4002 00 503 00 01146 0 0573 5.7296

In a je 2 3305 00 4005 00 499 00 01432 0 0859 -5.7296

Im a je 3 3610 00 3995 00 505 00 01719 0.4584 2 8648

Im a je 4 3613 00 4613 00 507 00 0 2865 -0 0573 185 6383

Im a je 5 3303 00 4617 00 493 00 -01432 0 4011 173 0333

Im a je 6 2997 00 4610 00 509 00 -01833 -0 2865 1816276

F jgure 4. S nx m age block.

■ Image 1 Image 2 Image 3 Image 4 Image 5 Image 6 Spline

Figure 5. Natural cubic spline.

0 < t < 1 tor two segments 0 < t < 1 for two segments

0 < t < 1 for two segments1

4 2. Experiments with Error Free EOPs

Spline param eters end spline location param eters are dependent upon various controls, and the unknow ns can be obtained by a com bined m odel of extended collinearity equations and the arc-length param eterization equations of splines. Splines in the object space are consider^ as tie lines in the sam e fashion as tie points in a conventional bundle block adjustm entt D ata on the exterior orientation param eters is consider^ as control inform ationin this experim entt A w el-know n fact inem ploying the last ^juares system is that good initial estim ates of true values m ake the ^st^ sw ifty convergent tow ards the correct axLtion.

Norm ally distrbut^l random noite is added to points in the im age space coordinate syst^ in all the experim ents. This has a zero m ean and a = ±5^m standarddeviaton. G enerally, the larger the noite level the m ore accurate are the approxim atons r^ui^ to achieve the ultim ate convergence of the results. A worstcase scenario forestimation is that the large noise level causes the proposed model not to converge tow ards the specific estim ates because the convergence radius is then proportional to the no;te level. The param eter estim atonis ^nsitive to the noite of the im age m easur^ ent Error propagation related to the noise in im age space observation is one of the m ost im portant elem ents in the estim atton theory. The proposed bundle block adjustm ent can be evaluate statistically using the variances and the covariances of param eters because a sm all variance indicates that the estim ated values have a sm all range and a large variance m ^ns thatthe estim ates are not properly calculated. The range of param eter variance is from zero in the case of error ftee param eters to infinity w ith completely unknown parameters. The result of one gpline segm ent is expre^d in Table 3 with cf as the initial values and % as the estim ates. The estim ated spline andspline locatonparam eters along with their standard deviations are established without the know Hedge of the poiht-to-pohtcoriespondence.

Table 3. Spine param eterand spline location param eter recovery.

Spline location param eters

Im a ge 1 Ima ge 2 Ima ge 3

t t t t t t

0.02 0 33 009 0.41 016 047

0.0415± 00046 0 3615± 00016 00917± 0 0017 0 4158 ± 0 0032 01412 ± 0 0043 0 4617± 00135

Im a ge 4 Ima ge 5 Ima ge 6

t tn t tï t t,

? 018 0 51 025 0.52 0 33 0 57

4 0 2174± 00098 04974± 00079 0 2647 ± 0 0817 0 5472 ± 0 0317 0 3133 ± 0 0127 0 6157 ± 01115

Spline param eters

a10 a11 a.2 a13

332217 7216 4514 27.15 4377 33 69 91

3335 0080 70.4660 48 88529 16.5634 4343.0712 63 0211

± 0 0004 ± 0 0585 ± 08310 ± 12083 ± 0 0004 ± 0 0258

th c10 c11 cL2 c1s

-17.49 13 68 48.82 1015 —27 63 2190

4 -28 .7770 98893 51.9897 81009 39 3702 13 3904

± 0 2193 ± 0 2067 ± 00006 ± 0 0589 ± 07139 ± 10103

If no random noise is added to im age points, the estim ates converge to the true values. The quality of initial estimates is important in the least squares system because it determ ines the iteration number of the system and the accuracy of the convergence. The assumption is that two points on one spline segmentare measured in each image so the totalnumberofequations is2 x 6 (the numberof images) x 2 (the num berof points) + 6 (the num berof the arc length), and the total num berof unknow ns is 12 (the num ber of spline param eters) + 12 (the num ber of spline location param eters). The redundancy (=the number of equations - the number of param eters), that is, the degrees of freedom , is six. W hile some of the geom etric constraints such as slope anddistance observations are dependent onthe extended collinearity equations using splines, other constraints such as slope and arc length increase the nonredundant inform ation in the adjustm ent to reduce the overall rank deficiency of the system .

The coplanarty approach is anotherm athem atical m odel of the perspective relationship betw een the im age and the object space features. The projection plane defined by the perspective center in the im age space and the plane including the straight line inthe object space are identical. B ecause the coplanarty condition is only for straight lines, the coplanarity approach cannot be extended to curves. Object space know ledge about the starting point of a spline can be em ployed in bundle block adjustm entt B ecause the control inform ation about a starting point is available for only three param eters of a total of 12 unknow n param eters to a spline, a spline with control inform ation about a starting point is called a partial control spline. Three spline param eters related to the starting point of a Spline are set to stochastic constraints and the result is seen in Table 4. The total number of equations is 2 x 6 (the numberof images) x 2 (the numberof points) + 6 (the numberof the arc length) = 30, and the total numberof unknowns is 9 (the numberof partial spline parameters) + 12 (the num berof spline

location parameters) = 21 so the redundancy is nine. A convergence of partial spline and spline location param eters has been archived with a partial control spline.

Table 4. Partial spline param eterand spline location param eter recovery.

Spline location param eters

Ima ge 1 Ima ge 2 Ima ge 3

tl t t> t t t

0.04 0 36 009 0.40 014 045

z 0.0525± 00067 0 3547± 0 0020 01128± 00047 0.4157± 0 0091 01575 ± 0 0028 04543± 0 0083

Ima ge 4 Ima ge 5 Ima ge 6

t tn t t11 t t2

021 0 50 0 27 0 54 031 061

01916± 00037 0 5128± 00087 0 2563± 00044 0 5319 ± 0 0056 0 2961 ± 0 0139 0 66239 ± 01147

Spline param eters

a2 ^3 b! b 2 b3

£ 0 7514 -52 87 3071 70 05 40 33 10 98

z 71.7099 -47 2220 -15 8814 62.3703 28 77260 7.1137

± 0 0795 ± 0.6872 ± 2 6439 ± 0 0579 ± 0 £473 ± 17699

£ 0 0 82 -30 72 10 51

71198 ± 0 9483 -35 3841 ± 13403 81557 ± 3 5852

In the next experm ent, spliine location param eters are estm ated w ith know n EO Ps and a full control spline. B ecause spline param eters and gpline location param eters are dependent upon other param eters, the unknowns can be obtained from the model of an observation equation with stochastic constraints. in this experm ent, spline param eters are set to stochastic constraints and the result is seen in Table 5.

Table 5. Spline location param eter recovery.

Spline location param eters

Image 1 Image 2

t t t3 t t t14

0 01 0 37 0 63 0 09 0 44 0 71

z 0 0589± 00015 0 3570± 00076 0 6712± 00197 01134± 00072 0 4175± 00054 0 7069± 00080

Image 3 Image 4

t t t5 tt t6

e 017 046 0 74 021 0 49 081

01757± 00031 04784± 00071 0 7631 ± 00095 0 2039± 00102 0 4869± 00030 0 8122 ± 0 0044

Image 5 Image 6

t t11 t7 tt t2

* 0 0 26 0 53 0 84 029 061 0 89

0 2544 0 5554 0 8597 0 3151 0.6284 0 9013

± 0 0050 ± 0 0069 ± 0 0089 ± 0 0095 ± 0 0052 ± 0 0086

The total number of equations is 2 x 6 (the number of images) x 3 (the number of points) = 36, and the total number of unknowns is 18 (the number of spline location parameters) so the redundancyis 18. Because spline Hocationparam eters are independent of eachother, the arc-length param etprizatdon is not required. The result indicates that a convergence of spline location param eters has been achieved w ith fixed spline param eters considered as stochastic constraints. The proposed m odel is robust w ith respect to the initial approxim atins of goline param eters. The uncertain inform ation related to the representation of a natural cubic goline is described in the dispersion m atrx.

4 3. Recovery of EOPs and Spline Param eters

The object space know ledge of golines is available to recover the exterior orientation param eters in a bundle block adjustm ent. C ontrol spline and partial control goline approaches are applied to verify the feasibility of using control inform ation w ith splines. In both cases, equations of the arc-length param eterizatjonare not necessaryifw e have enoughequations to solve the ^stem because spline param eters are independent of each other:. In the experim ent for a full control spline, the total number of equations is 2 x 6 (the numberof images) x 4 (the numberof points) + 3 (the numberof arc lengths) x 6 (the numberof images) = 66, and the totalnum berof unknowns is 36 (the numberof EOPs) + 24 (the number of goline location param eters) = 60. The redundancy is six. In the case of the partial control goline with one goline œgm ent, the total num berof equations is 2 x 6 (the num berof im ages) x 4 (the number of points) + 3 (the number of arc lengths) x 6 (the number of images) = 66, and the total number of unknowns is 36 (the number of EOPs) + 9 (the number of partial spline parameters) + 24 (the num ber of spline location param eters) = 69. Thus, one m ore segm ent is required to solve the underdetermined system . The total number of equations using two spline segments is 2 x 6 (the number of images) x 4 (the numberof points) x 2 (the number of spline œgments) + 3 (the numberof arc lengths) x 6 (the numberof images) x 2 (the numberof spline segments) = 132, and the total numberof unknowns is36 (the numberof EOPs) + 9 (the numberofpartialgoline param eters) x 2 (the num ber of goline œgm ents) + 24 (the num ber of goline location param eters) x 2 (the num ber of goline œgm ents) = 102. The redundancy is 30. A convergence of the EO Ps of an im age block and the goline param eters has been achieved in both experim ents.

Table 6 expresœs the convergence achievem ent of EO Ps and spline location param eters. The correlation coefficient between param eterXC and <p is high (p « 1) inthe dispersbnm atrix, that is,

two param eters are highly correlated am ong the EO Ps. The correlation coefficient betw een param eters YC and a> is approxim ately 0 85. In general, the correlation coefficient betw een param eters XC and q> is higherthan between param eters YC and a>.

B ecause a control spline provides the object space inform ation about the coordinate ^stem having datum defects of seven, tie spline parameters and EOPs can be recovered simultaneously. In the experim ent of com bined splines, the total num berof equations is 2 x 6 (the num berof im ages) x 3 (the numberof points) x 2 (the numberof splines) + 12 (the numberof arclengtths) x 2 (the numberof golines) = 96, and the total number of unknowns is 36 (the number of EOPs) + 12 (the number of tie goline param eters) + 18 (the num ber of tie spline Hocationparam eters) + 18 (the num ber of control goline location param eters) = 84.

Table 6. EOP and spline location param eterrecovery.

Param eter xc m ] YC [m] ZC [m] co deg] (p deg] k deg]

Image 1 3007 84 400117 50181 8.7090 -9 7976 -12 5845

30015852 ± 0 0154 40012238 ± 0 0215 503.2550 ± 01386 -0 8908 ± 0 3895 0 3252 ± 01351 6 0148 ± 0 8142

Image 2 3308 17 400117 497 52 10.23 8 3144 -5 5004

t 33051962 ± 0 3804 4004.9827 ± 01785 5012641 ± 0 2489 -01247 ± 0 0308 -0 5497 ± 0 0798 -5 2858 ± 0 4690

Image 3 i0 3612.68 3993 37 506 32 5.2731 7 2581 -10135

t 36118996 ± 01226 3995 77891 ± 0 0695 5051299 ± 0 0337 0.1486 ± 04467 01192 ± 0 0168 2 3372 ± 0 0794

Image 4 3619.75 4612 78 506 88 6.2571 -5 3482 183 66

t 3612.7128 ± 0 0258 4613.0145 ± 0 01895 507.0654 ± 0 0251 -0 0921 ± 0 7485 -0152 ± 0 4505 184 5016 ± 0 2289

Image 5 i0 330184 4618 63 49761 -61731 7 5182 187 7145

t 3302 8942 ± 0 0467 4617.0538 ± 0 0249 492.9424 ± 0 0704 -0 6347 ± 01413 0 2662 ± 0 8006 1719808 ± 0 6445

Image 6 ? 2999 59 4615 74 508 49 -71651 -4 8427 1851057

2997.9827 ± 0 0513 4610.1432 ± 0 0249 509.2952 ± 0 0401 -01360 ± 0 5659 -01279 ± 0 6225 1831789 ± 0 2271

Spline location param eters

Image 1 Image 2

tl t t3 to t t t* t2n

i 0 0.04 0 28 0 52 0 76 0 08 0 32 0 56 080

t 0 0432 ± 00033 0.2980 ± 0 0012 0.517 6 ± 0 0 039 0 7705 ± 0 0077 0.0813 ± 0 0082 0.3338 ± 0 0041 0.5715 ± 0 0039 0 8136 ± 0 0069

Image 3 Image 4

t t t15 t,! t t1n t16

f 0 012 0 36 060 0 84 016 0 40 0 64 0 88

001294 ± 0 0036 0.3578 ± 0 0092 0 602 4 ± 0 0 046 0 8437 ± 0 0079 0.1594 ± 0 0115 0.4112 ± 0 0057 0.6418 ± 0 0029 0 9783 ± 0 0037

Image 5 Image 6

t t1 t7 t23 t to t„ t*

¿¡0 020 0.44 068 0 92 0 24 0 48 0 72 0 96

02039 ± 0 0057 0.4461 ± 0 0125 0 671 3 ± 0 0 080 0.9264 ± 0 0061 0 2483 ± 0 0085 0 4860 ± 0 0073 0.7181 ± 0 0084 0 9613 ± 0 0079

Know ledge of objsct.space inform atton about a spline referred to as a full control spline is available prior to aerial triangulation.. A control spline is considered to be a stochastic constraint in the proposed

adjustm ent m odel andthe r^r^^ta^tdonof a control gpüne is the sam e as that of a tie splne. The resullt for com bined spünes that dem onstrates the feasibfHty of using tie spünes and control spünes for bundJe block adjustm ent is ilUusbated in Table 7.

leraton with an inaoreatspüne segment in which a spüne in the im age space does not lie on the projaation of a 3D spüne in the obj3ct space resulte in a divergence of the system . A control spüne is taken to be error free, but in reality this a^um pton is not correct. The accuracy of control spünes is propagat^ into the proposad bundlle block adjustm ent aügortbm , but initial data such as a G IS databan, m aps, ororthcphotos cannotbe witnout error.

Table 7. EO P, control, and tie gpüne param eterrecovery.

Param eter xc mi yc mi zc mi co degi cp [degi k degi

Image 1 3014 87 400718 500 79 0 9740 -8 6517 7 2155

3000 5917 ± 00011 40018935 ± 0 0059 503:2451 ± 01572 -0 0974 ± 01432 0:4297 ± 0 0974 6 6005 ± 02807

Image 2 3315 37 4008 57 503 31 -8 4225 -3 3232 7 2766

3305:1237 ± 00057 4005:0571 ± 0 0043 498 .8916 ± 0 0784 -0 5214 ± 0 3610 -01948 ± 01375 -61421 ± 0 5558

Image 3 £0 3613 85 399117 508 37 -13751 5:3783 4 3148

3609:5400 ± 01576 3995:1419 ± 0 0803 505:1791 ± 0 0428 4:5378 ± 5 4947 1:1746 ± 03610 2 2288 ± 04870

Image 4 £0 3618:46 4617 61 50318 8:5541 2:4287 182 7735

3613:1988 ± 00599 4612:8281 ± 0 0206 507:2056 ± 0 0472 1:1803 ± 0 2578 -0 4:068 ± 0 2979 185 7014 ± 01089

Image 5 £0 3305:71 4620 37 49117 -8 7148 -51487 1831114

3302:9716 ± 0 0718 4617:0808 ± 0 0592 492:9357 ± 0 0660 -06990 ± 01087 1:0485 ± 01437 172 8671 ± 0 2137

Image 6 3002.72 4613 63 49122 8 5475 5 0124 178 2353

2996 9737 ± 00315 4610:8773 ± 0 0672 509:3724 ± 0 0027 -3 2888 ± 0 0688 0:5672 ± 0 3837 182 2693 ± 0 2478

Control gpüne l:cato¡n param etErs

Image 1 Image 2 Image 3

t t t3 t t t* t t9 t15

0.04 0 36 0.66 011 0.41 0 71 016 046 0 76

0 0597 ± 0 0173 0:3495 ± 0 0085 0 6518 ± 0 0065 0 0982 ± 0 0074 0:4085 ± 0 0096 0:7087 ± 0 0067 0:1494 ± 0 0094 0:4499 ± 0 0089 0 7564 ± 0 0156

Image 4 Image 5 Image 6

t trn t16 t tT! t17 t t;, t„

019 0 51 079 0 24 0 54 0J36 029 0 59 0 91

0 2018 ± 0 0043 0:4984 ± 0 0078 0J3065 ± 0 0096 0 2573 ± 0 0086 0:5586 ± 0 0068 0:8553 ± 0 0110 0 3172 ± 0 0088 06137 ± 0 0078 0 8958 ± 0 0085

Tie gpline locatbn param etas

Image 1 Image 2 Image 3

t t t3 t t t4 t t ^5

0 03 0 34 0 67 0 09 039 0 71 014 047 0 73

0.0бВ0 0 3577 0 6694 01141 04032 0 £937 01495 04599 0761B

ç ± 0.0073 ± 0.00б7 ± 0.0033 ± 0.011б ± 0.0073 ± 0.0054 ± 0.0124 ± 0.0075 ± 0.0054

Table 7. Cont

Image 4 Image 5 Image 6

t t1n t16 t t11 tl7 t tl2 t8

021 049 081 026 0 56 0 83 031 0 58 0 92

1 01975 ± 0.0026 0 5109 ± 0 0019 08068 ± 00216 0 2488 ± 0 0773 0 5733 ± 0 0027 0 8527 ± 0 0138 0 3308 ± 0 0034 0 6142 ± 0 0115 0 9018 ± 0 0317

Tie spline param eters

ai0 a2 ^3 b0 bi

3341.44 7313 -48 32 -2072 4337 49 56 97

1 3335 0147 ± 0 0012 712914 ± 0 0478 -47 5124 ± 0.7959 -14 8527 ± 18668 4342 0369 ± 0 0009 62 4762 ± 0 0804

b 2 b3 Cn ct C, C3

36 55 2 57 4416 3 65 -28 22 7 99

1 28 0982 ± 04851 6 8679 ± 14219 515228 ± 0 0008 6.8220 ± 00421 -36 9681 ± 19215 13 0338 ± 2 0048

4 4. Tests with Real Data

In this action., actual experm ents w ith real data are undertaken to verify the feasibility of the proposed bundle block adjustm ent algorithm using splines for the recovery of EO Ps and spline parameters. M edium scale aerialimages covering the area of Jakobshavn Ibrae in W estGreenland are employed for thLis study. The aerial photographs were obtained by Kort and M atrkelstyrel&n (KM S: Danish National Survey and Cadastre) in 1985. KM S established aerial trangulation using GPS ground control points w ith a ± 1 pixel root m ean square error under favorable circum stances and im ages were oriented tothe W GS84 reference frame. Technical inform atononthe aerial im ages is described in Table 8.

Table 8. Inform atton about aerial im ages used in this study.

Vertical aerial photograph

Data 9 July 1985

Origin KM S

Focal length 87 75 mm

Photo scale 1150,000

Pixel size 12 pm

Scanning resolution 12

Ground sam pling 19 m

distance

The diapositive films were scanned with a RasterM aster phottogrammetric precision scanner, which has a m axim um im age resolution of 12 jam and a scan dim ension of 23 cm x 23 cm to obtain digital im ages for a ^ftcopy workstation as ^en in Figure 6.

Figure 6. Test images. (a) Image 762, (b) Image 764, (c) Image 766, and (d) Targetarea.

(c) (d)

The first experim ent is the recoveryof spline param eters w ithknownEO Ps obtainedbym anual operation using a softcopy workstation. A spline consists of four parts and the second segment param eters aie recovered. The total num ber of equations is 2 x 3 (the number of images) x 3 (the num ber of points) + 2 (the num ber of arc lengths) x 3 (the num ber of im ages) = 24, andthe total number of unknowns is 12 (the number of gpline param eters) + 9 (the number of spline location param eters) = 21 so the redundancy is three. Table 9 show s the convergence achievem entof spline and gpline location param eters.

Table 9. Spline param eterand gpline location param eter recovery.

Spline location param etErs

Im age 762 Im age 764

t t t tt t t

0.08 0 38 0.72 0 22 0 53 0 82

0.0844± 0 0046 0.4258± 00058 0£934± 00072 0 2224± 00175 0 5170± 00104 0 8272± 00156

Im age 766

0 32 0 59 0 88

0 3075 ± 0 0097 0.6176± 00148 0 9158± 00080

Table 9. Cont.

Spline param a11 sters a.2 ^3

535000.00 830 00 150 00 50 00

5353941732 ± 01273 867.6307± 0.7142 173 1357± 7 6540 243213± 213379

? 7671000 00 150 00 140 00 -300 00

7672048 3173 ± 0 2237 1431734 ± 16149 130 8147 ± 10 9058 -2901270± 26 7324

C.2 Cas

? 0 00 -10 00 50 00 50 00

4 21913 ± 0 0547 -3 7669± 01576 39 8003 ± 91572 27 7922 ± 19 6787

Estimation of spline parameters including their location parameters is established by the relationship betw een splines inthe object space and their projsctbninthe im age space w ithout the know ledge of the point-to-point correspondence. B ecause bundle block adjustm ent using splines does not require conjugate points generated by point-to-point correspondence know ledge, a m ore robust and flexible m atching algorithm can be adopted. Table 10 show s the available object space inform ation without know ledge of the point-to-point correspondence (the full control spline) • A H locations are assum ed as lying on the second spline segm entand the second spline segm ent as calculated from the ^ftcopy workstation is used as control inform ation.

Table 10. Spline location param eter recovery.

Spline location param eters

Im age 762 Imag e 764 Imag e 766

t t t t t t

015 0 60 030 0 75 0.45 0 90

4 01647± 00048 0 6177± 00091 0 2872 ± 0 0034 0 7481 ± 0 0093 0.4362± 0 0155 0 9249± 00087

The next experment is the recoveryof EO Ps w itth a control spline. The spline control points are (534415 91, 767199305, -18 97), (535394 52, 7672045.02, 2 127), (536110 66, 7672024 29, -13 £97), and (536654 0)4, 7671016 20, -2 51). Even though edge detectors are often used in digital photogram m etry and rem ote sensing software, the control points are extracted m anually because edge detectin is not our m ain goal. A m ong the three œgm ents, the second gpline segm ent is used for the EO P recovery. The inform ation of the control gpline is obtained by a m anual operation using the ^ftcopy workstation with an estim ated accuracy of ±1 pixel. The convergence radius of the proposed iterative algorithm is proportbnal to the estim ated accuracy level. The im age coordínate system is converti into the photo coordínate System using the interior orientation param eters from KM S. The association between a pointon a 3D spline segmentand a pointon a 2D image isnot establ&ned in tthis study. Of courra, 3D goline measurementin the stereo model using the ^ftcopy workstation cannot be without error ^ the accuracy of the control goline ispropagated into the recovery of EOPs. The result is illustrated in Table 11. The spline control inform ation is utilizó as sttochastic constraints in the adjustm entm odel. B ecause adding these constraints rem oves the rank deficiency of the Gauss-M arkov

m odel corresponding to spline param eters that are dependentupon spline location param eters, a bundle block adjustm ent can be m ade using only the extended collinearity equations for natural cubic splines.

Table 11. EO P and spline location param eter recovery.

Image xc mi yc mi zc mi co Icieg] (p dag] K meg]

762 r 547000.00 7659000.00 1400000 3.8472 2.1248 91.8101

1 547465 37 ± 15 0911 7658235.41 ± 13 0278 1370025 ± 5.4714 0.3622 ± 0 8148 0 5124 ± 01784 91.5124 ± 01717

764 r 546500.00 7670000 00 1350000 0.1125 0 6128 90.7015

1 546963 22 ± 12 5460 7672016.87 ± 171472 13708 82 ± 71872 -0 3258 ± 0 £913 -0 5217 ± 0 8632 91.1612 ± 11004

766 r 546000.00 768000.00 1370000 1.4871 5.9052 92.0975

1 546547 58 ± 13 8104 7685836.75 ± 12 1486 13712 20 ± 8 4854 1.2785 ± 14218 0.5468 ± 11957 92.9796 ± 0.6557

Spline location param etErs

Imag e 762 Imag e 764

t tt tt t1n t tt

008 0 32 0 56 080 016 0 40

00865± 0 0097 0 3192± 00159 0 5701 ± 0 0072 0 8167 ± 0 0088 01759 ± 0 0067 0 4167± 00085

Imag e 764 Imag e 766

t t! t3 tt tt t2

0 664 0 88 0 24 0 48 072 0 96

066557± 00131 0 8685± 00092 0 2471 ± 0 0086 0 4683 ± 0 0069 0 7251 ± 0 0141 0 9713± 00089

5. Conclusions

In this paper, traditional least squares of a bundle block adjustm ent process have been augm ented bysupport splines insteadof conventional point features. E stim ationof EO Ps andspline param eters including location param eters is established by the relationship betw een splines in the object space and their projection into the image space without any knowledge of the point-to-point correspondence. B ecause bundle block adjustm ent using splines does not require conjugate points generated by the point-to-point correspondence know ledge, a m ore reliable and flexible m atching algorithm can be adopted. Point-based aerial triangulation with experienced hum an operators is effective for traditional phottogram m etric activities but is not appropriate w ithin the autonom ous environm ent of digital phottogram m etry. Feature-based aerial triangulation is suitable for the developm en. of reliable and accurate auttom ation techniques. If linear features are em ployed as control features, they provide advantages over point features in aerial triangulation auttom ation. Point-based aerial triangulation based on m anual m easurem ent and the identification of conjugate points is less reliable than feature-based aerial triangulation because it has the lim itationLs of visibility (occllusionL), am biguity (repetitive patterns), and semantic information in the light of robust and appropriate automation. Auttom ation of aerial trangulation and pose estim ation is obstructed by the correspondence problem ,

but the em ploym ent of splines is one w ay to overcom e occlusion and am biguity issues. The m anual identification of corresponding entities in two im ages is crucial in the autom ation of photogram m etric tasks. A further problem of point-basedapproaches is theirw eak geom etric constraints as com pared with feature-based methods, so accurate initial values for the unknown parameters are required. Feature-based aerial triangulation can be im plem ented without con.jugate points because the m easured points in each im age are not the conjugate points in this proposed adjustm ent m odel. Thus, tie splines that do not appear in all the overlapped im ages together can be em ployed in feature-based aerial triangulation. Another advantage of employing splines is that the adoption of high level features increases the feasibility of geom etric inform ation and provides an appropriate analytical solution that em phasizes the redundancy of aerial triangulation.

3D linear features expreœed by 3D natural cubic gplines are em ployed as the m athem atical m odel of linear features in the object space and its counterpart in the projected im age space for bundle block adjustm entt To solve overparam eterization of 3D natural cubic splines, arc-length param eterization using Simpson's rule is developed, and in the case of straight lines and conic sections, gpline tangents can be additional equations to the overparam eterized ^stem . Photogram m etric triangulation by the proposed model, including the extended œllinearty and arc-length param eterization equations, is developed to show the feasibility of tie and control gplines for the estim ation of the exterior orientation of multiple im ages, splines, and gpline location param eters. A useful stochastic constraint for a gpline œgm entisexam ined for its utility to becom e a fullorpartial control gpline arch as known EO Ps with a tie, partial control, and full control spline, and unknow n EO Ps with a partial and full control spline. In addition, the inform ation content of an im age gpline is calculated and the feasibility of a tie spline and a control spline for a block adjustment is described. A simulation bundle block adjustment is im plem ented prior to the actual experiment with real data in order to evaluate the perform ance of the proposed algorithms. A simulation can control the measurement errors ^ that random noiœs m inim ally affect the overall geom etrry of a block:. The individual observations are generated based on the general situation of bundle block adjustment to estimate the properties of the proposed algorithm s. A sim ulation allow s adjustm ent for geom etric problem s or varying conditions w ithin individual experim ents.

Acknow ledgem ents

This research was supported by a grant(07KLSGC04) from Cutting-edge Urban Development— Korean Land Spatializatin Research Project funded by M inistry of Construction & Transportation of Korean governm ent.

References and Notes

1. Schernk, T. Digital Photogram m etry; Terra Science: Laurelville, OH , USA , 1999; Vol. 1.

2. Canny, JT. A Computational Approach to Edge Detection. IEEE Trans. Pattern Anal. Mach. Intell 1986,8,679-698.

3. Forrstner, W .; G ulch, E . A Fast O perator ior D etectionLand Precise Locationof D istinct Points, Corners and Centres of Circular Features. In ISPRS Intercom mission Workshop, Interlaken, Switzerland, Unterlaken, Sw itzerand, 1987; pp. 281-305.

4. Harric, C .; Stephens, M . A CombinedComer andEdge D etector. mProceedings of 4thAlvey Vision Conference, Cam bridge, UK, 1988; pp. 147-151.

5. M oravec, H . Towards Automatic Visual Obstacle Avoidance. In Proceedings of the International JointConfarence on Artificial Intelligence, Cambridge, M A , USA . 1977; p. 584.

6. Prewitt, J. Object Enhancement and Extraction in Picture Processing and Psychopictorics; Academic Press: New York,NY ,USA , 1970.

7. Sobel, I. Camera M odels and on Calibrating Computer Controlled Cameras for Perceiving 3D Scenes. Artif. Intell. 1974, 5, 185-198.

8. Sm ith, S ; Brady, J. SU SAN—A New Approach to Low LevelIm age Processing. Deface Research Agency,FamboxoughL, UK 1995. Tech.Rep. TR95SM S1c.

9. M arr, D . Vision. Freeman Publishers:New York, NY, USA, 1982.

10. Forrstner, W ; Gulch, E . A Feature Based Correspondence Algorithm for Im age M atching. Int. Arch. Phottogrramm . Rem ote Sens. 1986,26, 13-19.

11. Hannah, M . A System forDigial Stereo Image M atching.Phottogramm . Eng. Remote Sens. 1989, 55,1765-1770.

12. Schenk, T ; Totth, C .Li, J. Tow ards an Autonom ous System for O renting D igital Stereopairs. Phottoglamm . Eng. Remote Sens. 1991, 57, 1057-1064.

13. Schenk, T ; Toth, C . Towards a FullyAutom ated AelO-tr:angulat;bh System . InACSM /A.SPRS, Annual Convention TechnicalPapels, M inneapolis, MN,USA, 1993; Vol. 3, pp. 340-347.

14. Ebner, H ; Ohlhof, T . Utilization of GroundControlPoints ior Image Orentation withoutPoint Identification in Im age Space. Int.ABch. Photbglamm . Remote Sens. 1994, 30, 206-211.

15. Ackerman, F; Tsingas, V. Automatic Digital Aerial Trangulation. In Proceedings of the ACSM /ASPRS Annual Convention,Reno, NV,USA . 1994; pp. 1-12.

16. Haala, N ; Vosielman, G . RecogniionLof Road and River Patterns byRelational M atching. Mr Arch. Phottogramm . Remote Sens. 1992, 29, 969-975.

17. Drewniok, C ; Rohr, K . Automatic Exterior Orientation of Aerial Images in Urban Environm ents. Int. Arch.Phottogramm . Remote Sens. 1996, 31, 146-152.

18. Drewniok, C ; Rohr, K . Exterior Orientation—An Automatic Approach Based on Fitting Analytic Landmark M odels. ISPRS J". Photogramm . Remote Sens. 1997, 52, 132-145.

19. Zahran, M . Shape-based Multi-image M atching using the Principles of Generalized Hough Transform . PhD . Dissertation,Department of Civil and Environmental Engineering and Geodetic Science, Ohio State University: Columbus, OH, USA, 1997.

20. Schenk, T . D eterrm ining Transform ationParram eterrs between Surfaces w ithout Identical Points; Department of Civil and Environm ent Engineering and Geodetic Science, Ohio State University: Columbus, OH ,USA, 1998; Tech. Rep. No. 15.

21. Sarkar, S; Boyer, KL. Computing Perceptual Organization in Computer Vision; World Scientific: RiverEdge, NJ,USA , 1994.

22. Guy, G ; M edioni, G . Inferring Global Perceptual Contours from Local Features. Int. J. Com put. Vis.1996,20,113-133.

23. Casadei S.; M iter, S. Hierarchical Image Segmentation—Part 1: Detection of Regular Curves in a VectorGraph. Int. J. Com put. Vis. 1998, 27, 71-100.

24. Boyer, K L.; Sarkar, S. Perceptual O ranizaton in Computer Vision: Status, Challenges and Potential. Com put. Vis. Im age Underst. 1999, 76, 1-5.

25. M asry, S E.D igitalM apping using Entities- A New Concept. Phottogramm . Eng. Remote Sens. 1981,48,1561-1565.

26. Heikkila, J. Use of Linear Features inDigital Phottogram m etry. Phottogramm .J. Finl. 1991,12, 40-56.

27. Kubik, K . Phottogram m etric Restitution Based on Linear Features. lit. Arch. Phottogramm . Remote Sens. 1992, 29, 687-690.

28. Petsa, E .; Patias, P. Formulation and A sœs^ ent of Straight Line B asedA lgorithm s for D igital Phottogram m etry. Int Arch. Phottogramm . Remote Sens, M elbourne, Australia 1994, 30, 310-317.

29. Gulch, E . Line Photogramm etry: A Tool for Preciœ Localization of 3D Points and Lines in Automated Object Reconstruction. In Proceedings of Integrating Phottogram m etric Techniques with Scene Analysis and M achine Vision II SPIE, Orlando, FL, USA, 1995; pp. 2-12.

30. W iman, H .; Axels&on, P. Finding 3D -structures in M u iltiple Aerial Images using Lines and Regions. Int. Arch. Phottogramm . Remote Sens,Vienna, Austria 1996, 31, 953-959.

31. Chen, T ; Shiasaki, R S. Determination of C amera's Orientation Param eters Based on Line Features. Int. Arch. Phottogramm . Remote Sens. 1998, 32, 23-28.

32. Habib, AF.A erial Triangulaticnusing Point and Linear Features. mProoeedings of the ISPRS Conference on Auttomatic Extraction of G IS Objects from D igital Imagery, M unich, G erm any, 1999;pp.137-142.

33. H euvel, F. A L ine-phottogram m etric M athem atical M odel for the Recbrlstr]uotion of Polyhedral Objecte. In Videom etdcs VI, Proceedings of SP IE, D enver, CO, USA, 1999; Vol. 3641, pp. 60-71.

34. Tom m a sell i, a ; Poz, A . Line Based Orientation of Aerial Images. In Proceedings of the ISPRS Conference on Auttom atic Extraction of G IS O bjects from D igital Im agery, M unich, G erm any, 1999;pp.143-148.

35. Voœelman, G ; Veldhuis, H. M apping by Dragging and Fitting of W ire-frame M odels. Phottogramm . Eng. Remote Sens. 1999, 65, 769-776.

36. Forrstner, W . New O rdentationProcectuiies. Int. Arch. Phottogramm . Remote Sens. Spat. Inf. Sci. 2000,33,297-304.

37. Sm ith, M J; Park, DWG. Absolute and Exterior O rrientationushg Linear Features. Int. Arch. Phottogramm . Remote Sens. 2000, 33, 850-857.

38. Schenk:, T . Towards a Feature-based Phottogram m etry. Swed. Soc. Phottogramm . Remote Sens. 2002,1,143-150.

39. Tangelder, J; Erm es, P ; VosEEllm an, G ; H euvel, F. CAD -based Phottogram m etry for Reverse Engineering of Industrial Installation. Com put-Aided Civ. Infrastruc. Eng. 2003, 18, 264-274.

40. Parian, JA ; Gruen, A. Panoramic Camera Calibration using 3D StraightLines. In Proceedings of International Archives of. Phottogram m etry, Remote Sensing and Spatial Inform ation. Sciences, Berlin,Germ any, 2005.

41. Mulawa, D C ; M ikhail, EM . Phottogram m etric Treatment of Linear Features. Int. Arch. Phottogramm . Remote Sens. 1988, 27(B10),383-393.

42. M ulawa, D C . Estimation and Photogram metric Treatmentof Linear Features. PhD . Dissertation, Department of Civil Engineering, Purdue University": W est Lafayette, IN ,USA, 1989.

43. Mikhail, EM .; W eerawong, K. Feature-based Photogram m etric Object Construction. in Proceedingsofthe ASPRS/ACSM AnnualConventon, Reno, NV, USA . 1994; pp. 403-407.

44. Tommaseli A .; Tozzi, C. A Filtering-based Approach to Eye-in-hand Robot Vision. in Proceedings of IiternatinalArchives of Photogrammetry and Remote Sensing, W ashington, DC , USA, 1992; pp. 182-189.

45. Habib, A F.; Shin, S W .; M organ, M . Automatic Pose Estimation of Imagery using Free-form Control Linear Features. int. Arch. Photogramm . Rem ote Sens. Spat. inf. Sci. 2002b, 34, 150-155.

46. Habib, A ; M organ, M ; Kim, EM ; Cheng, R . Linear Features in Photogram m etric Activities. In XXth ISPRS Congress, Istanbul, Turkey, 2004; p. 610.

47. Mikhail, EM . Linear Features for Photogram m etric Restitution and Object Completion. In SPIE Proceedings-integrating Photogram metric Techniques with Scene Analysis and M achine Vision, Orlando, FL, USA, 1993.

48. Habib, AF; Asm am aw, A; Kelley, D . May, M . Linear Features in Photogrammetry; Department of Civil and Environmental Engineering and Geodetic Science, Ohio State University: Columbus, OH, USA 2000.

49. Schenk:, T. From Point-based to Feature-based Aerial Triangulation. ISPRS J. Photogramm . Rem ote Sens. 2004, 48, 315-329.

50. Zalm anson, G . H ierarchical Recovery of Exterior O rientation from Param etric and Natural 3D Curves. Ph D.D is^rtaton, D epartm ent of C ivil and Environm ental Engineering and G eodetc Science, Ohio State University: Columbus, OH, USA, 2000.

51. Besl P ; McKay, N. A Method for Registration of 3-D Shapes. IEEE Trans. Pattern Anal. Mach. IiteH. 1992, 14, 239-256.

52. Rabbani, T ; Dijkm an, S ; H euvel F ; Vo^elm an, G . An Integrated Approach for M odeling and GlobalRegistration of Point Clouds. ISPRS J. Photogramm . Remote Sens. 2007, 61,355-370.

53. Akav, A ; Zalm anson, G H ; Doytsher, Y . Linear Feature Based Aerial Triangulation. In XXth ISPRS Congress, Istanbul, Turkey Comm issbn m, 2004.

54. Lin, H T. Autonomous Recovery of Exterior Orientation of Imagery using Free-form Linear Features. PhD . Dissertation, Department of Civil and Environmental Engineering and Geodetic Science, Ohio State University: Columbus, OH, USA, 2002.

55. Gruen, A ; Akca, D . Least Squares 3D Surfface and Curve M atching. ISPRS J. Photogramm . Rem ote Sens. 2005, 59, 151-174.

56. W ang, H ; Kearney, J; Atkinson, K . Arc-length Parameterized Spline Curves for Realtime Simulation. In Proceedings of the 5th International Conference on Curves and Surfaces, San M alo,France, 2002; pp. 387-396.

57. Guenter, B . Parent, R. Computing the Arc Length of Parametric Curves. IEEE Comput. Graph. Appl. 1990, 10, 72-78.

58. Nasri, A H ; van Overveld, C W A M ; W yvil, B. A Recursive Subdivision Algorithm for Piecew is Circular Spline. C om put G raph. Forum 2001, 20, 35-45.

S9. Ji, Q .; Costa, S.; Haralick:, M .; Shapiro, L. A Robust Linear Least Squares Estimation of Camera Exterior O rientaton using M ultiple Geometric Features. ISPRS J. PhotDglaшm . Remoíe Sens. 2000,SS,7S-93.

© 2009 bythe authors; licen^see M olaCuLan D ielsityPresalvaЗih]htelnaЗ;bпal, Basel, Switzerland. This article is an open-access article distributed under the term s and aonditonLS of the C reative C om m ons A ttibutjon license (http //creatvecom m ons brз/liœпsas/foy/3 0/).

Copyright of Sensors (14248220) is the property of MDPI Publishing and its content may not be copied or emailed to multiple sites or posted to a listserv without the copyright holder's express written permission. However, users may print, download, or email articles for individual use.