Index: trunk/MagicSoft/Mars/Changelog =================================================================== --- trunk/MagicSoft/Mars/Changelog (revision 9238) +++ trunk/MagicSoft/Mars/Changelog (revision 9239) @@ -36,4 +36,8 @@ * msimreflector/MSimReflector.cc: - force sorting of the array + + * msimcamera/MSimPSF.[h,cc], msimcamera/MSimGeomCam.[h,cc], + msimcamera/MSimExcessNoise.[h,cc], msimcamera/MSimSignalCam.[h,cc]: + - added Index: trunk/MagicSoft/Mars/NEWS =================================================================== --- trunk/MagicSoft/Mars/NEWS (revision 9238) +++ trunk/MagicSoft/Mars/NEWS (revision 9239) @@ -11,4 +11,15 @@ So it is recommended although not necessary to reset your sequences in the database before further processing. + + ;NEW + + * A first implementation of a Monte Carlo simulation program + (ceres - Camera Electronics and REflector Simulation) has been + implemented. It can simulate a full telescope already including + a complete reflector simulation and a flexible trigger and readout + simulation. Currently, no perfect setup is implemented and + especially the output of the program (data structure) should + be considered alpha state. Also the user interface is not yet + ready. ;general Index: trunk/MagicSoft/Mars/mraw/MRawRunHeader.h =================================================================== --- trunk/MagicSoft/Mars/mraw/MRawRunHeader.h (revision 9238) +++ trunk/MagicSoft/Mars/mraw/MRawRunHeader.h (revision 9239) @@ -159,4 +159,7 @@ Bool_t IsValid() const { return fMagicNumber==0xc0c0 || fMagicNumber==0xc0c1; } Bool_t IsMonteCarloRun() const { return fRunType>0x00ff; } + Bool_t IsDataRun() const { return (fRunType&0xff)==kRTData; } + Bool_t IsPedestalRun() const { return (fRunType&0xff)==kRTPedestal; } + Bool_t IsCalibrationRun() const { return (fRunType&0xff)==kRTCalibration; } void Print(Option_t *t=NULL) const; Index: trunk/MagicSoft/Mars/msim/MPhotonData.cc =================================================================== --- trunk/MagicSoft/Mars/msim/MPhotonData.cc (revision 9238) +++ trunk/MagicSoft/Mars/msim/MPhotonData.cc (revision 9239) @@ -216,46 +216,2 @@ gLog << "Weight: " << fWeight << endl; } - -/* -#include -#include - -// -------------------------------------------------------------------------- -// -// Fill radial photon distance scaled by scale into 1D histogram. -// -void MPhotonData::FillRad(TH1 &hist, Float_t scale) const -{ - hist.Fill(TMath::Hypot(fPosX, fPosY)*scale); -} - -// -------------------------------------------------------------------------- -// -// Fill radial photon distance scaled by scale versus x into 2D histogram. -// -void MPhotonData::FillRad(TH2 &hist, Double_t x, Float_t scale) const -{ - hist.Fill(x, TMath::Hypot(fPosX, fPosY)*scale); -} - -// -------------------------------------------------------------------------- -// -// Fill photon position (x,y) scaled by scale into 2D histogram. -// (Note north in the camera plane is -y, and east -x) -// -void MPhotonData::Fill(TH2 &hist, Float_t scale) const -{ - hist.Fill(fPosY*scale, fPosX*scale); -} - -// -------------------------------------------------------------------------- -// -// Fill photon position (x,y) scaled by scale versus z into 3D histogram. -// (Note north in the camera plane is -y, and east -x) -// -void MPhotonData::Fill(TH3 &hist, Double_t z, Float_t scale) const -{ - hist.Fill(fPosY*scale, fPosX*scale); -} - -*/ Index: trunk/MagicSoft/Mars/msim/MPhotonEvent.cc =================================================================== --- trunk/MagicSoft/Mars/msim/MPhotonEvent.cc (revision 9238) +++ trunk/MagicSoft/Mars/msim/MPhotonEvent.cc (revision 9239) @@ -57,4 +57,6 @@ // Here is an example (how to remove every second entry) // +// --------------------------------------------------------------------- +// // Int_t cnt = 0; // @@ -91,4 +93,11 @@ // MPhotonEvent->Shrink(cnt); // +// --------------------------------------------------------------------- +// +// The flag for a sorted array is for speed reasons not in all conditions +// maintained automatically. Especially Add() doesn't reset it. +// +// So be sure that if you want to sort your array it is really sorted. +// // // Version 1: @@ -127,8 +136,10 @@ } +/* const char *MPhotonEvent::GetClassName() const { return static_cast(fData.GetClass())->GetName(); } +*/ // -------------------------------------------------------------------------- @@ -184,4 +195,26 @@ } +// -------------------------------------------------------------------------- +// +// Return a pointer to the first photon if available. +// +MPhotonData *MPhotonEvent::GetFirst() const +{ + return static_cast(fData.First()); +} + +// -------------------------------------------------------------------------- +// +// Return a pointer to the last photon if available. +// +MPhotonData *MPhotonEvent::GetLast() const +{ + return static_cast(fData.Last()); +} + +// -------------------------------------------------------------------------- +// +// Read the Event section from the file +// Int_t MPhotonEvent::ReadCorsikaEvt(istream &fin) { @@ -234,4 +267,5 @@ Shrink(n); + fData.UnSort(); SetReadyToSave(); @@ -241,4 +275,6 @@ } +// -------------------------------------------------------------------------- +// Int_t MPhotonEvent::ReadRflEvt(std::istream &fin) { @@ -291,113 +327,12 @@ } +// -------------------------------------------------------------------------- +// +// Print the array +// void MPhotonEvent::Print(Option_t *) const { fData.Print(); } - -/* -// -------------------------------------------------------------------------- -// -// Fills all photon distances scaled with scale (default=1) into a TH1 -// -void MPhotonEvent::FillRad(TH1 &hist, Float_t scale) const -{ - MPhotonData *ph=NULL; - - TIter Next(&fData); - while ((ph=(MPhotonData*)Next())) - ph->FillRad(hist, scale); -} - -// -------------------------------------------------------------------------- -// -// Fills all photon distances scaled with scale (default=1) versus x -// into a TH2 -// -void MPhotonEvent::FillRad(TH2 &hist, Double_t x, Float_t scale) const -{ - MPhotonData *ph=NULL; - - TIter Next(&fData); - while ((ph=(MPhotonData*)Next())) - ph->FillRad(hist, x, scale); -} - -// -------------------------------------------------------------------------- -// -// Fills all photons (x,y) scaled with scale (default=1) into a TH2. -// -void MPhotonEvent::Fill(TH2 &hist, Float_t scale) const -{ - MPhotonData *ph=NULL; - - TIter Next(&fData); - while ((ph=(MPhotonData*)Next())) - ph->Fill(hist, scale); -} - -// -------------------------------------------------------------------------- -// -// Fills all photons (x,y) scaled with scale (default=1) versus z into a TH3 -// -void MPhotonEvent::Fill(TH3 &hist, Double_t z, Float_t scale) const -{ - MPhotonData *ph=NULL; - - TIter Next(&fData); - while ((ph=(MPhotonData*)Next())) - ph->Fill(hist, z, scale); -} - -void MPhotonEvent::DrawPixelContent(Int_t num) const -{ -} - -#include "MHexagon.h" -#include "MGeomCam.h" -#include - -// ------------------------------------------------------------------------ -// -// Fill a reflector event. Sums all pixels in each pixel as the -// pixel contents. -// -// WARNING: Due to the estimation in DistanceToPrimitive and the -// calculation in pixels instead of x, y this is only a -// rough estimate. -// -Bool_t MPhotonEvent::GetPixelContent(Double_t &val, Int_t idx, const MGeomCam &cam, Int_t type) const -{ - // - // sum the photons content in each pixel - // - val = 0; - - MHexagon hex(cam[idx]); - - MPhotonData *ph=NULL; - - TIter Next(&fData); - - switch (type) - { - - case 1: // time - while ((ph=(MPhotonData*)Next())) - if (hex.DistanceToPrimitive(ph->GetPosY(), ph->GetPosX())<=0) - val += ph->GetTime(); - return val/fData.GetEntriesFast()>0; - - default: // signal - while ((ph=(MPhotonData*)Next())) - if (hex.DistanceToPrimitive(ph->GetPosY()*10, ph->GetPosX()*10)<=0) - val += cam.GetPixRatio(idx); - return val>0; - } - - return kFALSE; - -} -*/ // ------------------------------------------------------------------------ Index: trunk/MagicSoft/Mars/msim/MPhotonEvent.h =================================================================== --- trunk/MagicSoft/Mars/msim/MPhotonEvent.h (revision 9238) +++ trunk/MagicSoft/Mars/msim/MPhotonEvent.h (revision 9239) @@ -5,22 +5,10 @@ #include "MParContainer.h" #endif -/* -#ifndef MARS_MCamEvent -#include "MCamEvent.h" -#endif -*/ + #ifndef ROOT_TClonesArray #include #endif -// gcc 3.2 -//class ifstream; #include - -/* -class TH1; -class TH2; -class TH3; -*/ class MPhotonData; @@ -71,7 +59,9 @@ Changed(); } + + void SetSorted() { fSorted = kTRUE; } }; -class MPhotonEvent : public MParContainer//, public MCamEvent +class MPhotonEvent : public MParContainer { private: @@ -81,15 +71,20 @@ MPhotonEvent(const char *name=NULL, const char *title=NULL); - //void Clear(Option_t * = NULL); - void Print(Option_t * = NULL) const; + // TObjArray, FIXME: This could be improved if checking for IsSortable is omitted + void Sort(Bool_t force=kFALSE) { if (force) fData.UnSort(); fData.Sort(); } + Bool_t IsSorted() const { return fData.IsSorted(); } - const char *GetClassName() const; + + // Getter/Setter + Int_t GetNumPhotons() const { return fData.GetEntriesFast(); } TClonesArray &GetArray() { return fData; } const TClonesArray &GetArray() const { return fData; } - void Remove(TObject *o) { fData.Remove(o); } MPhotonData &Add(Int_t n); MPhotonData &Add(); + + MPhotonData *GetFirst() const; + MPhotonData *GetLast() const; MPhotonData &operator[](UInt_t idx); @@ -106,31 +101,25 @@ // for (int i=fData.GetSize()-1; i>=n; i--) // fData.RemoveAt(i); - static_cast(fData).FastRemove(n, fData.GetSize()-1); + const Bool_t sorted = fData.IsSorted(); + + MyClonesArray &loc = static_cast(fData); + + loc.FastRemove(n, fData.GetSize()-1); + + // If it was sorted before it is also sorted now. + if (sorted) + loc.SetSorted(); return fData.GetEntriesFast(); } -/* - void Expand(Int_t n) - { - for (int i=fData.GetSize(); i +! +! Copyright: Software Development, 2000-2009 +! +! +\* ======================================================================== */ + +////////////////////////////////////////////////////////////////////////////// +// +// MSimExcessNoise +// +// This task adds the noise (usually signal height, i.e. excess, dependent) +// to the photon signal. +// +// Input Containers: +// MCorsikaEvent +// +// Output Containers: +// MCorsikaEvent +// +////////////////////////////////////////////////////////////////////////////// +#include "MSimExcessNoise.h" + +#include + +#include "MLog.h" +#include "MLogManip.h" + +#include "MParList.h" + +#include "MPhotonEvent.h" +#include "MPhotonData.h" + +ClassImp(MSimExcessNoise); + +using namespace std; + +// -------------------------------------------------------------------------- +// +// Default Constructor. +// +MSimExcessNoise::MSimExcessNoise(const char* name, const char *title) +: fEvt(0) +{ + fName = name ? name : "MSimExcessNoise"; + fTitle = title ? title : "Task to simulate the excess dependant noise (conversion photon to signal height)"; +} + +// -------------------------------------------------------------------------- +// +// Check for the necessary parameter containers. +// +Int_t MSimExcessNoise::PreProcess(MParList *pList) +{ + fEvt = (MPhotonEvent*)pList->FindObject("MPhotonEvent"); + if (!fEvt) + { + *fLog << err << "MPhotonEvent not found... aborting." << endl; + return kFALSE; + } + + return kTRUE; +} + +// -------------------------------------------------------------------------- +// +// Change the weight of each signal according to the access noise +// +Int_t MSimExcessNoise::Process() +{ + const UInt_t num = fEvt->GetNumPhotons(); + for (UInt_t i=0; iGaus(oldw, 0.2*TMath::Sqrt(oldw)); + ph.SetWeight(neww); + } + + return kTRUE; +} Index: trunk/MagicSoft/Mars/msimcamera/MSimExcessNoise.h =================================================================== --- trunk/MagicSoft/Mars/msimcamera/MSimExcessNoise.h (revision 9239) +++ trunk/MagicSoft/Mars/msimcamera/MSimExcessNoise.h (revision 9239) @@ -0,0 +1,26 @@ +#ifndef MARS_MSimExcessNoise +#define MARS_MSimExcessNoise + +#ifndef MARS_MTask +#include "MTask.h" +#endif + +class MParList; +class MPhotonEvent; + +class MSimExcessNoise : public MTask +{ +private: + MPhotonEvent *fEvt; //! Event storing the photons + + // MTask + Int_t PreProcess(MParList *pList); + Int_t Process(); + +public: + MSimExcessNoise(const char *name=NULL, const char *title=NULL); + + ClassDef(MSimExcessNoise, 0) // Task to simulate the excess dependant noise (conversion photon to signal height) +}; + +#endif Index: trunk/MagicSoft/Mars/msimcamera/MSimGeomCam.cc =================================================================== --- trunk/MagicSoft/Mars/msimcamera/MSimGeomCam.cc (revision 9239) +++ trunk/MagicSoft/Mars/msimcamera/MSimGeomCam.cc (revision 9239) @@ -0,0 +1,197 @@ +/* ======================================================================== *\ +! +! * +! * 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 1/2009 +! +! Copyright: Software Development, 2000-2009 +! +! +\* ======================================================================== */ + +////////////////////////////////////////////////////////////////////////////// +// +// MSimGeomCam +// +// This task takes a photon list from a MPhotonEvent and checks which pixel +// from a MGeomCam is hit. The photons are tagged with the corresponding +// index. +// +// Additionally (and provisionally) it also calculates the photon event +// statistics. +// +// Input Containers: +// MPhotonEvent +// fNameGeomCam [MGeomCam] +// IntendedTrigPos [MParameterD] +// MRawRunHeader +// +// Output Containers: +// MPhotonStatistics +// -/- +// +////////////////////////////////////////////////////////////////////////////// +#include "MSimGeomCam.h" + +#include "MLog.h" +#include "MLogManip.h" + +#include "MParList.h" + +#include "MGeomCam.h" +#include "MGeomPix.h" + +#include "MPhotonEvent.h" +#include "MPhotonData.h" + +#include "MParameters.h" +#include "MRawRunHeader.h" + +ClassImp(MSimGeomCam); + +using namespace std; + +// -------------------------------------------------------------------------- +// +// Default Constructor. +// +MSimGeomCam::MSimGeomCam(const char* name, const char *title) + : fGeom(0), fEvt(0), fStat(0), fNameGeomCam("MGeomCam") +{ + fName = name ? name : "MSimGeomCam"; + fTitle = title ? title : "Task to tag each photon in a MPhotonEvent with a pixel index from a MGeomCam"; +} + +// -------------------------------------------------------------------------- +// +// Search for the needed parameter containers. +// +Int_t MSimGeomCam::PreProcess(MParList *pList) +{ + fGeom = (MGeomCam*)pList->FindObject(fNameGeomCam, "MGeomCam"); + if (!fGeom) + { + *fLog << inf << fNameGeomCam << " [MGeomCam] not found..." << endl; + + fGeom = (MGeomCam*)pList->FindCreateObj(fNameGeomCam); + if (!fGeom) + return kFALSE; + } + + fEvt = (MPhotonEvent*)pList->FindObject("MPhotonEvent"); + if (!fEvt) + { + *fLog << err << "MPhotonEvent not found... aborting." << endl; + return kFALSE; + } + + fPulsePos = (MParameterD*)pList->FindObject("IntendedPulsePos", "MParameterD"); + if (!fPulsePos) + { + *fLog << err << "IntendedPulsePos [MParameterD] not found... aborting." << endl; + return kFALSE; + } + + fHeader = (MRawRunHeader*)pList->FindObject("MRawRunHeader"); + if (!fHeader) + { + *fLog << err << "MRawRunHeader not found... aborting." << endl; + return kFALSE; + } + + fStat = (MPhotonStatistics*)pList->FindCreateObj("MPhotonStatistics"); + if (!fStat) + return kFALSE; + + return kTRUE; +} + + +// -------------------------------------------------------------------------- +// +Int_t MSimGeomCam::Process() +{ + const Int_t num = fEvt->GetNumPhotons(); + + Int_t cnt = 0; + for (Int_t i=0; iGetNumPixels(); idx++) + { + // FIXME: Improve search algorithm (2D Binary search?) + if (!(*fGeom)[idx].IsInside(ph.GetPosX()*10, ph.GetPosY()*10)) + continue; + + ph.SetTag(idx); + + (*fEvt)[cnt++] = ph; + + break; + } + } + + fEvt->Shrink(cnt); + //fEvt->Sort(); + + // ------ FIXME: Move somewhere else? MSimCalibrationSignal ------ + + if (!fEvt->IsSorted()) + { + *fLog << err << "ERROR - MPhotonEvent must be sorted!" << endl; + return kERROR; + } + + const Float_t first = cnt>0 ? fEvt->GetFirst()->GetTime() : 0; + const Float_t last = cnt>0 ? fEvt->GetLast()->GetTime() : -1; + + // Length (ns), Pulse position (Units ns) + const Float_t ns = fHeader->GetFreqSampling()/1000.*fHeader->GetNumSamplesHiGain(); + const Float_t pp = fPulsePos->GetVal(); + + // FIXME FIXME FIXME --- Make it depending on the trigger position? + + fStat->SetTime(first-pp, last+(ns-pp)); + fStat->SetMaxIndex(fGeom->GetNumPixels()-1); + fStat->SetReadyToSave(); + + // ---------------------------------------------------------------------- + + return kTRUE; +} + +// -------------------------------------------------------------------------- +// +// Read the parameters from the resource file. +// +// NameGeometry: MGeomCamDwarf +// +Int_t MSimGeomCam::ReadEnv(const TEnv &env, TString prefix, Bool_t print) +{ + Bool_t rc = kFALSE; + if (IsEnvDefined(env, prefix, "NameGeometry", print)) + { + rc = kTRUE; + SetNameGeomCam(GetEnvValue(env, prefix, "NameGeometry", fNameGeomCam)); + // FIXME: What about setup of this container? + } + + return rc; +} Index: trunk/MagicSoft/Mars/msimcamera/MSimGeomCam.h =================================================================== --- trunk/MagicSoft/Mars/msimcamera/MSimGeomCam.h (revision 9239) +++ trunk/MagicSoft/Mars/msimcamera/MSimGeomCam.h (revision 9239) @@ -0,0 +1,43 @@ +#ifndef MARS_MSimGeomCam +#define MARS_MSimGeomCam + +#ifndef MARS_MTask +#include "MTask.h" +#endif + +class MGeomCam; +class MParList; +class MPhotonEvent; +class MPhotonStatistics; +class MSignalCam; +class MParameterD; +class MRawRunHeader; + +class MSimGeomCam : public MTask +{ +private: + MGeomCam *fGeom; //! + MPhotonEvent *fEvt; //! Event stroing the photons + MPhotonStatistics *fStat; //! + MParameterD *fPulsePos; //! Intended pulse position in digitization window [ns] + MRawRunHeader *fHeader; //! Length of digitization window + + + TString fNameGeomCam; + + // MParContainer + Int_t ReadEnv(const TEnv &env, TString prefix, Bool_t print=kFALSE); + + // MTask + Int_t PreProcess(MParList *pList); + Int_t Process(); + +public: + MSimGeomCam(const char *name=NULL, const char *title=NULL); + + void SetNameGeomCam(const char *name="MGeomCam") { fNameGeomCam=name; } + + ClassDef(MSimGeomCam, 0) // Task to tag each photon in a MPhotonEvent with a pixel index from a MGeomCam +}; + +#endif Index: trunk/MagicSoft/Mars/msimcamera/MSimPSF.cc =================================================================== --- trunk/MagicSoft/Mars/msimcamera/MSimPSF.cc (revision 9239) +++ trunk/MagicSoft/Mars/msimcamera/MSimPSF.cc (revision 9239) @@ -0,0 +1,124 @@ +/* ======================================================================== *\ +! +! * +! * 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 1/2009 +! +! Copyright: Software Development, 2000-2009 +! +! +\* ======================================================================== */ + +////////////////////////////////////////////////////////////////////////////// +// +// MSimPSF +// +// This task makes a naiv simulation of the psf by smearing out the photons +// in a plane (camera plane) by a 2D-Gaussian with fSigma +// +// Input Containers: +// MPhotonEvent +// +// Output Containers: +// MPhotonEvent +// +////////////////////////////////////////////////////////////////////////////// +#include "MSimPSF.h" + +#include + +#include "MLog.h" +#include "MLogManip.h" + +#include "MParList.h" + +#include "MPhotonEvent.h" +#include "MPhotonData.h" + +ClassImp(MSimPSF); + +using namespace std; + +// -------------------------------------------------------------------------- +// +// Default Constructor. +// +MSimPSF::MSimPSF(const char* name, const char *title) + : fEvt(0), fSigma(1) +{ + fName = name ? name : "MSimPSF"; + fTitle = title ? title : "Task to do a naiv simulation of the psf by smearout in the camera plane"; +} + +// -------------------------------------------------------------------------- +// +// Search for MPhotonEvent +// +Int_t MSimPSF::PreProcess(MParList *pList) +{ + if (fSigma<=0) + return kSKIP; + + fEvt = (MPhotonEvent*)pList->FindObject("MPhotonEvent"); + if (!fEvt) + { + *fLog << err << "MPhotonEvent not found... aborting." << endl; + return kFALSE; + } + + return kTRUE; +} + +// -------------------------------------------------------------------------- +// +// Smear out all photons by a gaussian with fSigma +// +Int_t MSimPSF::Process() +{ + const UInt_t num = fEvt->GetNumPhotons(); + + // Loop over all mirrors + for (UInt_t i=0; iGaus(0, fSigma), gRandom->Gaus(0, fSigma)); + + // Add random smear out + ph.SetPosition(ph.GetPos2()+v); + } + + return kTRUE; +} + +// -------------------------------------------------------------------------- +// +// Sigma: 10 +// +Int_t MSimPSF::ReadEnv(const TEnv &env, TString prefix, Bool_t print) +{ + Bool_t rc = kFALSE; + + if (IsEnvDefined(env, prefix, "Sigma", print)) + { + rc = kTRUE; + fSigma = GetEnvValue(env, prefix, "Sigma", fSigma); + } + + return rc; +} Index: trunk/MagicSoft/Mars/msimcamera/MSimPSF.h =================================================================== --- trunk/MagicSoft/Mars/msimcamera/MSimPSF.h (revision 9239) +++ trunk/MagicSoft/Mars/msimcamera/MSimPSF.h (revision 9239) @@ -0,0 +1,33 @@ +#ifndef MARS_MSimPSF +#define MARS_MSimPSF + +#ifndef MARS_MTask +#include "MTask.h" +#endif + +class MParList; +class MPhotonEvent; + +class MSimPSF : public MTask +{ +private: + MPhotonEvent *fEvt; //! Event stroing the photons + + Double_t fSigma; // Gaussian sigma of the smearout + + // MParContainer + Int_t ReadEnv(const TEnv &env, TString prefix, Bool_t print=kFALSE); + + // MTask + Int_t PreProcess(MParList *pList); + Int_t Process(); + + +public: + MSimPSF(const char *name=NULL, const char *title=NULL); + + ClassDef(MSimPSF, 0) // Task to do a naiv simulation of the psf by smearout in the camera plane +}; + +#endif + Index: trunk/MagicSoft/Mars/msimcamera/MSimSignalCam.cc =================================================================== --- trunk/MagicSoft/Mars/msimcamera/MSimSignalCam.cc (revision 9239) +++ trunk/MagicSoft/Mars/msimcamera/MSimSignalCam.cc (revision 9239) @@ -0,0 +1,171 @@ +/* ======================================================================== *\ +! +! * +! * 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 1/2009 +! +! Copyright: Software Development, 2000-2009 +! +! +\* ======================================================================== */ + +////////////////////////////////////////////////////////////////////////////// +// +// MSimSignalCam +// +// This task takes a photon list (MPhotonEvent) and converts it into +// a MSignalCam container. Photons with kNightSky as source are not +// considered. The arrival time is just the average arrival time +// of the photons minus the one of the first. +// +// Input Containers: +// MPhotonEvent +// MPhotonStatistics +// +// Output Containers: +// MSignalCam +// +////////////////////////////////////////////////////////////////////////////// +#include "MSimSignalCam.h" + +#include "MLog.h" +#include "MLogManip.h" + +#include "MParList.h" + +#include "MSignalCam.h" +#include "MSignalPix.h" + +#include "MPhotonEvent.h" +#include "MPhotonData.h" + +ClassImp(MSimSignalCam); + +using namespace std; + +// -------------------------------------------------------------------------- +// +// Default Constructor. +// +MSimSignalCam::MSimSignalCam(const char* name, const char *title) + : fEvt(0), fStat(0), fSignal(0) +{ + fName = name ? name : "MSimSignalCam"; + fTitle = title ? title : "Task to convert a tagged MPhotonEvent list into MSignalCam"; +} + +// -------------------------------------------------------------------------- +// +// Search for the necessary parameter containers +// +Int_t MSimSignalCam::PreProcess(MParList *pList) +{ + fEvt = (MPhotonEvent*)pList->FindObject("MPhotonEvent"); + if (!fEvt) + { + *fLog << err << "MPhotonEvent not found... aborting." << endl; + return kFALSE; + } + + fStat = (MPhotonStatistics*)pList->FindObject("MPhotonStatistics"); + if (!fStat) + { + *fLog << err << "MPhotonStatistics not found... aborting." << endl; + return kFALSE; + } + + fSignal = (MSignalCam*)pList->FindCreateObj("MSignalCam"); + if (!fSignal) + return kFALSE; + + return kTRUE; +} + +// -------------------------------------------------------------------------- +// +// Set the size, i.e. the number of pixels, from MSignalCam. It has to be +// preset e.g. by MGeomApply in a preceeding ReInit. +// +Bool_t MSimSignalCam::ReInit(MParList *plist) +{ + const UInt_t npix = fSignal->GetNumPixels(); + + if (fCont.GetSize()!=npix) + { + fCont.Set(npix); + fTime.Set(npix); + } + + return kTRUE; +} + +// -------------------------------------------------------------------------- +// +// Loop over all photons to sum the photons and determine the average +// arrival time. Write the result into the MSignalCam. +// +// Photons from the NSB are ignored. +// +Int_t MSimSignalCam::Process() +{ + // FIXME: Check the maximum index in GetTag from the statistics container + + fCont.Reset(); + fTime.Reset(); + + // Loop over all photons in the event + const UInt_t num = fEvt->GetNumPhotons(); + for (UInt_t i=0; iGetTimeFirst(); + + // Loop over all pixels and set signal and arrival time. + // Set the pixels valid. + const UInt_t npix = fSignal->GetNumPixels(); + for (UInt_t idx=0; idxUniform(last-first) + pix.SetArrivalTime(fCont[idx]<=0 ? -1 : fTime[idx]/fCont[idx]-first); + pix.SetPixelUsed(); + pix.SetPixelCore(kFALSE); + } + fSignal->SetReadyToSave(); + + return kTRUE; +} Index: trunk/MagicSoft/Mars/msimcamera/MSimSignalCam.h =================================================================== --- trunk/MagicSoft/Mars/msimcamera/MSimSignalCam.h (revision 9239) +++ trunk/MagicSoft/Mars/msimcamera/MSimSignalCam.h (revision 9239) @@ -0,0 +1,45 @@ +#ifndef MARS_MSimSignalCam +#define MARS_MSimSignalCam + +#ifndef MARS_MTask +#include "MTask.h" +#endif + +#ifndef MARS_MArrayI +#include "MArrayI.h" +#endif + +#ifndef MARS_MArrayD +#include "MArrayD.h" +#endif + +class MParList; +class MPhotonEvent; +class MPhotonStatistics; +class MSignalCam; + +class MSimSignalCam : public MTask +{ +private: + MPhotonEvent *fEvt; //! Event containing the photons + MPhotonStatistics *fStat; //! Statistics about the event + MSignalCam *fSignal; //! Output container + + MArrayD fCont; //! Local buffer for contents + MArrayD fTime; //! local buffer for time + + // MParContainer + // Int_t ReadEnv(const TEnv &env, TString prefix, Bool_t print=kFALSE); + + // MTask + Int_t PreProcess(MParList *pList); + Bool_t ReInit(MParList *plist); + Int_t Process(); + +public: + MSimSignalCam(const char *name=NULL, const char *title=NULL); + + ClassDef(MSimSignalCam, 0) // Task to convert a tagged MPhotonEvent list into MSignalCam +}; + +#endif Index: trunk/MagicSoft/Mars/msimreflector/MSimReflector.cc =================================================================== --- trunk/MagicSoft/Mars/msimreflector/MSimReflector.cc (revision 9238) +++ trunk/MagicSoft/Mars/msimreflector/MSimReflector.cc (revision 9239) @@ -373,4 +373,5 @@ // This is to avoid time consuming calculation if there is no // chance of a coincidence. + // FIXME: Inmprove search algorithm (2D Binary search?) if (!mirror.CanHit(p)) continue; @@ -539,7 +540,5 @@ // Doesn't seem to be too time consuming. But we could also sort later! // (after cones, inside the camera) - fEvt->Sort(); - - // FIXME FIXME FIXME: Set maxindex, first and last time. + fEvt->Sort(kTRUE); return kTRUE;