Point::asin() - Code Metrics - Inspection of "Add extents to code coverage" - dvdoug/PHPCoord - Measure and Improve Code Quality continuously with Scrutinizer

Passed

Push — extents ( f35511...ba1d5f )

by Doug

created 2021-10-09 21:25 UTC

Point::asin() A

↳ Parent: Point

Complexity

Conditions	3
Paths	3

Size

Total Lines	9
Code Lines	5

Duplication

Lines	0
Ratio	0 %

Code Coverage

Tests	7
CRAP Score	3

Importance

Changes

Metric	Value
cc	3
eloc	5
nc	3
nop	1
dl	0
loc	9
ccs	7
cts	7
cp	1
crap	3
rs	10
c	0
b	0
f	0

<?php
/**
 * PHPCoord.
 *
 * @author Doug Wright
 */
declare(strict_types=1);

namespace PHPCoord;

use function abs;
use function acos;
use function asin;
use function atan;
use function atan2;
use function cos;
use DateTimeImmutable;
use const M_PI;
use PHPCoord\CoordinateOperation\CoordinateOperationMethods;
use PHPCoord\CoordinateOperation\CoordinateOperationParams;
filter:
    dependency_paths: ["lib/*"]
use PHPCoord\CoordinateOperation\CoordinateOperations;
use PHPCoord\CoordinateOperation\GeographicValue;
use PHPCoord\CoordinateReferenceSystem\CoordinateReferenceSystem;
use PHPCoord\Datum\Ellipsoid;
use PHPCoord\UnitOfMeasure\Angle\Angle;
use PHPCoord\UnitOfMeasure\Length\Length;
use PHPCoord\UnitOfMeasure\Length\Metre;
use PHPCoord\UnitOfMeasure\Scale\Coefficient;
use PHPCoord\UnitOfMeasure\Scale\Scale;
use PHPCoord\UnitOfMeasure\UnitOfMeasure;
use PHPCoord\UnitOfMeasure\UnitOfMeasureFactory;
use function sin;
use function sqrt;
use function sscanf;
use function str_starts_with;
use Stringable;
use function tan;

abstract class Point implements Stringable
{
    protected const ITERATION_CONVERGENCE_FORMULA = 1e-10;
    protected const ITERATION_CONVERGENCE_GRID = 0.0001;
    protected const METHODS_REQUIRING_HORIZONTAL_POINT = [
        CoordinateOperationMethods::EPSG_VERTICAL_OFFSET_AND_SLOPE => CoordinateOperationMethods::EPSG_VERTICAL_OFFSET_AND_SLOPE,
        CoordinateOperationMethods::EPSG_ZERO_TIDE_HEIGHT_TO_MEAN_TIDE_HEIGHT_EVRF2019 => CoordinateOperationMethods::EPSG_ZERO_TIDE_HEIGHT_TO_MEAN_TIDE_HEIGHT_EVRF2019,
        CoordinateOperationMethods::EPSG_VERTICAL_OFFSET_BY_GRID_INTERPOLATION_GTX => CoordinateOperationMethods::EPSG_VERTICAL_OFFSET_BY_GRID_INTERPOLATION_GTX,
    ];
    protected const METHODS_THAT_REQUIRE_DIRECTION = [
        CoordinateOperationMethods::EPSG_SIMILARITY_TRANSFORMATION => CoordinateOperationMethods::EPSG_SIMILARITY_TRANSFORMATION,
        CoordinateOperationMethods::EPSG_AFFINE_PARAMETRIC_TRANSFORMATION => CoordinateOperationMethods::EPSG_AFFINE_PARAMETRIC_TRANSFORMATION,
        CoordinateOperationMethods::EPSG_NADCON5_2D => CoordinateOperationMethods::EPSG_NADCON5_2D,
        CoordinateOperationMethods::EPSG_NADCON5_3D => CoordinateOperationMethods::EPSG_NADCON5_3D,
        CoordinateOperationMethods::EPSG_NTV2 => CoordinateOperationMethods::EPSG_NTV2,
        CoordinateOperationMethods::EPSG_ZERO_TIDE_HEIGHT_TO_MEAN_TIDE_HEIGHT_EVRF2019 => CoordinateOperationMethods::EPSG_ZERO_TIDE_HEIGHT_TO_MEAN_TIDE_HEIGHT_EVRF2019,
        CoordinateOperationMethods::EPSG_GEOCENTRIC_TRANSLATION_BY_GRID_INTERPOLATION_IGN => CoordinateOperationMethods::EPSG_GEOCENTRIC_TRANSLATION_BY_GRID_INTERPOLATION_IGN,
        CoordinateOperationMethods::EPSG_VERTICAL_OFFSET_BY_GRID_INTERPOLATION_GTX => CoordinateOperationMethods::EPSG_VERTICAL_OFFSET_BY_GRID_INTERPOLATION_GTX,
    ];

