-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathastrometry.go
5049 lines (3985 loc) · 174 KB
/
astrometry.go
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
586
587
588
589
590
591
592
593
594
595
596
597
598
599
600
601
602
603
604
605
606
607
608
609
610
611
612
613
614
615
616
617
618
619
620
621
622
623
624
625
626
627
628
629
630
631
632
633
634
635
636
637
638
639
640
641
642
643
644
645
646
647
648
649
650
651
652
653
654
655
656
657
658
659
660
661
662
663
664
665
666
667
668
669
670
671
672
673
674
675
676
677
678
679
680
681
682
683
684
685
686
687
688
689
690
691
692
693
694
695
696
697
698
699
700
701
702
703
704
705
706
707
708
709
710
711
712
713
714
715
716
717
718
719
720
721
722
723
724
725
726
727
728
729
730
731
732
733
734
735
736
737
738
739
740
741
742
743
744
745
746
747
748
749
750
751
752
753
754
755
756
757
758
759
760
761
762
763
764
765
766
767
768
769
770
771
772
773
774
775
776
777
778
779
780
781
782
783
784
785
786
787
788
789
790
791
792
793
794
795
796
797
798
799
800
801
802
803
804
805
806
807
808
809
810
811
812
813
814
815
816
817
818
819
820
821
822
823
824
825
826
827
828
829
830
831
832
833
834
835
836
837
838
839
840
841
842
843
844
845
846
847
848
849
850
851
852
853
854
855
856
857
858
859
860
861
862
863
864
865
866
867
868
869
870
871
872
873
874
875
876
877
878
879
880
881
882
883
884
885
886
887
888
889
890
891
892
893
894
895
896
897
898
899
900
901
902
903
904
905
906
907
908
909
910
911
912
913
914
915
916
917
918
919
920
921
922
923
924
925
926
927
928
929
930
931
932
933
934
935
936
937
938
939
940
941
942
943
944
945
946
947
948
949
950
951
952
953
954
955
956
957
958
959
960
961
962
963
964
965
966
967
968
969
970
971
972
973
974
975
976
977
978
979
980
981
982
983
984
985
986
987
988
989
990
991
992
993
994
995
996
997
998
999
1000
// Copyright 2022 HE Boliang
// All rights reserved.
package gofa
// SOFA Astrometry Tools
/*
Ab Apply stellar aberration
Apply aberration to transform natural direction into proper
direction.
Given:
pnat [3]float64 natural direction to the source (unit vector)
v [3]float64 observer barycentric velocity in units of c
s float64 distance between the Sun and the observer (au)
bm1 float64 sqrt(1-|v|^2): reciprocal of Lorenz factor
Returned:
ppr [3]float64 proper direction to source (unit vector)
Notes:
1. The algorithm is based on Expr. (7.40) in the Explanatory
Supplement (Urban & Seidelmann 2013), but with the following
changes:
o Rigorous rather than approximate normalization is applied.
o The gravitational potential term from Expr. (7) in
Klioner (2003) is added, taking into account only the Sun's
contribution. This has a maximum effect of about
0.4 microarcsecond.
2. In almost all cases, the maximum accuracy will be limited by the
supplied velocity. For example, if the SOFA Epv00 function is
used, errors of up to 5 microarcseconds could occur.
References:
Urban, S. & Seidelmann, P. K. (eds), Explanatory Supplement to
the Astronomical Almanac, 3rd ed., University Science Books
(2013).
Klioner, Sergei A., "A practical relativistic model for micro-
arcsecond astrometry in space", Astr. J. 125, 1580-1597 (2003).
Called:
Pdp scalar product of two p-vectors
*/
func Ab(pnat, v [3]float64, s, bm1 float64, ppr *[3]float64) {
var i int
var pdv, w1, w2, r2, w, r float64
var p [3]float64
pdv = Pdp(pnat, v)
w1 = 1.0 + pdv/(1.0+bm1)
w2 = SRS / s
r2 = 0.0
for i = 0; i < 3; i++ {
w = pnat[i]*bm1 + w1*v[i] + w2*(v[i]-pdv*pnat[i])
p[i] = w
r2 = r2 + w*w
}
r = sqrt(r2)
for i = 0; i < 3; i++ {
ppr[i] = p[i] / r
}
}
/*
Apcg Prepare for ICRS <-> GCRS, geocentric, special
For a geocentric observer, prepare star-independent astrometry
parameters for transformations between ICRS and GCRS coordinates.
The Earth ephemeris is supplied by the caller.
The parameters produced by this function are required in the
parallax, light deflection and aberration parts of the astrometric
transformation chain.
Given:
date1 float64 TDB as a 2-part...
date2 float64 ...Julian Date (Note 1)
ebpv [2][3]float64 Earth barycentric pos/vel (au, au/day)
ehp [3]float64 Earth heliocentric position (au)
Returned:
astrom ASTROM star-independent astrometry parameters:
pmt float64 PM time interval (SSB, Julian years)
eb [3]float64 SSB to observer (vector, au)
eh [3]float64 Sun to observer (unit vector)
em float64 distance from Sun to observer (au)
v [3]float64 barycentric observer velocity (vector, c)
bm1 float64 sqrt(1-|v|^2): reciprocal of Lorenz factor
bpn [3][3]float64 bias-precession-nutation matrix
along float64 unchanged
xpl float64 unchanged
ypl float64 unchanged
sphi float64 unchanged
cphi float64 unchanged
diurab float64 unchanged
eral float64 unchanged
refa float64 unchanged
refb float64 unchanged
Notes:
1. The TDB date date1+date2 is a Julian Date, apportioned in any
convenient way between the two arguments. For example,
JD(TDB)=2450123.7 could be expressed in any of these ways, among
others:
date1 date2
2450123.7 0.0 (JD method)
2451545.0 -1421.3 (J2000 method)
2400000.5 50123.2 (MJD method)
2450123.5 0.2 (date & time method)
The JD method is the most natural and convenient to use in cases
where the loss of several decimal digits of resolution is
acceptable. The J2000 method is best matched to the way the
argument is handled internally and will deliver the optimum
resolution. The MJD method and the date & time methods are both
good compromises between resolution and convenience. For most
applications of this function the choice will not be at all
critical.
TT can be used instead of TDB without any significant impact on
accuracy.
2. All the vectors are with respect to BCRS axes.
3. This is one of several functions that inserts into the astrom
structure star-independent parameters needed for the chain of
astrometric transformations ICRS <-> GCRS <-> CIRS <-> observed.
The various functions support different classes of observer and
portions of the transformation chain:
functions observer transformation
Apcg Apcg13 geocentric ICRS <-> GCRS
Apci Apci13 terrestrial ICRS <-> CIRS
Apco Apco13 terrestrial ICRS <-> observed
Apcs Apcs13 space ICRS <-> GCRS
Aper Aper13 terrestrial update Earth rotation
Apio Apio13 terrestrial CIRS <-> observed
Those with names ending in "13" use contemporary SOFA models to
compute the various ephemerides. The others accept ephemerides
supplied by the caller.
The transformation from ICRS to GCRS covers space motion,
parallax, light deflection, and aberration. From GCRS to CIRS
comprises frame bias and precession-nutation. From CIRS to
observed takes account of Earth rotation, polar motion, diurnal
aberration and parallax (unless subsumed into the ICRS <-> GCRS
transformation), and atmospheric refraction.
4. The context structure astrom produced by this function is used by
Atciq* and Aticq*.
Called:
Apcs astrometry parameters, ICRS-GCRS, space observer
*/
func Apcg(date1, date2 float64, ebpv [2][3]float64, ehp [3]float64, astrom *ASTROM) {
/* Geocentric observer */
pv := [2][3]float64{
{0.0, 0.0, 0.0},
{0.0, 0.0, 0.0},
}
/* Compute the star-independent astrometry parameters. */
Apcs(date1, date2, pv, ebpv, ehp, astrom)
}
/*
Apcg13 Prepare for ICRS <-> GCRS, geocentric
For a geocentric observer, prepare star-independent astrometry
parameters for transformations between ICRS and GCRS coordinates.
The caller supplies the date, and SOFA models are used to predict
the Earth ephemeris.
The parameters produced by this function are required in the
parallax, light deflection and aberration parts of the astrometric
transformation chain.
Given:
date1 float64 TDB as a 2-part...
date2 float64 ...Julian Date (Note 1)
Returned:
astrom ASTROM star-independent astrometry parameters:
pmt float64 PM time interval (SSB, Julian years)
eb [3]float64 SSB to observer (vector, au)
eh [3]float64 Sun to observer (unit vector)
em float64 distance from Sun to observer (au)
v [3]float64 barycentric observer velocity (vector, c)
bm1 float64 sqrt(1-|v|^2): reciprocal of Lorenz factor
bpn [3][3]float64 bias-precession-nutation matrix
along float64 unchanged
xpl float64 unchanged
ypl float64 unchanged
sphi float64 unchanged
cphi float64 unchanged
diurab float64 unchanged
eral float64 unchanged
refa float64 unchanged
refb float64 unchanged
Notes:
1. The TDB date date1+date2 is a Julian Date, apportioned in any
convenient way between the two arguments. For example,
JD(TDB)=2450123.7 could be expressed in any of these ways, among
others:
date1 date2
2450123.7 0.0 (JD method)
2451545.0 -1421.3 (J2000 method)
2400000.5 50123.2 (MJD method)
2450123.5 0.2 (date & time method)
The JD method is the most natural and convenient to use in cases
where the loss of several decimal digits of resolution is
acceptable. The J2000 method is best matched to the way the
argument is handled internally and will deliver the optimum
resolution. The MJD method and the date & time methods are both
good compromises between resolution and convenience. For most
applications of this function the choice will not be at all
critical.
TT can be used instead of TDB without any significant impact on
accuracy.
2. All the vectors are with respect to BCRS axes.
3. In cases where the caller wishes to supply his own Earth
ephemeris, the function Apcg can be used instead of the present
function.
4. This is one of several functions that inserts into the astrom
structure star-independent parameters needed for the chain of
astrometric transformations ICRS <-> GCRS <-> CIRS <-> observed.
The various functions support different classes of observer and
portions of the transformation chain:
functions observer transformation
Apcg Apcg13 geocentric ICRS <-> GCRS
Apci Apci13 terrestrial ICRS <-> CIRS
Apco Apco13 terrestrial ICRS <-> observed
Apcs Apcs13 space ICRS <-> GCRS
Aper Aper13 terrestrial update Earth rotation
Apio Apio13 terrestrial CIRS <-> observed
Those with names ending in "13" use contemporary SOFA models to
compute the various ephemerides. The others accept ephemerides
supplied by the caller.
The transformation from ICRS to GCRS covers space motion,
parallax, light deflection, and aberration. From GCRS to CIRS
comprises frame bias and precession-nutation. From CIRS to
observed takes account of Earth rotation, polar motion, diurnal
aberration and parallax (unless subsumed into the ICRS <-> GCRS
transformation), and atmospheric refraction.
5. The context structure astrom produced by this function is used by
Atciq* and Aticq*.
Called:
Epv00 Earth position and velocity
Apcg astrometry parameters, ICRS-GCRS, geocenter
*/
func Apcg13(date1, date2 float64, astrom *ASTROM) {
var ehpv, ebpv [2][3]float64
/* Earth barycentric & heliocentric position/velocity (au, au/d). */
Epv00(date1, date2, &ehpv, &ebpv)
/* Compute the star-independent astrometry parameters. */
Apcg(date1, date2, ebpv, ehpv[0], astrom)
}
/*
Apci Prepare for ICRS <-> CIRS, terrestrial, special
For a terrestrial observer, prepare star-independent astrometry
parameters for transformations between ICRS and geocentric CIRS
coordinates. The Earth ephemeris and CIP/CIO are supplied by the
caller.
The parameters produced by this function are required in the
parallax, light deflection, aberration, and bias-precession-nutation
parts of the astrometric transformation chain.
Given:
date1 float64 TDB as a 2-part...
date2 float64 ...Julian Date (Note 1)
ebpv [2][3]float64 Earth barycentric position/velocity (au, au/day)
ehp [3]float64 Earth heliocentric position (au)
x,y float64 CIP X,Y (components of unit vector)
s float64 the CIO locator s (radians)
Returned:
astrom ASTROM star-independent astrometry parameters:
pmt float64 PM time interval (SSB, Julian years)
eb [3]float64 SSB to observer (vector, au)
eh [3]float64 Sun to observer (unit vector)
em float64 distance from Sun to observer (au)
v [3]float64 barycentric observer velocity (vector, c)
bm1 float64 sqrt(1-|v|^2): reciprocal of Lorenz factor
bpn [3][3]float64 bias-precession-nutation matrix
along float64 unchanged
xpl float64 unchanged
ypl float64 unchanged
sphi float64 unchanged
cphi float64 unchanged
diurab float64 unchanged
eral float64 unchanged
refa float64 unchanged
refb float64 unchanged
Notes:
1. The TDB date date1+date2 is a Julian Date, apportioned in any
convenient way between the two arguments. For example,
JD(TDB)=2450123.7 could be expressed in any of these ways, among
others:
date1 date2
2450123.7 0.0 (JD method)
2451545.0 -1421.3 (J2000 method)
2400000.5 50123.2 (MJD method)
2450123.5 0.2 (date & time method)
The JD method is the most natural and convenient to use in cases
where the loss of several decimal digits of resolution is
acceptable. The J2000 method is best matched to the way the
argument is handled internally and will deliver the optimum
resolution. The MJD method and the date & time methods are both
good compromises between resolution and convenience. For most
applications of this function the choice will not be at all
critical.
TT can be used instead of TDB without any significant impact on
accuracy.
2. All the vectors are with respect to BCRS axes.
3. In cases where the caller does not wish to provide the Earth
ephemeris and CIP/CIO, the function Apci13 can be used instead
of the present function. This computes the required quantities
using other SOFA functions.
4. This is one of several functions that inserts into the astrom
structure star-independent parameters needed for the chain of
astrometric transformations ICRS <-> GCRS <-> CIRS <-> observed.
The various functions support different classes of observer and
portions of the transformation chain:
functions observer transformation
Apcg Apcg13 geocentric ICRS <-> GCRS
Apci Apci13 terrestrial ICRS <-> CIRS
Apco Apco13 terrestrial ICRS <-> observed
Apcs Apcs13 space ICRS <-> GCRS
Aper Aper13 terrestrial update Earth rotation
Apio Apio13 terrestrial CIRS <-> observed
Those with names ending in "13" use contemporary SOFA models to
compute the various ephemerides. The others accept ephemerides
supplied by the caller.
The transformation from ICRS to GCRS covers space motion,
parallax, light deflection, and aberration. From GCRS to CIRS
comprises frame bias and precession-nutation. From CIRS to
observed takes account of Earth rotation, polar motion, diurnal
aberration and parallax (unless subsumed into the ICRS <-> GCRS
transformation), and atmospheric refraction.
5. The context structure astrom produced by this function is used by
Atciq* and Aticq*.
Called:
Apcg astrometry parameters, ICRS-GCRS, geocenter
C2ixys celestial-to-intermediate matrix, given X,Y and s
*/
func Apci(date1, date2 float64, ebpv [2][3]float64, ehp [3]float64, x, y, s float64, astrom *ASTROM) {
/* Star-independent astrometry parameters for geocenter. */
Apcg(date1, date2, ebpv, ehp, astrom)
/* CIO based BPN matrix. */
C2ixys(x, y, s, &astrom.Bpn)
}
/*
Apci13 Prepare for ICRS <-> CIRS, terrestrial
For a terrestrial observer, prepare star-independent astrometry
parameters for transformations between ICRS and geocentric CIRS
coordinates. The caller supplies the date, and SOFA models are used
to predict the Earth ephemeris and CIP/CIO.
The parameters produced by this function are required in the
parallax, light deflection, aberration, and bias-precession-nutation
parts of the astrometric transformation chain.
Given:
date1 float64 TDB as a 2-part...
date2 float64 ...Julian Date (Note 1)
Returned:
astrom ASTROM star-independent astrometry parameters:
pmt float64 PM time interval (SSB, Julian years)
eb [3]float64 SSB to observer (vector, au)
eh [3]float64 Sun to observer (unit vector)
em float64 distance from Sun to observer (au)
v [3]float64 barycentric observer velocity (vector, c)
bm1 float64 sqrt(1-|v|^2): reciprocal of Lorenz factor
bpn [3][3]float64 bias-precession-nutation matrix
along float64 unchanged
xpl float64 unchanged
ypl float64 unchanged
sphi float64 unchanged
cphi float64 unchanged
diurab float64 unchanged
eral float64 unchanged
refa float64 unchanged
refb float64 unchanged
eo float64 equation of the origins (ERA-GST)
Notes:
1. The TDB date date1+date2 is a Julian Date, apportioned in any
convenient way between the two arguments. For example,
JD(TDB)=2450123.7 could be expressed in any of these ways, among
others:
date1 date2
2450123.7 0.0 (JD method)
2451545.0 -1421.3 (J2000 method)
2400000.5 50123.2 (MJD method)
2450123.5 0.2 (date & time method)
The JD method is the most natural and convenient to use in cases
where the loss of several decimal digits of resolution is
acceptable. The J2000 method is best matched to the way the
argument is handled internally and will deliver the optimum
resolution. The MJD method and the date & time methods are both
good compromises between resolution and convenience. For most
applications of this function the choice will not be at all
critical.
TT can be used instead of TDB without any significant impact on
accuracy.
2. All the vectors are with respect to BCRS axes.
3. In cases where the caller wishes to supply his own Earth
ephemeris and CIP/CIO, the function Apci can be used instead
of the present function.
4. This is one of several functions that inserts into the astrom
structure star-independent parameters needed for the chain of
astrometric transformations ICRS <-> GCRS <-> CIRS <-> observed.
The various functions support different classes of observer and
portions of the transformation chain:
functions observer transformation
Apcg Apcg13 geocentric ICRS <-> GCRS
Apci Apci13 terrestrial ICRS <-> CIRS
Apco Apco13 terrestrial ICRS <-> observed
Apcs Apcs13 space ICRS <-> GCRS
Aper Aper13 terrestrial update Earth rotation
Apio Apio13 terrestrial CIRS <-> observed
Those with names ending in "13" use contemporary SOFA models to
compute the various ephemerides. The others accept ephemerides
supplied by the caller.
The transformation from ICRS to GCRS covers space motion,
parallax, light deflection, and aberration. From GCRS to CIRS
comprises frame bias and precession-nutation. From CIRS to
observed takes account of Earth rotation, polar motion, diurnal
aberration and parallax (unless subsumed into the ICRS <-> GCRS
transformation), and atmospheric refraction.
5. The context structure astrom produced by this function is used by
Atciq* and Aticq*.
Called:
Epv00 Earth position and velocity
Pnm06a classical NPB matrix, IAU 2006/2000A
Bpn2xy extract CIP X,Y coordinates from NPB matrix
S06 the CIO locator s, given X,Y, IAU 2006
Apci astrometry parameters, ICRS-CIRS
Eors equation of the origins, given NPB matrix and s
*/
func Apci13(date1, date2 float64, astrom *ASTROM, eo *float64) {
var ehpv, ebpv [2][3]float64
var r [3][3]float64
var x, y, s float64
/* Earth barycentric & heliocentric position/velocity (au, au/d). */
Epv00(date1, date2, &ehpv, &ebpv)
/* Form the equinox based BPN matrix, IAU 2006/2000A. */
Pnm06a(date1, date2, &r)
/* Extract CIP X,Y. */
Bpn2xy(r, &x, &y)
/* Obtain CIO locator s. */
s = S06(date1, date2, x, y)
/* Compute the star-independent astrometry parameters. */
Apci(date1, date2, ebpv, ehpv[0], x, y, s, astrom)
/* Equation of the origins. */
*eo = Eors(r, s)
}
/*
Apco Prepare for ICRS <-> observed, terrestrial, special
For a terrestrial observer, prepare star-independent astrometry
parameters for transformations between ICRS and observed
coordinates. The caller supplies the Earth ephemeris, the Earth
rotation information and the refraction constants as well as the
site coordinates.
Given:
date1 float64 TDB as a 2-part...
date2 float64 ...Julian Date (Note 1)
ebpv [2][3]float64 Earth barycentric PV (au, au/day, Note 2)
ehp [3]float64 Earth heliocentric P (au, Note 2)
x,y float64 CIP X,Y (components of unit vector)
s float64 the CIO locator s (radians)
theta float64 Earth rotation angle (radians)
elong float64 longitude (radians, east +ve, Note 3)
phi float64 latitude (geodetic, radians, Note 3)
hm float64 height above ellipsoid (m, geodetic, Note 3)
xp,yp float64 polar motion coordinates (radians, Note 4)
sp float64 the TIO locator s' (radians, Note 4)
refa float64 refraction constant A (radians, Note 5)
refb float64 refraction constant B (radians, Note 5)
Returned:
astrom ASTROM star-independent astrometry parameters:
pmt float64 PM time interval (SSB, Julian years)
eb [3]float64 SSB to observer (vector, au)
eh [3]float64 Sun to observer (unit vector)
em float64 distance from Sun to observer (au)
v [3]float64 barycentric observer velocity (vector, c)
bm1 float64 sqrt(1-|v|^2): reciprocal of Lorenz factor
bpn [3][3]float64 bias-precession-nutation matrix
along float64 adjusted longitude (radians)
xpl float64 polar motion xp wrt local meridian (radians)
ypl float64 polar motion yp wrt local meridian (radians)
sphi float64 sine of geodetic latitude
cphi float64 cosine of geodetic latitude
diurab float64 magnitude of diurnal aberration vector
eral float64 "local" Earth rotation angle (radians)
refa float64 refraction constant A (radians)
refb float64 refraction constant B (radians)
Notes:
1. The TDB date date1+date2 is a Julian Date, apportioned in any
convenient way between the two arguments. For example,
JD(TDB)=2450123.7 could be expressed in any of these ways, among
others:
date1 date2
2450123.7 0.0 (JD method)
2451545.0 -1421.3 (J2000 method)
2400000.5 50123.2 (MJD method)
2450123.5 0.2 (date & time method)
The JD method is the most natural and convenient to use in cases
where the loss of several decimal digits of resolution is
acceptable. The J2000 method is best matched to the way the
argument is handled internally and will deliver the optimum
resolution. The MJD method and the date & time methods are both
good compromises between resolution and convenience. For most
applications of this function the choice will not be at all
critical.
TT can be used instead of TDB without any significant impact on
accuracy.
2. The vectors eb, eh, and all the astrom vectors, are with respect
to BCRS axes.
3. The geographical coordinates are with respect to the WGS84
reference ellipsoid. TAKE CARE WITH THE LONGITUDE SIGN
CONVENTION: the longitude required by the present function is
right-handed, i.e. east-positive, in accordance with geographical
convention.
The adjusted longitude stored in the astrom array takes into
account the TIO locator and polar motion.
4. xp and yp are the coordinates (in radians) of the Celestial
Intermediate Pole with respect to the International Terrestrial
Reference System (see IERS Conventions), measured along the
meridians 0 and 90 deg west respectively. sp is the TIO locator
s', in radians, which positions the Terrestrial Intermediate
Origin on the equator. For many applications, xp, yp and
(especially) sp can be set to zero.
Internally, the polar motion is stored in a form rotated onto the
local meridian.
5. The refraction constants refa and refb are for use in a
dZ = A*tan(Z)+B*tan^3(Z) model, where Z is the observed
(i.e. refracted) zenith distance and dZ is the amount of
refraction.
6. It is advisable to take great care with units, as even unlikely
values of the input parameters are accepted and processed in
accordance with the models used.
7. In cases where the caller does not wish to provide the Earth
Ephemeris, the Earth rotation information and refraction
constants, the function Apco13 can be used instead of the
present function. This starts from UTC and weather readings etc.
and computes suitable values using other SOFA functions.
8. This is one of several functions that inserts into the astrom
structure star-independent parameters needed for the chain of
astrometric transformations ICRS <-> GCRS <-> CIRS <-> observed.
The various functions support different classes of observer and
portions of the transformation chain:
functions observer transformation
Apcg Apcg13 geocentric ICRS <-> GCRS
Apci Apci13 terrestrial ICRS <-> CIRS
Apco Apco13 terrestrial ICRS <-> observed
Apcs Apcs13 space ICRS <-> GCRS
Aper Aper13 terrestrial update Earth rotation
Apio Apio13 terrestrial CIRS <-> observed
Those with names ending in "13" use contemporary SOFA models to
compute the various ephemerides. The others accept ephemerides
supplied by the caller.
The transformation from ICRS to GCRS covers space motion,
parallax, light deflection, and aberration. From GCRS to CIRS
comprises frame bias and precession-nutation. From CIRS to
observed takes account of Earth rotation, polar motion, diurnal
aberration and parallax (unless subsumed into the ICRS <-> GCRS
transformation), and atmospheric refraction.
9. The context structure astrom produced by this function is used by
Atioq, Atoiq, Atciq* and Aticq*.
Called:
Ir initialize r-matrix to identity
Rz rotate around Z-axis
Ry rotate around Y-axis
Rx rotate around X-axis
Anpm normalize angle into range +/- pi
C2ixys celestial-to-intermediate matrix, given X,Y and s
Pvtob position/velocity of terrestrial station
Trxpv product of transpose of r-matrix and pv-vector
Apcs astrometry parameters, ICRS-GCRS, space observer
Cr copy r-matrix
*/
func Apco(date1, date2 float64, ebpv [2][3]float64, ehp [3]float64, x, y, s float64, theta,
elong, phi float64, hm, xp, yp, sp float64, refa, refb float64, astrom *ASTROM) {
var r [3][3]float64
var a, b, eral, c float64
var pvc, pv [2][3]float64
/* Form the rotation matrix, CIRS to apparent [HA,Dec]. */
Ir(&r)
Rz(theta+sp, &r)
Ry(-xp, &r)
Rx(-yp, &r)
Rz(elong, &r)
/* Solve for local Earth rotation angle. */
a = r[0][0]
b = r[0][1]
// eral = ( a != 0.0 || b != 0.0 ) ? atan2(b, a) : 0.0;
if a != 0.0 || b != 0.0 {
eral = atan2(b, a)
} else {
eral = 0.0
}
astrom.Eral = eral
/* Solve for polar motion [X,Y] with respect to local meridian. */
a = r[0][0]
c = r[0][2]
astrom.Xpl = atan2(c, sqrt(a*a+b*b))
a = r[1][2]
b = r[2][2]
// astrom.Ypl = ( a != 0.0 || b != 0.0 ) ? -atan2(a, b) : 0.0;
if a != 0.0 || b != 0.0 {
astrom.Ypl = -atan2(a, b)
} else {
astrom.Ypl = 0.0
}
/* Adjusted longitude. */
astrom.Along = Anpm(eral - theta)
/* Functions of latitude. */
astrom.Sphi = sin(phi)
astrom.Cphi = cos(phi)
/* Refraction constants. */
astrom.Refa = refa
astrom.Refb = refb
/* Disable the (redundant) diurnal aberration step. */
astrom.Diurab = 0.0
/* CIO based BPN matrix. */
C2ixys(x, y, s, &r)
/* Observer's geocentric position and velocity (m, m/s, CIRS). */
Pvtob(elong, phi, hm, xp, yp, sp, theta, &pvc)
/* Rotate into GCRS. */
Trxpv(r, pvc, &pv)
/* ICRS <-> GCRS parameters. */
Apcs(date1, date2, pv, ebpv, ehp, astrom)
/* Store the CIO based BPN matrix. */
Cr(r, &astrom.Bpn)
}
/*
Apco13 Prepare for ICRS <-> observed, terrestrial
For a terrestrial observer, prepare star-independent astrometry
parameters for transformations between ICRS and observed
coordinates. The caller supplies UTC, site coordinates, ambient air
conditions and observing wavelength, and SOFA models are used to
obtain the Earth ephemeris, CIP/CIO and refraction constants.
The parameters produced by this function are required in the
parallax, light deflection, aberration, and bias-precession-nutation
parts of the ICRS/CIRS transformations.
Given:
utc1 float64 UTC as a 2-part...
utc2 float64 ...quasi Julian Date (Notes 1,2)
dut1 float64 UT1-UTC (seconds, Note 3)
elong float64 longitude (radians, east +ve, Note 4)
phi float64 latitude (geodetic, radians, Note 4)
hm float64 height above ellipsoid (m, geodetic, Notes 4,6)
xp,yp float64 polar motion coordinates (radians, Note 5)
phpa float64 pressure at the observer (hPa = mB, Note 6)
tc float64 ambient temperature at the observer (deg C)
rh float64 relative humidity at the observer (range 0-1)
wl float64 wavelength (micrometers, Note 7)
Returned:
astrom ASTROM star-independent astrometry parameters:
pmt float64 PM time interval (SSB, Julian years)
eb [3]float64 SSB to observer (vector, au)
eh [3]float64 Sun to observer (unit vector)
em float64 distance from Sun to observer (au)
v [3]float64 barycentric observer velocity (vector, c)
bm1 float64 sqrt(1-|v|^2): reciprocal of Lorenz factor
bpn [3][3]float64 bias-precession-nutation matrix
along float64 adjusted longitude (radians)
xpl float64 polar motion xp wrt local meridian (radians)
ypl float64 polar motion yp wrt local meridian (radians)
sphi float64 sine of geodetic latitude
cphi float64 cosine of geodetic latitude
diurab float64 magnitude of diurnal aberration vector
eral float64 "local" Earth rotation angle (radians)
refa float64 refraction constant A (radians)
refb float64 refraction constant B (radians)
eo float64 equation of the origins (ERA-GST)
Returned (function value):
int status: +1 = dubious year (Note 2)
0 = OK
-1 = unacceptable date
Notes:
1. utc1+utc2 is quasi Julian Date (see Note 2), apportioned in any
convenient way between the two arguments, for example where utc1
is the Julian Day Number and utc2 is the fraction of a day.
However, JD cannot unambiguously represent UTC during a leap
second unless special measures are taken. The convention in the
present function is that the JD day represents UTC days whether
the length is 86399, 86400 or 86401 SI seconds.
Applications should use the function Dtf2d to convert from
calendar date and time of day into 2-part quasi Julian Date, as
it implements the leap-second-ambiguity convention just
described.
2. The warning status "dubious year" flags UTCs that predate the
introduction of the time scale or that are too far in the
future to be trusted. See Dat for further details.
3. UT1-UTC is tabulated in IERS bulletins. It increases by exactly
one second at the end of each positive UTC leap second,
introduced in order to keep UT1-UTC within +/- 0.9s. n.b. This
practice is under review, and in the future UT1-UTC may grow
essentially without limit.
4. The geographical coordinates are with respect to the WGS84
reference ellipsoid. TAKE CARE WITH THE LONGITUDE SIGN: the
longitude required by the present function is east-positive
(i.e. right-handed), in accordance with geographical convention.
5. The polar motion xp,yp can be obtained from IERS bulletins. The
values are the coordinates (in radians) of the Celestial
Intermediate Pole with respect to the International Terrestrial
Reference System (see IERS Conventions 2003), measured along the
meridians 0 and 90 deg west respectively. For many
applications, xp and yp can be set to zero.
Internally, the polar motion is stored in a form rotated onto
the local meridian.
6. If hm, the height above the ellipsoid of the observing station
in meters, is not known but phpa, the pressure in hPa (=mB), is
available, an adequate estimate of hm can be obtained from the
expression
hm = -29.3 * tsl * log ( phpa / 1013.25 );
where tsl is the approximate sea-level air temperature in K
(See Astrophysical Quantities, C.W.Allen, 3rd edition, section
52). Similarly, if the pressure phpa is not known, it can be
estimated from the height of the observing station, hm, as
follows:
phpa = 1013.25 * exp ( -hm / ( 29.3 * tsl ) );
Note, however, that the refraction is nearly proportional to
the pressure and that an accurate phpa value is important for
precise work.
7. The argument wl specifies the observing wavelength in
micrometers. The transition from optical to radio is assumed to
occur at 100 micrometers (about 3000 GHz).
8. It is advisable to take great care with units, as even unlikely
values of the input parameters are accepted and processed in
accordance with the models used.
9. In cases where the caller wishes to supply his own Earth
ephemeris, Earth rotation information and refraction constants,
the function Apco can be used instead of the present function.
10)This is one of several functions that inserts into the astrom
structure star-independent parameters needed for the chain of
astrometric transformations ICRS <-> GCRS <-> CIRS <-> observed.
The various functions support different classes of observer and
portions of the transformation chain:
functions observer transformation
Apcg Apcg13 geocentric ICRS <-> GCRS
Apci Apci13 terrestrial ICRS <-> CIRS
Apco Apco13 terrestrial ICRS <-> observed
Apcs Apcs13 space ICRS <-> GCRS
Aper Aper13 terrestrial update Earth rotation
Apio Apio13 terrestrial CIRS <-> observed
Those with names ending in "13" use contemporary SOFA models to
compute the various ephemerides. The others accept ephemerides
supplied by the caller.
The transformation from ICRS to GCRS covers space motion,
parallax, light deflection, and aberration. From GCRS to CIRS
comprises frame bias and precession-nutation. From CIRS to
observed takes account of Earth rotation, polar motion, diurnal
aberration and parallax (unless subsumed into the ICRS <-> GCRS
transformation), and atmospheric refraction.
11)The context structure astrom produced by this function is used
by Atioq, Atoiq, Atciq* and Aticq*.
Called:
Utctai UTC to TAI
Taitt TAI to TT
Utcut1 UTC to UT1
Epv00 Earth position and velocity
Pnm06a classical NPB matrix, IAU 2006/2000A
Bpn2xy extract CIP X,Y coordinates from NPB matrix
S06 the CIO locator s, given X,Y, IAU 2006
Era00 Earth rotation angle, IAU 2000
Sp00 the TIO locator s', IERS 2000
Refco refraction constants for given ambient conditions
Apco astrometry parameters, ICRS-observed
Eors equation of the origins, given NPB matrix and s
*/
func Apco13(utc1, utc2, dut1 float64, elong, phi, hm, xp, yp float64, phpa, tc, rh, wl float64, astrom *ASTROM, eo *float64) int {
var j int
var tai1, tai2, tt1, tt2, ut11, ut12 float64
var ehpv, ebpv [2][3]float64
var r [3][3]float64
var x, y, s, theta, sp, refa, refb float64
/* UTC to other time scales. */
j = Utctai(utc1, utc2, &tai1, &tai2)
if j < 0 {
return -1
}
Taitt(tai1, tai2, &tt1, &tt2)
j = Utcut1(utc1, utc2, dut1, &ut11, &ut12)
if j < 0 {
return -1
}
/* Earth barycentric & heliocentric position/velocity (au, au/d). */
Epv00(tt1, tt2, &ehpv, &ebpv)
/* Form the equinox based BPN matrix, IAU 2006/2000A. */
Pnm06a(tt1, tt2, &r)
/* Extract CIP X,Y. */
Bpn2xy(r, &x, &y)
/* Obtain CIO locator s. */
s = S06(tt1, tt2, x, y)
/* Earth rotation angle. */
theta = Era00(ut11, ut12)
/* TIO locator s'. */
sp = Sp00(tt1, tt2)
/* Refraction constants A and B. */
Refco(phpa, tc, rh, wl, &refa, &refb)
/* Compute the star-independent astrometry parameters. */
Apco(tt1, tt2, ebpv, ehpv[0], x, y, s, theta,
elong, phi, hm, xp, yp, sp, refa, refb, astrom)
/* Equation of the origins. */
*eo = Eors(r, s)
/* Return any warning status. */
return j
}
/*
Apcs Prepare for ICRS <-> CIRS, space, special
For an observer whose geocentric position and velocity are known,
prepare star-independent astrometry parameters for transformations
between ICRS and GCRS. The Earth ephemeris is supplied by the