Changeset 9336 for trunk/MagicSoft/Mars/msim

Timestamp:

02/15/09 13:24:59 (16 years ago)

Author:

tbretz

Message:

*** empty log message ***

Location:

trunk/MagicSoft/Mars/msim

Files:

• MSimMMCS.cc (modified) (5 diffs)
• MSimMMCS.h (modified) (1 diff)
• MSimPointingPos.cc (modified) (5 diffs)
• MSimPointingPos.h (modified) (1 diff)

trunk/MagicSoft/Mars/msim/MSimMMCS.cc

@@ -84 +84 @@
         return kFALSE;
 
-    fPointingPos = (MPointingPos*)plist->FindObject("MPointingPos");
-    if (!fPointingPos)
+    fPointingTel = (MPointingPos*)plist->FindObject("MPointingPos");
+    if (!fPointingTel)
     {
         *fLog << err << "MPointingPos not found... aborting." << endl;
@@ -134 +134 @@
 
     //    fNumPheFromDNSB        MMcPedestalNSBAdd   // Number of phe/ns from diffuse NSB
+
+    // FIXME: Is there a way to write them as LAST entry in the file?
     mrh->SetNumSimulatedShowers(fRunHeader->GetNumEvents());
     mrh->SetImpactMax(fRunHeader->GetImpactMax());
+    mrh->SetCorsikaVersion(TMath::Nint(fRunHeader->GetProgramVersion()*100));
 
     // ----------------------------------------------------
@@ -148 +151 @@
     const UInt_t id = fRunHeader->GetParticleID();
 
+    // FIXME: Is there a way to write them as LAST entry in the file?
     rh->SetRunInfo(1, fRunHeader->GetRunNumber(), 0);
     rh->SetSourceInfo(MMcEvtBasic::GetParticleName(id));
@@ -166 +170 @@
     fMcEvtBasic->SetImpact(fEvtHeader->GetImpact());
 
-    fMcEvtBasic->SetTelescopeTheta(fPointingPos->GetZdRad());
-    fMcEvtBasic->SetTelescopePhi(fPointingPos->GetAzRad());
+    fMcEvtBasic->SetTelescopeTheta(fPointingTel->GetZdRad());
+    fMcEvtBasic->SetTelescopePhi(fPointingTel->GetAzRad());
 
     fMcEvtBasic->SetReadyToSave();
@@ -173 +177 @@
     static_cast<MMcEvtBasic&>(*fMcEvt) = *fMcEvtBasic;
 
+    // Convert from corsika franme to telescope frame, taking
+    // the magnetic field into account: tel = corsika+offset
     fMcEvt->SetTheta(fEvtHeader->GetZd());
-    fMcEvt->SetPhi(fEvtHeader->GetAz());
+    fMcEvt->SetPhi(fEvtHeader->GetAz()+fRunHeader->GetMagneticFieldAz());
 
     fMcEvt->SetEvtNumber(fEvtHeader->GetEvtNumber());

trunk/MagicSoft/Mars/msim/MSimMMCS.h

@@ -18 +18 @@
     MMcEvtBasic       *fMcEvtBasic;
     MMcEvt            *fMcEvt;
-    MPointingPos      *fPointingPos;
+    MPointingPos      *fPointingTel;  //! telescope poiting position in local (telescope) coordinate system
     MCorsikaEvtHeader *fEvtHeader;
     MCorsikaRunHeader *fRunHeader;

trunk/MagicSoft/Mars/msim/MSimPointingPos.cc

@@ -41 +41 @@
 //  Output Containers:
 //   MPointingPos
+//   PointingCorsika [MPointingPos]
 //
 //////////////////////////////////////////////////////////////////////////////
 #include "MSimPointingPos.h"
+
+#include <TRandom.h>
+#include <TVector3.h>
 
 #include "MLog.h"
@@ -64 +68 @@
 //
 MSimPointingPos::MSimPointingPos(const char* name, const char *title)
-    : fRunHeader(0), fEvtHeader(0), fPointing(0)
+    : fRunHeader(0), fEvtHeader(0), fPointingCorsika(0), fPointingLocal(0),
+    fOffTargetDistance(-1), fOffTargetPhi(-1)
+
 {
     fName  = name  ? name  : "MSimPointingPos";
@@ -77 +83 @@
 Int_t MSimPointingPos::PreProcess(MParList *pList)
 {
-    fPointing = (MPointingPos*)pList->FindCreateObj("MPointingPos");
-    if (!fPointing)
+    fPointingCorsika = (MPointingPos*)pList->FindCreateObj("MPointingPos", "PointingCorsika");
+    if (!fPointingCorsika)
+        return kFALSE;
+
+    fPointingLocal = (MPointingPos*)pList->FindCreateObj("MPointingPos");
+    if (!fPointingLocal)
         return kFALSE;
 
@@ -95 +105 @@
     }
 
+    // FIXED offset
+    // Diffuse?  ( FOV of camera folded with mirror diameter as Corsika input? )
+    // Hour angle!
+    // Angle to magnetic field!
+
+    if (IsOffTargetObservation())
+    {
+        *fLog << inf;
+        *fLog << "Off target observations switched on with" << endl;
+        *fLog <<"   a pointing distance of " << GetOffTargetDistance() << "deg ";
+        if (fOffTargetPhi<0)
+            *fLog << "randomly distributed in phi." << endl;
+        else
+            *fLog << "and phi=" << GetOffTargetPhi() << "deg." << endl;
+    }
+
     return kTRUE;
 }
+
+/*
+Bool_t MSimPointingPos::ReInit(MParList *pList)
+{
+    // FIXME: Check also the enlightened region on the ground!
+    return kTRUE;
+}
+*/
 
 // --------------------------------------------------------------------------
@@ -102 +136 @@
 Int_t MSimPointingPos::Process()
 {
+
     // If a view cone is given use the fixed telescope orientation
     const Bool_t fixed = fRunHeader->HasViewCone();
 
     // Local sky coordinates (direction of telescope axis)
-    const Double_t zd = fixed ? fRunHeader->GetZdMin() : fEvtHeader->GetZd()*TMath::RadToDeg();  // x==north
-    const Double_t az = fixed ? fRunHeader->GetAzMin() : fEvtHeader->GetAz()*TMath::RadToDeg();  // y==west
+    /*const*/ Double_t zd = fixed ? fRunHeader->GetZdMin() : fEvtHeader->GetZd()*TMath::RadToDeg();  // x==north
+    /*const*/ Double_t az = fixed ? fRunHeader->GetAzMin() : fEvtHeader->GetAz()*TMath::RadToDeg();  // y==west
+
+    if (!fixed && IsOffTargetObservation())
+    {
+        const Double_t theta = zd*TMath::DegToRad();
+        const Double_t phi   = az*TMath::DegToRad();
+
+        /*const*/ TVector3 source;
+        source.SetMagThetaPhi(1, theta, phi);
+
+        /*const*/ TVector3 point;
+        point.SetMagThetaPhi(1, theta+fOffTargetDistance, phi);
+
+        const Double_t delta = fOffTargetPhi>0 ? fOffTargetPhi :
+            gRandom->Uniform(TMath::TwoPi());
+
+        point.Rotate(delta, source);
+
+        zd = point.Theta()*TMath::RadToDeg();
+        az = point.Phi()*TMath::RadToDeg();
+    }
 
     // Setup the pointing position
-    fPointing->SetLocalPosition(zd, az);
+    fPointingCorsika->SetLocalPosition(zd, az);
+    fPointingLocal->SetLocalPosition(zd, az+fRunHeader->GetMagneticFieldAz());
+
+    // Calculate incident angle between magnetic field direction
+    // and pointing direction ( phi and theta? )
 
     return kTRUE;
 }
 
-/*
 Int_t MSimPointingPos::ReadEnv(const TEnv &env, TString prefix, Bool_t print)
 {
     Bool_t rc = kFALSE;
-    if (IsEnvDefined(env, prefix, "FileName", print))
+    if (IsEnvDefined(env, prefix, "OffTargetDistance", print))
     {
         rc = kTRUE;
-        SetFileName(GetEnvValue(env, prefix, "FileName", fFileName));
+        SetOffTargetDistance(GetEnvValue(env, prefix, "OffTargetDistance", GetOffTargetDistance()));
     }
 
-    if (IsEnvDefined(env, prefix, "UseTheta", print))
+    if (IsEnvDefined(env, prefix, "OffTargetPhi", print))
     {
         rc = kTRUE;
-        SetUseTheta(GetEnvValue(env, prefix, "UseTheta", fUseTheta));
+        SetOffTargetPhi(GetEnvValue(env, prefix, "OffTargetPhi", GetOffTargetPhi()));
     }
 
     return rc;
 }
-*/
+

trunk/MagicSoft/Mars/msim/MSimPointingPos.h

@@ -14 +14 @@
 {
 private:
-    MCorsikaRunHeader *fRunHeader;  //! Header storing event information
-    MCorsikaEvtHeader *fEvtHeader;  //! Header storing event information
-    MPointingPos      *fPointing;   //! Output storing telescope poiting position
+    MCorsikaRunHeader *fRunHeader;         //! Header storing event information
+    MCorsikaEvtHeader *fEvtHeader;         //! Header storing event information
+    MPointingPos      *fPointingCorsika;   //! Output storing telescope poiting position in corsika coordinate system (modulo magnetig field declination)
+    MPointingPos      *fPointingLocal;     //! Output storing telescope poiting position in local (telescope) coordinate system
+
+    Double_t fOffTargetDistance;  // [rad] Distance of the observed off-target position from the source
+    Double_t fOffTargetPhi;       // [rad] Rotation angle of the off-target position (phi==0 means south, phi=90 west) [0;2pi], phi<0 means random
 
     // MParContainer
-    //Int_t ReadEnv(const TEnv &env, TString prefix, Bool_t print=kFALSE);
+    Int_t ReadEnv(const TEnv &env, TString prefix, Bool_t print=kFALSE);
 
     // MTask
-    Int_t PreProcess(MParList *pList);
-    Int_t Process();
-
+    Int_t  PreProcess(MParList *pList);
+    //Bool_t ReInit(MParList *pList);
+    Int_t  Process();
 
 public:
     MSimPointingPos(const char *name=NULL, const char *title=NULL);
 
+    // Getter
+    Double_t GetOffTargetDistance() const { return fOffTargetDistance*TMath::RadToDeg(); }
+    Double_t GetOffTargetPhi() const { return fOffTargetPhi*TMath::RadToDeg(); }
+
+    // Setter
+    void SetOffTargetDistance(Double_t d=-1) { fOffTargetDistance = d*TMath::DegToRad(); }
+    void SetOffTargetPhi(Double_t p=-1) { fOffTargetPhi = p*TMath::DegToRad(); }
+
     // MSimPointingPos
+    Bool_t IsOffTargetObservation() const { return fOffTargetDistance>0; }
 
     // TObject