FUNCTION_BLOCK SUN_POS
TITLE = 'SUN_POS'
//
//this FUNCTION block calculates the sun position for a given date and time.
//the times are calculated in utc and have to be corrected for the given time zone.
//B is the angle from north and HR is the highth in degrees.
//
//uses: DT_DATE [FC6]
// DT_TOD [FC8]
// MODR [FC97]
// FLOOR2 [FC82]
// RAD [FC100]
// DEG [FC72]
// REFRACTION [FC206]
//
VERSION : '2.1'
AUTHOR : hugo
NAME : SUNPOS
FAMILY : TIDT
VAR_INPUT
latitude : REAL; (* latitude of geographical position *)
longitude : REAL; (* longitude of geographical position *)
utc : DT; (* world time *)
END_VAR
VAR_OUTPUT
B: REAL;
H: REAL;
HR: REAL;
END_VAR
VAR
g: REAL;
a: REAL;
d: REAL;
#t1: REAL;
n: REAL;
e: REAL;
c: REAL;
tau: REAL;
sin_d: REAL;
rlat: REAL;
sin_lat: REAL;
cos_lat: REAL;
cos_tau: REAL;
cos_d: REAL;
END_VAR
VAR_TEMP
temp_utc_DATE : DATE;
temp_utc_TOD : TOD;
temp : DINT;
END_VAR
CONST
PI := 3.14159265358979323846264338327950288; (* Kreiszahl PI *)
PI2 := 6.28318530717958647692528676655900576; (* PI * 2 *)
END_CONST
BEGIN
(* n is the julian date and the number of days since 1.1.2000-12:00 midday *)
(* be careful for step7 date startes 1.1.1990 instead of 1.1.1970 *)
(* for step7 this line must change *)
temp_utc_DATE := DT_DATE(IN := UTC); // DATE
temp_utc_TOD := DT_TOD(IN := UTC); // TIME_OF_DAY
temp := (DATE_TO_DINT(temp_utc_DATE) + 7305) * 86400;
temp := temp + TOD_TO_DINT(temp_utc_TOD) / 1000;
(* n is the julian date and the number of days since 1.1.2000-12:00 midday *)
(* be careful for step7 date startes 1.1.1990 instead of 1.1.1970 *)
(* for step7 this line must change *)
n := (temp - 946728000) * 0.000011574074074074;
g := MODR(IN:=6.240040768 + 0.01720197 * n, DIVI:=PI2);
d := MODR(IN:=4.89495042 + 0.017202792 * n, DIVI:=PI2) + 0.033423055 * SIN(g) + 0.000349066 * SIN(2.0*g);
e := 0.409087723 - 0.000000006981317008 * n;
cos_d := COS(d);
sin_d := SIN(d);
a := ATAN(COS(e) * sin_d / cos_d);
IF cos_d < 0.0 THEN a := a + PI; END_IF;
c := ASIN(SIN(e) * sin_d);
(* also here we must be very careful utc is from 1.1.1970 for step7 the formula must change *)
tau := RAD(MODR(IN:=6.697376 + (n - 0.25) * 0.0657098245037645 + DINT_TO_REAL(TOD_TO_DINT(DT_TOD(utc))) * 0.0000002785383333, DIVI:=24.0) * 15.0 + longitude) - a;
rlat := RAD(latitude);
sin_lat := SIN(rlat);
cos_lat := COS(rlat);
cos_tau := COS(tau);
#t1 := cos_tau * sin_lat - TAN(c) * cos_lat;
B := ATAN(SIN(tau) / #t1);
IF #t1< 0.0 THEN B := B + PI2; ELSE B := B + PI; END_IF;
B := DEG(MODR(IN:=b, DIVI:=PI2));
h := DEG(ASIN(COS(C) * cos_tau * cos_lat +SIN(c) * sin_lat));
IF h > 180.0 THEN h := h - 360.0; END_IF;
(* consider refraction *)
HR := h + REFRACTION(h);
END_FUNCTION_BLOCK