| 8043 | #endif |
| 8044 | |
| 8045 | int32 CALL_CONV swe_calc_pctr(double tjd, int32 ipl, int32 iplctr, int32 iflag, double *xxret, char *serr) |
| 8046 | { |
| 8047 | double t = 0, dt, daya[2], dtsave_for_defl = 0; |
| 8048 | double xx[6], xxctr[6], xxctr2[6], xx0[6], xxsv[24], xxsp[6], dx[6], xreturn[24]; |
| 8049 | double *xs; |
| 8050 | int i, j, niter; |
| 8051 | int32 iflag2, epheflag, retc; |
| 8052 | struct epsilon *oe; |
| 8053 | if (ipl == iplctr) { |
| 8054 | if (serr != NULL) |
| 8055 | sprintf(serr, "ipl and iplctr (= %d) must not be identical\n", ipl); |
| 8056 | return ERR; |
| 8057 | } |
| 8058 | iflag = plaus_iflag(iflag, ipl, tjd, serr); |
| 8059 | epheflag = iflag & SEFLG_EPHMASK; |
| 8060 | // this fills in obliquity and nutation values in swed |
| 8061 | swe_calc(tjd + swe_deltat_ex(tjd, epheflag, serr), SE_ECL_NUT, iflag, xx, serr); |
| 8062 | iflag &= ~(SEFLG_HELCTR|SEFLG_BARYCTR); |
| 8063 | iflag2 = epheflag; |
| 8064 | iflag2 |= (SEFLG_BARYCTR|SEFLG_J2000|SEFLG_ICRS|SEFLG_TRUEPOS|SEFLG_EQUATORIAL|SEFLG_XYZ|SEFLG_SPEED); |
| 8065 | iflag2 |= (SEFLG_NOABERR|SEFLG_NOGDEFL); |
| 8066 | retc = swe_calc(tjd, iplctr, iflag2, xxctr, serr); |
| 8067 | if (retc == ERR) |
| 8068 | return ERR; |
| 8069 | retc = swe_calc(tjd, ipl, iflag2, xx, serr); |
| 8070 | if (retc == ERR) |
| 8071 | return ERR; |
| 8072 | for (i = 0; i <= 5; i++) { |
| 8073 | xx0[i] = xx[i]; |
| 8074 | //xx[i] -= xxctr[i]; |
| 8075 | } |
| 8076 | /******************************* |
| 8077 | * light-time geocentric * |
| 8078 | *******************************/ |
| 8079 | if (!(iflag & SEFLG_TRUEPOS)) { |
| 8080 | /* number of iterations - 1 */ |
| 8081 | niter = 1; |
| 8082 | if (iflag & SEFLG_SPEED) { |
| 8083 | /* |
| 8084 | * Apparent speed is influenced by the fact that dt changes with |
| 8085 | * time. This makes a difference of several hundredths of an |
| 8086 | * arc second / day. To take this into account, we compute |
| 8087 | * 1. true position - apparent position at time t - 1. |
| 8088 | * 2. true position - apparent position at time t. |
| 8089 | * 3. the difference between the two is the part of the daily motion |
| 8090 | * that results from the change of dt. |
| 8091 | */ |
| 8092 | for (i = 0; i <= 2; i++) |
| 8093 | xxsv[i] = xxsp[i] = xx[i] - xx[i+3]; |
| 8094 | for (j = 0; j <= niter; j++) { |
| 8095 | for (i = 0; i <= 2; i++) { |
| 8096 | dx[i] = xxsp[i]; |
| 8097 | dx[i] -= (xxctr[i] - xxctr[i+3]); |
| 8098 | } |
| 8099 | /* new dt */ |
| 8100 | dt = sqrt(square_sum(dx)) * AUNIT / CLIGHT / 86400.0; |
| 8101 | for (i = 0; i <= 2; i++) /* rough apparent position at t-1 */ |
| 8102 | xxsp[i] = xxsv[i] - dt * xx0[i+3]; |
no test coverage detected