Fast Inverse Georectification

• Fast Inverse Georectification is a GPU-accelerated technique that maps hyperspectral imagery from sensor space to geospatial coordinates using mesh-based interpolation.
• It overcomes computational challenges and pixel artifacts inherent in direct georectification by employing perspective-correct barycentric interpolation.
• Benchmark results show over 100× speedup versus CPU-based methods, enabling real-time interactive analysis and batch export in remote sensing applications.

A fast inverse georectification technique enables hyperspectral remote sensing data to be rapidly and artifact-free mapped from sensor coordinates to geospatially consistent representations. Unlike conventional direct georectification, which is computationally intensive and prone to pixel coverage artifacts due to per-pixel ray tracing, the fast inverse method, as implemented in Spacecube, leverages GPU rasterization and mesh-based interpolation to achieve real-time, high-fidelity georectification suitable for interactive viewing, analysis, and batch export in hyperspectral imaging (Watson et al., 8 Jan 2026).

1. Mathematical Formulation

The georectification process involves transformations between three coordinate systems:

• Image pixel coordinates $(s, \ell)$: $s$ is the sample index (pixel column), $\ell$ is the line index (acquisition time/flight direction).
• Sensor (camera) coordinates $(X_c, Y_c, Z_c)$: with the origin at the optical center and $Z_c$ aligned with the optical axis.
• Geographic (NED) coordinates $(N, E, D)$: North-East-Down metric system derived from UTM and a specified ground height $H$.

Forward Georectification

Given $(s, \ell)$:

1. Normalize pixel positions: $u = \frac{s - c_x}{f_x}$, $v = \frac{\ell - c_y}{f_y}$.
2. Account for lens distortion (radial), derive direction $\mathbf{d}_c$ in camera space.
3. Apply extrinsic rotation $R(\ell)$ and translation $t(\ell)$ for each line.
4. Intersect the resulting world-space ray with the ground plane $D=0$ to obtain ground coordinates $(N, E)$ with $\lambda = -D_c / D_d$.

Inverse Georectification

Mapping $(N, E)$ to $(s, \ell)$ directly would require per-pixel inversion of the forward process, which is computationally prohibitive. Spacecube constructs a 2D mesh covering the imaged footprint for each input line. Each mesh vertex is attributed with both its ground position and its corresponding $(s, \ell)$. Rasterizing this mesh on the GPU yields, for every output pixel, the correct sample and line coordinates via perspective-correct barycentric interpolation, thus obviating the need for explicit per-pixel inversion.

Given triangle vertices $V_0, V_1, V_2$ at positions $X_i$ with sample/line $(s_i, \ell_i)$, each ground pixel $X$ is represented by

$X = \alpha_0 X_0 + \alpha_1 X_1 + \alpha_2 X_2$

and the sensor coordinates recover as

$$(s, \ell) = \alpha_0 (s_0, \ell_0) + \alpha_1 (s_1, \ell_1) + \alpha_2 (s_2, \ell_2),$$

with $(\alpha_0, \alpha_1, \alpha_2)$ the barycentric coordinates of $X$ within the triangle.

2. GPU/OpenGL Implementation

Data Structures

Each mesh vertex stores:

• Projected ground position $(N, E)$
• Depth (vertical distance for perspective correction)
• Sensor coordinates $(s, \ell)$

Input data for georectification includes:

• 2D textures of hyperspectral bands ($N_s \times N_\ell$ per band)
• Calibration textures: dark offsets, radiometric gains, line response

Shaders

Vertex Shader maps vertices from ground coordinates to device coordinates, passing $(s, \ell)$ through to the fragment shader using linear interpolation variables.

Fragment Shader performs:

1. Interpolation (with hardware perspective correction) of $(s, \ell)$.
2. Bilinear sampling of texture data per band.
3. Calibration adjustment: $$L_{b, \ell, s} = (DN_{b, \ell, s} - \text{dark}_{b, s}) \cdot \text{gain}_{b, s} \cdot \text{resp}_\ell$$.

This GPU pipeline enables all mesh triangles to be drawn in a single call per cube, with all output pixel computations executed in parallel and no explicit loops over input lines.

Example Shader Skeletons

Vertex Shader (GLSL-style pseudocode):

#version 330 core layout(location=0) in vec2 a_NE; // world coords layout(location=1) in float aDepth; // used for perspective layout(location=2) in float aSamp; // sample index layout(location=3) in float aLine; // line index out float vSamp; out float vLine; uniform vec2 uCenterNE; uniform vec2 uScaleNE; void main() { vec2 norm = (a_NE - uCenterNE) / uScaleNE; gl_Position = vec4(norm.xy, -aDepth, 1.0); vSamp = aSamp; vLine = aLine; }

Fragment Shader (GLSL-style pseudocode):

