Quadriviale Kuriositäten/ Die Stele vom Rocher des Doms/ SaekularerAufgang

Java-Programm zur Berechnung der säkularen Aufgänge von Gestirnen

/*
  Source file: SaekularerAufgang.java
  Programm: saekularer Aufgang von Gestirnen
  Urspruenglicher Autor: Uwe Pilz, Leipzig (Python)
  Ueberarbeitung und Ergaenzung: Markus Bautsch, Berlin (Java)
  Licence: public domain
  Date: 25. April 2022
  Version: 1.3 (15 May 2023)
  Programming language: Java
 */

/*
  Dieses Java-Programm berechnet unter der Beruecksichtigung der Eigenbewegung und der atmosphaerischen Refraktion den saekularen Aufgang eines Gestirns.
  Astronomische Parameter: aequatoriale Koordinaten und Eigenbewegung des Gestirns zur astronomischen Epoche
  Geographische Parameter: Beobachtungsjahr, geographische Breite, Hoehe ueber Normalnull
 */

public class SaekularerAufgang
{
	/* Deklaration der Klassenvariablen */

	/* Beobachtungszeitpunkt */
	private static final long firstObservationYear = -3000; // erstes Beobachtungsjahr
	private static final long yearIncrement = 500;          // Schrittweite fuer die Jahre
	private static final long lastObservationYear = 2000;   // letztes Beobachtungsjahr
	private static final long epoch = 2000;                 // astronomische Epoche J2000

	/* Beobachtungsort */
	private static double phi = java.lang.Math.toRadians (43.95); // geographischer Breitengrad in Radiant
	private static double REarth = 6371000;                       // Erdradius in Metern
	private static double elevation = 55;                         // Hoehe ueber Normalnull in Metern
	private static double chineInDegrees = java.lang.Math.toDegrees (java.lang.Math.acos (REarth / (REarth + elevation))); // Kimmwinkel in Bogengrad
	private static double altitudeIncrement = 1;                  // Schrittweite fuer die Hoehe ueber dem Horizont in Bogengrad

	/* Gestirn */
	private static java.lang.String name = "Phakt"; 
	private static double apparentMagnitude = 2.6; // Scheinbare Helligkeit in Groessenklassen
	/* Aequatoriale Koordinaten zur Epoche in Radiant */
	private static double alpha0 = java.lang.Math.toRadians (15 * (5 + 39/60.0 + 38.95/3600.0)); // Rektaszension
	private static double delta0 = java.lang.Math.toRadians (-(34 + 4/60.0 + 26.6/3600.0));      // Deklination
	/* Jaehrliche Veraenderung der aequatorialen Koordinaten zur Epoche in Millibogensekunden */
	private static double mu_RA =   4.6; // Eigenbewegung Rektaszension
	private static double mu_DE = -29.9; // Eigenbewegung Deklination

	/* Aequatoriale Koordinaten zum Beobachtungszeitpunkt in Radiant */
	private static double alpha; // Rektaszension zum Beobachtungszeitpunkt
	private static double delta; // Deklination zum Beobachtungszeitpunkt

	/*
	  Berechnung der Extinktion in Abhaengigkeit von der astronomischen Hoehe "heightInDegrees" als Reduktion der scheinbaren Helligkeit
	  Formel nach https://asterism.org/resources/atmospheric-extinction-and-refraction/
	 */
	private static double extinction (double heightInDegrees)
	{
		double z = java.lang.Math.toRadians (90 - heightInDegrees); // Zenithoehe
		double cosz = java.lang.Math.cos (z);
		double airMass = 1 / (cosz + 0.025 * java.lang.Math.exp (-11 * cosz)); // Luftmassen
		double extinction = 0.28 * airMass; // Extinktion
		return extinction;
	}

