Kalkulátory v astronomii

Diskuze o všem, co ještě patří do astrotechniky a jinam se nevešlo
Odpovědět
Uživatelský avatar
Psion
Příspěvky: 13039
Registrován: 02. 01. 2001, 05:03
Bydliště: Praha
Věk: 63
Kontaktovat uživatele:

Re: Kalkulátory v astronomii

#166

Příspěvek od Psion »

Upravil jsem program pro Commodore 128 pro zvýšení přesnosti výpočtů, protože v původním algoritmu je přesnost v RA +/- 1.5m a v DEC 1.5°. Nový algoritmus je výrazně přesnější. Výpočty jsou pro ekvinokcium 1950.

Původní algoritmus:
RA= 7h 37m 8.77s, DEC= 4° 5' 43.95''

Nový algoritmus:
RA= 7h 36m 42.26s, DEC= 4° 7' 52.47''

Správný výsledek:
RA= 7h 36m 43s DEC= +4°07' 50"

Kód: Vybrat vše

100 REM EFEMERIDY - UPRAVENO PRO COMMODORE 128 BASIC 7
105 REM VELKE PERIODICKE UHLY JSOU REDUKOVANY NA 0 AZ 360 STUPNU
110 DIM A(29)
120 DG=3.141592653/180
125 DEF FN ARCS(X1)=ATN(X1/SQR(1-X1*X1))*(180/3.141592653)
126 DEF FN MD(X1)=X1-360*INT(X1/360)

130 INPUT "TYP:";AS$
140 IF AS$="P" THEN INPUT "PERIAPSE=";A:GOTO 160
150 INPUT "AXIS=";A:INPUT "ECCENTRICITY=";B
160 INPUT "INCLINATION=";C
162 INPUT "LONGITUDE ASC.=";D
164 INPUT "ARGUMENT PERIAPSE=";O
170 PLAY "A"

180 GOSUB 540
190 N=Q
200 E=23.4457889
210 F=COS(D*DG)
220 G=SIN(D*DG)*COS(E*DG)
230 H=G*TAN(E*DG)
240 I=-SIN(D*DG)*COS(C*DG)
250 J=COS(D*DG)*COS(E*DG)*COS(C*DG)
252 J=J-SIN(C*DG)*SIN(E*DG)
260 K=COS(D*DG)*COS(C*DG)*SIN(E*DG)
262 K=K+SIN(C*DG)*COS(E*DG)

270 Y=F:X=I:GOSUB 500:F=U:I=V
295 Y=G:X=J:GOSUB 500:G=U:J=V
320 Y=H:X=K:GOSUB 500:H=U:K=V

340 PLAY "AB"
350 GOSUB 540
360 GOSUB 720
370 M=(Q-N)/A/SQR(A)
380 IF AS$="E" THEN GOSUB 580
390 IF AS$="P" THEN GOSUB 650

401 X=SIN((F+O+P)*DG)*I*L+R
402 Y=SIN((G+O+P)*DG)*J*L+S
403 Z=SIN((H+O+P)*DG)*K*L+T
420 GOSUB 500
430 D=SQR(V*V+Z*Z)
440 IF U<0 THEN U=U+360

450 REM PREVOD RA A DEC
460 E=U/15
470 D=FN ARCS(Z/D)
471 S1=INT(E)
472 M1=INT((E-S1)*60)
473 SE=((E-S1)*60-M1)*60
475 S2=INT(D)
476 M2=INT((D-S2)*60)
477 DE=((D-S2)*60-M2)*60
479 PRINT "RA=";S1;"HOD";M1;"MIN";SE;"SEK"
480 PRINT "DEC=";S2;"STU";ABS(M2);"MIN";ABS(DE);"VTE"
490 GOTO 340

500 U=ATN(Y/X)*(180/3.141592653)
510 U=U+90*(1-X/ABS(X))*Y/ABS(Y)
520 V=SQR(X*X+Y*Y)
530 RETURN

