Passed
Push — 4.x ( f9780c...a64886 )
by Doug
05:20
created

GeocentricPoint::positionVectorTransformation()   A

Complexity

Conditions 1
Paths 1

Size

Total Lines 22
Code Lines 12

Duplication

Lines 0
Ratio 0 %

Code Coverage

Tests 6
CRAP Score 1

Importance

Changes 0
Metric Value
eloc 12
c 0
b 0
f 0
dl 0
loc 22
ccs 6
cts 6
cp 1
rs 9.8666
cc 1
nc 1
nop 8
crap 1

How to fix   Many Parameters   

Many Parameters

Methods with many parameters are not only hard to understand, but their parameters also often become inconsistent when you need more, or different data.

There are several approaches to avoid long parameter lists:

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 DateTime;
13
use DateTimeImmutable;
14
use DateTimeInterface;
15
use InvalidArgumentException;
16
use PHPCoord\CoordinateOperation\AutoConversion;
17
use PHPCoord\CoordinateOperation\GeocentricValue;
18
use PHPCoord\CoordinateOperation\GeographicValue;
19
use PHPCoord\CoordinateReferenceSystem\Geocentric;
20
use PHPCoord\CoordinateReferenceSystem\Geographic;
21
use PHPCoord\CoordinateReferenceSystem\Geographic3D;
22
use PHPCoord\CoordinateSystem\Axis;
23
use PHPCoord\Exception\InvalidCoordinateException;
24
use PHPCoord\UnitOfMeasure\Angle\Angle;
25
use PHPCoord\UnitOfMeasure\Angle\Radian;
26
use PHPCoord\UnitOfMeasure\Length\Length;
27
use PHPCoord\UnitOfMeasure\Length\Metre;
28
use PHPCoord\UnitOfMeasure\Rate;
29
use PHPCoord\UnitOfMeasure\Scale\Scale;
30
use PHPCoord\UnitOfMeasure\Scale\Unity;
31
use PHPCoord\UnitOfMeasure\Time\Time;
32
use PHPCoord\UnitOfMeasure\Time\Year;
33
use function sprintf;
34
use function sqrt;
35
36
/**
37
 * Coordinate representing a point in ECEF geocentric form.
38
 */
