Index: trunk/FACT++/src/moon.cc
===================================================================
--- trunk/FACT++/src/moon.cc	(revision 14183)
+++ trunk/FACT++/src/moon.cc	(revision 14183)
@@ -0,0 +1,202 @@
+#include <libnova/solar.h>
+#include <libnova/lunar.h>
+#include <libnova/rise_set.h>
+#include <libnova/transform.h>
+
+#include "Time.h"
+#include "Configuration.h"
+
+using namespace std;
+
+// ------------------------------------------------------------------------
+
+class Moon
+{
+public:
+    Time time;
+
+    double ra;
+    double dec;
+
+    double zd;
+    double az;
+
+    double disk;
+
+    bool visible;
+
+    Time fRise;
+    Time fTransit;
+    Time fSet;
+
+    int state;
+
+    Moon() : time(Time::none)
+    {
+    }
+
+    // Could be done more efficient: Only recalcuate if
+    // the current time exceeds at least on of the stored times
+    Moon(double lon, double lat, const Time &t=Time()) : time(t)
+    {
+        const double JD = time.JD();
+
+        ln_lnlat_posn observer;
+        observer.lng = lon;
+        observer.lat = lat;
+
+        //observer.lng.degrees = -5;
+        //observer.lng.minutes = 36;
+        //observer.lng.seconds = 30;
+        //observer.lat.degrees = 42;
+        //observer.lat.minutes = 35;
+        //observer.lat.seconds = 40;
+
+        ln_rst_time moon;
+        ln_get_lunar_rst(JD-0.5, &observer, &moon);
+
+        fRise    = Time(moon.rise);
+        fTransit = Time(moon.transit);
+        fSet     = Time(moon.set);
+
+        //visible =
+        //    ((JD>moon.rise && JD<moon.set ) && moon.rise<moon.set) ||
+        //    ((JD<moon.set  || JD>moon.rise) && moon.rise>moon.set);
+
+        const bool is_up      = JD>moon.rise;
+        const bool is_sinking = JD>moon.transit;
+        const bool is_dn      = JD>moon.set;
+
+        ln_get_lunar_rst(JD+0.5, &observer, &moon);
+        if (is_up)
+            fRise = Time(moon.rise);
+        if (is_sinking)
+            fTransit = Time(moon.transit);
+        if (is_dn)
+            fSet = Time(moon.set);
+
+        ln_equ_posn pos;
+        ln_get_lunar_equ_coords(JD, &pos);
+
+        ln_hrz_posn hrz;
+        ln_get_hrz_from_equ (&pos, &observer, JD, &hrz);
+        az =    hrz.az;
+        zd = 90-hrz.alt;
+
+        ra  = pos.ra/15;
+        dec = pos.dec;
+
+        disk = ln_get_lunar_disk(JD)*100;
+        state = 0;
+        if (fRise   <fTransit && fRise   <fSet)     state = 0;  // not visible
+        if (fTransit<fSet     && fTransit<fRise)    state = 1;  // before culm
+        if (fSet    <fRise    && fSet    <fTransit) state = 2;  // after culm
+
+        visible = state!=0;
+
+        // 0: not visible
+        // 1: visible before cul
+        // 2: visible after cul
+    }
+
+    double Angle(double r, double d) const
+    {
+        const double theta0 = M_PI/2-d*M_PI/180;
+        const double phi0   = r*M_PI/12;
+
+        const double theta1 = M_PI/2-dec*M_PI/180;
+        const double phi1   = ra*M_PI/12;
+
+        const double x0 = sin(theta0) * cos(phi0);
+        const double y0 = sin(theta0) * sin(phi0);
+        const double z0 = cos(theta0);
+
+        const double x1 = sin(theta1) * cos(phi1);
+        const double y1 = sin(theta1) * sin(phi1);
+        const double z1 = cos(theta1);
+
+        double arg = x0*x1 + y0*y1 + z0*z1;
+        if(arg >  1.0) arg =  1.0;
+        if(arg < -1.0) arg = -1.0;
+
+        return acos(arg) * 180/M_PI;
+    }
+};
+
+// ========================================================================
+// ========================================================================
+// ========================================================================
+
+void SetupConfiguration(Configuration &conf)
+{
+    po::options_description control("Smart FACT");
+    control.add_options()
+        ("ra",        var<double>(), "Source right ascension")
+        ("dec",       var<double>(), "Source declination")
+        ("date-time", var<string>()
+#if BOOST_VERSION >= 104200
+         ->required()
+#endif                             
+                                   , "SQL time (UTC)")
+        ;
+
+    po::positional_options_description p;
+    p.add("date-time", 1); // The first positional options
+
+    conf.AddOptions(control);
+    conf.SetArgumentPositions(p);
+}
+
+void PrintUsage()
+{
+    cout <<
+        "moon - The moon calculator\n"
+        "\n"
+        "Usage: moon sql-datetime [--ra={ra} --dec={dec}]\n";
+    cout << endl;
+}
+
+int main(int argc, const char* argv[])
+{
+    Configuration conf(argv[0]);
+    conf.SetPrintUsage(PrintUsage);
+    SetupConfiguration(conf);
+
+    if (!conf.DoParse(argc, argv))
+        return 127;
+
+    if (conf.Has("ra")^conf.Has("dec"))
+    {
+        cout << "ERROR - Either --ra or --dec missing." << endl;
+        return 1;
+    }
+
+    Time time;
+    time.SetFromStr(conf.Get<string>("date-time"));
+
+    ln_lnlat_posn observer;
+    observer.lng = -(17.+53./60+26.525/3600);
+    observer.lat =   28.+45./60+42.462/3600; 
+
+    const Moon moon(observer.lng, observer.lat, time);
+
+    ln_equ_posn pos;
+    pos.ra  = conf.Get<double>("ra")*15;
+    pos.dec = conf.Get<double>("dec");
+
+    ln_hrz_posn hrz;
+    ln_get_hrz_from_equ(&pos, &observer, time.JD(), &hrz);
+
+    cout << setprecision(15);
+    cout << time.GetAsStr() << '\n';
+    cout << moon.visible    << '\n';
+    cout << moon.disk       << '\n';
+    cout << moon.zd         << '\n';
+
+    if (conf.Has("ra") && conf.Has("dec"))
+        cout << moon.Angle(pos.ra/15, pos.dec) << '\n';
+
+    cout << endl;
+
+    return 0;
+}
