Sonnenstandsberechnung für die intelligente Betätigung von Rolladen, Markisen und Lampen
Update 07.04.2016: Alternatives Skript mit weniger Variablennamen
Grundlagen
Zur Sonnenstandsberechnung hatte ich bereits ein Skript erstellt, das mit mehreren linearen Interpolationsverfahren die Sonnenposition mehr oder weniger genau berechnet. http://homematic-forum.de/forum/viewtopic.php?f=31&t=17300#p140282 Theorie dazu gibt´s hier: http://de.wikipedia.org/wiki/Sonnenstand
Im Hinblick auf die weitere Verbesserung der Berechnung wird hier nun ein umfangreicheres HM-Skript beschrieben, das die Sonnenposition relativ genau für die gesamten 24h eines Tages berechnet. Durch die Verfügbarkeit brauchbarer Daten der Sonnenposition auch in der Nacht (wenn die Sonne bereits unter dem Horizont ist, Elevation negativ), können auch nachtabhängige Steuerungsaufgaben stufenlos und dämmerungsabhängig einfach gelöst werden. Zeitsprünge durch Umstellung von Sommer- auf Winterzeit entfallen durch diese Steuerungsart vollständig.
Hier das Kochrezept
– sonne_elevation mit Variablentyp Zahl von -180 bis +180 grad
– sonne_azimut mit Variablentyp Zahl von-360 bis +360 grad
!berechg sonne_elevation, sonne_azimut; stand 14.03.2014
!Verfasser: funkleuchtturm
real phi = system.Latitude(); !Breitengrad holen
phi = 0.017453292 * phi; !umwandeln in bogenmass
!####### sin_phi und cos_phi mit taylorreihe berechnen
real temp = phi * phi;
real sin_phi =phi * ((temp * temp * 0.0083334) +1.0 - (temp * 0.1666667)); !sinus-naeherung
real cos_phi = (temp *temp *0.0416667) + 1.0 - (temp * 0.5); !cosinus-naeherung
!### berechnung sonnenzeit, alle zeiten in minuten ########
integer zeit_min = system.Date("%M").ToInteger() + 60*system.Date("%H").ToInteger();
integer tagesbeginn_min = system.SunriseTime("%M").ToInteger() + 60*system.SunriseTime("%H").ToInteger();
integer tagesende_min = system.SunsetTime("%M").ToInteger() + 60* system.SunsetTime("%H").ToInteger();
integer sonnenzeit =zeit_min + 720 - 0.5 *(tagesbeginn_min +tagesende_min);
if (sonnenzeit > 1440) {sonnenzeit = sonnenzeit -1440;}
if (sonnenzeit < 1) {sonnenzeit = 1440 + sonnenzeit;}
boolean nachmittag =false;
if (sonnenzeit > 720) {sonnenzeit =sonnenzeit - 720; nachmittag = true; }
else {sonnenzeit =720 -sonnenzeit;}
!##### berechnung sin_tau und cos_tau ############
real tau = 0.00436332313 * sonnenzeit; ! 15/60 * pi /180 * sonnenzeit [0 < tau < pi ]
if (tau < 1.570796327)
{temp = tau * tau;
real sin_tau =tau * ((temp * temp * 0.0083334) +1.0 - (temp *0.1666667));
tau= 1.570796327 - tau;
temp = tau * tau;
real cos_tau =tau * ((temp * temp * 0.0083334) +1.0 - (temp * 0.1666667));}
else
{real tau1 =3.141592654 - tau;
temp = tau1 * tau1;
real sin_tau =tau1 * ((temp * temp * 0.0083334) +1.0 - (temp * 0.1666667));
tau = tau - 1.570796327;
temp = tau * tau;
real cos_tau = (tau) *(-1.0)* ((temp * temp * 0.0083334) +1.0 - (temp * 0.1666667));}
!##### berechnung delta #######################
integer tageszahl = system.Date("%j").ToInteger();
tageszahl = tageszahl +10;
if (tageszahl > 365) {tageszahl = tageszahl - 365;}
if (tageszahl < 92) {real tag = 0.0172142 *tageszahl;temp = tag * tag;
real delta = (-0.410152) *((temp *temp *0.041666) + 1.0 - (temp * 0.5));}
if ((tageszahl >91) && (tageszahl < 184)) {tageszahl = 183 - tageszahl; real tag = 0.0172142 *tageszahl;
temp = tag * tag; real delta = (0.410152) *((temp *temp *0.041666) + 1.0 - (temp * 0.5));}
if ((tageszahl >183) && (tageszahl < 275)) {tageszahl = tageszahl - 183; real tag = 0.0172142 *tageszahl;
temp = tag * tag; real delta = (0.410152) *((temp *temp *0.041666) + 1.0 - (temp * 0.5));}
if ((tageszahl >274) && (tageszahl < 366)) {tageszahl = 365 - tageszahl; real tag = 0.0172142 *tageszahl;
temp = tag * tag; real delta = (-0.410152) *((temp *temp *0.041666) + 1.0 - (temp * 0.5));}
!##### berechnung sin_delta, cos_delta #######################
temp = delta * delta;
real sin_delta =delta * ((temp * temp * 0.0083334) +1.0 - (temp * 0.1666667)); !sinus-naeherung
real cos_delta = (temp *temp *0.0416667) + 1.0 - (temp * 0.5); !cosinus-naeherung
!##### berechnung tan_delta mit stueckweiser linearisierung des tan #######################
boolean vvorzeichen = true;
if (delta< 0.0) {vvorzeichen = false; delta = (-1.0) *delta;}
real tan_delta = 1.0233 * delta;
if (delta >=0.2618) {tan_delta = (1.1822*delta) - 0.0416;}
if (vvorzeichen == false) {tan_delta = (-1.0) * tan_delta;}
!##### berechnung sin_elevation und tan_azimut #######################
real sin_elevation = (sin_phi * sin_delta) +( cos_phi * cos_delta * cos_tau);
temp = sin_elevation * sin_elevation;
real sonne_elevation = sin_elevation * (1.0 + (0.1666667 * temp) + (0.075 * temp * temp));
sonne_elevation = 57.29577951 * sonne_elevation;
real nenner = (sin_phi*cos_tau) - (cos_phi * tan_delta);
if (nenner < 0.0) {boolean plus180 = true;}
real tan_azimut = sin_tau / nenner;
!##### berechnung sonne_azimut mit stueckweiser linearisierung des arctan ############
boolean vorzeichen = true;
if (tan_azimut < 0.0) {vorzeichen = false; tan_azimut = (-1.0)*tan_azimut;}
real sonne_azimut = 0.97723 * tan_azimut;
if ((tan_azimut >=0.2679)&&(tan_azimut < 0.5774)) {sonne_azimut = (0.84588* tan_azimut) + 0.035189;}
if ((tan_azimut >= 0.5774)&&(tan_azimut < 1.0)) {sonne_azimut = (0.6195* tan_azimut) + 0.1659;}
if ((tan_azimut >= 1.0)&&(tan_azimut < 1.3032)) {sonne_azimut = (0.43173* tan_azimut) + 0.3537;}
if ((tan_azimut >= 1.3032)&&(tan_azimut < 1.7321)) {sonne_azimut = (0.3052* tan_azimut) + 0.51856;}
if ((tan_azimut >= 1.7321)&&(tan_azimut < 2.4142)) {sonne_azimut = (0.1919* tan_azimut) + 0.7148;}
if ((tan_azimut >= 2.4142)&&(tan_azimut < 2.9459)) {sonne_azimut = (0.123* tan_azimut) + 0.88115;}
if ((tan_azimut >= 2.9459)&&(tan_azimut < 3.7321)) {sonne_azimut = (0.083312* tan_azimut) + 0.9981;}
if ((tan_azimut >= 3.7321)&&(tan_azimut < 5.0)) {sonne_azimut = (0.050792* tan_azimut) + 1.1194;}
if ((tan_azimut >= 5.0)&&(tan_azimut <7.0)) {sonne_azimut = (0.02775* tan_azimut) + 1.23465;}
if ((tan_azimut >= 7.0)&&(tan_azimut <12.0)) {sonne_azimut = (0.01175117* tan_azimut) + 1.346641;}
if ((tan_azimut >= 12.0)&&(tan_azimut <20.0)) {sonne_azimut = (0.004147854* tan_azimut) + 1.437881;}
if ((tan_azimut >= 20.0)&&(tan_azimut <50.0)) {sonne_azimut = (0.0009987* tan_azimut) + 1.5008639;}
if (tan_azimut >= 50.0) {sonne_azimut = (0.000099983* tan_azimut) + 1.54579974;}
if (sonne_azimut> 1.5707963278) {sonne_azimut = 1.5707963278;}
if (vorzeichen == false) {sonne_azimut = (-1.0) * sonne_azimut;}
sonne_azimut = 57.29577951 * sonne_azimut;
if (plus180 == true) {sonne_azimut = sonne_azimut + 180.0;}
!##### tageszeitliche korrektur und werte auf systemvariablen speichern ######
if (nachmittag == false)
{sonne_azimut = 180.0 - sonne_azimut; sonnenzeit = 720 - sonnenzeit;}
else
{sonne_azimut = sonne_azimut + 180.0;sonnenzeit = 720 + sonnenzeit;}
sonne_azimut = 0.1 *((sonne_azimut*10.0) .ToInteger());
sonne_elevation = 0.1 *((sonne_elevation*10.0) .ToInteger());
dom.GetObject("sonne_elevation").State(sonne_elevation);
dom.GetObject("sonne_azimut").State(sonne_azimut);Ein sehr einfaches beispielhaftes WebUI-Programm zum Betätigen einer Rolllade ist hier:
Mit der Systemvariablen roll_AutoManuell wird die automatische Steuerung ein- und ausgeschaltet. Diese Systemvariable kann mit einem weiteren WebUI-Programm und z.B. einer Fernbedienung ein- und ausgeschaltet werden. Die Möglichkeiten, mit Elevation und Azimut die Rollladen und Lampen eines Haus intelligent nach Licht und Schatten zu steuern, sind sehr vielfältig. Ich verweise hier auf mein Tutorial zur Rollladensteuerung.
Ach ja, zum Vergleich kann man die Daten mit einem „genauen“ Sonnenstandsrechner vergleichen: http://www.volker-quaschning.de/datserv/sunpos/index.php
Breiten und Längengrad eingeben und Berechnungsart „SUNAE“ wählen. Die Abweichungen sind bei mir meist unter einem Grad, was bei den verwendeten Näherungsverfahren für die trigonometrischen Funktionen ganz ordentlich ist. Für die Heimautomation ist das jedenfalls mehr als genau!
Ergänzung 14.03.2014
Wie man an den verwendeten Zahlenkolonnen im Skript erkennt, wurde bewusst weniger auf programmiertechnische „Schönheit“ ausgelegt, sondern vielmehr wurde auf möglichst wenig Rechenoperationen und damit Laufzeit optimiert. Dies ist deshalb wichtig, weil das Skript zyklisch häufig aufgerufen wird!
Geringe Unstetigkeiten im Verlauf des Elevations- und Azimutwinkels bedingt durch die verwendeten Näherungsverfahren wurden weiter optimiert, so dass jetzt für die Heimautomation mehr als ausreichende Genauigkeit der Sonnenkoordinaten erreicht wird. Insbesondere im Bereich des Nulldurchgangs der Elevation ist ein möglichst stetiger Winkelverlauf deshalb wichtig, weil die meisten Steuerungsaufgaben im Bereich der Dämmerung erfolgen und hier ein möglichst gleichmässiger stufenloser Verlauf von grossem Vorteil ist.
Das folgende Bild zeigt sehr schön den Azimut und die Elevation, aufgezeichnet mit Historian über 24h. Diese Kurve verschieben sich täglich um einen geringen Betrag, wie man an dem geringen Sprung in der Elevationskurve um ca 00:00h erkennt:
Genauigkeit
Der Vergleich mit Winkelergebnissen der Internet-Sonnenpositionsrechner ist abhängig von dem dort verwendeten Berechnungsverfahren:
Bei dem Berechnungsverfahren SUNAE war der Unterschied bei der Elevation meist weniger als ein halbes Grad. Der Azimutwinkel hatte Abweichungen bis zu etwa 2 Grad. Letztlich ist aber nicht die Genauigkeit entscheidend, sondern dass man im Bereich des Sonnenauf- und -untergangs mit dem Elevationswinkel eine gur aufgelöste stufenlose Zahl hat, von der man Schaltvorgänge reproduzierbar abhängig machen kann. Mit dem Azimut steuert man in der Regel sowieso nur Beschattungseinrichtungen, bei denen eine Winkelgenauigkeit von 2Grad mehr als genau sind.
Update 07.04.2016:
Alternatives Skript mit weniger verwendeter Systemvariablen
Vermutlich kann die CCU nur eine begrenzte Zahl von Namen in Skripten stabil verarbeiten. (man spricht von 200 ?!). Deshalb wurde das Sonnenstand-Skript diesbezüglich weiter optimiert, so daß beim nachfolgenden verbesserten Skript nur noch 11 Variablen verwendet werden. Allerdings ist die Lesbarkeit das Verständnis des Programmes dabei verloren gegangen. Aber man kann sich ja dafür das erste Listing ansehen, welches die gleiche Funktion allerdings mit deutlich mehr Variablennamen hat.
!berechnung sonne_elevation, sonne_azimut; stand 07.08.2016 verfasser: eugen stall, stall.biz
!Zwei Systemvariablen anlegen:
!sonne_elevation mit Variablentyp Zahl von -180 bis +180 grad
!sonne_azimut mit Variablentyp Zahl von-360 bis +360 grad
!dieses skript alle 4min aufrufen
real temp;
real sin_phi;
real cos_phi;
integer sonnenzeit;
boolean nachmittag;
real sin_tau;
real cos_tau;
real tau1;
real sin_delta;
real cos_delta;
boolean vorzeichen;
integer temp = system.Date("%M").ToInteger() + 60*system.Date("%H").ToInteger();
integer sonnenzeit =temp + 720 - 0.5 *((system.SunriseTime("%M").ToInteger() + 60*system.SunriseTime("%H").ToInteger()) +system.SunsetTime("%M").ToInteger() + 60* system.SunsetTime("%H").ToInteger());
if (sonnenzeit > 1440) {sonnenzeit = sonnenzeit -1440;}
if (sonnenzeit < 1) {sonnenzeit = 1440 + sonnenzeit;}
if (sonnenzeit > 720) {sonnenzeit =sonnenzeit - 720; nachmittag = true; }
else {sonnenzeit =720 -sonnenzeit;nachmittag =false;}
temp = 0.00436332313 * sonnenzeit;
if (temp < 1.570796327)
{real sin_tau =temp * ((temp * temp * temp * temp * 0.0083334) +1.0 - (temp * temp *0.1666667));
temp= 1.570796327 - temp;
cos_tau =temp * ((temp * temp * temp * temp * 0.0083334) +1.0 - (temp * temp * 0.1666667));}
else {real tau1 =3.141592654 - temp;
sin_tau =tau1 * ((tau1 * tau1 * tau1 * tau1 * 0.0083334) +1.0 - (tau1 * tau1 * 0.1666667));
temp = temp - 1.570796327;
real cos_tau = (temp) *(-1.0)* ((temp * temp * temp * temp * 0.0083334) +1.0 - (temp * temp * 0.1666667));}
temp = system.Date("%j").ToInteger() +10;
if (temp > 365) {temp = temp - 365;}
if (temp < 92) {real tag = 0.0172142 *temp;
tau1 = (-0.410152) *((tag * tag *tag * tag *0.041666) + 1.0 - (tag * tag * 0.5));}
if ((temp >91) && (temp < 184)) {temp = 183 - temp; real tag = 0.0172142 *temp;
tau1 = (0.410152) *((tag * tag *tag * tag *0.041666) + 1.0 - (tag * tag * 0.5));}
if ((temp >183) && (temp < 275)) {temp = temp - 183; real tag = 0.0172142 *temp;
tau1 = (0.410152) *((tag * tag *tag * tag *0.041666) + 1.0 - (tag * tag * 0.5));}
if ((temp >274) && (temp < 366)) {temp = 365 - temp; real tag = 0.0172142 *temp;
tau1 = (-0.410152) *((tag * tag *tag * tag *0.041666) + 1.0 - (tag * tag * 0.5));}
sin_delta =tau1 * ((tau1 * tau1 * tau1 * tau1 * 0.0083334) +1.0 - (tau1 * tau1 * 0.1666667)); !sinus-naeherung
cos_delta = (tau1 * tau1 * tau1 * tau1 *0.0416667) + 1.0 - (tau1 * tau1 * 0.5); !cosinus-naeherung
if (tau1 < 0.0) {vorzeichen = false; tau1 = (-1.0) * tau1;} else {vorzeichen = true;}
temp = 0.017453292*system.Latitude();
sin_phi =(temp * ((temp * temp * temp * temp * 0.0083334) +1.0 - (temp * temp * 0.1666667)));
cos_phi = ((temp * temp *temp * temp *0.0416667) + 1.0 - (temp * temp * 0.5));
temp = (sin_phi * sin_delta) +( cos_phi * cos_delta * cos_tau);
temp = temp * (1.0 + (0.1666667 * temp * temp) + (0.075 * temp * temp * temp * temp));
temp = 57.29577951 * temp;
temp= 0.1 *((temp *10.0) .ToInteger());
dom.GetObject("sonne_elevation").State(temp);
WriteLine(temp);
if (tau1 >=0.2618) {tau1 = (1.1822 * tau1) - 0.0416;}else {tau1 = 1.0233 * tau1;}
if (vorzeichen == false) {tau1 = (-1.0) * tau1;}
temp = (sin_phi*cos_tau) - (cos_phi * tau1);
if (temp < 0.0) {boolean plus180 = true;}
temp = sin_tau / temp;
vorzeichen = true;
if (temp < 0.0) {vorzeichen = false; temp = (-1.0)*temp;}
tau1 = 0.97723 * temp;
if ((temp >= 0.2679)&&(temp < 0.5774)) {tau1 = (0.84588* temp) + 0.035189;}
if ((temp >= 0.5774)&&(temp < 1.0)) {tau1 = (0.6195* temp) + 0.1659;}
if ((temp >= 1.0)&&(temp < 1.3032)) {tau1 = (0.43173* temp) + 0.3537;}
if ((temp >= 1.3032)&&(temp < 1.7321)) {tau1 = (0.3052* temp) + 0.51856;}
if ((temp >= 1.7321)&&(temp < 2.4142)) {tau1 = (0.1919* temp) + 0.7148;}
if ((temp >= 2.4142)&&(temp < 2.9459)) {tau1 = (0.123* temp) + 0.88115;}
if ((temp >= 2.9459)&&(temp < 3.7321)) {tau1 = (0.083312* temp) + 0.9981;}
if ((temp >= 3.7321)&&(temp < 5.0)) {tau1 = (0.050792* temp) + 1.1194;}
if ((temp >= 5.0)&&(temp <7.0)) {tau1 = (0.02775* temp) + 1.23465;}
if ((temp >= 7.0)&&(temp <12.0)) {tau1 = (0.01175117*temp) + 1.346641;}
if ((temp >= 12.0)&&(temp <20.0)) {tau1 = (0.004147854* temp) + 1.437881;}
if ((temp >= 20.0)&&(temp <50.0)) {tau1 = (0.0009987*temp) + 1.5008639;}
if (temp >= 50.0) {tau1 = (0.000099983* temp) + 1.54579974;}
if (tau1> 1.5707963278) {tau1 = 1.5707963278;}
if (vorzeichen == false) {tau1 = (-1.0) * tau1;}
tau1 = 57.29577951 * tau1;
if (plus180 == true) {tau1 = tau1 + 180.0;}
if (nachmittag == false)
{tau1 = 180.0 - tau1; sonnenzeit = 720 - sonnenzeit;}
else {tau1 = tau1 + 180.0;sonnenzeit = 720 + sonnenzeit;}
tau1 = 0.1 *((tau1*10.0) .ToInteger());
dom.GetObject("sonne_azimut").State(tau1);
WriteLine(tau1);
