Changeset 13778 for trunk

Timestamp:

05/18/12 17:45:26 (12 years ago)

Author:

tbretz

Message:

First draft of an implementation of sun-set, sun-rise, etc.

File:

• trunk/FACT++/src/smartfact.cc (modified) (3 diffs)

trunk/FACT++/src/smartfact.cc

@@ +1 @@
+#ifdef HAVE_LIBNOVA
+#include <libnova/solar.h>
+#include <libnova/lunar.h>
+#include <libnova/rise_set.h>
+#endif
+
 #include "Dim.h"
 #include "Event.h"
@@ -24 +30 @@
 
 #include "DimDescriptionService.h"
+
+// ------------------------------------------------------------------------
+
+#ifdef HAVE_LIBNOVA
+
+class Astro
+{
+    Time fSunRiseDayTime;
+    Time fSunRiseCivil;
+    Time fSunRiseAstronomical;
+    Time fSunRiseDarkTime;
+
+    Time fSunSetDayTime;
+    Time fSunSetCivil;
+    Time fSunSetAstronomical;
+    Time fSunSetDarkTime;
+
+    ln_rst_time moon;
+
+    ln_equ_posn moon_position;
+    double      moon_disk;
+
+    string description;
+
+    //observer.lng = -(17.+53./60+26.525/3600);
+    //observer.lat =   28.+45./60+42.462/3600;
+public:
+    Astro(double lon, double lat, /*const*/ Time/* &*/time=Time())
+    {
+        // observers location (Edinburgh), used to calc rst
+        struct ln_lnlat_posn observer;
+        observer.lng = -(17.+53./60+26.525/3600);
+        observer.lat =   28.+45./60+42.462/3600;
+
+        time = Time(2012, 5, 18, 23, 30);
+
+        // get Julian day from local time
+        const double JD = time.JD();
+
+        ln_rst_time sun_day;
+        ln_rst_time sun_civil;
+        ln_rst_time sun_astronomical;
+        ln_rst_time sun_dark;
+
+        // Warning: return code of 1 means circumpolar and is not checked!
+        ln_get_lunar_rst        (JD-0.5, &observer,      &moon);
+        ln_get_solar_rst        (JD-0.5, &observer,      &sun_day);
+        ln_get_solar_rst_horizon(JD-0.5, &observer, - 6, &sun_civil);
+        ln_get_solar_rst_horizon(JD-0.5, &observer, -12, &sun_astronomical);
+        ln_get_solar_rst_horizon(JD-0.5, &observer, -18, &sun_dark);
+
+        fSunSetDayTime       = Time(sun_day.set);
+        fSunSetCivil         = Time(sun_civil.set);
+        fSunSetAstronomical  = Time(sun_astronomical.set);
+        fSunSetDarkTime      = Time(sun_dark.set);
+
+        fSunRiseDayTime      = Time(sun_day.rise);
+        fSunRiseCivil        = Time(sun_civil.rise);
+        fSunRiseAstronomical = Time(sun_astronomical.rise);
+        fSunRiseDarkTime     = Time(sun_dark.rise);
+
+        const bool is_day   = JD>sun_day.rise;
+        const bool is_night = JD>sun_dark.set;
+
+        ln_get_lunar_rst        (JD+0.5, &observer,      &moon);
+        ln_get_solar_rst        (JD+0.5, &observer,      &sun_day);
+        ln_get_solar_rst_horizon(JD+0.5, &observer, - 6, &sun_civil);
+        ln_get_solar_rst_horizon(JD+0.5, &observer, -12, &sun_astronomical);
+        ln_get_solar_rst_horizon(JD+0.5, &observer, -18, &sun_dark);
+
+        if (is_day)
+        {
+            fSunRiseDayTime      = Time(sun_day.rise);
+            fSunRiseCivil        = Time(sun_civil.rise);
+            fSunRiseAstronomical = Time(sun_astronomical.rise);
+            fSunRiseDarkTime     = Time(sun_dark.rise);
+        }
+
+        if (is_night)
+        {
+            fSunSetDayTime      = Time(sun_day.set);
+            fSunSetCivil        = Time(sun_civil.set);
+            fSunSetAstronomical = Time(sun_astronomical.set);
+            fSunSetDarkTime     = Time(sun_dark.set);
+        }
+
+        // case 0: midnight to sun-rise | !is_day && !is_night | rise/set
+        // case 1: sun-rise to sun-set  |  is_day && !is_night | set /rise
+        // case 2: sun-set  to midnight |  is_day &&  is_night | rise/set
+
+        if (is_day^is_night)
+        {
+            cout << "SunSet:  " << fSunSetDayTime  << endl;
+            cout << "SunRise: " << fSunRiseDayTime << endl;
+        }
+        else
+        {
+            cout << "SunRise: " << fSunRiseDayTime << endl;
+            cout << "SunSet:  " << fSunSetDayTime  << endl;
+        }
+
+        //const double mjd = Time(JD-2400000.5).Mjd();
+
+        int state = is_day^is_night ? 4 : 0;
+        if (time>fSunSetDayTime)       state++;
+        if (time>fSunSetCivil)         state++;
+        if (time>fSunSetAstronomical)  state++;
+        if (time>fSunSetDarkTime)      state++;
+
+        if (time>fSunRiseDarkTime)     state++;
+        if (time>fSunRiseAstronomical) state++;
+        if (time>fSunRiseCivil)        state++;
+        if (time>fSunRiseDayTime)      state++;
+
+        string name[] =
+        {
+            "dark time",
+            "astronomical twilight",
+            "civil twilight",
+            "sunrise",
+            "day time",
+            "sunset",
+            "civil twilight",
+            "astronomical twilight",
+            "dark time"
+        };
+
+        description = state[name];
+
+
+        cout << "NOW: " << name[state] << " [" << state << "]" << endl << endl;
+
+/*
+        // RA, DEC
+        ln_get_lunar_equ_coords(JD, &moon_position);
+        //cout << "MoonRa:  " << equ.ra  << "\n";
+        //cout << "MoonDec: " << equ.dec << "\n";
+
+        // lunar disk, phase and bright limb
+        moon_disk = ln_get_lunar_disk(JD);
+
+        if (JD>sun_dark.set || JD<sun_normal.rise)
+        {
+            cout << "Sunrise will be at:" << endl;
+            cout << "Astronom. TwL: " << Time(sun_dark.rise        -2400000.5).GetAsStr() << endl;
+            cout << "Civil     TwL: " << Time(sun_astronomical.rise-2400000.5).GetAsStr() << endl;
+            cout << "Sunrise:       " << Time(sun_civil.rise       -2400000.5).GetAsStr() << endl;
+            cout << "Day time:      " << Time(sun_normal.rise      -2400000.5).GetAsStr() << endl;
+        }
+        else
+        {
+            cout << "Sunset will be at:" << endl;
+            cout << "Sunset:        " << Time(sun_normal.set      -2400000.5).GetAsStr() << endl;
+            cout << "Civil     TwL: " << Time(sun_civil.set       -2400000.5).GetAsStr() << endl;
+            cout << "Astronom. TwL: " << Time(sun_astronomical.set-2400000.5).GetAsStr() << endl;
+            cout << "Dark time:     " << Time(sun_dark.set        -2400000.5).GetAsStr() << endl;
+        }
+
+        if (JD>moon.rise && JD<moon.set)
+            cout << "Moon: VISIBLE (" << setprecision(2) <<  moon_disk*100 << "%)" << endl;
+
+        if (JD<moon.rise)
+            cout << "MoonRise: " << Time(moon.rise   -2400000.5) << endl;
+        if (JD<moon.transit)
+            cout << "MoonCulm: " << Time(moon.transit-2400000.5) << endl;
+        if (JD<moon.set)
+            cout << "MoonSet:  " << Time(moon.set    -2400000.5) << endl;
+            */
+        }
+
+    /*
+
+    Always show the next event:
+
+    What do we have now?
+    time between sunset and civil:        sunset
+    time between civial and astronominal: civil
+    time between astronominal and dark:   astronomical
+    time above dark:                      dark
+
+
+    */
+
+};
+
+#endif
 
 // ------------------------------------------------------------------------
@@ -408 +600 @@
         out << "#ffffff\t" << dir[idx] << '\n';
 
+        //Astro a(-(17.+53./60+26.525/3600), 28.+45./60+42.462/3600);
 
         ofstream(fPath+"/magicweather.txt") << out.str();