    /**
     * @internal
     */
    public function performOperation(string $srid, CoordinateReferenceSystem $to, bool $inReverse, array $additionalParams = []): self
    {
        $operation = CoordinateOperations::getOperationData($srid);
        $method = CoordinateOperationMethods::getFunctionName($operation['method']);
        $params = self::resolveParamsByOperation($srid, $operation['method'], $inReverse);

        if (isset(self::METHODS_REQUIRING_HORIZONTAL_POINT[$operation['method']])) {
            $params['horizontalPoint'] = $additionalParams['horizontalPoint'];
        }

        return $this->$method($to, ...$params);
    }

    protected static function resolveParamsByOperation(string $operationSrid, string $methodSrid, bool $inReverse): array
    {
        $params = [];
        $powerCoefficients = [];
        foreach (CoordinateOperationParams::getParamData($operationSrid) as $paramName => $paramData) {
            if (isset($paramData['fileProvider'])) {
                $params[$paramName] = (new $paramData['fileProvider']())->provideGrid();
            } else {
                if ($inReverse && $paramData['reverses']) {
                    $paramData['value'] *= -1;
                }
                if ($paramData['uom']) {
                    $param = UnitOfMeasureFactory::makeUnit($paramData['value'], $paramData['uom']);
                } else {
                    $param = $paramData['value'];
                }
                if (str_starts_with($paramName, 'Au') || str_starts_with($paramName, 'Bu')) {
                    $powerCoefficients[$paramName] = $param;
                } else {
                    $params[$paramName] = $param;
                }
            }
        }
        if ($powerCoefficients) {
            $params['powerCoefficients'] = $powerCoefficients;
        }
        if (isset(self::METHODS_THAT_REQUIRE_DIRECTION[$methodSrid])) {
            $params['inReverse'] = $inReverse;
        }

        return $params;
    }

    protected static function sign(float $number): int
    {
        if ($number < 0) {
            return -1;
        }

        return 1;
    }

