Changeset 14306 for trunk/FACT++/src/smartfact.cc

Timestamp:

08/02/12 13:22:38 (12 years ago)

Author:

tbretz

Message:

Fixed a performance issue with the calculation of the moon position, convert the light condition directly to currents, renamed the light-condition to current-prediction.

File:

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

trunk/FACT++/src/smartfact.cc

@@ -1855 +1855 @@
 
 #ifdef HAVE_NOVA
+
+    vector<pair<ln_equ_posn, double>> fMoonCoords;
+
+    void CalcMoonCoords(double jd)
+    {
+        jd = floor(jd);
+
+        fMoonCoords.clear();
+        for (double h=0; h<1; h+=1./(24*12))
+        {
+            ln_equ_posn  moon;
+            ln_get_lunar_equ_coords_prec(jd+h, &moon, 0.01);
+
+            const double disk = ln_get_lunar_disk(jd+h);
+
+            fMoonCoords.push_back(make_pair(moon, disk));
+        }
+    }
+
     pair<vector<float>, pair<Time, float>> GetVisibility(ln_equ_posn *src, ln_lnlat_posn *observer, double jd)
     {
@@ -1870 +1889 @@
         int cnt = 0;
 
+        int i=0;
+
         vector<float> alt;
-        for (double h=0; h<1; h+=1./(24*12))
+        for (double h=0; h<1; h+=1./(24*12), i++)
         {
             if (src==0)
-                ln_get_lunar_equ_coords(jd+h, &moon);
+                moon = fMoonCoords[i].first; //ln_get_lunar_equ_coords_prec(jd+h, &moon, 0.01);
 
             ln_hrz_posn hrz;
@@ -1905 +1926 @@
         const double jd1 = fmod(fSun.fRiseAstronomical.JD(), 1);
 
-        ln_equ_posn  moon;
-
         double max   = -1;
         double maxjd =  0;
 
         int cnt = 0;
-
-        vector<float> alt;
-        for (double h=0; h<1; h+=1./(24*12))
-        {
-            double lc = -1;
+        int i = 0;
+
+        vector<float> vec;
+        for (double h=0; h<1; h+=1./(24*12), i++)
+        {
+            double cur = -1;
 
             if (h>jd0 && h<jd1)
             {
-                const double disk = ln_get_lunar_disk(jd+h);
-                ln_get_lunar_equ_coords(jd+h, &moon);
+                ln_equ_posn moon  = fMoonCoords[i].first;//ln_get_lunar_equ_coords_prec(jd+h, &moon, 0.01);
+                const double disk = fMoonCoords[i].second;//ln_get_lunar_disk(jd+h);
 
                 ln_hrz_posn hrz;
@@ -1931 +1951 @@
                 const double angle = m.Angle(src->ra, src->dec);
 
-                lc = angle*hrz.alt*disk*M_PI*M_PI/180/180;
-
-                alt.push_back(lc);
+                const double lc = angle*hrz.alt*pow(disk, 6)/360/360;
+
+                cur = 7.7+4942*lc;
+
+                vec.push_back(cur); // Covert LC to pixel current in uA
             }
 
-            if (lc>max)
+            if (cur>max)
             {
-                max   = lc;
+                max   = cur;
                 maxjd = jd+h;
             }
 
-            if (h>jd0 && h<jd1 && lc>0)
+            if (h>jd0 && h<jd1 && cur>0)
                 cnt++;
         }
 
-        if (max<=15 || cnt==0)
+        if (max<0 || cnt==0)
             return make_pair(vector<float>(), make_pair(Time(), 0));
 
-        return make_pair(alt, make_pair(maxjd, maxjd>jd+jd0&&maxjd<jd+jd1?max:0));
+        return make_pair(vec, make_pair(maxjd, maxjd>jd+jd0&&maxjd<jd+jd1?max:0));
     }
 #endif
@@ -1959 +1981 @@
 
         Time now;
+
+        CalcMoonCoords(now.JD());
 
         fSun  = Sun (lon, lat, now);
@@ -2020 +2044 @@
         map<Time, pair<string, float>> lightcond;
         vector<vector<float>> alt;
-        vector<vector<float>> hlc;
+        vector<vector<float>> cur;
 
 #ifdef HAVE_NOVA
@@ -2069 +2093 @@
                 if (lc.first.size()>0)
                 {
-                    hlc.push_back(lc.first);
+                    cur.push_back(lc.first);
                     lightcond[lc.second.first] = make_pair(name, lc.second.second);
                 }
@@ -2111 +2135 @@
             }
 
+            out4 << setprecision(3);
             for (auto it=lightcond.begin(); it!=lightcond.end(); it++)
             {
@@ -2117 +2142 @@
                 out4 << "<B>" << it->second.first << "</B>";
                 if (it->second.second>0)
-                    out4 << " [" << nearbyint(90-it->second.second) << "&deg;]";
+                    out4 << " [" << nearbyint(it->second.second) << "]";
             }
 
@@ -2134 +2159 @@
             out2 << HTML::kWhite << '\t' << Time()-now << '\n';
 
-            WriteBinaryVec(now, "hist-visibility",      alt, 75, 15, "Alt "+title.str());
-            WriteBinaryVec(now, "hist-light-condition", hlc,  1,  0, "LC "+title.str());
+            WriteBinaryVec(now, "hist-visibility",         alt,  75, 15, "Alt "+title.str());
+            WriteBinaryVec(now, "hist-current-prediction", cur, 100,  0, "I "  +title.str());
         }
         catch (const exception &e)
@@ -2151 +2176 @@
         ofstream(fPath+"/source-list.data") << out2.str();
         ofstream(fPath+"/visibility.data") << out3.str();
-        ofstream(fPath+"/light-condition.data") << out4.str();
+        ofstream(fPath+"/current-prediction.data") << out4.str();
     }