Index: /trunk/Mars/fact/plots/quality.C
===================================================================
--- /trunk/Mars/fact/plots/quality.C	(revision 17520)
+++ /trunk/Mars/fact/plots/quality.C	(revision 17521)
@@ -47,4 +47,5 @@
 
     TGraph g;
+    g.SetName("Threshold");
 
     while (file.GetNextRow())
@@ -54,4 +55,5 @@
     }
 
+    g.SetMinimum(281);
     g.SetMarkerStyle(kFullDotMedium);
     g.GetXaxis()->SetTimeDisplay(true);
@@ -70,15 +72,20 @@
 #include <pair>
 
-vector<pair<double, ln_equ_posn>> vec;
-
-ln_equ_posn FindPointing(Double_t time)
-{
-    for (int i=0; i<vec.size(); i++)
-        if (time<vec[i].first)
-            return i==0 ? ln_equ_posn() : vec[i-1].second;
-
-    return vec[vec.size()-1].second;
-}
-
+vector<pair<double, Nova::EquPosn>> vecp;
+
+Nova::EquPosn FindPointing(Double_t time)
+{
+    for (int i=0; i<vecp.size(); i++)
+        if (time<vecp[i].first)
+        {
+            if (i==0)
+                return Nova::EquPosn();
+            else
+                return vecp[i-1].second;
+        }
+
+    return vecp[vecp.size()-1].second;
+}
+/*
 Float_t Prediction(Double_t time)
 {
@@ -93,9 +100,5 @@
 
     // Distance between source and moon
-    //const double angle =
-    //    MAstro::AngularDistance(TMath::DegToRad()*(90-moon.dec),
-    //                            TMath::DegToRad()*moon.ra,
-    //                            TMath::DegToRad()*(90-pos.dec),
-    //                            TMath::DegToRad()*pos.ra);
+    //const double angle = Nova::GetAngularSeparation(moon, hrzm)*TMath::DegToRad();
 
     // Current prediction
@@ -105,7 +108,72 @@
     double disk = Nova::GetLunarDisk(jd);
 
-    double lc = /*cang==0 ? 0 :*/ calt*pow(disk, 2.5);///sqrt(cang);
-
-    return 5+103*lc>4.5 ? 5+103*lc : 4.5;
+    //                                                                           semi-major axis
+    double lc = calt*pow(disk, 2.2)*pow(Nova::GetLunarEarthDist(jd)/384400, -2);///sqrt(cang);
+
+    return lc>0 ? 4+103*lc : 4;
+}*/
+
+Float_t Prediction(Double_t time)
+{
+    Double_t jd = time + 40587 + 2400000.5;
+    /*
+    Nova::EquPosn moon = Nova::GetLunarEquCoords(jd, 0.01);
+
+    Nova::EquPosn pos = FindPointing(time);
+
+    // get local position of moon
+    Nova::HrzPosn hrzm = Nova::GetHrzFromEqu(moon, jd);
+
+    // Distance between source and moon
+    const double angle = Nova::GetAngularSeparation(moon, pos);
+
+    // Distance between earth and moon relative to major semi-axis
+    const double dist  = Nova::GetLunarEarthDist(jd)/384400;
+
+    // Current prediction
+    double cang = 1-sin(angle*TMath::DegToRad());
+    double calt = sin(hrzm.alt*TMath::DegToRad());
+
+    double disk = Nova::GetLunarDisk(jd);
+
+    // light condition
+    double lc = sqrt(calt)*pow(disk, 2.3)*pow(cang+1, 1)*pow(dist, -2);
+
+    return lc>0 ? 7.2 + 69*lc : 7.2;
+    */
+
+    // Sun properties
+    Nova::EquPosn  sun  = Nova::GetSolarEquCoords(jd);
+    Nova::ZdAzPosn hrzs = Nova::GetHrzFromEqu(sun, jd);
+
+    // Get source position
+    Nova::EquPosn pos = FindPointing(time);
+
+    // Moon properties
+    Nova::EquPosn moon = Nova::GetLunarEquCoords(jd, 0.01);
+    Nova::HrzPosn hrzm = Nova::GetHrzFromEqu(moon, jd);
+    double        disk = Nova::GetLunarDisk(jd);
+
+    // Derived moon properties
+    double angle = Nova::GetAngularSeparation(moon, pos);
+    double edist = Nova::GetLunarEarthDist(jd)/384400;
+
+    // Current prediction
+    double sin_malt  = hrzm.alt<0 ? 0 : sin(hrzm.alt*TMath::DegToRad());
+    double sin_mdist = sin(angle*TMath::DegToRad());
+    double cos_salt  = cos(hrzs.zd*TMath::DegToRad());
+
+    double c0 = pow(disk,     2.52);
+    double c1 = pow(sin_malt, 0.72);
+    double c2 = pow(edist,   -2.00);
+    double c3 = exp(1.46*(1-sin_mdist)*(1-sin_mdist));
+    double c4 = exp(33.0)*exp(-20.1*(1-cos_salt)*(1-cos_salt));
+
+    double cur = 6.4 + 96.9*c0*c1*c2*c3 + c4;
+
+   // cout << cur << " " << hrzm.alt << " " << c1 << endl;
+
+    return cur;
+
 }
 
