5.3 Least-Squares Adjustment of Level Nets, Triangulation, and Plane Traverses

Key Takeaways

  • Parametric least squares expresses observations as direct functions of unknown parameters (V = AX - L), solving for unknown coordinate or elevation updates via X = (A^T P A)^(-1) A^T P L.
  • 1D leveling net adjustments construct design matrix A with entries +1 and -1 based on level line directional topology, weighting each line inversely to its length.
  • Non-linear observation equations in 2D triangulation/trilateration must be linearized using first-order Taylor series expansion around initial approximate station coordinates.
  • Least-squares plane traverse adjustment is superior to the empirical Compass (Bowditch) Rule because it incorporates rigorous instrument standard errors and computes coordinate error ellipses.
  • Standardized residuals (|v_i*| > 3.0) flag potential blunders in large surveying networks for rejection or re-observation.
Last updated: July 2026

Least-Squares Adjustment of Level Nets, Triangulation, and Plane Traverses

Modern geodetic data processing relies on rigorous matrix formulation of least squares to adjust complex 1D elevation networks, 2D triangulation/trilateration control figures, and plane traverses. This section details parametric matrix techniques ($A^T P A X = A^T P L$), non-linear model linearization, error ellipse computation, and residual quality control tests required for GELE mastery.


1. Parametric (Observation Equations) vs. Correlative (Condition Equations) Methods

Least-squares adjustments are executed through two primary formulation approaches:

FeatureParametric Method (Observation Equations)Correlative Method (Condition Equations)
Primary UnknownsParameter corrections to point coordinates or elevations ($X$)Lagrange multipliers (correlatives $K$)
FormulationExpresses each observation directly as a function of parameters ($V = AX - L$)Expresses exact geometric condition closure equations ($B V = W$)
Number of EquationsEquals number of parameters $u$Equals redundancy / degrees of freedom $r = n - u$
Matrix Normal System$(A^T P A)_{u \times u} X = A^T P L$$(B P^{-1} B^T)_{r \times r} K = W$
Best ApplicationComplex 2D/3D networks, GNSS vectors, hybrid observationsSimple closed leveling loops, geometric polygon triangles

2. Rigorous Matrix Formulation of Parametric Least Squares

For a system of $n$ observations and $u$ unknown parameters ($n > u$), the linear observation model in matrix notation is:

V=AXLV = A X - L

where:

  • $V_{n \times 1}$: Vector of observation residuals ($[v_1, v_2, \dots, v_n]^T$).
  • $A_{n \times u}$: Design Matrix containing partial derivatives of observation equations with respect to unknown parameters ($\partial f_i / \partial x_j$).
  • $X_{u \times 1}$: Vector of unknown parameter corrections ($[\delta x_1, \delta x_2, \dots, \delta x_u]^T$).
  • $L_{n \times 1}$: Vector of misclosures ($L = L_{observed} - L_{computed}^0$).
  • $P_{n \times n}$: Weight Matrix (diagonal matrix $P = \text{diag}(p_1, p_2, \dots, p_n)$ for independent observations).

2.1 Matrix Derivation of Normal Equations

Minimizing the weighted residual sum $V^T P V = (A X - L)^T P (A X - L)$:

VTPV=XTATPAX2XTATPL+LTPLV^T P V = X^T A^T P A X - 2 X^T A^T P L + L^T P L

Taking partial derivatives with respect to $X$ and setting to zero:

(VTPV)X=2ATPAX2ATPL=0    (ATPA)X=ATPL\frac{\partial(V^T P V)}{\partial X} = 2 A^T P A X - 2 A^T P L = 0 \implies (A^T P A) X = A^T P L

Solving for parameter corrections $X$:

X=(ATPA)1ATPLX = (A^T P A)^{-1} A^T P L

