var DtoR = Math.PI/180.0;        //Grad -> Bogenmass
var StoR = Math.PI/180.0/3600.0; //Bogensekunden -> Bogenmass
var RtoD = 180.0/Math.PI;        //Bogenmass -> Grad
var RtoH = (12 / Math.PI);       //Bogenmass -> Stunden

var pi2  = 2.0*Math.PI;          // 2 * Pi


/*****************************************************************************/
/* Name:    modpi                                                            */
/* Type:    function                                                         */
/* Purpose: reduce an angle to the interval (-pi, pi).                       */
/*****************************************************************************/

function modpi(x)
{
  return x - Math.floor(x / pi2 + 0.5) * pi2;
}   

/*****************************************************************************/
/* Name:    modpi2                                                           */
/* Typ:     Funktion                                                         */
/* Zweck:   Einen Winkel auf das Inervall [0,2*Pi] reduzieren                */
/* Argumente:                                                                */
/*          x: Winkel im Bogenmass                                           */
/* Rueckgabewert:                                                            */
/*          Winkel im Inervall [0,2*Pi]                                      */
/*****************************************************************************/
function modpi2(x)
{
    return (x - Math.floor(x / pi2) * pi2);
}                                            


/*****************************************************************************/
/* Name:    sgn                                                              */
/* Typ:     Funktion                                                         */
/* Zweck:   Implementation der Signum-Funktion                               */
/* Argumente:                                                                */
/*          a: Zahl                                                          */
/* Rueckgabewert:                                                            */
/*          0 wenn a=0, -1 wenn a < 0, 1 wenn a >0                           */
/*****************************************************************************/
function sgn(a)
{
    if (a==0)
      return 0;

    return (a>0) ? 1 : -1;
}


/*****************************************************************************/
/* Name:    SunEklPos                                                        */
/* Typ:     Funktion                                                         */
/* Zweck:   Berechnet geozentrische ekliptikale Laenge und Breite der Sonne  */
/*          fuer ein gegebenes Datum                                         */
/* Argumente:                                                                */
/*          JD: Datum                                                        */
/* Rueckgabewert:                                                            */
/*          Vektor mit ekliptikaler Laenge und Breite im Bogenmass           */
/*****************************************************************************/
function SunEklPos(JD)
{
  T = (JD - 2415020.0) / 36525.0;  //Juli. Jhd seit 1900

  //Bahnelemente der Erde, bezogen auf Aequinokitkum des Datums
  a = 1.0;
  ex= 0.016751 - 0.000042*T;  //Exzentizitaet
  i = 0;
  w = (101.2208 + 1.7192*T)*DtoR;
  M = (358.4758 + 35640*T + 359.0498*T)*DtoR;

  //Iteratives Loesen der Keplergleichung
  E = M;
  E1 = M + ex*Math.sin(E);
  
  while ( Math.abs(E1-E) > 1e-10 )
  {
    E = E1;
    E1 = M + ex*Math.sin(E);
  }

  r = a*(1-ex*Math.cos(E1));  // Entferung Erde-Sonne
  v = 2 * Math.atan( Math.sqrt((1+ex)/(1-ex)) * Math.tan(E1/2) );

  u = v+w;

  //Transformation Bahebene -> helioz. rechtwinklig ekliptikal
  x = r*(Math.cos(u) );
  y = r*( Math.sin(u)*Math.cos(i) );
  z = 0;

  //Transformation Heliozentrisch ekl. -> Geozentrisch ekl.
  x = -x;  y = -y;  

  lo = Math.atan(y/x);    //geoz. ekliptikale Laenge
  if(x < 0)
    lo = lo+Math.PI*sgn(y);
  lo = modpi2(lo);
    
  bo = 0;                //geoz. ekliptikale Breite

  return [lo, bo]; //geoz. ekliptikale Laenge und Breite
}


