Changeset 9845 for trunk/Mars

Timestamp:

08/11/10 17:05:32 (14 years ago)

Author:

tbretz

Message:

Implemented new parameters to describe the time spread in MHillasExt.

Location:

trunk/Mars

Files:

• Changelog (modified) (1 diff)
• NEWS (modified) (1 diff)
• mimage/MHillasExt.cc (modified) (7 diffs)
• mimage/MHillasExt.h (modified) (3 diffs)

trunk/Mars/Changelog

@@ -26 +26 @@
      - improved/fixed error calculation/display
 
+   * mimage/MHillasExt.[h,cc]:
+     - added new data members fTimeSpread and fTimeSpreadWeighted
+     - added new data members fSlopeSpread and fSlopeSpreadWeighted (not yet used)
+     - implemented code to calculate the time spreads
+     - increased class version id by 1
+
 
 

trunk/Mars/NEWS

@@ -51 +51 @@
      ns anymore, but in ns/deg. This means that the values for MAGIC
      have to be scaled by a factor of 10.
+
+   * Star now calculates new timing parameters, namely:
+
+   ** MHillasExt.fTimeSpread: 
+      This is the time spread (rms) around the mean of the arrival times.
+      It is calulated as 
+         <math>\sqrt{\sum{t_i}/n}</math>
+   ** MHillasExt.fTimeSpreadWeighted: 
+      This is the weighted time spread (rms) around the weighted mean of the arrival times.
+      It is calulated as 
+         <math>\frac{1}{\sum{n_i}}\sqrt{\sum{n_i^2(t_i-<wt>)^2}}}</math>
+            or
+         <math>\frac{1}{\sum{n_i}}\sqrt{\sum{n_i^2t_i^2}+<wt>^2\sum{n_i^2}-2<wt>\sum{t_in_i^2}}</math>
 
  ;trainenergy.C:

trunk/Mars/mimage/MHillasExt.cc

@@ -20 +20 @@
 !   Author(s): Wolfgang Wittek 06/2002 <mailto:wittek@mppmu.mpg.de>
 !
-!   Copyright: MAGIC Software Development, 2000-2007
+!   Copyright: MAGIC Software Development, 2000-2010
 !
 !
@@ -56 +56 @@
 // fSlopeLong    added
 // fSlopeTrans   added
+//
+// Version 5:
+// ----------
+// fTimeSpread           added
+// fTimeSpreadWeighted   added
+//
+// fSlopeSpread          added
+// fSlopeSpreadWeighted  added
 //
 //
@@ -138 +146 @@
     UInt_t cnt = 0;
 
-    Double_t sumx   = 0;
-    Double_t sumy   = 0;
-    Double_t sumt   = 0;
-    Double_t sumxy  = 0;
-    Double_t sumxt  = 0;
-    Double_t sumyt  = 0;
-    Double_t sumx2  = 0;
-    Double_t sumy2  = 0;
-    Double_t sumx3  = 0;
-    Double_t sumy3  = 0;
-    Double_t sumx2y = 0;
-    Double_t sumxy2 = 0;
+    Double_t sumx    = 0;
+    Double_t sumy    = 0;
+    Double_t sumt    = 0;
+    Double_t sumw    = 0;
+    Double_t sumxy   = 0;
+    Double_t sumxt   = 0;
+    Double_t sumyt   = 0;
+    Double_t sumtw   = 0;
+    Double_t sumx2   = 0;
+    Double_t sumy2   = 0;
+    Double_t sumw2   = 0;
+    Double_t sumt2   = 0;
+    Double_t sumx3   = 0;
+    Double_t sumy3   = 0;
+    Double_t sumtw2  = 0;
+    Double_t sumx2y  = 0;
+    Double_t sumxy2  = 0;
+    Double_t sumt2w2 = 0;
 
     const UInt_t npix = evt.GetNumPixels();
@@ -167 +181 @@
         const Double_t t = pix.GetArrivalTime();
 