2.2 Residuals and Variance-Covariance Matrices

  • Adjusted Residuals Vector: $V = A X - L$
  • A Posteriori Reference Variance ($\hat{\sigma}_0^2$): σ^02=VTPVnu=VTPVr\hat{\sigma}_0^2 = \frac{V^T P V}{n - u} = \frac{V^T P V}{r}
  • Parameter Variance-Covariance Matrix ($\Sigma_X$): ΣX=σ^02(ATPA)1=[σx12σx1x2σx2x1σx22]\Sigma_X = \hat{\sigma}_0^2 (A^T P A)^{-1} = \begin{bmatrix} \sigma_{x_1}^2 & \sigma_{x_1 x_2} & \dots \\ \sigma_{x_2 x_1} & \sigma_{x_2}^2 & \dots \\ \vdots & \vdots & \ddots \end{bmatrix}
  • Standard errors of adjusted parameters: $\sigma_{x_j} = \sqrt{\Sigma_{X, jj}}$.

3. 1D Level Network Adjustment: Step-by-Step Worked Matrix Example

Network Geometry

Consider a level net with two fixed control benchmarks BM-1 ($Elev = 10.000\text{ m}$) and BM-2 ($Elev = 25.000\text{ m}$) and two unknown junction points A ($H_A$) and B ($H_B$). Four elevation differences are measured:

Obs LineFromToObs $\Delta h$ (m)Distance $D$ (km)Weight $p_i = 1/D$
Line 1BM-1A$+5.100$$1.0$$1.0$
Line 2AB$+9.800$$2.0$$0.5$
Line 3BBM-2$+0.200$$1.0$$1.0$
Line 4BM-1B$+15.000$$2.0$$0.5$

Initial approximate elevations: $H_A^0 = 15.100\text{ m}$, $H_B^0 = 24.900\text{ m}$.

Step 1: Form Observation Equations ($v_i = a_{i1} \delta H_A + a_{i2} \delta H_B - l_i$)

  • Line 1: $H_A - 10.000 = 5.100 + v_1 \implies \delta H_A - (5.100 - 5.100) = \delta H_A - 0.000$
  • Line 2: $H_B - H_A = 9.800 + v_2 \implies \delta H_B - \delta H_A - (9.800 - (24.900 - 15.100)) = -\delta H_A + \delta H_B - 0.000$
  • Line 3: $25.000 - H_B = 0.200 + v_3 \implies -\delta H_B - (0.200 - (25.000 - 24.900)) = -\delta H_B - 0.100$
  • Line 4: $H_B - 10.000 = 15.000 + v_4 \implies \delta H_B - (15.000 - 14.900) = \delta H_B - 0.100$

Step 2: Assemble Matrices

A=[10110101],P=[1.000000.500001.000000.5],L=[0.0000.0000.1000.100]A = \begin{bmatrix} 1 & 0 \\ -1 & 1 \\ 0 & -1 \\ 0 & 1 \end{bmatrix}, \quad P = \begin{bmatrix} 1.0 & 0 & 0 & 0 \\ 0 & 0.5 & 0 & 0 \\ 0 & 0 & 1.0 & 0 \\ 0 & 0 & 0 & 0.5 \end{bmatrix}, \quad L = \begin{bmatrix} 0.000 \\ 0.000 \\ 0.100 \\ 0.100 \end{bmatrix}

Step 3: Compute Normal Matrix $N = A^T P A$

ATP=[10.50000.51.00.5]A^T P = \begin{bmatrix} 1 & -0.5 & 0 & 0 \\ 0 & 0.5 & -1.0 & 0.5 \end{bmatrix} N=ATPA=[(1+0.5)(0.5)(0.5)(0.5+1.0+0.5)]=[1.50.50.52.0]N = A^T P A = \begin{bmatrix} (1 + 0.5) & (-0.5) \\ (-0.5) & (0.5 + 1.0 + 0.5) \end{bmatrix} = \begin{bmatrix} 1.5 & -0.5 \\ -0.5 & 2.0 \end{bmatrix}

Step 4: Compute Right-Hand Side Vector $A^T P L$