	/* Berechnung der Praezession fuer das Jahr "year" nach IAU 2006 */
	private static void precession (long year)
	{
		/* Jahre in Bezug auf die Epoche */
		double jahre = year - epoch;

		/* Eigenbewegung beruecksichtigen */
		alpha = alpha0 + jahre * java.lang.Math.toRadians (mu_RA / 3600 / 1000);
		delta = delta0 + jahre * java.lang.Math.toRadians (mu_DE / 3600 / 1000);

		/* Praezessionsformel arbeitet mit Jahrhunderten */
		double jahrhunderte = jahre / 100.0;

		/* IAU 2006 */
		double zeta = (((((-0.0000003173 * jahrhunderte - 0.000005971) * jahrhunderte + 0.01801828)
				* jahrhunderte + 0.2988499) * jahrhunderte + 2306.083227) * jahrhunderte + 2.650545);
	    double z = (((((-0.0000002904 * jahrhunderte - 0.000028596) * jahrhunderte + 0.01826837)
	    		* jahrhunderte + 1.0927348) * jahrhunderte + 2306.077181) * jahrhunderte - 2.650545);
	    double theta = ((((-0.00000011274 * jahrhunderte - 0.000007089) * jahrhunderte - 0.04182264)
	    		* jahrhunderte - 0.4294934) * jahrhunderte + 2004.191903) * jahrhunderte;

		/* Umrechnung nach Bogenmass */
		zeta = java.lang.Math.toRadians (zeta / 3600);
		z = java.lang.Math.toRadians (z / 3600);
		theta = java.lang.Math.toRadians (theta / 3600);

		double A = java.lang.Math.cos (delta) * java.lang.Math.sin (alpha + zeta);
		double B = java.lang.Math.cos (theta) * java.lang.Math.cos (delta) * java.lang.Math.cos (alpha + zeta) - java.lang.Math.sin (theta) * java.lang.Math.sin (delta);
		double C = java.lang.Math.sin (theta) * java.lang.Math.cos (delta) * java.lang.Math.cos (alpha + zeta) + java.lang.Math.cos (theta) * java.lang.Math.sin (delta);

		alpha = java.lang.Math.atan2 (A, B) + z;
		if (alpha < 0)
		{
			alpha = alpha + 2 * java.lang.Math.PI;
		}
		delta = java.lang.Math.asin (C);
	}

	/* Unterprogramm fuer die Ausgabe einer Hoehe ueber dem Horizont */
	private static boolean printAltitudeLine (double altitude, double cosAscentAzimuth)
	{
		java.lang.System.out.printf ("%7.1f", altitude);
		java.lang.System.out.print ("° | ");

		boolean altitudeObservable = java.lang.Math.abs (cosAscentAzimuth) < 1; // Hoehe wird vom Gestirn erreicht
		if (altitudeObservable)
		{
			double azimuth = java.lang.Math.toDegrees (java.lang.Math.acos (cosAscentAzimuth));
			double visualMagnitude = apparentMagnitude + extinction (altitude);
			java.lang.System.out.printf ("%8.1f", azimuth);
			java.lang.System.out.printf ("° | %9.1f", 360 - azimuth);
			java.lang.System.out.printf ("° | %9.1f", visualMagnitude);
			java.lang.System.out.print ("m");
			java.lang.System.out.println ();
		}
		else
		{
			java.lang.System.out.println ("wird nicht erreicht");
		}
		return altitudeObservable;
	}

	/* Unterprogramm fuer die Berechnung der Hoehen ueber dem Horizont */
	private static void computeAltitudes ()
	{
		java.lang.System.out.println ("Aufgang und Untergang am südlichen Meridian:");
		java.lang.System.out.println ("         | östlicher | westlicher |  effektive");
		java.lang.System.out.println ("    Höhe |    Azimut |     Azimut | Helligkeit");

		double altitude = 0; // Aktuelle Hoehe ueber dem Horizont in Bogengrad
		boolean altitudeObservable = false;
		do // fuer alle Hoehen ueber dem Horizont
		{
			/* Aktuelle Hoehe currentAltitude ueber dem Horizont korrigiert fuer Kimmtiefe in Bogengrad */
			double currentAltitude = altitude - chineInDegrees;

			/* Refraktion R in Bogengrad */
			double arg = currentAltitude + (7.31 / (currentAltitude + 4.4));
			double rho = (1 / java.lang.Math.tan (java.lang.Math.toRadians (arg))) / 60; // Umrechnung von Bogenminuten nach Bogengrad

			/* Wirkliche Hoehe ohne Atmosphaere, dort ist der Aufgangsazimut */
			double trueAltitude = java.lang.Math.toRadians (currentAltitude - rho);

			/* Azimut des Aufgangs auf dem Horizont in Bogengrad */
			double cosAscentAzimuth = (java.lang.Math.sin (delta) - java.lang.Math.sin (phi) * java.lang.Math.sin (trueAltitude)) / (java.lang.Math.cos (phi) * java.lang.Math.cos (trueAltitude));

			/* Ergebnisse fuer Hoehe ausgeben */
			altitudeObservable = printAltitudeLine (altitude, cosAscentAzimuth);

			altitude = altitude + altitudeIncrement;
		} while (altitudeObservable);
		java.lang.System.out.println ();
	}

