Changeset 8154 for trunk/MagicSoft/Mars/mextralgo/MExtralgoDigitalFilter.h

Timestamp:

10/24/06 09:26:10 (18 years ago)

Author:

tbretz

Message:

*** empty log message ***

File:

• trunk/MagicSoft/Mars/mextralgo/MExtralgoDigitalFilter.h (modified) (6 diffs)

trunk/MagicSoft/Mars/mextralgo/MExtralgoDigitalFilter.h

@@ -2 +2 @@
 #define MARS_MExtralgoDigitalFilter
 
-#ifndef ROOT_TROOT
-#include <TROOT.h>
+#ifndef ROOT_TMatrix
+#include <TMatrix.h>
 #endif
 
@@ -12 +12 @@
 class TArrayF;
 
+//#include <TMatrix.h>
+
 class MExtralgoDigitalFilter
 {
 private:
     // Input
-    Float_t *fVal;
-    Int_t    fNum;
+    const Float_t *fVal;
+    Int_t fNum;
 
     Float_t const *fWeightsAmp;
     Float_t const *fWeightsTime;
+    Float_t const *fPulseShape;
+
+    const TMatrix *fAinv;
 
     const Int_t fWeightsPerBin; // Number of weights per data bin
@@ -30 +35 @@
     Float_t fSignal;
     Float_t fSignalDev;
+
+    Float_t GetChisq(const Int_t maxp, const Int_t frac, const Float_t sum) const;
+
+    inline Double_t ChiSq(const Double_t sum, const Int_t startv, const Int_t startw=0) const
+    {
+        //
+        // Slide with a window of size windowsize over the sample
+        // and multiply the entries with the corresponding weights
+        //
+        Double_t chisq = 0;
+
+        // Shift the start of the weight to the center of sample 0
+        Float_t const *w = fPulseShape + startw;
+
+        const Float_t *beg = fVal+startv;
+        for (Float_t const *pex=beg; pex<beg+fWindowSize; pex++)
+        {
+            const Double_t c = *w - *pex/sum;
+            chisq += c*c;
+            w += fWeightsPerBin;
+        }
+        return chisq;
+    }
 
     // Weights: Weights to evaluate
@@ -45 +73 @@
         Float_t const *w = weights + startw;
 
-        Float_t *const beg = fVal+startv;
+        const Float_t *beg = fVal+startv;
         for (Float_t const *pex=beg; pex<beg+fWindowSize; pex++)
         {
@@ -115 +143 @@
 
 public:
-    MExtralgoDigitalFilter(Int_t res, Int_t windowsize, Float_t *wa, Float_t *wt)
+    MExtralgoDigitalFilter(Int_t res, Int_t windowsize, Float_t *wa, Float_t *wt, Float_t *ps=0, TMatrix *ainv=0)
         : fVal(0), fNum(0), fWeightsAmp(wa+res/2), fWeightsTime(wt+res/2),
-        fWeightsPerBin(res), fWindowSize(windowsize),
+        fPulseShape(ps), fAinv(ainv), fWeightsPerBin(res), fWindowSize(windowsize),
         fTime(0), fTimeDev(-1), fSignal(0), fSignalDev(-1)
     {
     }
 
-    void SetData(Int_t n, Float_t *val) { fNum=n; fVal=val; }
+    void SetData(Int_t n, Float_t const *val) { fNum=n; fVal=val; }
 
     Float_t GetTime() const          { return fTime; }
@@ -134 +162 @@
 
     Float_t ExtractNoise(Int_t iter) const;
-    void Extract();
+    void Extract(Int_t maxpos=-1);
 
-    static Bool_t CalculateWeights(TH1 &shape, const TH2 &autocorr, TArrayF &wa, TArrayF &wt, Int_t wpb=-1);
-    static void CalculateWeights2(TH1F &shape, const TH2F &autocorr, TArrayF &wa, TArrayF &wt, Int_t wpb=-1);
+    static Int_t CalculateWeights(TH1 &shape, const TH2 &autocorr, TArrayF &wa, TArrayF &wt, Int_t wpb=-1);
+    static Int_t CalculateWeights2(TH1 &shape, const TH2 &autocorr, TArrayF &wa, TArrayF &wt, Int_t wpb=-1);
 };