Index: trunk/MagicSoft/Mars/Changelog =================================================================== --- trunk/MagicSoft/Mars/Changelog (revision 9331) +++ trunk/MagicSoft/Mars/Changelog (revision 9332) @@ -18,4 +18,55 @@ -*-*- END OF LINE -*-*- + 2009/02/13 Thomas Bretz + + * mcorsika/MCorsikaEvtHeader.h: + - added some comments on the coordinate system + + * ceres.rc: + - added some comments + - setup a psf for the dwarf mirror + + * mcorsika/MCorsikaRunHeader.h: + - added some comments + + * mjobs/MJSimulation.cc: + - added some histogram titles + - fixed the x, y coordinates of the impact histogram + - write "ElectronicNoise" + + * mjobs/MJStar.cc: + - don't force writing of containers which are not filles by ceres + - updated rule to support files without date + + * mpedestal/MPedestalPix.[h,cc]: + - added Print() function + + * mraw/MRawFileRead.cc: + - added a comment + + * msim/MPhotonData.h: + - added comments about coordinate system + + * msimcamera/MSimGeomCam.cc: + - when looking for the pixels convert from ceres-coordinates + to GeomCam camera-coordinates + + * msimreflector/MMirror.h: + - added Setter for SigmaPSF + + * msimreflector/MMirrorDisk.cc, msimreflector/MMirrorHex.cc, + msimreflector/MMirrorSquare.cc: + - now start counting in ReadM with 0 + + * msimreflector/MReflector.[h,cc]: + - added some documentation + - added function to set psf of all mirrors + - added psf optional to input file + - remove already evaluated tokens from TObjArray + - added reading of a default psf from the input file + + + + 2009/02/12 Thomas Bretz Index: trunk/MagicSoft/Mars/ceres.rc =================================================================== --- trunk/MagicSoft/Mars/ceres.rc (revision 9331) +++ trunk/MagicSoft/Mars/ceres.rc (revision 9332) @@ -74,6 +74,19 @@ # ------------------------------------------------------------------------- Reflector.Constructor: MReflector + +# For the file definition see MReflector::ReadFile + #Reflector.FileName: mreflector/reflector-magic.def +# ~0.020° +#Reflector.SetSigmaPSF: 6.0 +# ~0.025° +#Reflector.SetSigmaPSF: 7.5 + Reflector.FileName: mreflector/reflector-dwarf.def +# ~0.02° +#Reflector.SetSigmaPSF: 1.8 +# ~0.03° +Reflector.SetSigmaPSF: 2.6 + # --- DWARF --- @@ -86,4 +99,5 @@ #paraboloid_focal 488.87 +# Hex cones maximum D=8mm GeomCones.Constructor: MGeomCamDwarf(209.5, 7.5, 4.8887); MGeomCam.Constructor: MGeomCamDwarf(209.5, 15, 4.8887); Index: trunk/MagicSoft/Mars/mcorsika/MCorsikaRunHeader.h =================================================================== --- trunk/MagicSoft/Mars/mcorsika/MCorsikaRunHeader.h (revision 9331) +++ trunk/MagicSoft/Mars/mcorsika/MCorsikaRunHeader.h (revision 9332) @@ -33,6 +33,6 @@ Float_t fZdMin; // [rad] Zenith distance Float_t fZdMax; // [rad] Zenith distance - Float_t fAzMin; // [rad] Azimuth (north=0; west=90) - Float_t fAzMax; // [rad] Azimuth (north=0; west=90) + Float_t fAzMin; // [rad] Azimuth (north=0; east=90) + Float_t fAzMax; // [rad] Azimuth (north=0; east=90) (north denotes the magnet north which is defined to be in the geografic north!) Float_t fWavelengthMin; // [nm] Wavelength bandwidth lo edge Index: trunk/MagicSoft/Mars/mjobs/MJSimulation.cc =================================================================== --- trunk/MagicSoft/Mars/mjobs/MJSimulation.cc (revision 9331) +++ trunk/MagicSoft/Mars/mjobs/MJSimulation.cc (revision 9332) @@ -184,5 +184,5 @@ hist.SetLog(kTRUE, kTRUE, kFALSE); - hist.AddHist("MCorsikaEvtHeader.fX/100","MCorsikaEvtHeader.fY/100"); + hist.AddHist("-MCorsikaEvtHeader.fX/100","-MCorsikaEvtHeader.fY/100"); hist.SetDrawOption("colz"); hist.InitName("Impact;Impact;Impact"); @@ -361,14 +361,14 @@ mhew.SetTitle("Time between first and last photon hitting a detector"); - MFillH fillh1(&mhn1, "", "FillH1"); + MFillH fillh1(&mhn1, "", "FillCorsika"); MFillH fillh2(&mhn2, "", "FillH2"); MFillH fillh3(&mhn3, "", "FillH3"); MFillH filltp(&mhtp, "", "FillTriggerPos"); MFillH fillew(&mhew, "", "FillEvtWidth"); - fillh1.SetNameTab("H1", "Distribution of Muons as simulated"); - fillh2.SetNameTab("H2", "Distribution of Muons as available after all"); - fillh3.SetNameTab("H3", "Distribution after trigger"); - filltp.SetNameTab("TrigPos", "TriggerPosition w.r.t the first photon"); - fillew.SetNameTab("EvtWidth", "Time between first and last photon hitting a detector"); + fillh1.SetNameTab("Corsika", "Distribution as simulated by Corsika"); + fillh2.SetNameTab("Detectable", "Distribution of events hit the detector"); + fillh3.SetNameTab("Triggered", "Distribution of triggered events"); + filltp.SetNameTab("TrigPos", "TriggerPosition w.r.t the first photon"); + fillew.SetNameTab("EvtWidth", "Time between first and last photon hitting a detector"); MHPhotonEvent planeG(1); // Get from MaxImpact @@ -457,4 +457,5 @@ write3a.AddContainer("MMcRunHeader", "RunHeaders"); write3a.AddContainer("MMcCorsikaRunHeader", "RunHeaders", kFALSE); + write3a.AddContainer("ElectronicNoise", "RunHeaders", kFALSE); write3a.AddContainer("MMcEvt", "Events", kFALSE); @@ -511,4 +512,5 @@ // ------------------------------------------------------------------- + // FIXME: Remove isolated pixels MImgCleanStd clean(7, 4.5); clean.SetMethod(MImgCleanStd::kAbsolute); @@ -613,4 +615,5 @@ if (fCamera) { + // FIXME: MHCollectionArea Trigger Area! if (header.IsDataRun()) tasks.AddToList(&fillh3); Index: trunk/MagicSoft/Mars/mjobs/MJStar.cc =================================================================== --- trunk/MagicSoft/Mars/mjobs/MJStar.cc (revision 9331) +++ trunk/MagicSoft/Mars/mjobs/MJStar.cc (revision 9332) @@ -300,5 +300,5 @@ // Effective on-time need its own not to be skipped by (eg) image cleaning // Muons needs its own to have a unique SetReadyToSave - const TString rule(Form("s/([0-9]+(_M[12])?_[0-9.]+)_Y_(.*[.]root)$/%s\\/$1_I_$3/", Esc(fPathOut).Data())); + const TString rule(Form("s/(([0-9]+_)?(M[12]_)?[0-9.]+)_Y_(.*[.]root)$/%s\\/$1_I_$4/", Esc(fPathOut).Data())); MWriteRootFile write( 2, rule, fOverwrite?"RECREATE":"NEW", "Image parameters"); MWriteRootFile writet(2, rule, fOverwrite?"RECREATE":"NEW", "Image parameters"); // EffectiveOnTime @@ -333,10 +333,10 @@ // Monte Carlo Data write.AddContainer("MMcEvt", "Events"); - write.AddContainer("MMcTrig", "Events"); + write.AddContainer("MMcTrig", "Events", kFALSE); // Monte Carlo Run Headers write.AddContainer("MMcRunHeader", "RunHeaders"); - write.AddContainer("MMcTrigHeader", "RunHeaders"); - write.AddContainer("MMcFadcHeader", "RunHeaders"); - write.AddContainer("MMcConfigRunHeader", "RunHeaders"); + write.AddContainer("MMcTrigHeader", "RunHeaders", kFALSE); + write.AddContainer("MMcFadcHeader", "RunHeaders", kFALSE); + write.AddContainer("MMcConfigRunHeader", "RunHeaders", kFALSE); write.AddContainer("MMcCorsikaRunHeader", "RunHeaders"); } @@ -368,5 +368,5 @@ // Monte Carlo Data writem.AddContainer("MMcEvt", "Muons"); - writem.AddContainer("MMcTrig", "Muons"); + writem.AddContainer("MMcTrig", "Muons", kFALSE); } Index: trunk/MagicSoft/Mars/mpedestal/MPedestalPix.cc =================================================================== --- trunk/MagicSoft/Mars/mpedestal/MPedestalPix.cc (revision 9331) +++ trunk/MagicSoft/Mars/mpedestal/MPedestalPix.cc (revision 9332) @@ -95,4 +95,16 @@ } +void MPedestalPix::Print(Option_t *) const +{ + *fLog << all << "P=" << fPedestal; + if (fPedestalRms>=0) + *fLog << " RMS=" << fPedestalRms; + if (fPedestalABoffset>=0) + *fLog << " AB=" << fPedestalABoffset; + if (fNumEvents>0) + *fLog << " (N=" << fNumEvents << ")"; + *fLog << endl; +} + // ------------------------------------------------------------------------ // Index: trunk/MagicSoft/Mars/mpedestal/MPedestalPix.h =================================================================== --- trunk/MagicSoft/Mars/mpedestal/MPedestalPix.h (revision 9331) +++ trunk/MagicSoft/Mars/mpedestal/MPedestalPix.h (revision 9332) @@ -19,4 +19,5 @@ void Clear(Option_t *o=""); void Copy(TObject &object) const; + void Print(Option_t *o="") const; // Using histograms Index: trunk/MagicSoft/Mars/mraw/MRawFileRead.cc =================================================================== --- trunk/MagicSoft/Mars/mraw/MRawFileRead.cc (revision 9331) +++ trunk/MagicSoft/Mars/mraw/MRawFileRead.cc (revision 9332) @@ -202,4 +202,6 @@ if (fDisplay) { + // Show the new file name and the event number after which + // the new file has been opened TString txt = GetFileName(); txt += " @ "; Index: trunk/MagicSoft/Mars/msim/MPhotonData.h =================================================================== --- trunk/MagicSoft/Mars/msim/MPhotonData.h (revision 9331) +++ trunk/MagicSoft/Mars/msim/MPhotonData.h (revision 9332) @@ -24,9 +24,9 @@ class MCorsikaRunHeader; -class MPhotonData : public TObject //MParContainer +class MPhotonData : public TObject { private: - Float_t fPosX; // [cm] X (north) at observation level - Float_t fPosY; // [cm] Y (west) at observation level + Float_t fPosX; // [cm] "+west" "-east" (both at observation level) + Float_t fPosY; // [cm] "+south" "-north" (north denotes the magnet north which is defined to be in the geografic north!) Float_t fCosU; // [cos x] U direction cosine to x-axis Index: trunk/MagicSoft/Mars/msimcamera/MSimGeomCam.cc =================================================================== --- trunk/MagicSoft/Mars/msimcamera/MSimGeomCam.cc (revision 9331) +++ trunk/MagicSoft/Mars/msimcamera/MSimGeomCam.cc (revision 9332) @@ -148,5 +148,10 @@ { // FIXME: Improve search algorithm (2D Binary search?) - if (!(*fGeom)[idx].IsInside(ph.GetPosX()*10, ph.GetPosY()*10)) + // Here we convert the photons from the ceres-coordinate + // system which is viewed from the camera to the mirror + // into the camera coordinates which are viewed from + // the mirror to the camera. + // (x on the right, y upwards, right-handed) + if (!(*fGeom)[idx].IsInside(-ph.GetPosX()*10, ph.GetPosY()*10)) continue; Index: trunk/MagicSoft/Mars/msimreflector/MMirror.cc =================================================================== --- trunk/MagicSoft/Mars/msimreflector/MMirror.cc (revision 9331) +++ trunk/MagicSoft/Mars/msimreflector/MMirror.cc (revision 9332) @@ -126,6 +126,6 @@ psf /= 2; // The factor two because of the doubleing of the angle in the reflection psf /= F; // Scale the Gauss to the size of the PSF - psf *= n.Z(); // Normalize the addon vector to the normal vector - //psf *= n.Mag(); // Alternative! (Gaussian projected on the surface of a sphere) + //psf *= n.Z(); // + psf *= n.Mag(); // This means that the PSF is measured in the focal distance TVector3 dn(gx*psf, gy*psf, 0); // Instead of psf/F also atan(psf/F) might make sense @@ -136,4 +136,6 @@ } +// -------------------------------------------------------------------------- +// void MMirror::Print(Option_t *o) const { @@ -141,4 +143,6 @@ gLog << fNorm.X() << " " << fNorm.Y() << " " << fNorm.Z() << " "; gLog << fFocalLength << " "; + if (fSigmaPSF>0) + gLog << fSigmaPSF << " "; const TString n = ClassName(); Index: trunk/MagicSoft/Mars/msimreflector/MMirror.h =================================================================== --- trunk/MagicSoft/Mars/msimreflector/MMirror.h (revision 9331) +++ trunk/MagicSoft/Mars/msimreflector/MMirror.h (revision 9332) @@ -36,4 +36,6 @@ TVector2 operator-(const MMirror &m) const;// { return TVector2(X()-m.X(), Y()-m.Y()); } + void SetSigmaPSF(Double_t psf) { fSigmaPSF = psf; } + void SetFocalLength(Double_t f) { fFocalLength = f; } void SetPosition(const TVector3 &v) { fPos = v; } void SetNorm(const TVector3 &n) { @@ -54,4 +56,7 @@ const TVector3 &GetNorm() const { return fNorm; } + Double_t GetFocalLength() const { return fFocalLength; } + Double_t GetSigmaPSF() const { return fSigmaPSF; } + Double_t GetDist() const { return fPos.Perp(); } @@ -63,7 +68,4 @@ // ----- Basic function for parabolic mirror ----- - void SetFocalLength(Double_t f) { fFocalLength = f; } - Double_t GetFocalLength() const { return fFocalLength; } - Bool_t ExecuteReflection(MQuaternion &p, MQuaternion &u) const; Index: trunk/MagicSoft/Mars/msimreflector/MMirrorDisk.cc =================================================================== --- trunk/MagicSoft/Mars/msimreflector/MMirrorDisk.cc (revision 9331) +++ trunk/MagicSoft/Mars/msimreflector/MMirrorDisk.cc (revision 9332) @@ -144,8 +144,8 @@ Int_t MMirrorDisk::ReadM(const TObjArray &tok) { - if (tok.GetSize()<9) + if (tok.GetEntries()!=1) return -1; - Double_t r = atof(tok[8]->GetName()); + Double_t r = atof(tok[0]->GetName()); if (r<=0) Index: trunk/MagicSoft/Mars/msimreflector/MMirrorHex.cc =================================================================== --- trunk/MagicSoft/Mars/msimreflector/MMirrorHex.cc (revision 9331) +++ trunk/MagicSoft/Mars/msimreflector/MMirrorHex.cc (revision 9332) @@ -163,8 +163,8 @@ Int_t MMirrorHex::ReadM(const TObjArray &tok) { - if (tok.GetSize()<9) + if (tok.GetEntries()!=1) return -1; - const Double_t d = atof(tok[8]->GetName()); + const Double_t d = atof(tok[0]->GetName()); if (d<=0) Index: trunk/MagicSoft/Mars/msimreflector/MMirrorSquare.cc =================================================================== --- trunk/MagicSoft/Mars/msimreflector/MMirrorSquare.cc (revision 9331) +++ trunk/MagicSoft/Mars/msimreflector/MMirrorSquare.cc (revision 9332) @@ -160,8 +160,8 @@ Int_t MMirrorSquare::ReadM(const TObjArray &tok) { - if (tok.GetSize()<9) + if (tok.GetEntries()!=1) return -1; - Double_t l = atof(tok[8]->GetName()); + Double_t l = atof(tok[0]->GetName()); if (l<=0) Index: trunk/MagicSoft/Mars/msimreflector/MReflector.cc =================================================================== --- trunk/MagicSoft/Mars/msimreflector/MReflector.cc (revision 9331) +++ trunk/MagicSoft/Mars/msimreflector/MReflector.cc (revision 9332) @@ -65,4 +65,13 @@ // -------------------------------------------------------------------------- // +// Set the SigmaPSF of all mirrors currently stored. +// +void MReflector::SetSigmaPSF(Double_t psf) +{ + fMirrors.ForEach(MMirror, SetSigmaPSF)(psf); +} + +// -------------------------------------------------------------------------- +// // Calculate the maximum radius of th ereflector. This is not meant as // a precise number but as a rough estimate e.g. to bin a histogram. @@ -118,8 +127,105 @@ // -------------------------------------------------------------------------- // +// I/O helper for ReadFile to avoid calling "delete arr" before every return +// +MMirror *MReflector::EvalTokens(TObjArray &arr, Double_t defpsf) const +{ + if (arr.GetEntries()<9) + { + *fLog << err << "ERROR - Not enough arguments..." << endl; + return 0; + } + + const TVector3 pos(atof(arr[0]->GetName()), + atof(arr[1]->GetName()), + atof(arr[2]->GetName())); + + const TVector3 norm(atof(arr[3]->GetName()), + atof(arr[4]->GetName()), + atof(arr[5]->GetName())); + + const Double_t F = atof(arr[6]->GetName()); + + const TString val = arr[7]->GetName(); + + const Double_t psf = val.IsFloat() ? val.Atof() : -1; + + const UInt_t n = val.IsFloat() ? 9 : 8; + + TString type = arr[n-1]->GetName(); + type.Prepend("MMirror"); + + for (UInt_t i=0; iInheritsFrom(MMirror::Class())) + { + *fLog << err << "ERROR - Cannot create new instance of class " << type << ": " << endl; + *fLog << "Class doesn't inherit from MMirror." << endl; + return 0; + } + + MMirror *m = static_cast(cls->New()); + if (!m) + { + *fLog << err << "ERROR - Cannot create new instance of class " << type << ": " << endl; + *fLog << " - Class has no default constructor." << endl; + *fLog << " - An abstract member functions of a base class is not overwritten." << endl; + return 0; + } + + m->SetFocalLength(F); + m->SetPosition(pos); + m->SetNorm(norm); + m->SetSigmaPSF(psf>=0 ? psf : defpsf); + + if (m->ReadM(arr)>=0) + return m; + + *fLog << err << "ERROR - " << type << "::ReadM failed." << endl; + + delete m; + return 0; +} + +// -------------------------------------------------------------------------- +// // Read a reflector setup from a file. This needs improvemtn. -// FIXME: Documentation missing! -// -Bool_t MReflector::ReadFile(TString fname) +// +// The file structur is like: +// +// x y z nx ny nz F [psf] Type ... +// +// x: x-coordinate of the mirror's center +// y: y-coordinate of the mirror's center +// z: z-coordinate of the mirror's center +// nx: x-component of the normal vecor in the mirror center +// ny: y-component of the normal vecor in the mirror center +// nz: z-component of the normal vecor in the mirror center +// F: Focal distance of a spherical mirror +// [psf]: This number is the psf given in the units of x,y,z and +// defined at the focal distance F. It can be used to overwrite +// the second argument given in ReadFile for individual mirrors. +// Type: A instance of a mirrot of the class Type MMirrorType is created +// (Type can be, for example, Hex for for MMirrorHex). +// ...: Additional arguments as defined in MMirrorType::ReadM +// +// +// Coordinate System: +// The coordinate system is local in the reflectors frame. +// It is defined viewing from the back side of the reflector +// towards the camera. (x "right", y "up", z from reflector to camera) +// Note, that it is left-handed! +// +Bool_t MReflector::ReadFile(TString fname, Double_t defpsf) { SetTitle(fname); @@ -145,4 +251,6 @@ */ + Int_t cnt = 0; + while (1) { @@ -152,78 +260,36 @@ break; + // Count lines + cnt++; + + // Skip comments if (line.BeginsWith("#")) - { - //cout << line << endl; continue; - } - + + // Remove leading and trailing whitespaces line=line.Strip(TString::kBoth); + // Skip empty lines if (line.IsNull()) continue; + // Tokenize line TObjArray *arr = line.Tokenize(' '); - if (arr->GetSize()<8) + // Evaluate tokens + MMirror *m = EvalTokens(*arr, defpsf); + + // Delete now obsolete array + delete arr; + + // Check if a new mirror could be created successfully + if (!m) { - cout << "Skip3: " <GetName()), - atof((*arr)[1]->GetName()), - atof((*arr)[2]->GetName())); - - const TVector3 norm(atof((*arr)[3]->GetName()), - atof((*arr)[4]->GetName()), - atof((*arr)[5]->GetName())); - - const Double_t F = atof((*arr)[6]->GetName()); - - TString type = (*arr)[7]->GetName(); - type.Prepend("MMirror"); - - TString msg; - TClass *cls = MParContainer::GetClass(type); - if (!cls) - { - *fLog << err << dbginf << "ERROR - Class " << type << " not in dictionary." << endl; + *fLog << err << "Error in line " << cnt << ": " << line << endl; return kFALSE; } - if (!cls->InheritsFrom(MMirror::Class())) - { - *fLog << err << dbginf << "Cannot create new instance of class " << type << ": " << endl; - *fLog << "Class doesn't inherit from MMirror." << endl; - return kFALSE; - } - - MMirror *m = (MMirror*)cls->New(); - if (!m) - { - *fLog << err << dbginf << "Cannot create new instance of class " << type << ": " << endl; - *fLog << " - Class has no default constructor." << endl; - *fLog << " - An abstract member functions of a base class is not overwritten." << endl; - return kFALSE; - } - - m->SetFocalLength(F); - m->SetPosition(pos); - m->SetNorm(norm); - - Int_t n = m->ReadM(*arr); - if (n<=0) - { - *fLog << err << dbginf << "ERROR - ReadM failed." << endl; - return kFALSE; - } - + // Add new mirror to array fMirrors.Add(m); - - //maxr = TMath::Max(maxr, TMath::Hypot(pos[i].X()+24.75, pos[i].Y()+24.75)); - //maxr = TMath::Max(maxr, TMath::Hypot(pos.X()+24.75, pos.Y()+24.75)); - - delete arr; } @@ -285,13 +351,28 @@ // -------------------------------------------------------------------------- // +// SigmaPSF: -1 // FileName: reflector.txt // +// SigmaPSF can be used to set a default for the psf of the mirrors +// read from the file. Note, that this can be overwritten for individual +// mirrors in the file. +// +// For details on the file structure see MReflector::ReadFile +// Int_t MReflector::ReadEnv(const TEnv &env, TString prefix, Bool_t print) { Bool_t rc = kFALSE; + + Double_t psf = -1; + if (IsEnvDefined(env, prefix, "SetSigmaPSF", print)) + { + rc = kTRUE; + psf = GetEnvValue(env, prefix, "SetSigmaPSF", -1); + } + if (IsEnvDefined(env, prefix, "FileName", print)) { rc = kTRUE; - if (!ReadFile(GetEnvValue(env, prefix, "FileName", ""))) + if (!ReadFile(GetEnvValue(env, prefix, "FileName", ""), psf)) return kERROR; } Index: trunk/MagicSoft/Mars/msimreflector/MReflector.h =================================================================== --- trunk/MagicSoft/Mars/msimreflector/MReflector.h (revision 9331) +++ trunk/MagicSoft/Mars/msimreflector/MReflector.h (revision 9332) @@ -28,4 +28,7 @@ void InitMaxR(); + // Helper for I/O + MMirror *EvalTokens(TObjArray &arr, Double_t defpsf) const; + // MParContainer Int_t ReadEnv(const TEnv &env, TString prefix, Bool_t print); @@ -39,5 +42,5 @@ const MMirror *GetMirror(UInt_t idx) const { return idx>=GetNumMirrors()?0:*(GetFirstPtr()+idx); } - Bool_t ReadFile(TString fname); + Bool_t ReadFile(TString fname, Double_t defpsf=-1); Double_t GetMaxR() const { return fMaxR; } @@ -46,4 +49,6 @@ Int_t ExecuteReflector(MQuaternion &p, MQuaternion &u) const; + + void SetSigmaPSF(Double_t psf); // TObject