Return to nutate.c CVS log

Up to [Research Unix] / researchv10no / cmd / sky

Annotation of researchv10no/cmd/sky/nutate.c, revision 1.1

1.1     ! root        1: #include "sky.h"
        !             2: 
        !             3: extern struct nuttab
        !             4: {
        !             5:        float f[2];
        !             6:        char c[5];
        !             7: } nuttab[];
        !             8: 
        !             9: nutate()
        !            10: {
        !            11: 
        !            12: /*
        !            13:  *     uses radian, radsec, eday
        !            14:  *     sets psi, eps, dpsi, deps, obliq, gst, tobliq
        !            15:  */
        !            16: 
        !            17:        double msun, mnom, noded, dmoon;
        !            18:        struct nuttab *pp;
        !            19:        double arg;
        !            20:        double Eday;
        !            21:        double capT, capT2, capT3;
        !            22: 
        !            23: /*
        !            24:  *     nutation of the equinoxes is a irregular wobble of the
        !            25:  *     pole of the earths rotation whose magnitude is about
        !            26:  *     9 seconds of arc and whose period is about 18.6 years.
        !            27:  *
        !            28:  *     it depends upon the pull of the sun and moon on the
        !            29:  *     equatorial bulge of the earth.
        !            30:  *
        !            31:  *     psi and eps are the two angles which specify the
        !            32:  *     true pole with respect to the mean pole.
        !            33:  *
        !            34:  *     all coeffieients are from Supplement to A.E. 1984
        !            35:  */
        !            36: 
        !            37:        pp = &nuttab[0];
        !            38:        Eday = eday - 36525;
        !            39:        capT = capt-1;
        !            40:        capT2 = capT*capT;
        !            41:        capT3 = capT2*capT;
        !            42: 
        !            43:        mnom = 134.96298139 + 13.06499294724*Eday + 8.6972e-3*capT2
        !            44:                 + 17.8e-6*capT3;
        !            45:        msun = 357.52772333 + 0.98560028309*Eday - 0.16028e-3*capT2
        !            46:                 - 3.33e-6*capT3;
        !            47:        noded = 93.27191028 + 13.22935024060*Eday - 3.6825e-3*capT2
        !            48:                 + 3.05e-6*capT3;
        !            49:        dmoon = 297.85036306 + 12.19074911650*Eday - 1.91417e-3*capT2
        !            50:                 + 5.28e-6*capT3;
        !            51:        node = 125.04452222 - .05295376484*Eday + 2.07083e-3*capT2
        !            52:                 + 2.22e-6*capT3;
        !            53:        mnom *= radian;
        !            54:        msun *= radian;
        !            55:        noded *= radian;
        !            56:        dmoon *= radian;
        !            57:        node *= radian;
        !            58: 
        !            59: /*
        !            60:  *     long period terms
        !            61:  */
        !            62: 
        !            63: 
        !            64:        psi = -(17.1996+.01742*capT)*sin(node);
        !            65:        eps = (9.2025+0.00089*capT)*cos(node);
        !            66:        for(;;){
        !            67:                if(pp->f[0]==0.){
        !            68:                        pp++;
        !            69:                        break;
        !            70:                }
        !            71:                arg = pp->c[0]*mnom + pp->c[1]*msun + pp->c[2]*noded
        !            72:                        + pp->c[3]*dmoon + pp->c[4]*node;
        !            73:                psi += pp->f[0]*sin(arg);
        !            74:                eps += pp->f[1]*cos(arg);
        !            75:                pp++;
        !            76:        }
        !            77: 
        !            78: /*
        !            79:  *     short period terms
        !            80:  */
        !            81: 
        !            82:        dpsi = 0.;
        !            83:        deps = 0.;
        !            84:        if((flags&APPARENT)==0){
        !            85:                for(;;){
        !            86:                        if(pp->f[0]==0.){
        !            87:                                pp++;
        !            88:                                break;
        !            89:                        }
        !            90:                        arg = pp->c[0]*mnom + pp->c[1]*msun + pp->c[2]*noded
        !            91:                                + pp->c[3]*dmoon + pp->c[4]*node;
        !            92:                        dpsi += pp->f[0]*sin(arg);
        !            93:                        deps += pp->f[1]*cos(arg);
        !            94:                        pp++;
        !            95:                }
        !            96:        }
        !            97: 
        !            98: 
        !            99:        psi = (psi+dpsi)*radsec;
        !           100:        eps = (eps+deps)*radsec;
        !           101:        dpsi *= radsec;
        !           102:        deps *= radsec;
        !           103: 
        !           104:        obliq = 23.43929111 - .013004167*capT - 0.164e-6*capT2
        !           105:                 + 0.5036e-6*capT3;
        !           106:        obliq *= radian;
        !           107:        tobliq = obliq + eps;
        !           108: 
        !           109:        gst = 100.460618375 + 360.98564736629*Eday + 0.3882667e-3*capT2
        !           110:                 - 0.025e-6*capT3;
        !           111:        gst -= 180.;
        !           112:        gst = fmod(gst, 360.);
        !           113:        if(gst < 0.)
        !           114:                gst += 360.;
        !           115:        gst *= radian;
        !           116:        gst += psi*cos(tobliq);
        !           117: 
        !           118: /*
        !           119:  * print true obliquity, mean obliquity, mean GST, and true GST.
        !           120:  */
        !           121: 
        !           122:        if(flags&GOOBIE){
        !           123:                printf("Mean G. Sid. Time: ");
        !           124:                prhms(gst - psi*cos(obliq), 4);
        !           125:                printf("\n");
        !           126: 
        !           127:                printf("App. G. Sid. Time: ");
        !           128:                prhms(gst,4);
        !           129:                printf("\n");
        !           130: 
        !           131:                printf("Goobie: %6.3f %6.3f %6.3f %6.3f\n", psi/radsec,
        !           132:                        eps/radsec, dpsi/radsec, deps/radsec);
        !           133: 
        !           134:                printf("Obliquity: ");
        !           135:                prdms(tobliq,3);
        !           136:                printf("\n");
        !           137:        }
        !           138: 
        !           139: }

This archive runs on limited infrastructure. Preserving old code on modern bandwidth. Automated agents are requested to crawl responsibly.