/*****************************************************************************/
/* Name:    SunPos5                                                          */
/* Typ:     Funktion                                                         */
/* Zweck:   Berechnet RA und Dekl der Sonne fuer einen gegegenes Datum       */
/* Argumente:                                                                */
/*          JD: Datum                                                        */
/* Rueckgabewert:                                                            */
/*          Vektor mit RA und Dekl fuer die Epoche des Datums                */
/*****************************************************************************/
function SunPos5(JD)
{
  // Vektor mit geoz. ekl. Laenge und Breite der Sonne
  var lobo = SunEklPos(JD); 

  var lo =lobo[0];  //geoz. ekl. Laenge
  var bo =lobo[1];  //geoz. ekl. Breite

    //document.write("lo: "+lo*RtoD+"<P>");

  var T = (JD - 2451545.0) / 36525.0;     //Juli. Jhd seit 2000
  var se = (23.439291 - 0.013004*T)*DtoR; //Ekliptikschiefe

  //Deklination berechnen
  var s = Math.cos(se)*Math.sin(bo) + Math.sin(se)*Math.cos(bo)*Math.sin(lo);
  var Dekl= Math.atan(s/Math.sqrt(1-(s*s)));  

  //Rektaszension berechnen
  var c = Math.cos(lo)*Math.cos(bo);
  s = Math.sin(lo)*Math.cos(bo)*Math.cos(se) - Math.sin(bo)*Math.sin(se);

  var RA = Math.atan(s/c);
  if (c<0)
    RA = RA+Math.PI*sgn(s);
  RA = modpi2(RA);

  return [RA, Dekl]; //RA und Dekl bezogen auf Aequinokitkum des Datums
}




/*****************************************************************************/
/* Name:    EquationOfTime                                                   */
/* Typ:     Function                                                         */
/* Zweck:   Berechnet die Zeitgleichung fuer den gegebenen Zeitpunkt.        */
/* Argumente:                                                                */
/*   JD : Julianisches Datum                                                 */
/*   RA : Rektaszension der Sonne                                            */
/* Ausgabe:                                                                  */
/*   Zeitgleichung in Bogenmass                                              */
/*****************************************************************************/
function EquationOfTime(JD, RA)
{
  var T = (JD - 2451545.0) / 36525.0;
  var L0;

  var Obl = (23.439291 - 0.013004*T)*DtoR; //Ekliptikschiefe

  //Hilfsgroessen fuer Nutation
  var l = (218.316 + 481267.881*T) * DtoR; //mittlere Laenge des Mondes
  var m = (134.963 + 477198.867*T) * DtoR; //mittlere Anomalie des Mondes
  var KM= (125.045 -   1934.136*T) * DtoR; //mittlere Laenge des Mondknotens
  var Lo= (280.466 +  36000.770*T) * DtoR; //mittlere aenge der Sonne
  var Mo= (357.528 +  35999.050*T) * DtoR; //mittlere Anomalie der Sonne

  //Nutation in Laenge 
  var NutLon= (-17.200*Math.sin(KM)   + 0.206*Math.sin(2*KM)
               - 1.319*Math.sin(2*Lo) + 0.143*Math.sin(Mo)
               - 0.227*Math.sin(2*l)  + 0.071*Math.sin(m)   ) * StoR;

  T /= 10;
  L0 = T * (1/49931 - T * (1/15299 + T / 1988000));
  L0 = (280.4664567 + T * (360007.6982779 + T * (0.03032028 + L0))) * DtoR;

  return modpi(L0 - 0.0057183*DtoR - RA + NutLon*Math.cos(Obl));
}


