Passed
Push — master ( a784d9...14494d )
by Doug
60:49
created

OSTNOSGM15Grid::applyForwardAdjustment()   A

Complexity

Conditions 1
Paths 1

Size

Total Lines 8
Code Lines 4

Duplication

Lines 0
Ratio 0 %

Code Coverage

Tests 0
CRAP Score 2

Importance

Changes 1
Bugs 0 Features 0
Metric Value
cc 1
eloc 4
c 1
b 0
f 0
nc 1
nop 1
dl 0
loc 8
ccs 0
cts 5
cp 0
crap 2
rs 10
1
<?php
2
/**
3
 * PHPCoord.
4
 *
5
 * @author Doug Wright
6
 */
7
declare(strict_types=1);
8
9
namespace PHPCoord\CoordinateOperation;
10
11
use function abs;
12
use const PHP_MAJOR_VERSION;
13
use PHPCoord\CoordinateReferenceSystem\Projected;
0 ignored issues
show
Bug introduced by
The type PHPCoord\CoordinateReferenceSystem\Projected was not found. Maybe you did not declare it correctly or list all dependencies?

The issue could also be caused by a filter entry in the build configuration. If the path has been excluded in your configuration, e.g. excluded_paths: ["lib/*"], you can move it to the dependency path list as follows:

filter:
    dependency_paths: ["lib/*"]

For further information see https://scrutinizer-ci.com/docs/tools/php/php-scrutinizer/#list-dependency-paths

Loading history...
14
use PHPCoord\CoordinateSystem\Cartesian;
15
use PHPCoord\Datum\Datum;
16
use PHPCoord\ProjectedPoint;
17
use PHPCoord\UnitOfMeasure\Length\Metre;
18
use SplFileObject;
19
20
class OSTNOSGM15Grid extends SplFileObject
21
{
22
    private const GRID_SIZE = 1000;
23
24
    private const ITERATION_CONVERGENCE = 0.0001;
25
26
    public function applyForwardHorizontalAdjustment(ProjectedPoint $point): ProjectedPoint
27
    {
28
        $adjustment = $this->getAdjustment($point->getEasting()->asMetres(), $point->getNorthing()->asMetres());
29
30
        $easting = $point->getEasting()->add(new Metre($adjustment[0]));
31
        $northing = $point->getNorthing()->add(new Metre($adjustment[1]));
32
33
        return ProjectedPoint::createFromEastingNorthing($easting, $northing, Projected::fromSRID(Projected::EPSG_OSGB36_BRITISH_NATIONAL_GRID), $point->getCoordinateEpoch());
34
    }
35
36
    public function applyReverseHorizontalAdjustment(ProjectedPoint $point): ProjectedPoint
37
    {
38
        $osgb36NationalGrid = Projected::fromSRID(Projected::EPSG_OSGB36_BRITISH_NATIONAL_GRID);
39
        $etrs89NationalGrid = new Projected(
40
            'ETRS89 / National Grid',
41
            Cartesian::fromSRID(Cartesian::EPSG_2D_AXES_EASTING_NORTHING_E_N_ORIENTATIONS_EAST_NORTH_UOM_M),
42
            Datum::fromSRID(Datum::EPSG_EUROPEAN_TERRESTRIAL_REFERENCE_SYSTEM_1989_ENSEMBLE),
43
            $osgb36NationalGrid->getBoundingArea()
44
        );
45
46
        $adjustment = [0, 0];
47
        $easting = $point->getEasting();
48
        $northing = $point->getNorthing();
49
50
        do {
51
            $prevAdjustment = $adjustment;
52
            $adjustment = $this->getAdjustment($easting->asMetres(), $northing->asMetres());
53
            $easting = $point->getEasting()->subtract(new Metre($adjustment[0]));
54
            $northing = $point->getNorthing()->subtract(new Metre($adjustment[1]));
55
        } while (abs($adjustment[0] - $prevAdjustment[0]) > self::ITERATION_CONVERGENCE && abs($adjustment[1] - $prevAdjustment[1]) > self::ITERATION_CONVERGENCE);
56
57
        return ProjectedPoint::createFromEastingNorthing($easting, $northing, $etrs89NationalGrid, $point->getCoordinateEpoch());
58
    }
59
60
    public function getVerticalAdjustment(ProjectedPoint $point): Metre
61
    {
62
        $adjustment = $this->getAdjustment($point->getEasting()->asMetres(), $point->getNorthing()->asMetres());
63
64
        return new Metre($adjustment[2]);
65
    }
66
67
    private function getAdjustment(Metre $easting, Metre $northing): array
68
    {
69
        $eastIndex = (int) ($easting->getValue() / self::GRID_SIZE);
70
        $northIndex = (int) ($northing->getValue() / self::GRID_SIZE);
71
72
        $corner0 = $this->getRecord($eastIndex, $northIndex);
73
        $corner1 = $this->getRecord($eastIndex + 1, $northIndex);
74
        $corner2 = $this->getRecord($eastIndex + 1, $northIndex + 1);
75
        $corner3 = $this->getRecord($eastIndex, $northIndex + 1);
76
77
        $dx = $easting->getValue() - $corner0[1];
78
        $dy = $northing->getValue() - $corner0[2];
79
80
        $t = $dx / self::GRID_SIZE;
81
        $u = $dy / self::GRID_SIZE;
82
83
        $se = (1 - $t) * (1 - $u) * $corner0[3] + ($t) * (1 - $u) * $corner1[3] + ($t) * ($u) * $corner2[3] + (1 - $t) * ($u) * $corner3[3];
84
        $sn = (1 - $t) * (1 - $u) * $corner0[4] + ($t) * (1 - $u) * $corner1[4] + ($t) * ($u) * $corner2[4] + (1 - $t) * ($u) * $corner3[4];
85
        $sh = (1 - $t) * (1 - $u) * $corner0[5] + ($t) * (1 - $u) * $corner1[5] + ($t) * ($u) * $corner2[5] + (1 - $t) * ($u) * $corner3[5];
86
87
        return [$se, $sn, $sh];
88
    }
89
90
    private function getRecord(int $eastIndex, int $northIndex): array
91
    {
92
        $record = $northIndex * 701 + $eastIndex + 1;
93
94
        // https://bugs.php.net/bug.php?id=62004
95
        if (PHP_MAJOR_VERSION < 8) {
96
            --$record;
97
        }
98
99
        $this->seek($record);
100
101
        return $this->fgetcsv();
102
    }
103
}
104