    /**
     * Calculate surface distance between two points.
     */
    protected static function vincenty(GeographicValue $from, GeographicValue $to, Ellipsoid $ellipsoid): Length
    {
        $a = $ellipsoid->getSemiMajorAxis()->asMetres()->getValue();
        $b = $ellipsoid->getSemiMinorAxis()->asMetres()->getValue();
        $f = $ellipsoid->getInverseFlattening();
        $U1 = atan((1 - $f) * tan($from->getLatitude()->asRadians()->getValue()));
        $U2 = atan((1 - $f) * tan($to->getLatitude()->asRadians()->getValue()));
        $L = $to->getLongitude()->subtract($from->getLongitude())->asRadians()->getValue();

        $lambda = $L;
        do {
            $lambdaN = $lambda;

            $sinSigma = sqrt((cos($U2) * sin($lambda)) ** 2 + (cos($U1) * sin($U2) - sin($U1) * cos($U2) * cos($lambda)) ** 2);
            $cosSigma = sin($U1) * sin($U2) + cos($U1) * cos($U2) * cos($lambda);
            $sigma = atan2($sinSigma, $cosSigma);

            $sinAlpha = cos($U1) * cos($U2) * sin($lambda) / $sinSigma;
            $cosSqAlpha = (1 - $sinAlpha ** 2);
            $cos2SigmaM = $cosSqAlpha ? $cosSigma - (2 * sin($U1) * sin($U2) / $cosSqAlpha) : 0;
            $C = $f / 16 * $cosSqAlpha * (4 + $f * (4 - 3 * $cosSqAlpha));
            $lambda = $L + (1 - $C) * $f * $sinAlpha * ($sigma + $C * $sinSigma * ($cos2SigmaM + $C * $cosSigma * (-1 + 2 * $cos2SigmaM ** 2)));
        } while (abs($lambda - $lambdaN) >= static::ITERATION_CONVERGENCE_FORMULA && abs($lambda) < M_PI);

        // Antipodal case
        if (abs($lambda) >= M_PI) {
            if ($L >= 0) {
                $LPrime = M_PI - $L;
            } else {
                $LPrime = -M_PI - $L;
            }

            $lambdaPrime = 0;
            $sigma = M_PI - abs($U1 + $U2);
            $sinSigma = sin($sigma);
            $cosSqAlpha = 0.5;
            $sinAlpha = 0;

            do {
                $sinAlphaN = $sinAlpha;

                $C = $f / 16 * $cosSqAlpha * (4 + $f * (4 - 3 * $cosSqAlpha));
                $cos2SigmaM = cos($sigma) - 2 * sin($U1) * sin($U2) / $cosSqAlpha;
                $D = (1 - $C) * $f * ($sigma + $C * $sinSigma * ($cos2SigmaM + $C * cos($sigma) * (-1 + 2 * $cos2SigmaM ** 2)));
                $sinAlpha = ($LPrime - $lambdaPrime) / $D;
                $cosSqAlpha = (1 - $sinAlpha ** 2);
                $sinLambdaPrime = ($sinAlpha * $sinSigma) / (cos($U1) * cos($U2));
                $lambdaPrime = self::asin($sinLambdaPrime);
                $sinSqSigma = (cos($U2) * $sinLambdaPrime) ** 2 + (cos($U1) * sin($U2) + sin($U1) * cos($U2) * cos($lambdaPrime)) ** 2;
                $sinSigma = sqrt($sinSqSigma);
            } while (abs($sinAlpha - $sinAlphaN) >= static::ITERATION_CONVERGENCE_FORMULA);
        }

        $E = sqrt(1 + (($a ** 2 - $b ** 2) / $b ** 2) * $cosSqAlpha);
        $F = ($E - 1) / ($E + 1);

        $A = (1 + $F ** 2 / 4) / (1 - $F);
        $B = $F * (1 - 3 / 8 * $F ** 2);

        $deltaSigma = $B * $sinSigma * ($cos2SigmaM + $B / 4 * ($cosSigma * (-1 + 2 * $cos2SigmaM ** 2) - $B / 6 * $cos2SigmaM * (-3 + 4 * $sinSigma ** 2) * (-3 + 4 * $cos2SigmaM ** 2)));

        return new Metre($b * $A * ($sigma - $deltaSigma));
    }

    /**
     * General polynomial.
     * @param Coefficient[] $powerCoefficients
     */
    protected function generalPolynomialUnitless(
        float $xs,
        float $ys,
        UnitOfMeasure $ordinate1OfEvaluationPointInSourceCRS,
        UnitOfMeasure $ordinate2OfEvaluationPointInSourceCRS,
        UnitOfMeasure $ordinate1OfEvaluationPointInTargetCRS,
        UnitOfMeasure $ordinate2OfEvaluationPointInTargetCRS,
        Scale $scalingFactorForSourceCRSCoordDifferences,
        Scale $scalingFactorForTargetCRSCoordDifferences,
        Scale $A0,
        Scale $B0,
        array $powerCoefficients
    ): array {
        $xso = $ordinate1OfEvaluationPointInSourceCRS->getValue();
        $yso = $ordinate2OfEvaluationPointInSourceCRS->getValue();
        $xto = $ordinate1OfEvaluationPointInTargetCRS->getValue();
        $yto = $ordinate2OfEvaluationPointInTargetCRS->getValue();

        $U = $scalingFactorForSourceCRSCoordDifferences->asUnity()->getValue() * ($xs - $xso);
        $V = $scalingFactorForSourceCRSCoordDifferences->asUnity()->getValue() * ($ys - $yso);

        $mTdX = $A0->getValue();
        foreach ($powerCoefficients as $coefficientName => $coefficientValue) {
            if ($coefficientName[0] === 'A') {
                sscanf($coefficientName, 'Au%dv%d', $uPower, $vPower);

                $mTdX += $coefficientValue->getValue() * $U ** $uPower * $V ** $vPower;
            }
        }

        $mTdY = $B0->getValue();
        foreach ($powerCoefficients as $coefficientName => $coefficientValue) {
            if ($coefficientName[0] === 'B') {
                sscanf($coefficientName, 'Bu%dv%d', $uPower, $vPower);
                $mTdY += $coefficientValue->getValue() * $U ** $uPower * $V ** $vPower;
            }
        }

        $xt = $xs - $xso + $xto + $mTdX / $scalingFactorForTargetCRSCoordDifferences->asUnity()->getValue();
        $yt = $ys - $yso + $yto + $mTdY / $scalingFactorForTargetCRSCoordDifferences->asUnity()->getValue();

        return ['xt' => $xt, 'yt' => $yt];
    }

