Index: trunk/MagicSoft/Mars/datacenter/macros/plotoptical.C
===================================================================
--- trunk/MagicSoft/Mars/datacenter/macros/plotoptical.C	(revision 8134)
+++ trunk/MagicSoft/Mars/datacenter/macros/plotoptical.C	(revision 8134)
@@ -0,0 +1,536 @@
+/* ======================================================================== *\
+! $Name: not supported by cvs2svn $:$Id: plotoptical.C,v 1.1 2006-10-19 18:45:51 tbretz Exp $
+! --------------------------------------------------------------------------
+!
+! *
+! * This file is part of MARS, the MAGIC Analysis and Reconstruction
+! * Software. It is distributed to you in the hope that it can be a useful
+! * and timesaving tool in analysing Data of imaging Cerenkov telescopes.
+! * It is distributed WITHOUT ANY WARRANTY.
+! *
+! * Permission to use, copy, modify and distribute this software and its
+! * documentation for any purpose is hereby granted without fee,
+! * provided that the above copyright notice appear in all copies and
+! * that both that copyright notice and this permission notice appear
+! * in supporting documentation. It is provided "as is" without express
+! * or implied warranty.
+! *
+!
+!
+!   Author(s): Thomas Bretz, 05/2005 <mailto:tbretz@astro.uni-wuerzburg.de>
+!   Author(s): Daniela Dorner, 05/2005 <mailto:dorner@astro.uni-wuerzburg.de>
+!
+!   Copyright: MAGIC Software Development, 2000-2006
+!
+!
+\* ======================================================================== */
+
+/////////////////////////////////////////////////////////////////////////////
+//
+// plotdb.C
+// ========
+//
+// This macro is used to read quality parameters from the DB and plot them.
+// 
+// The parameters are from the following files:
+// calib*.root:mean conversion factor, mean arrival time, rms arrival time
+// (each parameter for inner and outer camera)
+// signal*.root: mean pedestal rms (for inner and outer camera)
+// star*.root: PSF, # of Muons, Effective OnTime, Muon rate,
+// Ratio MC/Data(MuonSize) and mean number of islands
+// 
+// In the DB these values are stored in the tables Calibration and Star.
+// 
+// Usage:
+//   .x plotdb.C   --> all values in the DB are plotted
+// You can chose are certain period:
+//   .x plotdb.C(25)   --> all values from period 25 are plotted
+// or a time period from a certain date to a certain date
+//   .x plotdb.C("2004-11-14 00:00:00", "2005-02-28 00:00:00")
+//  --> all values from 14.11.2004 0h to 28.2.2005 0h are plotted
+// or all data, but with dataset data highlighted
+//   .x plotdb.C("dataset.txt")
+//  --> the sequences defined in dataset.txt are highlighted (blue:on, red:off)
+//  --> You can also add a dataset-name as last argument to one of the
+//      calls above
+//
+// Make sure, that database and password are corretly set in a resource
+// file called sql.rc and the resource file is found.
+//
+/////////////////////////////////////////////////////////////////////////////
+#include <iostream>
+#include <iomanip>
+
+#include <TH1.h>
+#include <TEnv.h>
+#include <TPad.h>
+#include <TLine.h>
+#include <TText.h>
+#include <TFrame.h>
+#include <TStyle.h>
+#include <TGraph.h>
+#include <TCanvas.h>
+#include <TSQLRow.h>
+#include <TSQLResult.h>
+
+#include "MTime.h"
+#include "MAstro.h"
+#include "MDataSet.h"
+#include "MSQLMagic.h"
+#include "MStatusDisplay.h"
+
+class MPlot : public MParContainer
+{
+private:
+    MSQLMagic &fServer;
+
+    MDataSet *fDataSet;
+
+    TString fRequestFrom;
+    TString fRequestTo;
+    Int_t   fRequestPeriod;
+
+    Float_t fPlotMin;
+    Float_t fPlotMax;
+
+    Float_t fHistMin;
+    Float_t fHistMax;
+
+    TString fDescription;
+    TString fNameTab;
+
+    TString fCondition;
+
+    void PlotTable(TSQLResult &res, TString name, Float_t fmin, Float_t fmax, Float_t resolution)
+    {
+        gStyle->SetOptStat(111111);
+
+        TSQLRow *row;
+
+        TGraph gt;
+        gt.SetNameTitle(name, Form("%s vs Time", name.Data()));
+        gt.SetMarkerStyle(kFullDotMedium);
+
+        TGraph gz;
+        gz.SetNameTitle(name, Form("%s vs <Zd>", name.Data()));
+        gz.SetMarkerStyle(kFullDotMedium);
+
+        TGraph gt0, gt1;
+        gt0.SetMarkerColor(kRed);
+        gt1.SetMarkerColor(kBlue);
+        gt0.SetMarkerStyle(kFullDotLarge);
+        gt1.SetMarkerStyle(kFullDotLarge);
+
+        TGraph gz0, gz1;
+        gz0.SetMarkerColor(kRed);
+        gz1.SetMarkerColor(kBlue);
+        gz0.SetMarkerStyle(kFullDotLarge);
+        gz1.SetMarkerStyle(kFullDotLarge);
+
+        Int_t first = -1;
+        Int_t last  = -1;
+
+        while ((row=res.Next()))
+        {
+            const char *date = (*row)[0];
+            const char *zd   = (*row)[1];
+            const char *val  = (*row)[2];
+            const char *snum = res.GetFieldCount()>3 ? (*row)[3] : 0;
+            if (!date || !val || !zd)
+                continue;
+
+            MTime t(date);
+            if (!t.SetSqlDateTime(date))
+                continue;
+
+            if (fRequestPeriod>0 && MAstro::GetMagicPeriod(t.GetMjd())!=fRequestPeriod)
+                continue;
+
+            if (first<0)
+                first = TMath::Nint(TMath::Floor(t.GetMjd()));
+            last = TMath::Nint(TMath::Ceil(t.GetMjd()));
+
+            UInt_t seq = snum ? atoi(snum) : 0;
+    
+            Float_t value = atof(val);
+            Float_t zenith = atof(zd);
+
+            if (fDataSet)
+            {
+                if (fDataSet->HasOnSequence(seq))
+                {
+                    gt1.SetPoint(gt1.GetN(), t.GetAxisTime(), value);
+                    gz1.SetPoint(gz1.GetN(), zenith, value);
+                }
+
+                if (fDataSet->HasOffSequence(seq))
+                {
+                    gt0.SetPoint(gt0.GetN(), t.GetAxisTime(), value);
+                    gz0.SetPoint(gz0.GetN(), zenith, value);
+                }
+            }
+
+            gt.SetPoint(gt.GetN(), t.GetAxisTime(), value);
+            gz.SetPoint(gz.GetN(), zenith, value);
+        }
+
+        TString title = fNameTab.IsNull() ? name(name.First('.')+2, name.Length()) : fNameTab;
+        cerr << setprecision(4) << setw(10) << title << ":   ";
+        if (gt.GetN()==0)
+        {
+            cerr << "     <empty>" << endl;
+            return;
+        }
+        cerr << setw(8) << gt.GetMean(2) << "+-" << setw(8) << gt.GetRMS(2) << "   ";
+        if (gt0.GetN()>0 || gt1.GetN()>0)
+        {
+            cerr << setw(8) << gt1.GetMean(2) << "+-" << setw(8) << gt1.GetRMS(2) << "   ";
+            cerr << setw(8) << gt0.GetMean(2) << "+-" << setw(8) << gt0.GetRMS(2);
+        }
+        cerr << endl;
+
+        // If this is done earlier the plots remain empty since root 5.12/00
+        if (fmax>fmin)
+        {
+            gt.SetMinimum(fmin);
+            gt.SetMaximum(fmax);
+            gz.SetMinimum(fmin);
+            gz.SetMaximum(fmax);
+        }
+
+        gROOT->SetSelectedPad(0);
+
+        TCanvas &c = fDisplay ? fDisplay->AddTab(title) : *new TCanvas;
+        c.SetFillColor(kWhite);
+        c.SetBorderMode(0);
+        c.Divide(1,2);
+
+        TVirtualPad *pad = gPad;
+        pad->cd(2);
+        gPad->SetBorderMode(0);
+        gPad->SetFrameBorderMode(0);
+        gPad->SetGridy();
+
+        gPad->SetLeftMargin(0.06);
+        gPad->SetRightMargin(0.06);
+        gPad->SetBottomMargin(0.08);
+
+        TH1 *h = gt.GetHistogram();
+
+        h->SetXTitle("Time");
+        h->SetYTitle(name);
+        h->GetXaxis()->SetTimeDisplay(1);
+        h->GetYaxis()->SetTitleOffset(0.8);
+        h->GetXaxis()->SetTitleOffset(1.0);
+        h->GetXaxis()->SetLabelOffset(0.01);
+
+        gt.DrawClone("AP");
+        if (gt0.GetN()>0)
+            gt0.DrawClone("P");
+        if (gt1.GetN()>0)
+            gt1.DrawClone("P");
+
+        TLine l;
+        TText t;
+        Int_t num=0;
+        l.SetLineStyle(kDotted);
+        l.SetLineColor(kBlue);
+        t.SetTextColor(kBlue);
+        l.SetLineWidth(1);
+        t.SetTextSize(h->GetXaxis()->GetLabelSize());
+        t.SetTextAlign(21);
+        Int_t p0 = MAstro::GetMagicPeriod(first);
+        for (Int_t p = first; p<last; p++)
+        {
+            Int_t p1 = MAstro::GetMagicPeriod(p);
+            if (p1!=p0)
+            {
+                l.DrawLine(MTime(p).GetAxisTime(), h->GetMinimum(), MTime(p).GetAxisTime(), h->GetMaximum());
+                t.DrawText(MTime(p+15).GetAxisTime(), h->GetMaximum(), Form("%d", p1));
+                num++;
+            }
+            p0 = p1;
+        }
+        if (num<4)
+            gPad->SetGridx();
+
+        const Double_t min = fHistMin>fHistMax ? h->GetMinimum()-resolution/2 : fHistMin;
+        const Double_t max = fHistMin>fHistMax ? h->GetMaximum()+resolution/2 : fHistMax;
+
+        // Use this to save the pad with the time development to a file
+        //gPad->SaveAs(Form("plotdb-%s.eps", title.Data()));
+
+        pad->cd(1);
+        gPad->SetBorderMode(0);
+        gPad->SetFrameBorderMode(0);
+        gPad->Divide(2,1);
+
+        TVirtualPad *pad2 = gPad;
+        pad2->cd(1);
+        gPad->SetBorderMode(0);
+        gPad->SetFrameBorderMode(0);
+        gPad->SetGridx();
+        gPad->SetGridy();
+
+        const Int_t n = resolution>0 ? TMath::Nint((max-min)/resolution) : 50;
+
+        TH1F hist("Hist", Form("Distribution of %s", fDescription.IsNull() ? name.Data() : fDescription.Data()), n, min, max);
+        hist.SetDirectory(0);
+
+        for (int i=0; i<gt.GetN(); i++)
+            hist.Fill(gt.GetY()[i]);
+
+        if (fDescription.IsNull())
+            hist.SetXTitle(name);
+        hist.SetYTitle("Counts");
+
+        hist.DrawCopy("");
+
+        pad2->cd(2);
+        gPad->SetBorderMode(0);
+        gPad->SetFrameBorderMode(0);
+        gPad->SetGridy();
+
+        TH1 *h2 = gz.GetHistogram();
+
+        h2->SetXTitle("Zd");
+        h2->SetYTitle(name);
+
+        gz.DrawClone("AP");
+        if (gz0.GetN()>0)
+            gz0.DrawClone("P");
+        if (gz1.GetN()>0)
+            gz1.DrawClone("P");
+    }
+
+public:
+    MPlot(MSQLMagic &server) : fServer(server), fDataSet(NULL),
+        fRequestPeriod(-1), fPlotMin(0), fPlotMax(-1), fHistMin(0), fHistMax(-1)
+    {
+    }
+    ~MPlot()
+    {
+        if (fDataSet)
+            delete fDataSet;
+    }
+    void SetDataSet(const TString filename)
+    {
+        if (fDataSet)
+        {
+            delete fDataSet;
+            fDataSet = NULL;
+        }
+        if (!filename.IsNull())
+            fDataSet = new MDataSet(filename);
+    }
+    void SetPlotRange(Float_t min, Float_t max, Int_t n=5)
+    { fPlotMin = min; fPlotMax = max; }
+    void SetHistRange(Float_t min, Float_t max)
+    { fHistMin = min; fHistMax = max; }
+    void SetRequestRange(const char *from="", const char *to="")
+    { fRequestFrom = from; fRequestTo = to; }
+    void SetRequestPeriod(Int_t n=-1)
+    { fRequestPeriod = n; }
+    void SetCondition(const char *cond="")
+    { fCondition = cond; }
+    void SetDescription(const char *d, const char *t=0) { fDescription = d; fNameTab = t; }
+
+    Int_t QueryKeyOfSource(TString src)
+    {
+        return fServer.QueryKeyOfName("Object", src, kFALSE);
+    }
+
+    Bool_t Plot(const char *value, Float_t min=0, Float_t max=-1, Float_t resolution=0)
+    {
+        TString named  = "OpticalData.fTimestamp";
+        TString named2 = "fSkyLevel";
+        TString namev  = value;
+        TString join   = "fSequenceFirst";
+
+        TString tablev = namev(0, namev.First('.'));
+        TString valuev = namev(namev.First('.')+1, namev.Length());
+
+        TString tabled = named(0, named.First('.'));
+        TString valued = named(named.First('.')+1, named.Length());
+
+        TString query;
+        query  = Form("select %s, %s, %s ", valued.Data(), named2.Data(), valuev.Data());
+        query += Form("from %s ",           tabled.Data());
+
+        const Bool_t interval = !fRequestFrom.IsNull() && !fRequestTo.IsNull();
+
+        if (!fDataSet && !interval && tabled=="Star")
+        {
+            if (!query.Contains("Star.fSequenceFirst"))
+                query += "left join Star on Sequences.fSequenceFirst=Star.fSequenceFirst ";
+            query += "where Star.fEffOnTime>300 ";
+        }
+
+        if (interval)
+        {
+            query += query.Contains(" where ") ? "and " : "where ";
+            query += Form("fRunStart between '%s' and '%s' ",
+                          fRequestFrom.Data(), fRequestTo.Data());
+        }
+
+        if (!fCondition.IsNull())
+        {
+            query += "where ";
+            query += fCondition;
+            query += " ";
+        }
+
+        query += "order by fTimestamp";
+
+        TSQLResult *res = fServer.Query(query);
+        if (!res)
+            return kFALSE;
+
+        if (max>min)
+            PlotTable(*res, namev, min, max, resolution);
+        else
+            PlotTable(*res, namev, fPlotMin, fPlotMax, resolution);
+
+        delete res;
+        return kTRUE;
+    }
+};
+
+void plotall(MPlot &plot, TString source)
+{
+    TString cond = "fStatusKEY=13";
+    if (!source.IsNull())
+    {
+        const Int_t key = plot.QueryKeyOfSource(source);
+        if (key<0)
+            return;
+        cond += Form(" and fObjectKEY=%d", key);
+
+    }
+    plot.SetCondition(cond);
+
+    plot.SetDescription("Exposure;T_{E} [s]", "Expo");
+    plot.Plot("OpticalData.fExposure",  0, 900, 60);
+
+    plot.SetDescription("Sky Level;B [s^{-1}]", "SkyLvl");
+    plot.Plot("OpticalData.fSkyLevel/OpticalData.fExposure", 0, 5.0, 0.01);
+
+    plot.SetDescription("Full Width Half Maximum;FWHM [s^{-1}]", "Fwhm");
+    plot.Plot("OpticalData.fFWHM/OpticalData.fExposure",  0, 0.05, 0.001);
+
+    plot.SetDescription("Aperture Radius;R_{A}", "ApRad");
+    plot.Plot("OpticalData.fApertureRadius", 0, 10, 1);
+
+    plot.SetDescription("Instrumental Magnitude;M_{I}\\cdot s^{-1}", "InstMag/s");
+    plot.Plot("OpticalData.fInstrumentalMag/OpticalData.fExposure",  0, 0.2, 0.005);
+
+    plot.SetDescription("Instrumental Magnitude Error;\\sigma_{M}\\cdot s^{-1}", "MagErr/s");
+    plot.Plot("OpticalData.fInstrumentalMagErr/OpticalData.fExposure",  0, 0.01, 0.0002);
+
+    plot.SetDescription("Instrumental Magnitude;M_{I}", "InstMag");
+    plot.Plot("OpticalData.fInstrumentalMag",  0, 30, 0.5);
+
+    plot.SetDescription("Instrumental Magnitude Error;\\sigma_{M}", "MagErr");
+    plot.Plot("OpticalData.fInstrumentalMagErr",  0, 1, 0.01);
+}
+
+int plotoptical(TString from, TString to, const char *source=0)
+{
+    TEnv env("sql.rc");
+
+    MSQLMagic serv(env);
+    if (!serv.IsConnected())
+    {
+        cout << "ERROR - Connection to database failed." << endl;
+        return 0;
+    }
+
+    cout << "plotoptical" << endl;
+    cout << "-----------" << endl;
+    cout << endl;
+    cout << "Connected to " << serv.GetName() << endl;
+    cout << endl;
+
+    MStatusDisplay *d = new MStatusDisplay;
+    d->SetWindowName(serv.GetName());
+    d->SetTitle(serv.GetName());
+
+    MPlot plot(serv);
+//    plot.SetDataSet(dataset);
+    plot.SetDisplay(d);
+    plot.SetRequestRange(from, to);
+    plotall(plot, source);
+    d->SaveAsRoot("plotoptical.root");
+    d->SaveAsPS("plotoptical.ps");
+
+    return 1;
+}
+
+int plotoptical(const char *source)
+{
+    TEnv env("sql.rc");
+
+    MSQLMagic serv(env);
+    if (!serv.IsConnected())
+    {
+        cout << "ERROR - Connection to database failed." << endl;
+        return 0;
+    }
+
+    cout << "plotoptical" << endl;
+    cout << "-----------" << endl;
+    cout << endl;
+    cout << "Connected to " << serv.GetName() << endl;
+    cout << endl;
+
+    MStatusDisplay *d = new MStatusDisplay;
+    d->SetWindowName(serv.GetName());
+    d->SetTitle(serv.GetName());
+
+    MPlot plot(serv);
+//    plot.SetDataSet(ds);
+    plot.SetDisplay(d);
+    plot.SetRequestRange("", "");
+    plotall(plot, source);
+    d->SaveAsRoot("plotoptical.root");
+    d->SaveAsPS("plotoptical.ps");
+
+    return 1;
+}
+
+int plotoptical(Int_t period, const char *source="")
+{
+    TEnv env("sql.rc");
+
+    MSQLMagic serv(env);
+    if (!serv.IsConnected())
+    {
+        cout << "ERROR - Connection to database failed." << endl;
+        return 0;
+    }
+
+    cout << "plotoptical" << endl;
+    cout << "-----------" << endl;
+    cout << endl;
+    cout << "Connected to " << serv.GetName() << endl;
+    cout << endl;
+
+    MStatusDisplay *d = new MStatusDisplay;
+    d->SetWindowName(serv.GetName());
+    d->SetTitle(serv.GetName());
+
+    MPlot plot(serv);
+//    plot.SetDataSet(dataset);
+    plot.SetDisplay(d);
+    plot.SetRequestPeriod(period);
+    plotall(plot, source);
+    d->SaveAsRoot("plotoptical.root");
+    d->SaveAsPS("plotoptical.ps");
+
+    return 1;
+}
+
+int plotoptical()
+{
+    return plotoptical("", "");
+}