540 INPUT "DEN=";T:INPUT "MESIC=";S:INPUT "ROK=";R
550 IF S<=2 THEN R=R-1:S=S+12
560 Q=INT(365.25*R)+INT(30.6001*(S+1))+T-679003.5
562 Q=Q-INT(R/100)+INT(INT(R/100)/4)
570 RETURN

580 M=.985609*M
582 M=FN MD(M)
590 E=M
600 L=(M+(180*B/3.141592653)*SIN(E*DG)-E)
602 L=L/(1-B*COS(E*DG))
610 IF ABS(L)>1E-7 THEN E=E+L:GOTO 600
620 P=2*ATN(SQR((1+B)/(1-B))*TAN(E*DG/2))
622 P=P*180/3.141592653
630 L=A*(1-B*COS(E*DG))
640 RETURN

650 M=.0364911624*M
660 E=0
670 D=(2*E*E*E+M)/(E*E+1)/3
680 IF ABS(D-E)>1E-6 THEN E=D:GOTO 670
690 P=2*ATN(E)*(180/3.141592653)
700 L=A*(1+E*E)
710 RETURN

720 P=15020
730 L=36524.2199
740 R=(Q-P-.313)/L
750 S=(33282.423-Q)/L
760 L=3600
770 D=(((.018*S+.302)*S+2304.25+1.396*R)*S)/L
780 M=(((.42*S-.426)*S+2004.682-.853*R)*S)/L
790 E=((.001*S+.791)*S*S)/L+D

800 R=COS(D*DG)*COS(E*DG)*COS(M*DG)
802 R=R-SIN(E*DG)*SIN(D*DG)
810 S=SIN(D*DG)*COS(E*DG)
812 S=S+COS(D*DG)*SIN(E*DG)*COS(M*DG)
820 T=COS(D*DG)*SIN(M*DG)
830 U=-S
840 V=COS(D*DG)*COS(E*DG)
842 V=V-SIN(D*DG)*SIN(E*DG)*COS(M*DG)
850 W=-SIN(D*DG)*SIN(M*DG)
860 X=-T
870 Y=-SIN(E*DG)*SIN(M*DG)
880 Z=COS(M*DG)

890 D=(Q-P)/36525
900 M=((-33E-7*D-15E-5)*D+35999.04975)*D
902 M=(M+358.47583)/360
910 M=(M-INT(M))*360

920 A(28)=Q
930 Q=B
940 B=(126E-9*D-418E-7)*D+.01675104
950 A(27)=A
960 A=1+2E-7
970 GOSUB 590
975 X1=M

980 M=(3025E-7*D+36000.76892)*D
982 M=M+279.69668+P-X1
984 M=FN MD(M)

990 B=Q
1000 A=A(27)

1001 C=FN MD(153.23+22518.7541*D)
1002 E=FN MD(216.571+45037.5082*D)
1003 Q=FN MD(312.69+32964.3577*D)
1004 P=FN MD((-144E-5*D+445267.1142)*D+350.74)

1005 X1=543*SIN(C*DG)+1575*SIN(E*DG)
1006 X1=X1+1627*SIN(Q*DG)+3076*COS(P*DG)
1007 Y1=FN MD(65928.7155*D+353.4)
1008 X1=X1+927*SIN(Y1*DG)
1009 L=L+X1*1E-8

1010 X1=134*COS(C*DG)+154*COS(E*DG)
1011 X1=X1+200*COS(Q*DG)+179*SIN(P*DG)
1012 Y1=FN MD(231.19+20.2*D)
1013 X1=X1+178*SIN(Y1*DG)
1014 M=FN MD(M+X1*1E-5)

1020 Q=((503E-9*D-164E-8)*D-.0130125)*D
1022 Q=Q+23.452294

1030 C=L*COS(M*DG)
1040 D=L*SIN(M*DG)*COS(Q*DG)
1050 E=L*SIN(M*DG)*SIN(Q*DG)

1060 R=R*C+U*D+X*E
1062 S=S*C+V*D+Y*E
1064 T=T*C+W*D+Z*E

1080 Q=A(28)
1090 RETURN
Odpovědět