Index: /trunk/MagicSoft/Mars/Changelog =================================================================== --- /trunk/MagicSoft/Mars/Changelog (revision 4905) +++ /trunk/MagicSoft/Mars/Changelog (revision 4906) @@ -116,4 +116,8 @@ * mfileio/MReadRflFile.[h,cc]: - added missing abstract functions + + * star.cc, callisto.cc, mjobs/MJStar.[h,cc], + mjobs/MJCalibrateSignal.[h,cc]: + - added option to process monte carlo files Index: /trunk/MagicSoft/Mars/callisto.cc =================================================================== --- /trunk/MagicSoft/Mars/callisto.cc (revision 4905) +++ /trunk/MagicSoft/Mars/callisto.cc (revision 4906) @@ -50,4 +50,6 @@ gLog << " -c Calculate the calibration constants" << endl; gLog << " -y Extract and calibrate signal" << endl << endl; + gLog << " Input File Options:" << endl; + gLog << " -mc You must use this for MC files (PRELIMINARY)" << endl << endl; gLog << " Options:" << endl; gLog.Usage(); @@ -121,4 +123,5 @@ const Bool_t kOverwrite = arg.HasOnlyAndRemove("-f"); const Bool_t kForceExec = arg.HasOnlyAndRemove("-ff"); + const Bool_t kIsMC = arg.HasOnlyAndRemove("-mc"); const TString kInpathD = arg.GetStringAndRemove("--ind=", ""); @@ -374,5 +377,5 @@ // Where to search for calibration files - if (!job2.ProcessFile(job1.GetPedestalCam())) + if (!job2.ProcessFile(job1.GetPedestalCam(), kIsMC)) return -1; Index: /trunk/MagicSoft/Mars/mhist/MHEffectiveOnTime.cc =================================================================== --- /trunk/MagicSoft/Mars/mhist/MHEffectiveOnTime.cc (revision 4905) +++ /trunk/MagicSoft/Mars/mhist/MHEffectiveOnTime.cc (revision 4906) @@ -63,5 +63,6 @@ // MHEffectiveOnTime::MHEffectiveOnTime(const char *name, const char *title) - : fPointPos(0), fTime(0), fParam(0), fIsFinalized(kFALSE), fInterval(60) + : fPointPos(0), fTime(0), fParam(0), fIsFinalized(kFALSE), fInterval(60), + fNameProjDeltaT(Form("DeltaT_%p", this)), fNameProjTheta(Form("Theta_%p", this)) { // @@ -72,18 +73,32 @@ // Main histogram - fHTimeDiff.SetXTitle("\\Delta t [s]"); - fHTimeDiff.SetYTitle("\\Theta [\\circ]"); - fHTimeDiff.UseCurrentStyle(); - fHTimeDiff.SetDirectory(NULL); + fHDeltaT.SetName("DeltaT"); + fHDeltaT.SetXTitle("\\Delta t [s]"); + fHDeltaT.SetYTitle("\\Theta [\\circ]"); + fHDeltaT.UseCurrentStyle(); + fHDeltaT.SetDirectory(NULL); // effective on time versus theta - fHEffOn.SetName("EffOnTime"); - fHEffOn.SetTitle("Effective On Time T_{eff}"); - fHEffOn.SetXTitle("\\Theta [\\circ]"); - fHEffOn.SetYTitle("T_{eff} [s]"); - fHEffOn.UseCurrentStyle(); - fHEffOn.SetDirectory(NULL); - fHEffOn.GetYaxis()->SetTitleOffset(1.2); + fHEffOnTheta.SetName("EffOnTime"); + fHEffOnTheta.SetTitle("Effective On Time T_{eff}"); + fHEffOnTheta.SetXTitle("\\Theta [\\circ]"); + fHEffOnTheta.SetYTitle("T_{eff} [s]"); + fHEffOnTheta.UseCurrentStyle(); + fHEffOnTheta.SetDirectory(NULL); + fHEffOnTheta.GetYaxis()->SetTitleOffset(1.2); //fHEffOn.Sumw2(); + + // effective on time versus time + fHEffOnTime.SetName("EffOnTime"); + fHEffOnTime.SetTitle("Effective On Time T_{eff}"); + fHEffOnTime.SetXTitle("Time"); + fHEffOnTime.SetYTitle("T_{eff} [s]"); + fHEffOnTime.UseCurrentStyle(); + fHEffOnTime.SetDirectory(NULL); + fHEffOnTime.GetYaxis()->SetTitleOffset(1.2); + fHEffOnTime.GetXaxis()->SetLabelSize(0.033); + fHEffOnTime.GetXaxis()->SetTimeFormat("%H:%M:%S %F1995-01-01 00:00:00 GMT"); + fHEffOnTime.GetXaxis()->SetTimeDisplay(1); + fHEffOnTime.Sumw2(); // chi2/NDF versus theta @@ -128,12 +143,12 @@ // setup binning MBinning btheta("BinningTheta"); - btheta.SetEdgesCos(51, 0, 60); - - MBinning btime("BinningTimeDiff"); + btheta.SetEdgesCos(101, 0, 60); + + MBinning btime("BinningDeltaT"); btime.SetEdges(50, 0, 0.1); - MH::SetBinning(&fHTimeDiff, &btime, &btheta); - - btheta.Apply(fHEffOn); + MH::SetBinning(&fHDeltaT, &btime, &btheta); + + btheta.Apply(fHEffOnTheta); btheta.Apply(fHChi2); btheta.Apply(fHLambda); @@ -141,8 +156,4 @@ btheta.Apply(fHProb); } - -// FIXME: Just for a preliminary check -static Double_t testval = 0; -static Double_t testerr = 0; // -------------------------------------------------------------------------- @@ -167,5 +178,5 @@ return kFALSE; - const MBinning* binsdtime = (MBinning*)plist->FindObject("BinningTimeDiff"); + const MBinning* binsdtime = (MBinning*)plist->FindObject("BinningDeltaT"); const MBinning* binstheta = (MBinning*)plist->FindObject("BinningTheta"); if (!binstheta || !binsdtime) @@ -173,7 +184,7 @@ else { - SetBinning(&fHTimeDiff, binsdtime, binstheta); - - binstheta->Apply(fHEffOn); + SetBinning(&fHDeltaT, binsdtime, binstheta); + + binstheta->Apply(fHEffOnTheta); binstheta->Apply(fHChi2); binstheta->Apply(fHLambda); @@ -200,5 +211,5 @@ // nbins = number of Theta bins - const Int_t nbins = fHTimeDiff.GetNbinsY(); + const Int_t nbins = fHDeltaT.GetNbinsY(); TH1D *h=0; @@ -206,5 +217,5 @@ { // TH1D &h = *hist->ProjectionX("Calc-theta", i, i); - h = fHTimeDiff.ProjectionX(name, i, i, "E"); + h = fHDeltaT.ProjectionX(name, i, i, "E"); Double_t res[7]; @@ -213,6 +224,6 @@ // the effective on time is Nm/lambda - fHEffOn.SetBinContent(i, res[0]); - fHEffOn.SetBinError (i, res[1]); + fHEffOnTheta.SetBinContent(i, res[0]); + fHEffOnTheta.SetBinError (i, res[1]); // plot chi2-probability of fit @@ -346,26 +357,33 @@ TH1D *h0=0; - padsave->cd(1); - if ((h0 = (TH1D*)gPad->FindObject("TimeDiff"))) + padsave->GetPad(1)->cd(1); + if ((h0 = (TH1D*)gPad->FindObject(fNameProjDeltaT))) { - h0 = fHTimeDiff.ProjectionX("TimeDiff", -1, 9999, "E"); + h0 = fHDeltaT.ProjectionX(fNameProjDeltaT, -1, 9999, "E"); if (h0->GetEntries()>0) gPad->SetLogy(); } - padsave->cd(2); - if ((h0 = (TH1D*)gPad->FindObject("Theta"))) - fHTimeDiff.ProjectionY("Theta", -1, 9999, "E"); + padsave->GetPad(2)->cd(1); + if ((h0 = (TH1D*)gPad->FindObject(fNameProjTheta))) + fHDeltaT.ProjectionY(fNameProjTheta, -1, 9999, "E"); if (!fIsFinalized) FitThetaBins(); } - if (o==(TString)"paint") + + TH1D *h=0; + if (o==(TString)"theta") + h = &fHEffOnTheta; + if (o==(TString)"time") + h = &fHEffOnTime; + + if (h) { Double_t error = 0; - for (int i=0; iGetNbins(); i++) - error += fHEffOn.GetBinError(i); - - TLatex text(0.45, 0.94, Form("T_{eff} = %.1fs \\pm %.1fs", fHEffOn.Integral(), error)); + for (int i=0; iGetXaxis()->GetNbins(); i++) + error += h->GetBinError(i); + + TLatex text(0.45, 0.94, Form("T_{eff} = %.1fs \\pm %.1fs", h->Integral(), error)); text.SetBit(TLatex::kTextNDC); text.SetTextSize(0.04); @@ -386,5 +404,5 @@ AppendPad("fit"); - pad->Divide(2,2); + pad->Divide(2, 1, 0, 0); TH1 *h; @@ -392,5 +410,8 @@ pad->cd(1); gPad->SetBorderMode(0); - h = fHTimeDiff.ProjectionX("TimeDiff", -1, 9999, "E"); + gPad->Divide(1, 2, 0, 0); + pad->GetPad(1)->cd(1); + gPad->SetBorderMode(0); + h = fHDeltaT.ProjectionX(fNameProjDeltaT, -1, 9999, "E"); h->SetTitle("Distribution of \\Delta t [s]"); h->SetXTitle("\\Delta t [s]"); @@ -401,7 +422,16 @@ h->Draw(); + pad->GetPad(1)->cd(2); + gPad->SetBorderMode(0); + fHEffOnTime.Draw(); + AppendPad("time"); + pad->cd(2); gPad->SetBorderMode(0); - h = fHTimeDiff.ProjectionY("Theta", -1, 9999, "E"); + gPad->Divide(1, 3, 0, 0); + + pad->GetPad(2)->cd(1); + gPad->SetBorderMode(0); + h = fHDeltaT.ProjectionY(fNameProjTheta, -1, 9999, "E"); h->SetTitle("Distribution of \\Theta [\\circ]"); h->SetXTitle("\\Theta [\\circ]"); @@ -413,17 +443,44 @@ h->Draw(); - pad->cd(3); - gPad->SetBorderMode(0); - fHEffOn.Draw(); - AppendPad("paint"); - - pad->cd(4); + pad->GetPad(2)->cd(2); gPad->SetBorderMode(0); fHProb.Draw(); + + pad->GetPad(2)->cd(3); + gPad->SetBorderMode(0); + fHEffOnTheta.Draw(); + AppendPad("theta"); +} + +void MHEffectiveOnTime::FillTimeBin(Double_t val, Double_t err) +{ + // Get x-axis + TAxis &x = *fHEffOnTime.GetXaxis(); + + // Get number of bins + const Int_t n = x.GetNbins(); + + // Fill last bin with new values + fHEffOnTime.SetBinContent(n, fParam->GetVal()); + fHEffOnTime.SetBinError(n, fParam->GetErr()); + + // Get time range of histogram, get bin width + const Double_t lo = fHEffOnTime.GetXaxis()->GetXmin(); + const Double_t up = fHEffOnTime.GetXaxis()->GetXmax(); + const Double_t w = fHEffOnTime.GetXaxis()->GetBinWidth(1); + + // Enhance histogram by one bin + MBinning bins; + bins.SetEdges(n+2, lo, up+w); + bins.Apply(fHEffOnTime); + + // Transform overflow bin + // fHEffOnTime.SetBinContent(n+2, fHEffOnTime.GetBinContent(n+1)); + // fHEffOnTime.SetBinError(n+2, fHEffOnTime.GetBinError(n+1)); } void MHEffectiveOnTime::Calc() { - TH1D *h = fHTimeDiff.ProjectionX("", -1, 99999, "E"); + TH1D *h = fHDeltaT.ProjectionX("", -1, 99999, "E"); h->SetDirectory(0); @@ -438,9 +495,8 @@ const Double_t error = res[1]; - fParam->SetVal(val-fEffOnTime0, error-fEffOnErr0); + fParam->SetVal(val>fEffOnTime0?val-fEffOnTime0:0, error-fEffOnErr0); fParam->SetReadyToSave(); - testval += fParam->GetVal(); - testerr += fParam->GetErr(); + FillTimeBin(fParam->GetVal(), fParam->GetErr()); fEffOnTime0 = val; @@ -453,6 +509,6 @@ *fLog << Form("T_{eff} = %.1fs \\pm %.1fs", fParam->GetVal(), fParam->GetErr()); - *fLog << Form(" %.1f %.1f %.1f %.1f", - val, testval, error, testerr) << endl; + *fLog << Form(" %.1f %.1f", + val, error) << endl; fTime->AddMilliSeconds(fInterval*1000); @@ -468,5 +524,5 @@ if (!time) { - *fLog << err << "ERROR - MHEffectiveOnTimeTime::Fill without argument or container doesn't inherit from MTime... abort." << endl; + *fLog << err << "ERROR - MHEffectiveOnTime::Fill without argument or container doesn't inherit from MTime... abort." << endl; return kFALSE; } @@ -483,4 +539,8 @@ *fTime = *time; + MBinning bins; + bins.SetEdges(2, fTime->GetAxisTime(), fTime->GetAxisTime()+fInterval); + bins.Apply(fHEffOnTime); + // Make this just a ns before the first event fTime->Minus1ns(); @@ -493,5 +553,5 @@ } - fHTimeDiff.Fill(*time-fLastTime, fPointPos->GetZd(), w); + fHDeltaT.Fill(*time-fLastTime, fPointPos->GetZd(), w); fLastTime = *time; Index: /trunk/MagicSoft/Mars/mjobs/MJCalibrateSignal.cc =================================================================== --- /trunk/MagicSoft/Mars/mjobs/MJCalibrateSignal.cc (revision 4905) +++ /trunk/MagicSoft/Mars/mjobs/MJCalibrateSignal.cc (revision 4906) @@ -59,4 +59,5 @@ #include "MReadReports.h" +#include "MReadMarsFile.h" #include "MGeomApply.h" #include "MPointingPosCalc.h" @@ -180,5 +181,5 @@ } -Bool_t MJCalibrateSignal::ProcessFile(MPedestalCam &pedcam) +Bool_t MJCalibrateSignal::ProcessFile(MPedestalCam &pedcam, Bool_t ismc) { if (!fSequence.IsValid()) @@ -270,15 +271,16 @@ plist.AddToList(&tlist); - //MReadMarsFile read("Events"); - //read.DisableAutoScheme(); - MReadReports read; - read.AddTree("Events", "MTime.", kTRUE); - read.AddTree("Trigger"); - read.AddTree("Camera"); - read.AddTree("Drive"); - read.AddTree("CC"); - read.AddTree("Currents"); - read.AddFiles(iter); - //read.AddFiles(fnamein); + MReadReports readreal; + readreal.AddTree("Events", "MTime.", kTRUE); + readreal.AddTree("Trigger"); + readreal.AddTree("Camera"); + readreal.AddTree("Drive"); + readreal.AddTree("CC"); + readreal.AddTree("Currents"); + readreal.AddFiles(iter); + + MReadMarsFile readmc("Events"); + readmc.DisableAutoScheme(); + readmc.AddFiles(iter); MGeomApply apply; // Only necessary to create geometry @@ -332,5 +334,5 @@ write.AddContainer("MCerPhotEvt", "Events"); write.AddContainer("MPedPhotCam", "Events"); - write.AddContainer("MTime", "Events"); + write.AddContainer("MTime", "Events", kFALSE); write.AddContainer("MRawEvtHeader", "Events"); write.AddContainer("MArrivalTime", "Events"); @@ -383,5 +385,5 @@ // Now setup main tasklist - tlist.AddToList(&read); + tlist.AddToList(ismc ? (MTask*)&readmc : (MTask*)&readreal); tlist.AddToList(&tlist2, "Events"); tlist.AddToList(&pcalc, "Drive"); Index: /trunk/MagicSoft/Mars/mjobs/MJCalibrateSignal.h =================================================================== --- /trunk/MagicSoft/Mars/mjobs/MJCalibrateSignal.h (revision 4905) +++ /trunk/MagicSoft/Mars/mjobs/MJCalibrateSignal.h (revision 4906) @@ -31,5 +31,5 @@ MJCalibrateSignal(const char *name=NULL, const char *title=NULL); - Bool_t ProcessFile(MPedestalCam &cam); + Bool_t ProcessFile(MPedestalCam &cam, Bool_t ismc=kFALSE); ClassDef(MJCalibrateSignal, 0) // Tool to create a pedestal file (MPedestalCam) Index: /trunk/MagicSoft/Mars/mjobs/MJStar.cc =================================================================== --- /trunk/MagicSoft/Mars/mjobs/MJStar.cc (revision 4905) +++ /trunk/MagicSoft/Mars/mjobs/MJStar.cc (revision 4906) @@ -49,4 +49,5 @@ #include "MReadReports.h" +#include "MReadMarsFile.h" #include "MGeomApply.h" #include "MEventRateCalc.h" @@ -100,5 +101,5 @@ } -Bool_t MJStar::ProcessFile() +Bool_t MJStar::ProcessFile(Bool_t ismc) { if (!fSequence.IsValid()) @@ -145,12 +146,16 @@ plist.AddToList(&tlist); - MReadReports read; - read.AddTree("Events", "MTime.", kTRUE); - read.AddTree("Drive"); + MReadReports readreal; + readreal.AddTree("Events", "MTime.", kTRUE); + readreal.AddTree("Drive"); //read.AddTree("Trigger"); //read.AddTree("Camera"); //read.AddTree("CC"); //read.AddTree("Currents"); - read.AddFiles(iter); + readreal.AddFiles(iter); + + MReadMarsFile readmc("Events"); + readmc.DisableAutoScheme(); + readmc.AddFiles(iter); // ------------------ Setup general tasks ---------------- @@ -215,34 +220,44 @@ write.AddContainer("MNewImagePar2", "Events"); write.AddContainer("MImageParTime", "Events"); - write.AddContainer("MTime", "Events"); write.AddContainer("MRawEvtHeader", "Events"); - // Monte Carlo - write.AddContainer("MMcEvt", "Events", kFALSE); - write.AddContainer("MMcTrig", "Events", kFALSE); - // Run Header - write.AddContainer("MRawRunHeader", "RunHeaders"); - write.AddContainer("MBadPixelsCam", "RunHeaders"); - write.AddContainer("MGeomCam", "RunHeaders"); - //write.AddContainer("MObservatory", "RunHeaders"); - // Monte Carlo Headers - write.AddContainer("MMcTrigHeader", "RunHeaders", kFALSE); - write.AddContainer("MMcConfigRunHeader", "RunHeaders", kFALSE); - write.AddContainer("MMcCorsikaRunHeader", "RunHeaders", kFALSE); - // Drive - //write.AddContainer("MSrcPosCam", "Drive"); - write.AddContainer("MPointingPos", "Drive"); - write.AddContainer("MReportDrive", "Drive"); - write.AddContainer("MTimeDrive", "Drive"); - // Effective On Time - write.AddContainer("MEffectiveOnTime", "EffectiveOnTime"); - write.AddContainer("MTimeEffectiveOnTime", "EffectiveOnTime"); + if (ismc) + { + write.AddContainer("MPointingPos", "Events"); + // Monte Carlo + write.AddContainer("MMcEvt", "Events"); + write.AddContainer("MMcTrig", "Events"); + // Monte Carlo Headers + write.AddContainer("MMcTrigHeader", "RunHeaders"); + write.AddContainer("MMcConfigRunHeader", "RunHeaders"); + write.AddContainer("MMcCorsikaRunHeader", "RunHeaders"); + } + else + { + write.AddContainer("MTime", "Events"); + // Run Header + write.AddContainer("MRawRunHeader", "RunHeaders"); + write.AddContainer("MBadPixelsCam", "RunHeaders"); + write.AddContainer("MGeomCam", "RunHeaders"); + //write.AddContainer("MObservatory", "RunHeaders"); + // Drive + //write.AddContainer("MSrcPosCam", "Drive"); + write.AddContainer("MPointingPos", "Drive"); + write.AddContainer("MReportDrive", "Drive"); + write.AddContainer("MTimeDrive", "Drive"); + // Effective On Time + write.AddContainer("MEffectiveOnTime", "EffectiveOnTime"); + write.AddContainer("MTimeEffectiveOnTime", "EffectiveOnTime"); + } MTaskList tlist2; tlist2.AddToList(&apply); - tlist2.AddToList(&rate); - //tlist2.AddToList(&rate10000); - tlist2.AddToList(&fill8); - tlist2.AddToList(&fill9); - tlist2.AddToList(&filla); + if (!ismc) + { + tlist2.AddToList(&rate); + //tlist2.AddToList(&rate10000); + tlist2.AddToList(&fill8); + tlist2.AddToList(&fill9); + tlist2.AddToList(&filla); + } //tlist2.AddToList(&fillb); tlist2.AddToList(&clean); @@ -265,5 +280,5 @@ //tlist3.AddToList(&srcpos); - tlist.AddToList(&read); + tlist.AddToList(ismc ? (MTask*)&readmc : (MTask*)&readreal); tlist.AddToList(&tlist3, "Drive"); tlist.AddToList(&tlist2, "Events"); Index: /trunk/MagicSoft/Mars/mjobs/MJStar.h =================================================================== --- /trunk/MagicSoft/Mars/mjobs/MJStar.h (revision 4905) +++ /trunk/MagicSoft/Mars/mjobs/MJStar.h (revision 4906) @@ -14,5 +14,5 @@ MJStar(const char *name=NULL, const char *title=NULL); - Bool_t ProcessFile(); + Bool_t ProcessFile(Bool_t ismc=kFALSE); ClassDef(MJStar, 0) // Tool to create a pedestal file (MPedestalCam) Index: /trunk/MagicSoft/Mars/star.cc =================================================================== --- /trunk/MagicSoft/Mars/star.cc (revision 4905) +++ /trunk/MagicSoft/Mars/star.cc (revision 4906) @@ -49,4 +49,6 @@ gLog << " --debug-env Debug setting resources from file" << endl << endl; gLog << endl; + gLog << " Input Options:" << endl; + gLog << " -mc You must use this for MC files (PRELIMINARY)" << endl << endl; gLog << " Output options:" << endl; gLog << " -q Quit when job is finished" << endl; @@ -93,4 +95,5 @@ const Bool_t kOverwrite = arg.HasOnlyAndRemove("-f"); const Bool_t kForceExec = arg.HasOnlyAndRemove("-ff"); + const Bool_t kIsMC = arg.HasOnlyAndRemove("-mc"); const TString kInpath = arg.GetStringAndRemove("--ind=", ""); @@ -210,5 +213,5 @@ // job.SetPathIn(kInpath); // not yet needed - if (!job.ProcessFile()) + if (!job.ProcessFile(kIsMC)) { gLog << err << "Calculation of image parameters failed." << endl << endl;