+        Double_t nphot = pix.GetNumPhotons();        // [1]
+
         // --- time slope ----
-        sumx  += x;
-        sumy  += y;
-        sumt  += t;
+        sumx    += x;
+        sumy    += y;
+
+        sumt    += t;
+        sumt2   += t*t;
+
+        sumtw   += t*nphot;
+        sumtw2  += t*nphot*nphot;
+        sumt2w2 += t*t*nphot*nphot;
+
+        sumw    += nphot;
+        sumw2   += nphot*nphot;
 
         sumx2 += x*x;
@@ -181 +206 @@
         const Double_t dx = x - hil.GetMeanX();      // [mm]
         const Double_t dy = y - hil.GetMeanY();      // [mm]
-
-        Double_t nphot = pix.GetNumPhotons();        // [1]
 
         sumx3  += nphot * dx*dx*dx;
@@ -209 +232 @@
 
     //
+    // Time spread
+    //
+    fTimeSpread         = TMath::Sqrt(sumt2/cnt - sumt*sumt/cnt/cnt);
+    fTimeSpreadWeighted = TMath::Sqrt(sumt2w2 + sumw2*sumtw/sumw - 2*sumtw/sum2*sumtw2)/sumw;
+
+    //
     // Time slope
     //
@@ -225 +254 @@
     fSlopeLong  = detx==0 ? 0 : (cnt*dxt - sumt*dx)/detx;
     fSlopeTrans = dety==0 ? 0 : (cnt*dyt - sumt*dy)/dety;
+
+    //
+    // FIXME: Missing Time spread along slope
+    //
 
     //

trunk/Mars/mimage/MHillasExt.h

@@ -20 +20 @@
     Float_t fM3Trans; // [mm] 3rd moment (e-weighted) along minor axis
 
-    Float_t fSlopeLong;
-    Float_t fSlopeTrans;
+    Float_t fSlopeLong;   // [ns/mm] Slope of the arrival time along the major axis
+    Float_t fSlopeTrans;  // [ns/mm] Slope of the arrival time along the minor axis
+
+    // FIXME: We could try to use the Medians!
+    Float_t fTimeSpread;          // [ns] Spread (RMS) of the arrival times around the mean arrival time
+    Float_t fTimeSpreadWeighted;  // [ns] Weighted spread (RMS) of the arrival times around the weighted mean arrival time
+
+    Float_t fSlopeSpread;         // [ns] Spread (RMS) of the arrival time around the slope along the major axis
+    Float_t fSlopeSpreadWeighted; // [ns] Weighted spread (RMS) of the arrival time around the slope along the major axis
 
 public:
@@ -28 +35 @@
     void Reset();
 
-    Float_t GetAsym() const       { return fAsym; }
-    Float_t GetM3Long() const     { return fM3Long; }
-    Float_t GetM3Trans() const    { return fM3Trans; }
-    Float_t GetSlopeLong() const  { return fSlopeLong; }
-    Float_t GetSlopeTrans() const { return fSlopeTrans; }
+    Float_t GetAsym() const               { return fAsym; }
+    Float_t GetM3Long() const             { return fM3Long; }
+    Float_t GetM3Trans() const            { return fM3Trans; }
+    Float_t GetSlopeLong() const          { return fSlopeLong; }
+    Float_t GetSlopeTrans() const         { return fSlopeTrans; }
+    Float_t GetTimeSpread() const         { return fTimeSpread; }
+    Float_t GetTimeSpreadWeighted() const { return fTimeSpreadWeighted; }
 
     Int_t Calc(const MGeomCam &geom, const MSignalCam &pix,
@@ -44 +53 @@
     void Set(const TArrayF &arr);
 
-    ClassDef(MHillasExt, 4) // Storage Container for extended Hillas Parameter
+    ClassDef(MHillasExt, 5) // Storage Container for extended Hillas Parameter
 };
 #endif