Astronomy...
06-20-2015, 05:53 PM (This post was last modified: 06-21-2015 06:29 PM by salvomic.)
Post: #1
 salvomic Senior Member Posts: 1,396 Joined: Jan 2015
Astronomy...
hi everybody,
waiting for the evolution of the very interesting app by Marcel (I'm already using!), I'm developing also this program to get astronomical results for Sun, Moon and Planets...
Effemeridi (Astrolabio) is a program that I'm adapted from an old version I made many years ago for a TI pocket calculator.
The program has still some incorrect output (very few), so I'm glad to have your advice and hints to improve it, thank you.
The program calculates Sidereal Time, Longitude, Latitude (various coordinates), Rect Ascension and Declination, Azimuth and Height, phases of the Moon, Retrogradation of planets, Ascendent, after a book of S. Bouiges...
I'm working to add eclipses, and improving program also in "astrological" sense (I was asked for this) and Satellites of Jupiter and Saturn.
Also data for Pluto are still not so precise as those in the app of Marcel...

Some routines and data are still in Italian, sorry.

Salvo

The code (I'll post here the new versions for now):
Code:
 data(); calcN(); calcTS(); sunCalc(); planetCalc(); planets(); sun(); moon(); Mercurius(); Venus(); Mars(); Jupiter(); Saturn(); Uranus(); Neptunus(); Pluto(); ascendent(); zodiaco(); trunc_angle(); atan2(); thisday:= 2000.0101; lstd:=0; gmt:= 12; long:= 0; lat:=  0; m:= 1; dst:= 0; utc:=0; sunlist:= MAKELIST(X,X,1,25); planet_const:= MAKELIST(X,X,1,10); ε:= →HMS(23.4392911);  // approx 23.45 EXPORT Effemeridi() // from "Astrolabio ed Effemeridi", by Salvo Micciché, salvomic, 2015 // Credits: from a book, Calcule astronomique pour amateurs, by Serge Bouiges // Thanks Didier Lachieze, Marcel, Eried BEGIN   LOCAL ch;   HAngle:=1;   data(); // input data and calc Sun   CHOOSE(ch, "Effemeridi", "Sun", "Moon", "Planets", "Calc N", "Calc TS", "Ascendent", "Quit");   CASE   IF ch==1 THEN sun(); END;   IF ch==2 THEN moon(); END;   IF ch==3 THEN planets(); END;   IF ch==4 THEN calcN(); END;   IF ch==5 THEN calcTS(); END;   IF ch==6 THEN ascendent(); END;   IF ch==7 THEN HAngle:=0; RETURN; END;   DEFAULT   END; // case END; data() BEGIN   LOCAL dd, mm, yy, hh, mn, ss, loct;   yy:= IP(Date);   mm:= IP(FP(Date)*100);   dd:= FP(Date*100)*100;   hh:= IP(HMS→(Time));   mn:= IP(FP(HMS→(Time))*60);   ss:= IP(FP(FP(FP(HMS→(Time))*60))*60);   loct:= →HMS(hh+mn/60+ss/3600);   INPUT({ {m,[0],{15,15,1}}, {D,[0],{50,15,1}}, {Y,[0],{80,15,1}},   {lstd,[0],{65,30,2}},   {dst,0,{40,15,3}}, {utc,[0],{80,15,3}},    {long,[0],{20,25,4}},    {lat,[0],{70,25,4}} },   "Data: Use Shift+9 for H°M′S″",   {"Month:","Day :","Year :","Local Time (24):", "DST", "GMT",    "Long (-E):","Lat (+N)"},   {"Enter month", "Enter day", "Enter year (1901-2099)",    "Enter local time (decimal or H°M′S″)", "Daylight Save Time?", "UTC lag",   "Longitude", "Latitude"},    {1,1,1901,0, 0,0, 0, 0},{mm,dd,yy,loct, 0, 1, -14°43′41″, 36°55′15″});   // Ragusa (-14°43′41″, 36°55′15″) -    // Marina di Ragusa (-14°33′15″, 36°47′06″) - Scicli (-14°42'22″, 36°47' 22″)   D:=EVAL(D); m:=EVAL(m); Y:=EVAL(Y);   lstd:= lstd;    // Greenwich Mean Time   gmt:= lstd-utc-(1*dst); // UTC (DST, daylight save time)   thisday:= Y+m/100+D/10000;   sunCalc(); // calculate Sun data END; planets() BEGIN   LOCAL ch, nameplanet;   INPUT({{ch,{"Mercurius","Venus","Mars", "Jupiter", "Saturn", "Uranus", "Neptunus", "Pluto"}, {20,30,1}}},   "Choose planet",    "planet: ", "Choose the planet to calculate an press OK", 0,5 );   CASE   IF ch==1 THEN nameplanet:="Mercurius"; planetCalc(nameplanet); END;   IF ch==2 THEN nameplanet:="Venus"; planetCalc(nameplanet); END;   IF ch==3 THEN nameplanet:="Mars"; planetCalc(nameplanet); END;   IF ch==4 THEN nameplanet:="Jupiter"; planetCalc(nameplanet); END;   IF ch==5 THEN nameplanet:="Saturn"; planetCalc(nameplanet); END;   IF ch==6 THEN nameplanet:="Uranus"; planetCalc(nameplanet); END;   IF ch==7 THEN nameplanet:="Neptunus"; planetCalc(nameplanet); END;   IF ch==8 THEN nameplanet:="Pluto"; planetCalc(nameplanet); END;   DEFAULT   END; // case END; calcN() BEGIN   N:= DDAYS(1901.0101, thisday)+1+gmt/24; // N to GMT   RETURN N; END; calcTS() BEGIN   LOCAL TS, TSd, N0, TSL;   calcN();   N0:= DDAYS(1901.0101, thisday)+1+0; // N at 0 GMT (for TS)   TS:= ((23750.3 + 236.555362*(N0+gmt/24)) / 3600) MOD 24; // Sideral Time (GMT) in seconds   TS:= →HMS(TS); // in HMS   TSL:= (((23750.3 + 236.555362*N+lstd/24) / 3600) MOD 24); // local ST   TSL:= (TSL + gmt - 4*(HMS→(long))/60) MOD 24;   TSL:= →HMS(TSL); // in HMS   TSd:= →HMS((98.965 + 0.985647342*N0) MOD 360); // ST at 0 GMT in degrees   RETURN {TS, TSL, TSd, gmt}; END; ascendent() BEGIN   LOCAL ASC, ts, segno;   ts:= calcTS();   ASC:= atan2(-cos(ts(2)), sin(ts(2))*cos(ε)+tan(lat)*sin(ε));   IF ASC<0 THEN ASC:= ASC+360; END; IF ASC>360 THEN ASC:=ASC-360; END;   segno:= zodiaco(ASC); RETURN {→HMS(ASC), segno}; END; zodiaco(astro) BEGIN   LOCAL segno;   CASE  // Segni e case   IF astro>=0   AND astro <30  THEN segno:="Ariete"; END;   IF astro>=30  AND astro <60  THEN segno:="Toro"; END;   IF astro>=60  AND astro <90  THEN segno:="Gemelli"; END;   IF astro>=90  AND astro <120 THEN segno:="Cancro"; END;   IF astro>=120 AND astro <150 THEN segno:="Leone"; END;   IF astro>=150 AND astro <180 THEN segno:="Vergine"; END;   IF astro>=180 AND astro <210 THEN segno:="Bilancia"; END;   IF astro>=210 AND astro <240 THEN segno:="Scorpione"; END;   IF astro>=240 AND astro <270 THEN segno:="Sagittario"; END;   IF astro>=270 AND astro <300 THEN segno:="Capricorno"; END;   IF astro>=300 AND astro <330 THEN segno:="Acquario"; END;   IF astro>=330 AND astro <360 THEN segno:="Pesci"; END;   END; // case   RETURN segno; END; atan2(y,x) BEGIN // atan2 returns coordinates in correct angle for all four quadrants (DEG)     CASE     IF x==0 AND y==0 then return "undefined"; end;         // x=0, y=0       IF x>0 then return atan(y/x); end;                     // x>0     IF x<0 AND y>=0 then return atan(y/x)+180; end;         // x<0, y>=0     IF x==0 AND y>0 then return 90; end;                 // x=0, y>0     IF x<0 AND y<=0 then return atan(y/x)-180; end;         // x<0, y<0     IF x==0 AND y<0 then return -90; end;                // x=0, y<0    END;   RETURN;  END; trunc_angle(ang) // trunc decimal from a DEG form angle BEGIN →HMS((IP(HMS→(ang)*3600)/3600)); END; sunCalc() BEGIN   LOCAL ω, l, v, r, x, y, dz;   LOCAL α, δ, dayY, H, ts, az, z, h;   LOCAL SAD, DA, sorg, tram, culm, Hd, lt;   LOCAL as, ac, zc;   calcN();   L:= →HMS((278.965+0.985647342*N)) MOD 360; // long media   ω:= →HMS((281.235 + 0.0000469*N)) MOD 360; // long del perielio   M:= (L - ω); // anomalia vera   v:= M + 1.91934 * sin(M)+0.02*sin(2*M); // calcolo intermedio   l:= v + ω; // Longitudine Sole   lt:= l + 180; // Longitudine eliocentrica della Terra   r:= ((0.999721)/(1+0.01675*cos(v))); // Distanza Sole-Terra (UA)   x:= r * cos(l); y:= r * sin(l);   dayY:= DDAYS(IP(thisday)+1/100+1/10000, thisday); // calc day number   // δ:= →HMS((ε*sin((360*(284+dayY))/365)) MOD 360); // declinazione, formula di Cooper   δ:= →HMS((asin(sin(ε)*sin(l))) );   // α:= →HMS((acos(cos(l)/cos(δ))/15) MOD 24); // ascensione retta   α:= (acos(cos(IFTE(l>180, l-180, l))/cos(δ))); // ascensione retta   IF l>180 THEN α:= α +180; END;   α:= →HMS(α/15) MOD 24; // α in hms   ts:= calcTS();  // calcola tempo siderale   H:= (ts(2) - α) MOD 24; // Angolo orario   Hd:= →HMS(15*H); // Angolo orario in gradi   zc:= sin(lat)*sin(δ)+cos(lat)*cos(δ)*cos(Hd); // calc z, alt   dz:= →HMS(acos(zc)) ; // distanza zenitale   // alt:= →HMS(90 - z); // Altezza   h:= →HMS(ATAN(zc/sqrt(1-zc^2))); // Altezza   // as:= cos(δ)*sin(Hd)/sin(dz); // calc az   // ac:= (-cos(lat)*sin(δ)+sin(lat)*cos(δ)*cos(Hd))/sin(dz);   // az:= →HMS(ATAN(as/ac) ); // Azimuth   // IF (ac<0) THEN az:= (360 + az) MOD 360; END;   az:= atan2((-cos(δ) * sin(Hd)), cos(lat) *sin(δ)-sin(lat)*cos(δ)*cos(Hd));   IF az<0 THEN az:= az+360; END; IF az>360 THEN az:=az-360; END;   IF (az<0) THEN az:= 360+az; END;   E:= →HMS((460*sin(M)-592*sin(2*L))/60); // equazione del tempo in min (460s, 592s)   DA:= asin(tan(lat)*tan(δ)); // Differenza ascensionale   SAD:= (90 + DA); // Semiarco Diurno   // SAD:= →HMS(acos(-tan(lat)*tan(δ))); // Semiarco Diurno   culm:= (IFTE(α<ts(1),24+α,α) - ts(1)) ; // Culminazione   sorg:= (culm - SAD/15) MOD 24; // Alba (Dawn)   tram:= (culm + SAD/15) MOD 24; // Tramonto (Sunset)   sunlist:= {L,ω,l, M, α,δ, az, h, x,y,sorg, tram, culm, r, H, Hd, zc, dz, ts(1), DA, SAD, E, lt};   RETURN sunlist; END; sun() BEGIN   LOCAL ω, l, v, r, x, y, dz;   LOCAL α, δ, dayY, H, ts, az, z, h;   LOCAL SAD, DA, sorg, tram, culm, Hd, lt;   LOCAL as, ac, zc, segno;   L:= sunlist(1); ω:= sunlist(2);    l:= sunlist(3); M:= sunlist(4);    α:= sunlist(5); δ:= sunlist(6);    az:= sunlist(7); h:= sunlist(8);    x:= sunlist(9); y:= sunlist(10);   sorg:= sunlist(11); tram:= sunlist(12); culm:= sunlist(13);   r:= sunlist(14); H:= sunlist(15); Hd:= sunlist(16);   zc:= sunlist(17); dz:= sunlist(18);    ts:= calcTS(); // no 19   DA:= sunlist(20); SAD:= sunlist(21); E:= sunlist(22); lt:= sunlist(23);   segno:= zodiaco(l);   PRINT;   PRINT("Effemeridi del Sole a " + m+" "+D+" "+Y+" "+ lstd);   PRINT("");   PRINT("L longitudine media " + trunc_angle(L));   PRINT("ω long del perielio " + trunc_angle(ω));   PRINT("M anomalia vera " + trunc_angle(M));   PRINT("λ longitudine Sole " + trunc_angle(l) + " " + segno);   PRINT("r Distanza (UA) " + r);   PRINT("α Asc. r. "+ trunc_angle(α) + " δ Decl. "+ trunc_angle(δ));   PRINT("Az Azimuth " + trunc_angle(az) + " h Alt. " + trunc_angle(h));   PRINT(" ");   PRINT("Ang. orario "+ H + " " + trunc_angle(Hd));   PRINT("TS (UTC) " + trunc_angle(ts(1)));   PRINT("TS locale " + trunc_angle(ts(2)));   PRINT("Equazione del tempo (min) " + trunc_angle(E));   PRINT("Dif asc." + DA + " SAD " + trunc_angle(SAD));   PRINT("");   PRINT("Culmina a: " + trunc_angle(culm) + " UTC");   PRINT("Sorge a: " + trunc_angle(sorg) + " UTC");    PRINT("Tramonta a: " + trunc_angle(tram) + " UTC");   PRINT("Long elioc. Terra "+ trunc_angle(lt));   // RETURN {L, ω, M, v, l, r, x, y};   RETURN [[L, ω], [l, M], [α,δ], [az, h], [x,y], [sorg, tram], [culm, r], [H,ts(1)]]; END; moon() BEGIN   LOCAL ω,Ω, l, v, r, x, y, i;   LOCAL Ls, Ms, LongM, LatM, P, d, s;   LOCAL α,δ, ts, Hd, zc, dz, h, az;   LOCAL as, ac, DA, SAD, DD, culm, sorg, tram;   LOCAL b, dayY, xs, ys, ph, fase, eta;   LOCAL dd, hh, mm, segno;   i:= 5.13333333; // inclinazione 5°8' sull'eclittica   calcN();   // Sun   Ls:= sunlist(1); // Longitudine media   Ms:= sunlist(4); // Anomalia   // end Sun   L:= (33.231 + 13.17639653*N) MOD 360; // long. media   Ω:= (239.882 - 0.052953922*N) MOD 360; // Long nodo ascendente   M:= (18.294 + 13.06499245*N) MOD 360; // anomalia vera   ω:= (L - M) MOD 360; // perigeo medio   DD:= L - Ls; F:= L-Ω; // for the variations   LongM:= L + 6.28875 * sin(M); // calculus (h/ and above) of Moon longitude   LongM:= LongM + 0.2136 * sin(2*M); // equazione del centro   LongM:= LongM + 0.6583 * sin(2*DD); // variazione   LongM:= LongM - 0.1856 * sin(Ms); // equazione annuale   LongM:= LongM + 1.2740 * sin(2*DD-M); // evezione   LongM:= LongM - 0.1143 * sin(2*F); // riduzione all'eclittica   LongM:= LongM + 0.0588 * sin(2*DD-2*M);   LongM:= LongM + 0.0572 * sin(2*DD-M-Ms);   LongM:= LongM + 0.0533 * sin(2*DD+M);   LongM:= LongM + 0.0459 * sin(2*DD-Ms);   LongM:= LongM + 0.0410 * sin(M-Ms);   LongM:= LongM - 0.0305 * sin(M+Ms);   LongM:= LongM - 0.0348 * sin(DD);   LongM:= LongM MOD 360;   LatM := 5.1280 * sin(F); // termine principale   LatM := LatM + 0.2806 * sin(M+F);   LatM := LatM + 0.2777 * sin(M-F); // equazione del centro   LatM := LatM + 0.1732 * sin(2*DD-F); // grande ineguaglianza   LatM := LatM + 0.0554 * sin(2*DD-M+F); // evezione   LatM := LatM + 0.0462 * sin(2*DD-M-F); // evezione (2)   LatM := LatM + 0.0326 * sin(2*DD+F);   LatM := LatM MOD 360; b:= LatM; dayY:= DDAYS(IP(thisday)+1/100+1/10000, thisday); // calc day number   δ:= →HMS((asin(cos(ε)*sin(b)+sin(ε)*cos(b)*sin(LongM))) ); // declinazione (coord locali)   // α:= →HMS((acos(cos(b)*cos(LongM)/cos(δ))/15) MOD 24); // ascensione retta   α:= (acos(cos(b)*cos(IFTE(LongM>180, LongM-180, LongM))/cos(δ))); // ascensione retta   IF LongM>180 THEN α:= α +180; END;   α:= →HMS(α/15) MOD 24; // α in hms   ts:= calcTS();  // calcola tempo siderale   H:= (ts(2) - α) MOD 24; // Angolo orario   Hd:= →HMS(15*H); // Angolo orario in gradi   zc:= sin(lat)*sin(δ)+cos(lat)*cos(δ)*cos(Hd); // calc z, alt   dz:= →HMS(acos(zc)) ; // distanza zenitale :: h:= →HMS(90 - dz); // Altezza   h:= →HMS(ATAN(zc/sqrt(1-zc^2))); // Altezza (height)   as:= cos(δ)*sin(Hd)/sin(dz); // calc az   ac:= (-cos(lat)*sin(δ)+sin(lat)*cos(δ)*cos(Hd))/sin(dz);   // az:= →HMS(ATAN(as/ac) ); // Azimuth   // IF (ac<0) THEN az:= ((360 + az) MOD 360); END;   az:= atan2((-cos(δ) * sin(Hd)), cos(lat) *sin(δ)-sin(lat)*cos(δ)*cos(Hd));   IF az<0 THEN az:= az+360; END; IF az>360 THEN az:=az-360; END;   DA:= asin(tan(lat)*tan(δ)); // Differenza ascensionale   SAD:= (90 + DA); // Semiarco Diurno   // SAD:= acos(-tan(lat)*tan(δ)); // Semiarco Diurno   culm:= (IFTE(α<ts(1),24+α,α) - ts(1)) ; // Culminazione   sorg:= (culm - SAD/15) MOD 24; // Alba (Dawn)   tram:= (culm + SAD/15) MOD 24; // Tramonto (Set)   P:= →HMS((0.9508+0.0518*cos(M)) MOD 360); // parallasse   // Raggio quadratico medio della Terra: 6373,044737   d:= 1/sin(P); // distanza Terra-Luna in raggi terrestri   s:= asin((3/11)*sin(P)); // semidiametro (raggi terrestri)   ph:= abs(LongM - Ls); // Phase   CASE   // IF ph == 0 THEN fase:="New Moon"; END;   IF ph>280 AND ph<350 THEN fase:="Decrescente"; END;   IF ph>=260 AND ph<=280 THEN fase:="Last Quart"; END;   IF ph>190 AND ph<260 THEN fase:="Gibbosa Cal."; END;   IF ph>=170 AND ph<=190 THEN fase:="Full Moon"; END;   IF ph>100 AND ph<170 THEN fase:="Gibbosa Cresc."; END;   IF ph>=80 AND ph<=100 THEN fase:="1st Quart"; END;   IF ph>10 AND ph<80 THEN fase:="Crescente"; END;   DEFAULT fase:="New Moon";   END; // case   // età della Luna   dd:= IP(HMS→(ph/15));   hh:= IP(FP(HMS→(ph/15))*60);   IF (hh>24) THEN dd:= dd+1; hh:= hh MOD 24; END;   mm:= IP(FP(FP(FP(HMS→(ph/15))*60))*60);   eta:= →HMS(dd+hh/60+mm/3600); // Age of the Moon   segno:= zodiaco(LongM);   PRINT;   PRINT("Effemeridi della Luna a " + m+" "+D+" "+Y+" "+ lstd);   PRINT("L Long media " + trunc_angle(L));   PRINT("M anomalia vera " + trunc_angle(M));   PRINT("λ Longitudine " + trunc_angle(LongM) + " " + segno);   PRINT("β Latitudine " + trunc_angle(LatM));   PRINT("ω Perigeo medio " + trunc_angle(ω));   PRINT("Ω Nodo ascendente " + trunc_angle(Ω));   PRINT("α Asc. r. "+ trunc_angle(α) + " δ Decl. "+ trunc_angle(δ));   PRINT("α Ascensione retta (gradi) " + trunc_angle((α*15) MOD 360));   PRINT("Az Azimuth " + trunc_angle(az) + " h Alt. " + trunc_angle(h));   PRINT("");   PRINT("Ang. orario "+ trunc_angle(H) + " " + trunc_angle(Hd));   PRINT("TS (UTC) " + trunc_angle(ts(1)));   PRINT("TS (locale) " + trunc_angle(ts(2)));   PRINT("Diff. asc." + trunc_angle(DA) + " SAD " + trunc_angle(SAD));   PRINT("P Parallasse " + trunc_angle(P));   PRINT("d Distanza Terra - Luna "+ ROUND(d*6373.044737,3) + " km");   PRINT("Semidiametro "+  trunc_angle(→HMS(s)) + " :: " + ROUND(s*6373.044737,3) + " km");   PRINT(" ");   PRINT("Culmina a: " + trunc_angle(culm) + " UTC");   PRINT("Sorge a: " + trunc_angle(sorg) + " UTC");    PRINT("Tramonta a: " + trunc_angle(tram) + " UTC");   PRINT("Fase " + trunc_angle(ph) + " " + fase + " " + ROUND(100*→HMS(ph)/360,2)+"%");   PRINT("Età " + eta);   RETURN [[L, M],[LongM, LatM],[ω,Ω],[α,δ], [az, h], [P,0],[d,s], [sorg, tram], [culm, ts(1)],[ph, eta] ]; END; planetCalc(nameplanet) BEGIN   LOCAL ω,Ω, l, v, r, x, y, z;   LOCAL Ls, Ms, b, xs, ys, b_ec;   LOCAL x1, y1, z1, R, l_ec, H, Hd;   LOCAL ec, in, a, δ, α, ts, az, h;   LOCAL DA, SAD, tram, sorg, culm, dayY;   LOCAL as, ac, zc, dz, θ, mov, segno;   // calc of Sun (as basis)   xs:= sunlist(9); ys:= sunlist(10);   Ls:= sunlist(1); // Longitudine media   Ms:= sunlist(4); // Anomalia   // end Sun  EVAL(EXPR(nameplanet+"()")); // call planet-name function   L:= planet_const(1); // Longitudine media   ω:= planet_const(2); // perielio   Ω:= planet_const(3); // nodo   M:= planet_const(4); // anomalia   v:= planet_const(5);    ec:= planet_const(6); // eccentricità   in:= planet_const(7); // inclinazione   a:= planet_const(8);  // semiasse maggiore   θ:= planet_const(9);  // θs to calc if retrogade   b_ec:= asin(sin(in)*sin(v+ω-Ω)); // argomento di latitudine del perielio   l_ec:= acos(cos(v+ω-Ω)/cos(b_ec));  // long from eclittic   IF b_ec<0 THEN l_ec:= 360 - l_ec; END; // argomento secondario   l_ec:= →HMS((l_ec + Ω)) MOD 360; // aggiunge Nodo e riporta a 2PI   r:= a*(1-ec^2)/(1+ec*cos(v)); // 9.524899/(1+ 0.0559*cos(v)); raggio vettore   x1:= r*cos(b_ec)*cos(l_ec); y1:= r*cos(b_ec)*sin(l_ec); z1:=r*sin(b_ec); // cartesian   x:= x1+xs; y:= y1+ys; z:=z1; // cohordinates cartesian->spheric   R:= sqrt(x^2+y^2+z^2);    b:=asin(z/R); // latitudine  β   // l:=→HMS(atan(y/x)); // cart->spheric longtitudine λ   l:= atan2(y,x);   IF l<0 THEN l:= l+360; END; IF l>360 THEN l:=l-360; END;   dayY:= DDAYS(IP(thisday)+1/100+1/10000, thisday); // calc day number   // δ:= →HMS((ε*sin((360*(284+dayY))/365)) MOD 360); // declinazione, formula di Cooper   δ:= →HMS((asin(cos(ε)*sin(b)+sin(ε)*cos(b)*sin(l))) ); // declinazione (coord locali)   α:= (acos(cos(b)*cos(IFTE(l>180, l-180, l))/cos(δ))); // ascensione retta IF l>180 THEN α:= α +180; END;   α:= →HMS(α/15) MOD 24; // α in hms   ts:= calcTS();  // calcola tempo siderale   H:= (ts(2) - α) MOD 24; // Angolo orario   Hd:= →HMS(15*H); // Angolo orario in gradi   zc:= sin(lat)*sin(δ)+cos(lat)*cos(δ)*cos(Hd); // calc z, alt   dz:= →HMS(acos(zc)) ; // distanza zenitale :: h:= →HMS(90 - dz); // Altezza   h:= →HMS(ATAN(zc/sqrt(1-zc^2))); // Altezza (height)   as:= cos(δ)*sin(Hd)/sin(dz); // calc az   ac:= (-cos(lat)*sin(δ)+sin(lat)*cos(δ)*cos(Hd))/sin(dz);   // az:= →HMS(ATAN(as/ac) ); // Azimuth   az:= atan2((-cos(δ) * sin(Hd)), cos(lat) *sin(δ)-sin(lat)*cos(δ)*cos(Hd));   IF az<0 THEN az:= az+360; END; IF az>360 THEN az:=az-360; END;   DA:= →HMS(asin(tan(lat)*tan(δ))); // Differenza ascensionale   SAD:= →HMS(90 + DA); // Semiarco Diurno   // SAD:= acos(-tan(lat)*tan(δ)); // Semiarco Diurno   culm:= (IFTE(α<ts(1),24+α,α) - ts(1)) MOD 24 ; // Culminazione   sorg:= (culm - SAD/15) MOD 24; // Alba (Dawn)   tram:= (culm + SAD/15) MOD 24; // Tramonto (Set)   IF abs(l-220.4)>θ THEN mov:="diretto"; ELSE mov:="retrogrado"; END;   segno:= zodiaco(l);   PRINT;   PRINT("Effemeridi di " + nameplanet + " a " + m+" "+D+" "+Y+" "+ lstd);   PRINT("L Long media " + trunc_angle(L));   PRINT("ω Perielio " + trunc_angle(ω));   PRINT("Ω Nodo " + trunc_angle(Ω));   PRINT("M Anomalia vera " + trunc_angle(M));   PRINT("l_ec Long. eclittica " + trunc_angle(l_ec));   PRINT("b_ec Lat. eclittica " + trunc_angle(b_ec) + " v " + trunc_angle(v));   PRINT("r Raggio vettore " + r);   PRINT("Cartes. elioc.: x' " + ROUND(x1,3) + " y' " + ROUND(y1,3) + " z' " + ROUND(z1,3));   PRINT("Cartes. geoc.: x " + ROUND(x,3) + " y " + ROUND(y,3) + " z " + ROUND(z,3));   PRINT("Spher.: λ Longitudine " + trunc_angle(l) + " " + segno); // spheric   PRINT("Spher.: β Latitudine " + trunc_angle(b));   PRINT("Spher.: R Distanza " + ROUND(R,3));   PRINT("α Asc. r. "+ trunc_angle(α) + " δ Decl. "+ trunc_angle(δ));   PRINT("Az Azimuth " + trunc_angle(az) + " h Alt. " + trunc_angle(h));   PRINT("");   PRINT("Ang. orario "+ trunc_angle(H) + " " + trunc_angle(Hd));   PRINT("TS (UTC) " + trunc_angle(ts(1)));   PRINT("TS (locale) " + trunc_angle(ts(2)));   PRINT("Diff. ascensionale " + trunc_angle(DA));   PRINT("SAD semiarco diurno " + trunc_angle(SAD) + " :: " + trunc_angle(→HMS(SAD/15)));   PRINT("Culmina a: " + trunc_angle(culm) + " UTC");   PRINT("Sorge a: " + trunc_angle(sorg) + " UTC");    PRINT("Tramonta a: " + trunc_angle(tram) + " UTC");   PRINT("Movimento " + mov);   RETURN [[L,ω,Ω ],[M,l_ec,b_ec], [v,r,0], [x1,y1,z1],[x,y,z],[l,b,R], [α,δ,H], [az,h,ts(1)],[sorg, tram, culm]]; END; Mercurius() BEGIN   // Costanti per Mercurio   LOCAL ω,Ω, L, M, v;   LOCAL ec, in, a, θ;   ec:= 0.205615; // eccentricità   in:= 7.003; // inclinazione   a:= 0.387098; // semiasse maggiore   θ:= 35.6;  // θs to calc if retrogade   calcN();   // costanti planetarie   L:= →HMS((229.851 + 4.09237703*N) MOD 360); // Longitudine media   ω:= →HMS((75.913 + 0.00004253*N) MOD 360); // perielio   Ω:= →HMS((47.157 + 0.00003244*N) MOD 360); // nodo   M:= (L - ω); // Anomalia vera   v:= M + 23.4372*sin(M)+2.9810*sin(2*M)+0.5396*sin(3*M)+0.1098*sin(4*M);    planet_const:= {L, ω, Ω, M, v, ec, in, a, θ};   RETURN planet_const; END; Venus() BEGIN   // Costanti per Venere   LOCAL ω,Ω, L, M, v;   LOCAL ec, in, a, θ;   ec:= 0.006816; // eccentricità   in:= 3.393636; // inclinazione   a:= 0.72333; // semiasse maggiore   θ:= 13.0;  // θs to calc if retrogade   calcN();   // costanti planetarie   L:= →HMS((206.758 + 1.60216873*N) MOD 360); // Longitudine media   ω:= →HMS((130.154 + 0.00003757*N) MOD 360); // perielio   Ω:= →HMS((75.797 + 0.00002503*N) MOD 360); // nodo   M:= (L - ω); // Anomalia vera   v:= M + 0.7811*sin(M);    planet_const:= {L, ω, Ω, M, v, ec, in, a, θ};   RETURN planet_const; END; Mars() BEGIN   // Costanti per Marte   LOCAL ω,Ω, L, M, v;   LOCAL ec, in, a, θ;   ec:= 0.093309; // eccentricità   in:= 1.850303; // inclinazione   a:= 1.523678; // semiasse maggiore   θ:= 16.8;  // θs to calc if retrogade   calcN();   // costanti planetarie   L:= →HMS((124.765 + 0.52407108*N) MOD 360); // Longitudine media   ω:= →HMS((334.251 + 0.00005038*N) MOD 360); // perielio   Ω:= →HMS((48.794 + 0.00002127*N) MOD 360); // nodo   M:= (L - ω); // Anomalia vera   v:= M + 10.608*sin(M)+0.6216*sin(2*M)+0.0504*sin(3*M)+0.0047*sin(4*M);    planet_const:= {L, ω, Ω, M, v, ec, in, a, θ};   RETURN planet_const; END; Jupiter() BEGIN   // Costanti per Giove   LOCAL ω,Ω, L, M, v;   LOCAL ec, in, a, θ;   ec:= 0.048335; // eccentricità   in:= 1.308736; // inclinazione   a:= 5.202561; // semiasse maggiore   θ:= 54.43;  // θs to calc if retrogade   calcN();   // costanti planetarie   L:= →HMS((268.350 + 0.08312942*N) MOD 360); // Longitudine media   ω:= →HMS((12.737 + 0.00004408*N) MOD 360); // perielio   Ω:= →HMS((99.453 + 0.00002767*N) MOD 360); // nodo   M:= (L - ω); // Anomalia vera   v:= M + 5372*sin(M)+0.1662*sin(2*M)+0.007*sin(3*M);    planet_const:= {L, ω, Ω, M, v, ec, in, a, θ};   RETURN planet_const; END; Saturn() BEGIN   // Costanti per Saturno   LOCAL ω,Ω, L, M, v;   LOCAL ec, in, a, θ;   ec:= 0.055892; // eccentricità   in:= 2.492519; // inclinazione   a:= 9.554747; // semiasse maggiore   θ:= 65.53;  // θs to calc if retrogade   calcN();   // costanti planetarie   L:= →HMS((278.774 + 0.03349788*N) MOD 360); // Longitudine media   ω:= →HMS((91.117 + 0.00005362*N) MOD 360); // perielio   Ω:= →HMS((112.79 + 0.00002391*N) MOD 360); // nodo   M:= (L - ω); // Anomalia vera   v:= M + 6.4023*sin(M)+0.2235*sin(2*M)+0.011*sin(3*M);    planet_const:= {L, ω, Ω, M, v, ec, in, a, θ};   RETURN planet_const; END; Uranus() BEGIN   // Costanti per Urano   LOCAL ω,Ω, L, M, v;   LOCAL ec, in, a, θ;   ec:= 0.046344; // eccentricità   in:= 0.772464; // inclinazione   a:= 19.21814; // semiasse maggiore   θ:= 73.92;  // θs to calc if retrogade   calcN();   // costanti planetarie   L:= →HMS((248.487 + 0.01176902*N) MOD 360); // Longitudine media   ω:= →HMS((171.563 + 0.00004064*N) MOD 360); // perielio   Ω:= →HMS((73.482 + 0.00001365*N) MOD 360); // nodo   M:= (L - ω); // Anomalia vera   v:= M + 5.3092*sin(M)+0.1537*sin(2*M)+0.0062*sin(3*M);    planet_const:= {L, ω, Ω, M, v, ec, in, a, θ};   RETURN planet_const; END; Neptunus() BEGIN   // Costanti per Nettuno   LOCAL ω,Ω, L, M, v;   LOCAL ec, in, a, θ;   ec:= 0.008997; // eccentricità   in:= 1.779242; // inclinazione   a:= 30.10957; // semiasse maggiore   θ:= 77.63;  // θs to calc if retrogade   calcN();   // costanti planetarie   L:= →HMS((86.652 + 0.00602015*N) MOD 360); // Longitudine media   ω:= →HMS((46.742 + 0.000039*N) MOD 360); // perielio   Ω:= →HMS((130.693 + 0.00003009*N) MOD 360); // nodo   M:= (L - ω); // Anomalia vera   v:= M + 1.031*sin(M)+0.0058*sin(2*M);    planet_const:= {L, ω, Ω, M, v, ec, in, a, θ};   RETURN planet_const; END; Pluto() BEGIN   // Costanti per Plutone   LOCAL ω,Ω, L, M, v;   LOCAL ec, in, a, θ;   ec:= 0.250236; // eccentricità   in:= 17.17047; // inclinazione   a:= 39.438712; // semiasse maggiore   θ:= 79.41;  // θs to calc if retrogade   calcN();   // costanti planetarie   // Long, perielio, nodo, non precisi, considerati costanti   // L:= →HMS((229.851 + 4.09237703*N) MOD 360); // Longitudine media   M:= →HMS((230.67 + 0.00397943*N) MOD 360); // Anomalia vera (calcolata direttamente)   Ω:= →HMS((109.06 + 0.00003823*N) MOD 360); // nodo   ω:= →HMS((114.27 + Ω) MOD 360); // perielio   L:= →HMS((M + ω) MOD 360); // L (non sicuro)      v:= M + 28.45*sin(M)+4.38*sin(2*M)+0.97*sin(3*M)+0.24*sin(4*M);    planet_const:= {L, ω, Ω, M, v, ec, in, a, θ};   RETURN planet_const; END;

∫aL√0mic (IT9CLU) :: HP Prime 50g 41CX 71b 42s 39s 35s 12C 15C - DM42, DM41X - WP34s Prime Soft. Lib
 « Next Oldest | Next Newest »

 Messages In This Thread Astronomy... - salvomic - 06-20-2015 05:53 PM RE: Astronomy... - Marcel - 06-20-2015, 07:18 PM RE: Astronomy... - salvomic - 06-20-2015, 08:34 PM RE: Astronomy... - Marcel - 06-20-2015, 10:29 PM RE: Astronomy... - salvomic - 06-21-2015, 06:38 AM RE: Astronomy... - salvomic - 06-22-2015, 12:42 PM RE: Astronomy... - Marcel - 06-22-2015, 06:26 PM RE: Astronomy... - salvomic - 06-22-2015, 06:34 PM RE: Astronomy... - Marcel - 06-22-2015, 06:48 PM RE: Astronomy... - salvomic - 06-22-2015, 07:29 PM RE: Astronomy... - salvomic - 06-23-2015, 03:16 PM RE: Astronomy... - salvomic - 06-25-2015, 05:35 PM RE: Astronomy... - salvomic - 06-27-2015, 05:51 PM

User(s) browsing this thread: 1 Guest(s)