/*****************************************************************************/
/* Name:    Tierkreiszeichen                                                 */
/* Typ:     Function                                                         */
/* Zweck:   Gibt an, in welchem Tierkreiszeichen die Sonne steht.            */
/* Argumente:                                                                */
/*   JD : Julianisches Datum                                                 */
/* Ausgabe:                                                                  */
/*    Tierkreiszeichen, 0=Widder, 1=Stier, 2=Zwillinge, 3=Krebs, 4=Loewe,    */
/*    5=Jungfrau, 6=Waage, 7=Skorpion, 8=Schuetze, 9=Steinbock,              */
/*    10=Wassermann, 11=Fische                                               */
/*****************************************************************************/
function Tierkreiszeichen(JD)
{
  var lobo = SunEklPos(JD); 
  var lo   = lobo[0];       //ekl. Laenge der Sonne

  lo = modpi2(lo);  //lo auf [0,2*pi] reduzieren
  lo *= RtoD;       //lo im Gradmass

  var TKZ = Math.floor(lo/30.0); //Tierkeiszeichen bestimmen
  
  TKZ = (TKZ >=12) ? 0 : TKZ; //Dafuer sorgen, dass TKZ=12 zu TKZ=0 wird

  return TKZ;
}


/*****************************************************************************/
/* Name:    enter_Tierkreiszeichen                                           */
/* Typ:     Function                                                         */
/* Zweck:   Wann tritt die Sonne in ein anderes Tierkreiszeichen?            */
/* Argumente:                                                                */
/*   JD : Julianisches Datum, ab dem gesucht werden soll                     */
/* Ausgabe:                                                                  */
/*    Vektor mit Nr. des Tierkreiszeichens und JD des Eintrittsdatums        */
/*****************************************************************************/
function enter_Tierkreiszeichen(JD) 
{
  var dx=2.0;               //dx = Schrittweite in Tagen   
  var lobo = new Array(2);  //Vektor, enthaelt lo und bo
    
  do
  {
    //Im welchen Tierkreiszeichen steht die Sonne jetzt?
    n = Tierkreiszeichen(JD);

    //Im welchen Tierkreiszeichen steht die Sonne beim naechsten Schritt?
    n1= Tierkreiszeichen(JD+dx);

    if(n == n1)
      JD+=dx;
    else
      dx/=2.0;  //Schrittweite halbieren
  }
  while(dx >= 0.0001);

  return [n1, JD];
}


/*****************************************************************************/
/* Name:    get_Tierkreiszeichen                                             */
/* Typ:     Function                                                         */
/* Zweck:   gibt den Namen des Tierkreiszeichens aus                         */
/* Argumente:                                                                */
/*   n: Nr des Tierkreiszeichens, wobei 0<= n <=11                           */
/* Ausgabe:                                                                  */
/*   Name des Tierkreiszeichens                                              */
/*****************************************************************************/
function get_Tierkreiszeichen(n)
{
  var Zeichen = ["Widder", "Stier", "Zwillinge", "Krebs", "L&ouml;we",
                 "Jungfrau", "Waage", "Skorpion", "Schuetze", "Steinbock",
                 "Wassermann", "Fische"];

   return Zeichen[n];
} 


/*****************************************************************************/
/* Name:    enter_constel                                                    */
/* Typ:     Function                                                         */
/* Zweck:   Wann tritt die Sonne in ein anderes Sternbild?                   */
/* Argumente:                                                                */
/*   JD : Julianisches Datum, ab dem gesucht werden soll                     */
/* Ausgabe:                                                                  */
/*    Vektor mit Nr. des Sternbilds und JD des Eintrittsdatums               */
/* Anmerkung: Unbedingt constel.js mit ins HTML-File einbauen, wenn diese    */
/*    Funktion benutzt werden soll!                                          */
/*****************************************************************************/
function enter_constel(JD)
{
  var dx=2.0;                 //dx = Schrittweite in Tagen    
  var RADekl = new Array(2);  //Vektor, enthaelt RA und Dekl
    
  do
  {
    //Im welchen Sternbild steht die Sonne jetzt?
    RADekl=SunPos5(JD);
    n = whatconstel(RADekl[0], RADekl[1], JD);

    //Im welchen Sternbild steht die Sonne beim naechsten Schritt?
    RADekl=SunPos5(JD+dx);
    n1= whatconstel(RADekl[0], RADekl[1], JD+dx);

    if(n == n1)
      JD+=dx;
    else
      dx/=2.0;  //Schrittweite halbieren
  }
  while(dx >= 0.0001);

  return [n1, JD];
}
