Index: trunk/MagicSoft/Mars/manalysis/AnalysisLinkDef.h =================================================================== --- trunk/MagicSoft/Mars/manalysis/AnalysisLinkDef.h (revision 3100) +++ trunk/MagicSoft/Mars/manalysis/AnalysisLinkDef.h (revision 3101) @@ -76,4 +76,5 @@ #pragma link C++ class MArrivalTime+; #pragma link C++ class MArrivalTimeCalc+; +#pragma link C++ class MArrivalTimeCalc2+; #pragma link C++ class MMcCalibrationUpdate+; Index: trunk/MagicSoft/Mars/manalysis/MArrivalTimeCalc2.cc =================================================================== --- trunk/MagicSoft/Mars/manalysis/MArrivalTimeCalc2.cc (revision 3101) +++ trunk/MagicSoft/Mars/manalysis/MArrivalTimeCalc2.cc (revision 3101) @@ -0,0 +1,279 @@ +/* ======================================================================== *\ +! +! * +! * 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, 02/2004 +! Hendrik Bartko, 02/2004 +! +! Copyright: MAGIC Software Development, 2000-2004 +! +! +\* ======================================================================== */ + +////////////////////////////////////////////////////////////////////////////// +// +// MArrivalTimeCalc2 +// +// Calculates the arrival time as the mean time of the fWindowSize time slices +// which have the highest integral content. +// +// +////////////////////////////////////////////////////////////////////////////// +#include "MArrivalTimeCalc2.h" + +#include "MLog.h" +#include "MLogManip.h" + +#include "MParList.h" + +#include "MRawEvtData.h" +#include "MRawEvtPixelIter.h" + +#include "MPedestalCam.h" +#include "MPedestalPix.h" + +#include "MArrivalTime.h" + + + +ClassImp(MArrivalTimeCalc2); + +using namespace std; + +const Byte_t MArrivalTimeCalc2::fgSaturationLimit = 254; +const Byte_t MArrivalTimeCalc2::fgFirst = 0; +const Byte_t MArrivalTimeCalc2::fgLast = 14; +const Byte_t MArrivalTimeCalc2::fgWindowSize = 6; + +// -------------------------------------------------------------------------- +// +// Default constructor. +// +MArrivalTimeCalc2::MArrivalTimeCalc2(const char *name, const char *title) + : fSaturationLimit(fgSaturationLimit) +{ + + fName = name ? name : "MArrivalTimeCalc2"; + fTitle = title ? title : "Task to extract the signal from the FADC slices"; + + AddToBranchList("MRawEvtData.*"); + + SetRange(); +} + +void MArrivalTimeCalc2::SetRange(Byte_t hifirst, Byte_t hilast, Byte_t lofirst, Byte_t lolast, Byte_t windowsize) +{ + + fNumHiGainSamples = hilast-hifirst+1; + fNumLoGainSamples = lolast-lofirst+1; + + fHiGainFirst = hifirst; + fLoGainFirst = lofirst; + + + fWindowSize = windowsize & ~1; + + if (fWindowSize != windowsize) + *fLog << warn << "MArrivalTimeCalc2::SetRange - window size has to be even, set to: " << int(fWindowSize) << " samples " << endl; + + if (fWindowSize<2) + { + fWindowSize = 2; + *fLog << warn << "MArrivalTimeCalc2::SetRange - window size set to two samples" << endl; + } + + if (fWindowSize > fNumHiGainSamples) + { + fWindowSize = fNumLoGainSamples & ~1; + *fLog << warn << "MArrivalTimeCalc2::SetRange - window size set to " << int(fWindowSize) << " samples " << endl; + } + + if (fWindowSize > fNumLoGainSamples) + { + fWindowSize = fNumLoGainSamples & ~1; + *fLog << warn << "MArrivalTimeCalc2::SetRange - window size set to " << int(fWindowSize) << " samples " << endl; + } + + fWindowSizeSqrt = TMath::Sqrt((Float_t)fWindowSize); + +} + + + +// -------------------------------------------------------------------------- +// +// The PreProcess searches for the following input containers: +// - MRawEvtData +// - MPedestalCam +// +// The following output containers are also searched and created if +// they were not found: +// +// - MExtractedSignalCam +// +Int_t MArrivalTimeCalc2::PreProcess(MParList *pList) +{ + fRawEvt = (MRawEvtData*)pList->FindObject(AddSerialNumber("MRawEvtData")); + if (!fRawEvt) + { + *fLog << err << AddSerialNumber("MRawEvtData") << " not found... aborting." << endl; + return kFALSE; + } + + + + + fPedestals = (MPedestalCam*)pList->FindObject(AddSerialNumber("MPedestalCam")); + if (!fPedestals) + { + *fLog << err << AddSerialNumber("MPedestalCam") << " not found... aborting" << endl; + return kFALSE; + } + + fArrivalTime = (MArrivalTime*)pList->FindCreateObj(AddSerialNumber("MArrivalTime")); + if (!fArrivalTime) + return kFALSE; + + return kTRUE; +} + +void MArrivalTimeCalc2::FindSignalTime(Byte_t *ptr, Byte_t size, Float_t &time, Int_t &sat, Byte_t pedes) const +{ + const Byte_t *end = ptr + size; + + Int_t sum=0; // integral content of the actual window + Int_t max = 0; // highest integral content of all windows + + // + // Calculate the sum of the first fWindowSize slices + // + sat = 0; + Byte_t *p = ptr; + + + while (p= fSaturationLimit) + sat++; + } + + // + // Check for saturation in all other slices + // + while (p= fSaturationLimit) + sat++; + + // + // Calculate the i-th sum as + // sum_i+1 = sum_i + slice[i+8] - slice[i] + // This is fast and accurate (because we are using int's) + // + max=sum; + Byte_t *ptrmax=ptr; // pointer to the first slice of the maximum window + + for (p=ptr; p+fWindowSizemax) + { + max = sum; + ptrmax = p; + } + } + + + // now calculate the time for the maximum window + + + Int_t timesignalsum = 0; + Int_t timesum =0; + + + for (p=ptrmax; p < ptrmax + fWindowSize; p++) + { + timesignalsum += *p*(p-ptr); + timesum += p-ptr; + } + + + float pedsubsum = max - fWindowSize*pedes; + + float pedsubtimesignalsum = timesignalsum - timesum*pedes; + + time = pedsubsum != 0 ? pedsubtimesignalsum / pedsubsum : 1; + + +} + +// -------------------------------------------------------------------------- +// +// Calculate the integral of the FADC of fWindowSize time slices which have the +// highest signal +// +Int_t MArrivalTimeCalc2::Process() +{ + MRawEvtPixelIter pixel(fRawEvt); + + + Int_t sat=0; + while (pixel.Next()) + { + // + // Find signal in hi- and lo-gain + // + Float_t timehi, timelo; + Int_t sathi, satlo; + + // + // Take correspodning pedestal + // + const Int_t pixid = pixel.GetPixelId(); + + const MPedestalPix &ped = (*fPedestals)[pixid]; + + const Float_t pedes = ped.GetPedestal(); + // const Float_t pedrms = ped.GetPedestalRms(); + + + FindSignalTime(pixel.GetHiGainSamples()+fHiGainFirst, fNumHiGainSamples, timehi, sathi, pedes); + FindSignalTime(pixel.GetLoGainSamples()+fLoGainFirst, fNumLoGainSamples, timelo, satlo, pedes); + if (satlo) + sat++; + + if (!sathi) + { + fArrivalTime->SetTime(pixid,timehi); + } + else + fArrivalTime->SetTime(pixid,timelo); + + + } /* while (pixel.Next()) */ + + fArrivalTime->SetReadyToSave(); + + // + // Print a warning if event has saturationg lo-gains + // + if (sat) + *fLog << warn << "WARNING - Lo Gain saturated in " << sat << " pixels." << endl; + + + return kTRUE; +} Index: trunk/MagicSoft/Mars/manalysis/MArrivalTimeCalc2.h =================================================================== --- trunk/MagicSoft/Mars/manalysis/MArrivalTimeCalc2.h (revision 3101) +++ trunk/MagicSoft/Mars/manalysis/MArrivalTimeCalc2.h (revision 3101) @@ -0,0 +1,58 @@ +#ifndef MARS_MArrivalTimeCalc2 +#define MARS_MArrivalTimeCalc2 + +#ifndef MARS_MTask +#include "MTask.h" +#endif + +class MRawEvtData; +class MRawRunHeader; + +class MPedestalCam; +class MExtractedSignalCam; +class MArrivalTime; + +class MArrivalTimeCalc2 : public MTask +{ +private: + static const Byte_t fgSaturationLimit; + static const Byte_t fgFirst; + static const Byte_t fgLast; + static const Byte_t fgWindowSize; + + MPedestalCam *fPedestals; // Pedestals of all pixels in the camera + // MExtractedSignalCam *fSignals; // Extracted signal of all pixels in the camera + + MRawEvtData *fRawEvt; // raw event data (time slices) + MRawRunHeader *fRunHeader; // RunHeader information + + + MArrivalTime *fArrivalTime; // Arrival Time of FADC sample + + + Byte_t fHiGainFirst; // First hi gain to be used + Byte_t fLoGainFirst; // First lo gain to be used + + Byte_t fNumHiGainSamples; // Number of hi gain to be used + Byte_t fNumLoGainSamples; // Number of lo gain to be used + + Byte_t fWindowSize; // Number of gains in window + Float_t fWindowSizeSqrt; // Sqaure root of number of gains in window + + Byte_t fSaturationLimit; + + void FindSignalTime(Byte_t *ptr, Byte_t size, Float_t &time, Int_t &sat, Byte_t pedes) const; + + Int_t PreProcess(MParList *pList); + Int_t Process(); + +public: + MArrivalTimeCalc2(const char *name=NULL, const char *title=NULL); + + void SetRange(Byte_t hifirst=fgFirst, Byte_t hilast=fgLast, Byte_t lofirst=fgFirst, Byte_t lolast=fgLast, Byte_t windowsize=fgWindowSize); + void SetSaturationLimit(Byte_t lim) { fSaturationLimit = lim; } + + ClassDef(MArrivalTimeCalc2, 0) // Extracted Signal as highest integral content +}; + +#endif Index: trunk/MagicSoft/Mars/manalysis/Makefile =================================================================== --- trunk/MagicSoft/Mars/manalysis/Makefile (revision 3100) +++ trunk/MagicSoft/Mars/manalysis/Makefile (revision 3101) @@ -86,4 +86,5 @@ MArrivalTime.cc \ MArrivalTimeCalc.cc \ + MArrivalTimeCalc2.cc \ MMcCalibrationUpdate.cc