Index: /trunk/MagicSoft/Mars/Changelog =================================================================== --- /trunk/MagicSoft/Mars/Changelog (revision 8133) +++ /trunk/MagicSoft/Mars/Changelog (revision 8134) @@ -53,4 +53,6 @@ - added a fix to allow paths containing a tilde + * datacenter/macros/plotoptical.C: + - added 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 +! Author(s): Daniela Dorner, 05/2005 +! +! 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 +#include + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#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 ", 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 << " " << 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; pGetMinimum(), 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; icd(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("", ""); +}