Changeset 7173 for trunk/MagicSoft/Mars/mhflux

Timestamp:

06/28/05 10:22:44 (20 years ago)

Author:

tbretz

Message:

*** empty log message ***

Location:

trunk/MagicSoft/Mars/mhflux

Files:

• MHDisp.cc (modified) (9 diffs)
• MHDisp.h (modified) (2 diffs)

trunk/MagicSoft/Mars/mhflux/MHDisp.cc

@@ -64 +64 @@
 //
 MHDisp::MHDisp(const char *name, const char *title)
-    : fHilExt(0), fDisp(0)
+    : fHilExt(0), fDisp(0)//, fIsWobble(kFALSE)
 {
     //
@@ -78 +78 @@
     fHist.SetXTitle("x [\\circ]");
     fHist.SetYTitle("y [\\circ]");
-    fHist.SetZTitle("\\vartheta [deg^{2}]");
+    fHist.SetZTitle("Eq.cts");
 }
 
@@ -130 +130 @@
         return kFALSE;
     }
+
+    //const Double_t xmax = fHist.GetXaxis()->GetXmax();
 
     // Initialize all bins with a small (=0) value otherwise
     // these bins are not displayed
     for (int x=0; x<fHist.GetNbinsX(); x++)
-        for (int y=0; y<fHist.GetNbinsX(); y++)
-            fHist.Fill(fHist.GetXaxis()->GetBinCenter(x+1),
-                       fHist.GetYaxis()->GetBinCenter(y+1),
-                       0.0, 1e-10);
+        for (int y=0; y<fHist.GetNbinsY(); y++)
+        {
+            const Double_t cx = fHist.GetXaxis()->GetBinCenter(x+1);
+            const Double_t cy = fHist.GetYaxis()->GetBinCenter(y+1);
+            //if (TMath::Hypot(cx, cy)<xmax)
+                fHist.Fill(cx, cy, 0.0, 1e-10);
+        }
 
     return kTRUE;
@@ -160 +165 @@
         rho = fPointPos->RotationAngle(*fObservatory, *fTime);
 
+    // FIXME: Do wobble-rotation when subtracting?
+    if (!fHistOff/* && fIsWobble*/)
+        rho += TMath::Pi();
+
     // Get Disp from Parlist
     const Double_t disp = fDisp->GetVal();
 
-    // Calculate where disp is pointing
+    // Calculate both postiions where disp could point
     TVector2 pos1(hil->GetCosDelta(), hil->GetSinDelta());
     TVector2 pos2(hil->GetCosDelta(), hil->GetSinDelta());
-    pos1 *=  disp;
-    pos2 *= -disp;
-
+    pos1 *=  disp; // Vector of length  disp in direction of shower
+    pos2 *= -disp; // Vector of length -disp in direction of shower
+
+    // Move origin of vector to center-of-gravity of shower
     pos1 += hil->GetMean()*fMm2Deg;
     pos2 += hil->GetMean()*fMm2Deg;
@@ -207 +217 @@
         srcp = fSrcPos->GetXY();
 
+    // Derotate all position around camera center by -rho
     if (rho!=0)
     {
@@ -214 +225 @@
     }
 
+    // Shift the source position to 0/0
     pos1 -= srcp*fMm2Deg;
     pos2 -= srcp*fMm2Deg;
 
-    fHist.Fill(pos1.X(), pos1.Y(), 0.0, w*gweight);
-    fHist.Fill(pos2.X(), pos2.Y(), 0.0, w*(1-gweight));
+    //const Double_t xmax = fHist.GetXaxis()->GetXmax();
+
+    // Fill histograms
+    //if (pos1.Mod()<xmax)
+        fHist.Fill(pos1.X(), pos1.Y(), 0.0, w*gweight);
+    //if (pos2.Mod()<xmax)
+        fHist.Fill(pos2.X(), pos2.Y(), 0.0, w*(1-gweight));
 
     return kTRUE;
 }
-/*
-static Double_t FcnGauss2d(Double_t *x, Double_t *par)
-{
-    TVector2 v = TVector2(x[0], x[1]).Rotate(par[5]*TMath::DegToRad());
-
-    const Double_t g0 = TMath::Gaus(v.X(), par[1], par[2]);
-    const Double_t g1 = TMath::Gaus(v.Y(), par[3], par[4]);
-
-    //Gaus(Double_t x, Double_t mean=0, Double_t sigma=1, Bool_t norm=kFALSE);
-    return par[0]*g0*g1  + par[6]*(v.X()*v.X() + v.Y()*v.Y()) +  par[7];
-}*/
+
+Double_t MHDisp::GetOffSignal(TH1 &h) const
+{
+    const TAxis &axex = *h.GetXaxis();
+    const TAxis &axey = *h.GetYaxis();
+
+    Double_t sum = 0;
+    for (int x=0; x<h.GetNbinsX(); x++)
+        for (int y=0; y<h.GetNbinsY(); y++)
+        {
+            if (TMath::Hypot(axex.GetBinCenter(x+1), axey.GetBinCenter(y+1))>0.35)
+                sum += h.GetBinContent(x+1,y+1);
+        }
+
+    return sum;
+}
 
 void MHDisp::Paint(Option_t *o)
