Instantly share code, notes, and snippets.

# giuliano-macedo/reverse-jhs.md

Last active July 15, 2024 16:23
Show Gist options
• Save giuliano-macedo/006349278c4b6aaea13207291126182d to your computer and use it in GitHub Desktop.
More comprehensive reverse JHS formula for Transverse Mercator projection

# JHS Reverse formula

More comprehensive reverse JHS formula for computing latitude,longitude given northing, easting of a projection.

reference: Geomatics guidance Note number 7 part 2 - 3.5.3.1

(EPSG Dataset coordinate operation method code 9807)

## Constants

spheroid constants:

• $\text{f}$ flattening
• $a$ semi-major axis
• $e^2= 2\text{f} - \text{f}^2$ eccentricity

projection constants:

• $\varphi_0$ Latitude of orgin
• $\lambda_0$ Longitude of origin (sometimes called central meridian)
• $k_0$ scale factor
• $FE$ false easting
• $FN$ false northing

## Derived constants

$n=\text{f}/(2-\text{f})$

$B = [\frac{a}{1+n}](1+\frac{n^2}{4}+\frac{n^4}{64})$

• $h_1 = \frac{n}{2} - \frac{2}{3}n^2 + \frac{5}{16}n^3 + \frac{41}{180}n^4$

• $h_2 = \frac{13}{48}n^2 - \frac{3}{5}n^3 + \frac{557}{1440}n^4$

• $h_3 = \frac{61}{240}n^3 - \frac{103}{140}n^4$

• $h_4 = \frac{49561}{161280}n^4$

Obs: $h_1$ - $h_4$ are only used in the reverse formulas for calculating $M_0$ when $\varphi_o \notin \set{0, \frac{\pi}{2}, -\frac{\pi}{2}}$

• $h'_1 = \frac{n}{2} - \frac{2}{3} n^2 + \frac{37}{96} n^3 - \frac{1}{360} n^4$

• $h'_2 = \frac{1}{48} n^2 + \frac{1}{15} n^3 - \frac{437}{1440} n^4$

• $h'_3 = \frac{17}{480} n^3 - \frac{37}{840} n^4$

• $h'_4 = \frac{4397}{161280} n^4$

### $M_0$

meridional arc distance from equator to the projection origin $M_0$ is computed from:

if $\varphi_0=0$ then $M_0=0$

else if $\varphi_0=90^{\circ}\text{N}=\frac{\pi}{2}\text{ radians}$ then $M_0=B\frac{\pi}{2}$

else if $\varphi_0=90^{\circ}\text{S}=-\frac{\pi}{2}\text{ radians}$ then $M_0=-B\frac{\pi}{2}$

else:

$Q_0 = \text{asinh}(\tan{\varphi_0}) - [e\text{ atanh}(e\sin{\varphi_0})]$

$\beta_O=\text{atan}(\text{sinh } Q_0)$

$\xi_{O0}=\text{asin}(\sin B_0)$

• $\xi_{O1} = h_1 \sin(2\xi_{O0})$

• $\xi_{O2} = h_2 \sin(4\xi_{O0})$

• $\xi_{O3} = h_3 \sin(6\xi_{O0})$

• $\xi_{O4} = h_4 \sin(8\xi_{O0})$

$\xi_O= \xi_{O1}+\xi_{O2}+\xi_{O3}+\xi_{O4}$

$M_0=B\xi_O$

end

## Formula

Given easting $E$ and northing $N$ in meters, and $\varphi$, $\lambda$ as resulting longitute, latitude in radians:

$\eta' = \frac {E-FE}{B k_0}$

$\xi' = \frac { (N-FN) + k_0 M_0 }{B k_0}$

• $\xi'_1 = h'_1 \text{sin}(2\xi')\text{ cosh}(2\eta')$
• $\xi'_2 = h'_2 \text{sin}(4\xi')\text{ cosh}(4\eta')$
• $\xi'_3 = h'_3 \text{sin}(6\xi')\text{ cosh}(6\eta')$
• $\xi'_4 = h'_4 \text{sin}(8\xi')\text{ cosh}(8\eta')$
• $\xi'_0 = \xi' - (\xi'_1 + \xi'_2 + \xi'_3 + \xi'_4)$

• $\eta'_1 = h'_1 \text{cos}(2\xi')\text{ sinh}(2\eta')$
• $\eta'_2 = h'_2 \text{cos}(4\xi')\text{ sinh}(4\eta')$
• $\eta'_3 = h'_3 \text{cos}(6\xi')\text{ sinh}(6\eta')$
• $\eta'_4 = h'_4 \text{cos}(8\xi')\text{ sinh}(8\eta')$
• $\eta'_0 = \eta' - (\eta'_1 + \eta'_2 + \eta'_3 + \eta'_4)$

$\beta' = \text{asin}( \frac{\sin{\xi'_0}}{\text{cosh }\eta'_0} )$

$Q' = \text{asinh }(\tan{\beta'})$

$Q'' = Q' + [e \text{ atanh}(e \text{ tanh }Q')]$

$Q'' = Q' + [e \text{ atanh}(e \text{ tanh }Q'')]$ which should be iterated until the change in Q" is insignificant

Finally:

$\varphi = \text{atan }(\text{sinh } Q'')$

$\lambda = \lambda_0 + \text{asin }(\frac{\text{tanh}(\eta'_0)} {\cos \beta'} )$