c c Dieses File enthaelt Unterprogramme aus der alten Bibliothek libastr c Das File kann an ein Hauptprogramm anghaengt werden c C 1. Name: Subroutine BEJD U. Bastian, Okt. 1984 C C 2. Zweck: Besselsche Epochen in Julianische Daten umwandeln. C C 3. Aufruf: C CALL BEJD(BE,TJD) C C 4. Parameter: C C BE Input Besselsche Epoche (Jahre) Real*8 C C TJD Output Julianisches (Ephemeris) Datum (Tage) Real*8 C C SUBROUTINE BEJD(BE,TJD) IMPLICIT REAL*8 (A-Z) C TJD=(BE-1900.D0)*365.242198781D0+2415020.31352D0 RETURN END SUBROUTINE KSTPRZ(T0,T1,ZET,TETA,ZETA) C C BESTIMMUNG DER PRAEZESSINSWINKEL NACH DEN ANGABEN C BEI WOOLARD&CLEMENCE, S.262, UMGERECHNET UAF EPOCHE J2000.0 C UND BEZOGEN AUF 1 JULIANISCHES JAHRHUNDERT C C EINGABE: T0 = JULIANISCHES DATUM DER AUSGANGSEPOCHE C T1 = JULIANISCHES DATUM DER ENDEPOCHE C AUSGABE: PRAEZESSIONSWINKEL IM BOGENMASS C C IMPLICIT REAL*8 (A-H,O-Z) C PI = 3.141592653589793D0 C C EPOCHENDIFFERENZEN IN JULIANISCHE JAHRHUNDERTE UMRECHNEN C XJ2000 = 2451545.D0 DT0 = (T0-XJ2000)/36525.D0 DT1 = (T1-T0)/36525.D0 C ZETA = ( 2305.69970D0 + 1.3974397D0*DT0 + 0.000060D0*DT0*DT0)*DT1+ 1 ( 0.3020079D0 - 0.000270D0*DT0)*DT1*DT1 + 2 0.0179962D0*DT1*DT1*DT1 ZET = ( 2305.69970D0 + 1.3974397D0*DT0 + 0.000060D0*DT0*DT0)*DT1+ 1 ( 1.095432D0 + 0.00039D0*DT0)*DT1*DT1 + 2 0.018326D0*DT1*DT1*DT1 TETA = ( 2003.874600D0 - 0.8540465D0*DT0 - 0.000370D0*DT0*DT0)*DT1 1 +(-0.427073D0 - 0.000370D0*DT0)*DT1*DT1 - 2 0.041803D0*DT1*DT1*DT1 C C UMWANDLUNG INS BOGENMASS C X = (180.D0*3600.D0) ZETA = (ZETA*PI)/X TETA = (TETA*PI)/X ZET = (ZET *PI)/X RETURN END C SUBROUTINE PRZMTX (ZET,TETA,ZETA,PREC) C C BERECHNUNG DER PRAEZESSIONSMATRIX C C EINGABE DER PRAEZESSINSWINKEL IM BOGENMASS C IMPLICIT REAL*8 (A-H,O-Z) REAL*8 PREC(3,3) C SZET = DSIN(ZET) CZET = DCOS(ZET) SZETA = DSIN(ZETA) CZETA = DCOS(ZETA) STETA = DSIN(TETA) CTETA = DCOS(TETA) C PREC(1,1) = +CZETA*CTETA*CZET - SZETA*SZET PREC(1,2) = -SZETA*CTETA*CZET - CZETA*SZET PREC(1,3) = -STETA*CZET PREC(2,1) = +CZETA*CTETA*SZET + SZETA*CZET PREC(2,2) = -SZETA*CTETA*SZET + CZETA*CZET PREC(2,3) = -STETA*SZET PREC(3,1) = +CZETA*STETA PREC(3,2) = -SZETA*STETA PREC(3,3) = +CTETA RETURN END C SUBROUTINE DIRCOS (A,D,V) C C BERECHNUNG DES VEKTORS V DER RICHTUNGSKOSINUSSE C C EINGABE: A, D = ALPHA, DELTA IM BOGENMASS C AUSGABE: V = VEKTOR DER RICHTUNGSKOSINUSSE C IMPLICIT REAL*8 (A-H,O-Z) C REAL*8 V(3,1) C C V(1,1) = DCOS(A)*DCOS(D) V(2,1) = DSIN(A)*DCOS(D) V(3,1) = DSIN(D) RETURN END C SUBROUTINE GMPRD (A,B,R,N,M,L) C C MULTIPLIKATION ZWEIER ALLGEMEINER MATRIZEN, R = A.B C A = NAME DER 1. EINGABEMATRIX C B = NAME DER 2. EINGABEMATRIX C R = NAME DES PRODUKTES R = A.B C N = ZAHL DER ZEILEN VON A = ZAHL DER ZEILEN VON R C M = ZAHL DER SPALTEN VON A = ZAHL DER ZEILEN VON B C L = ZAHL DER SPALTEN VON B = ZAHL DER SPALTEN VON R C C IMPLICIT REAL*8 (A-H,O-Z) C DIMENSION A(N,M),B(M,L),R(N,L) C C DO 20 I = 1,N DO 20 J = 1,L R(I,J) = 0.D0 DO 20 K = 1,M 20 R(I,J) = R(I,J) + A(I,K)*B(K,J) RETURN END C SUBROUTINE ANGLE(V,A,D) C C BERECHNUNG VON ALPHA, DELTA AUS DEM VEKTOR DER RICHTUNGSKOSIN. C C EINGABE : V = VEKTOR DER RICHTUNGSKOSINUSSE C AUSGABE : A,D = ALPHA, DELTA C C IMPLICIT REAL*8 (A-H,O-Z) c REAL*8 V(3,1) integer i C C PI = 3.141592653589793D0 C C C UEBERGANG ZUM EINHEITSVEKTOR. C DIES IST DANN NOTWENDIG, WENN VERSEHENTLICH NICHT DIE C RICHTUNGSKOSINUSSE SELBST, SONDERN EIN VEKTOR, DER EIN C VIELFACHES DAVON IST, EINGEGEBEN WURDE C CALL SCLPRD (V,V,X,3) X = DSQRT(X) DO 10 I = 1,3 10 V(I,1) = V(I,1)/X C C Aufruf von d=darsin() auf d=dasin geaendert Ba. 1995 C D = DASIN(V(3,1)) A = DATAN2(V(2,1),V(1,1)) IF ( A.LT.0.D0 ) A = A + 2.D0*PI RETURN END SUBROUTINE REKON(IVERS,KATNR,INUM,YSTAR8,KENNY8,ISTAR2) C C C REKONSTRUKTION EINER BEOBACHTUNG AUS DER DIFFERENZ C C IVERS STEUERT DIE PRAEZEESION UND DAS LESEN DER DATEI C C IVERS = VERSION DES REFERENZKATALOGES, ZU DEM DIFFERENZEN C GEGEBEN SIND. C C IVERS IST EINE 2-STELLIGE ZAHL, C DEREN 1. ZIFFER AICH AUF DIE PRAEZESSION BEZIEHT, DIE C IM REFERENZKATALOG BENUTZT WURDE, UND C DIE 2. ZIFFER KENNZEICHNET ZUSAETZLICH DIE DATEI, C VON DER GELESEN WIRD. C DIE GESAMTE ZAHL = KANALNUMMER, UNTER DER DER REFERENZKATALOG C GELESEN WIRD. C C IAU(1976) - KATALOGE HABEN ALS 1. ZIFFER EINE 1 C NEWCOMB - KATALOGE HABEN ALS 1. ZIFFER EINE 2, 3 ODER 4 C STRUWE - KATALOGE HABEN ALS 1. ZIFFER EINE 5 C C FK5 - 2. ZIFFER 1 C FK4 - 1. ZIFFER 2; 2. ZIFFER 1 C FK3 - 1. ZIFFER 2; 2. ZIFFER 2 (TOTE STERNE OHNE 4000-ER NUM.) C NFK - 1. ZIFFER 2; 2. ZIFFER 3 C GC - 1. ZIFFER 2; 2. ZIFFER 4 C FK3R - 1. ZIFFER 2; 2. ZIFFER 5 C PFKSZ - 1. ZIFFER 2; 2. ZIFFER 6 C FK3 - 1. ZIFFER 2; 2. ZIFFER 7 (TOTE STERNE MIT 4000-ER NUM.) C EICHEL- 1. ZIFFER 2; 2. ZIFFER 8 C NFK-AV- 1. ZIFFER 2; 2. ZIFFER 9 C PGC - 1. ZIFFER 3; 2. ZIFFER 0 C NEWCOMB 1. ZIFFER 3; 2. ZIFFER 1 C EROS-ANH 1. ZIFFER 3; 2. ZIFFER 2 C = EROS-ANHALTSSTERNE AUS A.N. 241, P. 345 FF UND C A.N. 244, P. 385 FF. (KOPFF ET AL.) C C BEISPIEL: EIN BEOBACHTUNGSKATALOG GIBT DIFFERENZEN GEGEN DER FK4; C DANN IST IVERS = 21 UND DAS IST AUCH DIE KANALNUMMER C AUS DER DER FK4 ZU LESEN IST. C C IVERS = 21 C IMPLICIT REAL*8 (A-H,O-Z) REAL*8 YSTAR8(20) INTEGER*2 ISTAR2(30),n3(5),n6(5) INTEGER*4 KENNY8(20) CHARACTER STRALF*12,STRDEL*12 C C C open (99,file='/u/ari/s07/ari/neu') C C ICODE=ICODE+1 C WRITE(6,2449) KATNR 2449 FORMAT(/,' REKON; KATNR = ',I10) C C DIE ZIFFER N2 IN KENNY8(3) UND KENNY8(6) BERECHNEN. C FALLS DIESE ZIFFER = 2, DANN LIEGEN DIFFERENZEN GEGEN C EINEN REFERENZKATALOG VOR. C C WRITE(6,2249) KATNR 2249 FORMAT(//,' REKON; KATNR = ',I10,// ) C II = KENNY8(3)/10000 N2A = (KENNY8(3)-II*10000)/1000 C WRITE(6,2289) KENNY8(3),INUM,N2A 2289 FORMAT(3I7) IF ( (INUM.EQ.0) .AND. (N2A.EQ.2) ) WRITE(6,3000) C WRITE(6,3000) 3000 FORMAT(' ES LIEGEN DIFFERENZEN VOR ',/) II = KENNY8(6)/10000 N2D = (KENNY8(6)-II*10000)/1000 C WRITE(6,2289) KENNY8(6),INUM,N2D IF ( (INUM.EQ.0) .AND. (N2D.EQ.2) ) WRITE(6,3000) C WRITE(6,3000) C 3100 FORMAT(/,' FALSCHE KENNUNG FUER DIE DIFFERENZEN; ', F ' N2A, N2D = ', 2I5) 3110 FORMAT(' KENNUNG O.K. MINDESTENS EINE KOORD. NICHT BEOBACHTET;', F ' N2A, N2D = ',2I5,/) C IF ((N2A.NE.2).AND.(N2D.NE.2)) WRITE(6,3100) N2A,N2D C IF ((N2A.EQ.0).OR.(N2D.EQ.0)) WRITE(6,3110) N2A,N2D IF ((N2A.NE.2).AND.(N2D.NE.2)) WRITE(6,3100) N2A,N2D IF ((N2A.NE.2).AND.(N2D.NE.2)) GOTO 900 C c write (99,'(i5)')ivers C C DIFFERENZEN FESTHALTEN ZUM DRUCKEN C XDA = YSTAR8(3) XDD = YSTAR8(6) C C C C PRUEFEN OB ALPHA UND DELTA GEGEBEN SIND C C FALLS EINE KOORDINATE NICHT GEGEBEN, DANN AEQUINOX, EPOCHE, C UND KENNUNG DER ANDEREN KOORDINATE UEBERNEHMEN C IF (KENNY8(2).GE.0) GOTO 40 C ALPHA NICHT GEGEBEN, DEFINIEREN WIE DELTA: C YSTAR8( 2) = YSTAR8( 5) KENNY8( 2) = KENNY8( 5) YSTAR8(13) = YSTAR8( 5) KENNY8(13) = KENNY8( 5) 40 CONTINUE C IF (KENNY8(5).GE.0) GOTO 45 C DELTA NICHT GEGEBEN, DEFINIEREN WIE ALPHA: C YSTAR8( 5) = YSTAR8( 2) KENNY8( 5) = KENNY8( 2) YSTAR8(14) = YSTAR8( 2) KENNY8(14) = KENNY8( 2) 45 CONTINUE C IF (INUM.GT.1) GOTO 2600 IF (IVERS.EQ.11.AND.ICODE.EQ.1) WRITE(6,3011) IF (IVERS.EQ.21.AND.ICODE.EQ.1) WRITE(6,3021) IF (IVERS.EQ.22.AND.ICODE.EQ.1) WRITE(6,3022) IF (IVERS.EQ.23.AND.ICODE.EQ.1) WRITE(6,3023) IF (IVERS.EQ.24.AND.ICODE.EQ.1) WRITE(6,3024) IF (IVERS.EQ.25.AND.ICODE.EQ.1) WRITE(6,3025) IF (IVERS.EQ.26.AND.ICODE.EQ.1) WRITE(6,3026) IF (IVERS.EQ.27.AND.ICODE.EQ.1) WRITE(6,3027) IF (IVERS.EQ.28.AND.ICODE.EQ.1) WRITE(6,3028) IF (IVERS.EQ.29.AND.ICODE.EQ.1) WRITE(6,3029) IF (IVERS.EQ.30.AND.ICODE.EQ.1) WRITE(6,3030) IF (IVERS.EQ.31.AND.ICODE.EQ.1) WRITE(6,3031) IF (IVERS.EQ.32.AND.ICODE.EQ.1) WRITE(6,3032) IF (IVERS.EQ.33.AND.ICODE.EQ.1) WRITE(6,3033) C 3011 FORMAT(//,' ES WIRD DER FK5 ALS REFERENZKATALOG ANGENOMMEN') 3021 FORMAT(//,' ES WIRD DER FK4 ALS REFERENZKATALOG ANGENOMMEN') 3022 FORMAT(//,' ES WIRD DER FK3 ALS REFERENZKATALOG ANGENOMMEN') 3023 FORMAT(//,' ES WIRD DER NFK ALS REFERENZKATALOG ANGENOMMEN') 3024 FORMAT(//,' ES WIRD DER GC ALS REFERENZKATALOG ANGENOMMEN') 3025 FORMAT(//,' ES WIRD DER FK3R ALS REFERENZKATALOG ANGENOMMEN') 3026 FORMAT(//,' ES WIRD DER PFKSZ ALS REFERENZKATALOG ANGENOMMEN') 3027 FORMAT(//,' ES WIRD DER FK3 ALS REFERENZKATALOG ANGENOMMEN') 3028 FORMAT(//,' ES WIRD EICHELBERGER ALS', F ' REFERENZKATALOG ANGENOMMEN') 3029 FORMAT(//,' ES WIRD DER AV (NFK VERBESSERT) ALS', F ' REFERENZKATALOG ANGENOMMEN') 3030 FORMAT(//,' ES WIRD DER PGC ALS REFERENZKATALOG ANGENOMMEN') 3031 FORMAT(//,' ES WIRD NEWCOMB''S FC ALS REFERENZKAT. ANGENOMMEN') 3032 FORMAT(//,' ES WERDEN KOPFF''S EROS-ANHALTSSTERNE ANGENOMMEN') 3033 FORMAT(//,' ES WIRD DER KSV ALS REFERENZKATALOG ANGENOMMEN') C C if(icode.eq.1) *write(6,'(/,a,/)')' (Falls Differenzen zu verschiedenen Referenzka *talogen gegeben sind, siehe letzte Testseite)' C 2600 CONTINUE C C WRITE(6,2130 ) KATNR,(YSTAR8(I),I=1,2),XDA, C F (YSTAR8(I),I=4,5),XDD, C F (KENNY8(I),I=1,6) C C PRINT *,KENNY8(2) CALL REKPOS(IVERS,KATNR,INUM,YSTAR8,KENNY8,ISTAR2,ICODE, *n3,n6) istar2(26)=inum istar2(27)=n3(3) istar2(28)=n6(3) istar2(29)=n3(4) istar2(30)=n6(4) ystar8(19)=xda ystar8(20)=xdd C PRINT *,KENNY8(2) C C C DIE KENNUNGEN IERAS (= ISTAR2(10) ) UND IEDEC ( = ISTAR2(13) ) C ALS "222" DEFINIEREN, DENN DIE ZURUECKGEGEBENE POSITION C GILT PRINZIPIELL FUER AEQUINOX J2000.0, BASIERT AUF DER C IAU(1976) PRAEZESSION UND GILT FUER DIE JULIANISCHE EPOCHE C DER BEOBACHTUNG. C ISTAR2(10) = 222 ISTAR2(13) = 222 C C ALPHA = YSTAR8(3) DELTA = YSTAR8(6) C WRITE(6,2424) ALPHA,DELTA 2424 FORMAT(1X,2F20.8) CALL RADCHA(ALPHA,'HH MM SS.SSS',STRALF,IER,6) CALL RADCHA(DELTA,'VDD MM SS.SS',STRDEL,IER,6) C IF (ICODE.LT.2) WRITE(6,2000) 2000 FORMAT(///,' NR-REFK. EQUIN. EPOCHE DIFFER.', F ' EQUIN. EPOCHE DIFFER.', F ' KENNX8 (1 : 6)', F ' ALF(REK) DEL(REK)',//) C 2130 FORMAT(1X,I6,2X,2(2F10.1,F10.3,1X),2(2I3,I7),4X,A,4X,A,/) IF (ICODE.LT.6) FWRITE(6,2130 ) KATNR,(YSTAR8(I),I=1,2),XDA, F (YSTAR8(I),I=4,5),XDD, F (KENNY8(I),I=1,6), F STRALF,STRDEL C REWIND 21 C C INUM = INUM + 1 C RETURN C 900 STOP END C C REKONSTRUKTION EINER BEOBACHTUNG AUS DER DIFFERENZ C ================================================== C ZU EINEM ANDEREN KATALOG. C ========================= C C FILE: S.S25.ALL.ARIGFH.PROGFORT(REKPOS) C SUBROUTINE REKPOS(IVERS, KATNR,INUM,YSTAR8,KENNY8,ISTAR,ICODE, *n3,n6) C C................................................................ C C C IVERS = VERSION DES REFERENZKATALOGES, ZU DEM DIFFERENZEN C GEGEBEN SIND. C C IVERS IST EINE 2-STELLIGE ZAHL, C DEREN 1. ZIFFER AICH AUF DIE PRAEZESSION BEZIEHT, DIE C IM REFERENZKATALOG BENUTZT WURDE, UND C DIE 2. ZIFFER KENNZEICHNET ZUSAETZLICH DIE DATEI, C VON DER GELESEN WIRD. C DIE GESAMTE ZAHL = KANALNUMMER, UNTER DER DER REFERENZKATALOG C GELESEN WIRD. C C C IAU(1976) - KATALOGE HABEN ALS 1. ZIFFER EINE 1 C NEWCOMB - KATALOGE HABEN ALS 1. ZIFFER EINE 2, 3 ODER 4 C STRUWE - KATALOGE HABEN ALS 1. ZIFFER EINE 5 C C FK5 - 2. ZIFFER 1 C FK4 - 1. ZIFFER 2; 2. ZIFFER 1 C FK3 - 1. ZIFFER 2; 2. ZIFFER 2 (TOTE STERNE OHNE 4000-ER NUM.) C NFK - 1. ZIFFER 2; 2. ZIFFER 3 C GC - 1. ZIFFER 2; 2. ZIFFER 4 C FK3R - 1. ZIFFER 2; 2. ZIFFER 5 C PFKSZ - 1. ZIFFER 2; 2. ZIFFER 6 C FK3 - 1. ZIFFER 2; 2. ZIFFER 7 (TOTE STERNE MIT 4000-ER NUM.) C EICHEL- 1. ZIFFER 2; 2. ZIFFER 8 C NFK-AV- 1. ZIFFER 2; 2. ZIFFER 9 C PGC - 1. ZIFFER 3; 2. ZIFFER 0 C NEWCOMB 1. ZIFFER 3; 2. ZIFFER 1 C EROS-ANH 1. ZIFFER 3; 2. ZIFFER 2 C = EROS-ANHALTSSTERNE AUS A.N. 241, P. 345 FF UND C A.N. 244, P. 385 FF. (KOPFF ET AL.) C KSV - 1. ZIFFER 3; 2. ZIFFER 2 C C BEISPIEL: EIN BEOBACHTUNGSKATALOG GIBT DIFFERENZEN GEGEN DER FK4; C DANN IST IVERS = 21 UND DAS IST AUCH DIE KANALNUMMER C AUS DER DER FK4 ZU LESEN IST. C C KATNR : NUMMER DES BEOBACHTETEN STERNS IM REFERENZKATALOG C Z.B. FK4-NUMMER. C INUM : NUMMER DES REKORDS, DER GERADE BEARBEITET WIRD. C FUER INUM = 0 WERDEN EINIGE GROESSEN AUSGEGEBEN C YSTAR8 : WIE IM LESEPROGRAMM BESCHRIEBEN, DORT XSTAR8 C KENNY8 : WIE IM LESEPROGRAMM BESCHRIEBEN, DORT KENNX8 C C !!!!!!!! ACHTUNG !!!!!!!!!!! C C DIE REKONSTRUIERTEN POSITIONEN WERDEN IN DIE FELDER C XSTAR8(I) HINEINGESCHRIEBEN, EBENSO DIE EPOCHEN- UND C AEQUINOXANGABEN. C FALLS EINE KOORDINATE BEOBACHTET, ERFOLGT DIE AUSGABE C IN FOLGENDEN KONVENTIONEN: C C YSTAR8(1) : JULIANISCHE EPOCHE C YSTAR8(13): JULIANISCHE EPOCHE C YSTAR8(3) : BOGENMASS C C YSTAR8(4) : JULIANISCHE EPOCHE C YSTAR8(14): JULIANISCHE EPOCHE C YSTAR8(6) : BOGENMASS C C DIE FELDER "KENNY8" WERDEN ENTSPRECHEND DEFINIERT. C C................................................................ C C IMPLICIT REAL*8 (A-H,O-Z) REAL*8 YSTAR8(20) INTEGER*2 N3(5),N6(5) INTEGER*2 ISTAR(30) INTEGER*4 KENNY8(20) CHARACTER*53 TEXT C C TRANSFORMATIONS-PARAMETER ZUR UMRECHNUNG VON C 0.001 ZEITSEKUNDEN BZW. 0.01 BOGENSEKUNDEN INS C BOGENMASS DEFINIEREN. C C WRITE(6,2249) KATNR 2249 FORMAT(//, ' REKPOS',I10,//) TRA = 0.727220521664304D-7 TRD = 0.484813681109536D-7 C c TEXT = ' DER KATALOG WURDE IN BEIDEN KOORDINATEN BEOBACHTET' C C WRITE(6,2579) INUM 2579 FORMAT(//,' REKPOS; INUM = ',I10,//) C C AUS IVERS DIE ART DER PRAEZESSION BESTIMMEN, MIT C DER DER REFERENZKATALOG GERECHNET WURDE. C C KANALNUMMER (=LUN) FUERS LESEN HINSCHREIBEN C IP = IVERS/10 LUN = IVERS C C C STANDARDAEQUINOX J2000 DEFINIEREN C --------------------------------- C XJ2000 = 2451545.D0 C C DIE GROESSEN N1, N2, .... N5 IN DEN FELDERN KENNX8(3) UND C KENNX8(6) BERECHNEN UND IN N3(5), N6(5) SPEICHERN. C C WRITE(6,2244) KENNY8(3),KENNY8(6) N3(1) = KENNY8(3)/10000 N3(2) = (KENNY8(3)-N3(1)*10000)/1000 N3(3) = (KENNY8(3)-N3(1)*10000-N3(2)*1000)/100 N3(4) = (KENNY8(3)-N3(1)*10000-N3(2)*1000-N3(3)*100)/10 N3(5) = (KENNY8(3)-N3(1)*10000-N3(2)*1000-N3(3)*100-N3(4)*10) C C KENNX8(6) BERECHNEN UND IN N3(5), N6(5) SPEICHERN. C N6(1) = KENNY8(6)/10000 N6(2) = (KENNY8(6)-N6(1)*10000)/1000 N6(3) = (KENNY8(6)-N6(1)*10000-N6(2)*1000)/100 N6(4) = (KENNY8(6)-N6(1)*10000-N6(2)*1000-N6(3)*100)/10 N6(5) = (KENNY8(6)-N6(1)*10000-N6(2)*1000-N6(3)*100-N6(4)*10) C C WRITE(6,2244) N3 C WRITE(6,2244) N6 2244 FORMAT(1X,5I5) C C STERN IM REFERENZKATALOG SUCHEN, C DORT WIRD DIE REFERENZPOSITION AUF EPOCHE UND AEQUINOX DER C BEOBACHTUNG GERECHNET. C CALL GETREF(LUN,IVERS, KATNR,INUM, F YSTAR8(13),YSTAR8(14),KENNY8(13),KENNY8(14), F XEPJA,XEPJD,XAFEPA,XDFEPA,XAFEPD,XDFEPD,ICODE) C C DIFFERENZEN IN RADIANS UMRECHNEN C -------------------------------- C COS(DELTA) EVTL. ELIMINIEREN C UMRECHNUNG DER DIFFERENZEN IN 0.001 SEC UND 0.01 ARCSEC, C IF (ICODE.LT.2) WRITE(6,2500) N3(4),N6(4) 2500 FORMAT(///,' DIE EINHEITEN DER DIFFERNZEN SIND FESTGELEGT ', F ' DURCH N3(4) UND N6(4) = ',2I6 ,/, F ' VERGLEICHE HIERZU S.S25.ALL.ARIGFH.KENNUNG', F '(DIFFER)',//) C COSD = COS(XDF) C FAKTA = 1.D0 FAKTD = 1.D0 IF (N3(4).EQ.1) FAKTA = FAKTA*1000.D0 IF (N6(4).EQ.2) FAKTD = FAKTD*100.D0 IF ( (N3(4).LE.4) .AND. (N6(4).LE.4) ) GOTO 440 WRITE(6,2520) N3(4),N6(4) 2520 FORMAT(//,' EINHEITEN NICHT DEFINIERT', 2I5,//) GOTO 900 440 CONTINUE C IF(N3(5).EQ.1) FAKTA = 1.D0/COSD IF(N6(5).EQ.1) FAKTD = 1.D0/COSD C XDA = YSTAR8(3)*TRA*FAKTA XDD = YSTAR8(6)*TRD*FAKTD C WRITE(6,2329) XDA,XDD 2329 FORMAT(//,' DIFFERENZEN IM BOGENMASS: ',2F20.10,//) C C FALLS KOORDINATE NICHT BEOBACHTET, DANN DIFFERENZ 0 SETZEN C DANN STEHT IN XSTAR8(3) ODER XSTAR8(6) DIE POSITION DES C REFERENZKATALOGES. C AUSSERDEM WIRD DAS "BENUTZUNGSKENNZEICHEN" 0 GESETZT, C D.H.: DIE KOORDINATE IST NICHT BRAUCHBAR. C IF(N3(2).EQ.0) THEN XDA = 0.D0 ISTAR(11) = 0 C IF(YSTAR8(3).EQ.-8888.D0) C * TEXT=' NUR DEKLINATIONSKATALOG ' END IF IF(N6(2).EQ.0) THEN XDD = 0.D0 ISTAR(14) = 0 C IF(YSTAR8(6).EQ.-8888.D0) C * TEXT=' NUR RA KATALOG ' END IF C C IF(ICODE.LT.2)WRITE(6,'(1X,A)')TEXT C C BEOBACHTUNG REKONSTRUIEREN C -------------------------- C XAOBS = XAFEPA + XDA XDOBS = XDFEPD + XDD C YSTAR8(3) = XAOBS YSTAR8(6) = XDOBS C 2260 FORMAT(/,' RA, DEC VOR ANBRINGEN DER DIFFERENZEN',2F16.9,/, F ' DIFFERENZEN IN ALPHA UND DELTA: ',2F16.9,/, F ' RA, DEC NACH ANBRINGEN DER DIFFERENZEN',2F16.9,//) C C WRITE(6,2260) XAFEPA,XDFEPD,XDA,XDD,YSTAR8(3),YSTAR8(6) C C C TRANSFORMATION DER REKONSTRUIERTEN BEOBACHTUNG AUF J2000. C MIT DER IAU (1976) PRAEZEESION C IP = 1 C C C ALPHA TRANSFORMIEREN: C CALL TSPV2(IP,XEPJA ,XEPJA ,XJ2000,XEPJA , F XAOBS ,XDFEPA,0.D0 ,0.D0 ,0.D0,0.D0, F YA2000,Y1 ,Y2 ,Y3 ,YPAR,YVEL,IER) C C DELTA TRANSFORMIEREN: C CALL TSPV2(IP,XEPJD ,XEPJD ,XJ2000,XEPJD, F XAFEPD,XDOBS,0.D0 ,0.D0, 0.D0, 0.D0, 0 F Y1 ,YD2000,Y2 ,Y3 ,YPAR,YVEL,IER) C YSTAR8(3) = YA2000 YSTAR8(6) = YD2000 C C IF(INUM.EQ.0.AND.ICODE.EQ.1) WRITE(6,2220) XEPJA,XEPJD,XJ2000, C F XAOBS ,YA2000, XDOBS,YD2000 C WRITE(6,2220) XEPJA,XEPJD,XJ2000, C F XAOBS ,YA2000, XDOBS,YD2000 C 2220 FORMAT(//,' EPOCHENUEBERTRAGUNG AUF J2000 :'/, F ' AEQUIN = EPOCHE DER RA -BEOBACHTUNG.:',F16.3,/, F ' AEQUIN = EPOCHE DER DEC-BEOBACHTUNG.:',F16.3,/, F ' AEQUIN = EPOCHE J2000 :',F16.3,/, F ' ALPHA VOR UND NACH DER TRANSFORM.',2F16.9,/, F ' DELTA VOR UND NACH DER TRANSFORM.',2F16.9,///) C C C C ALLES IN J U L I A N I S C H E E P O C H E UMRECHNEN C =========== C KENNUNGEN ENTSPRECHEND DEFINIEREN: C ------------------------------------------------------------- C C AUFPASSEN: WENN KEINE DIFFERENZ GEGEBEN WAR, DANN IST HIER C ALS "AUS EINEM ANDEREN KATALOG" ZU DEFINIEREN. C YSTAR8( 1) = 2000.D0 YSTAR8( 2) = YSTAR8(13) C YSTAR8( 4) = 2000.D0 YSTAR8( 5) = YSTAR8(14) C KENNY8(1) = 2 KENNY8(2) = 2 KENNY8(13) = 2 IF (N3(2).EQ.2) N3(2) = 5 IF (N3(2).EQ.0) N3(2) = 3 KENNY8(3) = 10000*N3(1) + 1000*N3(2) KENNY8(4) = 2 KENNY8(5) = 2 KENNY8(14) = 2 IF (N6(2).EQ.2) N6(2) = 5 IF (N6(2).EQ.0) N6(2) = 3 KENNY8(6) = 10000*N6(1) + 1000*N6(2) 2300 FORMAT(' KENNY8( 1 : 6 ) ',/,1X,6I10) C WRITE(6,2300) (KENNY8(I),I=1,6) C RETURN 900 STOP END C C C REKONSTRUKTION EINER BEOBACHTUNG AUS DER DIFFERENZ C ================================================== C ZU EINEM ANDEREN KATALOG. C ========================= C C HIER: REFERENZKATALOG AUF EPOCHE DER BEOBACHTUNG TRANSFORMIEREN C EVTL. ELLIPTISCHE ABERRATION ELIMINIEREN. C C STEUERUNG VON PRAEZESSION UND ELL. ABERRATION C ERFOLGT UEBER IP (I-PRAEZESSION) C C DAS UNTERPROGRAMM WIRD VOM UNTERPROGRAMM REKPOS AUFGERUFEN C C FILE: S.S25.ALL.ARIGFH.PROGFORT(GETREF) C SUBROUTINE GETREF(LUN,IVERS, NRKAT,INUM, F XEPA,XEPD,KENEPA,KENEPD, F XEPJA,XEPJD,XAFEPA,XDFEPA,XAFEPD,XDFEPD,ICODE) C C................................................................ C EINGABE : C ========== C C LUN : KANALNUMMER ZUM LESEN DES REFERENZKATALOGS C IP : IM REFERENZKATALOG BENUTZTE PRAEZESSION. C IP = 1 : IAU(1976) C IP = 2, 3, 4 : NEWCOMB C IP = 5 : STRUWE C C NRKAT : NUMMER DES BEOBACHTETEN STERNS IM REFERENZKATALOG C Z.B. FK4-NUMMER. C INUM : NUMMER DES REKORDS, DER GERADE BEARBEITET WIRD. C FUER INUM = 0 WERDEN EINIGE GROESSEN AUSGEGEBEN C C XEPA : EPOCHE DER BEOBACHTUNG IN REKTASZENSION C KENEPA : ANGABE, WIE EPOCHE IM BEOB.KAT. GEGEBEN IST. C XEPD : EPOCHE DER BEOBACHTUNG IN DEKLINATION C KENEPD : ANGABE, WIE EPOCHE IM BEOB.KAT. GEGEBEN IST. C C AUSGABE: C ======== C XEPJA : EPOCHE DER BEOB. IN JULIANISCHEM DATUM C XEPJD : EPOCHE DER BEOB. IN JULIANISCHEM DATUM C XAF : ALPHA DES REF.KAT. ZU Epoche + Aequin. des Ref.Kat. C XDF : DELTA DES REF.KAT. ZU Epoche + Aequin. des Ref.Kat. C XAFEPA : ALPHA DES REF.KAT. ZU Alpha-Epoche + Aequin. der Beob. C XDFEPA : Delta DES REF.KAT. ZU Alpha-Epoche + Aequin. der Beob. C XAFEPD : ALPHA DES REF.KAT. ZU Delta-Epoche + Aequin. der Beob. C XDFEPD : Delta DES REF.KAT. ZU Delta-Epoche + Aequin. der Beob. C ........................................................... C C IMPLICIT REAL*8 (A-H,O-Z) C CHARACTER STRALF*9,STRDEL*9 C C INUM = 0 C WRITE(6,2579) INUM 2579 FORMAT(///,' GETREF; INUM = ',I10,///) C C TRANSFORMATIONS-PARAMETER ZUR UMRECHNUNG VON ZEITSEKUNDEN C BZW. BOGENSEKUNDEN INS BOGENMASS DEFINIEREN. C UMRECHNUNG VON TROPISCHEM IN JULIANISCHES JAHRHUNDERT. C DATA TRJUL /1.00002136/ DATA TRA /0.261799387799149D00/ DATA TRD /0.174532925199433D-1/ DATA TRASEC /0.727220521664304D-4/ DATA TRDSEC /0.484813681109536D-5/ C C WRITE(6,2323) XEPA,XEPD,KENEPA,KENEPD 2323 FORMAT(1X,2F20.9,2I15) C C C AEQUINOX UND EPOCHE DES REFERENZKATALOGS DEFINIEREN C --------------------------------------------------- C C IVERS = 11 = FK5 = JD 2451545.0 C IVERS = 21 = FK4 = B1950 C IVERS = 22 = FK3 = B1950 C IVERS = 23 = NFK = B1900 <======= 1900 C IVERS = 24 = GC = B1950 C IVERS = 25 = FK3R = B1950 C IVERS = 26 = PFKSZ= B1950 C IVERS = 27 = FK3 = B1950 (TOTE STERNE MIT 4000 ER NUMMERN) C IVERS = 28 = EICHELBERGER = B1925 <========= C IVERS = 29 = NFK - AV = B1900 <========= C IVERS = 30 = PGC = B1900 <========= C IVERS = 31 = NEWCOMB = B1900 <========= C IVERS = 32 = EROS-ANHALT = B1930 <========= C IVERS = 33 = KSV = B1950 <========= C C YEQREF = 1950.D0 IF (IVERS.EQ.23) YEQREF = 1900.D0 IF (IVERS.EQ.28) YEQREF = 1925.D0 IF (IVERS.EQ.29) YEQREF = 1900.D0 IF (IVERS.EQ.30) YEQREF = 1900.D0 IF (IVERS.EQ.31) YEQREF = 1900.D0 IF (IVERS.EQ.32) YEQREF = 1930.D0 IF (IVERS.EQ.33) YEQREF = 1950.D0 C CALL BEJD(YEQREF,XEQREF) C IF (IVERS.EQ.11) XEQREF = 2451545.D0 C C C WRITE(6,2323) YEQREF,XEQREF C C BENUTZTE PRAEZESSION IM REFERENZKATALOG. C ----------------------------------------- C IP RICHTIG UMDEFINIEREN, DAMIT NUR WERTE 1,2,3 ZUGELASSEN SIND C --------------------------------------------------------------- C IP = IVERS/10 IF (IP.EQ.3) IP = 2 IF (IP.EQ.4) IP = 2 IF (IP.EQ.5) IP = 3 C C PRAEZESSION IM REFERENZKATALOG C C IP = 1 : IAU(1976) C IP = 2, 3, 4 : NEWCOMB C IP = 5 : STRUWE C C 2001 FORMAT(' DER REFERENZKATALOG WIRD MIT IAU(1976) ', F ' PRAEZESSION REDUZIERT',//) 2002 FORMAT(' DER REFERENZKATALOG WIRD MIT NEWCOMB"S ', F ' PRAEZESSION REDUZIERT',//) 2005 FORMAT(' DER REFERENZKATALOG WIRD MIT STRUWE"S) ', F ' PRAEZESSION REDUZIERT',//) C C IF ((INUM.EQ.1).AND.(IP.EQ.1).AND.(ICODE.EQ.1)) WRITE(6,2001) C IF ((INUM.EQ.1).AND.(IP.EQ.2).AND.(ICODE.EQ.1)) WRITE(6,2002) C IF ((INUM.EQ.1).AND.(IP.EQ.3).AND.(ICODE.EQ.1)) WRITE(6,2002) C IF ((INUM.EQ.1).AND.(IP.EQ.4).AND.(ICODE.EQ.1)) WRITE(6,2002) C IF ((INUM.EQ.1).AND.(IP.EQ.5).AND.(ICODE.EQ.1)) WRITE(6,2005) C C C BEHANDLUNG DER ELLIPTISCHEN ABERRATION IM REFERENZKATALOG C --------------------------------------------------------- C C IEL = 1 : DER REFERENZKATALOG ENTHAELT ELLIPTISCHE C ABERRATION, DIESE WIRD DANN ELIMINIERT C IEL = 0 : DER REFERENZKATALOG ENTHAELT KEINE ELLIPTISHCE C ABERRATION, ES IST NICHTS ZU TUN. C C ES WIRD ANGENOMMEN, DASS ALLE KATALOGE, DIE MIT DER C IAU(1976)-PRAEZESSION GERECHNET SIND, KEINE ELLIPTISCHE ABERRATION C ENTHALTEN UND DASS ALLE KATALOGE, DIE NICHT MIT DER IAU(1976) C PRAEZESSION GERECHNET SIND, DIE ELLIPTISCHE ABERRRATION ENTHALTEN. C C IEL = 0 IF (IP.GT.1) IEL = 1 C 2100 FORMAT(//,' DER REFERENZKATALOG ENTHAELT KEINE ELLIPT. ABERR.',/, F ' DIESE BLEIBT AUCH DRAUSSEN. ',///) 2101 FORMAT(//,' DER REFERENZKATALOG ENTHAELT ELLIPTISCHE ', F ' ABERRATION; DIESE WIRD ELIMINIERT',//) C C IF ((INUM.EQ.1) .AND.(IEL.EQ.0).AND.(ICODE.EQ.1)) WRITE(6,2100) C IF ((INUM.EQ.1) .AND.(IEL.EQ.1).AND.(ICODE.EQ.1)) WRITE(6,2101) C C C BEOBACHTUNGSEPOCHE IN JULIANISCHEM DATUM UND ALS C BESSELSCHE EPOCHE RECHNEN C ------------------------------------------------ C C FALLS EINE EPOCHE NICHT GEGEBEN IST, DANN GLEICH DER C ANDEREN SETZEN. C C C WRITE(6,2257) KENEPA,KENEPD,XEPA,XEPD 2257 FORMAT(' GETREF: KENEPA, KENEPD, XEPA,XEPD:', 2I10,2F15.6) IF (KENEPA.EQ.0) XEPA = XEPD IF (KENEPD.EQ.0) XEPD = XEPA C WRITE(6,2258) KENEPA,KENEPD,XEPA,XEPD 2258 FORMAT(' GETREF: KENEPA, KENEPD, XEPA,XEPD:', 2I10,2F15.6) XEPJA = XEPA XEPJD = XEPD IF((KENEPA.EQ.1).OR.(KENEPD.EQ.1)) CALL BEJD(XEPA,XEPJA) IF((KENEPA.EQ.1).OR.(KENEPD.EQ.1)) CALL BEJD(XEPD,XEPJD) IF((KENEPA.EQ.2).OR.(KENEPD.EQ.2)) F XEPJA = 2451545.D0 + (XEPA-2000.D0)*365.25D0 IF((KENEPA.EQ.2).OR.(KENEPD.EQ.2)) F XEPJD = 2451545.D0 + (XEPD-2000.D0)*365.25D0 C C WRITE(6,2323) XEPA ,XEPD C WRITE(6,2323) XEPJA,XEPJD C C STERN IM REFERENZKATALOG AUFSUCHEN C ---------------------------------- C 500 CONTINUE C C FK4 AUS APFS-DATEI LESEN C C WRITE(6,2449) LUN 2449 FORMAT(' LUN =',I10) READ(LUN,1000,END=900) NRF,STRALF,XMYA,STRDEL,XMYD,XPAR,XVEL C PRINT *, NRF,XMYA,XMYD 1000 FORMAT(4X,I9,8X,A,F7.3,A,F6.2,F3.3,3X,F5.1) C 3377 FORMAT(/,' GETREF; NRKAT, NRF = ',2I10,/) C WRITE(6,3377) NRKAT,NRF IF (NRF.NE.NRKAT) GOTO 500 C CALL CHARAD(STRALF(1:9),'HHMMSSSSS',XAF,IER,6) CALL CHARAD(STRDEL(1:9),'VDDMMSSSS',XDF,IER,6) C C WRITE(6,2327) NRF,NRKAT 2327 FORMAT(1X,2I10) C IF (NRF.NE.NRKAT) GOTO 500 C C C XAF = (XAH + XAM/60.D0 + XAS/3600.D0)*TRA C XDF = (XDG + XDM/60.D0 + XDS/3600.D0)*TRD C IF(NV.EQ.MINUS) XDF = -XDF C WRITE(6,2323) XAF,XDF C C EIGENBEWEGUNGEN SIND GGF.AUF JULIANISCHES JAHRHUNDERT C UMZURECHNEN, DA DAS IM U.P. TSPV2 SO ERWARTET WIRD. C ES WIRD ANGENOMMEN, DASS IM FALLE IP UNGLEICH 1 DIE C EIGENBEWEGUNGEN AUF TROPISCHES JAHRHUNDERT BEZOGEN SIND. C XMYAF = XMYA*TRASEC XMYDF = XMYD*TRDSEC IF (IP.NE.1) XMYAF = XMYA*TRASEC*TRJUL IF (IP.NE.1) XMYDF = XMYD*TRDSEC*TRJUL C WRITE(6,2323) XMYAF,XMYDF C XPAR = XPAR*TRDSEC C WRITE(6,2323) XPAR C C C JETZT ELLIPTISCHE ABERRATION ELIMINIEREN C ---------------------------------------- C C WRITE(6,1717) XAF,XDF IF (IEL.EQ.0) GOTO 570 CALL KSTETR(YEQREF,C0,D0,EPS) CALL ETERM (XAF,XDF,C0,D0,EPS,CA,CD) XAF = XAF - CA XDF = XDF - CD C WRITE(6,1717) CA,CD C WRITE(6,1717) XAF,XDF C XXXA = CA/TRASEC C XXXD = CD/TRDSEC C WRITE(6,1717) XXXA,XXXD 1717 FORMAT(1X,2F20.10) 570 CONTINUE C C C TRANSFORMATION AUF EPOCHE UND AEQUINOX DER BEOBACHTUNG C ------------------------------------------------------ C C ALPHA TRANSFORMIEREN: C C WRITE(6,2345) IP 2345 FORMAT(' TRANSF. EP. BEOB.',I10) CALL TSPV2(IP,XEQREF,XEQREF,XEPJA,XEPJA, F XAF,XDF,XMYAF,XMYDF,XPAR,XVEL, F XAFEPA,XDFEPA,YMYAF,Y2 ,YPAR,YVEL,IER) C C DELTA TRANSFORMIEREN: C CALL TSPV2(IP,XEQREF,XEQREF,XEPJD,XEPJD, F XAF,XDF,XMYAF,XMYDF,XPAR,XVEL, F XAFEPD,XDFEPD,Y1 ,YMYDF,YPAR,YVEL,IER) C C IF(INUM.EQ.1.AND.ICODE.EQ.1) WRITE(6,2200) XEQREF,XEPJA,XEPJD, C F XAF,XAFEPA,XDF,XDFEPD,XMYAF,YMYAF, C F XMYDF,YMYDF,XPAR,YPAR,XVEL,YVEL C 2200 FORMAT(//,' TRANSFORMATION AUF EPOCHE DER BEOB. :',/, F ' AEQUIN = EPOCHE DES REFERENZKAT. :',F16.3,/, F ' AEQUIN = EPOCHE DER RA-BEOBACHTUNG. :',F16.3,/, F ' AEQUIN = EPOCHE DER DEC-BEOBACHTUNG.:',F16.3,/, F ' ALPHA VOR UND NACH DER TRANSFORM. :',2F16.9,/, F ' DELTA VOR UND NACH DER TRANSFORM. :',2F16.9,/, F ' MY-A VOR UND NACH DER TRANSFORM. :',2F16.9,/, F ' MY-D VOR UND NACH DER TRANSFORM. :',2F16.9,/, F ' PAR. VOR UND NACH DER TRANSFORM. :',2F16.9,/, F ' RV. VOR UND NACH DER TRANSFORM. :',2F16.9,///) C C XAF = YAF C XDF = YDF C C C BEOABACHTUNGSEPOCHE ALS JULIANISCHE EPOCHE AUSGEBEN C XEPA = 2000.D0 + (XEPJA-2451545.D0)/365.25D0 XEPD = 2000.D0 + (XEPJD-2451545.D0)/365.25D0 C WRITE(6,2323) XEPA,XEPD C RETURN C 900 CONTINUE c WRITE(6,1100) NRKAT 1100 FORMAT(////, ' DATEIENDE IM REFERENZKATALOG ERREICHT ',/, F ' FOLGENDE KATALOGNUMMER WURDE NICHT MEHR GEFUNDEN', F I10,///) STOP END SUBROUTINE TSPV2 (IP,XEQ,XEP,YEQ,YEP,XA,XD,XM,XMS,XP0,XVEL0, 1 YA,YD,YM,YMS,YP0,YVEL0,IER) C C EPOCHEN- UND AEQUINOX-UEBERTRAGUNG VON OERTERN U. EIGENBEWEGUNGEN C C EOCHEN-UEBERTRAGUNG MIT BENUTZUNG DER RAUMGESCHWINDIGKEIT C.............................................................. C ................... C C IP : STEUERT, WELCHE PRAEZESSIONSWERTE GENOMMEN WERDEN. C IP= 1 : NEUE IAU(1976) WERTE; C IP = 2 : NEWCOMBS PRAEZESSION; C IP = 3 : STRUVE'S PRAEZESSION; C C EINGABE: XEQ = JULIANISCHES DATUM DES AUSGANGS-AEQUINOXES C YEQ = JULIANISCHES DATUM DES END-AEQUINOXES C XEP = JULIANISCHES DATUM DER AUSGANGS-EPOCHE C YEP = JULIANISCHES DATUM DER AUSGANGS-EPOCHE C XA = REKTASZENSION FOER EPOCHE XEP, AEQUINOX XEQ C XD = DEKLINATION FUER EPOCHE XEP, AEQUINOX XEQ C XM = E.B. IN ALPHA FUER EPOCHE YEP, AEQUINOX XEQ C XMS = E.B. IN DELTA FUER EPOCHE XEP, AEQUINOX XEQ C XP0 = PARALLAXE ZUR EPOCHE XEP IM BOGENMASS C XVEL0 = RADIALGESCHWINDIGKEIT ZUR EPOCH XEP C C AUSGABE: YA = REKTASZENSION FUER EPOCHE YEP, AEQUINOX YEQ C YD = DEKLINATION FUER EPOCHE YEP, AEQUINOX YEQ C YM = E.B. IN ALPHA FUER EPOCHE YEP, AEQUINOX YEQ C YMS = E.B. IN DELTA FUER EPOCHE YEP, AEQUINOX YEQ C YP0 = PARALLAXE ZUR EPOCHE XEP IM BOGENMASS C YVEL0 = RADIALGESCHWINDIGKEIT ZUR EPOCH XEP C RETURN-CODE IER: = 0: RAD.VEL. UND PAR. UNGLEICH O C = 1: RAD.VEL. ODER PAR. GLEICH 0 C C ALLE WINKEL IM BOGENMASS C ZEITEINHEIT = JULIANISCHES JAHRHUNDERT C C IMPLICIT REAL*8 (A-H,O-Z) C REAL*8 PX(3,3) C C C C EPOCHENUEBERTRAGUNG DER OERTER UND EIGENBEWEGUNGEN C CALL TPMSPC (XEP,YEP,XA ,XD ,XM ,XMS ,XP0,XVEL0, 1 XA2,XD2,XM2,XMS2,YP0,YVEL0,IER) C C PRAEZESSIONSUEBERTRAGUNG VON OERTERN UND EIGENBEWEGUNGEN C CALL KPRZ (IP,XEQ,YEQ,ZET,TETA,ZETA) CALL PRZMTX( ZET,TETA,ZETA,PX) C CALL PRECES(XA2,XD2,XM2,XMS2,PX,XA3,XD3,XM3,XMS3) C YA = XA3 YD = XD3 YM = XM3 YMS=XMS3 RETURN END SUBROUTINE KSTETR(BY,C0,D0,EPS) C BESTIMMUNG DES "ELLIPTISCHEN TEILS DER ABERRATION DAY NUMBERS" C EPS WIRD BEZOGEN AUF DAS MITTLERE AEQUNOX FUER BEGINN DES C BESSELSCHEN JAHRES DER EPOCHE TA. DIE ZEIT IST IN C EPHEMERIDENZEIT (JULIANISCH) ZU RECHNEN (EXPL. SUPPL. P.98, 489) C BY = BESSELSCHES JAHR UND TEILE DAVON, FUER DAS GERECHNET WIRD C DT = DIFFERENZ ZWISCHEN 1900.0 UND BY IN JULIANISCHEN JAHRHUNDERTE C DTEPS = DIFFERENZ ZWISCHEN 1900.0 UND BY IN TROPISCHEN JAHRHUNDERT C KAPPA = ABERRATIONS-KONSTANTE = 20.496" (AB 1968) C C IMPLICIT REAL*8 (A-H,O-Z) REAL*8 T,BY,DT,DTEPS,PERL,EXC,EPS,C0,D0,KAPPA T = 3.1415926535897932D0/(180.D0*3600.D0) KAPPA = 20.496D0 DT = 365.2422D0*(BY-1900.D0)/36525.D0 DTEPS = (BY - 1900.D0)/100.D0 PERL = 281.22083 3D0 + 1.71917 5D0*DT + 0.00045 2778D0*DT*DT + 1 0.00000 33333D0*DT*DT*DT EXC = 0.01675 104D0 - 0.00004 180D0*DT - 0.00000 0126D0*DT*DT EPS = 23.45229 4D0 - 0.01301 25D0*DTEPS - 1 0.00000 164D0*DTEPS*DTEPS + 0.00000 0503D0*DTEPS*DTEPS*DTEPS PERL = PERL - 180.D0 PERL = PERL*3600.D0*T EPS = EPS*3600.D0*T KAPPA = KAPPA*T C0 = KAPPA*EXC*DCOS(PERL)*DCOS(EPS) D0 = KAPPA*EXC*DSIN(PERL) RETURN END SUBROUTINE ETERM(ALF,DEL,C0,D0,EPS,C1,C2) C C BERECHNUNG DER KORREKTIONEN, DIE VON DER ELLIPTIZITAET C DER ERDBAHN HERRUEHREN (LINEARER TEIL) C C LITERATUR: WOOLARD & CLEMENCE, SPHER. ASTR., P. 114 C EXPL. SUPPL. P. 48,151 C EPS = SCHIEFE DER EKLIPTIK, WIRD IN GESONDERTEM U.P. BERECHNET C C0,D0 = "ELLIPTISCHER TEIL DER ABERRATION DAY NUMBERS", WERDEN C KSTETR BERECHNET C C,CS,D,DS = STAR NUMBERS C A,D = ORT IM BOGENMASS C C1,C2 = E-TERM DER ELLIPTISCHEN ABERRATION (LINEARER TEIL) C IMPLICIT REAL*8 (A-H,O-Z) REAL*8 ALF,DEL,C1,C2,C0,D0,C,CS,D,DS,EPS C = DCOS(ALF)/DCOS(DEL) D = DSIN(ALF)/DCOS(DEL) CS = DTAN(EPS)*DCOS(DEL) - DSIN(ALF)*DSIN(DEL) DS = DCOS(ALF)*DSIN(DEL) C1 = C*C0 + D*D0 C2 = CS*C0 + DS*D0 RETURN END SUBROUTINE TPMSPC(T0,T1,A0,D0,XMY0,XMYS0,PAR0,VEL0, 1 A1,D1,XMY1,XMYS1,PAR1,VEL1,IER) C C C EPOCHENUEBERTRAGUNG VON OERTERN UND EIGENBEWEGUNGEN UNTER DER C ANNAHME KONSTANTER RAUMGESCHWINDIGKEIT. C C EINGABEDATEN: C.................. C T0, T1: ANFANGS-UND ENDEPOCHE IN JULIANISCHEM DATUM C A0, D0: ORT ZUR ANFANGSEPOCHE, IM BOGENMASS C XMY0,XMYS0: EIG.BEW.KOMP. IM BOGENMASS, PRO JUL. JAHRH. C PAR: PARALLAXE IM BOGENMASS ZUR ANFANGS-EPOCHE C VEL0: RADIALGSCHW. IN KM/SEC ZUR ANFANGSEPOCHE C C AUSGABEDATEN: C.................. C WIE EINGABEDATEN, MIT '1' STATT '0' C RETURN-CODE IER: = 0: RAD.VEL. UND PAR. UNGLEICH 0 C..................... C = 1: RAD.VEL. ODER PAR. GLEICH 0 C C IMPLICIT REAL*8 (A-H,O-Z) C DIMENSION R0(3),R0P(3),XA(3),XD(3),R(3) C PI = 3.141592653589793D0 TR = PI/180.D0 TR1 = 3600.D0/TR TR2 = 1.D0/TR1 IER = 0 F1 = 21.094953D0 F2 = 4.8481368D-6 1700 FORMAT(/,1X,3D25.10) C C WRITE(6,1700) T0,T1 C WRITE(6,1700) A0,D0 C WRITE(6,1700) XMY0,XMYS0 X1 = DABS(PAR0) X2 = DABS(VEL0) IF ( ( X1 .LT.1.D-15) .OR. ( X2 .LT.1.D-15) ) IER = 1 C C C EIGENBEWEGUNGEN UND PARALLAXE IN BOGENSEKUNDEN UMRECHNEN. C AA = A0 DA = D0 XMYA = XMY0 XMYSA = XMYS0 PARA = PAR0 VELA = VEL0 XMYA = XMYA*TR1 XMYSA=XMYSA*TR1 PARA = PARA*TR1 C WRITE(6,1700) XMYA,XMYSA,PARA C C C VEKTOR ZUM STERN UND GESCHWINDIGKEITSVEKTOR BERECHNEN, C REDUZIERT AUF DIE EINHEITSSPHAERE. C CALL DIRCOS(AA,DA,R0) C WRITE(6,1720) AA,DA C WRITE(6,1720) R0 1720 FORMAT(/,1X,3D25.13) CALL DLDA(AA,DA,XA) CALL DLDD(AA,DA,XD) C WRITE(6,1720) XA C WRITE(6,1720) XD DO 20 I = 1,3 20 R0P(I) = XMYA*XA(I)+XMYSA*XD(I) + F1*PARA*VELA*R0(I) C WRITE(6,1720) R0P C C C EPOCHENDIFFERENZ IN JULIANISCHEN JAHRHUNDERTEN: C DT = (T1-T0)/36525.D0 C WRITE(6,1720) DT C C C VEKTOR ZUM STERN ZUR ENDEPOCHE C DO 80 I = 1,3 80 R(I) = R0(I) + R0P(I)*DT*F2 C WRITE(6,1720) R C C C ENTFERNUNG DES STERNS IM BOGENMASS ZUR ENDEPOCHE: C CALL SCLPRD(R,R,E,3) E = DSQRT(E) C WRITE(6,1720) E PAR1 = PARA/E C WRITE(6,1720) PAR1 C C C UEBERGANG ZUM EINHEITSVEKTOR C ALPHA , DELTA AUR ENDEPOCHE BERECHNEN. C DO 100 I = 1,3 100 R(I) = R(I)/E C WRITE(6,1700) R C CALL ANGLE(R,A1,D1) C WRITE(6,1700) A1,D1 C C C NEUE RADIALGESCHWINDIGKEIT ZUR ENDEPOCHE: C CALL SCLPRD(R0P,R,VEL1,3) IF (IER.EQ.0) VEL1 = VEL1/(F1*PAR1) IF (IER.EQ.1) VEL1 = VELA C WRITE(6,1700) VEL1 C C C NEUE EIGENBEW.KOMP. ZUR ENDEPOCHE BERECHNEN: C CD = DCOS(D1) CALL DLDA(A1,D1,XA) CALL DLDD(A1,D1,XD) CALL SCLPRD(R0P,XA,XMY1,3) CALL SCLPRD(R0P,XD,XMYS1,3) C WRITE(6,1700) XMY1,XMYS1 XMY1 = XMY1/(E*CD*CD) XMYS1 = XMYS1/E C WRITE(6,1700) XMY1,XMYS1 C C EIGENBEW. UND PARALLAXE IN BOGENMASS UMWANDELN. C XMY1 = XMY1*TR2 XMYS1 = XMYS1*TR2 PAR1 = PAR1*TR2 C RETURN END SUBROUTINE KPRZ(IP,T0,T1,ZET,TETA,ZETA) C C BESTIMMUNG DER PRAEZESSIONSWINKEL NACH VERSCHIEDENEN KONSTANTEN C JE NACH WAHL DES PARAMETERS IP: C C IP = 1 : NEUE IAU(1976) PRAEZESSION, C IP = 2 : NEWCOMB'S PRAEZESSION, C IP = 3 : STRUVE'S PRAEZESSION. C C EINGABE: T0 = JULIANISCHES DATUM DER AUSGANGSEPOCHE C T1 = JULIANISCHES DATUM DER ENDEPOCHE C AUSGABE: PRAEZESSIONSWINKEL IM BOGENMASS C C IMPLICIT REAL*8 (A-H,O-Z) C PI = 3.141592653589793D0 C C EPOCHENDIFFERENZEN IN JULIANISCHE JAHRHUNDERTE UMRECHNEN C XJ2000 = 2451545.D0 DT0 = (T0-XJ2000)/36525.D0 DT1 = (T1-T0)/36525.D0 C C................................................................ IF (IP.NE.1) GOTO 200 C C BESTIMMUNG DER PRAEZESSINSWINKEL NACH DEN ANGABEN C BEI LIESKE, ASTRON.&ASTROPHYS., 73, 282-284,(1979). C ...................................................... C ZETA = ( 2306.2181D0 + 1.39656D0*DT0 - 0.000139D0*DT0*DT0)*DT1 + 1 ( 0.30188D0 - 0.000344D0*DT0)*DT1*DT1 + 2 0.017998D0*DT1*DT1*DT1 ZET = ( 2306.2181D0 + 1.39656D0*DT0 - 0.000139D0*DT0*DT0)*DT1 + 1 ( 1.09468D0 + 0.000066D0*DT0)*DT1*DT1 + 2 0.018203D0*DT1*DT1*DT1 TETA = ( 2004.3109D0 - 0.85330D0*DT0 - 0.000217D0*DT0*DT0)*DT1 + 1 (-0.42665D0 - 0.000217D0*DT0)*DT1*DT1 - 2 0.041833D0*DT1*DT1*DT1 C GOTO 700 C 200 CONTINUE C C .................................................... IF (IP.NE.2) GOTO 300 C C PRAEZESSIONSWINKEL NACH NEWCOMBS PRAEZESSIONSWERTEN C .................................................... C ZETA = ( 2305.69970D0 + 1.3974397D0*DT0 + 0.000060D0*DT0*DT0)*DT1+ 1 ( 0.3020079D0 - 0.000270D0*DT0)*DT1*DT1 + 2 0.0179962D0*DT1*DT1*DT1 ZET = ( 2305.69970D0 + 1.3974397D0*DT0 + 0.000060D0*DT0*DT0)*DT1+ 1 ( 1.095432D0 + 0.00039D0*DT0)*DT1*DT1 + 2 0.018326D0*DT1*DT1*DT1 TETA = ( 2003.874600D0 - 0.8540465D0*DT0 - 0.000370D0*DT0*DT0)*DT1 1 +(-0.427073D0 - 0.000370D0*DT0)*DT1*DT1 - 2 0.041803D0*DT1*DT1*DT1 C GOTO 700 C 300 CONTINUE C C ........................................................... IF (IP.NE.3) GOTO 400 C C BESTIMMUNG DER PRAEZESSIONSWINKEL NACH STRUVE'S PRAEZESSION C ........................................................... C ZETA = ( 2306.013272D0 + 1.424560853 D0*DT0)*DT1 + 1 0.316313512D0*DT1*DT1 ZET = ( 2306.013272D0 + 1.424560853 D0*DT0)*DT1 + 1 1.108347345D0*DT1*DT1 TETA = ( 2004.386800D0 - 0.8630368665D0*DT0)*DT1 + 1 (-0.4315184332D0)*DT1*DT1 C 400 CONTINUE C 700 CONTINUE C C UMWANDLUNG INS BOGENMASS C X = (180.D0*3600.D0) ZETA = (ZETA*PI)/X TETA = (TETA*PI)/X ZET = (ZET *PI)/X RETURN END C SUBROUTINE PRECES (XA,XD,XMUE,XMUES,PX,YA,YD,YMUE,YMUES) C C PRAEZESSIONSUEBERTRAGUNG VON EIGENBEWEGUNGEN UND OERTERN C EINE PRAEZESSION DER EIGENBEWEGUNGEN ERFORDERT STETS AUCH DIE C PRAEZESSION DER OERTER. C C EINGABE: XA = ALPHA FUER AUSGANGSAEQUINOX C XD = DELTA FUER AUSGANGSAEQUINOX C XMUE = E.B. IN ALPHA FUER AUSGANGSAEQUINOX C XMUES = E.B. IN DELTA FUER AUSGANGSAEQUINOX C PX = PRAEZESSIONSMATRIX C AUSGABE: YA = ALPHA FUER ENDAEQUINOX C YD = DELTA FUER ENDAEQUINOX C YMUE = E.B. IN ALPHA FUER ENDAEQUINOX C YMUES = E.B. IN DELTA FUER ENDAEQUINOX C ALLE WINKEL IM BOGENMASS C C IMPLICIT REAL*8 (A-H,O-Z) REAL*8 PX(3,3),VA(3,1),VD(3,1),XRES(3,1),YRES(3,1),ZRES(3,1) C C PRAEZESSION DES ORTES C CALL PRCPOS (XA,XD,YA,YD,PX) C C BERECHNUNG DER VEKTOREN DER ABGELEITETEN RICHTUNGSKOSINUSSE C CALL DLDA (XA,XD,VA) CALL DLDD (XA,XD,VD) C C PRAEZESSION DER EIGENBEWEGUNGEN C CALL GMPRD (PX,VA,XRES,3,3,1) CALL SMPY (XRES,XMUE,YRES,3,1) CALL GMPRD (PX,VD,XRES,3,3,1) CALL SMPY (XRES,XMUES,ZRES,3,1) CALL GMADD (YRES,ZRES,XRES,3,1) C C AUFLOESUNG NACH DEN EIGENBEWEGUNGS-KOMPONENETEN C CALL DLDA (YA,YD,VA) CALL SCLPRD (XRES,VA,YMUE,3) C YMUE = YMUE/(DCOS(YD)*DCOS(YD)) C CALL DLDD (YA,YD,VD) CALL SCLPRD (XRES,VD,YMUES,3) C YMUES = YMUES RETURN END C SUBROUTINE DLDA (A,D,V) C C BERECHNUNG DES NACH ALPHA ABGELEITETEN VEKTORS DER RICHTUNGSKOS. C C EINGABE : A, D = ALPHA, DELTA (IM BOGENMASS) C AUSGABE : V = VEKTOR DER ABGELETETEN RICHTUNGSKOSINUSSE C C IMPLICIT REAL*8 (A-H,O-Z) REAL*8 V(3,1) C V(1,1) = -DCOS(D)*DSIN(A) V(2,1) = +DCOS(D)*DCOS(A) V(3,1) = 0.D0 RETURN END C SUBROUTINE DLDD (A,D,V) C C BERECHNUNG DES NACH DELTA ABGELEITETEN VEKTORS DER C RICHTUNGSKOSINUSSE C C EINGABE: A, D = ALPHA,DELTA (IM BOGENMASS) C AUSGABE: V = VEKTOR DER RICHTUNGSKOSINUSSE C IMPLICIT REAL*8 (A-H,O-Z) REAL*8 V(3,1) C V(1,1) = -DSIN(D)*DCOS(A) V(2,1) = -DSIN(D)*DSIN(A) V(3,1) = +DCOS(D) RETURN END C SUBROUTINE PRCPOS (AOLD,DOLD,ANEW,DNEW,PX) C C PRAEZESSIONSUEBERTRAGUNG VON OERTERN C C EINGABE: AOLD = REKTASZENSION FUER AUSGANGSAEQUINOX C DOLD = DEKLINATION FUER AUSGANGSAEQUINOX C PX = PRAEZESSIONSMATRIX C AUSGABE: ANEW = REKTASZENSION FUER ENDAEQUINOX C DNEW = DEKLINATION FUER ENDAEQUINOX C ALLE WINKEL IM BOGENMASS C C IMPLICIT REAL*8 (A-H,O-Z) REAL*8 PX(3,3),RCOLD(3,1),RCNEW(3,1) C C C BERECHNUNG DER RICHTUNGSKOSINUSSE C CALL DIRCOS (AOLD,DOLD,RCOLD) C C PRAEZESSIONSUEBERTRAGUNG C C CALL GMPRD (PX,RCOLD,RCNEW,3,3,1) C C BERECHNUNG VON ALPHA, DELTA AUS DEN RICHTUNGSKOSINUSSEN C CALL ANGLE (RCNEW,ANEW,DNEW) C RETURN END C SUBROUTINE SMPY (A,C,R,N,M) C C MULTIPLIKATION EINER MATRIX MIT EINEM SKALAR C C A = NAME DER EINGABEMATRIX C C = SKALAR C R = NAME DER AUSGABEMATRIX, R = C.A C N = ZAHL DER ZEILEN VON A C M = ZAHL DER SPALTEN VON A C C IMPLICIT REAL*8 (A-H,O-Z) C DIMENSION A(N,M),R(N,M) C C DO 10 I = 1,N DO 10 J = 1,M 10 R(I,J) = C*A(I,J) RETURN END C SUBROUTINE GMADD (A,B,R,N,M) C C ADDITION ZWEIER MATRIZEN : R = A + B C C A = NAME DES 1. SUMMANDEN C B = NAME DES 2. SUMMANDEN C R = NAME DER SUMME: R = A + B C N = ZAHL DER ZEILEN VON A, B, R C M = ZAHL DER SPALTEN VON A, B, R C C IMPLICIT REAL*8 (A-H,O-Z) C DIMENSION A(N,M),B(N,M),R(N,M) C C DO 20 I = 1,N DO 20 J = 1,M 20 R(I,J) = A(I,J) + B(I,J) RETURN END SUBROUTINE JDBE(JD,BE) C C UMWANDLUNG DES JULIANISCHEN DATUMS IN BESSELSCHE EPOCHE C REAL*8 JD,BE C BE = (JD-2415020.31352D0)/365.242198781D0 + 1900.D0 C RETURN END SUBROUTINE GNOMAP(ALPHA0,DELTA0,ALPHA,DELTA,X,Y) IMPLICIT REAL*8 (A-H,O-Z) C C BERECHNUNG DER FENSTERGROESSE BZW.EINER STERNPOSITION BEI C GNOMONISCHER ABBILDUNG (SIMULATION EINER PHOTOPLATTE) C IN RECHTWINKLIGEN KOORDINATEN. C DIE STERNPOSITION ALPHA0,DELTA0 (ZU IDENTIFIZIERENDER STERN) IST C BERUEHRPUNKT DER TANGENTIALEBENE UND DAMIT URSPRUNG DES C KOORDINATENSYSTEMS MIT +X-ACHSE IN OESTLICHER ALPHA-RICHTUNG UND C +Y-ACHSE IN NOERDLICHER DELTA-RICHTUNG . C DIE STERNPOSITION ALPHA,DELTA (IDENTIFIKATIONSSTERN) WIRD IN X,Y C UEBERGEFUEHRT . C DIE EINHEITEN FUER X,Y ENTSPRECHEN DENEN AUF DER EINHEITSKUGEL, C SIND ALSO SOZUSAGEN IM BOGENMASS GEGEBEN. C DALFA=ALPHA-ALPHA0 SIND0=SIN(DELTA0) COSD0=COS(DELTA0) SIND =SIN(DELTA) COSD =COS(DELTA) COSALF=SIND0*SIND+COSD0*COSD*COS(DALFA) COSID0=COSD0*SIND SICOD0=SIND0*COSD X=COSD/COSALF*SIN(DALFA) Y=(COSID0-SICOD0*COS(DALFA))/COSALF RETURN END C SUBROUTINE TPMPOS (XA,XD,XMUE,XMUES,DMUFOR,T0,T1,YA,YD) C C EIGENBEWEGUNGSUEBERTRAGUNG VON MITTLEREN OERTERN C C EINGABE: XA = ALPHA ZUR EPOCHE T0 C XD = DELTA ZUR EPOCHE T0 C XMUE = E.B. IN ALPHA ZUR EPOCHE T0 C XMUES = E.B. IN DELTA ZUR EPOCHE T0 C DMUFOR = FORESHORTENING-EFFEKT C T0 = AUSGANGSEPOCHE C T1 = ENDEPOCHE C AUSGABE: YA = ALPHA ZUR ENDEPOCHE C YD = DELTA ZUR ENDEPOCHE C C ALLE WINKEL IM BOGENMASS C ZEITEINHEIT JULIANISCHES JAHRHUNDERT C C IMPLICIT REAL*8 (A-H,O-Z) REAL*8 XM(3,1) C PI = 3.141592653589793D0 C VORBEREITUNGEN DT = (T1-T0)/36525.D0 C SD = DSIN(XD) CD = DCOS(XD) SA = DSIN(XA) CA = DCOS(XA) C XMUE0 = DSQRT(XMUE*XMUE*CD*CD + XMUES*XMUES) XMUTOT = XMUE0 + 0.5D0*DMUFOR*DT C FALLS DER BETRAG DER GESAMTEN EIGENBEWEGUNG EXTREM KLEIN IST, C WIRD SIE ZU NULL GESETZT UND DAS PROGRAM OHNEWEITERE RECHNUNG C SOFORT VERLASSEN C IF(XMUE0.LT.1.D-20) *THEN YA=XA YD=XD RETURN END IF C BEGINN DER EIGENTLICHEN RECHNUNG C C SPSI = (XMUE*CD)/XMUE0 CPSI = XMUES/XMUE0 C SMT = DSIN(XMUTOT*DT) CMT = DCOS(XMUTOT*DT) C C BERECHNUNG DES VEKTOR "M.X" NACH MUELLER, S.115, FORMEL 4.94 C XM(1,1) = -SD*CA*CPSI*SMT - SA*SPSI*SMT + CD*CA*CMT XM(2,1) = -SD*SA*CPSI*SMT + CA*SPSI*SMT + CD*SA*CMT XM(3,1) = +CD*CPSI*SMT + SD*CMT C C BERECHNUNG DER WINKEL AUS DEN RICHTUNGSKOSINUSSEN C CALL ANGLE (XM,YA,YD) RETURN END C 1. Name: Subroutine EQUIN U.Bastian, Okt. 1984 C C 2. Zweck: Praezedieren von aequatorialen Koordin. und Eigenbew. C C 3. Aufruf: C CALL EQUIN(EQU,EQUNEU,ALP,DEL,PMA,PMD,ALPNEU,DELNEU,PMANEU,PMDNEU) C C 4. Parameter: C C EQU Input Ausgangsaequinox (Besselsche Epoche in Jahren) Real*8 C EQUNEU INPUT ZIELAEQUINOX (BESSELSCHE EPOCHE IN JAHREN) REAL*8 C ALP Input Rektaszension (Radians) Real*8 C DEL Input Deklination (Radians) Real*8 C PMA Input Eigenbewegg. (Proper Motion) in Rektasz. Real*8 C PMD Input Eigenbewegg. (Proper Motion) in Dekl. Real*8 C C ALPNEU Output transformierte Rektasz. (Radians) Real*8 C DELNEU Output transformierte Deklin. (Radians) Real*8 C PMANEU Output transformierte Eigenbew. in Rekt. Real*8 C PMDNEU Output transformierte Eigenbew. in Dekl. Real*8 C C 5. Notwendige Statements im rufenden Programm: Keine C C 6. Notwendige Steuerkarten: Keine C C 7. Weitere benoetigte Subroutines: C BEJD, KSTPRZ, PRZMTX, PRECES, C PRCPOS,DLDA,DLDD,SCLPRD,GMPRD,SMPY,GMADD C 8. Genaue Funktionsbeschreibung: C Die Subroutine praezediert einen Satz von aequatorialen Koordinaten C und zugehoerigen Eigenbewegungen. Ausgangs- und Zielaequinox werden C als Besselsche Epoche in Jahren gegeben. Es wird die Newcombsche C ("alte") Praezession verwendet. Exakte Koordinatendrehung, d.h. C keinerlei Naeherungen oder Reihenentwicklungen! C SUBROUTINE *EQUIN(EQU,EQUNEU,ALP,DEL,PMA,PMD,ALPNEU,DELNEU,PMANEU,PMDNEU) C IMPLICIT REAL*8 (A-H,O-Z) DIMENSION PX(3,3) C BESSELSCHE EPOCHEN IN JUL.DAT. VERWANDELN CALL BEJD(EQU,TJD) CALL BEJD(EQUNEU,TJDNEU) C PRAEZESSIONSWINKEL BERECHNEN CALL KSTPRZ(TJD,TJDNEU,ZET,TETA,ZETA) C PRAEZESSIONSMATRIX DARAUS CALL PRZMTX(ZET,TETA,ZETA,PX) C MATRIX ANBRINGEN CALL PRECES(ALP,DEL,PMA,PMD,PX, * ALPNEU,DELNEU,PMANEU,PMDNEU) C RETURN END C 1. Name: Subroutine EQUIN5 U.Bastian, Mai 1985 C C 2. Zweck: Praezedieren von aequatorialen Koordin. und Eigenbew. C ACHTUNG: Neue (IAU 1976) Praezession C C 3. Aufruf: C CALL EQUIN5(TJD,TJDNEU,ALP,DEL,PMA,PMD,ALPNEU,DELNEU,PMANEU,PMDNEU) C C 4. Parameter: C C TJD INPUT AUSGANGSAEQUINOX (JULIANSICHES DATUM) REAL*8 C TJDNEU INPUT ZIELAEQUINOX (JULIANISCHES DATUM) REAL*8 C ALP Input Rektaszension (Radians) Real*8 C DEL Input Deklination (Radians) Real*8 C PMA Input Eigenbewegg. (Proper Motion) in Rektasz. Real*8 C PMD Input Eigenbewegg. (Proper Motion) in Dekl. Real*8 C C ALPNEU Output transformierte Rektasz. (Radians) Real*8 C DELNEU Output transformierte Deklin. (Radians) Real*8 C PMANEU Output transformierte Eigenbew. in Rekt. Real*8 C PMDNEU Output transformierte Eigenbew. in Dekl. Real*8 C C 5. Notwendige Statements im rufenden Programm: Keine C C 6. Notwendige Steuerkarten: Keine C C 7. Weitere benoetigte Subroutines: C BEJD, KSTPRZ, PRZMT5, PRECES, C PRCPOS,DLDA,DLDD,SCLPRD,GMPRD,SMPY,GMADD C 8. Genaue Funktionsbeschreibung: C Die Subroutine praezediert einen Satz von aequatorialen Koordinaten C und zugehoerigen Eigenbewegungen. Ausgangs- und Zielaequinox werden C ALS JULIANISCHES DATUM GEGEBEN. ES WIRD DIE NEUE (IAU 1976) C PRAEZESSION VERWENDET. EXAKTE KOORDINATENDREHUNG, D.H. C keinerlei Naeherungen oder Reihenentwicklungen! C SUBROUTINE *EQUIN5(TJD,TJDNEU,ALP,DEL,PMA,PMD,ALPNEU,DELNEU,PMANEU,PMDNEU) C IMPLICIT REAL*8 (A-H,O-Z) DIMENSION PX(3,3) C BESSELSCHE EPOCHEN IN JUL.DAT. VERWANDELN C CALL BEJD(EQU,TJD) C CALL BEJD(EQUNEU,TJDNEU) C PRAEZESSIONSWINKEL BERECHNEN CALL KSTPR5(TJD,TJDNEU,ZET,TETA,ZETA) C PRAEZESSIONSMATRIX DARAUS CALL PRZMTX(ZET,TETA,ZETA,PX) C MATRIX ANBRINGEN CALL PRECES(ALP,DEL,PMA,PMD,PX, * ALPNEU,DELNEU,PMANEU,PMDNEU) C RETURN END SUBROUTINE KSTPR5(T0,T1,ZET,TETA,ZETA) C C BESTIMMUNG DER PRAEZESSINSWINKEL NACH DEN ANGABEN C BEI LIESKE, ASTRON.&ASTROPHYS., 73, 282-284,(1979). C C EINGABE: T0 = JULIANISCHES DATUM DER AUSGANGSEPOCHE C T1 = JULIANISCHES DATUM DER ENDEPOCHE C AUSGABE: PRAEZESSIONSWINKEL IM BOGENMASS C C IMPLICIT REAL*8 (A-H,O-Z) C PI = 3.141592653589793D0 C C EPOCHENDIFFERENZEN IN JULIANISCHE JAHRHUNDERTE UMRECHNEN C XJ2000 = 2451545.D0 DT0 = (T0-XJ2000)/36525.D0 DT1 = (T1-T0)/36525.D0 C ZETA = ( 2306.2181D0 + 1.39656D0*DT0 - 0.000139D0*DT0*DT0)*DT1 + 1 ( 0.30188D0 - 0.000344D0*DT0)*DT1*DT1 + 2 0.017998D0*DT1*DT1*DT1 ZET = ( 2306.2181D0 + 1.39656D0*DT0 - 0.000139D0*DT0*DT0)*DT1 + 1 ( 1.09468D0 + 0.000066D0*DT0)*DT1*DT1 + 2 0.018203D0*DT1*DT1*DT1 TETA = ( 2004.3109D0 - 0.85330D0*DT0 - 0.000217D0*DT0*DT0)*DT1 + 1 (-0.42665D0 - 0.000217D0*DT0)*DT1*DT1 - 2 0.041833D0*DT1*DT1*DT1 C C UMWANDLUNG INS BOGENMASS C X = (180.D0*3600.D0) ZETA = (ZETA*PI)/X TETA = (TETA*PI)/X ZET = (ZET *PI)/X RETURN END SUBROUTINE KALJUL (JAHR,MONAT,TAG,ST,MIN,SEC,JD) C C BERECHNUNG DES JULIANISCHEN DATUMS JD FUER DIE EPOCHE JAHR, C MONAT, TAG ST, MIN SEC IN BUERGERLICHER ZEIT. C C IMPLICIT REAL*8 (A-Z) INTEGER I,J,K,L,M,N C C I = JAHR J = MONAT K = TAG JD = K - 32075 + 1461*(I+4800+(J-14)/12)/4 + 367*(J-2-(J-14)/12*12 F )/12 - 3*((I+4900+(J-14)/12)/100)/4 + F ST/24.D0 + MIN/1440.D0 + SEC/86400.D0 - 0.5D0 RETURN END