ATPL=[1(0)0.5(0)+0(0.1)+0(0.1)0(0)+0.5(0)1.0(0.1)+0.5(0.1)]=[0.0000.050]A^T P L = \begin{bmatrix} 1(0) - 0.5(0) + 0(0.1) + 0(0.1) \\ 0(0) + 0.5(0) - 1.0(0.1) + 0.5(0.1) \end{bmatrix} = \begin{bmatrix} 0.000 \\ -0.050 \end{bmatrix}

Step 5: Invert Normal Matrix $N^{-1}$

det(N)=(1.5)(2.0)(0.5)(0.5)=3.00.25=2.75\det(N) = (1.5)(2.0) - (-0.5)(-0.5) = 3.0 - 0.25 = 2.75 N1=12.75[2.00.50.51.5]=[0.727270.181820.181820.54545]N^{-1} = \frac{1}{2.75} \begin{bmatrix} 2.0 & 0.5 \\ 0.5 & 1.5 \end{bmatrix} = \begin{bmatrix} 0.72727 & 0.18182 \\ 0.18182 & 0.54545 \end{bmatrix}

Step 6: Solve for Parameter Corrections $X = N^{-1} A^T P L$

X=[0.727270.181820.181820.54545][0.0000.050]=[0.00909 m0.02727 m]X = \begin{bmatrix} 0.72727 & 0.18182 \\ 0.18182 & 0.54545 \end{bmatrix} \begin{bmatrix} 0.000 \\ -0.050 \end{bmatrix} = \begin{bmatrix} -0.00909\text{ m} \\ -0.02727\text{ m} \end{bmatrix}

Step 7: Calculate Final Adjusted Elevations

  • $\hat{H}_A = H_A^0 + \delta H_A = 15.100 - 0.0091 = 15.0909\text{ m}$
  • $\hat{H}_B = H_B^0 + \delta H_B = 24.900 - 0.0273 = 24.8727\text{ m}$

4. Linearization of 2D Triangulation and Trilateration Equations

In 2D horizontal networks, distance and angle observation equations are non-linear functions of station coordinates $(X_A, Y_A, X_B, Y_B)$. They must be linearized using first-order Taylor series expansion around initial approximate coordinates $(X^0, Y^0)$.

Distance Observation Equation Linearization

SAB=(XBXA)2+(YBYA)2S_{AB} = \sqrt{(X_B - X_A)^2 + (Y_B - Y_A)^2}

δSAB=cos(αAB)δXAsin(αAB)δYA+cos(αAB)δXB+sin(αAB)δYB\delta S_{AB} = -\cos(\alpha_{AB}) \, \delta X_A - \sin(\alpha_{AB}) \, \delta Y_A + \cos(\alpha_{AB}) \, \delta X_B + \sin(\alpha_{AB}) \, \delta Y_B

where $\alpha_{AB}$ is the grid azimuth from station A to B.

Azimuth Observation Equation Linearization

αAB=arctan(XBXAYBYA)\alpha_{AB} = \arctan\left(\frac{X_B - X_A}{Y_B - Y_A}\right)

dαAB=cos(αAB)SABδXAsin(αAB)SABδYAcos(αAB)SABδXB+sin(αAB)SABδYBd\alpha_{AB} = \frac{\cos(\alpha_{AB})}{S_{AB}} \, \delta X_A - \frac{\sin(\alpha_{AB})}{S_{AB}} \, \delta Y_A - \frac{\cos(\alpha_{AB})}{S_{AB}} \, \delta X_B + \frac{\sin(\alpha_{AB})}{S_{AB}} \, \delta Y_B

The adjustment is iterated until parameter update vector $|X| < \epsilon$ (e.g. $< 0.1\text{ mm}$).


5. Plane Traverse Adjustment: Least Squares vs. Compass Rule

For decades, boundary surveys in the Philippines were adjusted using the empirical Compass (Bowditch) Rule. However, modern regulations under DENR Administrative Orders specify Least Squares for control networks.

