@@ -468,6 +468,14 @@ def julian_ephemeris_millennium(julian_ephemeris_century):
468
468
return jme
469
469
470
470
471
+ def mult_cos_add_mult (a , b , c , jme , scratch ):
472
+ np .multiply (c , jme , out = scratch )
473
+ np .add (b , scratch , out = scratch )
474
+ np .cos (scratch , out = scratch )
475
+ np .multiply (a , scratch , out = scratch )
476
+ return scratch
477
+
478
+
471
479
@jcompile ('float64(float64)' , nopython = True )
472
480
def heliocentric_longitude (jme ):
473
481
l0 = 0.0
@@ -477,36 +485,38 @@ def heliocentric_longitude(jme):
477
485
l4 = 0.0
478
486
l5 = 0.0
479
487
488
+ scratch = np .empty_like (jme )
489
+
480
490
for row in range (64 ):
481
- l0 += (HELIO_LONG_TABLE [0 , row , 0 ]
482
- * np . cos ( HELIO_LONG_TABLE [0 , row , 1 ]
483
- + HELIO_LONG_TABLE [0 , row , 2 ] * jme )
484
- )
491
+ l0 += mult_cos_add_mult (HELIO_LONG_TABLE [0 , row , 0 ],
492
+ HELIO_LONG_TABLE [0 , row , 1 ],
493
+ HELIO_LONG_TABLE [0 , row , 2 ],
494
+ jme , scratch )
485
495
for row in range (34 ):
486
- l1 += (HELIO_LONG_TABLE [1 , row , 0 ]
487
- * np . cos ( HELIO_LONG_TABLE [1 , row , 1 ]
488
- + HELIO_LONG_TABLE [1 , row , 2 ] * jme )
489
- )
496
+ l1 += mult_cos_add_mult (HELIO_LONG_TABLE [1 , row , 0 ],
497
+ HELIO_LONG_TABLE [1 , row , 1 ],
498
+ HELIO_LONG_TABLE [1 , row , 2 ],
499
+ jme , scratch )
490
500
for row in range (20 ):
491
- l2 += (HELIO_LONG_TABLE [2 , row , 0 ]
492
- * np . cos ( HELIO_LONG_TABLE [2 , row , 1 ]
493
- + HELIO_LONG_TABLE [2 , row , 2 ] * jme )
494
- )
501
+ l2 += mult_cos_add_mult (HELIO_LONG_TABLE [2 , row , 0 ],
502
+ HELIO_LONG_TABLE [2 , row , 1 ],
503
+ HELIO_LONG_TABLE [2 , row , 2 ],
504
+ jme , scratch )
495
505
for row in range (7 ):
496
- l3 += (HELIO_LONG_TABLE [3 , row , 0 ]
497
- * np . cos ( HELIO_LONG_TABLE [3 , row , 1 ]
498
- + HELIO_LONG_TABLE [3 , row , 2 ] * jme )
499
- )
506
+ l3 += mult_cos_add_mult (HELIO_LONG_TABLE [3 , row , 0 ],
507
+ HELIO_LONG_TABLE [3 , row , 1 ],
508
+ HELIO_LONG_TABLE [3 , row , 2 ],
509
+ jme , scratch )
500
510
for row in range (3 ):
501
- l4 += (HELIO_LONG_TABLE [4 , row , 0 ]
502
- * np . cos ( HELIO_LONG_TABLE [4 , row , 1 ]
503
- + HELIO_LONG_TABLE [4 , row , 2 ] * jme )
504
- )
511
+ l4 += mult_cos_add_mult (HELIO_LONG_TABLE [4 , row , 0 ],
512
+ HELIO_LONG_TABLE [4 , row , 1 ],
513
+ HELIO_LONG_TABLE [4 , row , 2 ],
514
+ jme , scratch )
505
515
for row in range (1 ):
506
- l5 += (HELIO_LONG_TABLE [5 , row , 0 ]
507
- * np . cos ( HELIO_LONG_TABLE [5 , row , 1 ]
508
- + HELIO_LONG_TABLE [5 , row , 2 ] * jme )
509
- )
516
+ l5 += mult_cos_add_mult (HELIO_LONG_TABLE [5 , row , 0 ],
517
+ HELIO_LONG_TABLE [5 , row , 1 ],
518
+ HELIO_LONG_TABLE [5 , row , 2 ],
519
+ jme , scratch )
510
520
511
521
l_rad = (l0 + l1 * jme + l2 * jme ** 2 + l3 * jme ** 3 + l4 * jme ** 4 +
512
522
l5 * jme ** 5 )/ 10 ** 8
@@ -518,16 +528,18 @@ def heliocentric_longitude(jme):
518
528
def heliocentric_latitude (jme ):
519
529
b0 = 0.0
520
530
b1 = 0.0
531
+ scratch = np .empty_like (jme )
532
+
521
533
for row in range (5 ):
522
- b0 += (HELIO_LAT_TABLE [0 , row , 0 ]
523
- * np . cos ( HELIO_LAT_TABLE [0 , row , 1 ]
524
- + HELIO_LAT_TABLE [0 , row , 2 ] * jme )
525
- )
534
+ b0 += mult_cos_add_mult (HELIO_LAT_TABLE [0 , row , 0 ],
535
+ HELIO_LAT_TABLE [0 , row , 1 ],
536
+ HELIO_LAT_TABLE [0 , row , 2 ],
537
+ jme , scratch )
526
538
for row in range (2 ):
527
- b1 += (HELIO_LAT_TABLE [1 , row , 0 ]
528
- * np . cos ( HELIO_LAT_TABLE [1 , row , 1 ]
529
- + HELIO_LAT_TABLE [1 , row , 2 ] * jme )
530
- )
539
+ b1 += mult_cos_add_mult (HELIO_LAT_TABLE [1 , row , 0 ],
540
+ HELIO_LAT_TABLE [1 , row , 1 ],
541
+ HELIO_LAT_TABLE [1 , row , 2 ],
542
+ jme , scratch )
531
543
532
544
b_rad = (b0 + b1 * jme )/ 10 ** 8
533
545
b = np .rad2deg (b_rad )
@@ -541,31 +553,33 @@ def heliocentric_radius_vector(jme):
541
553
r2 = 0.0
542
554
r3 = 0.0
543
555
r4 = 0.0
556
+ scratch = np .empty_like (jme )
557
+
544
558
for row in range (40 ):
545
- r0 += (HELIO_RADIUS_TABLE [0 , row , 0 ]
546
- * np . cos ( HELIO_RADIUS_TABLE [0 , row , 1 ]
547
- + HELIO_RADIUS_TABLE [0 , row , 2 ] * jme )
548
- )
559
+ r0 += mult_cos_add_mult (HELIO_RADIUS_TABLE [0 , row , 0 ],
560
+ HELIO_RADIUS_TABLE [0 , row , 1 ],
561
+ HELIO_RADIUS_TABLE [0 , row , 2 ],
562
+ jme , scratch )
549
563
for row in range (10 ):
550
- r1 += (HELIO_RADIUS_TABLE [1 , row , 0 ]
551
- * np . cos ( HELIO_RADIUS_TABLE [1 , row , 1 ]
552
- + HELIO_RADIUS_TABLE [1 , row , 2 ] * jme )
553
- )
564
+ r1 += mult_cos_add_mult (HELIO_RADIUS_TABLE [1 , row , 0 ],
565
+ HELIO_RADIUS_TABLE [1 , row , 1 ],
566
+ HELIO_RADIUS_TABLE [1 , row , 2 ],
567
+ jme , scratch )
554
568
for row in range (6 ):
555
- r2 += (HELIO_RADIUS_TABLE [2 , row , 0 ]
556
- * np . cos ( HELIO_RADIUS_TABLE [2 , row , 1 ]
557
- + HELIO_RADIUS_TABLE [2 , row , 2 ] * jme )
558
- )
569
+ r2 += mult_cos_add_mult (HELIO_RADIUS_TABLE [2 , row , 0 ],
570
+ HELIO_RADIUS_TABLE [2 , row , 1 ],
571
+ HELIO_RADIUS_TABLE [2 , row , 2 ],
572
+ jme , scratch )
559
573
for row in range (2 ):
560
- r3 += (HELIO_RADIUS_TABLE [3 , row , 0 ]
561
- * np . cos ( HELIO_RADIUS_TABLE [3 , row , 1 ]
562
- + HELIO_RADIUS_TABLE [3 , row , 2 ] * jme )
563
- )
574
+ r3 += mult_cos_add_mult (HELIO_RADIUS_TABLE [3 , row , 0 ],
575
+ HELIO_RADIUS_TABLE [3 , row , 1 ],
576
+ HELIO_RADIUS_TABLE [3 , row , 2 ],
577
+ jme , scratch )
564
578
for row in range (1 ):
565
- r4 += (HELIO_RADIUS_TABLE [4 , row , 0 ]
566
- * np . cos ( HELIO_RADIUS_TABLE [4 , row , 1 ]
567
- + HELIO_RADIUS_TABLE [4 , row , 2 ] * jme )
568
- )
579
+ r4 += mult_cos_add_mult (HELIO_RADIUS_TABLE [4 , row , 0 ],
580
+ HELIO_RADIUS_TABLE [4 , row , 1 ],
581
+ HELIO_RADIUS_TABLE [4 , row , 2 ],
582
+ jme , scratch )
569
583
570
584
r = (r0 + r1 * jme + r2 * jme ** 2 + r3 * jme ** 3 + r4 * jme ** 4 )/ 10 ** 8
571
585
return r
0 commit comments