	/*  Unterprogramm fuer die Berechnung eines Jahres */
	private static void computeYear (long currentYear)
	{
		precession (currentYear); // Berechnung der aequatorialen Koordinaten zum Zielzeitpunkt inklusive der Eigenbewegung

		java.lang.System.out.println ("Beobachtungsjahr = " + currentYear + " zur Epoche J" + epoch);
		java.lang.System.out.printf ("Rektaszension      = %5.1f", java.lang.Math.toDegrees (alpha));
		java.lang.System.out.println ("°");
		java.lang.System.out.printf ("Deklination        = %5.1f", java.lang.Math.toDegrees (delta));
		java.lang.System.out.println ("°");
		java.lang.System.out.printf ("Nördlichste Breite = %5.1f", 90 + java.lang.Math.toDegrees (delta));
		java.lang.System.out.println ("°");
		computeAltitudes ();
	}

	/* Schleife fuer die zu berechnenden Jahre */
	private static void computeYears ()
	{
		java.lang.System.out.println ("Berechnung der säkularen Aufgänge des Gestirns " + name);
		java.lang.System.out.println ("Scheinbare Helligkeit = " + apparentMagnitude + "m");
		java.lang.System.out.println ("Breitengrad der Beobachtung = " + java.lang.Math.toDegrees (phi) + "°");
		java.lang.System.out.println ("Höhe des Beobachtungspunkts über Normalnull = " + elevation + " Meter\n");
		java.lang.System.out.println ("Epoche = J" + epoch);
		java.lang.System.out.println ("Rektaszension zur Epoche = " + java.lang.Math.toDegrees (alpha0) + "°");
		java.lang.System.out.println ("Deklination   zur Epoche = " + java.lang.Math.toDegrees (delta0) + "°");
		java.lang.System.out.println ("Eigenbewegung Rektaszension pro Jahr, mu_RA = " + mu_RA + " mas/a");
		java.lang.System.out.println ("Eigenbewegung Deklination pro Jahr, mu_DE = " + mu_DE + " mas/a");
		java.lang.System.out.println ();

		long currentYear = firstObservationYear; // Aktuelles Jahr der Berechnung
		/* Die Schleife berechnet gegebenenfalls mehrere Jahre in der Vergangenheit */
		do // fuer alle Jahre "currentYear" bis zum Erreichen des Jahres der Epoche
		{
			computeYear (currentYear);
			currentYear = currentYear + yearIncrement;
		} while (currentYear <= lastObservationYear);
	}

	/* Aufrufbares Hauptprogramm "main" */
	public static void main (java.lang.String [] arguments)
	{
		computeYears ();
	}
}

Textausgabe

Für den Stern Phakt wird beispielsweise mit den oben angegebenen Parametern diese Textausgabe erzeugt:

Berechnung der säkularen Aufgänge des Gestirns Phakt

Scheinbare Helligkeit = 2.6m

Breitengrad der Beobachtung = 43.95°

Höhe des Beobachtungspunkts über Normalnull = 55.0 Meter

Epoche = J2000

Rektaszension zur Epoche = 84.91229166666668°

Deklination zur Epoche = -34.07405555555556°

Eigenbewegung Rektaszension pro Jahr, mu_RA = 4.6 mas/a

Eigenbewegung Deklination pro Jahr, mu_DE = -29.9 mas/a

Beobachtungsjahr = -3000 zur Epoche J2000

Rektaszension = 41,2°

Deklination = -46,6°

Nördlichste Breite = 43,4°