1 /* 2 3 -Procedure spkpos_c ( S/P Kernel, position ) 4 5 -Abstract 6 7 Return the position of a target body relative to an observing 8 body, optionally corrected for light time (planetary aberration) 9 and stellar aberration. 10 11 -Disclaimer 12 13 THIS SOFTWARE AND ANY RELATED MATERIALS WERE CREATED BY THE 14 CALIFORNIA INSTITUTE OF TECHNOLOGY (CALTECH) UNDER A U.S. 15 GOVERNMENT CONTRACT WITH THE NATIONAL AERONAUTICS AND SPACE 16 ADMINISTRATION (NASA). THE SOFTWARE IS TECHNOLOGY AND SOFTWARE 17 PUBLICLY AVAILABLE UNDER U.S. EXPORT LAWS AND IS PROVIDED "AS-IS" 18 TO THE RECIPIENT WITHOUT WARRANTY OF ANY KIND, INCLUDING ANY 19 WARRANTIES OF PERFORMANCE OR MERCHANTABILITY OR FITNESS FOR A 20 PARTICULAR USE OR PURPOSE (AS SET FORTH IN UNITED STATES UCC 21 SECTIONS 2312-2313) OR FOR ANY PURPOSE WHATSOEVER, FOR THE 22 SOFTWARE AND RELATED MATERIALS, HOWEVER USED. 23 24 IN NO EVENT SHALL CALTECH, ITS JET PROPULSION LABORATORY, OR NASA 25 BE LIABLE FOR ANY DAMAGES AND/OR COSTS, INCLUDING, BUT NOT 26 LIMITED TO, INCIDENTAL OR CONSEQUENTIAL DAMAGES OF ANY KIND, 27 INCLUDING ECONOMIC DAMAGE OR INJURY TO PROPERTY AND LOST PROFITS, 28 REGARDLESS OF WHETHER CALTECH, JPL, OR NASA BE ADVISED, HAVE 29 REASON TO KNOW, OR, IN FACT, SHALL KNOW OF THE POSSIBILITY. 30 31 RECIPIENT BEARS ALL RISK RELATING TO QUALITY AND PERFORMANCE OF 32 THE SOFTWARE AND ANY RELATED MATERIALS, AND AGREES TO INDEMNIFY 33 CALTECH AND NASA FOR ALL THIRD-PARTY CLAIMS RESULTING FROM THE 34 ACTIONS OF RECIPIENT IN THE USE OF THE SOFTWARE. 35 36 -Required_Reading 37 38 SPK 39 NAIF_IDS 40 FRAMES 41 TIME 42 43 -Keywords 44 45 EPHEMERIS 46 47 */ 48 49 #include "SpiceUsr.h" 50 #include "SpiceZfc.h" 51 #include "SpiceZmc.h" 52 53 spkpos_c(ConstSpiceChar * targ,SpiceDouble et,ConstSpiceChar * ref,ConstSpiceChar * abcorr,ConstSpiceChar * obs,SpiceDouble ptarg[3],SpiceDouble * lt)54 void spkpos_c ( ConstSpiceChar * targ, 55 SpiceDouble et, 56 ConstSpiceChar * ref, 57 ConstSpiceChar * abcorr, 58 ConstSpiceChar * obs, 59 SpiceDouble ptarg[3], 60 SpiceDouble * lt ) 61 /* 62 63 -Brief_I/O 64 65 Variable I/O Description 66 -------- --- -------------------------------------------------- 67 targ I Target body name. 68 et I Observer epoch. 69 ref I Reference frame of output position vector. 70 abcorr I Aberration correction flag. 71 obs I Observing body name. 72 ptarg O Position of target. 73 lt O One way light time between observer and target. 74 75 -Detailed_Input 76 77 targ is the name of a target body. Optionally, you may 78 supply the integer ID code for the object as 79 an integer string. For example both "MOON" and 80 "301" are legitimate strings that indicate the 81 moon is the target body. 82 83 The target and observer define a position vector 84 which points from the observer to the target. 85 86 et is the ephemeris time, expressed as seconds past 87 J2000 TDB, at which the position of the target body 88 relative to the observer is to be computed. `et' 89 refers to time at the observer's location. 90 91 ref is the name of the reference frame relative to which 92 the output position vector should be expressed. This 93 may be any frame supported by the SPICE system, 94 including built-in frames (documented in the Frames 95 Required Reading) and frames defined by a loaded 96 frame kernel (FK). 97 98 When `ref' designates a non-inertial frame, the 99 orientation of the frame is evaluated at an epoch 100 dependent on the selected aberration correction. See 101 the description of the output position vector `ptarg' 102 for details. 103 104 abcorr indicates the aberration corrections to be applied to 105 the position of the target body to account for 106 one-way light time and stellar aberration. See the 107 discussion in the Particulars section for 108 recommendations on how to choose aberration 109 corrections. 110 111 'abcorr' may be any of the following: 112 113 "NONE" Apply no correction. Return the 114 geometric position of the target body 115 relative to the observer. 116 117 The following values of 'abcorr' apply to the 118 "reception" case in which photons depart from the 119 target's location at the light-time corrected epoch 120 et-lt and *arrive* at the observer's location at `et': 121 122 "LT" Correct for one-way light time (also 123 called "planetary aberration") using a 124 Newtonian formulation. This correction 125 yields the position of the target at 126 the moment it emitted photons arriving 127 at the observer at `et'. 128 129 The light time correction uses an 130 iterative solution of the light time 131 equation (see Particulars for details). 132 The solution invoked by the "LT" option 133 uses one iteration. 134 135 "LT+S" Correct for one-way light time and 136 stellar aberration using a Newtonian 137 formulation. This option modifies the 138 position obtained with the "LT" option 139 to account for the observer's velocity 140 relative to the solar system 141 barycenter. The result is the apparent 142 position of the target---the position 143 as seen by the observer. 144 145 "CN" Converged Newtonian light time 146 correction. In solving the light time 147 equation, the "CN" correction iterates 148 until the solution converges (three 149 iterations on all supported platforms). 150 Whether the "CN+S" solution is 151 substantially more accurate than the 152 "LT" solution depends on the geometry 153 of the participating objects and on the 154 accuracy of the input data. In all 155 cases this routine will execute more 156 slowly when a converged solution is 157 computed. See the Particulars section 158 below for a discussion of precision of 159 light time corrections. 160 161 "CN+S" Converged Newtonian light time 162 correction and stellar aberration 163 correction. 164 165 166 The following values of 'abcorr' apply to the 167 "transmission" case in which photons *depart* from 168 the observer's location at `et' and arrive at the 169 target's location at the light-time corrected epoch 170 et+lt: 171 172 "XLT" "Transmission" case: correct for 173 one-way light time using a Newtonian 174 formulation. This correction yields the 175 position of the target at the moment it 176 receives photons emitted from the 177 observer's location at `et'. 178 179 "XLT+S" "Transmission" case: correct for one-way 180 light time and stellar aberration using a 181 Newtonian formulation. This option 182 modifies the position obtained with the 183 "XLT" option to account for the observer's 184 velocity relative to the solar system 185 barycenter. The computed target position 186 indicates the direction that photons 187 emitted from the observer's location must 188 be "aimed" to hit the target. 189 190 "XCN" "Transmission" case: converged 191 Newtonian light time correction. 192 193 "XCN+S" "Transmission" case: converged Newtonian 194 light time correction and stellar 195 aberration correction. 196 197 198 Neither special nor general relativistic effects are 199 accounted for in the aberration corrections applied 200 by this routine. 201 202 Case and blanks are not significant in the string 203 'abcorr'. 204 205 obs is the name of an observing body. Optionally, you may 206 supply the ID code of the object as an integer string. 207 For example, both "EARTH" and "399" are legitimate 208 strings to supply to indicate the observer is 209 Earth. 210 211 -Detailed_Output 212 213 ptarg is a Cartesian 3-vector representing the position of 214 the target body relative to the specified observer. 215 `ptarg' is corrected for the specified aberrations, and 216 is expressed with respect to the reference frame 217 specified by `ref'. The three components of `ptarg' 218 represent the x-, y- and z-components of the target's 219 position. 220 221 Units are always km. 222 223 `ptarg' points from the observer's location at `et' to 224 the aberration-corrected location of the target. 225 Note that the sense of this position vector is 226 independent of the direction of radiation travel 227 implied by the aberration correction. 228 229 Non-inertial frames are treated as follows: letting 230 ltcent be the one-way light time between the observer 231 and the central body associated with the frame, the 232 orientation of the frame is evaluated at et-ltcent, 233 et+ltcent, or `et' depending on whether the requested 234 aberration correction is, respectively, for received 235 radiation, transmitted radiation, or is omitted. ltcent 236 is computed using the method indicated by 'abcorr'. 237 238 lt is the one-way light time between the observer and 239 target in seconds. If the target position is 240 corrected for aberrations, then `lt' is the one-way 241 light time between the observer and the light time 242 corrected target location. 243 244 -Parameters 245 246 None. 247 248 -Exceptions 249 250 1) If name of target or observer cannot be translated to its 251 NAIF ID code, the error SPICE(IDCODENOTFOUND) is signaled. 252 253 2) If the reference frame `ref' is not a recognized reference 254 frame the error SPICE(UNKNOWNFRAME) is signaled. 255 256 3) If the loaded kernels provide insufficient data to 257 compute the requested position vector, the deficiency will 258 be diagnosed by a routine in the call tree of this routine. 259 260 4) If an error occurs while reading an SPK or other kernel file, 261 the error will be diagnosed by a routine in the call tree 262 of this routine. 263 264 -Files 265 266 This routine computes positions using SPK files that have been 267 loaded into the SPICE system, normally via the kernel loading 268 interface routine furnsh_c. See the routine furnsh_c and the SPK 269 and KERNEL Required Reading for further information on loading 270 (and unloading) kernels. 271 272 If the output position `ptarg' is to be expressed relative to a 273 non-inertial frame, or if any of the ephemeris data used to 274 compute `ptarg' are expressed relative to a non-inertial frame in 275 the SPK files providing those data, additional kernels may be 276 needed to enable the reference frame transformations required to 277 compute the position. These additional kernels may be C-kernels, PCK 278 files or frame kernels. Any such kernels must already be loaded 279 at the time this routine is called. 280 281 -Particulars 282 283 This routine is part of the user interface to the SPICE ephemeris 284 system. It allows you to retrieve position information for any 285 ephemeris object relative to any other in a reference frame that 286 is convenient for further computations. 287 288 This routine is identical in function to the routine SPKEZP 289 except that it allows you to refer to ephemeris objects by name 290 (via a character string). 291 292 293 Aberration corrections 294 ====================== 295 296 In space science or engineering applications one frequently 297 wishes to know where to point a remote sensing instrument, such 298 as an optical camera or radio antenna, in order to observe or 299 otherwise receive radiation from a target. This pointing problem 300 is complicated by the finite speed of light: one needs to point 301 to where the target appears to be as opposed to where it actually 302 is at the epoch of observation. We use the adjectives 303 "geometric," "uncorrected," or "true" to refer to an actual 304 position or state of a target at a specified epoch. When a 305 geometric position or state vector is modified to reflect how it 306 appears to an observer, we describe that vector by any of the 307 terms "apparent," "corrected," "aberration corrected," or "light 308 time and stellar aberration corrected." The SPICE Toolkit can 309 correct for two phenomena affecting the apparent location of an 310 object: one-way light time (also called "planetary aberration") and 311 stellar aberration. 312 313 One-way light time 314 ------------------ 315 316 Correcting for one-way light time is done by computing, given an 317 observer and observation epoch, where a target was when the observed 318 photons departed the target's location. The vector from the 319 observer to this computed target location is called a "light time 320 corrected" vector. The light time correction depends on the motion 321 of the target relative to the solar system barycenter, but it is 322 independent of the velocity of the observer relative to the solar 323 system barycenter. Relativistic effects such as light bending and 324 gravitational delay are not accounted for in the light time 325 correction performed by this routine. 326 327 Stellar aberration 328 ------------------ 329 330 The velocity of the observer also affects the apparent location 331 of a target: photons arriving at the observer are subject to a 332 "raindrop effect" whereby their velocity relative to the observer 333 is, using a Newtonian approximation, the photons' velocity 334 relative to the solar system barycenter minus the velocity of the 335 observer relative to the solar system barycenter. This effect is 336 called "stellar aberration." Stellar aberration is independent 337 of the velocity of the target. The stellar aberration formula 338 used by this routine does not include (the much smaller) 339 relativistic effects. 340 341 Stellar aberration corrections are applied after light time 342 corrections: the light time corrected target position vector is 343 used as an input to the stellar aberration correction. 344 345 When light time and stellar aberration corrections are both 346 applied to a geometric position vector, the resulting position 347 vector indicates where the target "appears to be" from the 348 observer's location. 349 350 As opposed to computing the apparent position of a target, one 351 may wish to compute the pointing direction required for 352 transmission of photons to the target. This also requires correction 353 of the geometric target position for the effects of light time 354 and stellar aberration, but in this case the corrections are 355 computed for radiation traveling *from* the observer to the target. 356 We will refer to this situation as the "transmission" case. 357 358 The "transmission" light time correction yields the target's 359 location as it will be when photons emitted from the observer's 360 location at `et' arrive at the target. The transmission stellar 361 aberration correction is the inverse of the traditional stellar 362 aberration correction: it indicates the direction in which 363 radiation should be emitted so that, using a Newtonian 364 approximation, the sum of the velocity of the radiation relative 365 to the observer and of the observer's velocity, relative to the 366 solar system barycenter, yields a velocity vector that points in 367 the direction of the light time corrected position of the target. 368 369 One may object to using the term "observer" in the transmission 370 case, in which radiation is emitted from the observer's location. 371 The terminology was retained for consistency with earlier 372 documentation. 373 374 Below, we indicate the aberration corrections to use for some 375 common applications: 376 377 1) Find the apparent direction of a target. This is 378 the most common case for a remote-sensing observation. 379 380 Use "LT+S" or "CN+S": apply both light time and stellar 381 aberration corrections. 382 383 Note that using light time corrections alone ("LT") is 384 generally not a good way to obtain an approximation to an 385 apparent target vector: since light time and stellar 386 aberration corrections often partially cancel each other, 387 it may be more accurate to use no correction at all than to 388 use light time alone. 389 390 391 2) Find the corrected pointing direction to radiate a signal 392 to a target. This computation is often applicable for 393 implementing communications sessions. 394 395 Use "XLT+S" or "XCN+S": apply both light time and stellar 396 aberration corrections for transmission. 397 398 399 3) Compute the apparent position of a target body relative 400 to a star or other distant object. 401 402 Use one of "LT", "CN", "LT+S", or "CN+S" as needed to match 403 the correction applied to the position of the distant 404 object. For example, if a star position is obtained from a 405 catalog, the position vector may not be corrected for 406 stellar aberration. In this case, to find the angular 407 separation of the star and the limb of a planet, the vector 408 from the observer to the planet should be corrected for 409 light time but not stellar aberration. 410 411 412 4) Obtain an uncorrected position vector derived directly from 413 data in an SPK file. 414 415 Use "NONE". 416 417 418 5) Use a geometric position vector as a low-accuracy estimate 419 of the apparent position for an application where execution 420 speed is critical: 421 422 Use "NONE". 423 424 425 6) While this routine cannot perform the relativistic 426 aberration corrections required to compute positions 427 with the highest possible accuracy, it can supply the 428 geometric positions required as inputs to these computations: 429 430 Use "NONE", then apply relativistic aberration 431 corrections (not available in the SPICE Toolkit). 432 433 434 Below, we discuss in more detail how the aberration corrections 435 applied by this routine are computed. 436 437 Geometric case 438 ============== 439 440 spkpos_c begins by computing the geometric position T(et) of the 441 target body relative to the solar system barycenter (SSB). 442 Subtracting the geometric position of the observer O(et) gives 443 the geometric position of the target body relative to the 444 observer. The one-way light time, 'lt', is given by 445 446 | T(et) - O(et) | 447 lt = ------------------- 448 c 449 450 The geometric relationship between the observer, target, and 451 solar system barycenter is as shown: 452 453 454 SSB ---> O(et) 455 | / 456 | / 457 | / 458 | / T(et) - O(et) 459 V V 460 T(et) 461 462 463 The returned position is 464 465 T(et) - O(et) 466 467 468 Reception case 469 ============== 470 471 When any of the options "LT", "CN", "LT+S", "CN+S" is selected 472 for `abcorr', spkpos_c computes the position of the target body at 473 epoch et-lt, where 'lt' is the one-way light time. Let T(t) and 474 O(t) represent the positions of the target and observer 475 relative to the solar system barycenter at time t; then 'lt' is 476 the solution of the light-time equation 477 478 | T(et-lt) - O(et) | 479 lt = ------------------------ (1) 480 c 481 482 The ratio 483 484 | T(et) - O(et) | 485 --------------------- (2) 486 c 487 488 is used as a first approximation to 'lt'; inserting (2) into the 489 right hand side of the light-time equation (1) yields the 490 "one-iteration" estimate of the one-way light time ("LT"). 491 Repeating the process until the estimates of 'lt' converge yields 492 the "converged Newtonian" light time estimate ("CN"). 493 494 Subtracting the geometric position of the observer O(et) gives 495 the position of the target body relative to the observer: 496 T(et-lt) - O(et). 497 498 SSB ---> O(et) 499 | \ | 500 | \ | 501 | \ | T(et-lt) - O(et) 502 | \ | 503 V V V 504 T(et) T(et-lt) 505 506 The light time corrected position vector is 507 508 T(et-lt) - O(et) 509 510 If correction for stellar aberration is requested, the target 511 position is rotated toward the solar system 512 barycenter-relative velocity vector of the observer. The 513 rotation is computed as follows: 514 515 Let r be the light time corrected vector from the observer 516 to the object, and v be the velocity of the observer with 517 respect to the solar system barycenter. Let w be the angle 518 between them. The aberration angle phi is given by 519 520 sin(phi) = v sin(w) / c 521 522 Let h be the vector given by the cross product 523 524 h = r X v 525 526 Rotate r by phi radians about h to obtain the apparent 527 position of the object. 528 529 530 Transmission case 531 ================== 532 533 When any of the options "XLT", "XCN", "XLT+S", "XCN+S" is 534 selected, spkpos_c computes the position of the target body T at 535 epoch et+lt, where 'lt' is the one-way light time. 'lt' is the 536 solution of the light-time equation 537 538 | T(et+lt) - O(et) | 539 lt = ------------------------ (3) 540 c 541 542 Subtracting the geometric position of the observer, O(et), 543 gives the position of the target body relative to the 544 observer: T(et-lt) - O(et). 545 546 SSB --> O(et) 547 / | * 548 / | * T(et+lt) - O(et) 549 / |* 550 / *| 551 V V V 552 T(et+lt) T(et) 553 554 The position component of the light-time corrected position 555 is the vector 556 557 T(et+lt) - O(et) 558 559 If correction for stellar aberration is requested, the target 560 position is rotated away from the solar system barycenter- 561 relative velocity vector of the observer. The rotation is 562 computed as in the reception case, but the sign of the 563 rotation angle is negated. 564 565 Precision of light time corrections 566 =================================== 567 568 Corrections using one iteration of the light time solution 569 ---------------------------------------------------------- 570 571 When the requested aberration correction is "LT", "LT+S", 572 "XLT", or "XLT+S", only one iteration is performed in the 573 algorithm used to compute 'lt'. 574 575 The relative error in this computation 576 577 | LT_ACTUAL - LT_COMPUTED | / LT_ACTUAL 578 579 is at most 580 581 (V/C)**2 582 ---------- 583 1 - (V/C) 584 585 which is well approximated by (V/C)**2, where V is the 586 velocity of the target relative to an inertial frame and C is 587 the speed of light. 588 589 For nearly all objects in the solar system V is less than 60 590 km/sec. The value of C is ~300000 km/sec. Thus the 591 one-iteration solution for LT has a potential relative error 592 of not more than 4e-8. This is a potential light time error of 593 approximately 2e-5 seconds per astronomical unit of distance 594 separating the observer and target. Given the bound on V cited 595 above: 596 597 As long as the observer and target are separated by less 598 than 50 astronomical units, the error in the light time 599 returned using the one-iteration light time corrections is 600 less than 1 millisecond. 601 602 The magnitude of the corresponding position error, given 603 the above assumptions, may be as large as (V/C)**2 * the 604 distance between the observer and the uncorrected target 605 position: 300 km or equivalently 6 km/AU. 606 607 In practice, the difference between positions obtained using 608 one-iteration and converged light time is usually much smaller 609 than the value computed above and can be insignificant. For 610 example, for the spacecraft Mars Reconnaissance Orbiter and 611 Mars Express, the position error for the one-iteration light 612 time correction, applied to the spacecraft-to-Mars center 613 vector, is at the 1 cm level. 614 615 Comparison of results obtained using the one-iteration and 616 converged light time solutions is recommended when adequacy of 617 the one-iteration solution is in doubt. 618 619 620 Converged corrections 621 --------------------- 622 623 When the requested aberration correction is "CN", "CN+S", 624 "XCN", or "XCN+S", as many iterations as are required for 625 convergence are performed in the computation of LT. Usually 626 the solution is found after three iterations. The relative 627 error present in this case is at most 628 629 (V/C)**4 630 ---------- 631 1 - (V/C) 632 633 which is well approximated by (V/C)**4. 634 635 The precision of this computation (ignoring round-off 636 error) is better than 4e-11 seconds for any pair of objects 637 less than 50 AU apart, and having speed relative to the 638 solar system barycenter less than 60 km/s. 639 640 The magnitude of the corresponding position error, given 641 the above assumptions, may be as large as (V/C)**4 * the 642 distance between the observer and the uncorrected target 643 position: 1.2 cm at 50 AU or equivalently 0.24 mm/AU. 644 645 However, to very accurately model the light time between 646 target and observer one must take into account effects due to 647 general relativity. These may be as high as a few hundredths 648 of a millisecond for some objects. 649 650 651 Relativistic Corrections 652 ========================= 653 654 This routine does not attempt to perform either general or 655 special relativistic corrections in computing the various 656 aberration corrections. For many applications relativistic 657 corrections are not worth the expense of added computation 658 cycles. If however, your application requires these additional 659 corrections we suggest you consult the astronomical almanac (page 660 B36) for a discussion of how to carry out these corrections. 661 662 663 -Examples 664 665 1) Load a planetary ephemeris SPK, then look up a series of 666 geometric positions of the moon relative to the earth, 667 referenced to the J2000 frame. 668 669 #include <stdio.h> 670 #include "SpiceUsr.h" 671 672 void main() 673 { 674 675 #define ABCORR "NONE" 676 #define FRAME "J2000" 677 678 /. 679 The name of the SPK file shown here is fictitious; 680 you must supply the name of an SPK file available 681 on your own computer system. 682 ./ 683 #define SPK "planetary_spk.bsp" 684 685 /. 686 ET0 represents the date 2000 Jan 1 12:00:00 TDB. 687 ./ 688 #define ET0 0.0 689 690 /. 691 Use a time step of 1 hour; look up 100 states. 692 ./ 693 #define STEP 3600.0 694 #define MAXITR 100 695 696 #define OBSERVER "earth" 697 #define TARGET "moon" 698 699 700 /. 701 Local variables 702 ./ 703 SpiceInt i; 704 705 SpiceDouble et; 706 SpiceDouble lt; 707 SpiceDouble pos [3]; 708 709 710 /. 711 Load the spk file. 712 ./ 713 furnsh_c ( SPK ); 714 715 /. 716 Step through a series of epochs, looking up a position vector 717 at each one. 718 ./ 719 for ( i = 0; i < MAXITR; i++ ) 720 { 721 et = ET0 + i*STEP; 722 723 spkpos_c ( TARGET, et, FRAME, ABCORR, 724 OBSERVER, pos, < ); 725 726 printf( "\net = %20.10f\n\n", et ); 727 printf( "J2000 x-position (km): %20.10f\n", pos[0] ); 728 printf( "J2000 y-position (km): %20.10f\n", pos[1] ); 729 printf( "J2000 z-position (km): %20.10f\n", pos[2] ); 730 } 731 } 732 733 734 -Restrictions 735 736 None. 737 738 -Literature_References 739 740 SPK Required Reading. 741 742 -Author_and_Institution 743 744 C.H. Acton (JPL) 745 B.V. Semenov (JPL) 746 N.J. Bachman (JPL) 747 W.L. Taber (JPL) 748 749 -Version 750 751 -CSPICE Version 3.0.1, 07-JUL-2014 (NJB) 752 753 Discussion of light time corrections was updated. Assertions 754 that converged light time corrections are unlikely to be 755 useful were removed. 756 757 -CSPICE Version 2.0.4, 04-APR-2008 (NJB) 758 759 Corrected minor error in description of XLT+S aberration 760 correction. 761 762 -CSPICE Version 2.0.3, 17-APR-2005 (NJB) 763 764 Error was corrected in example program: variable name `state' 765 was changed to `pos' in printf calls. 766 767 -CSPICE Version 2.0.2, 13-OCT-2003 (EDW) 768 769 Various minor header changes were made to improve clarity. 770 Added mention that 'lt' returns a value in seconds. 771 772 -CSPICE Version 2.0.1, 27-JUL-2003 (NJB) (CHA) 773 774 Various header corrections were made. 775 776 -CSPICE Version 2.0.0, 31-DEC-2001 (NJB) 777 778 Updated to handle aberration corrections for transmission 779 of radiation. Formerly, only the reception case was 780 supported. The header was revised and expanded to explain 781 the functionality of this routine in more detail. 782 783 -CSPICE Version 1.0.0, 29-MAY-1999 (NJB) (WLT) 784 785 -Index_Entries 786 787 using names get target position relative to an observer 788 position relative to observer corrected for aberrations 789 read ephemeris data 790 read trajectory data 791 792 -& 793 */ 794 795 { /* Begin spkpos_c */ 796 797 /* 798 Participate in error tracing. 799 */ 800 801 chkin_c ( "spkpos_c" ); 802 803 804 /* 805 Check the input strings to make sure the pointers are non-null 806 and the string lengths are non-zero. 807 */ 808 CHKFSTR ( CHK_STANDARD, "spkpos_c", targ ); 809 CHKFSTR ( CHK_STANDARD, "spkpos_c", ref ); 810 CHKFSTR ( CHK_STANDARD, "spkpos_c", abcorr ); 811 CHKFSTR ( CHK_STANDARD, "spkpos_c", obs ); 812 813 814 /* 815 Call the f2c'd Fortran routine. Use explicit type casts for every 816 type defined by f2c. 817 */ 818 spkpos_ ( ( char * ) targ, 819 ( doublereal * ) &et, 820 ( char * ) ref, 821 ( char * ) abcorr, 822 ( char * ) obs, 823 ( doublereal * ) ptarg, 824 ( doublereal * ) lt, 825 ( ftnlen ) strlen(targ), 826 ( ftnlen ) strlen(ref), 827 ( ftnlen ) strlen(abcorr), 828 ( ftnlen ) strlen(obs) ); 829 830 831 chkout_c ( "spkpos_c" ); 832 833 } /* End spkpos_c */ 834