Passed
Push — 4.0.x ( 414ad0...5c8c8d )
by Doug
05:04 queued 12s
created

Point::sign()   A

Complexity

Conditions 2
Paths 2

Size

Total Lines 7
Code Lines 3

Duplication

Lines 0
Ratio 0 %

Code Coverage

Tests 3
CRAP Score 2.0625

Importance

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