Comparison ParameterCompass (Bowditch) RuleLeast-Squares Adjustment
Theoretical RigorEmpirical rule; assumes distance & angle errors equalMathematically rigorous; minimizes $V^T P V$
Observation WeightingIgnores individual instrument accuraciesIncorporates individual distance/angle standard errors
Geometry HandlingLimited to simple closed loop traversesHandles complex interconnected networks, cross-ties, and GNSS
Error EstimationProvides no statistical variance for coordinatesComputes full covariance matrix $\Sigma_{XY}$ and confidence ellipses

6. Point Error Ellipses and Residual Analysis

6.1 Point Error Ellipses (Confidence Ellipses)

For a 2D station $i$, the coordinate covariance matrix is extracted from $\Sigma_X$:

ΣXY=[σx2σxyσxyσy2]\Sigma_{XY} = \begin{bmatrix} \sigma_x^2 & \sigma_{xy} \\ \sigma_{xy} & \sigma_y^2 \end{bmatrix}

The geometric parameters of the standard ($39.4%$ 2D probability) error ellipse are:

  1. Orientation angle of semi-major axis ($\theta$): tan(2θ)=2σxyσx2σy2\tan(2\theta) = \frac{2 \sigma_{xy}}{\sigma_x^2 - \sigma_y^2}

  2. Semi-major axis ($a$) and Semi-minor axis ($b$): a=12(σx2+σy2+(σx2σy2)2+4σxy2)a = \sqrt{ \frac{1}{2} \left( \sigma_x^2 + \sigma_y^2 + \sqrt{(\sigma_x^2 - \sigma_y^2)^2 + 4 \sigma_{xy}^2} \right) } b=12(σx2+σy2(σx2σy2)2+4σxy2)b = \sqrt{ \frac{1}{2} \left( \sigma_x^2 + \sigma_y^2 - \sqrt{(\sigma_x^2 - \sigma_y^2)^2 + 4 \sigma_{xy}^2} \right) }

To scale to a $95%$ confidence ellipse, multiply $a$ and $b$ by $K = \sqrt{\chi_{2, 0.05}^2} = \sqrt{5.991} \approx 2.4477$.

                 Y (North)
                     ^       / Semi-Major Axis (a)
                     |      / 
                     |     /---. 
                     |   .'  |  `.
                     |  /    |    \
                     | |     +-----|---> X (East)
                     |  \\   /     /
                     |   `.  \\  .'  
                     |     `---' \\ Semi-Minor Axis (b)
                     |            \\ Angle theta
                     +----------------------->

6.2 Residual Analysis and Outlier Detection

To detect blunders in large observations sets, Standardized Residuals ($v_i^*$) are evaluated:

vi=viσvi=viσ^0(P1A(ATPA)1AT)iiv_i^* = \frac{v_i}{\sigma_{v_i}} = \frac{v_i}{\hat{\sigma}_0 \sqrt{ (P^{-1} - A (A^T P A)^{-1} A^T)_{ii} }}

If $|v_i^*| > 3.0$, the observation is flagged as a statistical outlier (blunder) at the $99.7%$ confidence level and subjected to rejection or re-measurement.

Test Your Knowledge

In parametric least-squares matrix formulation, what is the dimensions of the normal equations matrix (A^T P A) for a network with n observations and u unknown parameters?

A
B
C
D
Test Your Knowledge

What matrix represents the full a posteriori covariance of adjusted parameter corrections in parametric least squares?

A
B
C
D
Test Your Knowledge

Which formula gives the orientation angle theta of the semi-major axis of a point error ellipse derived from coordinate variances sigma_x^2, sigma_y^2 and covariance sigma_xy?

A
B
C
D
Test Your Knowledge

What standardized residual threshold |v_i*| is standard in geodetic surveying data processing for identifying a statistical outlier (blunder)?

A
B
C
D