    /**
     * Reversible polynomial.
     */
    protected function reversiblePolynomialUnitless(
        float $xs,
        float $ys,
        Angle $ordinate1OfEvaluationPoint,
        Angle $ordinate2OfEvaluationPoint,
        Scale $scalingFactorForCoordDifferences,
        Scale $A0,
        Scale $B0,
        array $powerCoefficients
    ): array {
        $xo = $ordinate1OfEvaluationPoint->getValue();
        $yo = $ordinate2OfEvaluationPoint->getValue();

        $U = $scalingFactorForCoordDifferences->asUnity()->getValue() * ($xs - $xo);
        $V = $scalingFactorForCoordDifferences->asUnity()->getValue() * ($ys - $yo);

        $mTdX = $A0->getValue();
        foreach ($powerCoefficients as $coefficientName => $coefficientValue) {
            if ($coefficientName[0] === 'A') {
                sscanf($coefficientName, 'Au%dv%d', $uPower, $vPower);

                $mTdX += $coefficientValue->getValue() * $U ** $uPower * $V ** $vPower;
            }
        }

        $mTdY = $B0->getValue();
        foreach ($powerCoefficients as $coefficientName => $coefficientValue) {
            if ($coefficientName[0] === 'B') {
                sscanf($coefficientName, 'Bu%dv%d', $uPower, $vPower);
                $mTdY += $coefficientValue->getValue() * $U ** $uPower * $V ** $vPower;
            }
        }

        $xt = $xs + $mTdX * $scalingFactorForCoordDifferences->asUnity()->getValue();
        $yt = $ys + $mTdY * $scalingFactorForCoordDifferences->asUnity()->getValue();

        return ['xt' => $xt, 'yt' => $yt];
    }

    /**
     * Floating point vagaries mean that it's possible for inputs to be e.g. 1.00000000000001 which makes PHP give a
     * silent NaN as output so inputs need to be capped. atan/atan2 are not affected, they seem to cap internally.
     */
    protected static function acos(float $num): float
    {
        if ($num > 1.0) {
            $num = 1.0;
        } elseif ($num < -1) {
            $num = -1.0;
        }

        return acos($num);
    }

    /**
     * Floating point vagaries mean that it's possible for inputs to be e.g. 1.00000000000001 which makes PHP give a
     * silent NaN as output so inputs need to be capped. atan/atan2 are not affected, they seem to cap internally.
     */
    protected static function asin(float $num): float
    {
        if ($num > 1.0) {
            $num = 1.0;
        } elseif ($num < -1.0) {
            $num = -1.0;
        }

        return asin($num);
    }

    abstract public function getCRS(): CoordinateReferenceSystem;

    abstract public function getCoordinateEpoch(): ?DateTimeImmutable;

    abstract public function calculateDistance(self $to): Length;
}


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 asin;
14		use function atan;
15		use function atan2;
16		use function cos;
17		use DateTimeImmutable;
18		use const M_PI;
19		use PHPCoord\CoordinateOperation\CoordinateOperationMethods;
20		use PHPCoord\CoordinateOperation\CoordinateOperationParams;
		0 ignored issues – show Bug introduced 2021-01-03 21:47 UTC by Report Bug Copy Issue Report 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...
