Index: trunk/Mars/fact/processing/currents.C
===================================================================
--- trunk/Mars/fact/processing/currents.C	(revision 17075)
+++ trunk/Mars/fact/processing/currents.C	(revision 17100)
@@ -28,5 +28,4 @@
     double dev_avg  = 0;
     double dev_rms  = 0;
-    int    cnt  = 0;
     int    cnt_beg  = 0;
     int    cnt_end  = 0;
@@ -37,4 +36,7 @@
 
     double calc_dev_avg  = 0;
+    double calc_dev_rel  = 0;
+
+    TGraph g;
 
     while (file.GetNextRow())
@@ -76,21 +78,48 @@
 
         //calculate the light condition
-        double calt = sin(hrzm.alt*TMath::DegToRad());
         double disk = Nova::GetLunarDisk(jd);
-        double lc = calt*pow(disk, 2.5);
 
-        double Ical = lc>0 ? 4.5+103.0*lc : 4.5;
+        // Current prediction
+        //double cang = sin(angle);
+        const double calt = sin(hrzm.alt*TMath::DegToRad());
+
+        //                                                                           semi-major axis
+        const double lc = /*cang==0 ? 0 :*/ calt*pow(disk, 2.2)*pow(Nova::GetLunarEarthDist(jd)/384400, -2);///sqrt(cang);
+        //const double lc = cang==0 ? -1 : calt*pow(disk, 2.2)*pos(dist, -2);
+        const double Ical = lc>0 ? 4+103*lc : 4;
+
+        g.SetPoint(g.GetN(), time, med);
+
         calc_dev_avg += med - Ical;
+        calc_dev_rel += med/Ical;
 
-        cnt++;
     }
 
-    if (cnt==0)
+    if (g.GetN()==0)
     {
         if (diff<5./24/3600)
             return;
 
+        // return the last report before the run
+        //   (if not older than 5 minutes)
         cout << "result " << med_last << " 0 " << dev_last << " 0 0 0 0" << endl;
         return;
+    }
+
+    Double_t a0, a1;
+    Int_t ifail;
+    g.LeastSquareLinearFit(0, a0, a1, ifail);
+
+    double fluct_abs = 0;
+    double fluct_rel = 0;
+    if (ifail==0)
+    {
+        for (int i=0; i<g.GetN(); i++)
+            g.GetY()[i] -= a0+a1*g.GetX()[i];
+        fluct_abs = g.GetRMS(2);
+
+        for (int i=0; i<g.GetN(); i++)
+            g.GetY()[i] /= a0+a1*g.GetX()[i];
+        fluct_rel = g.GetRMS(2);
     }
 
@@ -98,13 +127,17 @@
     med_end = cnt_end>0 ? med_end/cnt_end : -1;
 
-    med_avg /= cnt;
-    med_rms /= cnt;
+    med_avg /= g.GetN();
+    med_rms /= g.GetN();
     med_rms = sqrt(med_rms-med_avg*med_avg);
 
-    dev_avg /= cnt;
-    dev_rms /= cnt;
+    dev_avg /= g.GetN();
+    dev_rms /= g.GetN();
     dev_rms = sqrt(dev_rms-dev_avg*dev_avg);
-    calc_dev_avg = calc_dev_avg / cnt;
+    calc_dev_avg /= g.GetN();
 
-    cout << "result " << med_avg << " " << med_rms << " " << dev_avg << " " << dev_rms << " " << med_beg << " " << med_end << " " << calc_dev_avg << endl;
+    calc_dev_rel -= g.GetN();
+    calc_dev_rel /= g.GetN();
+
+
+    cout << "result " << med_avg << " " << med_rms << " " << dev_avg << " " << dev_rms << " " << med_beg << " " << med_end << " " << calc_dev_avg << " " << calc_dev_rel << " " << fluct_abs << " " << fluct_rel << endl;
 }