39
class GeocentricPoint extends Point
40
{
41
    use AutoConversion;
42
43
    /**
44
     * X co-ordinate.
45
     */
46
    protected Length $x;
47
48
    /**
49
     * Y co-ordinate.
50
     */
51
    protected Length $y;
52
53
    /**
54
     * Z co-ordinate.
55
     */
56
    protected Length $z;
57
58
    /**
59
     * Coordinate reference system.
60
     */
61
    protected Geocentric $crs;
62
63
    /**
64
     * Coordinate epoch (date for which the specified coordinates represented this point).
65
     */
66
    protected ?DateTimeImmutable $epoch;
67
68 30
    protected function __construct(Length $x, Length $y, Length $z, Geocentric $crs, ?DateTimeInterface $epoch = null)
69
    {
70 30
        $this->crs = $crs;
71 30
        $this->x = Length::convert($x, $this->getAxisByName(Axis::GEOCENTRIC_X)->getUnitOfMeasureId());
72 30
        $this->y = Length::convert($y, $this->getAxisByName(Axis::GEOCENTRIC_Y)->getUnitOfMeasureId());
73 30
        $this->z = Length::convert($z, $this->getAxisByName(Axis::GEOCENTRIC_Z)->getUnitOfMeasureId());
74
75 30
        if ($epoch instanceof DateTime) {
76 6
            $epoch = DateTimeImmutable::createFromMutable($epoch);
77
        }
78 30
        $this->epoch = $epoch;
79 30
    }
80
81
    /**
82
     * @param Length $x refer to CRS for preferred unit of measure, but any length unit accepted
83
     * @param Length $y refer to CRS for preferred unit of measure, but any length unit accepted
84
     * @param Length $z refer to CRS for preferred unit of measure, but any length unit accepted
85
     */
86 30
    public static function create(Length $x, Length $y, Length $z, Geocentric $crs, ?DateTimeInterface $epoch = null): self
87
    {
88 30
        return new static($x, $y, $z, $crs, $epoch);
89
    }
90
91 20
    public function getX(): Length
92
    {
93 20
        return $this->x;
94
    }
95
96 20
    public function getY(): Length
97
    {
98 20
        return $this->y;
99
    }
100
101 20
    public function getZ(): Length
102
    {
103 20
        return $this->z;
104
    }
105
106 30
    public function getCRS(): Geocentric
107
    {
108 30
        return $this->crs;
109
    }
110
111 5
    public function getCoordinateEpoch(): ?DateTimeImmutable
112
    {
113 5
        return $this->epoch;
114
    }
115
116
    /**
117
     * Calculate surface distance between two points.
118
     * Note: this implementation is currently not accurate over long distances, it is the straight line distance, not
119
     * the surface distance.
120
     */
121 2
    public function calculateDistance(Point $to): Length
122
    {
123 2
        if ($to->getCRS()->getSRID() !== $this->crs->getSRID()) {
124 1
            throw new InvalidArgumentException('Can only calculate distances between two points in the same CRS');
125
        }
126
127
        /* @var GeocentricPoint $to */
128 1
        return new Metre(
129 1
            sqrt(
130 1
                ($to->getX()->getValue() - $this->x->getValue()) ** 2 +
131 1
                ($to->getY()->getValue() - $this->y->getValue()) ** 2 +
132 1
                ($to->getZ()->getValue() - $this->z->getValue()) ** 2
133
            )
134
        );
135
    }
136
137 3
    public function __toString(): string
138
    {
139 3
        return "({$this->x}, {$this->y}, {$this->z})";
140
    }
141
142
    /**
143
     * Geographic/geocentric conversions
144
     * In applications it is often concatenated with the 3- 7- or 10-parameter transformations 9603, 9606, 9607 or
145
     * 9636 to form a geographic to geographic transformation.
146
     */
147 2
    public function geographicGeocentric(
148
        Geographic $to
149
    ): GeographicPoint {
150 2
        $geocentricValue = new GeocentricValue($this->x, $this->y, $this->z, $to->getDatum());
151 2
        $asGeographic = $geocentricValue->asGeographicValue();
152
153 2
        return GeographicPoint::create($asGeographic->getLatitude(), $asGeographic->getLongitude(), $to instanceof Geographic3D ? $asGeographic->getHeight() : null, $to, $this->epoch);
154
    }
155
156
    /**
157
     * Coordinate Frame rotation (geocentric domain)
158
     * This method is a specific case of the Molodensky-Badekas (CF) method (code 1034) in which the evaluation point
159
     * is at the geocentre with coordinate values of zero. Note the analogy with the Position Vector method (code 1033)
160
     * but beware of the differences!
161
     */
162 6
    public function coordinateFrameRotation(
163
        Geocentric $to,
164
        Length $xAxisTranslation,
165
        Length $yAxisTranslation,
166
        Length $zAxisTranslation,
167
        Angle $xAxisRotation,
168
        Angle $yAxisRotation,
169
        Angle $zAxisRotation,
170
        Scale $scaleDifference
171
    ): self {
172 6
        return $this->coordinateFrameMolodenskyBadekas(
173 6
            $to,
174
            $xAxisTranslation,
175
            $yAxisTranslation,
176
            $zAxisTranslation,
177
            $xAxisRotation,
178
            $yAxisRotation,
179
            $zAxisRotation,
180
            $scaleDifference,
181 6
            new Metre(0),
182 6
            new Metre(0),
183 6
            new Metre(0)
184
        );
185
    }
186
187
    /**
188
     * Molodensky-Badekas (CF geocentric domain)
189
     * See method codes 1039 and 9636 for this operation in other coordinate domains and method code 1061 for opposite
190
     * rotation convention in geocentric domain.
191
     */
192 7
    public function coordinateFrameMolodenskyBadekas(
193
        Geocentric $to,
194
        Length $xAxisTranslation,
195
        Length $yAxisTranslation,
196
        Length $zAxisTranslation,
197
        Angle $xAxisRotation,
198
        Angle $yAxisRotation,
199
        Angle $zAxisRotation,
200
        Scale $scaleDifference,
201
        Length $ordinate1OfEvaluationPoint,
202
        Length $ordinate2OfEvaluationPoint,
203
        Length $ordinate3OfEvaluationPoint
204
    ): self {
205 7
        $xs = $this->x->asMetres()->getValue();
206 7
        $ys = $this->y->asMetres()->getValue();
207 7
        $zs = $this->z->asMetres()->getValue();
208 7
        $tx = $xAxisTranslation->asMetres()->getValue();
209 7
        $ty = $yAxisTranslation->asMetres()->getValue();
210 7
        $tz = $zAxisTranslation->asMetres()->getValue();
211 7
        $rx = $xAxisRotation->asRadians()->getValue();
212 7
        $ry = $yAxisRotation->asRadians()->getValue();
213 7
        $rz = $zAxisRotation->asRadians()->getValue();
214 7
        $M = 1 + $scaleDifference->asUnity()->getValue();
215 7
        $xp = $ordinate1OfEvaluationPoint->asMetres()->getValue();
216 7
        $yp = $ordinate2OfEvaluationPoint->asMetres()->getValue();
217 7
        $zp = $ordinate3OfEvaluationPoint->asMetres()->getValue();
218
219 7
        $xt = $M * ((($xs - $xp) * 1) + (($ys - $yp) * $rz) + (($zs - $zp) * -$ry)) + $tx + $xp;
220 7
        $yt = $M * ((($xs - $xp) * -$rz) + (($ys - $yp) * 1) + (($zs - $zp) * $rx)) + $ty + $yp;
221 7
        $zt = $M * ((($xs - $xp) * $ry) + (($ys - $yp) * -$rx) + (($zs - $zp) * 1)) + $tz + $zp;
222
223 7
        return static::create(new Metre($xt), new Metre($yt), new Metre($zt), $to, $this->epoch);
224
    }
225
226
    /**
227
     * Position Vector transformation (geocentric domain)
228
     * This method is a specific case of the Molodensky-Badekas (PV) method (code 1061) in which the evaluation point
229
     * is the geocentre with coordinate values of zero. Note the analogy with the Coordinate Frame method (code 1032)
230
     * but beware of the differences!
231
     */
232 6
    public function positionVectorTransformation(
233
        Geocentric $to,
234
        Length $xAxisTranslation,
235
        Length $yAxisTranslation,
236
        Length $zAxisTranslation,
237
        Angle $xAxisRotation,
238
        Angle $yAxisRotation,
239
        Angle $zAxisRotation,
240
        Scale $scaleDifference
241
    ): self {
242 6
        return $this->positionVectorMolodenskyBadekas(
243 6
            $to,
244
            $xAxisTranslation,
245
            $yAxisTranslation,
246
            $zAxisTranslation,
247
            $xAxisRotation,
248
            $yAxisRotation,
249
            $zAxisRotation,
250
            $scaleDifference,
251 6
            new Metre(0),
252 6
            new Metre(0),
253 6
            new Metre(0)
254
        );
255
    }
256
257
    /**
258
     * Molodensky-Badekas (PV geocentric domain)
259
     * See method codes 1062 and 1063 for this operation in other coordinate domains and method code 1034 for opposite
260
     * rotation convention in geocentric domain.
261
     */
262 7
    public function positionVectorMolodenskyBadekas(
263
        Geocentric $to,
264
        Length $xAxisTranslation,
265
        Length $yAxisTranslation,
266
        Length $zAxisTranslation,
267
        Angle $xAxisRotation,
268
        Angle $yAxisRotation,
269
        Angle $zAxisRotation,
270
        Scale $scaleDifference,
271
        Length $ordinate1OfEvaluationPoint,
272
        Length $ordinate2OfEvaluationPoint,
273
        Length $ordinate3OfEvaluationPoint
274
    ): self {
275 7
        $xs = $this->x->asMetres()->getValue();
276 7
        $ys = $this->y->asMetres()->getValue();
277 7
        $zs = $this->z->asMetres()->getValue();
278 7
        $tx = $xAxisTranslation->asMetres()->getValue();
279 7
        $ty = $yAxisTranslation->asMetres()->getValue();
280 7
        $tz = $zAxisTranslation->asMetres()->getValue();
281 7
        $rx = $xAxisRotation->asRadians()->getValue();
282 7
        $ry = $yAxisRotation->asRadians()->getValue();
283 7
        $rz = $zAxisRotation->asRadians()->getValue();
284 7
        $M = 1 + $scaleDifference->asUnity()->getValue();
285 7
        $xp = $ordinate1OfEvaluationPoint->asMetres()->getValue();
286 7
        $yp = $ordinate2OfEvaluationPoint->asMetres()->getValue();
287 7
        $zp = $ordinate3OfEvaluationPoint->asMetres()->getValue();
288
289 7
        $xt = $M * ((($xs - $xp) * 1) + (($ys - $yp) * -$rz) + (($zs - $zp) * $ry)) + $tx + $xp;
290 7
        $yt = $M * ((($xs - $xp) * $rz) + (($ys - $yp) * 1) + (($zs - $zp) * -$rx)) + $ty + $yp;
291 7
        $zt = $M * ((($xs - $xp) * -$ry) + (($ys - $yp) * $rx) + (($zs - $zp) * 1)) + $tz + $zp;
292
293 7
        return static::create(new Metre($xt), new Metre($yt), new Metre($zt), $to, $this->epoch);
294
    }
295
296
    /**
297
     * Geocentric translations
298
     * This method allows calculation of geocentric coords in the target system by adding the parameter values to the
299
     * corresponding coordinates of the point in the source system. See methods 1035 and 9603 for similar tfms
300
     * operating between other CRSs types.
301
     */
302 1
    public function geocentricTranslation(
303
        Geocentric $to,
304
        Length $xAxisTranslation,
305
        Length $yAxisTranslation,
306
        Length $zAxisTranslation
307
    ): self {
308 1
        return $this->positionVectorTransformation(
309 1
            $to,
310
            $xAxisTranslation,
311
            $yAxisTranslation,
312
            $zAxisTranslation,
313 1
            new Radian(0),
314 1
            new Radian(0),
315 1
            new Radian(0),
316 1
            new Unity(0)
317
        );
318
    }
319
320
    /**
321
     * Time-dependent Coordinate Frame rotation (geocen)
322
     * Note the analogy with the Time-dependent Position Vector transformation (code 1053) but beware of the
323
     * differences!  The Position Vector convention is used by IAG. See method codes 1057 and 1058 for similar methods
324
     * operating between other CRS types.
325
     */
326 2
    public function timeDependentCoordinateFrameRotation(
327
        Geocentric $to,
328
        Length $xAxisTranslation,
329
        Length $yAxisTranslation,
330
        Length $zAxisTranslation,
331
        Angle $xAxisRotation,
332
        Angle $yAxisRotation,
333
        Angle $zAxisRotation,
334
        Scale $scaleDifference,
335
        Rate $rateOfChangeOfXAxisTranslation,
336
        Rate $rateOfChangeOfYAxisTranslation,
337
        Rate $rateOfChangeOfZAxisTranslation,
338
        Rate $rateOfChangeOfXAxisRotation,
339
        Rate $rateOfChangeOfYAxisRotation,
340
        Rate $rateOfChangeOfZAxisRotation,
341
        Rate $rateOfChangeOfScaleDifference,
342
        Time $parameterReferenceEpoch
343
    ): self {
344 2
        if ($this->epoch === null) {
345 1
            throw new InvalidCoordinateException('This transformation requires an epoch, none given');
346
        }
347
348
        // Points use PHP DateTimes for ease of use, but transformations use decimal years...
349 1
        $pointEpoch = Year::fromDateTime($this->epoch);
350 1
        $yearsToAdjust = $pointEpoch->subtract($parameterReferenceEpoch)->getValue();
351 1
        $xAxisTranslation = $xAxisTranslation->add($rateOfChangeOfXAxisTranslation->getChangePerYear()->multiply($yearsToAdjust));
0 ignored issues
show
Bug introduced by
The method multiply() does not exist on PHPCoord\UnitOfMeasure\UnitOfMeasure. It seems like you code against a sub-type of said class. However, the method does not exist in PHPCoord\UnitOfMeasure\Rate. Are you sure you never get one of those? ( Ignorable by Annotation )

If this is a false-positive, you can also ignore this issue in your code via the ignore-call  annotation

351
        $xAxisTranslation = $xAxisTranslation->add($rateOfChangeOfXAxisTranslation->getChangePerYear()->/** @scrutinizer ignore-call */ multiply($yearsToAdjust));
Loading history...
352 1
        $yAxisTranslation = $yAxisTranslation->add($rateOfChangeOfYAxisTranslation->getChangePerYear()->multiply($yearsToAdjust));
353 1
        $zAxisTranslation = $zAxisTranslation->add($rateOfChangeOfZAxisTranslation->getChangePerYear()->multiply($yearsToAdjust));
354 1
        $xAxisRotation = $xAxisRotation->add($rateOfChangeOfXAxisRotation->getChangePerYear()->multiply($yearsToAdjust));
355 1
        $yAxisRotation = $yAxisRotation->add($rateOfChangeOfYAxisRotation->getChangePerYear()->multiply($yearsToAdjust));
356 1
        $zAxisRotation = $zAxisRotation->add($rateOfChangeOfZAxisRotation->getChangePerYear()->multiply($yearsToAdjust));
357 1
        $scaleDifference = $scaleDifference->add($rateOfChangeOfScaleDifference->getChangePerYear()->multiply($yearsToAdjust));
358
359 1
        return $this->coordinateFrameRotation($to, $xAxisTranslation, $yAxisTranslation, $zAxisTranslation, $xAxisRotation, $yAxisRotation, $zAxisRotation, $scaleDifference);
360
    }
361
362
    /**
363
     * Time-dependent Position Vector tfm (geocentric)
364
     * Note the analogy with the Time-dependent Coordinate Frame rotation (code 1056) but beware of the differences!
365
     * The Position Vector convention is used by IAG. See method codes 1054 and 1055 for similar methods operating
366
     * between other CRS types.
367
     */
368 2
    public function timeDependentPositionVectorTransformation(
369
        Geocentric $to,
370
        Length $xAxisTranslation,
371
        Length $yAxisTranslation,
372
        Length $zAxisTranslation,
373
        Angle $xAxisRotation,
374
        Angle $yAxisRotation,
375
        Angle $zAxisRotation,
376
        Scale $scaleDifference,
377
        Rate $rateOfChangeOfXAxisTranslation,
378
        Rate $rateOfChangeOfYAxisTranslation,
379
        Rate $rateOfChangeOfZAxisTranslation,
380
        Rate $rateOfChangeOfXAxisRotation,
381
        Rate $rateOfChangeOfYAxisRotation,
382
        Rate $rateOfChangeOfZAxisRotation,
383
        Rate $rateOfChangeOfScaleDifference,
384
        Time $parameterReferenceEpoch
385
    ): self {
386 2
        if ($this->epoch === null) {
387 1
            throw new InvalidCoordinateException('This transformation requires an epoch, none given');
388
        }
389
390
        // Points use PHP DateTimes for ease of use, but transformations use decimal years...
391 1
        $pointEpoch = Year::fromDateTime($this->epoch);
392 1
        $yearsToAdjust = $pointEpoch->subtract($parameterReferenceEpoch)->getValue();
393 1
        $xAxisTranslation = $xAxisTranslation->add($rateOfChangeOfXAxisTranslation->getChangePerYear()->multiply($yearsToAdjust));
394 1
        $yAxisTranslation = $yAxisTranslation->add($rateOfChangeOfYAxisTranslation->getChangePerYear()->multiply($yearsToAdjust));
395 1
        $zAxisTranslation = $zAxisTranslation->add($rateOfChangeOfZAxisTranslation->getChangePerYear()->multiply($yearsToAdjust));
396 1
        $xAxisRotation = $xAxisRotation->add($rateOfChangeOfXAxisRotation->getChangePerYear()->multiply($yearsToAdjust));
397 1
        $yAxisRotation = $yAxisRotation->add($rateOfChangeOfYAxisRotation->getChangePerYear()->multiply($yearsToAdjust));
398 1
        $zAxisRotation = $zAxisRotation->add($rateOfChangeOfZAxisRotation->getChangePerYear()->multiply($yearsToAdjust));
399 1
        $scaleDifference = $scaleDifference->add($rateOfChangeOfScaleDifference->getChangePerYear()->multiply($yearsToAdjust));
400
401 1
        return $this->positionVectorTransformation($to, $xAxisTranslation, $yAxisTranslation, $zAxisTranslation, $xAxisRotation, $yAxisRotation, $zAxisRotation, $scaleDifference);
402
    }
403
404
    /**
405
     * Time-specific Coordinate Frame rotation (geocen)
406
     * Note the analogy with the Time-specific Position Vector method (code 1065) but beware of the differences!
407
     */
408 4
    public function timeSpecificCoordinateFrameRotation(
409
        Geocentric $to,
410
        Length $xAxisTranslation,
411
        Length $yAxisTranslation,
412
        Length $zAxisTranslation,
413
        Angle $xAxisRotation,
414
        Angle $yAxisRotation,
415
        Angle $zAxisRotation,
416
        Scale $scaleDifference,
417
        Time $transformationReferenceEpoch
418
    ): self {
419 4
        if ($this->epoch === null) {
420 1
            throw new InvalidCoordinateException(sprintf('This transformation is only valid for epoch %s, none given', $transformationReferenceEpoch->getValue()));
421
        }
422
423
        // Points use PHP DateTimes for ease of use, but transformations use decimal years...
424 3
        $pointEpoch = Year::fromDateTime($this->epoch);
425
426 3
        if (abs($pointEpoch->getValue() - $transformationReferenceEpoch->getValue()) > 0.001) {
427 1
            throw new InvalidCoordinateException(sprintf('This transformation is only valid for epoch %s, got %s', $transformationReferenceEpoch, $pointEpoch));
428
        }
429
430 2
        return $this->coordinateFrameRotation($to, $xAxisTranslation, $yAxisTranslation, $zAxisTranslation, $xAxisRotation, $yAxisRotation, $zAxisRotation, $scaleDifference);
431
    }
432
433
    /**
434
     * Time-specific Position Vector transform (geocen)
435
     * Note the analogy with the Time-specifc Coordinate Frame method (code 1066) but beware of the differences!
436
     */
437 3
    public function timeSpecificPositionVectorTransformation(
438
        Geocentric $to,
439
        Length $xAxisTranslation,
440
        Length $yAxisTranslation,
441
        Length $zAxisTranslation,
442
        Angle $xAxisRotation,
443
        Angle $yAxisRotation,
444
        Angle $zAxisRotation,
445
        Scale $scaleDifference,
446
        Time $transformationReferenceEpoch
447
    ): self {
448 3
        if ($this->epoch === null) {
449 1
            throw new InvalidCoordinateException(sprintf('This transformation is only valid for epoch %s, none given', $transformationReferenceEpoch->getValue()));
450
        }
451
452
        // Points use PHP DateTimes for ease of use, but transformations use decimal years...
453 2
        $pointEpoch = Year::fromDateTime($this->epoch);
454
455 2
        if (abs($pointEpoch->getValue() - $transformationReferenceEpoch->getValue()) > 0.001) {
456 1
            throw new InvalidCoordinateException(sprintf('This transformation is only valid for epoch %s, got %s', $transformationReferenceEpoch, $pointEpoch));
457
        }
458
459 1
        return $this->positionVectorTransformation($to, $xAxisTranslation, $yAxisTranslation, $zAxisTranslation, $xAxisRotation, $yAxisRotation, $zAxisRotation, $scaleDifference);
460
    }
461
462 2
    public function asGeographicValue(): GeographicValue
463
    {
464 2
        return (new GeocentricValue($this->x, $this->y, $this->z, $this->getCRS()->getDatum()))->asGeographicValue();
465
    }
466
}
467