Index: /trunk/MagicSoft/Mars/mjobs/MJCalibrateSignal.cc =================================================================== --- /trunk/MagicSoft/Mars/mjobs/MJCalibrateSignal.cc (revision 6034) +++ /trunk/MagicSoft/Mars/mjobs/MJCalibrateSignal.cc (revision 6035) @@ -34,6 +34,4 @@ // but this is NOT recommended! // -// Resource file entries are case sensitive! -// ///////////////////////////////////////////////////////////////////////////// #include "MJCalibrateSignal.h" @@ -62,4 +60,22 @@ #include "MCalibrationRelTimeCam.h" #include "MCalibrationChargePINDiode.h" +#include "MCalibrationHiLoCam.h" + +#include "MCalibrationChargeCalc.h" +#include "MCalibrationRelTimeCalc.h" + +#include "MCalibrationIntensityChargeCam.h" +#include "MCalibrationIntensityBlindCam.h" +#include "MCalibrationIntensityRelTimeCam.h" +#include "MCalibrationIntensityQECam.h" +#include "MBadPixelsIntensityCam.h" + +#include "MHCalibrationChargeCam.h" +#include "MHCalibrationChargeBlindCam.h" +#include "MHCalibrationChargePINDiode.h" +#include "MHCalibrationRelTimeCam.h" +#include "MHCamera.h" + +#include "MCalibCalcFromPast.h" #include "MReadReports.h" @@ -72,7 +88,10 @@ #include "MMcPedestalCopy.h" #include "MPointingPosCalc.h" +#include "MInterpolatePointingPos.h" #include "MPedCalcFromLoGain.h" #include "MExtractor.h" #include "MExtractTimeAndCharge.h" +#include "MExtractPINDiode.h" +#include "MExtractBlindPixel.h" #include "MFCosmics.h" #include "MTaskEnv.h" @@ -96,5 +115,6 @@ // MJCalibrateSignal::MJCalibrateSignal(const char *name, const char *title) - : fIsDataCheck(kFALSE) + : fIsDataCheck(kFALSE), fIsInterlaced(kTRUE), fIsRelTimesUpdate(kFALSE), + fIsHiLoCalibration(kTRUE) { fName = name ? name : "MJCalibrateSignal"; @@ -102,5 +122,5 @@ } -Bool_t MJCalibrateSignal::WriteResult() +Bool_t MJCalibrateSignal::WriteResult(TObjArray &cont) { if (fPathOut.IsNull()) @@ -115,9 +135,4 @@ TFile file(oname, "RECREATE"); - if (!file.IsOpen()) - { - *fLog << err << "ERROR - Couldn't open file " << oname << " for writing..." << endl; - return kFALSE; - } *fLog << inf << " - MStatusDisplay..." << flush; @@ -129,5 +144,5 @@ *fLog << inf << "ok." << endl; - return kTRUE; + return fIsInterlaced ? WriteContainer(cont) : kTRUE; } @@ -188,9 +203,15 @@ Bool_t MJCalibrateSignal::CheckEnvLocal() { + SetDataCheck(GetEnv("DataCheck", IsDataCheck())); + SetInterlaced(GetEnv("Interlaced", fIsInterlaced)); + SetRelTimesUpdate(GetEnv("RelTimesUpdate", fIsRelTimesUpdate)); + SetHiLoCalibration(GetEnv("HiLoCalibration", fIsHiLoCalibration)); + return kTRUE; } -Bool_t MJCalibrateSignal::ProcessFile(MPedestalCam &pedcamab, MPedestalCam &pedcam2, MPedestalCam &pedcam) +Bool_t MJCalibrateSignal::ProcessFile(MPedestalCam &pedcamab, MPedestalCam &pedcambias, + MPedestalCam &pedcamextr) { if (!fSequence.IsValid()) @@ -230,10 +251,32 @@ // Read File - /* - MCalibrationIntensityChargeCam calcam; - MCalibrationIntensityQECam qecam; - MCalibrationIntensityBlindCam bndcam; - MCalibrationIntensityRelTimeCam tmcam; - */ + MCalibrationIntensityChargeCam ichcam; + MCalibrationIntensityQECam iqecam; + MCalibrationIntensityBlindCam iblcam; + MCalibrationIntensityRelTimeCam itmcam; + MBadPixelsIntensityCam ibdcam; + + TObjArray interlacedcont; + interlacedcont.Add(&ichcam); + interlacedcont.Add(&iqecam); + interlacedcont.Add(&iblcam); + interlacedcont.Add(&itmcam); + interlacedcont.Add(&ibdcam); + + MHCalibrationChargeCam hchacam; + MHCalibrationChargeBlindCam hbndcam; + MHCalibrationChargePINDiode hpndiod; + MHCalibrationRelTimeCam hrelcam; + + hchacam.SetLoGain(kTRUE); + hchacam.SetOscillations(kFALSE); + hrelcam.SetOscillations(kFALSE); + + interlacedcont.Add(&hchacam); + interlacedcont.Add(&hbndcam); + interlacedcont.Add(&hpndiod); + if (fIsRelTimesUpdate) + interlacedcont.Add(&hrelcam); + MCalibrationChargeCam calcam; MCalibrationQECam qecam; @@ -241,8 +284,14 @@ MCalibrationChargePINDiode pind; MCalibrationRelTimeCam tmcam; + MCalibrationHiLoCam hilcam; + + if (fIsHiLoCalibration) + interlacedcont.Add(&hilcam); + MBadPixelsCam badpix; MExtractor *extractor1=0; MExtractor *extractor2=0; + MExtractor *extractor3=0; TString geom; @@ -257,4 +306,6 @@ return kFALSE; + badpix.Print(); + *fLog << all; if (extractor1) @@ -263,4 +314,5 @@ extractor1->Print(); *fLog << endl; + extractor3 = (MExtractor*)extractor1->Clone(); } else @@ -289,4 +341,8 @@ plist.AddToList(&badcam); plist.AddToList(&calibcont); + plist.AddToList(&hilcam); + + if (fIsInterlaced) + plist.AddToList(&interlacedcont); // Setup Tasklist @@ -326,4 +382,8 @@ readreal.AddFiles(iter); + MInterpolatePointingPos pextr; +// pextr.SetDebug(); // take it out + pextr.AddFiles(&iter); + MGeomApply apply; // Only necessary to create geometry if (!geom.IsNull()) @@ -333,8 +393,8 @@ // Make sure that pedcamab has the correct name pedcamab.SetName("MPedestalFundamental"); - pedcam.SetName("MPedestalFromExtractorRndm"); - pedcam2.SetName("MPedestalFromExtractor"); - plist.AddToList(&pedcam); - plist.AddToList(&pedcam2); + pedcamextr.SetName("MPedestalFromExtractorRndm"); + pedcambias.SetName("MPedestalFromExtractor"); + plist.AddToList(&pedcamextr); + plist.AddToList(&pedcambias); plist.AddToList(&pedcamab); @@ -344,5 +404,9 @@ MFTriggerPattern ftp; ftp.RequireCalibration(); - + // ftp.RequireCalibration(kTRUE); + // ftp.DenyCalibration(); + // ftp.DenyPedestal(); + // ftp.DenyPinDiode(); + // ftp.SetInverted(); MContinue conttp(&ftp, "ContTrigPattern"); // --> tlist2 @@ -389,5 +453,8 @@ } if (extractor2) - extractor2->SetPedestals(&pedcamab); + extractor2->SetPedestals(&pedcamab); + + if (extractor3) + extractor3->SetPedestals(&pedcamab); MFCosmics fcosmics; @@ -399,6 +466,9 @@ MTaskEnv taskenv1("ExtractSignal"); MTaskEnv taskenv2("ExtractTime"); + MTaskEnv taskenv3("ExtractInterlaced"); taskenv1.SetDefault(extractor1); taskenv2.SetDefault(extractor2); + taskenv3.SetDefault(extractor3); + MCalibrateData calib; if (filetype==3) // MC file @@ -419,4 +489,39 @@ } + MExtractPINDiode pinext; + MExtractBlindPixel bldext; + + MFTriggerPattern fcalib("CalibFilter"); + MCalibrationChargeCalc chcalc; + MCalibrationRelTimeCalc recalc; + MCalibCalcFromPast pacalc; + + fcalib.DenyCalibration(); + fcalib.SetInverted(); + + pinext.SetPedestals(&pedcamab); + bldext.SetPedestals(&pedcamab); + chcalc.SetPedestals(&pedcamextr); + + pacalc.SetChargeCalc(&chcalc); + if (fIsRelTimesUpdate) + pacalc.SetRelTimeCalc(&recalc); + pacalc.SetCalibrate(&calib); + + // + // Intensity Calibration histogramming + // + MFillH filpin("MHCalibrationChargePINDiode", "MExtractedSignalPINDiode", "FillPINDiode"); + MFillH filbnd("MHCalibrationChargeBlindCam", "MExtractedSignalBlindPixel", "FillBlindCam"); + MFillH filcam("MHCalibrationChargeCam", "MExtractedSignalCam", "FillChargeCam"); + MFillH filtme("MHCalibrationRelTimeCam", "MArrivalTimeCam", "FillRelTime"); + MFillH filhil("MHCalibrationHiLoCam", "MExtractedSignalCam", "FillHiLoRatio"); + + filpin.SetBit(MFillH::kDoNotDisplay); + filbnd.SetBit(MFillH::kDoNotDisplay); + filcam.SetBit(MFillH::kDoNotDisplay); + filtme.SetBit(MFillH::kDoNotDisplay); + filhil.SetBit(MFillH::kDoNotDisplay); + MCalibrateRelTimes caltm; MBadPixelsCalc bpcal; @@ -434,25 +539,32 @@ MHCamEvent evt1(2, "PedRmsFLG", "Pedestal RMS from Lo Gain;;\\sigma_{p} [fadc/sl]"); MHCamEvent evt2(0, "Extra'd", "Extracted Signal;;S [fadc/sl]"); - //MHCamEvent evt3(4, "PedPhot", "Calibrated Pedestal;;P [\\gamma]"); + MHCamEvent evt3(4, "PedPhot", "Calibrated Pedestal;;P [\\gamma]"); MHCamEvent evt4(5, "PedRMS", "Calibrated Pedestal RMS;;\\sigma_{p} [\\gamma]"); MHCamEvent evt5(0, "Interp'd", "Interpolated Signal;;S [\\gamma]"); - MHCamEvent evt6(2, "Unsuitable", "Fraction of unsuitable events per Pixel;;[1]"); + MHCamEvent evt6(2, "Unsuitable", "Unsuitable event ratio;;%"); MHCamEvent evt7(0, "Times", "Arrival Time;;T [slice]"); + // MHCamEvent evt8(0, "HiLoConv", "Ratio Hi-Lo Gain Signal;;Ratio [1]"); + // MHCamEvent evt9(5, "HiLoDiff", "Time Diff. Hi-Lo Gain Signal;;\\Delta T [3.33ns]"); + evt0.EnableVariance(); evt1.EnableVariance(); evt2.EnableVariance(); - //evt3.EnableVariance(); + evt3.EnableVariance(); evt4.EnableVariance(); evt5.EnableVariance(); evt7.EnableVariance(); + // evt8.EnableVariance(); + // evt9.EnableVariance(); MFillH fill0(&evt0, "MPedestalFundamental", "FillPedFLG"); MFillH fill1(&evt1, "MPedestalFromExtractor", "FillPedRmsFLG"); MFillH fill2(&evt2, "MExtractedSignalCam", "FillExtracted"); - //MFillH fill3(&evt3, "MPedPhotFundamental", "FillPedPhot"); + MFillH fill3(&evt3, "MPedPhotFundamental", "FillPedPhot"); MFillH fill4(&evt4, "MPedPhotFromExtractor", "FillPedRMS"); MFillH fill5(&evt5, "MCerPhotEvt", "FillInterpolated"); MFillH fill6(&evt6, "MBadPixelsCam", "FillUnsuitable"); MFillH fill7(&evt7, "MArrivalTime", "FillTimes"); + // MFillH fill8(&evt8, "MCalibrationHiLoCam", "FillRatio"); + // MFillH fill9(&evt9, "MCalibrationHiLoCam", "FillDiff"); // The second rule is for the case reading raw-files! @@ -494,5 +606,5 @@ write.AddContainer("MTimeTrigger", "Trigger", kFALSE); // Slow-Control: Drive - write.AddContainer("MPointingPos", "Drive", kFALSE); + write.AddContainer("MPointingPos", "Events", kFALSE); write.AddContainer("MReportDrive", "Drive", kFALSE); write.AddContainer("MTimeDrive", "Drive", kFALSE); @@ -503,6 +615,6 @@ // Now setup tasklist for events MTaskList tlist2; + tlist2.AddToList(&decode); - tlist2.AddToList(&conttp); tlist2.AddToList(&apply); tlist2.AddToList(&merge); @@ -515,19 +627,51 @@ tlist2.AddToList(&pedlo3); } + if (fIsInterlaced) + { + pinext.SetFilter(&fcalib); + bldext.SetFilter(&fcalib); + taskenv3.SetFilter(&fcalib); + pacalc.SetFilter(&fcalib); + filcam.SetFilter(&fcalib); + filtme.SetFilter(&fcalib); + filbnd.SetFilter(&fcalib); + filpin.SetFilter(&fcalib); + chcalc.SetFilter(&fcalib); + recalc.SetFilter(&fcalib); + + tlist2.AddToList(&fcalib); + tlist2.AddToList(&pinext); + tlist2.AddToList(&bldext); + tlist2.AddToList(&taskenv3); + tlist2.AddToList(&pacalc); + tlist2.AddToList(&filcam); + if (fIsRelTimesUpdate) + tlist2.AddToList(&filtme); + tlist2.AddToList(&filbnd); + tlist2.AddToList(&filpin); + tlist2.AddToList(&chcalc); + if (fIsRelTimesUpdate) + tlist2.AddToList(&recalc); + } + + tlist2.AddToList(&conttp); + if (extractor1) + tlist2.AddToList(&taskenv1); + if (extractor2) + tlist2.AddToList(&taskenv2); tlist2.AddToList(&fill0); tlist2.AddToList(&fill1); - if (extractor1) - tlist2.AddToList(&taskenv1); - if (extractor2) - tlist2.AddToList(&taskenv2); tlist2.AddToList(&contcos); + if (fIsHiLoCalibration) + tlist2.AddToList(&filhil); tlist2.AddToList(&fill2); tlist2.AddToList(&calib); if (extractor2 || extractor1->InheritsFrom("MExtractTimeAndCharge")) tlist2.AddToList(&caltm); + tlist2.AddToList(&bpcal); tlist2.AddToList(&treat); tlist2.AddToList(&fill6); - //tlist2.AddToList(&fill3); + tlist2.AddToList(&fill3); tlist2.AddToList(&fill4); tlist2.AddToList(&fill5); @@ -535,4 +679,9 @@ tlist2.AddToList(&fill7); + if (fIsHiLoCalibration) + { + // tlist2.AddToList(&fill8); + // tlist2.AddToList(&fill9); + } // Setup List for Drive-tree MPointingPosCalc pcalc; @@ -540,7 +689,14 @@ // Now setup main tasklist tlist.AddToList(IsDataCheck() ? (MTask*)&rawread : (MTask*)&readreal); + tlist2.AddToList(&pextr); tlist.AddToList(&tlist2, IsDataCheck()?"All":"Events"); - if (!IsDataCheck()) - tlist.AddToList(&pcalc, "Drive"); + + // + // !!! Replaced by MInterpolatePointingPos !!! + // + // if (!IsDataCheck()) + // tlist.AddToList(&pcalc, "Drive"); + + tlist.AddToList(&write); @@ -562,5 +718,7 @@ tlist.PrintStatistics(); - if (!WriteResult()) + DisplayResult(plist); + + if (!WriteResult(interlacedcont)) return kFALSE; @@ -570,2 +728,39 @@ return kTRUE; } + + +void MJCalibrateSignal::DisplayResult(MParList &plist) +{ + + if (!fDisplay) + return; + + if (!fIsHiLoCalibration) + return; + + // + // Get container from list + // + MGeomCam &geomcam = *(MGeomCam*)plist.FindObject("MGeomCam"); + MCalibrationHiLoCam *cam = (MCalibrationHiLoCam*)plist.FindObject("MCalibrationHiLoCam"); + + // Create histograms to display + MHCamera disp1 (geomcam, "HiLoConv", "Ratio Amplification HiGain vs. LoGain (Charges)"); + MHCamera disp2 (geomcam, "HiLoDiff", "Arrival Time Diff. HiGain vs. LoGain (Times)"); + + disp1.SetCamContent(*cam, 0); + disp1.SetCamError( *cam, 1); + disp2.SetCamContent(*cam, 5); + disp2.SetCamError( *cam, 6); + + disp1.SetYTitle("R [1]"); + disp2.SetYTitle("\\Delta T [FADC sl.]"); + + + TCanvas &c1 = fDisplay->AddTab("HiLoConv"); + c1.Divide(2,3); + + disp1.CamDraw(c1, 1, 2, 1); + disp2.CamDraw(c1, 2, 2, 1); +} + Index: /trunk/MagicSoft/Mars/mjobs/MJCalibrateSignal.h =================================================================== --- /trunk/MagicSoft/Mars/mjobs/MJCalibrateSignal.h (revision 6034) +++ /trunk/MagicSoft/Mars/mjobs/MJCalibrateSignal.h (revision 6035) @@ -20,9 +20,15 @@ { private: + Bool_t fIsDataCheck; - + Bool_t fIsInterlaced; // Distinguish interlaced from other calibration + Bool_t fIsRelTimesUpdate; // Choose to update relative times from interlaced + Bool_t fIsHiLoCalibration; // Choose to calibrate the high-gain vs. low-gains + Bool_t CheckEnvLocal(); - Bool_t WriteResult(); + void DisplayResult(MParList &plist); + + Bool_t WriteResult(TObjArray &cont); Bool_t ReadCalibration(TObjArray &o, MBadPixelsCam &bpix, MExtractor* &ext1, MExtractor* &ext2, TString &geom) const; @@ -36,4 +42,8 @@ Bool_t ProcessFile(MPedestalCam &camab, MPedestalCam &cam1, MPedestalCam &cam2); + void SetInterlaced ( const Bool_t b=kTRUE ) { fIsInterlaced = b; } + void SetRelTimesUpdate ( const Bool_t b=kTRUE ) { fIsRelTimesUpdate = b; } + void SetHiLoCalibration( const Bool_t b=kTRUE ) { fIsHiLoCalibration = b; } + ClassDef(MJCalibrateSignal, 0) // Tool to create a pedestal file (MPedestalCam) }; Index: /trunk/MagicSoft/Mars/mpointing/MInterpolatePointingPos.cc =================================================================== --- /trunk/MagicSoft/Mars/mpointing/MInterpolatePointingPos.cc (revision 6035) +++ /trunk/MagicSoft/Mars/mpointing/MInterpolatePointingPos.cc (revision 6035) @@ -0,0 +1,408 @@ +/* ======================================================================== *\ +! +! * +! * 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): Marcos Lopez 03/2004 +! +! Copyright: MAGIC Software Development, 2000-2004 +! +! +\* ======================================================================== */ + +///////////////////////////////////////////////////////////////////////////// +// +// MExtrapolatePointingPos +// +// In the PreProcess, read the drive report and store it in an TSpline. +// In the Process, use the TSpline to calculate the PointingPos for the +// time of the current event. +// +// Input Containers: +// MRawEvtData +// +// Output Containers: +// MPointingPos +// +// +///////////////////////////////////////////////////////////////////////////// + +#include "MInterpolatePointingPos.h" + +#include "MLog.h" +#include "MLogManip.h" + +#include + +#include "MTaskList.h" +#include "MParList.h" +#include "MEvtLoop.h" +#include "MReadReports.h" +#include "MReportDrive.h" +#include "MPointingPos.h" +#include "MTime.h" +#include "MRawRunHeader.h" +#include "MDirIter.h" + +#include + +ClassImp(MInterpolatePointingPos); + +using namespace std; + +const Int_t MInterpolatePointingPos::fgNumStartEvents = 1; +// -------------------------------------------------------------------------- +// +// Constructor +// +MInterpolatePointingPos::MInterpolatePointingPos(const char *name, const char *title) + : fEvtTime(NULL), fPointingPos(NULL), fRunHeader(NULL), fDirIter(NULL), + fSplineZd(NULL), fSplineAz(NULL), fSplineRa(NULL), fSplineDec(NULL), + fTimeMode(MInterpolatePointingPos::kEventTime) + +{ + fName = name ? name : "MInterpolatePointingPos"; + fTitle = title ? title : "Task to interpolate the pointing positons from drive reports"; + + fFilename = ""; + fDebug = kFALSE; + + SetNumStartEvents(); +} + + +MInterpolatePointingPos::~MInterpolatePointingPos() +{ + Clear(); +} + +void MInterpolatePointingPos::Clear(Option_t *o) +{ + if(fSplineZd) + delete fSplineZd; + if(fSplineAz) + delete fSplineAz; + if(fSplineRa) + delete fSplineRa; + if(fSplineDec) + delete fSplineDec; +} + +// --------------------------------------------------------------------------- +// +// Read the drive report file for the whole night, a build from it the splines +// +Bool_t MInterpolatePointingPos::ReadDriveReport() +{ + + *fLog << endl << "["<< GetName() + << "]: Loading report file \"" << fFilename << "\" into TSpline..." << endl; + + // + // ParList + // + MParList plist; + MTaskList tlist; + plist.AddToList(&tlist); + + + // + // TaskList + // + MReadReports read; + read.AddTree("Drive"); + if (fDirIter) + read.AddFiles(*fDirIter); + else + read.AddFile(fFilename); // after the reading of the trees!!! + + read.AddToBranchList("MReportDrive.*"); + + tlist.AddToList(&read); + + // + // EventLoop + // + MEvtLoop evtloop; + evtloop.SetParList(&plist); + + if (!evtloop.PreProcess()) + return kFALSE; + + TArrayD reportTime(fNumStartEvents); + TArrayD currentZd(fNumStartEvents); + TArrayD currentAz(fNumStartEvents); + TArrayD nominalZd(fNumStartEvents); + TArrayD nominalAz(fNumStartEvents); + TArrayD ra(fNumStartEvents); + TArrayD dec(fNumStartEvents); + + Int_t n=1; + while (tlist.Process()) + { + MReportDrive* report = (MReportDrive*)plist.FindObject("MReportDrive"); + MTime* reporttime = (MTime*)plist.FindObject("MTimeDrive"); + + if(n==1) + fFirstDriveTime = *reporttime; + else + fLastDriveTime = *reporttime; + + // + // Update the number of entries + // + if (n>fNumStartEvents) + { + reportTime.Set(n); + currentZd.Set(n); + currentAz.Set(n); + nominalZd.Set(n); + nominalAz.Set(n); + ra.Set(n); + dec.Set(n); + } + // + // Sometimes there are two reports with the same time + // + if (reporttime->GetTime() == reportTime[n-1]) + { + cout << warn <<"["<< GetName() + << "]: Warning: this report has the same time that the previous one...skipping it " << endl; + continue; + } + + reportTime[n-1] = reporttime->GetTime(); + currentZd [n-1] = report->GetCurrentZd(); + currentAz [n-1] = report->GetCurrentAz(); + nominalZd [n-1] = report->GetNominalZd(); + nominalAz [n-1] = report->GetNominalAz(); + ra [n-1] = report->GetRa(); + dec [n-1] = report->GetDec(); + + if (fDebug) + { + *fLog << " GetTime(): " << reporttime->GetTime() << endl; + *fLog << " GetCurrentZd(): " << report->GetCurrentZd() << endl; + } + n++; + } + + tlist.PrintStatistics(); + + *fLog << "["<< GetName() << "]: loaded " << n-1 << " ReportDrive from " + << fFirstDriveTime << " to " << fLastDriveTime << endl << endl; + + if (fDebug) + { + for (int i=0;iDivide(2,2); + c->cd(1); + fSplineZd->Draw(); + c->cd(2); + fSplineAz->Draw(); + c->cd(3); + fSplineRa->Draw(); + c->cd(4); + fSplineDec->Draw(); + c->Modified(); + c->Update(); + c->SaveAs("pointing.root"); + } + return kTRUE; +} + + +// -------------------------------------------------------------------------- +// +// Input: +// - MTime +// +// Output: +// - MPointingPos +// +Int_t MInterpolatePointingPos::PreProcess( MParList *pList ) +{ + + Clear(); + + fRunHeader = (MRawRunHeader*)pList->FindObject(AddSerialNumber("MRawRunHeader")); + if (!fRunHeader) + { + *fLog << err << "MRunHeader not found... aborting." << endl; + return kFALSE; + } + + fEvtTime = (MTime*)pList->FindObject("MTime"); + if (!fEvtTime) + { + *fLog << err << "MTime not found... aborting." << endl; + return kFALSE; + } + + fPointingPos = (MPointingPos*)pList->FindCreateObj("MPointingPos"); + if (!fPointingPos) + return kFALSE; + + if (fFilename.IsNull() && !fDirIter) + { + *fLog << err << "File name is empty or no MDirIter has been handed over... aborting" << endl; + return kFALSE; + } + + if( !ReadDriveReport() ) + return kFALSE; + + return kTRUE; +} + + +// -------------------------------------------------------------------------- +// +// Get the run start time, and get the pointing position for that time +// +Int_t MInterpolatePointingPos::Process() +{ + + //const Int_t run = fRunHeader->GetRunNumber(); + const MTime* StartRunTime = &fRunHeader->GetRunStart(); + Int_t time = StartRunTime->GetTime(); + + switch(fTimeMode) + { + case kRunTime: + time = StartRunTime->GetTime(); + + // + // Check that we have drive report for this time + // + if( *StartRunTimefLastDriveTime) + { + *fLog << err << GetDescriptor() << ": Run time " << *StartRunTime + << " outside range of drive reports (" << fFirstDriveTime + << ", " << fLastDriveTime << ")" << endl; + + if ( *fEvtTimefLastDriveTime ) time = fLastDriveTime.GetTime(); + + *fLog << " PointingPos: time = " << time << " (" << *fEvtTime << ") (zd, az, ra, dec) = (" << fSplineZd->Eval( time )<< ", " <Eval( time )<< ", " <Eval( time ) << ", " <Eval( time )<< ")" << endl; + + //fError = kTRUE; + //return kFALSE; + } + break; + + case kEventTime: + time = fEvtTime->GetTime(); + + if (fDebug) + { + cout << " real time : " << time + << " first time: " << fFirstDriveTime.GetTime() + << " last time: " << fLastDriveTime.GetTime() << endl; + } + // + // Check that we have drive report for this time + // + if( *fEvtTimefLastDriveTime) + { + if (fDebug) + { + *fLog << err << GetDescriptor() << ": Run time = " + << *fEvtTime << " outside range of drive reports (" + << fFirstDriveTime << ", "<< fLastDriveTime << ")" << endl; + + if ( *fEvtTimefLastDriveTime ) time = fLastDriveTime.GetTime(); + + *fLog << " PointingPos: time = " << time << " (" << *fEvtTime << ") (zd, az, ra, dec) = (" << fSplineZd->Eval( time )<< ", " <Eval( time )<< ", " <Eval( time ) << ", " <Eval( time )<< ")" << endl; + } + + //fError = kTRUE; + //return kFALSE; + //return kCONTINUE; + } + break; + } + + + + + // + // Check that we have drive report for this time + // + //if( *StartRunTimefLastDriveTime) + if( *StartRunTime>fLastDriveTime) + { + *fLog << err << GetDescriptor() << ": Run time " << *StartRunTime + << " outside range of drive reports (" << fFirstDriveTime + << ", " << fLastDriveTime << ")" << endl; + return kERROR; + } + + //if(run < 20000) + // time = fEvtTime->GetTime(); + //else + // time = StartRunTime->GetTime(); + + + const Double_t zd = fSplineZd->Eval( time ); + const Double_t az = fSplineAz->Eval( time ); + const Double_t ra = fSplineRa->Eval( time ); + const Double_t dec = fSplineDec->Eval( time ); + + + if(TMath::Abs(zd)>90 || TMath::Abs(az)>360 || ra>24 || ra<0 || TMath::Abs(dec)>90){ + *fLog << err << GetDescriptor() << ": Wrong Interpolated Pointing Position." << endl; + *fLog << " PointingPos: time = " << time << " (" << *fEvtTime << ") (zd, az, ra, dec) = (" << zd << ", " << az << ", " << ra << ", " << dec << ")" << endl; + } + + fPointingPos->SetLocalPosition( zd, az ); + fPointingPos->SetSkyPosition( ra*TMath::DegToRad()/15, dec*TMath::DegToRad()); + + return kTRUE; +} + + + +Int_t MInterpolatePointingPos::ReadEnv(const TEnv &env, TString prefix, Bool_t print) +{ + Bool_t rc = kFALSE; + + if (IsEnvDefined(env, prefix, "NumStartEvents", print)) + { + SetNumStartEvents(GetEnvValue(env, prefix, "NumStartEvents", fNumStartEvents)); + rc = kTRUE; + } + return rc; +} Index: /trunk/MagicSoft/Mars/mpointing/MInterpolatePointingPos.h =================================================================== --- /trunk/MagicSoft/Mars/mpointing/MInterpolatePointingPos.h (revision 6035) +++ /trunk/MagicSoft/Mars/mpointing/MInterpolatePointingPos.h (revision 6035) @@ -0,0 +1,85 @@ +#ifndef MARS_MInterpolatePointingPos +#define MARS_MInterpolatePointingPos + +#ifndef MARS_MTask +#include "MTask.h" +#endif + +#ifndef ROOT_TSpline +#include +#endif + +#ifndef MARS_MTime +#include "MTime.h" +#endif + +class MTime; +class MPointingPos; +class MRawRunHeader; +class MDirIter; + +class MInterpolatePointingPos : public MTask +{ +public: + + enum TimeMode_t { + kRunTime, + kEventTime + }; + +private: + + static const Int_t fgNumStartEvents; //! Default for fNumStartEvents (now set to: 10000) + + Int_t fNumStartEvents; // Start number of allowed events + + Bool_t fDebug; + + TString fFilename; + MTime fFirstDriveTime; + MTime fLastDriveTime; + + MTime *fEvtTime; //! Raw event time + MPointingPos *fPointingPos; //! Telescope pointing postion + MRawRunHeader *fRunHeader; //! Run Header + MDirIter *fDirIter; //! Dir Iter + + TSpline3* fSplineZd; //! Zd vs. time + TSpline3* fSplineAz; //! Az vs. time + TSpline3* fSplineRa; //! Ra vs. time + TSpline3* fSplineDec; //! Dec vs. time + + Int_t PreProcess(MParList *pList); + Int_t Process(); + Bool_t ReadDriveReport(); + + TimeMode_t fTimeMode; + + + Int_t ReadEnv(const TEnv &env, TString prefix, Bool_t print); + +public: + + MInterpolatePointingPos(const char *name=NULL, const char *title=NULL); + + ~MInterpolatePointingPos(); + + void AddFiles(MDirIter *dir) { fDirIter = dir; } + void AddFile(const char *name) { fFilename = name; } + + void SetTimeMode( TimeMode_t mode) { fTimeMode = mode; } + void SetDebug( const Bool_t b=kTRUE) { fDebug = b; } + + void Clear(Option_t *o=""); + + Int_t GetNumStartEvents() const { return fNumStartEvents; } + void SetNumStartEvents ( const Int_t i=fgNumStartEvents ) { fNumStartEvents = i; } + + ClassDef(MInterpolatePointingPos, 1) // Interpolate the drive pointing positions +}; + +#endif + + + + Index: unk/MagicSoft/Mars/mtemp/mmpi/MInterpolatePointingPos.cc =================================================================== --- /trunk/MagicSoft/Mars/mtemp/mmpi/MInterpolatePointingPos.cc (revision 6034) +++ (revision ) @@ -1,408 +1,0 @@ -/* ======================================================================== *\ -! -! * -! * 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): Marcos Lopez 03/2004 -! -! Copyright: MAGIC Software Development, 2000-2004 -! -! -\* ======================================================================== */ - -///////////////////////////////////////////////////////////////////////////// -// -// MExtrapolatePointingPos -// -// In the PreProcess, read the drive report and store it in an TSpline. -// In the Process, use the TSpline to calculate the PointingPos for the -// time of the current event. -// -// Input Containers: -// MRawEvtData -// -// Output Containers: -// MPointingPos -// -// -///////////////////////////////////////////////////////////////////////////// - -#include "MInterpolatePointingPos.h" - -#include "MLog.h" -#include "MLogManip.h" - -#include - -#include "MTaskList.h" -#include "MParList.h" -#include "MEvtLoop.h" -#include "MReadReports.h" -#include "MReportDrive.h" -#include "MPointingPos.h" -#include "MTime.h" -#include "MRawRunHeader.h" -#include "MDirIter.h" - -#include - -ClassImp(MInterpolatePointingPos); - -using namespace std; - -const Int_t MInterpolatePointingPos::fgNumStartEvents = 1; -// -------------------------------------------------------------------------- -// -// Constructor -// -MInterpolatePointingPos::MInterpolatePointingPos(const char *name, const char *title) - : fEvtTime(NULL), fPointingPos(NULL), fRunHeader(NULL), fDirIter(NULL), - fSplineZd(NULL), fSplineAz(NULL), fSplineRa(NULL), fSplineDec(NULL), - fTimeMode(MInterpolatePointingPos::kEventTime) - -{ - fName = name ? name : "MInterpolatePointingPos"; - fTitle = title ? title : "Task to interpolate the pointing positons from drive reports"; - - fFilename = ""; - fDebug = kFALSE; - - SetNumStartEvents(); -} - - -MInterpolatePointingPos::~MInterpolatePointingPos() -{ - Clear(); -} - -void MInterpolatePointingPos::Clear(Option_t *o) -{ - if(fSplineZd) - delete fSplineZd; - if(fSplineAz) - delete fSplineAz; - if(fSplineRa) - delete fSplineRa; - if(fSplineDec) - delete fSplineDec; -} - -// --------------------------------------------------------------------------- -// -// Read the drive report file for the whole night, a build from it the splines -// -Bool_t MInterpolatePointingPos::ReadDriveReport() -{ - - *fLog << endl << "["<< GetName() - << "]: Loading report file \"" << fFilename << "\" into TSpline..." << endl; - - // - // ParList - // - MParList plist; - MTaskList tlist; - plist.AddToList(&tlist); - - - // - // TaskList - // - MReadReports read; - read.AddTree("Drive"); - if (fDirIter) - read.AddFiles(*fDirIter); - else - read.AddFile(fFilename); // after the reading of the trees!!! - - read.AddToBranchList("MReportDrive.*"); - - tlist.AddToList(&read); - - // - // EventLoop - // - MEvtLoop evtloop; - evtloop.SetParList(&plist); - - if (!evtloop.PreProcess()) - return kFALSE; - - TArrayD reportTime(fNumStartEvents); - TArrayD currentZd(fNumStartEvents); - TArrayD currentAz(fNumStartEvents); - TArrayD nominalZd(fNumStartEvents); - TArrayD nominalAz(fNumStartEvents); - TArrayD ra(fNumStartEvents); - TArrayD dec(fNumStartEvents); - - Int_t n=1; - while (tlist.Process()) - { - MReportDrive* report = (MReportDrive*)plist.FindObject("MReportDrive"); - MTime* reporttime = (MTime*)plist.FindObject("MTimeDrive"); - - if(n==1) - fFirstDriveTime = *reporttime; - else - fLastDriveTime = *reporttime; - - // - // Update the number of entries - // - if (n>fNumStartEvents) - { - reportTime.Set(n); - currentZd.Set(n); - currentAz.Set(n); - nominalZd.Set(n); - nominalAz.Set(n); - ra.Set(n); - dec.Set(n); - } - // - // Sometimes there are two reports with the same time - // - if (reporttime->GetTime() == reportTime[n-1]) - { - cout << warn <<"["<< GetName() - << "]: Warning: this report has the same time that the previous one...skipping it " << endl; - continue; - } - - reportTime[n-1] = reporttime->GetTime(); - currentZd [n-1] = report->GetCurrentZd(); - currentAz [n-1] = report->GetCurrentAz(); - nominalZd [n-1] = report->GetNominalZd(); - nominalAz [n-1] = report->GetNominalAz(); - ra [n-1] = report->GetRa(); - dec [n-1] = report->GetDec(); - - if (fDebug) - { - *fLog << " GetTime(): " << reporttime->GetTime() << endl; - *fLog << " GetCurrentZd(): " << report->GetCurrentZd() << endl; - } - n++; - } - - tlist.PrintStatistics(); - - *fLog << "["<< GetName() << "]: loaded " << n-1 << " ReportDrive from " - << fFirstDriveTime << " to " << fLastDriveTime << endl << endl; - - if (fDebug) - { - for (int i=0;iDivide(2,2); - c->cd(1); - fSplineZd->Draw(); - c->cd(2); - fSplineAz->Draw(); - c->cd(3); - fSplineRa->Draw(); - c->cd(4); - fSplineDec->Draw(); - c->Modified(); - c->Update(); - c->SaveAs("pointing.root"); - } - return kTRUE; -} - - -// -------------------------------------------------------------------------- -// -// Input: -// - MTime -// -// Output: -// - MPointingPos -// -Int_t MInterpolatePointingPos::PreProcess( MParList *pList ) -{ - - Clear(); - - fRunHeader = (MRawRunHeader*)pList->FindObject(AddSerialNumber("MRawRunHeader")); - if (!fRunHeader) - { - *fLog << err << "MRunHeader not found... aborting." << endl; - return kFALSE; - } - - fEvtTime = (MTime*)pList->FindObject("MTime"); - if (!fEvtTime) - { - *fLog << err << "MTime not found... aborting." << endl; - return kFALSE; - } - - fPointingPos = (MPointingPos*)pList->FindCreateObj("MPointingPos"); - if (!fPointingPos) - return kFALSE; - - if (fFilename.IsNull() && !fDirIter) - { - *fLog << err << "File name is empty or no MDirIter has been handed over... aborting" << endl; - return kFALSE; - } - - if( !ReadDriveReport() ) - return kFALSE; - - return kTRUE; -} - - -// -------------------------------------------------------------------------- -// -// Get the run start time, and get the pointing position for that time -// -Int_t MInterpolatePointingPos::Process() -{ - - //const Int_t run = fRunHeader->GetRunNumber(); - const MTime* StartRunTime = &fRunHeader->GetRunStart(); - Int_t time = StartRunTime->GetTime(); - - switch(fTimeMode) - { - case kRunTime: - time = StartRunTime->GetTime(); - - // - // Check that we have drive report for this time - // - if( *StartRunTimefLastDriveTime) - { - *fLog << err << GetDescriptor() << ": Run time " << *StartRunTime - << " outside range of drive reports (" << fFirstDriveTime - << ", " << fLastDriveTime << ")" << endl; - - if ( *fEvtTimefLastDriveTime ) time = fLastDriveTime.GetTime(); - - *fLog << " PointingPos: time = " << time << " (" << *fEvtTime << ") (zd, az, ra, dec) = (" << fSplineZd->Eval( time )<< ", " <Eval( time )<< ", " <Eval( time ) << ", " <Eval( time )<< ")" << endl; - - //fError = kTRUE; - //return kFALSE; - } - break; - - case kEventTime: - time = fEvtTime->GetTime(); - - if (fDebug) - { - cout << " real time : " << time - << " first time: " << fFirstDriveTime.GetTime() - << " last time: " << fLastDriveTime.GetTime() << endl; - } - // - // Check that we have drive report for this time - // - if( *fEvtTimefLastDriveTime) - { - if (fDebug) - { - *fLog << err << GetDescriptor() << ": Run time = " - << *fEvtTime << " outside range of drive reports (" - << fFirstDriveTime << ", "<< fLastDriveTime << ")" << endl; - - if ( *fEvtTimefLastDriveTime ) time = fLastDriveTime.GetTime(); - - *fLog << " PointingPos: time = " << time << " (" << *fEvtTime << ") (zd, az, ra, dec) = (" << fSplineZd->Eval( time )<< ", " <Eval( time )<< ", " <Eval( time ) << ", " <Eval( time )<< ")" << endl; - } - - //fError = kTRUE; - //return kFALSE; - //return kCONTINUE; - } - break; - } - - - - - // - // Check that we have drive report for this time - // - //if( *StartRunTimefLastDriveTime) - if( *StartRunTime>fLastDriveTime) - { - *fLog << err << GetDescriptor() << ": Run time " << *StartRunTime - << " outside range of drive reports (" << fFirstDriveTime - << ", " << fLastDriveTime << ")" << endl; - return kERROR; - } - - //if(run < 20000) - // time = fEvtTime->GetTime(); - //else - // time = StartRunTime->GetTime(); - - - const Double_t zd = fSplineZd->Eval( time ); - const Double_t az = fSplineAz->Eval( time ); - const Double_t ra = fSplineRa->Eval( time ); - const Double_t dec = fSplineDec->Eval( time ); - - - if(TMath::Abs(zd)>90 || TMath::Abs(az)>360 || ra>24 || ra<0 || TMath::Abs(dec)>90){ - *fLog << err << GetDescriptor() << ": Wrong Interpolated Pointing Position." << endl; - *fLog << " PointingPos: time = " << time << " (" << *fEvtTime << ") (zd, az, ra, dec) = (" << zd << ", " << az << ", " << ra << ", " << dec << ")" << endl; - } - - fPointingPos->SetLocalPosition( zd, az ); - fPointingPos->SetSkyPosition( ra*TMath::DegToRad()/15, dec*TMath::DegToRad()); - - return kTRUE; -} - - - -Int_t MInterpolatePointingPos::ReadEnv(const TEnv &env, TString prefix, Bool_t print) -{ - Bool_t rc = kFALSE; - - if (IsEnvDefined(env, prefix, "NumStartEvents", print)) - { - SetNumStartEvents(GetEnvValue(env, prefix, "NumStartEvents", fNumStartEvents)); - rc = kTRUE; - } - return rc; -} Index: unk/MagicSoft/Mars/mtemp/mmpi/MInterpolatePointingPos.h =================================================================== --- /trunk/MagicSoft/Mars/mtemp/mmpi/MInterpolatePointingPos.h (revision 6034) +++ (revision ) @@ -1,85 +1,0 @@ -#ifndef MARS_MInterpolatePointingPos -#define MARS_MInterpolatePointingPos - -#ifndef MARS_MTask -#include "MTask.h" -#endif - -#ifndef ROOT_TSpline -#include -#endif - -#ifndef MARS_MTime -#include "MTime.h" -#endif - -class MTime; -class MPointingPos; -class MRawRunHeader; -class MDirIter; - -class MInterpolatePointingPos : public MTask -{ -public: - - enum TimeMode_t { - kRunTime, - kEventTime - }; - -private: - - static const Int_t fgNumStartEvents; //! Default for fNumStartEvents (now set to: 10000) - - Int_t fNumStartEvents; // Start number of allowed events - - Bool_t fDebug; - - TString fFilename; - MTime fFirstDriveTime; - MTime fLastDriveTime; - - MTime *fEvtTime; //! Raw event time - MPointingPos *fPointingPos; //! Telescope pointing postion - MRawRunHeader *fRunHeader; //! Run Header - MDirIter *fDirIter; //! Dir Iter - - TSpline3* fSplineZd; //! Zd vs. time - TSpline3* fSplineAz; //! Az vs. time - TSpline3* fSplineRa; //! Ra vs. time - TSpline3* fSplineDec; //! Dec vs. time - - Int_t PreProcess(MParList *pList); - Int_t Process(); - Bool_t ReadDriveReport(); - - TimeMode_t fTimeMode; - - - Int_t ReadEnv(const TEnv &env, TString prefix, Bool_t print); - -public: - - MInterpolatePointingPos(const char *name=NULL, const char *title=NULL); - - ~MInterpolatePointingPos(); - - void AddFiles(MDirIter *dir) { fDirIter = dir; } - void AddFile(const char *name) { fFilename = name; } - - void SetTimeMode( TimeMode_t mode) { fTimeMode = mode; } - void SetDebug( const Bool_t b=kTRUE) { fDebug = b; } - - void Clear(Option_t *o=""); - - Int_t GetNumStartEvents() const { return fNumStartEvents; } - void SetNumStartEvents ( const Int_t i=fgNumStartEvents ) { fNumStartEvents = i; } - - ClassDef(MInterpolatePointingPos, 1) // Interpolate the drive pointing positions -}; - -#endif - - - -