@@ -131,9 +199,9 @@
     while (file.GetNextRow())
     {
-        ln_equ_posn p;
+        Nova::EquPosn p;
         p.ra  = ra*15;
         p.dec = dec;
 
-        vec.push_back(make_pair(time, p));
+        vecp.push_back(make_pair(time, p));
     }
 
@@ -169,4 +237,6 @@
     TGraph g1;
     TGraph g2;
+    g1.SetName("Currents");
+    g2.SetName("Prediction");
 
     while (file.GetNextRow())
@@ -234,4 +304,6 @@
 
     TGraph g1, g2;
+    g1.SetName("TriggerRate");
+    g2.SetName("RelOnTime");
 
     while (file.GetNextRow())
@@ -246,5 +318,5 @@
 
     g1.SetMinimum(0);
-    g1.SetMaximum(200);
+    g1.SetMaximum(269);
     g1.SetMarkerStyle(kFullDotMedium);
     g1.GetXaxis()->SetTimeDisplay(true);
@@ -260,5 +332,5 @@
     gPad->GetCanvas()->cd(4);
 
-    gPad->SetGridx();
+    gPad->SetGrid();
     gPad->SetTopMargin(0);
     gPad->SetBottomMargin(0);
@@ -278,7 +350,45 @@
     g2.DrawClone("AP");
 
-
-    return 0;
-}
+    return 0;
+}
+
+void PlotRateQC(UInt_t night, MSQLServer &serv)
+{
+    TString query = 
+        "LEFT JOIN AnalysisResultsRunLP USING(fNight, fRunID) "
+        "WHERE fRunTypeKey=1 AND NOT ISNULL (AnalysisResultsRunLP.fNumEvtsAfterCleaning) AND fNight=";
+    query += night;
+
+    TTree *t = serv.GetTree("RunInfo", query);
+    if (!t)
+        return;
+
+    int save = gErrorIgnoreLevel;
+    gErrorIgnoreLevel = kFatal;
+
+    gROOT->SetSelectedPad(0);
+    gPad->GetCanvas()->cd(3);
+
+    t->Draw("AnalysisResultsRunLP.fNumEvtsAfterCleaning/AnalysisResultsRunLP.fOnTimeAfterCuts:(RunInfo.fRunStart+RunInfo.fRunStop)/2+9131*24*3600", "", "same");
+    TGraph *g = (TGraph*)gPad->GetPrimitive("Graph");
+    if (g)
+    {
+        g->SetName("CleaningRate");
+        g->SetMarkerColor(kRed);
+        g->SetMarkerStyle(29);//kFullDotMedium);
+    }
+
+    t->Draw("AnalysisResultsRunLP.fNumEvtsAfterQualCuts/AnalysisResultsRunLP.fOnTimeAfterCuts:(RunInfo.fRunStart+RunInfo.fRunStop)/2+9131*24*3600", "", "same");
+    g = (TGraph*)gPad->GetPrimitive("Graph");
+    if (g)
+    {
+        g->SetName("RateAfterQC");
+        g->SetMarkerColor(kBlue);
+        g->SetMarkerStyle(29);//kFullDotMedium);
+    }
+
+    gErrorIgnoreLevel = save;
+}
+
 
 Int_t PlotPointing(TArrayD **vec, TString fname)