21		use PHPCoord\CoordinateOperation\CoordinateOperations;
22		use PHPCoord\CoordinateOperation\GeographicValue;
23		use PHPCoord\CoordinateReferenceSystem\CoordinateReferenceSystem;
24		use PHPCoord\Datum\Ellipsoid;
25		use PHPCoord\UnitOfMeasure\Angle\Angle;
26		use PHPCoord\UnitOfMeasure\Length\Length;
27		use PHPCoord\UnitOfMeasure\Length\Metre;
28		use PHPCoord\UnitOfMeasure\Scale\Coefficient;
29		use PHPCoord\UnitOfMeasure\Scale\Scale;
30		use PHPCoord\UnitOfMeasure\UnitOfMeasure;
31		use PHPCoord\UnitOfMeasure\UnitOfMeasureFactory;
32		use function sin;
33		use function sqrt;
34		use function sscanf;
35		use function str_starts_with;
36		use Stringable;
37		use function tan;
38
39		abstract class Point implements Stringable
40		{
41		protected const ITERATION_CONVERGENCE_FORMULA = 1e-10;
42		protected const ITERATION_CONVERGENCE_GRID = 0.0001;
43		protected const METHODS_REQUIRING_HORIZONTAL_POINT = [
44		CoordinateOperationMethods::EPSG_VERTICAL_OFFSET_AND_SLOPE => CoordinateOperationMethods::EPSG_VERTICAL_OFFSET_AND_SLOPE,
45		CoordinateOperationMethods::EPSG_ZERO_TIDE_HEIGHT_TO_MEAN_TIDE_HEIGHT_EVRF2019 => CoordinateOperationMethods::EPSG_ZERO_TIDE_HEIGHT_TO_MEAN_TIDE_HEIGHT_EVRF2019,
46		CoordinateOperationMethods::EPSG_VERTICAL_OFFSET_BY_GRID_INTERPOLATION_GTX => CoordinateOperationMethods::EPSG_VERTICAL_OFFSET_BY_GRID_INTERPOLATION_GTX,
47		];
48		protected const METHODS_THAT_REQUIRE_DIRECTION = [
49		CoordinateOperationMethods::EPSG_SIMILARITY_TRANSFORMATION => CoordinateOperationMethods::EPSG_SIMILARITY_TRANSFORMATION,
50		CoordinateOperationMethods::EPSG_AFFINE_PARAMETRIC_TRANSFORMATION => CoordinateOperationMethods::EPSG_AFFINE_PARAMETRIC_TRANSFORMATION,
51		CoordinateOperationMethods::EPSG_NADCON5_2D => CoordinateOperationMethods::EPSG_NADCON5_2D,
52		CoordinateOperationMethods::EPSG_NADCON5_3D => CoordinateOperationMethods::EPSG_NADCON5_3D,
53		CoordinateOperationMethods::EPSG_NTV2 => CoordinateOperationMethods::EPSG_NTV2,
54		CoordinateOperationMethods::EPSG_ZERO_TIDE_HEIGHT_TO_MEAN_TIDE_HEIGHT_EVRF2019 => CoordinateOperationMethods::EPSG_ZERO_TIDE_HEIGHT_TO_MEAN_TIDE_HEIGHT_EVRF2019,
55		CoordinateOperationMethods::EPSG_GEOCENTRIC_TRANSLATION_BY_GRID_INTERPOLATION_IGN => CoordinateOperationMethods::EPSG_GEOCENTRIC_TRANSLATION_BY_GRID_INTERPOLATION_IGN,
56		CoordinateOperationMethods::EPSG_VERTICAL_OFFSET_BY_GRID_INTERPOLATION_GTX => CoordinateOperationMethods::EPSG_VERTICAL_OFFSET_BY_GRID_INTERPOLATION_GTX,
57		];
58
59		/**
60		* @internal
61		*/
62		public function performOperation(string $srid, CoordinateReferenceSystem $to, bool $inReverse, array $additionalParams = []): self
63		{
64		$operation = CoordinateOperations::getOperationData($srid);
65		$method = CoordinateOperationMethods::getFunctionName($operation['method']);
66		$params = self::resolveParamsByOperation($srid, $operation['method'], $inReverse);
67
68		if (isset(self::METHODS_REQUIRING_HORIZONTAL_POINT[$operation['method']])) {
69		$params['horizontalPoint'] = $additionalParams['horizontalPoint'];
70		}
71	63200
72		return $this->$method($to, ...$params);
73	63200	}
74
75	63200	protected static function resolveParamsByOperation(string $operationSrid, string $methodSrid, bool $inReverse): array
76	63200	{
77	63200	$params = [];
78	63200	$powerCoefficients = [];
79	504	foreach (CoordinateOperationParams::getParamData($operationSrid) as $paramName => $paramData) {
80	504	if (isset($paramData['fileProvider'])) {
81		$params[$paramName] = (new $paramData['fileProvider']())->provideGrid();
82		} else {
83	504	if ($inReverse && $paramData['reverses']) {
84	504	$paramData['value'] *= -1;
85		}
86		if ($paramData['uom']) {
87	62763	$param = UnitOfMeasureFactory::makeUnit($paramData['value'], $paramData['uom']);
88		} else {
89		$param = $paramData['value'];
90	63200	}
91		if (str_starts_with($paramName, 'Au') \|\| str_starts_with($paramName, 'Bu')) {
92	63200	$powerCoefficients[$paramName] = $param;
93	409	} else {
94		$params[$paramName] = $param;
95		}
96	63200	}
97	63200	}
98		if ($powerCoefficients) {
99		$params['powerCoefficients'] = $powerCoefficients;
100		}
101		if (isset(self::METHODS_THAT_REQUIRE_DIRECTION[$methodSrid])) {
102		$params['inReverse'] = $inReverse;
103		}
104	63200
105	3	return $params;
106		}
107	63200
108		protected static function sign(float $number): int
109	63200	{
110		if ($number < 0) {
111		return -1;
112		}
113
114		return 1;
115		}
116
117		/**
118		* Calculate surface distance between two points.
119		*/
120		protected static function vincenty(GeographicValue $from, GeographicValue $to, Ellipsoid $ellipsoid): Length
121		{
122		$a = $ellipsoid->getSemiMajorAxis()->asMetres()->getValue();
123		$b = $ellipsoid->getSemiMinorAxis()->asMetres()->getValue();
124		$f = $ellipsoid->getInverseFlattening();
125	63209	$U1 = atan((1 - $f) * tan($from->getLatitude()->asRadians()->getValue()));
126		$U2 = atan((1 - $f) * tan($to->getLatitude()->asRadians()->getValue()));
127	63209	$L = $to->getLongitude()->subtract($from->getLongitude())->asRadians()->getValue();
128	63209
129	63209	$lambda = $L;
130	520	do {
131	520	$lambdaN = $lambda;
132	520
133	520	$sinSigma = sqrt((cos($U2) * sin($lambda)) ** 2 + (cos($U1) * sin($U2) - sin($U1) * cos($U2) * cos($lambda)) ** 2);
134	520	$cosSigma = sin($U1) * sin($U2) + cos($U1) * cos($U2) * cos($lambda);
135		$sigma = atan2($sinSigma, $cosSigma);
136
137	62727	$sinAlpha = cos($U1) * cos($U2) * sin($lambda) / $sinSigma;
138		$cosSqAlpha = (1 - $sinAlpha ** 2);
139		$cos2SigmaM = $cosSqAlpha ? $cosSigma - (2 * sin($U1) * sin($U2) / $cosSqAlpha) : 0;
140	63209	$C = $f / 16 * $cosSqAlpha * (4 + $f * (4 - 3 * $cosSqAlpha));
141	62910	$lambda = $L + (1 - $C) * $f * $sinAlpha * ($sigma + $C * $sinSigma * ($cos2SigmaM + $C * $cosSigma * (-1 + 2 * $cos2SigmaM ** 2)));
142		} while (abs($lambda - $lambdaN) >= static::ITERATION_CONVERGENCE_FORMULA && abs($lambda) < M_PI);
143
144	63209	// Antipodal case
145		if (abs($lambda) >= M_PI) {
146		if ($L >= 0) {
147	63200	$LPrime = M_PI - $L;
148		} else {
149	63200	$LPrime = -M_PI - $L;
150	63200	}
151	63200
152	59474	$lambdaPrime = 0;
153	173	$sigma = M_PI - abs($U1 + $U2);
154		$sinSigma = sin($sigma);
155	59331	$cosSqAlpha = 0.5;
156	14267	$sinAlpha = 0;
157
158	59331	do {
159	59318	$sinAlphaN = $sinAlpha;
160
161	402	$C = $f / 16 * $cosSqAlpha * (4 + $f * (4 - 3 * $cosSqAlpha));
162		$cos2SigmaM = cos($sigma) - 2 * sin($U1) * sin($U2) / $cosSqAlpha;
163	59331	$D = (1 - $C) * $f * ($sigma + $C * $sinSigma * ($cos2SigmaM + $C * cos($sigma) * (-1 + 2 * $cos2SigmaM ** 2)));
164	99	$sinAlpha = ($LPrime - $lambdaPrime) / $D;
165		$cosSqAlpha = (1 - $sinAlpha ** 2);
166	59331	$sinLambdaPrime = ($sinAlpha * $sinSigma) / (cos($U1) * cos($U2));
167		$lambdaPrime = self::asin($sinLambdaPrime);
168		$sinSqSigma = (cos($U2) * $sinLambdaPrime) ** 2 + (cos($U1) * sin($U2) + sin($U1) * cos($U2) * cos($lambdaPrime)) ** 2;
169		$sinSigma = sqrt($sinSqSigma);
170	63200	} while (abs($sinAlpha - $sinAlphaN) >= static::ITERATION_CONVERGENCE_FORMULA);
171	99	}
172
173	63200	$E = sqrt(1 + (($a 2 - $b 2) / $b ** 2) * $cosSqAlpha);
174	430	$F = ($E - 1) / ($E + 1);
175
176		$A = (1 + $F ** 2 / 4) / (1 - $F);
177	63200	$B = $F * (1 - 3 / 8 * $F ** 2);
178
179		$deltaSigma = $B * $sinSigma * ($cos2SigmaM + $B / 4 * ($cosSigma * (-1 + 2 * $cos2SigmaM ** 2) - $B / 6 * $cos2SigmaM * (-3 + 4 * $sinSigma ** 2) * (-3 + 4 * $cos2SigmaM ** 2)));
180
181		return new Metre($b * $A * ($sigma - $deltaSigma));
182		}
183
184		/**
185		* General polynomial.
186		* @param Coefficient[] $powerCoefficients
187		*/
188		protected function generalPolynomialUnitless(
189		float $xs,
190		float $ys,
191		UnitOfMeasure $ordinate1OfEvaluationPointInSourceCRS,
192		UnitOfMeasure $ordinate2OfEvaluationPointInSourceCRS,
193		UnitOfMeasure $ordinate1OfEvaluationPointInTargetCRS,
194	8685	UnitOfMeasure $ordinate2OfEvaluationPointInTargetCRS,
195		Scale $scalingFactorForSourceCRSCoordDifferences,
196	8685	Scale $scalingFactorForTargetCRSCoordDifferences,
197	891	Scale $A0,
198		Scale $B0,
199		array $powerCoefficients
200	7839	): array {
201		$xso = $ordinate1OfEvaluationPointInSourceCRS->getValue();
202		$yso = $ordinate2OfEvaluationPointInSourceCRS->getValue();
203		$xto = $ordinate1OfEvaluationPointInTargetCRS->getValue();
204		$yto = $ordinate2OfEvaluationPointInTargetCRS->getValue();
205
206	153	$U = $scalingFactorForSourceCRSCoordDifferences->asUnity()->getValue() * ($xs - $xso);
207		$V = $scalingFactorForSourceCRSCoordDifferences->asUnity()->getValue() * ($ys - $yso);
208	153
209	153	$mTdX = $A0->getValue();
210	153	foreach ($powerCoefficients as $coefficientName => $coefficientValue) {
211	153	if ($coefficientName[0] === 'A') {
212	153	sscanf($coefficientName, 'Au%dv%d', $uPower, $vPower);
		1 ignored issue – show Comprehensibility Best Practice introduced 2020-12-24 20:47 UTC by Report Bug Copy Issue Report The variable `$vPower` seems to be never defined. Loading history...
213	153	$mTdX += $coefficientValue->getValue() * $U ** $uPower * $V ** $vPower;
214		}
215	153	}
216
217	153	$mTdY = $B0->getValue();
218		foreach ($powerCoefficients as $coefficientName => $coefficientValue) {
219	153	if ($coefficientName[0] === 'B') {
220	153	sscanf($coefficientName, 'Bu%dv%d', $uPower, $vPower);
221	153	$mTdY += $coefficientValue->getValue() * $U ** $uPower * $V ** $vPower;
222		}
223	153	}
224	153
225	153	$xt = $xs - $xso + $xto + $mTdX / $scalingFactorForTargetCRSCoordDifferences->asUnity()->getValue();
226	153	$yt = $ys - $yso + $yto + $mTdY / $scalingFactorForTargetCRSCoordDifferences->asUnity()->getValue();
227	153
228	153	return ['xt' => $xt, 'yt' => $yt];
229		}
230
231	153	/**
232	27	* Reversible polynomial.
233	27	*/
234		protected function reversiblePolynomialUnitless(
235		float $xs,
236		float $ys,
237		Angle $ordinate1OfEvaluationPoint,
238	27	Angle $ordinate2OfEvaluationPoint,
239	27	Scale $scalingFactorForCoordDifferences,
240	27	Scale $A0,
241	27	Scale $B0,
242	27	array $powerCoefficients
243		): array {
244		$xo = $ordinate1OfEvaluationPoint->getValue();
245	27	$yo = $ordinate2OfEvaluationPoint->getValue();
246
247	27	$U = $scalingFactorForCoordDifferences->asUnity()->getValue() * ($xs - $xo);
248	27	$V = $scalingFactorForCoordDifferences->asUnity()->getValue() * ($ys - $yo);
249	27
250	27	$mTdX = $A0->getValue();
251	27	foreach ($powerCoefficients as $coefficientName => $coefficientValue) {
252	27	if ($coefficientName[0] === 'A') {
253	27	sscanf($coefficientName, 'Au%dv%d', $uPower, $vPower);
		1 ignored issue – show Comprehensibility Best Practice introduced 2020-12-24 20:47 UTC by Report Bug Copy Issue Report The variable `$vPower` seems to be never defined. Loading history...
254	27	$mTdX += $coefficientValue->getValue() * $U ** $uPower * $V ** $vPower;
255	27	}
256	27	}
257
258		$mTdY = $B0->getValue();
259	153	foreach ($powerCoefficients as $coefficientName => $coefficientValue) {
260	153	if ($coefficientName[0] === 'B') {
261		sscanf($coefficientName, 'Bu%dv%d', $uPower, $vPower);
262	153	$mTdY += $coefficientValue->getValue() * $U ** $uPower * $V ** $vPower;
263	153	}
264		}
265	153
266		$xt = $xs + $mTdX * $scalingFactorForCoordDifferences->asUnity()->getValue();
267	153	$yt = $ys + $mTdY * $scalingFactorForCoordDifferences->asUnity()->getValue();
268
269		return ['xt' => $xt, 'yt' => $yt];
270		}
271
272		/**
273		* Floating point vagaries mean that it's possible for inputs to be e.g. 1.00000000000001 which makes PHP give a
274	54	* silent NaN as output so inputs need to be capped. atan/atan2 are not affected, they seem to cap internally.
275		*/
276		protected static function acos(float $num): float
277		{
278		if ($num > 1.0) {
279		$num = 1.0;
280		} elseif ($num < -1) {
281		$num = -1.0;
282		}
283
284		return acos($num);
285		}
286
287	54	/**
288	54	* Floating point vagaries mean that it's possible for inputs to be e.g. 1.00000000000001 which makes PHP give a
289	54	* silent NaN as output so inputs need to be capped. atan/atan2 are not affected, they seem to cap internally.
290	54	*/
291		protected static function asin(float $num): float
292	54	{
293	54	if ($num > 1.0) {
294		$num = 1.0;
295	54	} elseif ($num < -1.0) {
296	54	$num = -1.0;
297	54	}
298	54
299	54	return asin($num);
300		}
301
302		abstract public function getCRS(): CoordinateReferenceSystem;
303	54
304	54	abstract public function getCoordinateEpoch(): ?DateTimeImmutable;
305	54
306	54	abstract public function calculateDistance(self $to): Length;
307		}
308

dvdoug / PHPCoord

Push — extents ( f35511...ba1d5f )

Point::asin() A

Complexity

Size

Duplication

Code Coverage

Importance

Duplication Side-by-Side

Filter issues like