#version 330 core in float vSamp; in float vLine; uniform sampler2D uBandTex[3]; uniform sampler2D uDarkTex; uniform sampler2D uRadTex; uniform sampler2D uRespTex; out vec3 outColor; void main() { float t_s = (vSamp + 0.5) / Nsamples; float t_l = (vLine + 0.5) / Nlines; float val[3]; for(int c=0; c<3; ++c) val[c] = texture(uBandTex[c], vec2(t_s,t_l)).r; vec2 dc = texture(uDarkTex, vec2(t_s,0.5)).rg; vec2 rg = texture(uRadTex, vec2(t_s,0.5)).rg; float resp = texture(uRespTex, vec2(t_l,0.5)).r; for(int c=0; c<3; ++c) val[c] = (val[c] - dc[c]) * rg[c] * resp; outColor = vec3(val[0], val[1], val[2]); }

3. Interpolation, Artifacts, and Error Properties

The GPU pipeline relies on two levels of hardware interpolation:

1. Perspective-correct barycentric interpolation of $(s, \ell)$ within mesh triangles ensures that every covered pixel receives exactly one projected sample coordinate, eliminating undersampling and gaps.
2. Bilinear interpolation within texture fetch units for band data provides smoothing and mitigates artifacts across sensor coordinates.

Due to the mesh’s construction and the triangle coverage, every output pixel mapped to the surveyed ground area is hit by exactly one triangle, guaranteeing that the output image is contiguous and free from holes—an advantage over forward-mapping methods that require post-hoc interpolation to fill pixel gaps, or which may introduce over-write artifacts. The only significant geometric error arises from the assumption of local planarity in the ground model, which empirical measurements show to remain below 1% of a sample width even with significant roll/pitch variations (Watson et al., 8 Jan 2026).

4. Performance Evaluation

Benchmarking on a workstation with a Xeon E-2176M and a Quadro P2000 GPU (64GB RAM, data on RAM-disk):

Method	Wall-Clock Time	Real-Time Factor
Spacecube (GPU inv.)	11 s	×1400 real time
Spectronon (CPU dir.)	21 m 29 s	×0.12 real time
PARGE (CPU dir.)	~35 m 48 s	×0.075 real time

For a dataset of 40 cubes ($900 \times 1000 \times 300$ bands), Spacecube processes at approximately 200 Mpixels/s (3 bands) on a midrange GPU, achieving over 100× speedup and more than 10× real-time throughput compared to typical CPU-based direct georectification pipelines (Watson et al., 8 Jan 2026). Doubling the number of bands doubles runtime linearly; quadrupling output pixel count (i.e., halving GSD) increases load as expected but remains GPU-bound. Even on consumer GPUs (e.g., GTX 1660 Ti), throughput sufficient for >1000% of real-time is observed for standard cube sizes.

5. Workflow Pipeline and Practical Usage

The complete pipeline includes:

A. Data ingestion and calibration: - Reading ENVI cubes. - Loading sensor calibration (dark offsets, radiometric gains, line exposure) from metadata. - Converting raw digital numbers to radiance (or reflectance) via $$L_{b, \ell, s} = (DN_{b, \ell, s} - \text{dark}_{b, s}) \cdot \text{gain}_{b, s} \cdot \text{resp}_\ell$$.

B. Georectification setup: - Converting INS/GNSS sensor timestamps and poses to per-line extrinsics $(R(\ell), t(\ell))$ in UTM/NED. - Interpolating these to line start times for precise ground mapping.

C. Ground mesh generation: - For each line $\ell$, project edge sensor rays to ground to define line endpoints. - Connect endpoints as a triangle strip mesh covering the imaged ground region, with each vertex tagged by $(s, \ell)$ and $(N, E)$.

D. Rasterization (Rendering): - Uploading mesh and texture data to GPU once per cube. - Adjusting view, scale, or selected bands interactively or in batch. - Drawing the mesh to georectified output; results are either displayed or saved.

User-Configurable Parameters

Parameter	Effect	Notes
Ground height	Shifts mesh, affects intersection $\lambda$	Errors create parallax in overlaps
GSD	Reduces/increases output resolution	Coarser GSD increases speed
Calibration	Toggles dark/gain multiplication	Small GPU cost for switching calibration
Band selection	Controls number of bands exported	More bands = more passes (linear scaling)
Mipmapping	Enhances performance when zoomed out	Can cause seams at cube borders

6. Algorithmic Advantages and Limitations

This fast inverse georectification pipeline recasts the problem into a memory-bandwidth-limited, fully parallel process. It delivers gap-free georectified outputs, deterministic coverage, and high throughput by exploiting hardware rasterization and interpolation. The primary limitation is the flat-ground assumption, but the geometric errors introduced remain acceptably low for most low-altitude remote sensing operations. A plausible implication is that this paradigm directly supports real-time hyperspectral data acquisition and interactive analytics in agricultural and small-UAS contexts.

7. Conclusion

Spacecube's fast inverse georectification technique demonstrates that recasting georectification as a GPU mesh rasterization problem delivers major improvements over CPU-based direct approaches. Mathematical mappings are unchanged in substance, yet inversion is achieved without brute-force computation, resulting in an operational system that is routinely over 100× faster, runs at >10× real time, and robustly prevents pixel coverage artifacts. The pipeline relies only on a local planar-ground model assumption, whose empirical errors are negligible in targeted applications (Watson et al., 8 Jan 2026).