@@ -304,4 +414,5 @@
 
     TGraph g;
+    g.SetName("Zd");
 
     while (file.GetNextRow())
@@ -309,10 +420,10 @@
             g.SetPoint(g.GetN(), time*24*3600, 90-zd);
 
-    g.SetMinimum(0);
+    g.SetMinimum(1);
     g.SetMaximum(90);
     g.SetMarkerStyle(kFullDotMedium);
     g.GetXaxis()->SetTimeDisplay(true);
     g.GetXaxis()->SetTimeFormat("%H:%M %F1995-01-01 00:00:00 GMT");
-    g.GetXaxis()->SetLabelSize(0.1);
+    g.GetXaxis()->SetLabelSize(0.12);
     g.GetYaxis()->SetLabelSize(0.1);
     g.GetYaxis()->SetTitle("ELEVATION");
@@ -323,8 +434,8 @@
     return 0;
 }
-/*
-Int_t PlotMem(TArrayD **vec, TString fname)
-{
-    fname += ".FAD_CONTROL_STATISTICS1.fits";
+
+Int_t PlotTemperature1(TArrayD **vec, TString fname)
+{
+    fname += ".TEMPERATURE_DATA.fits";
 
     fits file(fname.Data());
@@ -338,28 +449,27 @@
 
     Double_t time;
-    UInt_t buf[5];
-    ULong64_t mem[4];
-    UInt_t rate[2];
+    Float_t  temp;
 
     if (!file.SetPtrAddress("Time",  &time))
         return -1;
-    if (!file.SetPtrAddress("bufferInfo", buf) ||
-        !file.SetPtrAddress("memInfo",    mem) ||
-        !file.SetPtrAddress("rateNew",    rate))
+    if (!file.SetPtrAddress("T", &temp))
         return -1;
 
     TGraph g;
+    g.SetName("ContainerTemp");
 
     while (file.GetNextRow())
         if (Contains(vec, time))
-            g.SetPoint(g.GetN(), time*24*3600, mem[1]);
-
+            g.SetPoint(g.GetN(), time*24*3600, temp);
+
+    g.SetMinimum(-5);
+    g.SetMaximum(49);
     g.SetMarkerStyle(kFullDotMedium);
+    g.SetMarkerColor(kRed);
     g.GetXaxis()->SetTimeDisplay(true);
-    g.GetXaxis()->SetTimeFormat("%H:%M");
-    g.GetXaxis()->SetTimeOffset(0, "gmt");
-    g.GetXaxis()->SetLabelSize(0.12);
+    g.GetXaxis()->SetTimeFormat("%H:%M %F1995-01-01 00:00:00 GMT");
+    g.GetXaxis()->SetLabelSize(0.1);
     g.GetYaxis()->SetLabelSize(0.1);
-    //g.GetYaxis()->SetTitle("ELEVATION");
+    g.GetYaxis()->SetTitle("TEMP");
     g.GetYaxis()->SetTitleOffset(0.2);
     g.GetYaxis()->SetTitleSize(0.1);
@@ -368,6 +478,174 @@
     return 0;
 }
-*/
-void quality(UInt_t y=0, UInt_t m=0, UInt_t d=0)
+
+Int_t PlotTemperature2(TArrayD **vec, TString fname)
+{
+    fname += ".MAGIC_WEATHER_DATA.fits";
+
+    fits file(fname.Data());
+    if (!file)
+    {
+        cerr << fname << ": " << gSystem->GetError() << endl;
+        return -2;
+    }
+
+    //cout << fname << endl;
+
+    Double_t time;
+    Float_t  temp;
+
+    if (!file.SetPtrAddress("Time",  &time))
+        return -1;
+    if (!file.SetPtrAddress("T", &temp))
+        return -1;
+
+    TGraph g;
+    g.SetName("OutsideTemp");
+
+    while (file.GetNextRow())
+        if (Contains(vec, time))
+            g.SetPoint(g.GetN(), time*24*3600, temp);
+
+    g.SetMarkerStyle(kFullDotMedium);
+    g.SetMarkerColor(kBlue);
+    g.DrawClone("P");
+
+    return 0;
+}
+
+Int_t PlotTemperature3(TArrayD **vec, TString fname)
+{
+    fname += ".FSC_CONTROL_TEMPERATURE.fits";
+
+    fits file(fname.Data());
+    if (!file)
+    {
+        cerr << fname << ": " << gSystem->GetError() << endl;
+        return -2;
+    }
+
+    //cout << fname << endl;
+
+    Double_t time;
+    Float_t  temp[31];
+
+    if (!file.SetPtrAddress("Time",  &time))
+        return -1;
+    if (!file.SetPtrAddress("T_sens", temp))
+        return -1;
+
+    TGraph g, g1, g2;
+    g.SetName("SensorTempAvg");
+    g1.SetName("SensorTempMin");
+    g2.SetName("SensorTempMax");
+
+    while (file.GetNextRow())
+        if (Contains(vec, time))
+        {
+            float min =  100;
+            float max = -100;
+            double avg = 0;
+            int num = 0;
+            for (int i=0; i<31; i++)
+                if (temp[i]!=0)
+                {
+                    avg += temp[i];
+                    num++;
+
+                    min = TMath::Min(min, temp[i]);
+                    max = TMath::Max(max, temp[i]);
+                }
+
+            g.SetPoint(g.GetN(), time*24*3600, avg/num);
+            g1.SetPoint(g1.GetN(), time*24*3600, min);
+            g2.SetPoint(g2.GetN(), time*24*3600, max);
+        }
+
+    g.SetMarkerStyle(kFullDotMedium);
+    g.DrawClone("P");
+
+    /*
+     g1.SetLineWidth(1);
+     g1.DrawClone("L");
+
+     g2.SetLineWidth(1);
+     g2.DrawClone("L");
+    */
+    return 0;
+}
+
+Int_t PlotTemperature4(TArrayD **vec, TString fname)
+{
+    fname += ".FAD_CONTROL_TEMPERATURE.fits";
+
+    fits file(fname.Data());
+    if (!file)
+    {
+        cerr << fname << ": " << gSystem->GetError() << endl;
+        return -2;
+    }
+
+    //cout << fname << endl;
+
+    Double_t time;
+    Float_t  temp[160];
+
+    if (!file.SetPtrAddress("Time",  &time))
+        return -1;
+    if (!file.SetPtrAddress("Data1", temp) &&
+        !file.SetPtrAddress("temp", temp))
+        return -1;
+
+    Int_t num = file.GetN("temp")==0 ? file.GetN("Data1") : file.GetN("temp");
+    Int_t beg = num==82 ? 2 : 0;
+
+    TGraphErrors g1;
+    TGraph g2,g3;
+
+    g1.SetName("FadTempAvg");
+    g2.SetName("FadTempMin");
+    g3.SetName("FadTempMax");
+
+    while (file.GetNextRow())
+        if (Contains(vec, time))
+        {
+            double avg = 0;
+            double rms = 0;
+            float min =  100;
+            float max = -100;
+            for (int i=beg; i<num; i++)
+            {
+                avg += temp[i];
+                rms += temp[i]*temp[i];
+
+                min = TMath::Min(min, temp[i]);
+                max = TMath::Max(max, temp[i]);
+            }
+
+            avg /= num-beg;
+            rms /= num-beg;
+
+            g1.SetPoint(g1.GetN(), time*24*3600, avg);
+            g1.SetPointError(g1.GetN()-1, 0, sqrt(rms-avg*avg));
+
+            g2.SetPoint(g2.GetN(), time*24*3600, min);
+            g3.SetPoint(g3.GetN(), time*24*3600, max);
+        }
+
+    g1.SetLineColor(kGreen);
+    g1.DrawClone("[]");
+
+    g2.SetLineColor(kGreen);
+    g2.SetLineWidth(1);
+    g2.DrawClone("L");
+
+    g3.SetLineColor(kGreen);
+    g3.SetLineWidth(1);
+    g3.DrawClone("L");
+
+    return 0;
+}
+
+void quality(UInt_t y=0, UInt_t m=0, UInt_t d=0, TString outpath="./")
 {
     // To get correct dates in the histogram you have to add
@@ -377,7 +655,7 @@
     {
         UInt_t nt = MTime(MTime(-1).GetMjd()-1.5).GetNightAsInt();
-        y = nt/10000;
+        y =  nt/10000;
         m = (nt/100)%100;
-        d = nt%100;
+        d =  nt%100;
 
         cout << y << "/" << m << "/" << d << endl;
@@ -409,8 +687,8 @@
     {
         TString query;
-        query += "SELECT fRunID, fRunStart, fRunStop FROM runinfo";
+        query += "SELECT fRunID, fRunStart, fRunStop FROM RunInfo";
         query += " WHERE fNight=";
         query += night;
-        query += " AND fRunTypeKEY=1 ORDER BY fRunStart";
+        query += " AND fRunTypeKey=1 ORDER BY fRunStart";
 
         TSQLResult *res = serv.Query(query);
@@ -438,12 +716,15 @@
         if (n==0)
             cout << "WARNING - No data runs in db, displaying all data." << endl;
+        else
+            cout << "Num: " << n << "\n" << endl;
     }
 
     TCanvas *c = new TCanvas("quality", Form("Quality %04d/%02d/%02d", y, m, d), 1280, 960);
-    c->Divide(1, 5, 1e-5, 1e-5);
+    c->Divide(1, 6, 1e-5, 1e-5);
 
     gROOT->SetSelectedPad(0);
+
     c->cd(1);
-    gPad->SetGridx();
+    gPad->SetGrid();
     gPad->SetTopMargin(0);
     gPad->SetRightMargin(0.001);
@@ -454,5 +735,5 @@
     gROOT->SetSelectedPad(0);
     c->cd(2);
-    gPad->SetGridx();
+    gPad->SetGrid();
     gPad->SetTopMargin(0);
     gPad->SetRightMargin(0.001);
@@ -463,5 +744,5 @@
     gROOT->SetSelectedPad(0);
     c->cd(3);
-    gPad->SetGridx();
+    gPad->SetGrid();
     gPad->SetTopMargin(0);
     gPad->SetBottomMargin(0);
@@ -469,22 +750,27 @@
     gPad->SetLeftMargin(0.04);
     cout << PlotRate(runs, fname) << endl;
+    cout << PlotRateQC(night, serv) << endl;
 
     gROOT->SetSelectedPad(0);
     c->cd(5);
-    gPad->SetGridx();
+    gPad->SetGrid();
+    gPad->SetTopMargin(0);
+    gPad->SetBottomMargin(0);
+    gPad->SetRightMargin(0.001);
+    gPad->SetLeftMargin(0.04);
+    cout << PlotPointing(runs, fname) << endl;
+
+    gROOT->SetSelectedPad(0);
+    c->cd(6);
+    gPad->SetGrid();
     gPad->SetTopMargin(0);
     gPad->SetRightMargin(0.001);
     gPad->SetLeftMargin(0.04);
-    cout << PlotPointing(runs, fname) << endl;
-/*
-    gROOT->SetSelectedPad(0);
-    c->cd(6);
-    gPad->SetGridx();
-    gPad->SetTopMargin(0);
-    gPad->SetRightMargin(0.001);
-    gPad->SetLeftMargin(0.04);
-    cout << PlotMem(runs, fname) << endl;
-*/
-    c->SaveAs(Form("%04d%02d%02d-quality.pdf", y, m, d));
+    cout << PlotTemperature1(runs, fname) << endl;
+    cout << PlotTemperature2(runs, fname) << endl;
+    cout << PlotTemperature3(runs, fname) << endl;
+    cout << PlotTemperature4(runs, fname) << endl;
+
+    c->SaveAs(Form("%s/%04d%02d%02d-quality.png", outpath.Data(), y, m, d));
 }