Changeset 9849 for trunk/Mars/mimage

Timestamp:

08/11/10 17:59:53 (14 years ago)

Author:

tbretz

Message:

Fixed MHillasExt::fTimeSpread and implemented calculation of fSlopeSpread.

File:

• trunk/Mars/mimage/MHillasExt.cc (modified) (5 diffs)

trunk/Mars/mimage/MHillasExt.cc

@@ -153 +153 @@
     Double_t sumxt   = 0;
     Double_t sumyt   = 0;
+    Double_t sumxw   = 0;
+    Double_t sumyw   = 0;
     Double_t sumtw   = 0;
     Double_t sumx2   = 0;
@@ -160 +162 @@
     Double_t sumx3   = 0;
     Double_t sumy3   = 0;
+    Double_t sumxw2  = 0;
+    Double_t sumyw2  = 0;
+    Double_t sumtx2  = 0;
+    Double_t sumty2  = 0;
     Double_t sumtw2  = 0;
     Double_t sumx2y  = 0;
     Double_t sumxy2  = 0;
+    Double_t sumtxw2 = 0;
+    Double_t sumtyw2 = 0;
     Double_t sumt2w2 = 0;
 
@@ -190 +198 @@
         sumt2   += t*t;
 
-        sumtw   += t*nphot;
-        sumtw2  += t*nphot*nphot;
-        sumt2w2 += t*t*nphot*nphot;
-
         sumw    += nphot;
         sumw2   += nphot*nphot;
+
+        sumxw   += x*nphot;
+        sumyw   += x*nphot;
+        sumtw   += t*nphot;
+
+        sumtx2  += x*nphot*nphot;
+        sumty2  += y*nphot*nphot;
+        sumtw2  += t*nphot*nphot;
+
+        sumtxw2 += t*x*nphot*nphot;
+        sumtyw2 += t*y*nphot*nphot;
+        sumt2w2 += t*t*nphot*nphot;
 
         sumx2 += x*x;
@@ -235 +251 @@
     //
     fTimeSpread         = TMath::Sqrt(sumt2/cnt - sumt*sumt/cnt/cnt);
-    fTimeSpreadWeighted = TMath::Sqrt(sumt2w2 + sumw2*sumtw/sumw - 2*sumtw/sum2*sumtw2)/sumw;
+    fTimeSpreadWeighted = TMath::Sqrt(sumt2w2 + (sumw2 - 2*sumtw2)*sumtw/sumw)/sumw;
 
     //
@@ -254 +270 @@
     fSlopeLong  = detx==0 ? 0 : (cnt*dxt - sumt*dx)/detx;
     fSlopeTrans = dety==0 ? 0 : (cnt*dyt - sumt*dy)/dety;
+
+    //
+    // Calculate time spread around longitudinal time slope
+    //
+    const Double_t sumdt     = sumt    -   fSlopeLong*(c*sumx    + s*sumy);
+    const Double_t sumdtw    = sumtw   -   fSlopeLong*(c*sumxw   + s*sumyw);
+    const Double_t sumdtw2   = sumtw2  -   fSlopeLong*(c*sumxw2  + s*sumyw2);
+    const Double_t sumdt2    = sumt2   - 2*fSlopeLong*(c*sumxt   + s*sumyt)   + fSlopeLong*fSlopeLong * (c*c*sumx2 + s*s*sumy2 + 2*c*s*sumxy);
+    const Double_t sumdt2w2  = sumt2w2 - 2*fSlopeLong*(c*sumtxw2 + s*sumtyw2) + fSlopeLong*fSlopeLong * (c*c*sumx2 + s*s*sumy2 + 2*c*s*sumxy) * sumw2;
+
+    fSlopeSpread         = TMath::Sqrt(sumdt2/cnt - sumdt*sumdt/cnt/cnt);
+    fSlopeSpreadWeighted = TMath::Sqrt(sumdt2w2 + (sumw2 - 2*sumdtw2)*sumdtw/sumw)/sumw;
 
     //