Passed
Push — master ( 1132b4...5d5536 )
by Doug
09:03
created

Point::vincenty()   B

Complexity

Conditions 7
Paths 6

Size

Total Lines 62
Code Lines 46

Duplication

Lines 0
Ratio 0 %

Code Coverage

Tests 43
CRAP Score 7.0005

Importance

Changes 0
Metric Value
cc 7
eloc 46
c 0
b 0
f 0
nc 6
nop 3
dl 0
loc 62
ccs 43
cts 44
cp 0.9773
crap 7.0005
rs 8.2448

How to fix   Long Method   

Long Method

Small methods make your code easier to understand, in particular if combined with a good name. Besides, if your method is small, finding a good name is usually much easier.

For example, if you find yourself adding comments to a method's body, this is usually a good sign to extract the commented part to a new method, and use the comment as a starting point when coming up with a good name for this new method.

Commonly applied refactorings include:

1
<?php
2
/**
3
 * PHPCoord.
4
 *
5
 * @author Doug Wright
6
 */
7
declare(strict_types=1);
8
9
namespace PHPCoord;
10
11
use function abs;
12
use function acos;
13
use function array_reverse;
14
use function array_values;
15
use function asin;
16
use function atan;
17
use function atan2;
18
use function cos;
19
use DateTimeImmutable;
20
use function in_array;
21
use function lcfirst;
22
use const M_PI;
23
use const PHP_MAJOR_VERSION;
24
use PHPCoord\CoordinateOperation\CoordinateOperationMethods;
25
use PHPCoord\CoordinateOperation\CoordinateOperationParams;
0 ignored issues
show
Bug introduced by
The type PHPCoord\CoordinateOpera...ordinateOperationParams 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...
26
use PHPCoord\CoordinateOperation\CoordinateOperations;
0 ignored issues
show
Bug introduced by
The type PHPCoord\CoordinateOperation\CoordinateOperations 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...
27
use PHPCoord\CoordinateOperation\GeographicValue;
28
use PHPCoord\CoordinateReferenceSystem\CoordinateReferenceSystem;
29
use PHPCoord\CoordinateSystem\Axis;
30
use PHPCoord\Datum\Ellipsoid;
31
use PHPCoord\UnitOfMeasure\Angle\Angle;
32
use PHPCoord\UnitOfMeasure\Length\Length;
33
use PHPCoord\UnitOfMeasure\Length\Metre;
34
use PHPCoord\UnitOfMeasure\Scale\Coefficient;
35
use PHPCoord\UnitOfMeasure\Scale\Scale;
36
use PHPCoord\UnitOfMeasure\UnitOfMeasure;
37
use PHPCoord\UnitOfMeasure\UnitOfMeasureFactory;
38
use function preg_match;
39
use function sin;
40
use function sqrt;
41
use function sscanf;
42
use function str_replace;
43
use Stringable;
44
use function tan;
45
use function ucwords;
46
47
abstract class Point implements Stringable
48
{
49
    protected const ITERATION_CONVERGENCE = 1e-12;
50
51
    /**
52
     * @internal
53
     */
54 11
    public function performOperation(string $srid, CoordinateReferenceSystem $to, bool $inReverse): self
55
    {
56 11
        $operations = self::resolveConcatenatedOperations($srid, $inReverse);
57
58 11
        $point = $this;
59 11
        foreach ($operations as $operationSrid => $operation) {
60 11
            $method = CoordinateOperationMethods::getFunctionName($operation['method']);
61 11
            if (isset($operation['source_crs'])) {
62
                $destCRS = CoordinateReferenceSystem::fromSRID($inReverse ? $operation['source_crs'] : $operation['target_crs']);
63
            } else {
64 11
                $destCRS = $to;
65
            }
66
67 11
            $params = self::resolveParamsByOperation($operationSrid, $operation['method'], $inReverse);
68
69 11
            if (PHP_MAJOR_VERSION >= 8) {
70
                $point = $point->$method($destCRS, ...$params);
71
            } else {
72 11
                $point = $point->$method($destCRS, ...array_values($params));
73
            }
74
        }
75
76 11
        $point->crs = $to; //some operations are reused across CRSses (e.g. ETRS89 and WGS84), so the $destCRS of the final suboperation might not be the intended target
77
78 11
        return $point;
79
    }
80
81 11
    protected static function camelCase(string $string): string
82
    {
83 11
        $string = str_replace([' ', '-'], '', ucwords($string, ' -'));
84 11
        if (!preg_match('/[ABC][uv\d]/', $string)) {
85 11
            $string = lcfirst($string);
86
        }
87
88 11
        return $string;
89
    }
90
91 11
    protected static function resolveConcatenatedOperations(string $operationSrid, bool $inReverse): array
92
    {
93 11
        $operations = [];
94 11
        $operation = CoordinateOperations::getOperationData($operationSrid);
95 11
        if (isset($operation['operations'])) {
96
            foreach ($operation['operations'] as $subOperation) {
97
                $subOperationData = CoordinateOperations::getOperationData($subOperation['operation']);
98
                $subOperationData['source_crs'] = $subOperation['source_crs'];
99
                $subOperationData['target_crs'] = $subOperation['target_crs'];
100
                $operations[$subOperation['operation']] = $subOperationData;
101
            }
102
        } else {
103 11
            $operations[$operationSrid] = $operation;
104
        }
105
106 11
        if ($inReverse) {
107 2
            $operations = array_reverse($operations, true);
108
        }
109
110 11
        return $operations;
111
    }
112
113 11
    protected static function resolveParamsByOperation(string $operationSrid, string $methodSrid, bool $inReverse): array
114
    {
115 11
        $params = [];
116 11
        $powerCoefficients = [];
117 11
        foreach (CoordinateOperationParams::getParamData($operationSrid) as $paramName => $paramData) {
118 11
            $value = $paramData['value'];
119 11
            if ($inReverse && $paramData['reverses']) {
120 1
                $value *= -1;
121
            }
122 11
            if ($paramData['uom']) {
123 11
                $param = UnitOfMeasureFactory::makeUnit($value, $paramData['uom']);
124
            } else {
125
                $param = $paramData['value'];
126
            }
127 11
            $paramName = static::camelCase($paramName);
128 11
            if (preg_match('/^(Au|Bu)/', $paramName)) {
129 3
                $powerCoefficients[$paramName] = $param;
130
            } else {
131 11
                $params[$paramName] = $param;
132
            }
133
        }
134 11
        if ($powerCoefficients) {
135 3
            $params['powerCoefficients'] = $powerCoefficients;
136
        }
137 11
        if (in_array($methodSrid, [
138
            CoordinateOperationMethods::EPSG_SIMILARITY_TRANSFORMATION,
139
            CoordinateOperationMethods::EPSG_AFFINE_PARAMETRIC_TRANSFORMATION,
140
        ], true)) {
141
            $params['inReverse'] = $inReverse;
142
        }
143
144 11
        return $params;
145
    }
146
147 213
    protected function getAxisByName(string $name): ?Axis
148
    {
149 213
        foreach ($this->getCRS()->getCoordinateSystem()->getAxes() as $axis) {
150 213
            if ($axis->getName() === $name) {
151 213
                return $axis;
152
            }
153
        }
154
155 135
        return null;
156
    }
157
158 10
    protected static function sign(float $number): int
159
    {
160 10
        if ($number < 0) {
161
            return -1;
162
        }
163
164 10
        return 1;
165
    }
166
167
    /**
168
     * Calculate surface distance between two points.
169
     */
170 14
    protected static function vincenty(GeographicValue $from, GeographicValue $to, Ellipsoid $ellipsoid): Length
171
    {
172 14
        $a = $ellipsoid->getSemiMajorAxis()->asMetres()->getValue();
173 14
        $b = $ellipsoid->getSemiMinorAxis()->asMetres()->getValue();
174 14
        $f = $ellipsoid->getInverseFlattening();
175 14
        $U1 = atan((1 - $f) * tan($from->getLatitude()->asRadians()->getValue()));
176 14
        $U2 = atan((1 - $f) * tan($to->getLatitude()->asRadians()->getValue()));
177 14
        $L = $to->getLongitude()->subtract($from->getLongitude())->asRadians()->getValue();
178
179 14
        $lambda = $L;
180
        do {
181 14
            $lambdaN = $lambda;
182
183 14
            $sinSigma = sqrt((cos($U2) * sin($lambda)) ** 2 + (cos($U1) * sin($U2) - sin($U1) * cos($U2) * cos($lambda)) ** 2);
184 14
            $cosSigma = sin($U1) * sin($U2) + cos($U1) * cos($U2) * cos($lambda);
185 14
            $sigma = atan2($sinSigma, $cosSigma);
186
187 14
            $sinAlpha = cos($U1) * cos($U2) * sin($lambda) / $sinSigma;
188 14
            $cosSqAlpha = (1 - $sinAlpha ** 2);
189 14
            $cos2SigmaM = $cosSqAlpha ? $cosSigma - (2 * sin($U1) * sin($U2) / $cosSqAlpha) : 0;
190 14
            $C = $f / 16 * $cosSqAlpha * (4 + $f * (4 - 3 * $cosSqAlpha));
191 14
            $lambda = $L + (1 - $C) * $f * $sinAlpha * ($sigma + $C * $sinSigma * ($cos2SigmaM + $C * $cosSigma * (-1 + 2 * $cos2SigmaM ** 2)));
192 14
        } while (abs($lambda - $lambdaN) >= static::ITERATION_CONVERGENCE && abs($lambda) < M_PI);
193
194
        // Antipodal case
195 14
        if (abs($lambda) >= M_PI) {
196 3
            if ($L >= 0) {
197 3
                $LPrime = M_PI - $L;
198
            } else {
199
                $LPrime = -M_PI - $L;
200
            }
201
202 3
            $lambdaPrime = 0;
203 3
            $sigma = M_PI - abs($U1 + $U2);
204 3
            $sinSigma = sin($sigma);
205 3
            $cosSqAlpha = 0.5;
206 3
            $sinAlpha = 0;
207
208
            do {
209 3
                $sinAlphaN = $sinAlpha;
210
211 3
                $C = $f / 16 * $cosSqAlpha * (4 + $f * (4 - 3 * $cosSqAlpha));
212 3
                $cos2SigmaM = cos($sigma) - 2 * sin($U1) * sin($U2) / $cosSqAlpha;
213 3
                $D = (1 - $C) * $f * ($sigma + $C * $sinSigma * ($cos2SigmaM + $C * cos($sigma) * (-1 + 2 * $cos2SigmaM ** 2)));
214 3
                $sinAlpha = ($LPrime - $lambdaPrime) / $D;
215 3
                $cosSqAlpha = (1 - $sinAlpha ** 2);
216 3
                $sinLambdaPrime = ($sinAlpha * $sinSigma) / (cos($U1) * cos($U2));
217 3
                $lambdaPrime = self::asin($sinLambdaPrime);
218 3
                $sinSqSigma = (cos($U2) * $sinLambdaPrime) ** 2 + (cos($U1) * sin($U2) + sin($U1) * cos($U2) * cos($lambdaPrime)) ** 2;
219 3
                $sinSigma = sqrt($sinSqSigma);
220 3
            } while (abs($sinAlpha - $sinAlphaN) >= static::ITERATION_CONVERGENCE);
221
        }
222
223 14
        $E = sqrt(1 + (($a ** 2 - $b ** 2) / $b ** 2) * $cosSqAlpha);
224 14
        $F = ($E - 1) / ($E + 1);
225
226 14
        $A = (1 + $F ** 2 / 4) / (1 - $F);
227 14
        $B = $F * (1 - 3 / 8 * $F ** 2);
228
229 14
        $deltaSigma = $B * $sinSigma * ($cos2SigmaM + $B / 4 * ($cosSigma * (-1 + 2 * $cos2SigmaM ** 2) - $B / 6 * $cos2SigmaM * (-3 + 4 * $sinSigma ** 2) * (-3 + 4 * $cos2SigmaM ** 2)));
230
231 14
        return new Metre($b * $A * ($sigma - $deltaSigma));
232
    }
233
234
    /**
235
     * General polynomial.
236
     * @param Coefficient[] $powerCoefficients
237
     */
238 2
    protected function generalPolynomialUnitless(
239
        float $xs,
240
        float $ys,
241
        UnitOfMeasure $ordinate1OfEvaluationPointInSourceCRS,
242
        UnitOfMeasure $ordinate2OfEvaluationPointInSourceCRS,
243
        UnitOfMeasure $ordinate1OfEvaluationPointInTargetCRS,
244
        UnitOfMeasure $ordinate2OfEvaluationPointInTargetCRS,
245
        Scale $scalingFactorForSourceCRSCoordDifferences,
246
        Scale $scalingFactorForTargetCRSCoordDifferences,
247
        Scale $A0,
248
        Scale $B0,
249
        array $powerCoefficients
250
    ): array {
251 2
        $xso = $ordinate1OfEvaluationPointInSourceCRS->getValue();
252 2
        $yso = $ordinate2OfEvaluationPointInSourceCRS->getValue();
253 2
        $xto = $ordinate1OfEvaluationPointInTargetCRS->getValue();
254 2
        $yto = $ordinate2OfEvaluationPointInTargetCRS->getValue();
255
256 2
        $U = $scalingFactorForSourceCRSCoordDifferences->asUnity()->getValue() * ($xs - $xso);
257 2
        $V = $scalingFactorForSourceCRSCoordDifferences->asUnity()->getValue() * ($ys - $yso);
258
259 2
        $mTdX = $A0->getValue();
260 2
        foreach ($powerCoefficients as $coefficientName => $coefficientValue) {
261 2
            if ($coefficientName[0] === 'A') {
262 2
                sscanf($coefficientName, 'Au%dv%d', $uPower, $vPower);
1 ignored issue
show
Comprehensibility Best Practice introduced by
The variable $vPower seems to be never defined.
Loading history...
263 2
                $mTdX += $coefficientValue->getValue() * $U ** $uPower * $V ** $vPower;
264
            }
265
        }
266
267 2
        $mTdY = $B0->getValue();
268 2
        foreach ($powerCoefficients as $coefficientName => $coefficientValue) {
269 2
            if ($coefficientName[0] === 'B') {
270 2
                sscanf($coefficientName, 'Bu%dv%d', $uPower, $vPower);
271 2
                $mTdY += $coefficientValue->getValue() * $U ** $uPower * $V ** $vPower;
272
            }
273
        }
274
275 2
        $xt = $xs - $xso + $xto + $mTdX / $scalingFactorForTargetCRSCoordDifferences->asUnity()->getValue();
276 2
        $yt = $ys - $yso + $yto + $mTdY / $scalingFactorForTargetCRSCoordDifferences->asUnity()->getValue();
277
278 2
        return ['xt' => $xt, 'yt' => $yt];
279
    }
280
281
    /**
282
     * Reversible polynomial.
283
     */
284 4
    protected function reversiblePolynomialUnitless(
285
        float $xs,
286
        float $ys,
287
        Angle $ordinate1OfEvaluationPoint,
288
        Angle $ordinate2OfEvaluationPoint,
289
        Scale $scalingFactorForCoordDifferences,
290
        Scale $A0,
291
        Scale $B0,
292
        array $powerCoefficients
293
    ): array {
294 4
        $xo = $ordinate1OfEvaluationPoint->getValue();
295 4
        $yo = $ordinate2OfEvaluationPoint->getValue();
296
297 4
        $U = $scalingFactorForCoordDifferences->asUnity()->getValue() * ($xs - $xo);
298 4
        $V = $scalingFactorForCoordDifferences->asUnity()->getValue() * ($ys - $yo);
299
300 4
        $mTdX = $A0->getValue();
301 4
        foreach ($powerCoefficients as $coefficientName => $coefficientValue) {
302 4
            if ($coefficientName[0] === 'A') {
303 4
                sscanf($coefficientName, 'Au%dv%d', $uPower, $vPower);
1 ignored issue
show
Comprehensibility Best Practice introduced by
The variable $vPower seems to be never defined.
Loading history...
304 4
                $mTdX += $coefficientValue->getValue() * $U ** $uPower * $V ** $vPower;
305
            }
306
        }
307
308 4
        $mTdY = $B0->getValue();
309 4
        foreach ($powerCoefficients as $coefficientName => $coefficientValue) {
310 4
            if ($coefficientName[0] === 'B') {
311 4
                sscanf($coefficientName, 'Bu%dv%d', $uPower, $vPower);
312 4
                $mTdY += $coefficientValue->getValue() * $U ** $uPower * $V ** $vPower;
313
            }
314
        }
315
316 4
        $xt = $xs + $mTdX * $scalingFactorForCoordDifferences->asUnity()->getValue();
317 4
        $yt = $ys + $mTdY * $scalingFactorForCoordDifferences->asUnity()->getValue();
318
319 4
        return ['xt' => $xt, 'yt' => $yt];
320
    }
321
322
    /**
323
     * Floating point vagaries mean that it's possible for inputs to be e.g. 1.00000000000001 which makes PHP give a
324
     * silent NaN as output so inputs need to be capped. atan/atan2 are not affected, they seem to cap internally.
325
     */
326
    protected static function acos(float $num): float
327
    {
328
        if ($num > 1.0) {
329
            $num = 1.0;
330
        } elseif ($num < -1) {
331
            $num = -1.0;
332
        }
333
334
        return acos($num);
335
    }
336
337
    /**
338
     * Floating point vagaries mean that it's possible for inputs to be e.g. 1.00000000000001 which makes PHP give a
339
     * silent NaN as output so inputs need to be capped. atan/atan2 are not affected, they seem to cap internally.
340
     */
341 42
    protected static function asin(float $num): float
342
    {
343 42
        if ($num > 1.0) {
344
            $num = 1.0;
345 42
        } elseif ($num < -1.0) {
346
            $num = -1.0;
347
        }
348
349 42
        return asin($num);
350
    }
351
352
    abstract public function getCRS(): CoordinateReferenceSystem;
353
354
    abstract public function getCoordinateEpoch(): ?DateTimeImmutable;
355
356
    abstract public function calculateDistance(self $to): Length;
357
}
358