@@ -240 +262 @@
     pad->cd(1);
 
+    // Project on data onto yx-plane
     fHist.GetZaxis()->SetRange(0,9999);
     TH1 *h1=fHist.Project3D("yx_on");
-    gStyle->SetPalette(1, 0);
+
+    // Set Glow-palette (PaintSurface doesn't allow more than 99 colors)
+    MH::SetPalette(fHistOff?"glowsym":"glow1", 99);
+    h1->SetContour(99);
+
+    Double_t scale = 1;
+    if (fHistOff)
+    {
+        // Project off data onto yx-plane and subtract it from on-data
+        fHistOff->GetZaxis()->SetRange(0,9999);
+        TH1 *h=fHistOff->Project3D("yx_off");
+
+        scale = -1;
+
+        //if (!fIsWobble)
+        {
+            const Double_t h1off = GetOffSignal(*h1);
+            const Double_t hoff  = GetOffSignal(*h);
+            scale = hoff==0?-1:-h1off/hoff;
+        }
+
+        h1->Add(h, scale);
+        delete h;
+
+        // Force calculation of minimum, maximum
+        h1->SetMinimum();
+        h1->SetMaximum();
+
+        // Find maximum
+        const Double_t max = TMath::Max(TMath::Abs(h1->GetMinimum()),
+                                        TMath::Abs(h1->GetMaximum()));
+
+        // Set new minimum, maximum
+        h1->SetMinimum(-max);
+        h1->SetMaximum( max);
+    }
 
     Int_t ix, iy, iz;
@@ -304 +362 @@
         const Double_t s  = MMath::SignificanceLiMa(e, b);
 
-        h2->SetTitle(Form("P=(%.2f\\circ/%.2f\\circ) \\omega=%.2f\\circ \\sigma=%.1f E=%.0f B=%.0f", x0, y0, func.GetParameter(2), s, e-b, b));
+        h2->SetTitle(Form("P=(%.2f\\circ/%.2f\\circ) \\omega=%.2f\\circ \\sigma=%.1f E=%.0f B=%.0f f=%.2f", x0, y0, func.GetParameter(2), s, e-b, b, scale));
     }
    /*
@@ -386 +444 @@
     h->Draw();
 }
+
+// --------------------------------------------------------------------------
+//
+// The following resources are available:
+//
+//    MHDisp.Wobble: on/off
+//
+/*
+Int_t MHFalseSource::ReadEnv(const TEnv &env, TString prefix, Bool_t print)
+{
+    Bool_t rc = kFALSE;
+    if (IsEnvDefined(env, prefix, "DistMin", print))
+    {
+        rc = kTRUE;
+        SetMinDist(GetEnvValue(env, prefix, "DistMin", fMinDist));
+    }
+    if (IsEnvDefined(env, prefix, "DistMax", print))
+    {
+        rc = kTRUE;
+        SetMaxDist(GetEnvValue(env, prefix, "DistMax", fMaxDist));
+    }
+
+    if (IsEnvDefined(env, prefix, "DWMin", print))
+    {
+        rc = kTRUE;
+        SetMinDW(GetEnvValue(env, prefix, "DWMin", fMinDist));
+    }
+    if (IsEnvDefined(env, prefix, "DWMax", print))
+    {
+        rc = kTRUE;
+        SetMaxDW(GetEnvValue(env, prefix, "DWMax", fMaxDist));
+    }
+
+    return rc;
+}
+*/

trunk/MagicSoft/Mars/mhflux/MHDisp.h

@@ -13 +13 @@
 {
 private:
-    MHillasExt  *fHilExt; //!
-    MParameterD *fDisp;   //!
+    MHillasExt  *fHilExt;  //!
+    MParameterD *fDisp;    //!
 
-    Double_t fM3lCut;     //!
-    Double_t fXi;         //!
-    Double_t fXiTheta;    //!
+    Double_t     fM3lCut;  //!
+    Double_t     fXi;      //!
+    Double_t     fXiTheta; //!
+
+    //Bool_t       fIsWobble;
+
+    // MHDisp
+    Double_t GetOffSignal(TH1 &h) const;
 
 public:
@@ -27 +32 @@
     Bool_t Fill(const MParContainer *par, const Stat_t w=1);
 
+    // TObject
     void Paint(Option_t *o="");
     void Draw(Option_t *o="");
+
+    // MParContainer
+    //Int_t ReadEnv(const TEnv &env, TString prefix, Bool_t print=kFALSE);
 
     ClassDef(MHDisp, 1) //3D-histogram in alpha, x and y